Code covered by the BSD License  

Highlights from
Huygens-Fresnel integral approximation, free-form apertures and rough surfaces

image thumbnail

Huygens-Fresnel integral approximation, free-form apertures and rough surfaces

by

 

Simulate wave propagation through free-form apertures, or off rough surfaces.

planeWave2Plane(lightsource, surface, options)
function surface = planeWave2Plane(lightsource, surface, options)
% planeWave2Plane - calculates the complex E field vector for every point
% on a grid if it is illuminated by a plane wavefront.
%
% http://dx.doi.org/10.1016/j.measurement.2013.05.003
% Section 3.1, Fig. 2
%
% A ......... surface grid, [xi; yi; zi]
% As ........ surface grid spacing, [dx dy]
% An ........ main surface normal
% Ad ........ minimal distance from the surface to the origin
% iMag ...... irradiance
% ir ........ displacement of the laser aperture
% in ........ laser beam direction
% iSigma .... beam width
% iLambda ... laser wavelength

% map to more convenient variables
A		= surface.grid;
As		= surface.resolution;
An		= surface.normal;
Ad		= surface.normalDistance;
iMag	= lightsource.irradiance;
ir		= lightsource.translations;
in		= lightsource.normal;
iLambda = lightsource.wavelength;
lx		= lightsource.coordinateSystem(:,1);
ly		= lightsource.coordinateSystem(:,2);

% find intersection laser beam/surface
% solve Hesse normal form: Anx*x + Any*y + Anz*z - Ad == 0
% where ir + lambda * in = [x y z]
% intersection point = aperture + lambda * propagation vector
lambda = (Ad - An' * ir) / (An' * in);
p = ir + lambda .* in;

% distance grid points - intersection point
d = A - p(:,ones(1,size(A,2)));

% project distances 'd' on the direction of propagation 'in'
dn = sum(d.*in(:,ones(1,size(d,2)))) / (in' * in);

% component of 'd' parallel to aperture, shortest distance to laser beam
da = d - dn(ones(3,1),:).*in(:,ones(1,size(d,2)));

% irradiance * planeWave(1) * exp(1i * dphi); pixel area is contained in
% huygens-fresnel integral
surface.E = sqrt(iMag) * 1 .* exp(1i * 2 * pi .* dn ./ iLambda); % (As(1)*As(2))

% compute vectors of ideal reflection for every point, mirror light
% propagation vector on surface normals
if isempty(surface.mainReflectionVectors)
	displayMessage(2, [' ' surface.name ': main reflection vectors']);
	surface.mainReflectionVectors = getReflection(in, surface.gridNormals);
end

Contact us