jiebaoguo.rar_四步相移

clear all close all ph=peaks(256).*3; %预设原始相位 mu=0.3; %噪声系数 noise=rands(256,256).*mu.*pi; %噪声 ph=ph+noise; %叠加有噪声的预设原始相位 figure,imshow(ph,[]); %显示原始含噪声连续相位 colormap(jet); U=exp(j*ph); wrapped_phs =angle(U); %对应的包裹相位 figure,imshow(wrapped_phs,[]); %显示包裹相位 colormap(jet); phi=wrapped_phs; %将包裹相位wrapped_phs赋值给a [m,n]=size(phi); %计算二维包裹相位的大小（行、列数） M=m-1; N=n-1; coosai=angle(U); %cosi为缠绕相位矩阵 X=zeros(M,N+1); for P=1:N+1 for i=1:M X(i,P)=mod(coosai(i+1,P)-coosai(i,P)+pi,2*pi)-pi; end end XX=zeros(M,N+2); XX(:,1)=X(:,2); XX(:,2:N+2)=X; derta_X=zeros(M+2,N+2);%derta_X 为缠绕相位X方向上的差分矩阵 derta_X(1,:)=-XX(1,:); derta_X(2:M+1,:)=XX; derta_X(M+2,:)=-XX(M,:); for jj=1:N for ii=1:M+1 Y(ii,jj)=mod(coosai(ii,jj+1)-coosai(ii,jj)+pi,2*pi)-pi; end end YY=zeros(M+1,N+2); YY(:,1)=-Y(:,1); YY(:,2:N+1)=Y; YY(:,N+2)=-Y(:,N); derta_Y=zeros(M+2,N+2);%derta_Y 为缠绕相位Y方向上的差分矩阵 derta_Y(1,:)=YY(2,:); derta_Y(2:M+2,:)=YY; %%%%%%%%%%%%%%以上是求出x和y方向的差分矩阵%%%%%%%%%%%%%%%%%%%%%% Lou=zeros(M+1,N+1); for iii=1:M+1 for jjj=1:N+1 Lou(iii,jjj)=derta_X(iii+1,jjj+1)-derta_X(iii,jjj+1)+derta_Y(iii+1,jjj+1)-derta_Y(iii+1,jjj); end end f=size(2*M,2*N); f(1:(M+1),1:(N+1))=Lou; f(1:(M+1),(N+2):2*N)=Lou(1:(M+1),N:-1:2); f((M+2):2*M,1:(N+1))=Lou(M:-1:2,1:(N+1)); f((M+2):2*M,(N+2):2*N)=Lou(M:-1:2,N:-1:2); P=fft2(f); %%%%%%%%%以上是求出P(m,n)，即lou的二维离散傅里叶变换%%%%%%%%%% for p=1:2*M for q=1:2*N k1=2*cos((p-1)*pi/(M)); k2=2*cos((q-1)*pi/(N)); PH(p,q)=P(p,q)/(k1+k2-4); end end PH(1,1)=0; phi1=ifft2(PH); phi=phi1(1:M+1,1:N+1); figure,imshow(real(phi),[]) %显示解包裹相位 colormap(jet); %为检验解包裹效果，对解包裹相位作再包裹运算 % U=exp(j.*real(phi3)); %用解包裹相位构建光场 % rewrap_ph=atan(imag(U)./real(U)); %作再包裹运算 % figure,imshow(rewrap_ph,[]) %显示再包裹结果 % colormap(jet); RMSE2=sqrt((sum(sum((real(phi)-ph).^2)))/(m*n)); RMSE2