How to compute Fresnel diffraction from rectangular aperture

I have used a certain code to compute Fresnel diffraction pattern from a rectangular aperture in matlab but i couldn't get the results. Can anyone please rectify the errors for me. Given below is the code used.

1. function D=rectan(amm,bmm,Wmm,smm,q0mm,lnm,t) 2. l=1nm*1e-6; 3. a=amm*sqrt(2/(1*q0mm)); b=bmm*sqrt(2/(1*q0mm)); 4. W=Wmm*sqrt(2/(1*q0mm)); s=smm*sqrt(2/(1*q0mm)); 5. Cu2=mfun('FresnelC',a:s:W+a); 6. Cu1=mfun('FresnelC',-a:s:W-a); 7. Su2=mfun('FresnelS',a:s:W+a); 8. Su1=mfun('FresnelS',-a:s:W-a); 9. A=(Cu1-Cu2).^2+(Su1-Su2).^2; 10. Cv2=mfun('FresnelC',b:s:W+b); 11. Cv1=mfun('FresnelC',-b:s:W-b); 12. Sv2=mfun('FresnelS',b:s:W+b); 13. Sv1=mfun('FresnelS',-b:s:W-b); 14. B=(Cv1-Cv2).^2+(Sv1-Sv2).^2; 15. B=B';n=size(B);B=repmat(B(:,1),1,n); 16. A=A';A=repmat(A(:,1),1,n);A=(A)'; 17. D=B*A; 18. D=t*D/max(max(D)); 19. m=2*fix(((2*W)/s)/2); 20. E=zeros(m+1,m+1); 21. for q=1:1:m/2+1; 22. for p=1:1:m/2+1; E(m/2+2-p,q+m/2)=D(p,q); end 23. end 24. for q=m+1:-1:m/2+2; 25. for p=1:1:m/2+1; E(p,-q+2+m)=E(p,q); end 26. end 27. for q=1:1:m+1; 28. for p=1:1:m/2; E(-p+2+m,q)=E(p,q); end 29. end 30. y=-W:s:W;ymm= y*sqrt((1*q0mm)/2); 31. imagesc(ymm,ymm,E,[0 1]); colormap(gray);

Answer 1

[Ahmet Cinar]

download this function:

https://www.mathworks.com/matlabcentral/fileexchange/34134-fresnel-integral

replacement for fresnel integral.

Try this:

function E = rectan(amm,bmm,Wmm,smm,Mask2Sensor,lambda,t)

a = amm*sqrt(2/(lambda*Mask2Sensor));

b = bmm *sqrt(2/(lambda* Mask2Sensor));

W = Wmm*sqrt(2/(lambda*Mask2Sensor));

s = smm*sqrt(2/(lambda*Mask2Sensor));

[Cu2, Su2] = fcs(a:s:W+a);

[Cu1, Su1] = fcs(-a:s:W-a);

[Cv2, Sv2] = fcs(b:s:W+b);

[Cv1, Sv1] = fcs(-b:s:W-b);

% Cu2 = mfun('FresnelC',a:s:W+a);

% Cu1 = mfun('FresnelC',-a:s:W-a);

% Su2 = mfun('FresnelS',a:s:W+a);

% Su1 = mfun('FresnelS',-a:s:W-a);

% Cv2 = mfun('FresnelC',b:s:W+b);

% Cv1 = mfun('FresnelC',-b:s:W-b);

% Sv2 = mfun('FresnelS',b:s:W+b);

% Sv1 = mfun('FresnelS',-b:s:W-b);

A = (Cu1-Cu2).^2+(Su1-Su2).^2;

B = (Cv1-Cv2).^2+(Sv1-Sv2).^2;

B = B';A = A';

n = length(B);B = repmat(B(:,1),1,n);

A = repmat(A(:,1),1,n);A = (A)';

D = B*A;

D = t*D/max(D(:));

m = 2*fix(((2*W)/s)/2);

E = zeros(m+1,m+1);

for q = 1:(m/2+1)

for p = 1:(m/2+1)

E(((m/2)+2-p),(q+(m/2))) = D(p,q);

end

end

for q = (m+1):-1:((m/2)+2)

for p = 1:((m/2)+1)

E(p,(-q+2+m)) = E(p, q);

end

end

for q = 1:(m+1)

for p = 1:(m/2)

E((-p+2+m),q) = E(p,q);

end

end

y = -W:s:W;

ymm = y*sqrt((lambda*Mask2Sensor)/2);

imagesc(ymm,ymm,E,[0 1]);colormap(gray);

[Ahmet Cinar]

keep everything in mm, including lambda. Run this -

rectan(0.055,0.055,11.2640/2,Sp,40,632e-6,1);

[huy dinh]

sorry bro. I'm a newbie and i can't use your code because my matlab can't found function 'Fcs'. Pls help me. Tks u so much

[Harsh Raj]

No problem, bro!