clear all
%close all
clc
clear
I=imread('1.bmp');
imshow(I)
fftimage=fftshift(fft2(double(I)));%对条纹图像I进行二维傅里叶变换
[p,q]=size(fftimage);
figure,plot(mat2gray(abs(fftimage(257,:))));%显示频谱图
arrinput=ginput;%鼠标选取高低频之间的波谷位置，点一下，按回车
arrtemp=fftimage(257,:);
[FringeMax,FringeMaxIndex]=max(arrtemp(arrinput(1,1):512));
%得到从选取点开始向右的最大值(波峰)，并记录位置
FringeMaxIndex=round(arrinput(1,1)+FringeMaxIndex-1);
%得到1级谱(基频)的波峰位置
maxvalue=max(fftimage(257,FringeMaxIndex:512));
%取出1级谱的位置。257表示第257行(中间行，谱的位置)
[m,n]=find(fftimage==maxvalue);%取出1级谱最大值位置
fringe_num=abs(n-m);%确定条纹个数
sidelobe=zeros(512);
w=8;%设定频域滤波窗口宽度，窗口为6个像素(根据情况改变)
sidelobe(:,n-w:n+w)=fftimage(:,n-w:n+w);%取出窗口提取的基频分量
sidelobe_ifft=(ifft2(ifftshift(sidelobe)));%变换基频到中心，并傅里叶反变换
sidelobe_phase=angle(sidelobe_ifft);%得到条纹图的折叠相位(相位主值)
phi=sidelobe_phase;
[m,n]=size(phi);
%phi0=peaks(200)*3;
%phi0=rgb2gray(phi0);
%figure;imshow(phi0)
%figure
%mesh(phi0);
%imagesc(phi0)
%phi=angle(exp(j*phi0));
%[m,n]=size(phi);figure;imshow(phi)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
phi1=zeros(m,n);
phi1(1:m,1:n)=phi;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
phidx=zeros(m,n);
phidy=zeros(m,n);
phidx(1:m-1,:)=phi1(2:m,:)-phi1(1:m-1,:);
%phidx=phidx-2*pi*(phidx>pi); %缠绕
%phidx=phidx+2*pi*(phidx<-pi);
phidx=mod(phidx+pi,2*pi)-pi;
phidy(:,1:n-1)=phi1(:,2:n)-phi1(:,1:n-1);
%phidy=phidy-2*pi*(phidy>pi); %缠绕
%phidy=phidy+2*pi*(phidy<-pi);
phidy=mod(phidy+pi,2*pi)-pi;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
lou=zeros(m,n);
loudx=zeros(m,n);
loudy=zeros(m,n);
loudx(2:m,:)=phidx(2:m,:)-phidx(1:m-1,:);
loudx(1,:)=phidx(1,:);
loudy(:,2:n)=phidy(:,2:n)-phidy(:,1:n-1);
loudy(:,1)=phidy(:,1);
lou=loudx+loudy;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
PP=dct2(lou);
for ii=1:m
for jj=1:n
k1=2*cos((ii-1)*pi/(m));
k2=2*cos((jj-1)*pi/(n));
PH(ii,jj)=PP(ii,jj)/(k1+k2-4);
end
end
PH(1,1)=0;
phi3=idct2(PH);
phi3=phi3(1:m,1:n);
%phi3=real(phi3);
figure;imshow(phi3)
figure
%mesh(phi3)
imagesc(phi3)
评论0