g.m

这是用来仿真四通道计算全息的程序

A=imread('1.bmp');B=imread('2.bmp');C=imread('3.bmp');D=imread('4.bmp');
a=rand(64);b=rand(64);c=rand(64);d=rand(64);
aa=exp(i*2*pi.*a);bb=exp(i*2*pi.*b);cc=exp(i*2*pi.*c);dd=exp(i*2*pi.*d);
A=double(A).*aa;B=double(B).*bb;C=double(C).*cc;D=double(D).*dd;
A0=fft2(A);B0=fft2(B);C0=fft2(C);D0=fft2(D);
Afft2=fftshift(A0);Bfft2=fftshift(B0);Cfft2=fftshift(C0);Dfft2=fftshift(D0);
A1=abs(Afft2);B1=abs(Bfft2);C1=abs(Cfft2);D1=abs(Cfft2);
p01=angle(Afft2);p02=angle(Bfft2);p03=angle(Cfft2);p04=angle(Dfft2);
P1=p01/(2*pi);P2=p02/(2*pi);P3=p03/(2*pi);P4=p04/(2*pi);
A1max=max(max(A1));B1max=max(max(B1));C1max=max(max(C1));D1max=max(max(D1));
A1=A1/A1max;B1=B1/B1max;C1=C1/C1max;D1=D1/D1max;
% W1=asin(A1)/pi;W2=asin(B1)/pi;W3=asin(C1)/pi;W4=asin(D1)/pi;
nm=1e-9;
H=632.8*nm;
z=5e+5.*nm;
for m=1:64
    for n=1:64
        A10(m,n)=A1(m,n).*exp(j*p01(m,n)).*exp(-j*pi*(m^2+n^2)/(H*z))+B1(m,n).*exp(-j*p02(m,n)).*exp(j*pi*(m^2+n^2)/(H*z));
        A11(m,n)=C1(m,n).*exp(j*p03(m,n)).*exp(-j*pi*(m^2+n^2)/(H*z))+D1(m,n).*exp(-j*p04(m,n)).*exp(j*pi*(m^2+n^2)/(H*z));
    end
end
A10=double(A10);A11=double(A11);
T01=fftshift(fft2(A10));T02=fftshift(fft2(A11));
T1=abs(T01);T2=abs(T02);
R1=angle(T01)/(2*pi);R2=angle(T02)/(2*pi);
T1max=max(max(T1));T2max=max(max(T2));
TW1=T1/T1max;TW2=T2/T2max;
TW1=asin(TW1)/pi;TW2=asin(TW1)/pi;
for h=1:64
    a=1/4;b=1/4;
    for k=1:64
        if  h+R1(h,k)+TW1(h,k)/2-a/2>h+1/2
            if  k+R2(h,k)+TW2(h,k)/2-b/2>k+1/2
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2-1,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2-1,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
            elseif   k+R2(h,k)-TW2(h,k)/2+b/2<k-1/2 
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            else  
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2-1,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2-1,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            end
        elseif h+R1(h,k)-TW1(h,k)/2+a/2<h-1/2
            if  k+R2(h,k)+TW2(h,k)/2-b/2>k+1/2
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
            elseif   k+R2(h,k)-TW2(h,k)/2+b/2<k-1/2 
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on 
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            else
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2+1,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            end
        else  h+R1(h,k)+TW1(h,k)/2-a/2<h+1/2& k+R2(h,k)-TW2(h,k)/2+b/2>k-1/2
            if  k+R2(h,k)+TW2(h,k)/2-b/2>k+1/2
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2-1,a,b],'Facecolor','k');
                hold on
            elseif  k+R2(h,k)-TW2(h,k)/2+b/2<k-1/2
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2+1,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            else
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)-TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)-TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
                rectangle('Position',[h+R1(h,k)+TW1(h,k)/2-a/2,k+R2(h,k)+TW2(h,k)/2-b/2,a,b],'Facecolor','k');
                hold on
            end
        end
    end
end