Digital Image Processing Course Assignment implemented over Python
Posted on :: Digital Image Processing Python Software

We were given assignment for the Digital Image processing elective at my undergrad in semester 7. Instead of solving these numericals manually, I've written the following script to compute the required image matrices that were asked.

  • Contains histogram equalization for contrast enhancement (1st ques) and DCT, Hadamard transform codes.

  • Finds the global threshold using iterative mean method

  • To run the program, you can pass in the registration number like

    python .\main.py <reg_no>

  • This code needs numpy package for execution.

Code

import sys
import numpy as np
import matplotlib.pyplot as plt
IMG_DIM = (4,4)

def print_dict(d):
    for key, value in d.items():
        print(f"{key} \t {value}")

def hist_equalize(img):
    """
        Takes in an image numpy array and performs histogram equalization
        Returns the histogram equalized image array
        - Calculates the required bit for equalizing by itself
    """
    max_val = np.max(img)
    if max_val>=8: no_of_bits = 4# find the number of bits required to represent
    else: no_of_bits = 3
    max_val = int(2**no_of_bits)
    print(f"Using {int(no_of_bits)} bits for histogram equalization.\nHence the Max value is {max_val-1}.\n")
    img_dict = dict() # For calculating the number of occurences of each pixel values
    for i in range(max_val): # in order to account for pixel values that don't occur in the image
        img_dict[i] = 0
    for i in range(IMG_DIM[0]*IMG_DIM[1]): # update the img_dict with pixel count in given image
        img_dict[img[i//4,i%4]] = np.count_nonzero(img==img[i//4,i%4])
    print(f"Image and its distribution before constrast enhancement: \nPixel value\tCount")
    print_dict(img_dict)
    tab = dict() # Histogram equalization table dictionary
    cumm_sum = 0
    for key, value in img_dict.items():
        tab[key] = [value,0,0] # [occurence of that pixel intensity, cummulative sum, ceil(cummulative_sum*(max-1)/(no_of_pixels))
    for key, value in tab.items():
        cumm_sum+=value[0]
        tab[key][1] = cumm_sum
        tab[key][2] = np.ceil((cumm_sum*(max_val-1))/(IMG_DIM[0]*IMG_DIM[1]))
    print(f"\nHistogram equalization mapping function is:")
    print_dict(tab)
    hist_eq_img = np.zeros(img.shape)
    for key, value in tab.items(): # Mapping r to s (histogram equalization)
        hist_eq_img[img == np.float32(key)] = value[2]
    print(f"The histogram equalized image is: \n{hist_eq_img}")
    return hist_eq_img

def dct(img): # assuming a square matrix img
    """
        Performs Discrete Cosine transform on the input image array. Constructs kernel by itself
        - Assumes the input image array is a square matrix
        - Returns the DCT transform applied image
    """
    def a(u): # coefficient
        if u==0: return (1/IMG_DIM[0])**(0.5)
        else: return (2/IMG_DIM[0])**(0.5)
    kernel = np.zeros(IMG_DIM)
    for u in range(IMG_DIM[0]): # kernel formation
        for x in range(IMG_DIM[1]):
            kernel[u,x] = np.float32(a(u)*np.cos(((2*x+1)*u*np.pi)/(2*IMG_DIM[0])))
    print(f"The DCT kernel formed is: \n{kernel}\n")
    dct_img = np.matmul(kernel,np.matmul(img,np.transpose(kernel))) # applying the transform
    print(f"The DCT transformed matrix for the given image is:\n G = \n{dct_img}")
    return dct_img

def hdt(img):
    """
        Performs Discrete Cosine transform on the input image array, constructs kernel by itself
        - Assumes the input image array is a square matrix
        - Returns the DCT transform applied image
    """
    h0_kernel = (0.5**(0.5))*np.array([[1,1],[1,-1]]) # H1 kernel of hadamard transform
    count = img.shape[0]/h0_kernel.shape[0] # used to find the number of kronekar product operations to be performed to obtain the kernel
    kernel = h0_kernel # intial kernel
    while count:
        kernel = np.kron(h0_kernel,kernel)
        count-=2
    print(f"\nThe Hadamard Transform kernel is:\n{kernel}\n")
    hdt_img = np.matmul(kernel,np.matmul(img,np.transpose(kernel))) # technically you dont need to take the transpose
    print(f"The Hadamard transformed matrix for the given image is:\n G = \n{hdt_img}")
    return hdt_img

def compute_thresh(img):
    """
        Computes the global threshold for an image for segmentation
        - Computes the global threshold by finding the average of mean pixel values of two threshold levels
            - threshold = (u_l+u_g)/2 => threshold is updated over 40 iterations
        - Assumes 40 iterations of threshold value updation
        - Returns the global threshold value
    """
    intial_threshold = np.mean(img)
    print(f"\nInitial threshold: {intial_threshold}")
    for iter in range(40):
        u_l, u_g = 0, 0
        count_l, count_g = 0, 0
        for i in range(img.shape[0]):
            for j in range(img.shape[1]):
                if img[i,j] > intial_threshold:
                    count_g +=1
                    u_g = (u_g*(count_g-1) + img[i,j])/count_g
                else:
                    count_l +=1
                    u_l = (u_l*(count_l-1) + img[i,j])/count_l
            intial_threshold = (u_l+u_g)/2
    print(f"Thresholding value after 40 iterations of mean updation {intial_threshold}")
    return intial_threshold

def image_segment(img, threshold):
    """
        Segments the image for a given threshold
        - Returns a boolean ndarray of the thresholded binary image
    """
    print(f"The segmented image for the given global thresold of {threshold} is\n{1*(img>=threshold)}")
    return 1*(img>=threshold)

if __name__ == "__main__":
    if(len(sys.argv)==1):
        reg_no = input("Enter your registration number: ")
    else:
        reg_no = str(sys.argv[1])
    
    img_no = reg_no+str(2*int(reg_no))
    img_no = img_no[len(img_no)%IMG_DIM[0]:]
    print(f"Your matrix is formed from: {img_no}")
    img = list()
    for i in range(IMG_DIM[0]):
        temp = list()
        for j in range(IMG_DIM[1]):
            temp.append(float(img_no[IMG_DIM[0]*i+j]))
        img.append(temp)
    img = np.array(img)
    print(f"Your image is: \n{np.array2string(img)}")

    # 1. constrast enhancement
    print("\n\nQ1 Histogram equalization<explains the reasons why yourself>\n\n")
    hist_eq_img = hist_equalize(img)

    # 2. Bit plane slicing... do it yourself
    print("\n\nQ2 Bit plane slicing <do it yourself>\n")

    # 3.1 Discrete Cosine Transform
    print("\n\nQ3.1 DCT transform of equalized image <show internal matrices, this just prints the kernel and DCT op>\n")
    dct_img = dct(hist_eq_img)

    # 3.2 Discrete Hadamard Transform
    print("\n\nQ3.2 Hadamard transform of equalized image <show internal matrices, this just prints the kernel and DCT op>\n")
    had_img = hdt(hist_eq_img)

    # 4. Global threshold value and segmentation
    print("\n\nQ4 Thresholding and segmentation of equalized image <Write the theory for thresholding and segmentation (refer slides - global thresholding topic)>\n")
    # Global threshold computation
    threshold = compute_thresh(hist_eq_img)
    
    # Segment the image!
    seg_img = image_segment(hist_eq_img,threshold) # returns an binary numpy array
    
    # 5. Morphological process
    print("\n\nQ5 Morphological process (dilation and erosion, use a 3x3 structing element) <do it yourself>\n")

Thanks to the Code reviewers

  • Dhasharadh
  • Ananthakrishnan
  • Anusha
Table of Contents