空間周波数の操作２（２次元フーリエスペクトル）

空間周波数の操作１では、肝心のフーリエ変換を行っておりませんでした。

以下のサンプルプログラムでは、フーリエ変換を行い、その結果得られた（パワー）スペクトルを２次元の画像として表示しています。

ハニング窓をかけていない左端の画像では、スペクトル画像に十字の輝線が描かれていることに注意してください。これは枠の効果と思われます。ハニング窓をかけた、中央、右端の画像では輝線が描かれていません。

スペクトル画像を描くとき、サンプルではフーリエ変換した結果（複素数）の絶対値を取り、さらに対数を取っています。対数を取るのは、画像の周辺部分（高空間周波数帯域）での変化を見やすくするためです。また絶対値ではなく、複素共役（conj）と掛け合わせたパワーの値を表示する方法もあります。

補足ですが、Visual Psychopysics という本（Lu & Dosher, 2014）のP51にフーリエ変換、およびフィルタリングに関する説明があり大変参考になります。

function spectrum2

close all;

clear all;

imgFileName = 'grayRabbit.jpg'; % 画像のファイル名

im1 = imread(imgFileName); % 画像のデータを読み込む。

% 画像のファイル名を分解

[pathstr, fileName, ext] = fileparts(imgFileName);

%% 必要に応じて正規化

im2 = double(im1);

imMean=mean(im2(:)); % 平均輝度

imSD=std(im2(:)); % 輝度の標準偏差

im2=(im2-imMean)./imSD; % 平均を０、標準偏差を１に。

im2 = im2.*25; % 標準偏差の変更（25は適当）

im2 = im2+128; % 平均値の変更（128は256階調の半分）

% 正規化しないときは次の行をコメントアウト

%im2 = im1;

% 正規化した画像をjpgファイルとして出力。

imwrite(uint8(im2), [fileName '-norm.jpg'], 'jpeg');

%% 座標の設定

imSize = size(im1);

halfImgSize = imSize(2)/2; % 画像が正方形、つまり、imSize(1) = imSize(2)を前提としている。

[X, Y] = meshgrid(-halfImgSize:(halfImgSize-1), halfImgSize:-1:-(halfImgSize-1)); % 行がY方向、列がX方向であることに注意。Y軸は上が正であることに注意。

[TH, R] = cart2pol(X,Y); % X,Yを極座標へ変換

%% 出力結果を表示するfigureの設定

ScreenSize = get(0, 'ScreenSize');

figure('Position', [ScreenSize(3)./4, ScreenSize(4)./4, 1000, 800]); % 画像の位置

%% 正規化した画像を表示

subplot(2,3,1);

ShowGrayScale(im2);

%% 画像データにハニング窓（四角形）をかける

absXY = max(abs(X), abs(Y));

% absXYは0からhalfImgSizeを取り得る。

% 画像の中心部では、absXYは0で、hanningの値は1となる。

% absXYが最大(halfImgSize）のとき、hanningの値は0となる。

hanning1 = 0.5 - 0.5 * cos(2*pi*(halfImgSize/2+absXY./2)/halfImgSize);

im3 = double(im2) .* hanning1;

subplot(2,3,2);

ShowGrayScale(im3);

%% 画像データにハニング窓（円形）をかける

R2 = min(R, halfImgSize); % halfImgSizeを超えない半径のデータを作る。

% R2は0からhalfImgSizeを取り得る。

% 画像の中心部では、R2は0で、hanningの値は1となる。

% R2が最大(halfImgSize）のとき、hanningの値は0となる。

hanning2 = 0.5 - 0.5 * cos(2*pi*(halfImgSize/2+R2./2)/halfImgSize);

im4 = double(im2) .* hanning2;

subplot(2,3,3);

ShowGrayScale(im4);

%% im2のFFT

Fim2 = fftshift(fft2(im2));

% スペクトル画像

subplot(2,3,4);

ShowGrayScale(log(abs(Fim2)));

%ShowGrayScale(log(Fim2.*conj(Fim2))); % こういう方法もある。

%% im3のFFT

Fim3 = fftshift(fft2(im3));

% スペクトル画像

subplot(2,3,5);

ShowGrayScale(log(abs(Fim3)));

%% im4のFFT

Fim4 = fftshift(fft2(im4));

% スペクトル画像

subplot(2,3,6);

ShowGrayScale(log(abs(Fim4)));

function ShowGrayScale(imageData)

% 画像をグレースケールで表示する。

imagesc(imageData);

axis('image');

axis('off');

colormap(gray);