diff options
Diffstat (limited to 'beam_intensity_at_point_from_image.m')
| -rw-r--r-- | beam_intensity_at_point_from_image.m | 35 |
1 files changed, 35 insertions, 0 deletions
diff --git a/beam_intensity_at_point_from_image.m b/beam_intensity_at_point_from_image.m new file mode 100644 index 0000000..f305d16 --- /dev/null +++ b/beam_intensity_at_point_from_image.m @@ -0,0 +1,35 @@ +function intensity=beam_intensity_at_point_from_image(x,y,z, img, xpos, ypos, lambda) + +[Ny,Nx] = size(img); + +xim=linspace(xpos(1), xpos(end), Nx); +yim=linspace(ypos(1), ypos(end), Ny); + +[ind_x, ind_y] = meshgrid(1:Nx, 1:Ny); +[rvec_x, rvec_y] = meshgrid(xim-x, yim-y); % rvec_x and rvec_y are matrices of Ny x Nx +rvec_z=z+rvec_x*0; % to get same matrix size as rvec_x and rvec_y + +dist=sqrt(rvec_x.^2+rvec_y.^2+rvec_z.^2); + + +E_ampl = img ./dist .* exp(1i*2*pi.*dist./lambda); + +% simple case: far field - Fraunhofer limit + +E_tot= (sum(E_ampl(:))); +intensity = E_tot*conj(E_tot); + +%%%% Below is exact case: Fresnel diffraction +%kvec_x= rvec_x./dist; +%kvec_y= rvec_y./dist; +%kvec_z= rvec_z./dist; +% +%Ex= E_ampl * kvec_x; +%Ey= E_ampl * kvec_y; +%Ez= E_ampl * kvec_z; +%% summing all E field components +%Ex_tot = sum (Ex(:)); +%Ey_tot = sum (Ey(:)); +%Ez_tot = sum (Ez(:)); +%intensity = Ex_tot*conj(Ex_tot) + Ey_tot*conj(Ey_tot) + Ez_tot*conj(Ez_tot); + |