四步移相位解包程序

clc; clear all clf X1=imread('a.bmp'); X2=imread('b.bmp'); X3=imread('c.bmp'); X4=imread('d.bmp'); figure(1); I1=imresize(X1,1,'bilinear'); imshow(I1); figure(2); I2=imresize(X2,1,'bilinear'); imshow(I2); figure(3); I3=imresize(X3,1,'bilinear'); imshow(I3); figure(4); I4=imresize(X4,1,'bilinear'); imshow(I4); [M N]=size(I1); I1=double(I1); I2=double(I2); I3=double(I3); I4=double(I4); % [M N]=size(I1); P1=rgb2gray(I1); P2=rgb2gray(I2); P3=rgb2gray(I3); P4=rgb2gray(I4); for j=1:615 for i=2:M phase(i,j)=atan2(I2(i,j)-I4(i,j),I1(i,j)-I3(i,j)); end end %四步相移法计算出相位 figure(5); imshow(mat2gray(phase)); figure(6); imshow(phase); n=zeros(M,615); %解包 n(1,1)=0; for i=2:M if abs(phase(1,i)-phase(1,i-1))<pi n(1,i)=n(1,i-1); elseif phase(1,i)-phase(1,i-1)<=-pi n(1,i)=n(1,i-1)+1; elseif phase(1,i)-phase(1,i-1)>=pi n(1,i)=n(1,i-1)-1; end end for i=2:461 for j=1:615 if abs(phase(i,j)-phase(i-1,j))<pi n(i,j)=n(i-1,j); elseif phase(i,j)-phase(i-1,j)<=-pi n(i,j)=n(i-1,j)+1; elseif phase(i,j)-phase(i-1,j)>=pi n(i,j)=n(i-1,j)-1; end end end pphase=phase+2*pi.*n; figure(7); imshow(mat2gray(pphase)); figure(8); surf(pphase(2:end-1,2:end-1));