二维傅里叶变换相位提取

二维傅里叶变换相位提取的matlab模拟程序。

%%% 二维傅里叶变换相位提取

clear all

N=512;      % 抽样点数
L=1;        % x，y方向长度
a=1;        % 背景强度
b=0.8;      % 调制度
r0=1/16;
f_x=r0*N;    % x方向载频
f_y=r0*N;    % y方向载频
delta=L/N;   % 空域采样间隔
delta_f=1/（N*delta）; % 频域采样间隔

[x y]=meshgrid（（-N/2:N/2-1）*delta）;
[fx fy]=meshgrid（（-N/2:N/2-1）*delta_f）;    

[XY]=meshgrid（1:N）;

W=0.65e-6.*（x.^2+y.^2）.^2;
PV=max（max（W））-min（min（W））;
RMS=sqrt（sum（sum（W.^2）））/N;
lamda=532e-9;
k=2*pi/lamda;
phi_real=k.*W;
i=a+b.*cos（2*pi*f_x.*x+2*pi*f_y.*y+phi_real）;
bg=sum（sum（i））/（N*N）;
  
 i1=uint8（255*i）;
 figureimagesc（i1）;
 colormap gray
 axis square
 
flag_1=（（sqrt（（X-N/2）.^2+（Y-N/2）.^2））>N/2）;
flag_2=（（sqrt（（X-N/2）.^2+（Y-N/2）.^2））flag_3=（X>（N/2-15））&（X<（N/2+15））;
flag_4=（Y>（N/2-15））&（Y<（N/2+15））;

flag=flag_1|flag_2;%|flag_3|flag_4;

flag1=not（flag）;

 i_P=i.*flag1;
 
 i1=uint8（255/2*i_P）;
 figureimagesc（i1）;
 colormap gray
 axis square

 %%%%%%%%%%%%%%%%%%%%%%%%% 三 求出载频值  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
 I_P=ft2（i_Pdelta）;
 
 [j k]=find（abs（I_P）==max（max（abs（I_P））））; 
 s0=j（1）;  t0=k（1）;                       % 找出零级频谱坐标
 
 [j k]=find（abs（I_P（1:（N/2-N/128）1:（N/2-N/128）））...
 ==max（max（abs（I_P（1:（N/2-N/128）1:（N/2-N/128））））））;  % 找出负一级频谱坐标
 s1=j（1）;t1=k（1）;

fx0=s0-s1;                                   % 求出载频值      
fy0=t0-t1;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 四 条纹延拓 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    i_E=（i-bg）.*flag1;
%   i1=uint8（255*abs（i_E））;
%   figureimagesc（i1）;
%   colormap gray
%   axis square

for m=1:10
 
    I_E=ft2（i_Edelta）;
%  
%  I1=uint8（255*abs（I_E））;
%  figureimagesc（I1）;
%  colormap gray
%  axis square
%  
 %%%%%%%%%%%%%%%%%%% 找出一级和负一频谱坐标 %%%%%%%%%%%%%%%%%%%%%%
 
[j k]=find（abs（I_E（1:（N/2-N/128）1:（N/2-N/128）））...
==max（max（abs（I_E（1:（N/2-N/128）1:（N/2-N/128））））））;  % 在第三象限内找出一级频谱坐标
s11=j（1）;  t11=k（1）; 

[j k]=find（abs（I_E（（N/2+N/128）:N（N/2+N/128）:N））...
==max（max（abs（I_E（（N/2+N/128）:N（N/2+N/128）:N）））））;  % 在第二象限内找出负一级频谱坐标
s22=j（1）+N/2+N/128-1;  t22=k（1）+N/2+N/128-1; 

%%%%%%%%%%%%%%%%%%%%%  滤波器设计 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

r1=1/16;
H1=（（sqrt（（X-s11）.^2+（Y-t11）.^2））

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  滤波  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

I_E=I_E.*H1;

 %%%%%%%%%%%%%%%%%%%%%%%%% 反傅里叶变换 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
i_E=ift2（I_Edelta_f）;


 i_E=i_E.*flag+（i-bg）.*flag1;
 
end
 
 i1=uint8（255*i_E）;
 figureimagesc（i1）;
 colormap gray
 axis square

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       2668  2012-10-17 10:48  explam-fft.m

----------- ---------  ---------- -----  ----

                 2668                    1