proj1_tif_blurry code

import numpy as np
import cv2
import matplotlib.pyplot as plt

# Function to manually apply a Gaussian filter
def apply_gaussian(image, kernel_size=5, sigma=1.0):
    kernel = np.fromfunction(lambda x, y: (1/2*np.pi*sigma**2) * np.exp(-((x-(kernel_size-1)//2)**2 + (y-(kernel_size-1)//2)**2) / (2*sigma**2)), (kernel_size, kernel_size))
    kernel /= np.sum(kernel)
    
    # Manually apply convolution using the Gaussian kernel
    image_padded = np.pad(image, (kernel_size // 2,), mode='reflect')
    height, width = image.shape
    filtered_image = np.zeros_like(image)
    
    for i in range(height):
        for j in range(width):
            filtered_image[i, j] = np.sum(kernel * image_padded[i:i+kernel_size, j:j+kernel_size])
    
    return filtered_image

# Function to downsample (reduce size by a factor of 2)
def downsample(image):
    return image[::2, ::2]

# Function to upsample (increase size by a factor of 2)
def upsample(image):
    upsampled = np.zeros((image.shape[0]*2, image.shape[1]*2), dtype=image.dtype)
    upsampled[::2, ::2] = image
    return upsampled

# Function to crop the inner 90% of the image
def crop_inner_90(image):
    height, width = image.shape
    h_crop = int(height * 0.075)  # Crop 7.5% on top and bottom
    w_crop = int(width * 0.075)   # Crop 7.5% on left and right
    return image[h_crop:height - h_crop, w_crop:width - w_crop]

# Function to create Gaussian pyramid
def gaussian_pyramid(image, levels):
    pyramid = [image]
    for i in range(1, levels):
        image = downsample(apply_gaussian(image))
        pyramid.append(image)
    return pyramid

# Function to create Laplacian pyramid
def laplacian_pyramid(gaussian_pyramid):
    laplacian_pyramid = []
    for i in range(len(gaussian_pyramid) - 1):
        upsampled = upsample(gaussian_pyramid[i+1])
        # Ensure the sizes match after upsampling (sometimes slight mismatches)
        upsampled = cv2.resize(upsampled, (gaussian_pyramid[i].shape[1], gaussian_pyramid[i].shape[0]))
        laplacian = gaussian_pyramid[i] - upsampled
        laplacian_pyramid.append(laplacian)
    laplacian_pyramid.append(gaussian_pyramid[-1])  # The last level is kept as is
    return laplacian_pyramid

# Function to compute L2 norm (SSD) with inner 90% of the image
def compute_l2_norm(image1, image2):
    height = min(image1.shape[0], image2.shape[0])
    width = min(image1.shape[1], image2.shape[1])
    image1_cropped = image1[:height, :width]
    image2_cropped = image2[:height, :width]
    image1_inner = crop_inner_90(image1_cropped)
    image2_inner = crop_inner_90(image2_cropped)
    return np.sqrt(np.sum((image1_inner - image2_inner) ** 2))

# Function to compute NCC (Normalized Cross-Correlation) with inner 90% of the image
def compute_ncc(image1, image2):
    height = min(image1.shape[0], image2.shape[0])
    width = min(image1.shape[1], image2.shape[1])
    image1_cropped = image1[:height, :width]
    image2_cropped = image2[:height, :width]
    image1_inner = crop_inner_90(image1_cropped)
    image2_inner = crop_inner_90(image2_cropped)
    image1_norm = (image1_inner - np.mean(image1_inner)) / np.std(image1_inner)
    image2_norm = (image2_inner - np.mean(image2_inner)) / np.std(image2_inner)
    return np.sum(image1_norm * image2_norm) / (np.sqrt(np.sum(image1_norm ** 2)) * np.sqrt(np.sum(image2_norm ** 2)))

# Function to shift an image by a given offset
def shift_image(image, dx, dy):
    return np.roll(np.roll(image, dx, axis=1), dy, axis=0)

# Function to manually align images using either L2 or NCC metric
def align_images(reference_image, target_image, metric='l2', displacement_range=15):
    best_offset = (0, 0)
    best_score = float('inf') if metric == 'l2' else -1  # Initialize best score

    # Search over all possible displacements
    for dx in range(-displacement_range, displacement_range + 1):
        for dy in range(-displacement_range, displacement_range + 1):
            # Shift the target image by dx and dy
            shifted_target = np.roll(np.roll(target_image, dx, axis=1), dy, axis=0)

            # Calculate the score
            if metric == 'l2':
                score = compute_l2_norm(reference_image, shifted_target)
                if score < best_score:
                    best_score = score
                    best_offset = (dx, dy)
            elif metric == 'ncc':
                score = compute_ncc(reference_image, shifted_target)
                if score > best_score:  # Maximize NCC
                    best_score = score
                    best_offset = (dx, dy)

    return best_offset

# Function to process the image (input is a .tif, output will be a .jpg)
def process_image(image_path, metric='l2', displacement_range=15, pyramid_levels=3):
    # Load the .tif image
    image = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    if image is None:
        raise ValueError(f"Error loading image {image_path}. Ensure it's a valid image file.")

    height = image.shape[0] // 3
    blue = image[:height, :]
    green = image[height:2 * height, :]
    red = image[2 * height:, :]

    # Build Gaussian pyramids for each channel
    blue_pyramid = gaussian_pyramid(blue, pyramid_levels)
    green_pyramid = gaussian_pyramid(green, pyramid_levels)
    red_pyramid = gaussian_pyramid(red, pyramid_levels)

    # Align the Green and Red channels to the Blue channel at the top of the pyramid
    offset_green = align_images(blue_pyramid[-1], green_pyramid[-1], metric=metric, displacement_range=displacement_range)
    offset_red = align_images(blue_pyramid[-1], red_pyramid[-1], metric=metric, displacement_range=displacement_range)

    # Shift the Green and Red channels based on the calculated offsets
    green_aligned = shift_image(green, offset_green[0], offset_green[1])
    red_aligned = shift_image(red, offset_red[0], offset_red[1])

    # Crop all channels to the same size
    min_height = min(blue.shape[0], green_aligned.shape[0], red_aligned.shape[0])
    min_width = min(blue.shape[1], green_aligned.shape[1], red_aligned.shape[1])
    blue_cropped = blue[:min_height, :min_width]
    green_cropped = green_aligned[:min_height, :min_width]
    red_cropped = red_aligned[:min_height, :min_width]

    # Merge the aligned channels into a color image
    aligned_image = cv2.merge([blue_cropped, green_cropped, red_cropped])

    # Display the aligned image
    plt.imshow(cv2.cvtColor(aligned_image, cv2.COLOR_BGR2RGB))
    plt.title('Aligned Color Image')
    plt.axis('off')
    plt.show()

    # Save the aligned image as .jpg
    output_path = image_path.replace('.tif', '_aligned.jpg')
    cv2.imwrite(output_path, aligned_image)

# Example usage
image_path = 'emir.tif'  # Provide the path to your .tif file
process_image(image_path, metric='ncc', displacement_range=15, pyramid_levels=5)