diff --git a/datasets/raise/split_and_compare_prnus_of_subgroups.py b/datasets/raise/split_and_compare_prnus_of_subgroups.py
index 50f6c22..af55cf3 100644
--- a/datasets/raise/split_and_compare_prnus_of_subgroups.py
+++ b/datasets/raise/split_and_compare_prnus_of_subgroups.py
@@ -1,12 +1,34 @@
 from PIL import Image
 import numpy as np
+import matplotlib.pyplot as plt
+from tqdm import tqdm
+from skimage.restoration import denoise_tv_chambolle
 import sys
 
 sys.path.insert(0, '../../algorithms/image_utils/')
 
 from image_utils import showImageWithMatplotlib
 
-imagePath = 'flat-field/TIF/flat_001.tif'
-imagePil = Image.open(imagePath).convert('F')
-imageNpArray = np.array(imagePil)
-showImageWithMatplotlib(imageNpArray)
\ No newline at end of file
+NUMBER_OF_SUBGROUPS = 2
+fig, axes = plt.subplots(NUMBER_OF_SUBGROUPS, 1)
+fig.suptitle(f'PRNU estimation and comparison for {NUMBER_OF_SUBGROUPS} subgroups with different number of flat-images')
+
+IMAGES_FOLDER = 'flat-field/TIF'
+imagesFileNames = os.listdir(IMAGES_FOLDER)
+numberOfImagesPerSubgroup = len(imagesFileNames) // NUMBER_OF_SUBGROUPS
+
+# Assume random image order to not introduce a bias.
+for subgroupIndex in range(NUMBER_OF_SUBGROUPS):
+    axis = axes[subgroupIndex]
+    imagesPrnuEstimateNpArray = []
+    for imageFileName in tqdm(imagesFileNames[numberOfImagesPerSubgroup * subgroupIndex : numberOfImagesPerSubgroup * (subgroupIndex + 1)]):
+        imagePath = f'{IMAGES_FOLDER}/{imageFileName}'
+        imagePil = Image.open(imagePath).convert('F')
+        imageNpArray = np.array(imagePil)
+        imagePrnuEstimateNpArray = imageNpArray - denoise_tv_chambolle(imageNpArray, weight=0.2, channel_axis=-1)
+        imagesPrnuEstimateNpArray += [imagePrnuEstimateNpArray]
+    cameraPrnuEstimateNpArray = np.array(imagesPrnuEstimateNpArray).mean(axis = 0)
+    axis.set_title(f'PRNU estimate for subgroup {subgroupIndex} with ')
+    axis.imshow(cameraPrnuEstimateNpArray)
+
+plt.show()
\ No newline at end of file