關于傅裡葉譜和相角的進一步說明：
在頻域下，原始圖像經過dft變換後，可得到兩個矩陣，分别是複數的實部和虛部；Z=A+Bi
在複數坐标系下，某一點像素可由值和角度表示；
如果隻保留相角資訊，則隻需将r=1，即實部和虛部同時縮小r倍；
如果隻保留幅值資訊，則隻需将a值和b值同時設定為
都設定為45,135,225,315度
如果使用圖像a的幅值圖像b的相角，則按下圖計算
distance1/distance2 = r/r2 = i/i2
同理，用a的相位和b的幅值亦可得到；
按照書中婦女圖像4.27和白色矩形4.24（上篇文章實驗圖檔）實驗
4.27婦女原圖和相位圖
使用幅值和相角分别重建圖像：
使用婦女的相角和矩形的譜 ------- 矩形的相角和婦女的譜 分别重建圖像：
由此可見：
相角是描述圖像中關鍵的形狀特性的；譜（幅值）是描述圖像灰階資訊的
代碼實作：
#include "opencv2/opencv.hpp"
using namespace cv;
//将幅度歸一，相角保持不變
void one_amplitude(Mat &complex_r, Mat &complex_i, Mat &dst )
{
Mat temp[] = {Mat::zeros(complex_r.size(),CV_32FC1), Mat::zeros(complex_r.size(),CV_32FC1)};
float realv=0.0,imaginv=0.0;
for(int i=0;i<complex_r.cols;i++){
for( int j = 0; j < complex_r.rows; j++ ){
realv = complex_r.at<float>(i,j);
imaginv = complex_i.at<float>(i,j);
float distance=sqrt(realv*realv+imaginv*imaginv);
temp[0].at<float>(i,j) =realv/distance;
temp[1].at<float>(i,j) =imaginv/distance;
}
}
merge(temp, 2, dst);
}
//将相角歸一，幅值保持不變
void one_angel(Mat &complex_r, Mat &complex_i, Mat &dst)
{
Mat temp[] = {Mat::zeros(complex_r.size(),CV_32FC1), Mat::zeros(complex_r.size(),CV_32FC1)};
float realv=0.0,imaginv=0.0;
for(int i=0;i<complex_r.cols;i++){
for( int j = 0; j < complex_r.rows; j++ ){
realv = complex_r.at<float>(i,j);
imaginv = complex_i.at<float>(i,j);
float distance=sqrt(realv*realv+imaginv*imaginv);
temp[0].at<float>(i,j) =distance/sqrt(2);
temp[1].at<float>(i,j) =distance/sqrt(2);
}
}
merge(temp, 2, dst);
}
//使用1的幅值和2的相位合并
void mixed_amplitude_with_phase(Mat &real1, Mat &imag1, Mat &real2, Mat &imag2, Mat &dst)
{
if(real1.size() != real2.size() ){
std::cerr << "image don't ==" << std::endl;
return;
}
Mat temp[] = {Mat::zeros(real1.size(),CV_32FC1), Mat::zeros(real1.size(),CV_32FC1)};
float realv1=0.0, imaginv1=0.0, realv2=0.0, imaginv2=0.0;
for(int i=0;i<real1.cols;i++){
for( int j = 0; j < real1.rows; j++ ){
realv1 = real1.at<float>(i,j);
imaginv1 = imag1.at<float>(i,j);
realv2 = real2.at<float>(i,j);
imaginv2 = imag2.at<float>(i,j);
float distance1=sqrt(realv1*realv1+imaginv1*imaginv1);
float distance2=sqrt(realv2*realv2+imaginv2*imaginv2);
temp[0].at<float>(i,j) =(realv2*distance1)/distance2;
temp[1].at<float>(i,j) =(imaginv2*distance1)/distance2;
}
}
merge(temp, 2, dst);
}
//使用1的相位和2的幅值合并
void mixed_phase_with_amplitude(Mat &real1, Mat &imag1, Mat &real2, Mat &imag2, Mat &dst)
{
if(real1.size() != real2.size() ){
std::cerr << "image don't ==" << std::endl;
return;
}
Mat temp[] = {Mat::zeros(real1.size(),CV_32FC1), Mat::zeros(real1.size(),CV_32FC1)};
float realv1=0.0, imaginv1=0.0, realv2=0.0, imaginv2=0.0;
for(int i=0;i<real1.cols;i++){
for( int j = 0; j < real1.rows; j++ ){
realv1 = real1.at<float>(i,j);
imaginv1 = imag1.at<float>(i,j);
realv2 = real2.at<float>(i,j);
imaginv2 = imag2.at<float>(i,j);
float distance1=sqrt(realv1*realv1+imaginv1*imaginv1);
float distance2=sqrt(realv2*realv2+imaginv2*imaginv2);
temp[0].at<float>(i,j) =(realv1*distance2)/distance1;
temp[1].at<float>(i,j) =(imaginv1*distance2)/distance1;
}
}
merge(temp, 2, dst);
}
cv::Mat fourior_inverser( Mat &_complexim )
{
Mat dst;
Mat iDft[]={Mat::zeros(_complexim.size(),CV_32F),Mat::zeros(_complexim.size(),CV_32F)};//建立兩個通道，類型為float，大小為填充後的尺寸
idft(_complexim,_complexim);//傅立葉逆變換
split(_complexim,iDft);//結果貌似也是複數
magnitude(iDft[0],iDft[1],dst);//分離通道，主要擷取0通道
// dst += Scalar::all(1); // switch to logarithmic scale
// log(dst, dst);
//歸一化處理，float類型的顯示範圍為0-255,255為白色，0為黑色
normalize(dst, dst, 0, 255, NORM_MINMAX);
dst.convertTo(dst, CV_8U);
return dst;
}
void move_to_center( Mat ¢er_img)
{
int cx = center_img.cols/2;
int cy = center_img.rows/2;
Mat q0(center_img, Rect(0, 0, cx, cy)); // Top-Left - Create a ROI per quadrant
Mat q1(center_img, Rect(cx, 0, cx, cy)); // Top-Right
Mat q2(center_img, Rect(0, cy, cx, cy)); // Bottom-Left
Mat q3(center_img, Rect(cx, cy, cx, cy)); // Bottom-Right
Mat tmp; // swap quadrants (Top-Left with Bottom-Right)
q0.copyTo(tmp);
q3.copyTo(q0);
tmp.copyTo(q3);
q1.copyTo(tmp); // swap quadrant (Top-Right with Bottom-Left)
q2.copyTo(q1);
tmp.copyTo(q2);
}
void fast_dft( cv::Mat &src_img, cv::Mat &real_img, cv::Mat &ima_img )
{
src_img.convertTo(src_img,CV_32FC1);
///快速傅裡葉變換/
int oph = getOptimalDFTSize(src_img.rows);
int opw = getOptimalDFTSize(src_img.cols);
Mat padded;
copyMakeBorder(src_img, padded, 0, oph-src_img.rows, 0, opw-src_img.cols,
BORDER_CONSTANT, Scalar::all(0));
Mat temp[] = {padded, Mat::zeros(padded.size(),CV_32FC1)};
Mat complexI;
merge(temp, 2, complexI);
dft(complexI, complexI); //傅裡葉變換
split(complexI, temp); //顯示頻譜圖
temp[0].copyTo(real_img);
temp[1].copyTo(ima_img);
}
int main( int argc, char *argv[] )
{
const char* filename = argc >=2 ? argv[1] : "427.tif";
Mat image = imread( filename, IMREAD_GRAYSCALE);
Mat image2 = imread( "424.tif", IMREAD_GRAYSCALE);
if( image.empty() || image2.empty())
return -1;
imshow("woman_src",image);
resize( image2, image2, cv::Size(), 0.5, 0.5);
imshow("sqrt_src",image2);
Mat woman_real, woman_imag;
Mat sqrt_real, sqrt_imag;
fast_dft( image, woman_real, woman_imag );
fast_dft( image2, sqrt_real, sqrt_imag );
Mat img_range, img_angle;
one_amplitude(woman_real, woman_imag, img_range ); //原圖像的幅值歸一化
one_angel(woman_real, woman_imag, img_angle ); //原圖像的相位歸一化
Mat woman_amp2sqrt_angle, sqrt_amp2woman_angle;
mixed_amplitude_with_phase(woman_real,woman_imag,
sqrt_real, sqrt_imag, woman_amp2sqrt_angle);
mixed_phase_with_amplitude(woman_real,woman_imag,
sqrt_real, sqrt_imag, sqrt_amp2woman_angle);
Mat amplitude, angle, amplitude_src;
magnitude(woman_real, woman_imag, amplitude); //計算原圖像幅值
phase(woman_real, woman_imag, angle); //計算原圖像相位
// cartToPolar(temp[0], temp[1],amplitude, angle);
move_to_center(amplitude); //幅值亮度集合到中心
divide(amplitude, amplitude.cols*amplitude.rows, amplitude_src );
imshow("amplitude_src",amplitude_src);
amplitude += Scalar::all(1); // switch to logarithmic scale
log(amplitude, amplitude);
normalize(amplitude, amplitude, 0, 255, NORM_MINMAX); //歸一化 友善顯示，和實際資料沒有關系
amplitude.convertTo(amplitude, CV_8U);
imshow("amplitude",amplitude);
normalize(angle, angle, 0, 255, NORM_MINMAX); //歸一化 友善顯示，和實際資料沒有關系
angle.convertTo(angle, CV_8U);
imshow("angle",angle);
//*******************************************************************
Mat inverse_amp = fourior_inverser(img_range); //傅裡葉逆變換
Mat inverse_angle= fourior_inverser(img_angle);
Mat inverse_dst1 = fourior_inverser(woman_amp2sqrt_angle);
Mat inverse_dst2 = fourior_inverser(sqrt_amp2woman_angle);
move_to_center(inverse_angle);
imshow("inverse_angle", inverse_amp);
imshow("inverse_amp", inverse_angle);
imshow("woman_amp2sqrt_angle", inverse_dst1);
imshow("sqrt_amp2woman_angle", inverse_dst2);
imwrite("phase.jpg", angle);
waitKey(0);
return 1;
}