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.

plane2Plane(objA, objB, options)
function [objB options] = plane2Plane(objA, objB, options)
% plane2plane - use the Huygens-Fresnel principle to superpose phasors from
% a source plane on a destination plane
%
% objA ... input
% objB ... output

% check inputs
if isempty(objA) || isempty(objA.E)
	displayMessage(-1, 'ERROR: ''objA'' or ''objA.E'' are missing');
end
if isempty(options)
	displayMessage(-1, 'ERROR: ''options'' is missing');
end

% set wavelength
objB.wavelength = objA.wavelength;

% choose array stride, 14+2 arrays have to fit in GPU memory
if isempty(objA.stride)
	% available memory on the graphics card/computer
	if options.useGpu
		availableMemory = GPUmem();
	else
		mem = memory();
		availableMemory = mem.MemAvailableAllArrays;
	end
	
	% number of arrays used by 'huygensFresnel.m' + safety
	nr1 = options.gpuMemoryUsage(1) + options.gpuMemoryOverhead;
	nr2 = options.gpuMemoryUsage(2);
	
	% number of elements per array
	if isempty(objA.aperture)
		N1 = objA.nGrid(1) * objA.nGrid(2) * 8; % size in bytes
	else
		N1 = numel(objA.aperture) * 8;
	end
	
	if isempty(objB.aperture)
		N2 = objB.nGrid(1) * objB.nGrid(2) * 8; % size in bytes
	else
		N2 = numel(objB.aperture) * 8;
	end
	
	objA.stride = 1;
	while(availableMemory < nr1*N1/objA.stride + nr2*N2)
		objA.stride = objA.stride + 1;
	end
end

displayMessage(2, ['objA.stride = ' num2str(objA.stride)]);

% empty cell array to store filenames
options.fileList = cell(objA.stride);

% check filename
if isempty(objB.fileName)
	objB.fileName = [options.folder filesep options.title ' - ' objA.name '-' objB.name];
end

displayMessage(2, 'create object folder');
mkdir(objB.fileName);      	% objects folder

for n = 1:objA.stride
	% assemble filename
	options.fileList{n} = [objB.fileName filesep 'E' num2str(n) '.mat'];

	% execute only if no file exists
	if exist(options.fileList{n},'file') == 0 || options.overwrite == true  
		% split input surface in partitions (objA.stride)
		% convert to matrices, calculate and sava a sub-image
		displayMessage(2, ['starting sub-image (' options.fileList{n} ')']);
		
		% calculate and save fields
		E = huygensFresnel(objA, objB, options, n);
		save(options.fileList{n},'E');
	else
		displayMessage(2, ['skipping sub-image (' options.fileList{n} ')']);
	end	
end

% empty results field
objB.E = zeros([1,size(objB.grid,2)]);

% compose complete Image
displayMessage(1, 'start sub-image assembly');
for n = 1:objA.stride
	temp = load(options.fileList{n},'E');
	objB.E = objB.E + temp.E;
end

Contact us