Commit 9aa13e8c authored by Marco Berzborn's avatar Marco Berzborn

added some more documentation to sph mimo simulation

parent 50502600
function varargout = ita_sph_mimo_error_simulation(varargin) function varargout = ita_sph_mimo_error_simulation(varargin)
%ITA_SPH_MIMO_ERROR_SIMULATION - Simulate aliasing and noise errors in %ITA_SPH_MIMO_ERROR_SIMULATION - Simulate aliasing and noise errors in
% spherical MIMO systems. % spherical MIMO systems.
% This function % This function simulates the mismatch errors from noise, sampling dispacement and aliasing
% for a MIMO system comprised of a spherical loudspeaker array and a spherical microphone
% array.
% %
% Syntax: % Syntax:
% itaResult = ita_sph_mimo_error_simulation(options) % itaResult = ita_sph_mimo_error_simulation(options)
% %
% Options (default): % Options (default):
% 'opt1' (defaultopt1) : description % 'SNR' (60) : SNR at each transducer
% 'opt2' (defaultopt1) : description % 'nRuns' (5) : calculate average of nRuns number of source-receiver orientations
% 'opt3' (defaultopt1) : description % 'dirMeasurementFile' ( [] ) : Filename for a directivity file that is to be included
%
% Directivity files need to be in hdf5 format. The name of data fields need to be 'fullDirRe' and
% 'fullDirIm' for the real and imaginary part respectively.
% %
% Example: % Example:
% audioObjOut = ita_sph_mimo_error_simulation(audioObjIn) % [individualTerms, allTerms, totalError] = ita_sph_mimo_error_simulation(sourceParams, receiverParams, opts)
% %
% See also: % See also:
% ita_toolbox_gui, ita_read, ita_write, ita_generate % ita_toolbox_gui, ita_read, ita_write, ita_generate
...@@ -26,7 +31,7 @@ function varargout = ita_sph_mimo_error_simulation(varargin) ...@@ -26,7 +31,7 @@ function varargout = ita_sph_mimo_error_simulation(varargin)
% </ITA-Toolbox> % </ITA-Toolbox>
% Author: Marco Berzborn -- Email: marco.berzborn@rwth-aachen.de % Author: Marco Berzborn -- Email: marco.berzborn@akustik.rwth-aachen.de
% Created: 29-Mar-2016 % Created: 29-Mar-2016
...@@ -47,6 +52,7 @@ sArgs = struct('pos1_source','struct',... ...@@ -47,6 +52,7 @@ sArgs = struct('pos1_source','struct',...
'dirMeasurementFile',[]); 'dirMeasurementFile',[]);
[source,receiver,sArgs] = ita_parse_arguments(sArgs,varargin); [source,receiver,sArgs] = ita_parse_arguments(sArgs,varargin);
% save all struct fields as variables for speed improvements in the parfor loop
receiverNmax = receiver.Nmax; receiverNmax = receiver.Nmax;
receiverSampling = receiver.sampling; receiverSampling = receiver.sampling;
...@@ -74,7 +80,7 @@ simWNG = sArgs.simWNG; ...@@ -74,7 +80,7 @@ simWNG = sArgs.simWNG;
yMIMOgroundTruth = zeros(sArgs.nRuns,1,numel(kVec)); yMIMOgroundTruth = zeros(sArgs.nRuns,1,numel(kVec));
% check if source and receiver sampling have a unique radius within a certain tolerance
if numel(unique(sourceSampling.r)) > 1 if numel(unique(sourceSampling.r)) > 1
uniRSort = sort(unique(sourceSampling.r)); uniRSort = sort(unique(sourceSampling.r));
if uniRSort(1) > uniRSort(end)*(1-2*eps) && uniRSort(end) < uniRSort(1)*(1+2*eps) if uniRSort(1) > uniRSort(end)*(1-2*eps) && uniRSort(end) < uniRSort(1)*(1+2*eps)
...@@ -90,8 +96,10 @@ if numel(unique(receiverSampling.r)) > 1 ...@@ -90,8 +96,10 @@ if numel(unique(receiverSampling.r)) > 1
end end
end end
% include measured source directivity
dirMeasurementFile = sArgs.dirMeasurementFile; dirMeasurementFile = sArgs.dirMeasurementFile;
%% separate error bounds
% initialize constant variables
receiverYalias = ita_sph_base(receiverSampling,receiverNmaxAlias); receiverYalias = ita_sph_base(receiverSampling,receiverNmaxAlias);
sourceYalias = ita_sph_base(sourceSampling,sourceNmaxAlias); sourceYalias = ita_sph_base(sourceSampling,sourceNmaxAlias);
sourceGalias = ita_sph_aperture_function_sla(sourceSampling,sourceNmaxAlias,sourcerMem); sourceGalias = ita_sph_aperture_function_sla(sourceSampling,sourceNmaxAlias,sourcerMem);
...@@ -100,11 +108,9 @@ sourceG = sourceGalias(:,1:(sourceNmax+1)^2); ...@@ -100,11 +108,9 @@ sourceG = sourceGalias(:,1:(sourceNmax+1)^2);
sourceY = sourceYalias(:,1:(sourceNmax+1)^2); sourceY = sourceYalias(:,1:(sourceNmax+1)^2);
receiverY = receiverYalias(:,1:(receiverNmax+1)^2); receiverY = receiverYalias(:,1:(receiverNmax+1)^2);
%% init % declare non constant variables
% sma terms
eMismatchReceiver = zeros(sArgs.nRuns,1,numel(freqVec)); eMismatchReceiver = zeros(sArgs.nRuns,1,numel(freqVec));
eAliasReceiver = zeros(sArgs.nRuns,1,numel(freqVec)); eAliasReceiver = zeros(sArgs.nRuns,1,numel(freqVec));
% sla terms
eMismatchSource = zeros(sArgs.nRuns,1,numel(freqVec)); eMismatchSource = zeros(sArgs.nRuns,1,numel(freqVec));
eAliasSource = zeros(sArgs.nRuns,1,numel(freqVec)); eAliasSource = zeros(sArgs.nRuns,1,numel(freqVec));
elapsedTime = zeros(sArgs.nRuns,1); elapsedTime = zeros(sArgs.nRuns,1);
...@@ -185,6 +191,7 @@ for idxRun = 1:sArgs.nRuns ...@@ -185,6 +191,7 @@ for idxRun = 1:sArgs.nRuns
Msource = sourceB .* (sourceG.'.*sourceY'); Msource = sourceB .* (sourceG.'.*sourceY');
end end
else else
% include measured source directivity
dirMat = h5read(dirMeasurementFile,'/dir/fullRe',[1,1,idxFreq],[sourceSampling.nPoints,(sourceNmaxAlias+1)^2,1]) +... dirMat = h5read(dirMeasurementFile,'/dir/fullRe',[1,1,idxFreq],[sourceSampling.nPoints,(sourceNmaxAlias+1)^2,1]) +...
1i* h5read(dirMeasurementFile,'/dir/fullIm',[1,1,idxFreq],[sourceSampling.nPoints,(sourceNmaxAlias+1)^2,1]); 1i* h5read(dirMeasurementFile,'/dir/fullIm',[1,1,idxFreq],[sourceSampling.nPoints,(sourceNmaxAlias+1)^2,1]);
Msource = dirMat(:,1:(sourceNmax+1)^2).'; Msource = dirMat(:,1:(sourceNmax+1)^2).';
...@@ -235,8 +242,6 @@ for idxRun = 1:sArgs.nRuns ...@@ -235,8 +242,6 @@ for idxRun = 1:sArgs.nRuns
end end
end end
% alias error needed here since antialias bf needs them and parfor
% requires the variables to be initialized
if simSMA if simSMA
if numel(unique(receiverSampling.r)) == 1 if numel(unique(receiverSampling.r)) == 1
EreceiverAlias = receiverYalias*receiverBalias; EreceiverAlias = receiverYalias*receiverBalias;
...@@ -256,7 +261,7 @@ for idxRun = 1:sArgs.nRuns ...@@ -256,7 +261,7 @@ for idxRun = 1:sArgs.nRuns
EsourceAlias = []; EsourceAlias = [];
end end
% just use random dnm for simulation % use random weighting coefficients for simulation
sourceLambda = sourceRandPattern / norm(sourceRandPattern) .* ita_sph_base(sourceLookDir,sourceNmax)' * 4*pi/(sourceNmax+1)^2; sourceLambda = sourceRandPattern / norm(sourceRandPattern) .* ita_sph_base(sourceLookDir,sourceNmax)' * 4*pi/(sourceNmax+1)^2;
receiverLambda = receiverRandPattern / norm(receiverRandPattern) .* ita_sph_base(receiverLookDir,receiverNmax)' * 4*pi/(receiverNmax+1)^2; receiverLambda = receiverRandPattern / norm(receiverRandPattern) .* ita_sph_base(receiverLookDir,receiverNmax)' * 4*pi/(receiverNmax+1)^2;
......
function varargout = ita_sph_modal_strength(varargin) function varargout = ita_sph_modal_strength(varargin)
%ITA_SPH_MODAL_STRENGTH - +++ Short Description here +++ %ITA_SPH_MODAL_STRENGTH - Calculate the modal strength function for a spherical array
% This function ++++ FILL IN INFO HERE +++ % This function calculates the modal strength for a plane wave incident
% onto a spherical microphone aray. The output is a diagonal matrix if the sampling
% has a unique radius. If not, the output is a matrix with full rank. In case the wave
% number k is given as a vector the output is of size N x M x nBins.
% %
% Syntax: % Syntax:
% audioObjOut = ita_sph_modal_strength(audioObjIn, options) % B = ita_sph_modal_strength(sampling, Nmax, k, 'rigid')
% %
% Options (default): % Options (default):
% 'opt1' (defaultopt1) : description % 'transducer' (microphone) : microphone or loudspeaker array
% 'opt2' (defaultopt1) : description % 'hankelKind' (2) : choose kind of Hankel function
% 'opt3' (defaultopt1) : description
% %
% Example: % Example:
% audioObjOut = ita_sph_modal_strength(audioObjIn) % B = ita_sph_modal_strength(sampling, Nmax, k, 'rigid')
% %
% See also: % See also:
% ita_toolbox_gui, ita_read, ita_write, ita_generate % ita_toolbox_gui, ita_read, ita_write, ita_generate
...@@ -25,7 +27,7 @@ function varargout = ita_sph_modal_strength(varargin) ...@@ -25,7 +27,7 @@ function varargout = ita_sph_modal_strength(varargin)
% </ITA-Toolbox> % </ITA-Toolbox>
% Author: Marco Berzborn -- Email: marco.berzborn@rwth-aachen.de % Author: Marco Berzborn -- Email: marco.berzborn@akustik.rwth-aachen.de
% Created: 26-Feb-2016 % Created: 26-Feb-2016
%% Initialization and Input Parsing %% Initialization and Input Parsing
......
...@@ -4,15 +4,13 @@ function varargout = ita_sph_sampling_displacement(varargin) ...@@ -4,15 +4,13 @@ function varargout = ita_sph_sampling_displacement(varargin)
% positions taken from a given sampling grid % positions taken from a given sampling grid
% %
% Syntax: % Syntax:
% audioObjOut = ita_sph_sampling_displacement(audioObjIn, options) % samplingDisplaced = ita_sph_sampling_displacement(sampling, opts)
% %
% Options (default): % Options (default):
% 'opt1' (defaultopt1) : description % 'relativeError' ([0.1,0.1,0.1]) : relative displacement in percent/100
% 'opt2' (defaultopt1) : description
% 'opt3' (defaultopt1) : description
% %
% Example: % Example:
% audioObjOut = ita_sph_sampling_displacement(audioObjIn) % samplingDisplaced = ita_sph_sampling_displacement(sampling, [0.01,0.01,0.01])
% %
% See also: % See also:
% ita_toolbox_gui, ita_read, ita_write, ita_generate % ita_toolbox_gui, ita_read, ita_write, ita_generate
...@@ -26,7 +24,7 @@ function varargout = ita_sph_sampling_displacement(varargin) ...@@ -26,7 +24,7 @@ function varargout = ita_sph_sampling_displacement(varargin)
% </ITA-Toolbox> % </ITA-Toolbox>
% Author: Marco Berzborn -- Email: marco.berzborn@rwth-aachen.de % Author: Marco Berzborn -- Email: marco.berzborn@akustik.rwth-aachen.de
% Created: 29-Mar-2016 % Created: 29-Mar-2016
...@@ -66,4 +64,4 @@ if sArgs.relativeError(:,1) == 0 ...@@ -66,4 +64,4 @@ if sArgs.relativeError(:,1) == 0
end end
varargout{1} = samplingError; varargout{1} = samplingError;
end end
\ No newline at end of file
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment