Skip to content
Snippets Groups Projects

Fix time-domain calculations of the optimal observation filter

Merged Imported Fix time-domain calculations of the optimal observation filter
ObsFiltTDFix into main
Merged Imported Achilles Kappis requested to merge ObsFiltTDFix into main
Viewing commit 09f272de
Prev Next
Show latest version
1 file
+ 116
− 132
Compare changes
  • Side-by-side
  • Inline
Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
+ 116
− 132
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
% Date: 16/07/2024 (DD/MM/YYYY)
% Date: 18/09/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -13,17 +13,35 @@
% Input
%
% gType [char/string]: The type of the source geometry. At the moment the
% available options are "Single", "Array", "Grid",
% "Cube" and "Dual". The last is an arrangement with 2
% sensors "xLen" apart that can be translated and
% rotated.
% available options are "Single", "Array", "Cube" and
% "Dual". The last is an arrangement a distance apart
% that can be translated and rotated.
%
% geoDim [numeric] (Optional): The dimensions of the geometry. This is a
% vector with three elements holding the
% dimensions of the setup in each Cartesian
% direction. Based on the "gType" some or all
% of them are used. For "gType" set either as
% "Array" or "Dual", the first value is used
% as the length of the "Array" or distance
% between the two position in "Dual". When
% "gType" is "Single" this argument is
% ignored. [Default: ones(3, 1)].
%
% xLen [numeric] (Optional): The width of the "Array" geometry, the length
% of the side of the "Grid" and "Cube geometries
% and the distance between positions in "Dual"
% geometry.
% [Default: Array/Grid/Cube - 1, Dual - 0.2].
% nSens [numeric] (Optional): This can be either a real integer value
% denoting the number of sensors in the array
% or a vector with three elements denoting
% the number of sensor along each Cartesian
% dimension. If a dimension is equal to zero
% the corresponding value is ignored. If this
% is a scalar value the used must ensure that
% it is correctly divisible over the non-zero
% dimensions of the geometry uniformly (see
% Notes for info). For the "Array" geometry,
% only the first value (if "nSens" is a
% vector) is used and for the "Single"
% geometry the argument is ignored.
% [Default: 10 * ones(3, 1)].
%
% xOff [numeric] (Optional): X-axis offset. [Default: 0].
%
@@ -31,28 +49,10 @@
%
% zOff [numeric] (Optional): Z-axis offset. [Default: 0].
%
% nSens [numeric] (Optional): The number of sensors in the array.
% For the "Grid" geometry the square root of
% "nSens" must be an integer value. For the
% "Cube" geometry the cubic root of "nSens"
% must be an integer value. This parameter is
% ignored for the "Single" geometry.
% [Default: Array/Grid - 100, Cube - 1000].
%
% orient [char/string/numeric] (Optional): The orientation of the
% geometry. For the "Single" and
% "Cube" geometries this value is
% ignored. For the "Dual"
% configuration only numerical
% values are used to declare the
% rotation of the configuration
% with respect to the local z-axis
% in degrees. For the rest of the
% configurations, the string
% values allowed are "XY", "XZ"
% and "YZ" and declare the plane
% on which the configuration will
% be placed. [Default: "XY" | 0].
% orient [numeric] (Optional): This is a real vector holding the rotation
% for the geometry in degrees, along each
% Cartesian axis. The rotations are performed
% clockwise. [Default: zeros(3, 1)].
%
% --------------------------------------------------
% Output
@@ -62,24 +62,26 @@
% number of sources and the other dimension correspond to
% the x, y and z coodrinates.
%
% vPosMesh [numeric]: The matrix has dimensions sqrt(N)xsqrt(N)x3 where N
% is the number of sources. It contains the x, y and z
% Cartesian coordinates for each position in a
% rectangular grid with dimensions sqrt(N)xsqrt(N).
% This is used only with the "Grid" geometry, otherwise
% it is empty.
% vPosMesh [numeric]: The matrix has dimensions nSens(1)xnSens(2)x3. It
% holds the x, y and z Cartesian coordinates for each
% position of the "Grid" geometry. If the geometry is
% otherwise not "Grid", this is an empty array.
%
% --------------------------------------------------
% Notes
%
% - Dependencies: rotMat3d(): To rotate sources.
% - Dependencies: rotMat3d(): To rotate the geometry.
%
% - If the number of sources are provided as a single integer, the user
% must make sure that they can be uniformly distributed in the non-zero
% dimensions of the geometry. For example, for the "Cube" geometry with
% two non-zero dimensions, the square root of the number of sensors must
% be an integer and for three non-zero dimensions its third root must be
% an integer. For the "Array" configuration, there is not an issue as the
% setup is one dimensional.
%
% - The position of the geometry is calculated like this: First the
% geometry is placed on the x-y plane, then rotated to match the
% orientation defined by the "orient" argument and finally, the offsets
% on each axis are applied.
% --------------------------------------------------
function [vPos, vPosMesh] = virtMicGeo(gType, xLen, xOff, yOff, zOff, nSens, orient)
function [vPos, vPosMesh] = virtMicGeo(gType, geoDim, nSens, xOff, yOff, zOff, orient)
% ====================================================
% Check for number of arguments
% ====================================================
@@ -91,81 +93,67 @@ function [vPos, vPosMesh] = virtMicGeo(gType, xLen, xOff, yOff, zOff, nSens, ori
% ====================================================
% Validate mandatory arguments
validateattributes(gType, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "Geometry type", 1);
validatestring(gType, ["Single", "Dual", "Array", "Grid", "Cube"], mfilename, "Geometry type", 1);
validatestring(gType, ["Single", "Dual", "Array", "Cube"], mfilename, "Geometry type", 1);
% Validate optional arguments
if nargin > 1 && ~isempty(xLen) && ~strcmpi(gType, "Single")
validateattributes(xLen, "numeric", {'scalar', 'real', 'nonnan', 'finite', 'positive'}, mfilename, "Width of the geometry, or distance between positions in Dual geometry", 2);
if nargin > 1 && ~isempty(geoDim) && ~strcmpi(gType, "Single")
validateattributes(geoDim, "numeric", {'vector', 'real', 'nonnan', 'finite', 'nonempty', 'positive'}, mfilename, "Geometry dimensions", 2);
if numel(geoDim) > 3 || numel(geoDim) == 2
error("The dimensions argument must be either a scalar or a vector with three elements");
end
if isscalar(geoDim)
geoDim = geoDim * ones(3, 1);
end
else
if strcmpi(gType, "Dual")
xLen = 0.2;
else
xLen = 1;
geoDim = ones(3, 1);
end
if nargin > 2 && ~isempty(nSens) && sum(strcmpi(gType, ["Single", "Dual"])) == 0
validateattributes(nSens, "numeric", {'vector', 'real', 'nonempty', 'nonnan', 'nonnegative', 'integer', 'finite'}, mfilename, "Number of sensors in the geometry", 3);
if numel(nSens) > 3 || numel(nSens) == 2
error("The number of sensors must be either a scalar or a vector with three elements");
end
if isscalar(nSens) && strcmpi(gType, "Cube")
if mod(nthroot(nSens, sum(geoDim ~= 0)), 1) ~= 0
error("The number of sources is cannot be divided uniformly over the non-zero dimensions");
else
nSens = nthroot(nSens, sum(geoDim ~= 0)) * double(geoDim ~= 0);
end
end
elseif strcmpi(gType, "Dual")
nSens = [2, 0, 0];
else
nSens = 10 * ones(3, 1);
end
if nargin > 2 && ~isempty(xOff)
validateattributes(xOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "X-offset", 3);
if nargin > 3 && ~isempty(xOff)
validateattributes(xOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "X-offset", 4);
else
xOff = 0;
end
if nargin > 3 && ~isempty(yOff)
validateattributes(yOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "Y-offset", 4);
if nargin > 4 && ~isempty(yOff)
validateattributes(yOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "Y-offset", 5);
else
yOff = 0;
end
if nargin > 4 && ~isempty(zOff)
validateattributes(zOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "Z-offset", 5);
if nargin > 5 && ~isempty(zOff)
validateattributes(zOff, "numeric", {'scalar', 'real', 'nonnan', 'finite'}, mfilename, "Z-offset", 6);
else
zOff = 0;
end
if nargin > 5 && ~isempty(nSens)
if ~strcmpi(gType, "Single")
validateattributes(nSens, "numeric", {'scalar', 'real', 'nonnegative', 'finite', 'nonnan'}, mfilename, "Number of sensors", 6);
end
else
if sum(strcmpi(gType, ["Array", "Grid"])) > 0
nSens = 1e2;
elseif strcmpi(gType, "Cube")
nSens = 1e3;
else
nSens = 0;
end
end
if nargin > 6 && ~isempty(orient)
if ischar(orient) || isstring(orient)
validateattributes(orient, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "Geometry orientation", 7);
validatestring(orient, ["XY", "XZ", "YZ"], mfilename, "Geometry orientation", 7);
else
if ~strcmpi(gType, "Dual")
error("Planes can be provided for the orientation of only the Array and Grid geometries.");
end
validateattributes(orient, {'numeric'}, {'scalar', 'nonempty', 'nonnan', 'finite', 'real'}, mfilename, "Geometry orientation", 7);
end
if nargin > 6 && ~isempty(orient) && ~strcmpi(gType, "Single")
validateattributes(orient, {'numeric'}, {'vector', 'nonempty', 'nonnan', 'finite', 'real', 'numel', 3}, mfilename, "Geometry rotation around the Cartesian axes", 7);
else
if strcmpi(gType, "Dual")
orient = 0;
else
orient = "XY";
end
orient = zeros(3, 1);
end
% ====================================================
% Check we have all needed parameters and they have
% acceptable values
% ====================================================
% For "Grid" geometry the square root of the number of sensors must be an integral value
if strcmpi(gType, "Grid") && mod(sqrt(nSens), 1) ~= 0
error("For the Grid geometry, the square root of the number of sensors must be an integer.");
elseif strcmpi(gType, "Cube") && mod(nthroot(nSens, 3), 1) ~= 0
error("For the Cube geometry, the cubic root of the number of sensors must be an integer (a tolerance of 1e-6 has been used for the calculation).");
end
% ====================================================
% Calculate virtual microphone positions
@@ -173,53 +161,49 @@ function [vPos, vPosMesh] = virtMicGeo(gType, xLen, xOff, yOff, zOff, nSens, ori
switch lower(gType)
case "single"
% Create a single point at offset coordinates and return
vPos = [xOff, yOff, zOff];
vPosMesh = [];
return;
vPos = zeros(1, 3);
case "dual"
vPos = [-xLen/2, 0, 0; ...
xLen/2, 0, 0];
vPos = [-geoDim(1)/2, 0, 0; ...
geoDim(1)/2, 0, 0];
case "array"
vPos = linspace(-xLen/2, xLen/2, nSens);
vPos = linspace(-geoDim(1)/2, geoDim(1)/2, nSens(1));
vPos = [vPos(:), zeros(numel(vPos), 2)];
case "grid"
vPos = linspace(-xLen/2, xLen/2, sqrt(nSens));
[x, y] = meshgrid(vPos, vPos);
vPos = [x(:), y(:), zeros(numel(x), 1)];
case "cube"
vPos = linspace(-xLen/2, xLen/2, round(nSens^(1/3)));
[x, y, z] = meshgrid(vPos, vPos, vPos);
% Make we don't get empty arrays when dimensions are 0
if nSens(1) == 0
x = 0;
else
x = linspace(-geoDim(1)/2, geoDim(1)/2, nSens(1));
end
if nSens(2) == 0
y = 0;
else
y = linspace(-geoDim(2)/2, geoDim(2)/2, nSens(2));
end
if nSens(3) == 0
z = 0;
else
z = linspace(-geoDim(3)/2, geoDim(3)/2, nSens(3));
end
[x, y, z] = meshgrid(x, y, z);
vPos = [x(:), y(:), z(:)];
otherwise
error("Oops... something went wrong here... not known geometry...!!!");
end
% Rotate
if ~strcmpi(gType, "Cube")
if isnumeric(orient)
vPos = vPos * rotMat3d(0, 0, -orient, "Degs");
elseif strcmpi(orient, "XZ")
if strcmpi(gType, "array")
vPos = vPos * rotMat3d(0, 0, 90, "Degs");
else
vPos = vPos * rotMat3d(90, 0, 0, "Degs");
end
elseif strcmpi(orient, "YZ")
vPos = vPos * rotMat3d(0, 90, 0, "Degs");
end
end
vPos = vPos * rotMat3d(-orient(1), -orient(2), -orient(3), "Degs");
% Translate
if ~strcmpi(gType, "Cube")
vPos = vPos + [xOff, yOff, zOff];
end
vPos = vPos + [xOff, yOff, zOff];
% For "Grid" and "Cube" geometry provide the coordinates in a "mesh format"
if strcmpi(gType, "Grid")
vPosMesh = reshape(vPos, sqrt(nSens), [], 3);
elseif strcmpi(gType, "Cube")
vPosMesh = reshape(vPos, round(nSens^(1/3)), round(nSens^(1/3)), [], 3);
else
% For "Cube" geometry provide the coordinates in a "mesh format"
if nargout > 1 && strcmpi(gType, "Cube")
vPosMesh = cat(3, x, y, z);
elseif nargout > 1
vPosMesh = [];
end
end
\ No newline at end of file