Spatial High-pass Filters.

Posted on Edited on In Views: Waline: Word count in article: 11k Reading time ≈ 10 mins.

This is an automatically translated post by LLM. The original post is in Chinese. If you find any translation errors, please leave a comment to help me improve the translation. Thanks!

High-pass filters are often used to enhance images and extract edge information. They have a wide range of applications in image processing and image recognition. The main types of high-pass filters are as follows: unsharp mask, Sobel operator, Laplacian operator, and Canny algorithm.

Unsharp Mask

This algorithm is often used for image enhancement and is not commonly used for edge extraction. The implementation method is as follows:

  • For the original image \(f(x,y)\), apply Gaussian blur to obtain the smoothed image \(\overline{f(x,y)}\).
  • Take the difference between the original image and the smoothed image to obtain the mask \(g_{mask}(x,y)=f(x,y)-\overline{f(x,y)}\).
  • Add the original image and k times the mask to obtain the enhanced image \(g(x,y)=f(x,y)+k*g_{mask}(x,y)\).
  • Generally, the value of k is set to 1 for image enhancement. If k is greater than 1, it is a high-boost image.

The image processing effect of this algorithm is as follows:

From left to right: original image, obtained mask (normalized to 0-255), enhanced image

From the results in the above figure, it can be seen that the mask image has high brightness at the edges and low brightness elsewhere. After overlaying with the original image, the brightness of the original image's edges becomes larger. In this case, the method used here is to normalize the overlaid image as a whole and then multiply it by 255 to obtain an integer value. Because the brightness of the original image's edges is large after overlaying with the mask, the bright parts in the original image become darker after normalization. If another method is used, such as treating values greater than 255 as 255 and values less than 0 as 0, it can effectively solve this problem, but the edges of the image will not be as obvious as the previous method. The specific method to be used depends on the actual situation.

Sobel Operator

The Sobel operator is often used for edge extraction in images. It uses the first-order derivative of the image to extract edges, which is denoted as \(\nabla f\).

We know that the gradient of a binary function consists of the directional derivatives in the x and y directions. The Sobel operator is the same, with two operators for the x and y directions, respectively, as follows: \[ \begin{bmatrix}-1&0&1\\-1&0&2\\-1&0&1\end{bmatrix} ~~~~~~~~~~\begin{bmatrix}-1&-2&-1\\0&0&0\\1&2&1\end{bmatrix} \]

  • Convolve the image with the above two operators separately to obtain the x-direction and y-direction edges \(G_x\) and \(G_y\).
  • The first-order derivative of the image \(G\) can be calculated using \(G_x\) and \(G_y\), which is \(G=\sqrt{G_x^2+G_y^2}\). Sometimes, for fast calculation, \(G=|G_x|+|G_y|\) or \(G=max\{|G_x|,|G_y|\}\) can be used.
  • Normalize G to the range of 0-255 to obtain the Sobel edges of the image.

From left to right: original image, x-direction edges, y-direction edges, Sobel edges

The Sobel operator performs well in extracting obvious edges in images and also performs well in extracting detailed edges in images. This can be illustrated by the processing results of another image:

It can be seen that even the small edge contours on the chimney are well represented by the edges extracted by the Sobel operator.

Laplacian Operator

\[ [f(x+1)-f(x)]-[f(x)-f(x-1)]=f(x+1)+f(x-1)-2f(x) \] From this, the Laplacian operator can be derived as: \[ \begin{bmatrix}0&1&0\\1&-4&1\\0&1&0\end{bmatrix} \] Sometimes, the second-order derivative in the diagonal direction is also added, and the operator becomes: \[ \begin{bmatrix}1&1&1\\1&-8&1\\1&1&1\end{bmatrix} \] In practice, the following two operators are often used: \[ \begin{bmatrix}0&-1&0\\-1&4&-1\\0&-1&0\end{bmatrix} \begin{bmatrix}-1&-1&-1\\-1&8&-1\\-1&-1&-1\end{bmatrix} \]

Assuming that the extracted image edges are \(L(x,y)\), the algorithm for image enhancement is \(g(x,y)=f(x,y)+c*L(x,y)\)

If the above two operators are used to extract the edges, c=-1; if the following two operators are used, c=1.

The image edge information extracted using the Laplacian operator is as follows:

From left to right: original image, Laplacian edges without diagonal, Laplacian edges with diagonal

The edges extracted by the Laplacian operator have two edge lines at each edge, which is determined by the properties of its second-order derivative. Compared with the edges extracted by the Sobel operator, the Laplacian edges are more detailed and capture the "edges of edges". This can be seen more clearly from the comparison in the following figure:

The left side is the result of the Sobel operator, and the right side is the result of the Laplacian operator

For slightly more complex images, the edges extracted by the Laplacian operator are too detailed, making it difficult to see many areas clearly, which poses some difficulties for human visual perception. However, in some image recognition fields, such as using satellite images to identify ground vehicles, the vehicles on the ground are often small color blocks, and the Laplacian operator can well outline the edge contours and some details inside these vehicles. The Sobel operator is not as effective in capturing these internal details. At the same time, the property of "edges of edges" makes the edges extracted by the Laplacian operator suitable for image enhancement.

Canny Algorithm

The Canny algorithm is an optimization of the Sobel edge extraction. The Sobel operator represents all edges in the final image, regardless of their strength. This results in many invalid edges being extracted and displayed in the resulting image. The basic idea of the Canny algorithm is to filter out these edge information and only keep the pixels that are most likely to be edges. The implementation method is as follows:

  • Similar to the Sobel operator, calculate \(G_x\) and \(G_y\).
  • Calculate the weight \(weight=\sqrt{G_x^2+G_y^2}\) and the angle \(angle=atan\frac{G_y}{G_x}\) for each pixel based on \(G_x\) and \(G_y\).
  • Discretize the angle to the nearest multiple of \(45^o\).
  • For each pixel \((x,y)\) in the image, compare its weight \(weight(x,y)\) with the weights of the two pixels in the \(angle(x,y)\) direction and \(-angle(x,y)\) direction. If the weight of the pixel is not the largest, set it to 0.
  • Double threshold detection: set upper and lower thresholds for the brightness of the edges, and perform another filtering on the edges of the image. Finally, normalize the edge information to the range of 0-255. The thresholds can be manually set.

From left to right: original image, edges without double threshold detection, edges with lower threshold of 150, edges with lower threshold of 200

Appendix

References

[1] Rafael C. Gonzalez and Richard E. Woods, "Digital Image Processing", 3rd Edition, Beijing: Publishing House of Electronics Industry, 2017.

[2] Brook_icv, "Image Processing Basics (4): Gaussian Filter Detailed Explanation" [Online]. Available: https://www.cnblogs.com/wangguchangqing/p/6407717.html#autoid-4-1-0

[3] Naughty Stone 7788121, "Image Edge Detection: Canny Operator, Prewitt Operator, and Sobel Operator" [Online]. Available: https://www.jianshu.com/p/bed4ffe996a1

Source Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
import cv2 as cv
import numpy as np
import math

sigma=1.5   # Parameter of the Gaussian filter

def add_zeros(img,edge):    # Add zeros to the edges of the image
    shape=img.shape
    temp=np.zeros((shape[0]+2*edge,shape[1]+2*edge))
    for i in range(shape[0]):
        for j in range(shape[1]):
            temp[i+edge][j+edge]=img[i][j][0]
    return temp

def f(x,y):     # Define the 2D Gaussian distribution function
    return 1/(math.pi*sigma**2)*math.exp(-(x**2+y**2)/(2*sigma**2))

def gauss(n):   # Generate an n*n Gaussian filter
    mid=n//2
    filt=np.zeros((n,n))
    for i in range(n):
        for j in range(n):
            filt[i][j]=f(i-mid,j-mid)/f(-mid,-mid)
    return filt.astype(np.uint8)

def gauss_filter(img,n):    # Apply n*n convolutional blocks of Gaussian filtering to the image
    filt=gauss(n)
    con=1/np.sum(filt)
    shape=img.shape
    temp=add_zeros(img,n//2)
    result=np.zeros((shape[0],shape[1],1))
    for i in range(shape[0]):
        for j in range(shape[1]):
            tmp=0
            for k in range(n):
                for l in range(n):
                    tmp+=filt[k][l]*temp[i+k][j+l]
            result[i][j][0]=con*tmp
    return result.astype(np.uint8)

def unsharp_mask(img,n,is_mask=0):        # Unsharp mask using n*n Gaussian blur, return the mask if is_mask=1
    shape=img.shape
    new_img=np.zeros((shape[0],shape[1],1))
    for i in range(shape[0]):
        for j in range(shape[1]):
            new_img[i][j][0]=img[i][j][0]
    mask=new_img-gauss_filter(img,n)
    for i in range(shape[0]):
        for j in range(shape[1]):
            if i==0 or j==0 or i==shape[0]-1 or j==shape[1]-1:
                mask[i][j][0]=0
    result=new_img+mask
    result=result-np.min(result)
    result=result/np.max(result)*255
    mask=mask-np.min(mask)
    mask=mask/np.max(mask)*255
    if is_mask:
        return mask.astype(np.uint8)
    return result.astype(np.uint8)

sobelx=[[-1,0,1],[-2,0,2],[-1,0,1]]
sobely=[[-1,-2,-1],[0,0,0],[1,2,1]]
laplace4=[[0,-1,0],[-1,4,-1],[0,-1,0]]
laplace8=[[-1,-1,-1],[-1,8,-1],[-1,-1,-1]]

def filt_3(img,filt):          # Arbitrary 3*3 filter (image, operator)
    shape=img.shape
    temp=add_zeros(img,1)
    result=np.zeros((shape[0],shape[1],1))
    for i in range(shape[0]):
        for j in range(shape[1]):
            tmp=0
            for k in range(3):
                for l in range(3):
                    tmp+=filt[k][l]*temp[i+k][j+l]
            result[i][j][0]=tmp
    return result

def laplace_edge(img,filt): # Laplacian edges after normalization
    tmp=filt_3(img,filt)
    tmp=tmp-np.min(tmp)
    shape=tmp.shape
    for i in range(shape[0]):
        for j in range (shape[1]):
            if i==0 or j==0 or i==shape[0]-1 or j==shape[1]-1:
                tmp[i][j][0]=0
    tmp=tmp/np.max(tmp)*255
    return tmp.astype(np.uint8)

def laplace(img,filt):  # Overlay of the original image and the Laplacian edges
    tmp=filt_3(img,filt)
    shape=img.shape
    result=np.zeros((shape[0],shape[1]))
    for i in range(shape[0]):
        for j in range(shape[1]):
            result[i][j]=tmp[i][j][0]+img[i][j][0]
            if i==0 or j==0 or i==shape[0]-1 or j==shape[1]-1:
                result[i][j]=0
    result-=np.min(result)
    result=result/np.max(result)*255
    return result.astype(np.uint8)

def sobel(img):     # Extract Sobel edges of the image
    shape=img.shape
    sobx=filt_3(img,sobelx)
    soby=filt_3(img,sobely)
    result=np.zeros((shape[0],shape[1]))
    for i in range(shape[0]):
        for j in range(shape[1]):
            if i==0 or j==0 or i==shape[0]-1 or j==shape[1]-1:
                result[i][j]=0
            else:
                result[i][j]=math.sqrt(sobx[i][j][0]**2+soby[i][j][0]**2)
    result=result/np.max(result)*255
    return result.astype(np.uint8)

def canny(img,n=3):     # Extract image edges using the Canny algorithm (blur operation using an n*n Gaussian filter)
    de=[[1,0,-1,0],[1,1,-1,-1],[0,1,0,-1],[-1,1,1,-1]]
    shape=img.shape
    tmp=gauss_filter(img,n)
    sobx=filt_3(tmp,sobelx)
    soby=filt_3(tmp,sobely)
    weight,angle,result=np.zeros((shape[0],shape[1])),np.zeros((shape[0],shape[1])),np.zeros((shape[0],shape[1]))
    angle=angle.astype(np.int)
    for i in range(shape[0]):
        for j in range(shape[1]):
            weight[i][j]=math.sqrt(sobx[i][j][0]**2+soby[i][j][0]**2)
            if sobx[i][j][0]:
                angle[i][j]=round((math.atan(soby[i][j][0]/sobx[i][j][0])/(math.pi/4)-0.5))%4
    for i in range(shape[0]-2):
        for j in range(shape[1]-2):
            tmp_i,tmp_j=i+1,j+1
            if weight[tmp_i][tmp_j]<=weight[tmp_i+de[angle[tmp_i][tmp_j]][0]][tmp_j+de[angle[tmp_i][tmp_j]][1]] and weight[tmp_i][tmp_j]<=weight[tmp_i+de[angle[tmp_i][tmp_j]][2]][tmp_j+de[angle[tmp_i][tmp_j]][3]]:
                result[tmp_i][tmp_j]=0
            else:
                result[tmp_i][tmp_j]=weight[tmp_i][tmp_j]
    result=result/np.max(result)*255
    mean=np.mean(img)
    for i in range(shape[0]):
        for j in range(shape[1]):
            if result[i][j]<100:
                result[i][j]=0
    return result.astype(np.uint8)


filename=["test3_corrupt.pgm","test4.tif"]
for i in filename:
    img=cv.imread(i)
    cv.imwrite(i+"_mask.bmp",unsharp_mask(img,3,1))
    cv.imwrite(i+"_unsharp_mask.bmp",unsharp_mask(img,3))
    cv.imwrite(i+"_sobel.bmp",sobel(img))
    cv.imwrite(i+"_canny.bmp",canny(img,3))
    cv.imwrite(i+"laplace4_edge.bmp",laplace_edge(img,laplace4))
    cv.imwrite(i+"laplace8_edge.bmp",laplace_edge(img,laplace8))