求解四步相移法的光栅相位的matlab程序

求解四步相移法的光栅相位的matlab程序

clear all
clc
h=fspecial（‘gaussian‘91）;
%度如图 并求平均值
I1=imread（‘a01.bmp‘‘bmp‘）;
I1=filter2（hdouble（I1（::1）））;
I2=imread（‘a02.bmp‘‘bmp‘）;
I2=filter2（hdouble（I2（::1）））;
I3=imread（‘a03.bmp‘‘bmp‘）;
I3=filter2（hdouble（I3（::1）））;
I4=imread（‘a04.bmp‘‘bmp‘）;
I4=filter2（hdouble（I4（::1）））;

% I1=double（I1（201:400201:500））;
% I2=double（I2（201:400201:500））;
% I3=double（I3（201:400201:500））;
% I4=double（I4（201:400201:500））;
I1=double（I1）;
I2=double（I2）;
I3=double（I3）;
I4=double（I4）;
[widthheight]=size（I1）;
% width=200;
% height=300;
%解包
for i=1:width
    for j=1:height
        if （I1（ij）==I3（ij））
            if（I2（ij）>I4（ij））
                phase（ij）=pi/2;
            elseif（I2（ij）                phase（ij）=3*pi/2;
            else
                phase（ij）=0;
            end
        elseif （I1（ij）>I3（ij））
            if I2（ij）>=I4（ij）
                phase（ij）=atan（（I2（ij）-I4（ij））/（I1（ij）-I3（ij）））;
            else
                phase（ij）=2*pi-atan（（I4（ij）-I2（ij））/（I1（ij）-I3（ij）））; 
            end
        else
            if I2（ij）>=I4（ij）
                phase（ij）=pi-atan（（I2（ij）-I4（ij））/（I3（ij）-I1（ij）））;
            else
                phase（ij）=pi+atan（（I4（ij）-I2（ij））/（I3（ij）-I1（ij）））; 
            end
        end
    end
end
phase=unwrap（phase）;
%消除2pi突变
% midx=fix（width/2）;     
% for j=1:width
%     for m=height-1:-1:1
%         a=phase（jm）-phase（jm+1）;
%         if a>pi
%            for n=m:-1:1
%                phase（jn）=phase（jn）-2*pi;
%            end
%         end
%         if a<-pi
%            for n=m:-1:1
%             phase（jn）=phase（jn）+2*pi;
%            end
%         end
%     
%         if a==pi
%            if phase（jm）-phase（jm-1）>0
%               for n=m:-1:1
%                   phase（jn）=phase（jn）-2*pi;
%               end
%            end
%         end
%         if a==-pi
%            if phase（jm）-phase（jm-1）<0
%               for n=m:-1:1
%                   phase（jn）=phase（jn）+2*pi;
%               end
%            end
%         end  
%     end
% end
% 
% for k=1:height    
%     for m=midx:-1:1
%         c=phase（mk）-phase（m+1k）;
%         if c>pi
%            for n=m:-1:1
%                phase（nk）=phase（nk）-2*pi;
%            end
%         end
%         if c<-pi
%            for n=m:-1:1
%                phase（nk）=phase（nk）+2*pi;
%            end
%         end
%     
%         if c==pi
%            if phase（mk）-phase（mk-1）>0
%               for n=m:-1:1
%                phase（nk）=phase（nk）-2*pi;;
%               end
%            end
%         end
%         if c==-pi
%            if phase（mk）-phase（mk-1）<0
%               for n=m:-1:1
%                  phase（nk）=phase（nk）+2*pi;
%               end
%            end
%         end 
%     end
%     for m=midx:width
%         b=phase（mk）-phase（m-1k）;
%         if b>pi
%     

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

     文件       4068  2008-01-28 21:45  sibuyixiang.m

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

                 4068                    1