最後一天了，大家知道傅立葉轉換適合用來處理光啊，頻率啊。

你們知道傅立葉轉換也很適合用來處理圖形嗎？

今天繼續延伸前面的系列，跟著助教來寫圖片的快速傅立葉轉換程式。

程式碼請見下方：

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

#定義計算二進制的位元反轉函數
def bit_reversal(number, bitsize):
    binary = bin(number) #將整數轉換為二進制表示
    reverse = binary[-1:1:-1] #反轉二進位字符串
    reverse = str(reverse) + "0" * (bitsize - len(reverse)) #根據位元大小 (bitsize)，對反轉後的二進位字符串補0
    return int(reverse, 2) #將反轉後的二進制數轉為整數

#定義根據蝴蝶圖進行數據合併的FFT_element函數   
def FFT_element(array_e, array_o): #array_e和array_o 分別是原始數據中所有偶數和奇數索引的位置
    N = len(array_e) * 2  #將一個長度為 N 的陣列拆分成兩個長度為 N/2 的陣列進行處理
    W_array = np.array([np.exp(complex(0, 2*np.pi/N*i)) for i in range(int(N/2))]) #製作傅立葉變換中的旋轉因子
    output_array = np.hstack((array_e + array_o*W_array, array_e - array_o*W_array)) #將兩組數據合併
    return output_array

# 一維的快速傅立葉轉換
def FFT_1D(input):
    length = len(input) #取input的長度，得輸入數據的總數量
    bit_size = int(np.log2(length)) #計算輸入數據長度的對數（以 2 為底），即 log2(length)，用來確定位反轉(bit reversal)所需的位數。
    output = [] #初始化一個空列表，稍後將用來存放經過位反轉後重新排列的數據
    for i in range(length): 
        bit_reverse = bit_reversal(i, bit_size) #計算i 的位反轉結果。例: i 是二進位的 001（十進位為 1），在三位的情況下，它的位反轉結果是 100（十進位為 4）
        output.append(np.array([input[bit_reverse]])) #將重新排列的array加入到 rearrange 列表
    for s in range(bit_size): 
        tmp = [] #初始化一個空列表
        for s2 in range(int(len(output)/2)): #遍歷目前數據的1/2
            tmp.append(FFT_element(output[2 * s2], output[2 * s2 + 1])) #帶入FFT_unit對數據執行蝴蝶操作，產生新的頻域數據
        output = tmp
    return output[0] #得最終的計算結果

# 二維快速傅立葉轉換
def FFT_2D(input):
    len_x, len_y = input.shape
    output = np.zeros((len_x, len_y), dtype=complex)
    for i in range(len_x):
        output[i, :] = FFT_1D(input[i, :])
    tmp = output
    for j in range(len_y):
        output[:, j] = FFT_1D(tmp[:, j])
    return output


def FT_shift_2D(input):
    len_x, len_y = input.shape
    shift_x = int(len_x / 2)
    shift_y = int(len_y / 2)
    return np.roll(input, (shift_x, shift_y), (0, 1))

def fill_zeros(array):
    len_x, len_y = array.shape
    bit_x = int(np.log2(len_x))+1  #計算新陣列的行數所需的位元長度，+1用更大的次方來容納擴展後的陣列
    bit_y = int(np.log2(len_y))+1  #同上
    output = np.zeros((2**bit_x, 2**bit_y))  #建立擴展後的零填充陣列
    output[0:len_x, 0:len_y] = array[:,:]  # 將原始陣列填充進去
    return output

#讀取影像並進行翻轉和填充操作
img = cv2.imread("Frieren potion.jpeg", cv2.IMREAD_REDUCED_GRAYSCALE_2)
img2 = cv2.flip(img, -1)
img3 = fill_zeros(img2)

far_field = FFT_2D(img3)
far_field = FT_shift_2D(far_field)
len_x, len_y = far_field.shape #取得len_x和len_y的長度
intensity1 = [[abs(far_field[j,i])**2
                  for i in range(len_y)]
                  for j in range(len_x)] #計算光強度

#對頻譜進行逆傅立葉變換後的結果，轉換回空間域的影像
reverse = FFT_2D(far_field) #將遠場繞射結果作FFT_2D轉換
len_x, len_y = reverse.shape #取得len_x和len_y的長度
intensity2 = [[abs(reverse[j,i])**2
                  for i in range(len_y)]
                  for j in range(len_x)] #計算光強度


plt.subplot(3,1,1)
plt.contourf(img3, levels = 90, cmap = "grey")
plt.subplot(3,1,2)
plt.contourf(intensity1, levels = 90, cmap = "grey")
plt.subplot(3,1,3)
plt.contourf(intensity2, levels = 90, cmap = "grey")
plt.show()

原始阿嬤圖片為：

(原本想用我婆莉央，但怕大家只顧著看水溝，還是收斂點好惹)

程式碼跑出的圖如下：

第一張圖img2是FFT轉換前，預先做灰階處理和翻轉。
第二張圖是遠場繞射傅立葉轉換的結果。
第三張圖是以第二張圖進行二次傅立葉轉換的結果。

⋯⋯

怎麼好像哪裡怪怪的＠＠

img2 = cv2.flip(img, -1)應該要把啊罵上下左右顛倒，可是img2卻只有左右翻轉。
會設計這一行也是因為傅立葉轉換會上下左右翻轉，所以為了讓轉換回來的圖是正的，我們在轉換前先翻轉一次。

如果這邊改成將影像沿Y軸翻轉(1)
img2 = cv2.flip(img, 1)

這樣就如同設計那樣，最後出來的圖是正的。

問了一下chatGPT，為何我的cv2.flip，1和-1出來的結果這麼詭異？

chatGPT：應該是你的opencv版本或設定和一般定義的不同。

我：⋯⋯好喔⋯⋯

第30篇真的漫長，感謝大家陪我這文組摸索物理世界到現在，明天還會有一篇總結，希望之後還會繼續學習寫作！
再次感恩讚嘆！謝謝謝助教！謝謝無料ChatGPT！