diff --git a/scripts/convert_images.py b/scripts/convert_images.py
index 690dbc3..9e3b805 100644
--- a/scripts/convert_images.py
+++ b/scripts/convert_images.py
@@ -2,33 +2,54 @@
 import cv2
 import argparse
 import numpy as np
-
+import rasterio
 from PIL import Image
-from deep_image_matching import logger
 
 
 def load_img(
     path_to_img: str,
+    reader: str,
 ):
-    #img = cv2.imread(path_to_img, cv2.IMREAD_ANYDEPTH)
-    logger.debug(
-        "Image loaded with pillow .."
-    )
-    img = Image.open(path_to_img)
-    img = img.convert("L")
-    logger.info(f"Image mode: {img.mode}")
-    img = np.array(img)
-    logger.info(f"Image mode: {img.dtype}")
+    """
+    Load an image from the given path using the specified reader.
+
+    Args:
+        path_to_img (str): Path to the image file.
+        reader (str): The image reader to use. Supported options are "rasterio", "opencv", and "pillow".
+
+    Returns:
+        numpy.ndarray: The loaded image as a NumPy array.
+    """
+    if reader == "rasterio":
+        with rasterio.open(path_to_img) as dataset:
+            img = dataset.read(1)
+
+    elif reader == "opencv":
+        img = cv2.imread(path_to_img, cv2.IMREAD_ANYDEPTH)
+
+    elif reader == "pillow":
+        img = Image.open(path_to_img)
+        img = img.convert("L")
+        img = np.array(img)
 
     return img
 
 
-def convert_images(input_folder, output_folder, target_format, normalize, method):
-    # Create the output folder if it doesn't exist
+def convert_images(input_folder, output_folder, target_format, normalize, method, reader):
+    """
+    Convert images in a folder to a target image format using OpenCV.
+
+    Args:
+        input_folder (str): Path to the folder containing input images.
+        output_folder (str): Path to the folder where converted images will be saved.
+        target_format (str): Target image format (e.g., 'png', 'jpg', 'jpeg', 'tiff').
+        normalize (bool): Flag indicating whether to normalize pixel values to cover the whole range 0-255.
+        method (str): Normalization method. Supported options are "single_image" and "image_set".
+        reader (str): Library to read the images. Supported options are "opencv", "pillow", and "rasterio".
+    """
     if not os.path.exists(output_folder):
         os.makedirs(output_folder)
 
-    # Get a list of all image files in the input folder
     image_files = [
         f
         for f in os.listdir(input_folder)
@@ -55,7 +76,7 @@ def convert_images(input_folder, output_folder, target_format, normalize, method
         )
     ]
     if len(image_files) == 0:
-        logger.error(
+        print(
             "No images in the folder. Check that the path to the folder is correct or that the format is supported."
         )
         quit()
@@ -66,7 +87,7 @@ def convert_images(input_folder, output_folder, target_format, normalize, method
 
     for image_file in image_files:
         image_path = os.path.join(input_folder, image_file)
-        img = load_img(image_path)
+        img = load_img(image_path, reader)
 
         MIN.append(np.min(img))
         MAX.append(np.max(img))
@@ -75,15 +96,14 @@ def convert_images(input_folder, output_folder, target_format, normalize, method
     abs_max = np.max(MAX)
 
     for i, image_file in enumerate(image_files):
-        # Read the image
-        logger.info(f"[IMAGE{i}]")
+        print(f"[IMAGE{i}]")
         image_path = os.path.join(input_folder, image_file)
         output_path = os.path.join(
             output_folder, os.path.splitext(image_file)[0] + "." + target_format
         )
-        img = load_img(image_path)
+        img = load_img(image_path, reader)
         min_original, max_original = np.min(img), np.max(img)
-        logger.info(f"Range: min_original {min_original}  max_original {max_original}")
+        print(f"Range: min_original {min_original}  max_original {max_original}")
 
         if not normalize:
             if img.dtype == np.uint8:
@@ -116,7 +136,6 @@ def convert_images(input_folder, output_folder, target_format, normalize, method
 
 
 if __name__ == "__main__":
-    # Parse command-line arguments
     parser = argparse.ArgumentParser(
         description="Convert images in a folder to a target image format using OpenCV."
     )
@@ -135,18 +154,23 @@ def convert_images(input_folder, output_folder, target_format, normalize, method
         help="Normalize pixel values if set to cover the whole range 0-255.",
     )
     parser.add_argument(
-        "-m",
         "--method",
         choices=["single_image", "image_set"],
         help="Normalize respect max and min of each single image or respect max and min of the entire image set.",
     )
+    parser.add_argument(
+        "--reader",
+        choices=["opencv", "pillow", "rasterio"],
+        default="pillow",
+        help="Library to read the images. Default is 'pillow'.",
+    )
     args = parser.parse_args()
 
-    # Call the function to convert images
     convert_images(
         args.input_folder,
         args.output_folder,
         args.target_format,
         args.normalize,
         args.method,
+        args.reader,
     )
diff --git a/scripts/normalize_images.py b/scripts/normalize_images.py
index 81e381f..d8ffdca 100644
--- a/scripts/normalize_images.py
+++ b/scripts/normalize_images.py
@@ -5,96 +5,168 @@
 
 from pathlib import Path
 
-def img_gradient(image):
-    blurred_image = cv2.GaussianBlur(image, (3, 3), 0)
-    grad_x = cv2.Sobel(blurred_image, cv2.CV_64F, 1, 0, ksize=3)
-    grad_y = cv2.Sobel(blurred_image, cv2.CV_64F, 0, 1, ksize=3)
-    gradient_magnitude = np.sqrt(grad_x**2 + grad_y**2)
-    gradient_magnitude = cv2.convertScaleAbs(gradient_magnitude)
-    normalized_gradient = cv2.normalize(gradient_magnitude, None, 0, 255, cv2.NORM_MINMAX)
-    return normalized_gradient
-
-def normalize_image(image_path, output_path, patch_size=300, method="minmax"):
-    # Read the image in grayscale mode
-    image = cv2.imread(str(image_path), cv2.IMREAD_GRAYSCALE)
-    image = cv2.GaussianBlur(image, (5, 5), 0)
-    if image is None:
-        raise ValueError(f"Image at {image_path} could not be read.")
-
-    if method == "gradient":
-        gradient_magnitude = img_gradient(image)
-        cv2.imwrite(str(output_path), gradient_magnitude)
-        return
-
-    # Get the dimensions of the image
-    rows, cols = image.shape
-    
-    # Create an empty array to store the normalized image
-    normalized_image = np.zeros_like(image, dtype=np.float32)
-    
-    # Define the half patch size
-    half_patch_size = patch_size // 2
-    
-    # Iterate over each pixel in the image
-    for i in range(rows):
-        for j in range(cols):
-            # Define the boundaries of the patch
-            r_min = max(0, i - half_patch_size)
-            r_max = min(rows, i + half_patch_size + 1)
-            c_min = max(0, j - half_patch_size)
-            c_max = min(cols, j + half_patch_size + 1)
-            
-            # Extract the patch
-            patch = image[r_min:r_max, c_min:c_max]
-            
-            # Calculate the mean and standard deviation of the patch
-            patch_mean = np.mean(patch)
-            patch_std = np.std(patch)
-            #masked_patch = np.ma.masked_outside(patch, patch_mean - 2 * patch_std, patch_mean + 2 * patch_std)
-            #robust_mean = np.mean(masked_patch)
-            #robust_std = np.std(masked_patch)
-            min_value = np.min(patch)
-            max_value = np.max(patch)
-            
-            ## Normalize the pixel with respect to the patch
-            #if patch_std > 0:
-            #    normalized_image[i, j] = (image[i, j] - patch_mean) / patch_std
-            #else:
-            #    normalized_image[i, j] = 0  # If the standard deviation is 0, set the normalized value to 0
-
-            #norm_value = int((max_value-min_value)*image[i,j]/255)
-
-            if method == "minmax":
-                norm_value = (image[i,j]-min_value)*255/(max_value-min_value)
-                if norm_value > 0 and norm_value < 256:
-                    normalized_image[i, j] = norm_value
-                else:
-                    normalized_image[i, j] = 0
-            elif method == "robust":
-                norm_value = (image[i,j]-patch_mean-patch_std)*255/(2*patch_std)
-                if norm_value > 0 and norm_value < 256:
-                    normalized_image[i, j] = norm_value
-                else:
-                    normalized_image[i, j] = 0
-    
-    # Scale the normalized image to the range [0, 255]
-    #normalized_image = cv2.normalize(normalized_image, None, 0, 255, cv2.NORM_MINMAX)
-    
-    # Convert the normalized image to uint8 type
-    normalized_image = normalized_image.astype(np.uint8)
+
+class ImageNormalizer:
+    """
+    A class for normalizing images using different methods.
+    """
+
+    def __init__(
+            self, 
+            reduce_noise: bool = True,
+            noise_kernel: tuple = (5, 5),
+            ):
+        """
+        Initialize the ImageNormalizer object.
+
+        Args:
+            reduce_noise (bool): Whether to apply noise reduction. Default is True.
+            noise_kernel (tuple): The size of the noise reduction kernel. Default is (5, 5).
+        """
+        self.reduce_noise = reduce_noise
+        self.noise_kernel = noise_kernel
+
+    def img_gradient(
+            self,
+            image_path: Path,
+            output_path: Path,
+            ksize: int = 3,
+            ):
+        """
+        Export the gradient of the image.
+
+        Args:
+            image_path (Path): The path to the input image.
+            output_path (Path): The path to save the normalized image.
+            ksize (int): The size of the Sobel kernel. Default is 3.
+        """
+        image = cv2.imread(str(image_path), cv2.IMREAD_GRAYSCALE)
+        if image is None:
+            raise ValueError(f"Image at {image_path} could not be read.")
+        if self.reduce_noise:
+            image = cv2.GaussianBlur(image, self.noise_kernel, 0)
+        grad_x = cv2.Sobel(image, cv2.CV_64F, 1, 0, ksize=ksize)
+        grad_y = cv2.Sobel(image, cv2.CV_64F, 0, 1, ksize=ksize)
+        gradient_magnitude = np.sqrt(grad_x**2 + grad_y**2)
+        gradient_magnitude = cv2.convertScaleAbs(gradient_magnitude)
+        normalized_gradient = cv2.normalize(
+            gradient_magnitude, None, 0, 255, cv2.NORM_MINMAX
+        )
+        cv2.imwrite(str(output_path), normalized_gradient)
     
-    # Save the normalized image
-    cv2.imwrite(str(output_path), normalized_image)
+    def normalize_kernel(
+            self,
+            image_path: Path,
+            output_path: Path,
+            patch_size: int = 30,
+            method: str = "minmax",
+            ):
+        """
+        Normalize an image using kernel method.
+
+        Args:
+            image_path (Path): The path to the input image.
+            output_path (Path): The path to save the normalized image.
+            patch_size (int): The size of the patch. Default is 30.
+            method (str): The normalization method to use. Default is "minmax".
+        """
+        image = cv2.imread(str(image_path), cv2.IMREAD_GRAYSCALE)
+        if image is None:
+            raise ValueError(f"Image at {image_path} could not be read.")
+        if self.reduce_noise:
+            image = cv2.GaussianBlur(image, self.noise_kernel, 0)
+
+        rows, cols = image.shape
+        normalized_image = np.zeros_like(image, dtype=np.float32)
+        half_patch_size = patch_size // 2
+
+        for i in range(rows):
+            for j in range(cols):
+                # Define the boundaries of the patch
+                r_min = max(0, i - half_patch_size)
+                r_max = min(rows, i + half_patch_size + 1)
+                c_min = max(0, j - half_patch_size)
+                c_max = min(cols, j + half_patch_size + 1)
+
+                # Extract the patch
+                patch = image[r_min:r_max, c_min:c_max]
+
+                if method == "minmax":
+                    min_value = np.min(patch)
+                    max_value = np.max(patch)
+                    norm_value = (image[i, j] - min_value) * 255 / (max_value - min_value)
+                    if norm_value > 0 and norm_value < 256:
+                        normalized_image[i, j] = norm_value
+                    else:
+                        normalized_image[i, j] = 0
+
+                elif method == "meanstd":
+                    patch_mean = np.mean(patch)
+                    patch_std = np.std(patch)
+                    # masked_patch = np.ma.masked_outside(patch, patch_mean - 2 * patch_std, patch_mean + 2 * patch_std)
+                    # robust_mean = np.mean(masked_patch)
+                    # robust_std = np.std(masked_patch)
+                    norm_value = (
+                        (image[i, j] - patch_mean - patch_std) * 255 / (2 * patch_std)
+                    )
+                    if norm_value > 0 and norm_value < 256:
+                        normalized_image[i, j] = norm_value
+                    else:
+                        normalized_image[i, j] = 0
+
+        normalized_image = normalized_image.astype(np.uint8)
+        cv2.imwrite(str(output_path), normalized_image)
 
 
 if __name__ == "__main__":
     parser = argparse.ArgumentParser(description="Normalize an image")
     parser.add_argument("images", type=Path, help="Path to image folder to normalize")
-    parser.add_argument("output", type=Path, help="Path to output folder to save normalized images")
-    parser.add_argument("--patch_size", type=int, default=30, help="Size of the patch (default: 30)")
+    parser.add_argument(
+        "output", type=Path, help="Path to output folder to save normalized images"
+    )
+    parser.add_argument(
+        "--patch_size", type=int, default=30, help="Size of the patch (default: 30)"
+    )
+    parser.add_argument(
+        "--reduce_noise", type=bool, default=True, help="Apply noise reduction (default: True)"
+    )
+    parser.add_argument(
+        "--noise_kernel", type=tuple, default=(5, 5), help="Apply noise reduction (default: True)"
+    )
+    parser.add_argument(
+        "--method", type=str, choices=["norm_kernel_minmax", "norm_kernel_meanstd", "gradient"], default=(5, 5), help="Apply noise reduction (default: True)"
+    )
     args = parser.parse_args()
-    
+
     images = os.listdir(args.images)
+    if not os.path.exists(args.output):
+        os.makedirs(args.output)
+
+    image_normalizer = ImageNormalizer(
+        reduce_noise = args.reduce_noise,
+        noise_kernel = args.noise_kernel,
+        )
 
     for img in images:
-        normalize_image(args.images / img, args.output / img, args.patch_size)
\ No newline at end of file
+        if args.method == "norm_kernel_minmax":
+            image_normalizer.normalize_kernel(
+                args.images / img, 
+                args.output / img, 
+                patch_size=args.patch_size,
+                method="minmax",
+                )
+    
+        elif args.method == "norm_kernel_meanstd":
+            image_normalizer.normalize_kernel(
+                args.images / img, 
+                args.output / img, 
+                patch_size=args.patch_size,
+                method="meanstd",
+                )
+        
+        elif args.method == "gradient":
+            image_normalizer.img_gradient(
+                args.images / img, 
+                args.output / img, 
+                ksize=3,
+                )
\ No newline at end of file