2024-06-06碼農

1. 二維離散傅立葉變換

DFT 是 Discrete Fourier Transform 即離散傅立葉變換的簡稱。二維離散傅立葉變換（2D Discrete Fourier Transform，簡稱 2D DFT）是將二維離散訊號（例如數位影像）從空間域變換到頻率域的一種數學工具。

1.1 定義

二維離散傅立葉變換的定義如下：

設 f(x,y) 是一個 M×N 的影像，其中 x=0,1,…,M−1 和 y=0,1,…,N−1。則其二維離散傅立葉變換 F(u,v) 定義為：

其中，u=0,1,…,M−1 和 v=0,1,…,N−1。

二維離散傅立葉逆變換：

1.2 矩陣乘法表示

二維離散傅立葉變換可以表示為矩陣乘法，這是一種更直觀、更易於理解的表示方式。

其中，

F 是 M×N 的傅立葉變換矩陣，代表變換後的影像在頻率域的數據。

f 是 M×N 的影像矩陣。

W 是 M×N 的傅立葉變換基矩陣，由傅立葉變換系數構成。

傅立葉變換基矩陣 W 的元素定義為：

矩陣乘法表示傅立葉變換的過程，是將原始影像每個像素值與基矩陣中對應元素進行乘積並求和，得到變換後的每個頻率分量。

二維離散傅立葉逆變換使用矩陣乘法表示：

其中，是 W 的共軛轉置。

矩陣乘法表示的優勢:

簡潔直觀：用矩陣乘法表示二維離散傅立葉變換，可以直觀地理解變換過程，並將計算過程轉化為矩陣運算，方便電腦實作。

易於編程：矩陣乘法是編程中常用的操作，可以用各種程式語言方便地實作二維離散傅立葉變換。

⾼效計算：基於矩陣乘法的快速傅立葉變換（FFT）演算法，可以顯著提高計算效率，是影像處理中常用的傅立葉變換實作方法。

便於擴充套件：矩陣乘法可以擴充套件到更高維度的傅立葉變換。

二維離散傅立葉變換按照定義的公式來計算，其計算量非常大。在實際套用中，通常使用快速傅立葉變換（FFT）演算法來計算二維離散傅立葉變換。 FFT 演算法可以將二維離散傅立葉變換的計算量從降低到。

2. 二維離散傅立葉變換的性質

2.1 尺度變換

如果影像在空間域中進行尺度變換，則其傅立葉變換在頻率域中會發生相應的尺度反變換。

f(x,y) 將其沿著 x 方向縮放 s 倍，或者沿著 y 方向縮放 t 倍，則其傅立葉變換 F(u,v) 沿著 u 方向縮放 1/s 倍，或者沿著 v 方向縮放 1/t 倍。

2.2 平移性

如果 f(x,y) 沿著 x 方向平移 p 個單位，或者沿著 y 方向平移 q 個單位，則其傅立葉變換 F(u,v) 沿著 u 方向平移 −p 個單位，或者沿著 v 方向平移 −q 個單位。

該性質意味著影像在空間域中平移，其傅立葉變換在頻率域中對應方向上發生平移。

2.3 旋轉不變性

如果 f(x,y) 繞原點旋轉 θ 角度，則其傅立葉變換 F(u,v) 繞原點旋轉 −θ 角度。

該性質意味著影像在空間域中旋轉，其傅立葉變換在頻率域中對應方向上旋轉。

2.4 周期性

二維離散傅立葉變換的周期性體現在以下兩個方面：

水平方向周期性：將二維離散傅立葉變換的結果沿水平方向平移 M 個單位（即影像寬度），則結果將與原始結果相同。即：

垂直方向周期性：將二維離散傅立葉變換的結果沿垂直方向平移 N 個單位（即影像高度），則結果將與原始結果相同。即：

二維離散傅立葉變換的周期性可以用以下公式來表示：

其中，k 和 l 是任意整數。

2.5 共軛對稱性

二維離散傅立葉變換的實部和虛部具有共軛對稱性。即

該性質意味著二維離散傅立葉變換的頻譜在實數軸和虛軸對稱。

2.6 傅立葉頻譜和相角

二維離散傅立葉變換通常是復函式，因此可以用極座標形式表示：

振幅：

稱為傅立葉頻譜(頻譜)。

而稱為相角(相位譜)。

傅立葉頻譜它通常用影像表示，其中亮度表示振幅的大小。傅立葉頻譜可以用來分析影像的頻率特征。相角是二維離散傅立葉變換的相位值，它表示影像中各個頻率分量的相位資訊。相角可以用來分析影像的相位特征。

2.7 二維離散摺積定理

兩個函式在空間域中摺積的結果的傅立葉變換等於這兩個函式的傅立葉變換在頻率域中的乘積。

設 f(x,y) 和 h(x,y) 是兩個 M×N 的影像，則其摺積 g(x,y) 定義為：

其中，x=0,1,…,M−1 和 y=0,1,…,N−1。

根據摺積定理，有：

其中，G(u,v)、F(u,v) 和 H(u,v) 分別是 g(x,y)、f(x,y) 和 h(x,y) 的二維離散傅立葉變換。

3. 範例

3.1 將頻率域影像轉換成空間域的影像

下面的程式碼，展示了影像從空間域轉換到頻域，再轉換回空間域的過程。

#include<opencv2/opencv.hpp> #include<opencv2/core.hpp> #include<opencv2/highgui.hpp> usingnamespacestd; usingnamespace cv; // fft 變換後進行頻譜中心化 voidfftshift(cv::Mat &plane0, cv::Mat &plane1) { int cx = plane0.cols / 2; int cy = plane0.rows / 2; cv::Mat q0_r(plane0, cv::Rect(0, 0, cx, cy)); // 元素座標表示為(cx, cy) cv::Mat q1_r(plane0, cv::Rect(cx, 0, cx, cy)); cv::Mat q2_r(plane0, cv::Rect(0, cy, cx, cy)); cv::Mat q3_r(plane0, cv::Rect(cx, cy, cx, cy)); cv::Mat temp; q0_r.copyTo(temp); //左上與右下交換位置(實部) q3_r.copyTo(q0_r); temp.copyTo(q3_r); q1_r.copyTo(temp); //右上與左下交換位置(實部) q2_r.copyTo(q1_r); temp.copyTo(q2_r); cv::Mat q0_i(plane1, cv::Rect(0, 0, cx, cy)); //元素座標(cx,cy) cv::Mat q1_i(plane1, cv::Rect(cx, 0, cx, cy)); cv::Mat q2_i(plane1, cv::Rect(0, cy, cx, cy)); cv::Mat q3_i(plane1, cv::Rect(cx, cy, cx, cy)); q0_i.copyTo(temp); //左上與右下交換位置(虛部) q3_i.copyTo(q0_i); temp.copyTo(q3_i); q1_i.copyTo(temp); //右上與左下交換位置(虛部) q2_i.copyTo(q1_i); temp.copyTo(q2_i); } intmain() { Mat src = imread(".../girl.jpg"); imshow("src", src); Mat gray; cvtColor(src, gray, COLOR_BGR2GRAY); gray.convertTo(gray, CV_32FC1); cv::Mat planes[] = {Mat_<float>(gray), cv::Mat::zeros(src.size() , CV_32FC1) }; cv::Mat complexI; cv::merge(planes, 2, complexI); // 合並通道（把兩個矩陣合並為一個2通道的Mat類容器） cv::dft(complexI, complexI); // 進行傅立葉變換 // 分離通道（陣列分離） cv::split(complexI, planes); // 計算幅值 cv::Mat mag,mag_log,mag_nor; cv::magnitude(planes[0], planes[1], mag); // 幅值對數化：log（1+m） mag += Scalar::all(1); cv::log(mag, mag_log); cv::normalize(mag_log, mag_nor, 1,0, NORM_MINMAX); // 頻譜中心化 fftshift(planes[0], planes[1]); cv::Mat dst; // 再次頻譜中心化 fftshift(planes[0], planes[1]); cv::merge(planes, 2, dst); // 實部與虛部合並 cv::idft(dst, dst); // idft 結果也為復數 dst = dst / dst.rows / dst.cols; cv::split(dst, planes);//分離通道，主要獲取通道 convertScaleAbs(gray,gray); convertScaleAbs(planes[0],planes[0]); imshow("gray", gray); imshow("result", planes[0]); waitKey(0); return0; }

3.2 添加盲浮水印

盲浮水印 功能是指將浮水印以不可見的形式添加到圖片中，這樣的浮水印不會圖片觀感產生影響，同時也保證了圖片的原創性。當發現圖片被盜取後，可以透過提取盲浮水印，將不可見的浮水印提取出來，驗證圖片的歸屬。

下面的程式碼展示了添加盲浮水印的流程:

原圖 -> 透過傅立葉變換，在頻域上添加浮水印 -> 最佳化由 dft 操作產生的影像 -> 頻域圖 -> 傅立葉逆變換 -> 最終在空間域生成含有盲浮水印的影像

#include<opencv2/opencv.hpp> #include<opencv2/core.hpp> #include<opencv2/highgui.hpp> usingnamespacestd; usingnamespace cv; // 頻譜中心化 voidshiftDFT(cv::Mat &magnitudeImage) { // 如果影像的尺寸是奇數的話對影像進行裁剪並重新排列（減去補充部份） magnitudeImage = magnitudeImage(cv::Rect(0, 0, magnitudeImage.cols & -2, magnitudeImage.rows & -2)); // 重新排列影像的象限，使得影像的中心在象限的原點 int cx = magnitudeImage.cols / 2; int cy = magnitudeImage.rows / 2; cv::Mat q0 = cv::Mat(magnitudeImage, cv::Rect(0, 0, cx, cy)); // 左上 cv::Mat q1 = cv::Mat(magnitudeImage, cv::Rect(cx, 0, cx, cy)); // 右上 cv::Mat q2 = cv::Mat(magnitudeImage, cv::Rect(0, cy, cx, cy)); // 左下 cv::Mat q3 = cv::Mat(magnitudeImage, cv::Rect(cx, cy, cx, cy)); // 右下 // 交換象限 cv::Mat tmp = cv::Mat(); // 左上與右下交換 q0.copyTo(tmp); q3.copyTo(q0); tmp.copyTo(q3); // 右上與左下交換 q1.copyTo(tmp); q2.copyTo(q1); tmp.copyTo(q2); } // 將頻域的影像轉換為空間域的影像 Mat transformImage(cv::Mat complexImage, vector<cv::Mat> allPlanes) { cv::Mat invDFT; cv::idft(complexImage, invDFT, cv::DFT_SCALE | cv::DFT_REAL_OUTPUT, 0); cv::Mat restoredImage; invDFT.convertTo(restoredImage, CV_8U); // 合並多通道 allPlanes.erase(allPlanes.begin()); allPlanes.insert(allPlanes.begin(), restoredImage); cv::Mat dst; cv::merge(allPlanes, dst); return dst; } intmain() { Mat src = imread(".../girl.jpg"); imshow("src", src); Mat complexI; // 傅立葉變換結果，復數 vector<Mat> planes; vector<Mat> allPlanes; // 將多通道分為單鍊結（因為讀入的是彩色圖） cv::split(src, allPlanes); // 只獲取 B 通道 Mat padded = cv::Mat(); if (allPlanes.size() > 1) { padded = allPlanes[0]; } else { padded = src; } padded.convertTo(padded, CV_32F); // 將單鍊結擴充套件至雙鍊結，以接收 DFT 的復數結果 planes.push_back(padded); planes.push_back(cv::Mat::zeros(padded.size(), CV_32F)); // 將 planes 陣列組合合並成一個多通道 Mat cv::merge(planes, complexI); // 進行離散傅立葉變換 cv::dft(complexI, complexI); // 添加文本浮水印 string text = "tony test"; Scalar scalar = cv::Scalar(0,0,0,0); Point point = cv::Point(300, 300); int fontFace = cv::FONT_HERSHEY_TRIPLEX; double fontScale = 8.0; putText(complexI, text, point, fontFace, fontScale, scalar); flip(complexI, complexI, -1); putText(complexI, text, point, fontFace, fontScale, scalar); flip(complexI, complexI, -1); planes.clear(); cv::split(complexI, planes); // 計算幅值矩陣 Mat mag; cv::magnitude(planes[0], planes[1], mag); mag += Scalar::all(1); // 轉換到對數尺度 cv::log(mag, mag); // 剪下和重分布振幅圖像限 shiftDFT(mag); // 歸一化，用 0 到 255 之間的浮點值將矩陣變換為視覺化的影像格式 mag.convertTo(mag, CV_8UC1); normalize(mag, mag, 0, 255, cv::NORM_MINMAX, CV_8UC1); imshow("Frequency Domain Image", mag); Mat result = transformImage(complexI, allPlanes); planes.clear(); imshow("Spatial Domain Image", result); waitKey(0); return0; }