數字影像處理（一）之灰度轉換和卷積python實現

使用Python實現數字影像處理中如下功能：

1. 彩色影像轉成灰度影像
2. 實現影像的相關&卷積操作
3. 實現影像的高斯核卷積

使用的庫和python版本如下：

• imageio：2.9.0 用於讀取磁碟中的圖片檔案
• numpy：1.20.3 用於矩陣等操作
• matplotlib：3.4.2 用於畫圖
• python：3.8.11

讀取影像

在進行影像處理操作前，首先需要對影像進行讀取。這裡使用imageio庫對圖片進行讀取，並將其轉成numpy陣列。

下面定義一個covert_img_to_array函式，用於讀取圖片。

def covert_img_to_array(self, path:str) -> np.array:
    """[將圖片轉成Array便於處理]

        Args:
            path (str): [圖片儲存位置]

        Returns:
            np.array: [返回numpy陣列，陣列元素uint8]
        """
    return np.array(imageio.imread(path))

展示圖片

使用matplotlib庫用於展示圖片，為了更高的展示如片，定義下show_img函式，當不指定col或者row時儘量以方正的形式去展示圖片。

def show_img(self,title:str, imgs:list, cmaps:list,row:int = 0,col:int = 0):
    """展示圖片 len(imgs) must equal to the len of cmaps

    Args:
        title (str): [影像標題]
        imgs (list): [圖片元組]
        cmaps (list): [mask,plt以何種形式展示圖片，可參考官方文件使用：'gray'表示灰度圖，None表示彩色圖]
        row (int, optional): [指令row]. Defaults to 0.
        col (int, optional): [指令col]. Defaults to 0.
    """
    if len(imgs) != len(cmaps):
        print("圖片和mask的len必須相同")
    else:
        if row == 0 and col !=0:
            row = np.ceil(len(imgs)/col).astype("uint8")
        elif row!=0 and col == 0:
            col = np.ceil(len(imgs)/row).astype("uint8")
        elif row*col < len(imgs):
            # 儘量以方正的形式去展示圖片
            row = np.ceil(np.sqrt(len(imgs))).astype("uint8")
            col = np.ceil(len(imgs)/row).astype("uint8")

        for index,img in enumerate(imgs):
            plt.subplot(row,col,index+1)
            plt.imshow(img,cmap=cmaps[index])
        plt.suptitle(title)
        plt.show()

彩色影像轉成灰度影像

彩色影像一般來說RGB表示的。也就是說，如果有一張64*64大小的圖片，那麼它在numpy中便是以64*64*3的shape進行儲存的。將RGB圖片轉成灰度圖有兩種方式：

1. \(gray=\frac{R+G+B}{3}\)
2. \(gray=R*0.2989 + G*0.5870 + B*0.1140\) 這種灰度轉換稱之為NTSC標準，考慮了人類的彩色感知體驗。

下面定義covert_rgb_to_gray函式，其中method如果為average，則使用第一種方式灰度轉換方式；預設為NTSC，使用第二種方式轉換。

def covert_rgb_to_gray(self, image:np.array, method:str = 'NTSC') -> np.array:
    """將RGB影像轉成gray影像

    Args:
        image (np.array): [rgb影像]
        method (str, optional): [轉換模式]. Defaults to 'NTSC'.

    Returns:
        Array: [返回的灰度影像]
    """
    if method == 'average':
        gray_img = image[:,:,0]/3+image[:,:,1]/3+image[:,:,2]/3
    
    else:
        gray_img = image[:,:,0]*0.2989 + image[:,:,1]*0.5870 + image[:,:,2]*0.1140
    return gray_img

影像卷積

影像卷積的公式如下所示，\(g\)代表輸入的畫素矩陣，\(w\)代表的是權重係數矩陣也就是所謂的卷積核kernel。

\[h(x,y) =\sum_{s=-a}^{a} \sum_{t=-b}^{b} w(s,t)g(x-s,y-t) \]

這裡有一個很需要值得注意的點，那就是相關操作。相關操作和卷積很類似，相關操作的公式如下：

\[h(x,y) =\sum_{s=-a}^{a} \sum_{t=-b}^{b} w(s,t)g(x+s,y+t) \]

在網路有一些部落格文章，在解釋卷積的時候，使用的是第一個公式，但是在做計算或者實現程式碼的時候卻用的是第二個公式，這樣做是不對的。因為卷積的kernel與相關的kernel相差了\(180^{\circ}\)。

但是值得注意的是，在卷積神經網路中，實際上使用的數學公式是相關相關運算，如下圖所示。因為在CNN中，kernel的引數是學習過來的，kernel是否翻轉並不會影響結果。

理解卷積

前置知識：

卷積定理指出，函式卷積的傅立葉變換是函式傅立葉變換的乘積。至於推導，可以查一下資料。

\[\mathcal{F}\{f * g\}=\mathcal{F}\{f\} \cdot \mathcal{F}\{g\} \]

提一下影像卷積的含義。如果一個如下的均值濾波器對影像進行卷積，從人類的直覺進行出發，可以去除噪聲和平滑影像。（在影像中，一般影像噪聲的頻率比較大，影像邊緣部分的頻率也比較大。 因此使用均值濾波器可以去除噪聲和平滑影像。）

\[1 / 9\left[\begin{array}{lll} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{array}\right] \]

那麼為什麼會造成這種現象呢？如何從數學的角度來解釋均值濾波器的作用呢？

如下所示，圖左邊是一個一維均值濾波器的函式影像，圖右邊是均值函式在頻域上面的影像。在右邊影像上，可以發現一個很明顯的特點：頻率越高，\(F(\mu)\)越小。

那麼如果將\(F(\mu)\)與某另外一個頻域上面的函式（比如影像）相乘，顯而易見，如果影像的頻率越高，則\(F(\mu)\)與之相乘被拖下水的的程度就越大。也就是說，相乘之後，頻率低的就被抬上去了，頻率高的被拉下去了。

說的細一點，其實從上圖可以看到，隨著頻率的增大，\(F(\mu)\)並不是嚴格的下降，中間有一個波浪的起伏，這樣會在邊緣造成一些不好的現象。但是高斯濾波不會有這種情況。後面會介紹高斯濾波。

均值濾波器的二維頻域圖如下所示：

矩陣點積

下面定義矩陣點積函式。

def __matrix_dot_product(self,matrix,kernel):
    """矩陣點乘 [1,2,3]*[4,5,6] = 1*4 + 2*5 + 3*6 = 32
    
    Args:
        matrix ([type]): [部分影像]
        kernel ([type]): [kernel]

    Returns:
        [type]: [點乘結果]
    """
    if len(matrix) != len(kernel):
        print("點積失敗，大小不一致")
    else:
        # 速度快
        return (np.multiply(matrix,kernel)).sum()

        # result = 0
        # for i, row_nums in enumerate(matrix):
        #     for j,num in enumerate(row_nums):
        #         result += num * kernel[i][j]
        # return result 

影像padding

如果不對影像進行padding的話，會造成一個現象，影像越卷越小。在卷積的時候，我們希望卷積後的影像大小與原影像保持一致（CNN網路可能會越卷越小），因此需要對影像進行padding。padding有兩種方式，一種在填充0，一種是填充與其距離最近的元素。下圖中影像周圍虛線部分就是padding的元素。

下面是實現padding操作的具體函式。實際上，可以直接使用np.pad操作實現。（但是我的作業要求不能使用pad操作，只能自己實現）

    def __padding(self, padding_type:str, image:np.array, padding_w:int, padding_h:int):
        """對圖片進行padding

        Args:
            padding_type (str): [padding方式]
            image (np.array): [圖片]
            padding_w (int): [寬度pdding]
            padding_h (int): [高度padding，一般來說padding_w = padding_h]

        Returns:
            [type]: [返回padding之後的結果]
        """
        image_w = image.shape[0]
        image_h = image.shape[1]

        padding_image = np.zeros((image_w+padding_w*2,image_h+padding_h*2))
        padding_image[padding_w:padding_w+image_w,padding_h:padding_h+image_h] = image

        if padding_type == 'zero':
            return padding_image

        if padding_type == "replicate": 
            # 補充四個角
            padding_image[0:padding_w+1,0:padding_h+1] = image[0,0]
            padding_image[image_w+padding_w-1:,0:padding_h+1] = image[image_w-1,0]
            padding_image[0:padding_w+1,image_h+padding_h-1:] = image[0,image_h-1]
            padding_image[image_w+padding_w-1:,image_h+padding_h-1:] = image[image_w-1,image_h-1]

            # 補充旁邊的元素
            for i in range(padding_w+1,image_w+padding_w-1):
                padding_image[i,0:padding_h] = image[i-padding_w,0]
                padding_image[i,image_h+padding_h:] = image[i-padding_w,image_h-1]
                
            
            for i in range(padding_h+1,image_h+padding_h-1):
                padding_image[0:padding_w,i] = image[0,i-padding_h]
                padding_image[image_w+padding_w:,i] = image[image_w-1,i-padding_h]
            return padding_image

如果想使得卷積之後的結果與原影像一致，padding_w,padding_h為卷積核大小的一半（向下取整，卷積核大小一般是奇數）。比如核的大小是\(5 \times 5\)，那麼padding的長寬便是\(2\)。

影像相關操作

前面說過影像的卷積實際上就是將kernel進行翻轉\(180^{\circ}\)，然後進行相關運算，因此可以先定義相關操作函式：

def corr2D(self, image:np.array, kernel:np.array, padding:str = 'zero') -> np.array:
    """對圖片進行相關運算。

    Args:
        image (np.array): [(*,*)shape的圖片]
        kernel (np.array): [kernel，kernel為奇數]
        padding (str, optional): [zero以零填充，replicate以鄰近的填充]. Defaults to 'zero'.

    Returns:
        [type]: [description]
    """
    kernel_size_w = kernel.shape[0]
    kernel_size_h = kernel.shape[1]

    image_w,image_h = image.shape
    
    
    padding_w = kernel_size_w // 2
    padding_h = kernel_size_h // 2

    # 將圖片padding起來
    padding_image = self.__padding(padding,image,padding_w,padding_h)

    new_image = np.zeros((image_w,image_h))
    for i in range(image_w):
        for j in range(image_h):
            new_image[i][j] = self.__matrix_dot_product(padding_image[i:i+kernel_size_w,j:j+kernel_size_h],kernel)

    return new_image.clip(0,255).astype("uint8")

卷積操作

旋轉kernel

旋轉kernel的程式碼很簡單，如下所示，通過以下操作可以將行和列翻轉（相當於反轉了\(180^{\circ}\)）。

def flip_180(self, arr: np.array) -> np.array:
    return arr[::-1,::-1]

卷積

將kernel繼續寧翻轉，然後進行相關運算便是卷積了。

def conv2D(self, image:np.array, kernel:np.array, padding:str = 'zero') -> np.array:
    """二維卷積

    Args:
        image (np.array): [(*,*)shape的圖片]
        kernel (np.array): [kernel，kernel為奇數]
        padding (str, optional): [zero以零填充，replicate以鄰近的填充]. Defaults to 'zero'.

    Returns:
        [type]: [卷積好的結果]
    """
    return self.corr2D(image,self.flip_180(kernel),padding)

高斯核

二維高斯核的公式如下所示：

\[G(x, y,\sigma_x,\sigma_y)=\frac{1}{2 \pi \sigma_{x}\sigma_{y}} e^{-\left(\frac{x^{2}}{2{\sigma_x}^2} + \frac{y^{2}}{2{\sigma_y}^2}\right)} \]

二維高斯核的頻域圖如下所示。

下面是二維高斯濾波函式的定義，其中\(\sigma_x=\sigma_y=sig\)。並對卷積核進行歸一化，使得所有元素加起來和為1。

    def gauss_2d_kernel(self,sig,m=0):
        """產生高斯核

        Args:
            sig ([type]): [高斯核引數 sigx = sigy]
            m (int, optional): [高斯kernel的大小]. Defaults to 0. if m=0，then m = ceil(3*sig)*2 +1

        Returns:
            [type]: [m*m大小的高斯核]
        """
        fit_m = math.ceil(3 * sig)*2+1

        if m == 0:
            m = fit_m
        if m < fit_m:
            print("你的核的size應該大一點")
        
        # 中心點
        center = m //2
        kernel = np.zeros(shape=(m,m))
        for i in range(m):
            for j in range(m):
                kernel[i][j] = (1/(2*math.pi*sig**2))*math.e**(-((i-center)**2+(j-center)**2)/(2*sig**2))
		# 歸一化
        return kernel/(kernel.sum())

結果

灰度轉換結果

高斯核卷積

參考

• 數字影像處理（第三版）