diff --git a/Deblur_Algorithm.py b/Deblur_Algorithm.py
new file mode 100644
index 0000000000000000000000000000000000000000..adef6d9e0f3f2cb84925889357a42f7d0bb2348b
--- /dev/null
+++ b/Deblur_Algorithm.py
@@ -0,0 +1,151 @@
+"""
+Deblur Algorithm
+
+input image : F16.tiff
+filters : Mean Filter, Gaussian Filter
+paddings : zero, replicate, mirror
+"""
+
+import cv2
+import numpy as np
+import matplotlib.pyplot as plt
+from scipy.sparse import csr_matrix
+from scipy.sparse.linalg import spsolve
+from datetime import datetime
+from pytz import timezone
+
+# calculate global matrix M in CSR format
+def M(input_data, filter_data, stride):
+ num_elements= n*m*filter_data.size
+ # allocate memory for data, row indices, column indices
+ data=np.zeros(num_elements)
+ row_indices= np.zeros(num_elements, dtype=int)
+ col_indices= np.zeros(num_elements, dtype=int)
+ # initialize counters for data storage
+ count_i = 0 # for row indices
+ count_data=0 # for the current position
+ # operation for convolution
+ for i in range(0, input_data.shape[0] - filter_data.shape[0] + 1, stride):
+ for j in range(0, input_data.shape[1] - filter_data.shape[1] + 1, stride):
+ tmp = np.array(input_data[i:i + filter_data.shape[0], j:j + filter_data.shape[1]]).ravel()
+ tmp2 = np.nonzero(tmp)[0]
+ for k in range(len(tmp2)):
+ data[count_data]= filter_data.ravel()[tmp2[k]]
+ row_indices[count_data]= count_i
+ col_indices[count_data]= tmp[tmp2[k]] - 1
+ count_data +=1
+ count_i += 1
+ data= data[:count_data]
+ row_indices=row_indices[:count_data]
+ col_indices=col_indices[:count_data]
+ result_sparse = csr_matrix((data, (row_indices, col_indices)), shape=(n * m, n * m)) #store in CSR format
+ return result_sparse
+
+# solve MX=Y using spsolve
+def solve(M_sparse, img_blur):
+ try:
+ Y = np.array(img_blur)
+ img_deblur = spsolve(M_sparse, Y.ravel())
+ img_deblur = np.reshape(img_deblur, (n, m))
+ # if inverse matrix doesn't exists, exit the program
+ except (np.linalg.LinAlgError, RuntimeError):
+ img_deblur = np.zeros((n, m))
+ print("Inverse of M does not exist.")
+ exit()
+ return img_deblur
+
+# calculate errors
+def InfNorm(A, B):
+ if A.shape != B.shape:
+ raise ValueError("The shapes of two matrixes must be the same")
+ diff = np.abs(A - B)
+ infnorm = np.max((diff))
+ return infnorm
+
+def apply_filter_and_deblur(img, K):
+ # blur using openCV
+ if border == 0:
+ img_blur = cv2.filter2D(img, -1, K, borderType=cv2.BORDER_CONSTANT)
+ elif border == 1:
+ img_blur = cv2.filter2D(img, -1, K, borderType=cv2.BORDER_REPLICATE)
+ elif border == 2:
+ img_blur = cv2.filter2D(img, -1, K, borderType=cv2.BORDER_REFLECT_101)
+ # calculate M
+ M_sparse = M(border_X, K, 1)
+ # solve MX=Y using spsolve
+ img_deblur = solve(M_sparse, img_blur)
+ # calculate errors between original image and deblurred image
+ error = InfNorm(img128, img_deblur)
+ return img_blur, img_deblur, error
+
+# input image
+img = cv2.imread("images/F16.tiff", cv2.IMREAD_GRAYSCALE)
+img128 = cv2.resize(img, (128, 128))
+img128 = img128.astype(np.float64)/255.0 # image normalization
+
+# show results using matplotlib
+fig = plt.figure()
+rows = 1
+cols = 3
+
+# result 1. original image
+ax1 = fig.add_subplot(rows, cols, 1)
+ax1.imshow(img128, cmap='gray')
+ax1.set_title('F16')
+ax1.axis("off")
+
+n, m = 128, 128
+X = (np.arange(1, n * m + 1)).reshape(n, m)
+
+# choose the padding type
+border_type = ['zero', 'replicate', 'mirror']
+border = int(input(f"{border_type[:]}(Choose from 0-2) : " ))
+if border == 0: # zero padding
+ border_X = np.pad(X, ((1, 1), (1, 1)), mode='constant', constant_values=0)
+elif border == 1: # replicate padding
+ border_X = np.pad(X, ((1, 1), (1, 1)), mode='edge')
+elif border == 2: # mirror padding
+ border_X = np.pad(X, ((1, 1), (1, 1)), mode='reflect')
+
+# choose the filter type
+filter_type = ['Mean', 'Gaussian']
+kernel_size = 3
+kernel = int(input(f"{filter_type[:]}(Choose from 0-1) : " ))
+if kernel == 0: # Mean filter
+ K = np.ones((kernel_size,kernel_size),np.float64)/(kernel_size*kernel_size)
+ img_blur, img_deblur, min_error = apply_filter_and_deblur(img128, K)
+elif kernel == 1: # Gaussian filter
+ min_error = float('inf')
+ best_sigma = None
+ sigma_values = np.arange(1, 101, 1)
+ # find best sigma values
+ for sigma in sigma_values:
+ tmp3 = cv2.getGaussianKernel(kernel_size, sigma)
+ K = np.outer(tmp3, tmp3.T)
+ img_blur, img_deblur, error = apply_filter_and_deblur(img128, K)
+ if error < min_error:
+ min_error = error
+ best_sigma = sigma
+ best_deblur = img_deblur
+ img_deblur = best_deblur
+ print("Best sigma value: ", best_sigma)
+
+# result 2. blurred image using openCV
+ax2 = fig.add_subplot(rows, cols, 2)
+ax2.imshow(img_blur, cmap='gray')
+ax2.set_title(f'{filter_type[kernel]}')
+ax2.axis("off")
+
+#result 3. deblurred image using M()
+ax3 = fig.add_subplot(rows, cols, 3)
+ax3.imshow(img_deblur, cmap='gray')
+ax3.set_title('Deblur')
+ax3.axis("off")
+
+# print error
+print("Error:", min_error)
+
+# save results
+date_today = datetime.now(timezone('Asia/Seoul')).strftime("%Y%m%d_%H%M")
+plt.savefig('Results/{}.jpg'.format(date_today), format='jpg')
+plt.close()
\ No newline at end of file
diff --git a/F16.tiff b/F16.tiff
new file mode 100644
index 0000000000000000000000000000000000000000..8d0769987d66cc7c66d8e157d26403eca57b4b62
Binary files /dev/null and b/F16.tiff differ