src/methods/Elmehdi_lahmar/malvar.py

+import numpy as np
+from scipy.signal import convolve2d
+from src.forward_model import CFA
+
+def malvar(y: np.ndarray, op: CFA) -> np.ndarray:
+    """
+    Malvar-He-Cutler demosaicing algorithm for Bayer pattern.
+    """
+    # Convert the mosaicked image to the initial estimated channels.
+    z = op.adjoint(y)
+
+    # Define convolution kernels
+    kernel_G_at_RB = np.array([[0, 0, -1, 0, 0], [0, 0, 2, 0, 0], [-1, 2, 4, 2, -1], [0, 0, 2, 0, 0], [0, 0, -1, 0, 0]]) / 8
+    kernel_RB_at_G = np.array([[0, 0, 0.5, 0, 0], [0, -1, 0, -1, 0], [-1, 4, 5, 4, -1], [0, -1, 0, -1, 0], [0, 0, 0.5, 0, 0]]) / 8
+    kernel_RB_at_RB = np.array([[0, 0, -1.5, 0, 0], [0, 2, 0, 2, 0], [-1.5, 0, 6, 0, -1.5], [0, 2, 0, 2, 0], [0, 0, -1.5, 0, 0]]) / 8
+
+    # Interpolate each channel
+    R = z[:, :, 0]
+    G = z[:, :, 1]
+    B = z[:, :, 2]
+
+    mask = op.mask
+    R_m, G_m, B_m = mask[:, :, 0], mask[:, :, 1], mask[:, :, 2]
+
+    # Interpolate G at R and B locations
+    G = np.where(np.logical_or(R_m == 1, B_m == 1), convolve2d(y, kernel_G_at_RB, mode='same'), G)
+
+    # Interpolate R at G and B locations, B at R and G locations
+    R = np.where(np.logical_or(G_m == 1, B_m == 1), convolve2d(y, kernel_RB_at_G, mode='same'), R)
+    B = np.where(np.logical_or(R_m == 1, G_m == 1), convolve2d(y, kernel_RB_at_G, mode='same'), B)
+
+    # Interpolate R at B locations and B at R locations
+    R = np.where(B_m == 1, convolve2d(y, kernel_RB_at_RB, mode='same'), R)
+    B = np.where(R_m == 1, convolve2d(y, kernel_RB_at_RB, mode='same'), B)
+
+    # Combine channels
+    return np.clip(np.stack((R, G, B), axis=-1), 0, 1)
+
+
+
+
+
+####
+####
+####
+
+####      ####                ####        #############
+####      ######              ####      ##################
+####      ########            ####      ####################
+####      ##########          ####      ####        ########
+####      ############        ####      ####            ####
+####      ####  ########      ####      ####            ####
+####      ####    ########    ####      ####            ####
+####      ####      ########  ####      ####            ####
+####      ####  ##    ######  ####      ####          ######
+####      ####  ####      ##  ####      ####    ############
+####      ####  ######        ####      ####    ##########
+####      ####  ##########    ####      ####    ########
+####      ####      ########  ####      ####
+####      ####        ############      ####
+####      ####          ##########      ####
+####      ####            ########      ####
+####      ####              ######      ####
+
+# 2023
+# Authors: Mauro Dalla Mura and Matthieu Muller
\ No newline at end of file

src/methods/Elmehdi_lahmar/reconstruct.py

+"""The main file for the reconstruction.
+This file should NOT be modified except the body of the 'run_reconstruction' function.
+Students can call their functions (declared in others files of src/methods/your_name).
+"""
+
+
+import numpy as np
+from src.methods.Elmehdi_lahmar.malvar import malvar
+from src.forward_model import CFA
+
+def run_reconstruction(y: np.ndarray, cfa: str) -> np.ndarray:
+    """
+    Run the demosaicing process using Malvar-He-Cutler algorithm for Bayer pattern.
+
+    Args:
+    y (np.ndarray): The mosaicked image to be reconstructed.
+    cfa (str): Name of the CFA, expected to be 'bayer'.
+
+    Returns:
+    np.ndarray: The demosaicked image.
+    """
+    if cfa != 'bayer':
+        raise ValueError("Malvar-He-Cutler demosaicing only supports Bayer CFA pattern.")
+
+    input_shape = (y.shape[0], y.shape[1], 3)
+    op = CFA(cfa, input_shape)
+
+    return malvar(y, op)
+
+
+
+
+
+
+####
+####
+####
+
+####      ####                ####        #############
+####      ######              ####      ##################
+####      ########            ####      ####################
+####      ##########          ####      ####        ########
+####      ############        ####      ####            ####
+####      ####  ########      ####      ####            ####
+####      ####    ########    ####      ####            ####
+####      ####      ########  ####      ####            ####
+####      ####  ##    ######  ####      ####          ######
+####      ####  ####      ##  ####      ####    ############
+####      ####  ######        ####      ####    ##########
+####      ####  ##########    ####      ####    ########
+####      ####      ########  ####      ####
+####      ####        ############      ####
+####      ####          ##########      ####
+####      ####            ########      ####
+####      ####              ######      ####
+
+# 2023
+# Authors: Mauro Dalla Mura and Matthieu Muller

src/methods/GERAYELI_Kourosh/GERAYELI_Kourosh.py

+
+# Importing libraries
+import os
+import colour
+from colour_demosaicing import (
+    demosaicing_CFA_Bayer_bilinear,
+    demosaicing_CFA_Bayer_Malvar2004,
+    demosaicing_CFA_Bayer_Menon2007,
+    mosaicing_CFA_Bayer)
+
+from src.utils import psnr,ssim
+
+# Image path
+image_pathes = ['images/img_1.png','images/img_2.png','images/img_3.png','images/img_4.png']
+
+for i in image_pathes:
+    LIGHTHOUSE_IMAGE = colour.io.read_image(i)
+    img = LIGHTHOUSE_IMAGE
+    colour.plotting.plot_image(
+        colour.cctf_encoding(LIGHTHOUSE_IMAGE))
+
+
+    # Mosaicing
+    CFA = mosaicing_CFA_Bayer(LIGHTHOUSE_IMAGE)
+
+    colour.plotting.plot_image(
+        colour.cctf_encoding(CFA),
+        text_kwargs={'text': 'Lighthouse - CFA - RGGB'})
+
+    colour.plotting.plot_image(
+        colour.cctf_encoding(mosaicing_CFA_Bayer(LIGHTHOUSE_IMAGE, 'BGGR')), 
+        text_kwargs={'text': 'Lighthouse - CFA - BGGR'});
+
+
+    # Demosaicing bilinear
+    colour.plotting.plot_image(
+        colour.cctf_encoding(demosaicing_CFA_Bayer_bilinear(CFA)), 
+        text_kwargs={'text': 'Demosaicing - Bilinear'});
+
+    recons_bilinear = colour.cctf_encoding(demosaicing_CFA_Bayer_bilinear(CFA))
+
+
+    # demosaicing Malvar
+    recons_malvar = colour.cctf_encoding(demosaicing_CFA_Bayer_Malvar2004(CFA))
+    colour.plotting.plot_image(
+        colour.cctf_encoding(demosaicing_CFA_Bayer_Malvar2004(CFA)), 
+        text_kwargs={'text': 'Demosaicing - Malvar (2004)'});
+
+    # demosaicing Menon
+    recons_menon = colour.cctf_encoding(demosaicing_CFA_Bayer_Menon2007(CFA))
+    colour.plotting.plot_image(
+        colour.cctf_encoding(demosaicing_CFA_Bayer_Menon2007(CFA)), 
+        text_kwargs={'text': 'Demosaicing - Menon (2007)'});
+
+
+    print('bilinear : ')
+    print(f'PSNR: {psnr(img, recons_bilinear):.2f}')
+    print(f'SSIM: {ssim(img, recons_bilinear):.4f}')
+
+    print('Malvar : ')
+    print(f'PSNR: {psnr(img, recons_malvar):.2f}')
+    print(f'SSIM: {ssim(img, recons_malvar):.4f}')
+
+    print('Menon')
+    print(f'PSNR: {psnr(img, recons_menon):.2f}')
+    print(f'SSIM: {ssim(img, recons_menon):.4f}')
\ No newline at end of file

src/methods/Mirabito Jules/reconstruct.py

+"""The main file for the reconstruction.
+This file should NOT be modified except the body of the 'run_reconstruction' function.
+Students can call their functions (declared in others files of src/methods/your_name).
+"""
+
+
+import numpy as np
+
+
+from colour_demosaicing import (
+    demosaicing_CFA_Bayer_bilinear,
+    demosaicing_CFA_Bayer_Malvar2004,
+    demosaicing_CFA_Bayer_Menon2007,
+)
+
+
+def run_reconstruction(y: np.ndarray, cfa: str) -> np.ndarray:
+    """Performs demosaicking on y.
+
+    Args:
+        y (np.ndarray): Mosaicked image to be reconstructed.
+        cfa (str): Name of the CFA. Can be bayer or quad_bayer.
+
+    Returns:
+        np.ndarray: Demosaicked image.
+    """
+
+    #The code implemented here comes from an open-source python library that has optimized the two methods presented.
+    #This library is based on the "color" library, specific to color image processing.
+    #The source code for the following functions can be found in the "demosaicing" folder in the following git hub: https://github.com/colour-science/colour-demosaicing
+    #Copyright 2015 Colour Developers -
+
+    if cfa == 'bayer' : 
+        #res = demosaicing_CFA_Bayer_bilinear(y,'GRBG')
+        res = demosaicing_CFA_Bayer_Menon2007(y,'GRBG')
+        #res = demosaicing_CFA_Bayer_Malvar2004(y,'GRBG')
+    else : 
+        print('Error - Not implemented')
+        return
+
+    return res
+
+####
+####
+####
+
+####      ####                ####        #############
+####      ######              ####      ##################
+####      ########            ####      ####################
+####      ##########          ####      ####        ########
+####      ############        ####      ####            ####
+####      ####  ########      ####      ####            ####
+####      ####    ########    ####      ####            ####
+####      ####      ########  ####      ####            ####
+####      ####  ##    ######  ####      ####          ######
+####      ####  ####      ##  ####      ####    ############
+####      ####  ######        ####      ####    ##########
+####      ####  ##########    ####      ####    ########
+####      ####      ########  ####      ####
+####      ####        ############      ####
+####      ####          ##########      ####
+####      ####            ########      ####
+####      ####              ######      ####
+
+# 2023
+# Authors: Mauro Dalla Mura and Matthieu Muller

src/methods/RHOUCH_Oussama/clip_values.py

+import numpy as np
+
+def clip_values(img: np.ndarray) -> np.ndarray:
+    """Clip values of an image between 0 and 1
+    
+    Args:
+        img (np.ndarray): image to clip
+        
+    Returns:
+        np.ndarray: clipped image
+    """
+    # Clip values between 0 and 1
+    img[img < 0] = 0
+    img[img > 1] = 1 
+    return img

src/methods/RHOUCH_Oussama/find_direct_neighbors.py

+import numpy as np  
+
+def find_direct_neighbors(neighbors: list) -> list:
+    """
+    Calculate the gradients in horizontal, vertical, and diagonal directions for a pixel
+    based on its neighboring pixels.
+
+    Args:
+        neighbors (list): The neighbors of a pixel in a 3x3 grid. The order of neighbors
+                          is assumed to be [top-left, top, top-right, right, center, left,
+                          bottom-left, bottom, bottom-right].
+
+    Returns:
+        list: The gradients [horizontal, vertical, diagonal-x, diagonal-y].
+    """
+    
+    # Horizontal gradient: Difference between the right and left neighbors
+    horizontal_gradient = (neighbors[3] - neighbors[5]) / 2
+
+    # Vertical gradient: Difference between the top and bottom neighbors
+    vertical_gradient = (neighbors[1] - neighbors[7]) / 2
+
+    # Diagonal gradient along x-axis: Difference between top-right and bottom-left neighbors
+    diagonal_gradient_x = (neighbors[2] - neighbors[6]) / (2 * np.sqrt(2))
+
+    # Diagonal gradient along y-axis: Difference between top-left and bottom-right neighbors
+    diagonal_gradient_y = (neighbors[0] - neighbors[8]) / (2 * np.sqrt(2))
+    
+    return [horizontal_gradient, vertical_gradient, diagonal_gradient_x, diagonal_gradient_y]

src/methods/RHOUCH_Oussama/find_neighbors.py

+import numpy as np
+
+def find_neighbors(img: np.ndarray, channel: int, i: int, j: int, N: int, M: int) -> list:
+    """
+    Retrieve the 3x3 neighborhood of a pixel in a specified channel of an image.
+    The function wraps around the image edges to handle border pixels.
+
+    Args:
+        img (np.ndarray): The image to process.
+        channel (int): The index of the channel to process.
+        i (int): The row index of the pixel.
+        j (int): The column index of the pixel.
+        N (int): Height of the image.
+        M (int): Width of the image.
+        
+    Returns:
+        np.ndarray: An array of neighbor pixel values in the specified channel.
+    """
+    
+    neighbors = []
+
+    for di in [-1, 0, 1]:
+        for dj in [-1, 0, 1]:
+            neighbor_i = (i + di) % N
+            neighbor_j = (j + dj) % M
+
+            neighbors.append(img[neighbor_i, neighbor_j, channel])
+
+    return np.array(neighbors)

src/methods/RHOUCH_Oussama/find_weights.py

+import numpy as np
+from .find_neighbors import find_neighbors
+from .find_direct_neighbors import find_direct_neighbors
+
+def find_weights(img: np.ndarray, direct_neighbors: list, channel: int, i: int, j: int, N: int, M: int) -> list:
+    """
+    Find the weights of the neighbors of a pixel in the image.
+    
+    Args:
+        img (np.ndarray): The image to process.
+        direct_neighbors (list): The list of direct neighbors of the pixel.
+        channel (int): The index of the channel to process.
+        i (int): The row index of the pixel.
+        j (int): The column index of the pixel.
+        N (int): Height of the image.
+        M (int): Width of the image.
+        
+    Returns:
+        list: The list of weights of the neighbors of the pixel.
+    """
+    
+    [Dx, Dy, Dxx, Dyy] = direct_neighbors
+    E = []
+    c = 1
+    
+    for k in [-1, 0, 1]:
+        for l in [-1, 0, 1]:
+            n = find_neighbors(img, channel, i + k, j + l, N, M)
+            dd = find_direct_neighbors(n)
+
+            sqrt_arguments = {
+                1: 1 + Dyy * 2 + dd[3] * 2,
+                3: 1 + Dxx * 2 + dd[2] * 2,
+                2: 1 + Dy * 2 + dd[1] * 2,
+                4: 1 + Dx * 2 + dd[0] * 2
+            }
+
+            value = sqrt_arguments.get(c, 1)
+            if value < 0:
+                E.append(0)
+            else:
+                E.append(1 / np.sqrt(value))
+
+            c += 1
+            
+    return E

src/methods/RHOUCH_Oussama/interpolate.py

+def interpolate(neighbors: list, weights: list) -> float:
+    """
+    Interpolate the pixel value.
+    
+    Args:
+        neighbors (list): The neighbors of the pixel.
+        weights (list): The weights of the direct neighbors.
+        
+    Returns:
+        float: The interpolated value.
+    """
+    
+    return (weights[1] * neighbors[1] + weights[3] * neighbors[3] + weights[4] * neighbors[5] + weights[6] * neighbors[7]) / (weights[1] + weights[1] + weights[4] + weights[6])
+  
\ No newline at end of file

src/methods/RHOUCH_Oussama/interpolate_neighboring.py

+def interpolate_neighboring(neighbors_0: list, neighbors_1: list, weights: list) -> float:
+    """
+    Interpolate the neighboring pixels.
+    
+    Args:
+        neighbors_0 (list): The neighbors of the pixel in the channel to process.
+        neighbors_1 (list): The neighbors of the pixel in the green channel.
+        weights (list): The weights of the direct neighbors.
+        
+    Returns:
+        float: The interpolated value.
+    """
+    
+    return neighbors_1[4] * (((weights[0] * neighbors_0[0]) / neighbors_1[0]) + ((weights[2] * neighbors_0[2]) / neighbors_1[2]) + ((weights[5] * neighbors_0[6]) / neighbors_1[6]) + ((weights[7] * neighbors_0[8]) / neighbors_1[8])) / (weights[0] + weights[2] + weights[5] + weights[7])
+    
\ No newline at end of file

src/methods/RHOUCH_Oussama/main.ipynb

+%% Cell type:markdown id: tags:
+
+## Main notebook for experimenting with the algorithms
+
+Here is a simple notebook to test your code. You can modify it as you please.
+
+Remember to restart the jupyter kernel each time you modify a file.
+
+%% Cell type:code id: tags:
+
+``` python
+%%capture
+
+%pip install -r requirements.txt
+```
+
+%% Cell type:code id: tags:
+
+``` python
+import matplotlib.pyplot as plt
+
+from src.utils import load_image, save_image, psnr, ssim
+from src.forward_model import CFA
+from src.methods.RHOUCH_Oussama.reconstruct import run_reconstruction
+```
+
+%% Cell type:markdown id: tags:
+
+### Load the input image
+
+%% Cell type:code id: tags:
+
+``` python
+image_path = 'images/img_4.png'
+
+img = load_image(image_path)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Shows some information on the image
+plt.imshow(img)
+plt.show()
+print(f'Shape of the image: {img.shape}.')
+```
+
+%% Cell type:markdown id: tags:
+
+### Definition of the forward model
+
+To setup the forward operator we just need to instanciate the `CFA` class. This class needs two arguments: `cfa_name` being the kind of pattern (bayer or kodak), and `input_shape` the shape of the inputs of the operator.
+
+This operation is linear and can be represented by a matrix $A$ but no direct access to this matrix is given (one can create it if needed). However the method `direct` allows to perform $A$'s operation. Likewise the method `adjoint` will perform the operation of $A^T$.
+
+For example let $X \in \mathbb R^{M \times N \times 3}$ the input RGB image in natural shape. Then we got $x \in \mathbb R^{3MN}$ (vectorized version of $X$) and $A \in \mathbb R^{MN \times 3MN}$, leading to:
+
+\begin{equation*}
+    y = Ax \in \mathbb R^{MN} \quad \text{and} \quad z = A^Ty  \in \mathbb R^{3MN}
+\end{equation*}
+
+However thanks to `direct` and `adjoint` there is no need to work with vectorized images, except if it is interesting to create the matrix $A$ explicitly.
+
+%% Cell type:code id: tags:
+
+``` python
+cfa_name = 'bayer' # bayer or quad_bayer
+op = CFA(cfa_name, img.shape)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Shows the mask
+plt.imshow(op.mask[:10, :10])
+plt.show()
+print(f'Shape of the mask: {op.mask.shape}.')
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Applies the mask to the image
+y = op.direct(img)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Applies the adjoint operation to y
+z = op.adjoint(y)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+fig, axs = plt.subplots(2, 3, figsize=(15, 10))
+axs[0, 0].imshow(img)
+axs[0, 0].set_title('Input image')
+axs[0, 1].imshow(y, cmap='gray')
+axs[0, 1].set_title('Output image')
+axs[0, 2].imshow(z)
+axs[0, 2].set_title('Adjoint image')
+axs[1, 0].imshow(img[800:864, 450:514])
+axs[1, 0].set_title('Zoomed input image')
+axs[1, 1].imshow(y[800:864, 450:514], cmap='gray')
+axs[1, 1].set_title('Zoomed output image')
+axs[1, 2].imshow(z[800:864, 450:514])
+axs[1, 2].set_title('Zoomed adjoint image')
+plt.show()
+```
+
+%% Cell type:markdown id: tags:
+
+### Run the reconstruction
+
+Here the goal is to reconstruct the image `img` using only `y` and `op` (using `img` is forbidden).
+
+To run the reconstruction we simply call the function `run_reconstruction`. This function takes in argument the image to reconstruct and the kind of CFA used (bayer or kodak). All the parameters related to the reconstruction itself must be written inside `run_reconstruction`.
+
+%% Cell type:code id: tags:
+
+``` python
+res = run_reconstruction(y, cfa_name)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Prints some information on the reconstruction
+print(f'Size of the reconstruction: {res.shape}.')
+```
+
+%% Cell type:markdown id: tags:
+
+### Quantitative and qualitative results
+
+%% Cell type:code id: tags:
+
+``` python
+fig, axs = plt.subplots(1, 2, figsize=(15, 15))
+axs[0].imshow(img)
+axs[0].set_title('Original image')
+axs[1].imshow(res)
+axs[1].set_title('Reconstructed image')
+plt.show()
+```
+
+%% Cell type:code id: tags:
+
+``` python
+# Computes some metrics
+print(f'PSNR: {psnr(img, res):.2f}')
+print(f'SSIM: {ssim(img, res):.4f}')
+```
+
+%% Cell type:code id: tags:
+
+``` python
+reconstructed_path = 'output/reconstructed_image.png'
+
+save_image(reconstructed_path, res)
+```
+
+%% Cell type:code id: tags:
+
+``` python
+```

src/methods/RHOUCH_Oussama/process_blue_channel.py

+import numpy as np
+from .find_neighbors import find_neighbors
+from .find_direct_neighbors import find_direct_neighbors
+from .find_weights import find_weights
+from .interpolate import interpolate
+from .interpolate_neighboring import interpolate_neighboring
+from .clip_values import clip_values
+from .process_channel import process_channel
+from .process_interpolation import process_interpolation
+
+def process_blue_channel(img: np.ndarray, N: int, M: int) -> np.ndarray:
+    """
+    Process the blue channel of an image
+    
+    Args:
+        img (np.ndarray): image to process
+        N (int): height of the image
+        M (int): width of the image
+        
+    Returns:
+        np.ndarray: processed blue channel of the image
+    """
+    
+    img = process_interpolation(img, 2, N, M)
+    img = process_channel(img, 2, N, M)
+    img = clip_values(img)
+    
+    return img
\ No newline at end of file

src/methods/RHOUCH_Oussama/process_channel.py

+import numpy as np
+from .find_neighbors import find_neighbors
+from .find_direct_neighbors import find_direct_neighbors
+from .find_weights import find_weights
+from .interpolate import interpolate
+
+def process_channel(img: np.ndarray, channel_index: int, N: int, M: int) -> np.ndarray:
+    """
+    Generic function to process a specific channel in the image.
+
+    Args:
+        img (np.ndarray): The image to process.
+        channel_index (int): The index of the channel to process.
+        N (int): Height of the image.
+        M (int): Width of the image.
+
+    Returns:
+        np.ndarray: The processed image.
+    """
+
+    for i in range(N):
+        for j in range(M):
+            if(img[i, j, channel_index] == 0):
+                neighbors = find_neighbors(img, channel_index, i, j, N, M)
+                direct_neighbors = find_direct_neighbors(neighbors)
+                weights = find_weights(img, direct_neighbors, channel_index, i, j, N, M)
+                img[i, j, channel_index] = interpolate(neighbors, weights)
+    
+    return img
+                

src/methods/RHOUCH_Oussama/process_green_channel.py

+import numpy as np
+from .clip_values import clip_values
+from .process_channel import process_channel
+
+def process_green_channel(img: np.ndarray, N: int, M: int) -> np.ndarray:
+    """Process the green channel of an image
+    
+    Args:
+        img (np.ndarray): image to process
+        N (int): height of the image
+        M (int): width of the image
+        
+    Returns:
+        np.ndarray: processed green channel of the image
+    """
+    
+    img = process_channel(img, 1, N, M)                
+    img = clip_values(img)
+    
+    return img
\ No newline at end of file