From be1583d6dac42fd91d5f97220a6af67238bfe465 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 14 Nov 2024 16:35:16 +0000
Subject: [PATCH 01/60] Update CHANGELOG.md after rebasing with new version
---
CHANGELOG.md | 17 +++++++++++++++++
1 file changed, 17 insertions(+)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 9777709..9bbebfc 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,3 +1,20 @@
+### v0.3.4 ###
+**Signal Processing - Generic**\
+\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
+
+**Signal Processing - Array Processing**\
+\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
+\* Modify `firstOrderDMA()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
+\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
+\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
+\* Correct erroneous output calculations in `firstOrderDma()`.
+
+** Signal Processing - Measurements**\
+\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+
+-------------------------------------------
+
+
### v0.3.3 ###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
--
GitLab
From edc8aaf24282e4191db0851fd793e730177c634a Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 14 Nov 2024 16:35:32 +0000
Subject: [PATCH 02/60] Update the version in project README.md
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 1586033..6482412 100644
--- a/README.md
+++ b/README.md
@@ -87,7 +87,7 @@ The project is licensed under the ***MIT License***, which is a rather unrestric
## Versioning ##
-The project uses [semantic versioning](https://semver.org/) and the current version is **v0.3.3**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
+The project uses [semantic versioning](https://semver.org/) and the current version is **v0.3.4**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
#### **Important**
--
GitLab
From c35e82eb6a24ec19c3ef792a623820d23e3c0aa1 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:53:23 +0000
Subject: [PATCH 03/60] Change the version from v0.3.3 to v0.4 to reflect the
backward-incompatibility
---
CHANGELOG.md | 19 +------------------
1 file changed, 1 insertion(+), 18 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 9bbebfc..f27e5f7 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,21 +1,4 @@
-### v0.3.4 ###
-**Signal Processing - Generic**\
-\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
-
-**Signal Processing - Array Processing**\
-\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
-\* Modify `firstOrderDMA()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
-\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
-\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
-\* Correct erroneous output calculations in `firstOrderDma()`.
-
-** Signal Processing - Measurements**\
-\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
-
--------------------------------------------
-
-
-### v0.3.3 ###
+### v0.4 ###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
--
GitLab
From 4f7526f65939ebe80c2f0a08bf70024261e4434e Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:54:56 +0000
Subject: [PATCH 04/60] =?UTF-8?q?Change=20the=20version=20in=20the=20proje?=
=?UTF-8?q?ct=E2=80=99s=20README.md=20to=20v0.4?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 6482412..035f5b3 100644
--- a/README.md
+++ b/README.md
@@ -87,7 +87,7 @@ The project is licensed under the ***MIT License***, which is a rather unrestric
## Versioning ##
-The project uses [semantic versioning](https://semver.org/) and the current version is **v0.3.4**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
+The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.0**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
#### **Important**
--
GitLab
From 577fca9ea5f896e69063ef3e9e72938ef788df54 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:47:52 +0000
Subject: [PATCH 05/60] Fix conflict in CHANGELOG.md to rebase
---
CHANGELOG.md | 19 ++++++++++++++++++-
1 file changed, 18 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index f27e5f7..670d1a9 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,4 +1,21 @@
-### v0.4 ###
+### v0.4.1 ###
+**Signal Processing - Generic**\
+\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
+
+**Signal Processing - Array Processing**\
+\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
+\* Modify `firstOrderDMA()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
+\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
+\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
+\* Correct erroneous output calculations in `firstOrderDma()`.
+
+** Signal Processing - Measurements**\
+\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+
+-------------------------------------------
+
+
+### v0.4.0 ###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
--
GitLab
From b3760c272469fc7cad212432e047333a8ecf8f23 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:49:35 +0000
Subject: [PATCH 06/60] Update README.md to fix conflict for rebasing
---
README.md | 1 +
1 file changed, 1 insertion(+)
diff --git a/README.md b/README.md
index 035f5b3..b7d1f76 100644
--- a/README.md
+++ b/README.md
@@ -1,3 +1,4 @@
+<<<<<<< HEAD
# IN-NOVA
This repository holds the open-access codebase used within the [IN-NOVA](https://in-nova-horizon.eu/) project. The codebase is divided into various sections, each relevant to a different topic/problem. The sections will be updated and more sections will be added on an as-needed basis, since the purpose of the project is not to create a software library to tackle the scientific problems.
--
GitLab
From bee8811643f1ab5bbfe3b5f0f3d537350f438afa Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:58:14 +0000
Subject: [PATCH 07/60] Remove leftovers after conflict fixes and update
version to the correct v0.4.1 value
---
README.md | 3 +--
1 file changed, 1 insertion(+), 2 deletions(-)
diff --git a/README.md b/README.md
index b7d1f76..5eb9af2 100644
--- a/README.md
+++ b/README.md
@@ -1,4 +1,3 @@
-<<<<<<< HEAD
# IN-NOVA
This repository holds the open-access codebase used within the [IN-NOVA](https://in-nova-horizon.eu/) project. The codebase is divided into various sections, each relevant to a different topic/problem. The sections will be updated and more sections will be added on an as-needed basis, since the purpose of the project is not to create a software library to tackle the scientific problems.
@@ -88,7 +87,7 @@ The project is licensed under the ***MIT License***, which is a rather unrestric
## Versioning ##
-The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.0**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
+The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.1**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
#### **Important**
--
GitLab
From 57c9d363fbb3c43252fbb22bec92ef888245a1e1 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:16:46 +0000
Subject: [PATCH 08/60] Fix bug in rcvGeo.m returning the microphone positions
of all 'positive' and then 'negative' microphones instead of microphone pairs
---
Utilities/Geometries/MATLAB/Functions/rcvGeo.m | 12 ++++++------
1 file changed, 6 insertions(+), 6 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
index 514a470..b6a16de 100644
--- a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 14/11/2024 (DD/MM/YYYY)
+% Date: 17/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -391,27 +391,27 @@ function [omniPos, fig8Pos, triPos, boxPos, box2DPos, tetPos, fig8Vec, triVec, b
% ====================================================
% Provide a Nx3 "Figure-of-Eight" matrix
if nargout > 6
- fig8Vec = reshape(fig8Pos, 3, []);
+ fig8Vec = reshape(permute(fig8Pos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "Triangle" matrix
if nargout > 7
- triVec = reshape(triPos, 3, []);
+ triVec = reshape(permute(triPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "Box" matrix
if nargout > 8
- boxVec = reshape(boxPos, 3, []);
+ boxVec = reshape(permute(boxPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "2DBox" matrix
if nargout > 9
- box2DVec = reshape(box2DPos, 3, []);
+ box2DVec = reshape(permute(box2DPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 Tetrahedron matrix
if nargout > 10
- tetVec = reshape(tetPos, 3, []);
+ tetVec = reshape(permute(tetPos, [1, 3, 2]), 3, []);
end
end
end
\ No newline at end of file
--
GitLab
From b9ef65a7e32dcbcbed09c89c097aa85bb6d60a6c Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:24:16 +0000
Subject: [PATCH 09/60] Fix input check in firstOrderDma.m to check first
dimension of input argument instead of second
---
.../Array Processing/MATLAB/Functions/firstOrderDma.m | 6 +++---
1 file changed, 3 insertions(+), 3 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
index f31cdb6..82e4168 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 05/11/2024 (DD/MM/YYYY)
+% Date: 17/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -77,8 +77,8 @@ function [h, output] = firstOrderDma(input, freq, d, pPattern, beta)
% ====================================================
validateattributes(input, {'numeric'}, {'3d', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Input', 1);
- if mod(size(input, 2), 2) ~= 0
- error("Second dimension of 'input' parameter must have even length.")
+ if mod(size(input, 1), 2) ~= 0
+ error("First dimension of 'input' parameter must have even length.")
end
validateattributes(freq, {'numeric'}, {'real', 'nonnan', 'finite', 'nonempty', 'vector', 'numel', size(input, 3)}, mfilename, 'Frequencies', 2);
--
GitLab
From b9f73abc22cd11aa37aec2b436be4150c1b99281 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:24:51 +0000
Subject: [PATCH 10/60] Update CHANGELOG.md with fixes in rcvGeo.m and
firstOrderDma.m input argument check
---
CHANGELOG.md | 6 +++++-
1 file changed, 5 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 670d1a9..e3f2cba 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -8,10 +8,14 @@
\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
\* Correct erroneous output calculations in `firstOrderDma()`.
+\* Correct input check in `firstOrderDma()` which checked the second dimension to be even, while the first should be checked (the number of microphones).
-** Signal Processing - Measurements**\
+**Signal Processing - Measurements**\
\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+**Utilities - Geometries**\
+\* Fix bug where `rcvGeo()` returned the microphone positions of the "positive" microphones first and then the "negative". Now the microphone positions are returned in "microphone pairs" (they can be used correctly with `firstOrderDma()`).
+
-------------------------------------------
--
GitLab
From ad31eefcedb2ef3db9a7d4c018933b91d2e04449 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:42:55 +0000
Subject: [PATCH 11/60] Fix bug in rcvGeo.m throwing error when asking for box
and box2d vector positions, after the last fix
---
Utilities/Geometries/MATLAB/Functions/rcvGeo.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
index b6a16de..ff86c5d 100644
--- a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
@@ -401,12 +401,12 @@ function [omniPos, fig8Pos, triPos, boxPos, box2DPos, tetPos, fig8Vec, triVec, b
% Provide a Nx3 "Box" matrix
if nargout > 8
- boxVec = reshape(permute(boxPos, [1, 3, 2]), 3, []);
+ boxVec = reshape(permute(boxPos, [1, 3, 4, 2]), 3, []);
end
% Provide a Nx3 "2DBox" matrix
if nargout > 9
- box2DVec = reshape(permute(box2DPos, [1, 3, 2]), 3, []);
+ box2DVec = reshape(permute(box2DPos, [1, 3, 4, 2]), 3, []);
end
% Provide a Nx3 Tetrahedron matrix
--
GitLab
From 447a317f81177eef71acc62a529dd6992a714b5e Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:53:23 +0000
Subject: [PATCH 12/60] Change the version from v0.3.3 to v0.4 to reflect the
backward-incompatibility
---
CHANGELOG.md | 1 +
1 file changed, 1 insertion(+)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 09a9cee..972202f 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,3 +1,4 @@
+<<<<<<< HEAD
### v0.4.1 ###
** Optimisation - Memetic Algorithm**\
\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
--
GitLab
From 1774996c1745d0eee1b5888fdd440ab4d6c6289b Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 15 Nov 2024 21:47:52 +0000
Subject: [PATCH 13/60] Fix conflict in CHANGELOG.md to rebase
---
CHANGELOG.md | 18 +++++++++++++++++-
1 file changed, 17 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 972202f..c84bc50 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,4 +1,20 @@
-<<<<<<< HEAD
+### v0.4.2 ###
+**Signal Processing - Generic**\
+\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
+
+**Signal Processing - Array Processing**\
+\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
+\* Modify `firstOrderDMA()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
+\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
+\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
+\* Correct erroneous output calculations in `firstOrderDma()`.
+
+** Signal Processing - Measurements**\
+\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+
+-------------------------------------------
+
+
### v0.4.1 ###
** Optimisation - Memetic Algorithm**\
\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
--
GitLab
From 725e4b00d41c9189e60c4140fa5ba7444722acb8 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 21 Nov 2024 18:55:00 +0000
Subject: [PATCH 14/60] Update README.md with new version after rebase
---
README.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/README.md b/README.md
index 4aac359..2dc6e7f 100644
--- a/README.md
+++ b/README.md
@@ -87,7 +87,7 @@ The project is licensed under the ***MIT License***, which is a rather unrestric
## Versioning ##
-The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.1**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
+The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.2**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
#### **Important**
--
GitLab
From 2533c7cbf8c1ebefdfddb3b11dc3c3c112f82753 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:16:46 +0000
Subject: [PATCH 15/60] Fix bug in rcvGeo.m returning the microphone positions
of all 'positive' and then 'negative' microphones instead of microphone pairs
---
Utilities/Geometries/MATLAB/Functions/rcvGeo.m | 12 ++++++------
1 file changed, 6 insertions(+), 6 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
index 514a470..b6a16de 100644
--- a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 14/11/2024 (DD/MM/YYYY)
+% Date: 17/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -391,27 +391,27 @@ function [omniPos, fig8Pos, triPos, boxPos, box2DPos, tetPos, fig8Vec, triVec, b
% ====================================================
% Provide a Nx3 "Figure-of-Eight" matrix
if nargout > 6
- fig8Vec = reshape(fig8Pos, 3, []);
+ fig8Vec = reshape(permute(fig8Pos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "Triangle" matrix
if nargout > 7
- triVec = reshape(triPos, 3, []);
+ triVec = reshape(permute(triPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "Box" matrix
if nargout > 8
- boxVec = reshape(boxPos, 3, []);
+ boxVec = reshape(permute(boxPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 "2DBox" matrix
if nargout > 9
- box2DVec = reshape(box2DPos, 3, []);
+ box2DVec = reshape(permute(box2DPos, [1, 3, 2]), 3, []);
end
% Provide a Nx3 Tetrahedron matrix
if nargout > 10
- tetVec = reshape(tetPos, 3, []);
+ tetVec = reshape(permute(tetPos, [1, 3, 2]), 3, []);
end
end
end
\ No newline at end of file
--
GitLab
From 7b5f19055081dcb400be5817a01f795aaf4b8372 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:24:16 +0000
Subject: [PATCH 16/60] Fix input check in firstOrderDma.m to check first
dimension of input argument instead of second
---
.../Array Processing/MATLAB/Functions/firstOrderDma.m | 6 +++---
1 file changed, 3 insertions(+), 3 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
index f31cdb6..82e4168 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 05/11/2024 (DD/MM/YYYY)
+% Date: 17/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -77,8 +77,8 @@ function [h, output] = firstOrderDma(input, freq, d, pPattern, beta)
% ====================================================
validateattributes(input, {'numeric'}, {'3d', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Input', 1);
- if mod(size(input, 2), 2) ~= 0
- error("Second dimension of 'input' parameter must have even length.")
+ if mod(size(input, 1), 2) ~= 0
+ error("First dimension of 'input' parameter must have even length.")
end
validateattributes(freq, {'numeric'}, {'real', 'nonnan', 'finite', 'nonempty', 'vector', 'numel', size(input, 3)}, mfilename, 'Frequencies', 2);
--
GitLab
From e8c6e69a1bf7a6524397e7ab3a1e75bb80b42f33 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 17 Nov 2024 17:42:55 +0000
Subject: [PATCH 17/60] Fix bug in rcvGeo.m throwing error when asking for box
and box2d vector positions, after the last fix
---
Utilities/Geometries/MATLAB/Functions/rcvGeo.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
index b6a16de..ff86c5d 100644
--- a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
@@ -401,12 +401,12 @@ function [omniPos, fig8Pos, triPos, boxPos, box2DPos, tetPos, fig8Vec, triVec, b
% Provide a Nx3 "Box" matrix
if nargout > 8
- boxVec = reshape(permute(boxPos, [1, 3, 2]), 3, []);
+ boxVec = reshape(permute(boxPos, [1, 3, 4, 2]), 3, []);
end
% Provide a Nx3 "2DBox" matrix
if nargout > 9
- box2DVec = reshape(permute(box2DPos, [1, 3, 2]), 3, []);
+ box2DVec = reshape(permute(box2DPos, [1, 3, 4, 2]), 3, []);
end
% Provide a Nx3 Tetrahedron matrix
--
GitLab
From 90130129a652cd068f228be81cc46256cb16ca1c Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Mon, 25 Nov 2024 14:14:20 +0000
Subject: [PATCH 18/60] Remove uncessecary checks for parallel for-loop in
locSearch.m
---
.../Memetic Algorithm/MATLAB/Functions/locSearch.m | 14 ++++----------
1 file changed, 4 insertions(+), 10 deletions(-)
diff --git a/Optimisation/Memetic Algorithm/MATLAB/Functions/locSearch.m b/Optimisation/Memetic Algorithm/MATLAB/Functions/locSearch.m
index 9d7a8b7..9bf65a3 100644
--- a/Optimisation/Memetic Algorithm/MATLAB/Functions/locSearch.m
+++ b/Optimisation/Memetic Algorithm/MATLAB/Functions/locSearch.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 30/11/2023 (DD/MM/YYYY)
+% Date: 25/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -161,13 +161,6 @@ function optimPop = locSearch(pop, fitVals, nNeighbourPts, selMeth, stepSize, nI
% ====================================================
% Perform the local search
% ====================================================
- % Handle parallel execution
- if locOptParExec
- nWorkers = gcp().NumWorkers;
- else
- nWorkers = 0;
- end
-
% Parameters
rngFactor = stepSize(2, :) - stepSize(1, :); % Factor determining the range of the values
sftFactor = stepSize(1, :); % Factor shifting the mean of the values
@@ -177,7 +170,7 @@ function optimPop = locSearch(pop, fitVals, nNeighbourPts, selMeth, stepSize, nI
optimPop = zeros(size(pop));
% Go through the population
- parfor (popIdx = 1:size(pop, 1), nWorkers)
+ parfor popIdx = 1:size(pop, 1)
optimPop(popIdx, :) = pop(popIdx, :);
% Loop conditions
@@ -230,4 +223,5 @@ function optimPop = locSearch(pop, fitVals, nNeighbourPts, selMeth, stepSize, nI
end
end
end
- end
\ No newline at end of file
+ end
+end
\ No newline at end of file
--
GitLab
From 165517f0f7671a2bb69892d2c39cdf51098573a9 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Mon, 25 Nov 2024 14:19:52 +0000
Subject: [PATCH 19/60] Update CHANGELOG.md with changes in the memetic
algorithm
---
CHANGELOG.md | 3 +++
1 file changed, 3 insertions(+)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index c84bc50..c8d9087 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -12,6 +12,9 @@
** Signal Processing - Measurements**\
\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+** Optimisation - Memetic Algorithm**\
+\* Removed redundant code to handle serial for-loop if parfor is not available.
+
-------------------------------------------
--
GitLab
From 27d19adb89d5548a5601fd43d7d56071216c9574 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Tue, 26 Nov 2024 18:59:02 +0000
Subject: [PATCH 20/60] Add function to calculate the sound field due to line
sources in the frequency domain (tonal)
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 132 +++++++++++++++++++
1 file changed, 132 insertions(+)
create mode 100644 Sound Fields/MATLAB/Functions/lineSrcField.m
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
new file mode 100644
index 0000000..a1373c0
--- /dev/null
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -0,0 +1,132 @@
+%% Calculate wavefield generated by line sources
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 26/11/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the wavefield generated by line sources at
+% specified locations.
+% --------------------------------------------------
+% Input
+%
+% sPos [numeric]: The position of the sources in Cartesian coordinates. The
+% variable should be a matrix with dimensions 3xN where N
+% is the number of sources and the rows represent the x, y
+% and z coordinates respectively. IMPORTANT: The z
+% components are set to 0 internally to ensure a 2D
+% problem.
+%
+% rPos [numeric]: The position of the receivers in Cartesian coordinates.
+% The variable should be a matrix with dimensions 3xM where
+% M is the number of receivers and the rows represent the
+% x, y and z coordinates respectively. IMPORTANT: The z
+% components are set to 0 internally to ensure a 2D
+% problem.
+%
+% freqs [numeric]: The frequencies for which to calculate the wave field.
+% It must be a scalar or vector.
+%
+% Q [numeric] (Optional): The source strength of the line sources. This
+% variable must be either a scalar corresponding to
+% all sources having the same strength, or a matrix
+% with dimensions NxK, where K is the number of
+% frequencies of interest. The strengths can be
+% complex numbers. [Default: 1].
+%
+% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
+%
+% rho [numeric (Optional): The density of the medium in kg/(m^3).
+% [Default: 1.204].
+%
+% --------------------------------------------------
+% Output
+%
+% waveField [numeric]: The complex pressures at the position of the
+% receivers. The variable is a matrix of dimensions
+% MxNxF where M is the number of receivers, N the
+% number of sources and F the number of frequencies.
+%
+% --------------------------------------------------
+% Notes
+%
+% Dependencies: - twoPtDist(): To calculate the distances between sources
+% and receivers.
+% --------------------------------------------------
+function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 6);
+ nargoutchk(0, 1);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(sPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Source positions", 1);
+ validateattributes(rPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Receiver positions", 2);
+ validateattributes(f, "numeric", {'vector', 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Frequencies", 3);
+
+ % Validate optional arguments (and set their values)
+ if nargin > 5 && ~isempty(rho)
+ validateattributes(rho, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Medium density", 6);
+ else
+ rho = 1.204;
+ end
+
+ if nargin > 4 && ~isempty(c)
+ validateattributes(c, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Speed of propagation", 5);
+ else
+ c = 343;
+ end
+
+ if nargin > 3 && ~isempty(Q)
+ validateattributes(Q, "numeric", {'2d', 'finite', 'nonempty', 'nonnan'}, mfilename, "Source strength(s)", 4);
+
+ % Make sure Q has the correct dimensions
+ if ~isscalar(Q) && sum((size(Q) ~= [size(sPos, 2), length(f)])) > 0
+ error("Source strength(s) must be either a scalar or its length must match the number of sources.");
+ end
+ else
+ Q = 1;
+ end
+
+
+ % ====================================================
+ % "Condition" arguments
+ % ====================================================
+ % Make sure frequencies is a column vector
+ f = f(:);
+
+ % Make sure we work in 2D (z components are zero)
+ sPos(3, :) = 0;
+ rPos(3, :) = 0;
+
+
+ % =============================================
+ % Calculate variables
+ % ============================================
+ s2rDist = twoPtDist(sPos, rPos).';
+
+ % Create a vector with appropriate dimensions for the source strength(s)
+ if isscalar(Q)
+ Q = ones(size(sPos, 2), length(f)) * Q;
+ end
+
+
+ % =============================================
+ % Calculate variables
+ % ============================================
+ % Go through the frequencies
+ for fIdx = numel(f):-1:1
+ % Calculate parameters
+ w = 2 * pi * f(fIdx); % Radial frequency
+ k = w/c; % Wavenumber
+
+ % Calculate complex pressure
+ waveField(:, :, fIdx) = (-1j * w * rho * Q(:, fIdx).')/(2 * pi) .* besselh(0, 1, k * s2rDist);
+ end
+end
\ No newline at end of file
--
GitLab
From ce9c9e1e962344785028847d156ad2770c293135 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Tue, 26 Nov 2024 19:02:47 +0000
Subject: [PATCH 21/60] Update CHANGELOG.md with the addition of the
lineSrcField() function to calculate sound field due to line sources
---
CHANGELOG.md | 39 +++++++++++++++++++++------------------
1 file changed, 21 insertions(+), 18 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index c8d9087..3c9e3a7 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,4 +1,4 @@
-### v0.4.2 ###
+###v0.4.2###
**Signal Processing - Generic**\
\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
@@ -9,16 +9,19 @@
\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
\* Correct erroneous output calculations in `firstOrderDma()`.
-** Signal Processing - Measurements**\
+**Signal Processing - Measurements**\
\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
-** Optimisation - Memetic Algorithm**\
+**Optimisation - Memetic Algorithm**\
\* Removed redundant code to handle serial for-loop if parfor is not available.
+**Sound Fields**\
+\+ Add function (`lineSrcField()`) to calculate the sound field due to line sources in the frequency domain (tonal).
+
-------------------------------------------
-### v0.4.1 ###
+###v0.4.1###
** Optimisation - Memetic Algorithm**\
\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
\* Update `selParents()` to use the parallel for-loop (`parfor`) whenever it is available.
@@ -26,7 +29,7 @@
-------------------------------------------
-### v0.4 ###
+###v0.4###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
@@ -37,7 +40,7 @@
-------------------------------------------
-### v0.3.2 ###
+###v0.3.2###
**Signal Processing - Generic**\
\+ Add option to use negative delay values with the `winSincFracDel()` function.
@@ -56,14 +59,14 @@
-------------------------------------------
-### v0.3.1 ###
+###v0.3.1###
**Virtual Sensing - Remote Microphone Technique**\
\+ Add the option to add "noise" to the monitoring microphone signals in the time domain estimation of the optimal observation filters in `obsFiltTD()` to replicate the feature of `obsFilt()`. This change **breaks** backwards compatibility.
-------------------------------------------
-### v0.3.0 ###
+###v0.3.0###
**General**\
\* Convert the project to accept column position and direction vectors as input. The affected functions are:
- `srcGeo()`
@@ -84,7 +87,7 @@
-------------------------------------------
-### v0.2.8 ###
+###v0.2.8###
**Utilities - Generic**\
\+ Add the option to pick the order the rotations will be performed in the `rotMat3d()` function. This results is a **backward incompatible** change.
@@ -94,7 +97,7 @@
-------------------------------------------
-### v0.2.7 ###
+###v0.2.7###
**Signal Processing - Measurements**\
\+ Introduction of the *Measurements* "topic" in the Signal Processing part of the codebase.\
\+ Function to estimate the impulse responses from sweep measurements.
@@ -102,14 +105,14 @@
-------------------------------------------
-### v0.2.6 ###
+###v0.2.6###
**Signal Processing - Generic**\
\+ Add function to perform fractional octave smoothing of spectra.
-------------------------------------------
-### v0.2.5 ###
+###v0.2.5###
**Utilities - Geometries**\
\* Combine translation offsets to a single vector input argument for the `rcvGeo()` MATLAB function.\
\* Combine translation offsets to a single vector input argument for th `virtMicGeo()` MATLAB function.\
@@ -118,7 +121,7 @@
-------------------------------------------
-### v0.2.4 ###
+###v0.2.4###
**Virtual Sensing**\
\* Corrected time-domain observation filter calculations.\
\* Removed scaling from the cross-correlation calculations in the time-domain optimal observation filters.\
@@ -133,7 +136,7 @@
-------------------------------------------
-### v0.2.3 ###
+###v0.2.3###
**Virtual Sensing**\
\+ Added function to estimate the observation filters in the time-domain.\
\+ Added function to perform estimation with observation filters in the time-domain.
@@ -148,21 +151,21 @@
-------------------------------------------
-### v0.2.2 ###
+###v0.2.2###
**Virtual Sensing**\
\* Update the observation filter and estimation with observation filter functions with better input argument checking and changed the name of variables to match those in the literature.
-------------------------------------------
-### v0.2.1 ###
+###v0.2.1###
**Virtual Sensing**\
\* Fix a bug where noise was added to the power spectral density matrix for the optimal observation filter calculations in the noiseless case.
-------------------------------------------
-### v0.2.0 ###
+###v0.2.0###
**Utilities**\
\+ Added a bisection method for single-valued functions.\
\* Improved rotation matrix calculation function with degrees and radian calculations available.\
@@ -213,5 +216,5 @@
-------------------------------------------
-### v0.1.0 ###
+###v0.1.0###
This is the initial version of the project.
--
GitLab
From 84d723a8d19bce7b85f944e32861d6e7c4dc4397 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Wed, 27 Nov 2024 00:02:00 +0000
Subject: [PATCH 22/60] Fix a typo in ptSrcField.m header
---
Sound Fields/MATLAB/Functions/ptSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/ptSrcField.m b/Sound Fields/MATLAB/Functions/ptSrcField.m
index fd110c0..67aca00 100644
--- a/Sound Fields/MATLAB/Functions/ptSrcField.m
+++ b/Sound Fields/MATLAB/Functions/ptSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 29/09/2024 (DD/MM/YYYY)
+% Date: 26/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -25,7 +25,7 @@
% freqs [numeric]: The frequencies for which to calculate the wave field.
% It must be a scalar or vector.
%
-% Q [numeric] (Optional): The source strength of the plane waves. This
+% Q [numeric] (Optional): The source strength of the point sources. This
% variable must be either a scalar corresponding to
% all sources having the same strength, or a matrix
% with dimensions NxK, where K is the number of
--
GitLab
From 3428746b12cb614f6d43c8e12d4ca68e9140330a Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Wed, 27 Nov 2024 17:45:49 +0000
Subject: [PATCH 23/60] Fix the phase term in lineSrcField.m
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 8 ++++----
1 file changed, 4 insertions(+), 4 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
index a1373c0..c6edbab 100644
--- a/Sound Fields/MATLAB/Functions/lineSrcField.m
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 26/11/2024 (DD/MM/YYYY)
+% Date: 27/11/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -123,10 +123,10 @@ function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
% Go through the frequencies
for fIdx = numel(f):-1:1
% Calculate parameters
- w = 2 * pi * f(fIdx); % Radial frequency
- k = w/c; % Wavenumber
+ w = 2 * pi * f(fIdx);
+ kr = w * s2rDist/c;
% Calculate complex pressure
- waveField(:, :, fIdx) = (-1j * w * rho * Q(:, fIdx).')/(2 * pi) .* besselh(0, 1, k * s2rDist);
+ waveField(:, :, fIdx) = (1j * w * rho * Q(:, fIdx).')/(2 * pi) .* besselh(0, 1, kr);
end
end
\ No newline at end of file
--
GitLab
From 34e5e24f78fec52a90851812f72aa2bf888a7650 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Tue, 3 Dec 2024 15:20:23 +0000
Subject: [PATCH 24/60] Correct equation giving the complex pressure in
lineSrcField.m
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
index c6edbab..83c2ac5 100644
--- a/Sound Fields/MATLAB/Functions/lineSrcField.m
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 27/11/2024 (DD/MM/YYYY)
+% Date: 03/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -127,6 +127,6 @@ function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
kr = w * s2rDist/c;
% Calculate complex pressure
- waveField(:, :, fIdx) = (1j * w * rho * Q(:, fIdx).')/(2 * pi) .* besselh(0, 1, kr);
+ waveField(:, :, fIdx) = 1i * rho * c * Q(:, fIdx).' .* besselh(0, 2, kr)/4;
end
end
\ No newline at end of file
--
GitLab
From 5b4ce8a54ca88b54908c904e52f5a17281e452b5 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 6 Dec 2024 10:43:47 +0000
Subject: [PATCH 25/60] Update parameter list and checks in firstOrderDma.m
---
.../MATLAB/Functions/firstOrderDma.m | 37 +++++++++++--------
1 file changed, 22 insertions(+), 15 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
index 82e4168..6377e7b 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 17/11/2024 (DD/MM/YYYY)
+% Date: 04/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -64,39 +64,44 @@
% Notes
%
% --------------------------------------------------
-function [h, output] = firstOrderDma(input, freq, d, pPattern, beta)
+function [h, output] = firstOrderDma(freq, d, pPattern, beta, input)
% ====================================================
% Check for number of arguments
% ====================================================
- narginchk(3, 5);
+ narginchk(2, 5);
nargoutchk(0, 2);
% ====================================================
% Validate input arguments
% ====================================================
- validateattributes(input, {'numeric'}, {'3d', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Input', 1);
+ validateattributes(freq, {'numeric'}, {'real', 'nonnan', 'finite', 'nonempty', 'vector'}, mfilename, 'Frequencies', 1);
+ validateattributes(d, {'numeric'}, {'scalar', 'real', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Inter-element distance', 2);
- if mod(size(input, 1), 2) ~= 0
- error("First dimension of 'input' parameter must have even length.")
- end
+ % Check input
+ if nargin > 4 && ~isempty(input)
+ validateattributes(input, {'numeric'}, {'3d', 'nonnan', 'finite', 'nonempty', 'size', [NaN, NaN, length(freq)]}, mfilename, 'Input', 5);
- validateattributes(freq, {'numeric'}, {'real', 'nonnan', 'finite', 'nonempty', 'vector', 'numel', size(input, 3)}, mfilename, 'Frequencies', 2);
- validateattributes(d, {'numeric'}, {'scalar', 'real', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Inter-element distance', 3);
+ if mod(size(input, 1), 2) ~= 0
+ error("First dimension of 'input' parameter must have even length.")
+ end
+ else
+ input = [];
+ end
% Check beta
- if nargin > 4
- validateattributes(beta, {'numeric'}, {'scalar', 'real', 'finite', 'nonempty', 'nonnan', '<=', 1, '>=', 0}, mfilename, "Beta", 5)
+ if nargin > 3 && ~isempty(beta)
+ validateattributes(beta, {'numeric'}, {'scalar', 'real', 'finite', 'nonempty', 'nonnan', '<=', 1, '>=', 0}, mfilename, "Beta", 4)
else
beta = 0;
end
% Check alpha (angle of null)
- if nargin > 3
+ if nargin > 2 && ~isempty(pPattern);
if isstring(pPattern) || ischar(pPattern)
- validateattributes(pPattern, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, 'Polar pattern', 4);
+ validateattributes(pPattern, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, 'Polar pattern', 3);
elseif isnumeric(pPattern)
- validateattributes(pPattern, {'numeric'}, {'scalar', 'real', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Angle of null', 4);
+ validateattributes(pPattern, {'numeric'}, {'scalar', 'real', 'nonnan', 'finite', 'nonempty'}, mfilename, 'Angle of null', 3);
end
else
pPattern = "dipole";
@@ -151,7 +156,7 @@ function [h, output] = firstOrderDma(input, freq, d, pPattern, beta)
% ====================================================
% Calculate array output
% ====================================================
- if nargout > 1
+ if nargout > 1 && ~isempty(input)
% Go through the frequencies
for freqIdx = length(freq):-1:1
for pairIdx = size(input, 1)/2:-1:1
@@ -159,6 +164,8 @@ function [h, output] = firstOrderDma(input, freq, d, pPattern, beta)
output(pairIdx, :, freqIdx) = h(:, freqIdx)' * input(pairIdx * 2 - 1:pairIdx * 2, :, freqIdx);
end
end
+ else
+ output = [];
end
end
--
GitLab
From 1ff66dbfee2d890edfb5dc476354a22b363c4cbe Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 6 Dec 2024 10:47:47 +0000
Subject: [PATCH 26/60] Add backward-incompatible changes of firstOrderDma.m
and change the version number to 0.5.0 in CHANGELOG.md
---
CHANGELOG.md | 7 ++++---
1 file changed, 4 insertions(+), 3 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 3c9e3a7..1bc8f7d 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,13 +1,14 @@
-###v0.4.2###
+###v0.5.0###
**Signal Processing - Generic**\
\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
**Signal Processing - Array Processing**\
\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
-\* Modify `firstOrderDMA()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
+\* Modify `firstOrderDma()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
\* Correct erroneous output calculations in `firstOrderDma()`.
+\* Moved the `input` input argument to the last of the input arguments list and updated the input argument checks. This change **breaks backward compatibility**.
**Signal Processing - Measurements**\
\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
@@ -29,7 +30,7 @@
-------------------------------------------
-###v0.4###
+###v0.4.0###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
--
GitLab
From c2475bd4bd2fe23170c6540ec0b158c252126e47 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 6 Dec 2024 11:01:17 +0000
Subject: [PATCH 27/60] Update README.md with changes in the version, project
information, drop of Julia implementation, thought of porting to Python and
minor changes to wording and section titles
---
README.md | 18 ++++++++++--------
1 file changed, 10 insertions(+), 8 deletions(-)
diff --git a/README.md b/README.md
index 2dc6e7f..6b81538 100644
--- a/README.md
+++ b/README.md
@@ -1,5 +1,5 @@
# IN-NOVA
-This repository holds the open-access codebase used within the [IN-NOVA](https://in-nova-horizon.eu/) project. The codebase is divided into various sections, each relevant to a different topic/problem. The sections will be updated and more sections will be added on an as-needed basis, since the purpose of the project is not to create a software library to tackle the scientific problems.
+This repository holds the open-access codebase used within the "Directional microphone arrays for remote microphone virtual sensing" sub-project, under the [IN-NOVA](https://in-nova-horizon.eu/) project . The codebase is divided into various sections, each relevant to a different topic/problem. The sections will be updated and more sections will be added on an as-needed basis, since the purpose of the project is not to create a software library to tackle the scientific problems.
The codebase under each section may or may not include code implementations in different languages but no plan exists to port each implementation in all available programming languages. More information about the codebase, the available implementations and topics/fields can be found in the [Documentation](https://gitlab.com/in-nova/in-nova/-/tree/main/Documentation/latex?ref_type=heads) of the project and its [Wiki](https://gitlab.com/in-nova/in-nova/-/wikis/home).
@@ -24,6 +24,8 @@ The current structure of the codebase is shown below. More information can be fo
- Memetic Algorithm
- [Signal Processing](https://gitlab.com/in-nova/in-nova/-/tree/main/Signal%20Processing?ref_type=heads)
- Array Processing
+ - Frequency domain
+ - Time domain
- Generic
- Measurements
- [Sound fields](https://gitlab.com/in-nova/in-nova/-/tree/main/Sound%20Fields?ref_type=heads)
@@ -44,11 +46,11 @@ If you find any bugs, other issues and/or inconsistencies you can [e-mail](mailt
## Roadmap
-The repository will be updated on an as-needed basis.
+There is no development roadmap for this project. The codebase is developed to cover the needs of the project as they show up, so the repository will be updated on an as-needed basis.
Although the repository is not considered to be *under construction* anymore, there are significant gaps in its structure. At the moment, the main focus is placed on updating code with new features and bug fixes and documenting ideas and useful features for the future (mainly as issues).
-As a parallel slow-moving project, porting the codebase to [Julia](https://julialang.org/) has also started, aiming to provide full functionality in the future. Although Julia and MATLAB were designed to tackle similar (if not the same) problems, they are quite different. The facilities and intricacies of Julia dictate the need to dig deeper into the language in order to provide an efficient, maintainable and intuitive codebase that will parallel the current MATLAB implementation in functionality. More time must be spent in the design phase for the Julia implementations, which will, most probably, keep the progress to a slow pace. Although it is a long shot, the target is for the Julia codebase to reach the MATLAB implementation's state and then move forward in tandem.
+The port of the codebase to [Julia](https://julialang.org/) started a while ago, but it has currently been abandonded. The development requires considerable effort to provide a consistent and performant codebase. Instead of spending time to learn Julia to efficiently port the code, the effort will be put to improve the current codebase and possibly port it to [Python](https://www.python.org/), as it is one of the main languages used in computational sciences and signal processing.
## Contributing
@@ -70,16 +72,16 @@ To try and maintain a coherent style in the codebase you are encouraged to use t
Irrespective of style, it is highly advisable to name your variables in a meaningful manner and keep them as short as possible while still communicating their purpose. Comments are highly recommended in your code, especially where they can facilitate a better understanding of the code and the intent of the developer. Over-commented code is better than under-commented code.
-## Authors and acknowledgment
+## Contributors and acknowledgment
-#### Authors
+#### Contributors
The main contributor is Achilles Kappis.
#### Acknowledgments
-Part of the `ptsOnSphere()` MATLAB/Octave function is an adaptation of a function provided by Prof. Jordan Cheer of the Institute of Sound and Vibration Research (ISVR), University of Southampton.
+Part of the `ptsOnSphere()` MATLAB/Octave function is an adaptation of a function provided by Prof. Jordan Cheer of the Institute of Sound & Vibration Research (ISVR), University of Southampton.
## License
@@ -87,9 +89,9 @@ The project is licensed under the ***MIT License***, which is a rather unrestric
## Versioning ##
-The project uses [semantic versioning](https://semver.org/) and the current version is **v0.4.2**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
+The project uses [semantic versioning](https://semver.org/) and the current version is **v0.5.0**. The project hasn't reached v1.0.0 yet so changes that break backwards compatibility may (and have been) introduced in any version update. For information about the changes in each version consult the [CHANGELOG](https://gitlab.com/in-nova/in-nova/-/blob/main/CHANGELOG.md?ref_type=heads).
#### **Important**
-If you are unable, or do not want to contribute under this license please do **not** contribute to this project at all. All content within the project must be under the same license.
+If you are unable, or do not want to contribute under the current license please do **not** contribute to this project at all. All content within the project must be under the same license.
--
GitLab
From c6203b565b72a7ee3f7302c3a0427ff77664447a Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 8 Dec 2024 23:13:46 +0000
Subject: [PATCH 28/60] Fix the line source pressure equation to use Hankel
function of first kind and replace speed with frequency in the scale factor
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
index 83c2ac5..ff4e6ca 100644
--- a/Sound Fields/MATLAB/Functions/lineSrcField.m
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 03/12/2024 (DD/MM/YYYY)
+% Date: 08/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -127,6 +127,6 @@ function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
kr = w * s2rDist/c;
% Calculate complex pressure
- waveField(:, :, fIdx) = 1i * rho * c * Q(:, fIdx).' .* besselh(0, 2, kr)/4;
+ waveField(:, :, fIdx) = 1i * rho * w * Q(:, fIdx).' .* besselh(0, 1, kr)/4;
end
end
\ No newline at end of file
--
GitLab
From cc7e69e0cf967d05192f75ca8a5fcad4120225fc Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 13 Dec 2024 12:31:43 +0000
Subject: [PATCH 29/60] Correct the Green's function for the line source sound
field
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
index ff4e6ca..c63ce09 100644
--- a/Sound Fields/MATLAB/Functions/lineSrcField.m
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 08/12/2024 (DD/MM/YYYY)
+% Date: 13/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -127,6 +127,6 @@ function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
kr = w * s2rDist/c;
% Calculate complex pressure
- waveField(:, :, fIdx) = 1i * rho * w * Q(:, fIdx).' .* besselh(0, 1, kr)/4;
+ waveField(:, :, fIdx) = -rho * w * Q(:, fIdx).' .* besselh(0, 1, kr)/4;
end
end
\ No newline at end of file
--
GitLab
From 73f9f815b47df90eaeffa44ec858a4594efc8dff Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 13 Dec 2024 12:33:15 +0000
Subject: [PATCH 30/60] Just removed a space in CHANGELOG.md
---
CHANGELOG.md | 442 +++++++++++++++++++++++++--------------------------
1 file changed, 221 insertions(+), 221 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 1bc8f7d..2f46b30 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,221 +1,221 @@
-###v0.5.0###
-**Signal Processing - Generic**\
-\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
-
-**Signal Processing - Array Processing**\
-\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
-\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
-\* Modify `firstOrderDma()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
-\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
-\* Correct erroneous output calculations in `firstOrderDma()`.
-\* Moved the `input` input argument to the last of the input arguments list and updated the input argument checks. This change **breaks backward compatibility**.
-
-**Signal Processing - Measurements**\
-\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
-
-**Optimisation - Memetic Algorithm**\
-\* Removed redundant code to handle serial for-loop if parfor is not available.
-
-**Sound Fields**\
-\+ Add function (`lineSrcField()`) to calculate the sound field due to line sources in the frequency domain (tonal).
-
--------------------------------------------
-
-
-###v0.4.1###
-** Optimisation - Memetic Algorithm**\
-\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
-\* Update `selParents()` to use the parallel for-loop (`parfor`) whenever it is available.
-
--------------------------------------------
-
-
-###v0.4.0###
-**Utilities - Generic**\
-\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
-\* The following functions are changed internally to apply the changes in `rotMat3d()`:
- - `rcvGeo()`
- - `srcGeo()`
- - `virtMicGeo()`
-
--------------------------------------------
-
-
-###v0.3.2###
-**Signal Processing - Generic**\
-\+ Add option to use negative delay values with the `winSincFracDel()` function.
-
-**Signal Processing - Array Processing**\
-\* Fix error in `arrMan()` in the inner product of the source directions matrix and the sensor position vectors.\
-\* Fix error in `arrMan()` where the returned array manifold had incorrect dimensions. It was $N_{\textrm{d}} \times N_{\textrm{m}} \times N_{\textrm{k}}$, where $N_{\textrm{d}}$ is the number of source directions, $N_{\textrm{m}}$ the number of sensors and $N_{\textrm{k}}$ the number of frequencies. Now it is $N_{\textrm{m}} \times N_{\textrm{d}} \times N_{\textrm{k}}$.
-
-**Utilities - Geometries**\
-\* Fix bug in `virtMicGeo()` where translation and rotation were not applied to the geometry mesh coordinates array.
-
-**Virtual Sensing - Remote Microphone Technique**\
-\* Update `obsFiltEstTD()`, removing the observatoin filter replication when a single filter is provided for many trials/sound field realisations to reduce the memory footprint of the function.\
-\* Update `obsFiltTD()`, performing the delay operation for the virtual microphone signals only when needed and added check for $SNR = \infty$ to calculate the SNR matrix appropriately.
-\* Fix bug in `obsFiltTD()`, returning the condition number without taking into account the regularisation and SNR values.
-
--------------------------------------------
-
-
-###v0.3.1###
-**Virtual Sensing - Remote Microphone Technique**\
-\+ Add the option to add "noise" to the monitoring microphone signals in the time domain estimation of the optimal observation filters in `obsFiltTD()` to replicate the feature of `obsFilt()`. This change **breaks** backwards compatibility.
-
--------------------------------------------
-
-
-###v0.3.0###
-**General**\
-\* Convert the project to accept column position and direction vectors as input. The affected functions are:
-- `srcGeo()`
-- `rcvGeo()`
-- `virtMicGeo()`
-- `arrMan()`
-- `firstOrderDMA()`
-- `ptsOnSphere()`
-- `twoPtDist()`
-- `circHarm()`
-- `sphHarm()`
-- `planeWave()`
-- `planeWaveSH()`
-- `ptSrcField()`
-- `ptSrcFieldSH()`
-- `ptSrcFieldTD()`
-
--------------------------------------------
-
-
-###v0.2.8###
-**Utilities - Generic**\
-\+ Add the option to pick the order the rotations will be performed in the `rotMat3d()` function. This results is a **backward incompatible** change.
-
-**Utilities - Geometries**\
-\* Update the geometry generating functions `rcvGeo()`, `srcGeo()` and `virtMicGeo` to accommodate the changes in `rotMat3d()`.
-
--------------------------------------------
-
-
-###v0.2.7###
-**Signal Processing - Measurements**\
-\+ Introduction of the *Measurements* "topic" in the Signal Processing part of the codebase.\
-\+ Function to estimate the impulse responses from sweep measurements.
-
--------------------------------------------
-
-
-###v0.2.6###
-**Signal Processing - Generic**\
-\+ Add function to perform fractional octave smoothing of spectra.
-
--------------------------------------------
-
-
-###v0.2.5###
-**Utilities - Geometries**\
-\* Combine translation offsets to a single vector input argument for the `rcvGeo()` MATLAB function.\
-\* Combine translation offsets to a single vector input argument for th `virtMicGeo()` MATLAB function.\
-\* Add the ability to rotate the receiver geometries in `rcvGeo()` MATLAB function.
-
--------------------------------------------
-
-
-###v0.2.4###
-**Virtual Sensing**\
-\* Corrected time-domain observation filter calculations.\
-\* Removed scaling from the cross-correlation calculations in the time-domain optimal observation filters.\
-\* Fixed bug in time-domain observation filter that screwed up dimensions of the filters when only one source was present.
-
-**Sound Fields**\
-\* Fixed bug in time-domain point source signal calculations that screwed up dimensions of the resulting arrays when only one source was present.
-
-**Utilities - Geometries**\
-\* The virtual microphone geometry generation function has been heavily modified and its interface has changed in a **backward incompatible** way. It is now more generic and allows for "arbitrary" size of each Cartesian dimension of the geometries as well as different number of sensors along each dimension.
-
--------------------------------------------
-
-
-###v0.2.3###
-**Virtual Sensing**\
-\+ Added function to estimate the observation filters in the time-domain.\
-\+ Added function to perform estimation with observation filters in the time-domain.
-
-**Signal Processing - Generic**\
-\* Fixed bugs, improved calculations and removed output arguments for the windowed sinc FIR fractional delay filters
-\+ Added four window functions in the calculation of windowed sinc FIR fractional delay filters.\
-
-**Sound Fields**\
-\+ Added function to calculate signals generated by ideal point sources in the time-domain.
-
--------------------------------------------
-
-
-###v0.2.2###
-**Virtual Sensing**\
-\* Update the observation filter and estimation with observation filter functions with better input argument checking and changed the name of variables to match those in the literature.
-
--------------------------------------------
-
-
-###v0.2.1###
-**Virtual Sensing**\
-\* Fix a bug where noise was added to the power spectral density matrix for the optimal observation filter calculations in the noiseless case.
-
--------------------------------------------
-
-
-###v0.2.0###
-**Utilities**\
-\+ Added a bisection method for single-valued functions.\
-\* Improved rotation matrix calculation function with degrees and radian calculations available.\
-\+ Added function to calculate a double-sided spectrum from a single-sided spectrum.\
-\+ Added a function to calculate the intersection of three vectors.\
-\+ Added a pair of functions to check for even or odd elements.\
-\+ Added a function to pick values from a vector based on given probabilities.\
-\+ Added a function to pick randomly unique rows from a matrix.\
-\+ Added a Heaviside step function.\
-\+ Added function to swap argument values.
-
-**Optimisation**\
-\+ Added a MATLAB Memetic Algorithm (MA) implementation.\
-\* The algorithm calls provided functions and can solve generic problems.
-
-**Control**\
-\+ Added tonal control in the frequency domain.\
-\* Control is contained in a single function.\
-\* Can be used with or without virtual sensing.\
-\* Implementations of optimal control and FxLMS (still frequency domain) calculations are available.
-
-**Signal Processing - Array Processing**\
-\+ Added array manifold (steering vector) calculation function.
-
-**Signal Processing - Generic**\
-\+ Added frequency band calculation function with 1/1 octave and 1/3 octave bands.\
-\+ Add fractional delay filter impulse response generation function.
-
-**Virtual Sensing - Generic**\
-\+ Added a multiple coherence calculation function.
-
-**Virtual Sensing - Remote Microphone Technique**\
-\* Divide the observation filter and estimation with the observation filter in two functions. The return arguments are split appropriately to the function they relate.\
-\+ Added option to include noise in the monitoring microphone power spectral density matrix with specified SNR value.
-
-**Sound Fields**\
-\+ Added plane wave calculation function.\
-\+ Added Circular Harmonics calculation function.\
-\+ Added Spherical Harmonics calculation function.\
-\+ Added Discrete Circular Fourier Transform (DCFT) calculation function.\
-\+ Added Inverse Discrete Circular Fourier Transform (IDCFT) calculation function.\
-\+ Added Discrete Spherical Fourier Transform (DSFT) calculation function.\
-\+ Added Inverse Discrete Spherical Fourier Transform (IDSFT) calculation function.\
-\+ Added function to extrapolate a sound field in the Spherical Harmonics domain.\
-\+ Added function to calculate sound field generated by a point source in the Spherical Harmonics domain (truncated order).\
-\+ Added function to calculate sound field generated by a plane wave in the Spherical Harmonics domain (truncated order).
-
--------------------------------------------
-
-
-###v0.1.0###
-This is the initial version of the project.
+###v0.5.0###
+**Signal Processing - Generic**\
+\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
+
+**Signal Processing - Array Processing**\
+\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
+\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
+\* Modify `firstOrderDma()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
+\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
+\* Correct erroneous output calculations in `firstOrderDma()`.
+\* Moved the `input` input argument to the last of the input arguments list and updated the input argument checks. This change **breaks backward compatibility**.
+
+**Signal Processing - Measurements**\
+\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+
+**Optimisation - Memetic Algorithm**\
+\* Removed redundant code to handle serial for-loop if parfor is not available.
+
+**Sound Fields**\
+\+ Add function (`lineSrcField()`) to calculate the sound field due to line sources in the frequency domain (tonal).
+
+-------------------------------------------
+
+
+###v0.4.1###
+** Optimisation - Memetic Algorithm**\
+\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
+\* Update `selParents()` to use the parallel for-loop (`parfor`) whenever it is available.
+
+-------------------------------------------
+
+
+###v0.4.0###
+**Utilities - Generic**\
+\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
+\* The following functions are changed internally to apply the changes in `rotMat3d()`:
+ - `rcvGeo()`
+ - `srcGeo()`
+ - `virtMicGeo()`
+
+-------------------------------------------
+
+
+###v0.3.2###
+**Signal Processing - Generic**\
+\+ Add option to use negative delay values with the `winSincFracDel()` function.
+
+**Signal Processing - Array Processing**\
+\* Fix error in `arrMan()` in the inner product of the source directions matrix and the sensor position vectors.\
+\* Fix error in `arrMan()` where the returned array manifold had incorrect dimensions. It was $N_{\textrm{d}} \times N_{\textrm{m}} \times N_{\textrm{k}}$, where $N_{\textrm{d}}$ is the number of source directions, $N_{\textrm{m}}$ the number of sensors and $N_{\textrm{k}}$ the number of frequencies. Now it is $N_{\textrm{m}} \times N_{\textrm{d}} \times N_{\textrm{k}}$.
+
+**Utilities - Geometries**\
+\* Fix bug in `virtMicGeo()` where translation and rotation were not applied to the geometry mesh coordinates array.
+
+**Virtual Sensing - Remote Microphone Technique**\
+\* Update `obsFiltEstTD()`, removing the observatoin filter replication when a single filter is provided for many trials/sound field realisations to reduce the memory footprint of the function.\
+\* Update `obsFiltTD()`, performing the delay operation for the virtual microphone signals only when needed and added check for $SNR = \infty$ to calculate the SNR matrix appropriately.
+\* Fix bug in `obsFiltTD()`, returning the condition number without taking into account the regularisation and SNR values.
+
+-------------------------------------------
+
+
+###v0.3.1###
+**Virtual Sensing - Remote Microphone Technique**\
+\+ Add the option to add "noise" to the monitoring microphone signals in the time domain estimation of the optimal observation filters in `obsFiltTD()` to replicate the feature of `obsFilt()`. This change **breaks** backwards compatibility.
+
+-------------------------------------------
+
+
+###v0.3.0###
+**General**\
+\* Convert the project to accept column position and direction vectors as input. The affected functions are:
+- `srcGeo()`
+- `rcvGeo()`
+- `virtMicGeo()`
+- `arrMan()`
+- `firstOrderDMA()`
+- `ptsOnSphere()`
+- `twoPtDist()`
+- `circHarm()`
+- `sphHarm()`
+- `planeWave()`
+- `planeWaveSH()`
+- `ptSrcField()`
+- `ptSrcFieldSH()`
+- `ptSrcFieldTD()`
+
+-------------------------------------------
+
+
+###v0.2.8###
+**Utilities - Generic**\
+\+ Add the option to pick the order the rotations will be performed in the `rotMat3d()` function. This results is a **backward incompatible** change.
+
+**Utilities - Geometries**\
+\* Update the geometry generating functions `rcvGeo()`, `srcGeo()` and `virtMicGeo` to accommodate the changes in `rotMat3d()`.
+
+-------------------------------------------
+
+
+###v0.2.7###
+**Signal Processing - Measurements**\
+\+ Introduction of the *Measurements* "topic" in the Signal Processing part of the codebase.\
+\+ Function to estimate the impulse responses from sweep measurements.
+
+-------------------------------------------
+
+
+###v0.2.6###
+**Signal Processing - Generic**\
+\+ Add function to perform fractional octave smoothing of spectra.
+
+-------------------------------------------
+
+
+###v0.2.5###
+**Utilities - Geometries**\
+\* Combine translation offsets to a single vector input argument for the `rcvGeo()` MATLAB function.\
+\* Combine translation offsets to a single vector input argument for th `virtMicGeo()` MATLAB function.\
+\* Add the ability to rotate the receiver geometries in `rcvGeo()` MATLAB function.
+
+-------------------------------------------
+
+
+###v0.2.4###
+**Virtual Sensing**\
+\* Corrected time-domain observation filter calculations.\
+\* Removed scaling from the cross-correlation calculations in the time-domain optimal observation filters.\
+\* Fixed bug in time-domain observation filter that screwed up dimensions of the filters when only one source was present.
+
+**Sound Fields**\
+\* Fixed bug in time-domain point source signal calculations that screwed up dimensions of the resulting arrays when only one source was present.
+
+**Utilities - Geometries**\
+\* The virtual microphone geometry generation function has been heavily modified and its interface has changed in a **backward incompatible** way. It is now more generic and allows for "arbitrary" size of each Cartesian dimension of the geometries as well as different number of sensors along each dimension.
+
+-------------------------------------------
+
+
+###v0.2.3###
+**Virtual Sensing**\
+\+ Added function to estimate the observation filters in the time-domain.\
+\+ Added function to perform estimation with observation filters in the time-domain.
+
+**Signal Processing - Generic**\
+\* Fixed bugs, improved calculations and removed output arguments for the windowed sinc FIR fractional delay filters
+\+ Added four window functions in the calculation of windowed sinc FIR fractional delay filters.\
+
+**Sound Fields**\
+\+ Added function to calculate signals generated by ideal point sources in the time-domain.
+
+-------------------------------------------
+
+
+###v0.2.2###
+**Virtual Sensing**\
+\* Update the observation filter and estimation with observation filter functions with better input argument checking and changed the name of variables to match those in the literature.
+
+-------------------------------------------
+
+
+###v0.2.1###
+**Virtual Sensing**\
+\* Fix a bug where noise was added to the power spectral density matrix for the optimal observation filter calculations in the noiseless case.
+
+-------------------------------------------
+
+
+###v0.2.0###
+**Utilities**\
+\+ Added a bisection method for single-valued functions.\
+\* Improved rotation matrix calculation function with degrees and radian calculations available.\
+\+ Added function to calculate a double-sided spectrum from a single-sided spectrum.\
+\+ Added a function to calculate the intersection of three vectors.\
+\+ Added a pair of functions to check for even or odd elements.\
+\+ Added a function to pick values from a vector based on given probabilities.\
+\+ Added a function to pick randomly unique rows from a matrix.\
+\+ Added a Heaviside step function.\
+\+ Added function to swap argument values.
+
+**Optimisation**\
+\+ Added a MATLAB Memetic Algorithm (MA) implementation.\
+\* The algorithm calls provided functions and can solve generic problems.
+
+**Control**\
+\+ Added tonal control in the frequency domain.\
+\* Control is contained in a single function.\
+\* Can be used with or without virtual sensing.\
+\* Implementations of optimal control and FxLMS (still frequency domain) calculations are available.
+
+**Signal Processing - Array Processing**\
+\+ Added array manifold (steering vector) calculation function.
+
+**Signal Processing - Generic**\
+\+ Added frequency band calculation function with 1/1 octave and 1/3 octave bands.\
+\+ Add fractional delay filter impulse response generation function.
+
+**Virtual Sensing - Generic**\
+\+ Added a multiple coherence calculation function.
+
+**Virtual Sensing - Remote Microphone Technique**\
+\* Divide the observation filter and estimation with the observation filter in two functions. The return arguments are split appropriately to the function they relate.\
+\+ Added option to include noise in the monitoring microphone power spectral density matrix with specified SNR value.
+
+**Sound Fields**\
+\+ Added plane wave calculation function.\
+\+ Added Circular Harmonics calculation function.\
+\+ Added Spherical Harmonics calculation function.\
+\+ Added Discrete Circular Fourier Transform (DCFT) calculation function.\
+\+ Added Inverse Discrete Circular Fourier Transform (IDCFT) calculation function.\
+\+ Added Discrete Spherical Fourier Transform (DSFT) calculation function.\
+\+ Added Inverse Discrete Spherical Fourier Transform (IDSFT) calculation function.\
+\+ Added function to extrapolate a sound field in the Spherical Harmonics domain.\
+\+ Added function to calculate sound field generated by a point source in the Spherical Harmonics domain (truncated order).\
+\+ Added function to calculate sound field generated by a plane wave in the Spherical Harmonics domain (truncated order).
+
+-------------------------------------------
+
+
+###v0.1.0###
+This is the initial version of the project.
--
GitLab
From f38336954038f11d29dddcab13c75c9d99d8a960 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 15 Dec 2024 18:58:46 +0000
Subject: [PATCH 31/60] Fix error in rcvGeo.m applying incorrect rotation to
the circular array elements
---
Utilities/Geometries/MATLAB/Functions/rcvGeo.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
index ff86c5d..cea25d6 100644
--- a/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/rcvGeo.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 17/11/2024 (DD/MM/YYYY)
+% Date: 15/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -343,7 +343,7 @@ function [omniPos, fig8Pos, triPos, boxPos, box2DPos, tetPos, fig8Vec, triVec, b
end
case "uca"
% Get the angle of the positions on the circle
- az = cart2sph(omniPos(1, :), omniPos(2, :), omniPos(3, :));
+ az = rad2deg(cart2sph(omniPos(1, :), omniPos(2, :), omniPos(3, :)));
% Rotate and position tetrahedrals
for measPosIdx = numel(az):-1:1
--
GitLab
From d433d8f5b0a2162a0c57330e25efe4198eab2eea Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 15 Dec 2024 18:59:55 +0000
Subject: [PATCH 32/60] Fix bug in virtMicGeo.m not returning the z coordinates
in the mesh matrix
---
Utilities/Geometries/MATLAB/Functions/virtMicGeo.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m b/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
index 3159cfe..3b18c17 100644
--- a/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 14/11/2024 (DD/MM/YYYY)
+% Date: 15/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -198,7 +198,7 @@ function [vPos, vPosMesh] = virtMicGeo(gType, geoDim, nSens, trans, rot, rotOrd)
% For "Cube" geometry provide the coordinates in a "mesh format"
if nargout > 1 && strcmpi(gType, "Cube")
- vPosMesh = reshape(vPos.', size(x, 1), size(x, 2), 3);
+ vPosMesh = reshape(vPos.', size(x, 1), size(x, 2), [], 3);
elseif nargout > 1
vPosMesh = [];
end
--
GitLab
From 147023d9d0ec24070342a3bd2220201ea2e4513c Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 15 Dec 2024 19:06:07 +0000
Subject: [PATCH 33/60] Add fixes of rcvGeo and virtMicGeo to the CHANGELOG.md
---
CHANGELOG.md | 6 +++++-
1 file changed, 5 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 2f46b30..8c34389 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -17,7 +17,11 @@
\* Removed redundant code to handle serial for-loop if parfor is not available.
**Sound Fields**\
-\+ Add function (`lineSrcField()`) to calculate the sound field due to line sources in the frequency domain (tonal).
+\+ Add `lineSrcField()` to calculate the sound field due to line sources in the frequency domain (tonal).
+
+**Utilities - Geometries**\
+\* Fix error in `rcvGeo()` applying incorrect rotation to the elements in the Uniform Circular Array.\
+\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.
-------------------------------------------
--
GitLab
From 7ebb503f72a5165f0647534fb23c8c196a2a2c3f Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 19 Dec 2024 13:39:14 +0000
Subject: [PATCH 34/60] Parallelise for-loops in obsFiltEst.m
---
.../MATLAB/Functions/obsFiltEst.m | 6 +++---
1 file changed, 3 insertions(+), 3 deletions(-)
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEst.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEst.m
index 99787e1..1b0946d 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEst.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEst.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 01/10/2024 (DD/MM/YYYY)
+% Date: 19/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -109,7 +109,7 @@ function [est, err, errSqr, normErrSqr, See] = obsFiltEst(Pm, O, Pe, Svv)
% Sum of squared estimation errors
if nargout > 2
- for eIdx = size(Pe, 1):-1:1
+ parfor eIdx = 1:size(Pe, 1)
errSqr(eIdx, 1) = err(eIdx, :) * Svv * err(eIdx, :)';
end
end
@@ -117,7 +117,7 @@ function [est, err, errSqr, normErrSqr, See] = obsFiltEst(Pm, O, Pe, Svv)
% Normalised squared errors
if nargout > 3
% Calculate power true spectral density at virtual microphones
- for eIdx = length(errSqr):-1:1
+ parfor eIdx = 1:length(errSqr)
See(eIdx, 1) = Pe(eIdx, :) * Svv * Pe(eIdx, :)';
end
--
GitLab
From 02ca71c6ef2d6cd0f337c9b32c3d4fdc79fb20e6 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 19 Dec 2024 13:40:17 +0000
Subject: [PATCH 35/60] Parallelise for-loop in ptSrcField.m
---
Sound Fields/MATLAB/Functions/ptSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/ptSrcField.m b/Sound Fields/MATLAB/Functions/ptSrcField.m
index 67aca00..b6fcaa1 100644
--- a/Sound Fields/MATLAB/Functions/ptSrcField.m
+++ b/Sound Fields/MATLAB/Functions/ptSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 26/11/2024 (DD/MM/YYYY)
+% Date: 19/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -113,7 +113,7 @@ function waveField = ptSrcField(sPos, rPos, f, Q, c, rho)
% Calculate variables
% ============================================
% Go through the frequencies
- for fIdx = numel(f):-1:1
+ parfor fIdx = 1:numel(f)
% Calculate parameters
w = 2 * pi * f(fIdx); % Radial frequency
k = w/c; % Wavenumber
--
GitLab
From 6165602b528ee97efbb3916506718b54147a44a9 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 19 Dec 2024 13:40:52 +0000
Subject: [PATCH 36/60] Parallelise for-loop in lineSrcField.m
---
Sound Fields/MATLAB/Functions/lineSrcField.m | 4 ++--
1 file changed, 2 insertions(+), 2 deletions(-)
diff --git a/Sound Fields/MATLAB/Functions/lineSrcField.m b/Sound Fields/MATLAB/Functions/lineSrcField.m
index c63ce09..be28339 100644
--- a/Sound Fields/MATLAB/Functions/lineSrcField.m
+++ b/Sound Fields/MATLAB/Functions/lineSrcField.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 13/12/2024 (DD/MM/YYYY)
+% Date: 19/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -121,7 +121,7 @@ function waveField = lineSrcField(sPos, rPos, f, Q, c, rho)
% Calculate variables
% ============================================
% Go through the frequencies
- for fIdx = numel(f):-1:1
+ parfor fIdx = 1:numel(f)
% Calculate parameters
w = 2 * pi * f(fIdx);
kr = w * s2rDist/c;
--
GitLab
From 853a99999993b2747371d13dae57624aa895fcc4 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Thu, 19 Dec 2024 13:58:09 +0000
Subject: [PATCH 37/60] Add parallelised for-loop in ptSrcField(),
lineSrcField() and obsFiltEst() update information in CHANGELOG.md
---
CHANGELOG.md | 5 +++++
1 file changed, 5 insertions(+)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 8c34389..913a4ce 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -18,11 +18,16 @@
**Sound Fields**\
\+ Add `lineSrcField()` to calculate the sound field due to line sources in the frequency domain (tonal).
+\* Parallelize the for-loops in `ptSrcField()` and `lineSrcField()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.
+
**Utilities - Geometries**\
\* Fix error in `rcvGeo()` applying incorrect rotation to the elements in the Uniform Circular Array.\
\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.
+** Virtual Sensing - Remote Microphone Technique**\
+\* Parallelize the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.
+
-------------------------------------------
--
GitLab
From c1f8942e9b4cc0004abfa0bd6a43d234afe9632c Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 21 Dec 2024 21:48:03 +0000
Subject: [PATCH 38/60] Parallelise one for-loop in obsFiltTD.m
---
.../MATLAB/Functions/obsFiltTD.m | 8 +++++---
1 file changed, 5 insertions(+), 3 deletions(-)
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m
index 8559c8c..cfeafab 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 14/10/2024 (DD/MM/YYYY)
+% Date: 21/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -241,9 +241,11 @@ function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, condNum, mMtx, Omean, RmeMean, RmmM
for jIdx = size(m, 3):-1:1
% Go through the monitoring microphones
for mIdx = size(m, 2):-1:1
+ tmp = m(:, mIdx, jIdx);
+
% Calculate the cross-correlations between virtual and monitoring microphones
- for eIdx = size(e, 2):-1:1
- corr = xcorr(m(:, mIdx, jIdx), e(:, eIdx, jIdx), filtLen);
+ parfor eIdx = 1:size(e, 2)
+ corr = xcorr(tmp, e(:, eIdx, jIdx), filtLen);
Rme(:, mIdx, eIdx, jIdx) = corr(filtLen + 1:-1:2);
end
--
GitLab
From d60afc0e8cee503af061bd11eeb22031d9801ec6 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 21 Dec 2024 21:54:31 +0000
Subject: [PATCH 39/60] Change name of obsFiltEstTD.m to obsFiltEstTd.m
---
.../Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m
index 52b957a..ee9417f 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 13/10/2024 (DD/MM/YYYY)
+% Date: 21/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
--
GitLab
From 2e2d88acf8bea0c40aa1a4a779b5a64f993454d2 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 21 Dec 2024 21:55:26 +0000
Subject: [PATCH 40/60] Update CHANGELOG.md with parallelised for-loop in
obsFiltTD.m and change of names for time domain observation filter functions
---
CHANGELOG.md | 6 ++++--
1 file changed, 4 insertions(+), 2 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 913a4ce..92a0b69 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -20,13 +20,15 @@
\+ Add `lineSrcField()` to calculate the sound field due to line sources in the frequency domain (tonal).
\* Parallelize the for-loops in `ptSrcField()` and `lineSrcField()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.
-
**Utilities - Geometries**\
\* Fix error in `rcvGeo()` applying incorrect rotation to the elements in the Uniform Circular Array.\
\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.
** Virtual Sensing - Remote Microphone Technique**\
-\* Parallelize the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.
+\* Parallelise the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
+\* Parallelise on of the for-loops in `obsFiltTD()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
+\* Change function name from `obsFiltTD()` to `obsFiltTd()`.\
+\* Change function name from `obsFiltEstTD()` to `obsFiltEstTd()`.
-------------------------------------------
--
GitLab
From 9fde22b3969bd8b074c12f721af72eb06adb3c78 Mon Sep 17 00:00:00 2001
From: Achilles Kappis <axilleaz@protonmail.com>
Date: Sat, 21 Dec 2024 21:58:20 +0000
Subject: [PATCH 41/60] Rename time domain observation filter functions online
---
.../MATLAB/Functions/{obsFiltEstTD.m => obsFiltEstTd.m} | 0
.../MATLAB/Functions/{obsFiltTD.m => obsFiltTd.m} | 0
2 files changed, 0 insertions(+), 0 deletions(-)
rename Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/{obsFiltEstTD.m => obsFiltEstTd.m} (100%)
rename Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/{obsFiltTD.m => obsFiltTd.m} (100%)
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m
similarity index 100%
rename from Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTD.m
rename to Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m
similarity index 100%
rename from Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m
rename to Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m
--
GitLab
From a15c831f83620bd0b90ce1e429d70c4dc1cf1fc6 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 27 Dec 2024 10:00:26 +0000
Subject: [PATCH 42/60] Update obsFilt.m to return the condition of the
perturbed power spectral density matrix
---
.../MATLAB/Functions/obsFilt.m | 10 ++++------
1 file changed, 4 insertions(+), 6 deletions(-)
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFilt.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFilt.m
index 46a7368..71cb92b 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFilt.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFilt.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 04/10/2024 (DD/MM/YYYY)
+% Date: 25/12/2024 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -154,13 +154,11 @@ function [oOpt, Sme, Smm, condNum] = obsFilt(Pe, Pm, Svv, regFacs, snrVal)
snrMtx = 0;
end
snrMtx = snrMtx .* eye(size(Pm, 1));
-
- % Calcualte observation filters
- invQty = Smm + snrMtx + regMat;
- oOpt = Sme/invQty;
+ Smm = Smm + snrMtx + regMat;
+ oOpt = Sme/Smm;
% Condition number of auto-correlation (power spectrum) matrix
if nargout > 3
- condNum = cond(invQty);
+ condNum = cond(Smm);
end
end
\ No newline at end of file
--
GitLab
From 921032ee750e9a64627f359df6052055ea9d3a8b Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Fri, 27 Dec 2024 10:09:22 +0000
Subject: [PATCH 43/60] Update CHANGELOG.md with changes in obsFilt.m and
rearranged the updates a bit
---
CHANGELOG.md | 9 +++++----
1 file changed, 5 insertions(+), 4 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 92a0b69..d629cdb 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -2,6 +2,9 @@
**Signal Processing - Generic**\
\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
+**Signal Processing - Measurements**\
+\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
+
**Signal Processing - Array Processing**\
\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
@@ -10,9 +13,6 @@
\* Correct erroneous output calculations in `firstOrderDma()`.
\* Moved the `input` input argument to the last of the input arguments list and updated the input argument checks. This change **breaks backward compatibility**.
-**Signal Processing - Measurements**\
-\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
-
**Optimisation - Memetic Algorithm**\
\* Removed redundant code to handle serial for-loop if parfor is not available.
@@ -28,7 +28,8 @@
\* Parallelise the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
\* Parallelise on of the for-loops in `obsFiltTD()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
\* Change function name from `obsFiltTD()` to `obsFiltTd()`.\
-\* Change function name from `obsFiltEstTD()` to `obsFiltEstTd()`.
+\* Change function name from `obsFiltEstTD()` to `obsFiltEstTd()`.\
+\* Update `obsFilt()` to return the condition of the perturbed power spectral density matrix.
-------------------------------------------
--
GitLab
From 2b847cd11e29a8eaf2827d978880f9d28fb919f4 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 13:27:07 +0000
Subject: [PATCH 44/60] Update header documentation of firstOrderDma.m to
resemble the changes in the order of input arguments
---
.../MATLAB/Functions/firstOrderDma.m | 18 +++++++++---------
1 file changed, 9 insertions(+), 9 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
index 6377e7b..b7b0440 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDma.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 04/12/2024 (DD/MM/YYYY)
+% Date: 04/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -12,13 +12,6 @@
% --------------------------------------------------
% Input
%
-% input [numeric]: The input to the array. 3D array/matrix with dimensions
-% MxIxF, where M represents the number of microphones and
-% must be even, I is the number of measurments (or number
-% of sources), and F is the number of frequencies of
-% interest. Each microphone pair is treated as a first
-% order DMA.
-%
% freq [numeric]: The frequencies of interest. A vector with number of
% elements matching the third dimension of the "input"
% parameter.
@@ -50,6 +43,13 @@
% appropriate value is automatically set, and
% this argument is ignored. [Default: 0]
%
+% input [numeric]: The input to the array. 3D array/matrix with dimensions
+% MxIxF, where M represents the number of microphones and
+% must be even, I is the number of measurments (or number
+% of sources), and F is the number of frequencies of
+% interest. Each microphone pair is treated as a first
+% order DMA.
+%
% --------------------------------------------------
% Output
%
@@ -97,7 +97,7 @@ function [h, output] = firstOrderDma(freq, d, pPattern, beta, input)
end
% Check alpha (angle of null)
- if nargin > 2 && ~isempty(pPattern);
+ if nargin > 2 && ~isempty(pPattern)
if isstring(pPattern) || ischar(pPattern)
validateattributes(pPattern, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, 'Polar pattern', 3);
elseif isnumeric(pPattern)
--
GitLab
From a8f90b96ccc1016a20353acde77652be2670c74f Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:07:21 +0000
Subject: [PATCH 45/60] Update the reference in arrManTD.m
---
.../Array Processing/MATLAB/Functions/arrManTD.m | 5 ++---
1 file changed, 2 insertions(+), 3 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m b/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
index 81d5499..1639b2a 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
@@ -84,9 +84,8 @@
% calculated internally with dimensions equal to the mean of each
% dimension over all the sensor positions.
%
-% - The implementation is based on "First-Order Differential Microphone
-% Arrays from a Time-Domain Broadband Perspective" by Buchris, Cohen and
-% Benesty.
+% - The implementation is based on "On the Design of Time-Domain
+% Differential Microphone Arrays" by Buchris, Cohen and Benesty.
%
% - Dependencies: * winSincFracDel() to calculate the fractional delay
% filters.
--
GitLab
From 5f36b0abfc5f079ce1b7973228277e2543561113 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:08:05 +0000
Subject: [PATCH 46/60] Add function to calculate filter and array output of a
first order DMA in the time domain
---
.../MATLAB/Functions/firstOrderDmaTd.m | 149 ++++++++++++++++++
1 file changed, 149 insertions(+)
create mode 100644 Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
new file mode 100644
index 0000000..97530c4
--- /dev/null
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
@@ -0,0 +1,149 @@
+%% 1st order Time-Domain Directional Microphone Array
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 04/01/2025 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the filters and output of a first order
+% Differential Microphone Array (DMA) in the time-domain.
+% --------------------------------------------------
+% Input
+%
+% sMtx [numeric]: The steering matrix of the differential array. This must
+% be an NxMx2 array, where N is 2 * L with, L the length of
+% the desired FIR beamforming filter. M is the length of
+% time-domain the steering vectors (see notes). The third
+% dimension corresponds to the number of microphones in the
+% array; at the moment only first order arrays with two
+% microphones are supported.
+%
+% del [numeric] (Optional): This is a positive real integer corresponding
+% to the delay index at which the right-hand side
+% vector will have its non-zero sample. For an
+% explanation look at the literature referenced
+% in the notes, where ways to calculate the
+% optimal value are briefly discussed.
+% [Default: 14].
+%
+% regFac [numeric] (Optional): This is a positive scalar that acts to
+% regularise the inversion of the constraints
+% matrix. [Default: 0].
+%
+% sig [numeric] (Optional): This must be an KxPx2 array, corresponding to
+% the signals to be filtered by the array, with
+% real valued elements. The dimensions are as
+% follows: K is the length of the signals in
+% samples, P are the number of signals to be
+% filtered and the last dimension corresponds to
+% the fact that the signals must be given
+% separately for each microphone. If no signals
+% are given, the output arguments "out" and
+% "outSum" will be empty. [Default: []].
+%
+% --------------------------------------------------
+% Output
+%
+% hMtx [numeric]: The beamforming FIR filter in vector format. This will be
+% an Nx1 vector, where N is again 2 * L, with L the length
+% of the beamforming FIR filter. The filters corresponding
+% to the microphones are stacked and not interleaved.
+%
+% h [numeric]: The beamforming FIR filters given for each microphone
+% separately. This is an Lx2 matrix, with L the length of the
+% filters.
+%
+% out [numeric]: The output of each microphone (after filtering) for each
+% signal provided in the input argument "sigs". This is a
+% KxPx2 array, with K the length of the signals and P the
+% number of signals.
+%
+% outSum [numeric]: This is the sum of the output signals (after filtering)
+% of each microphone. This is a KxP matrix and
+% corresponds to "sum(out, 3)".
+%
+% --------------------------------------------------
+% Notes
+%
+% - The implementation is based on: "On the Design of Time-Domain
+% Differential Microphone Arrays" by Buchris, Cohen and Benesty.
+%
+% --------------------------------------------------
+function [hMtx, h, out, outSum] = firstOrderDmaTD(sMtx, del, regFac, sig)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(1, 4);
+ nargoutchk(0, 4);
+
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ validateattributes(sMtx, "numeric", {'3d', 'real', 'nonempty', 'nonnan', 'finite', 'size', [NaN, NaN, 2]}, mfilename, "Angle of null response", 1);
+ validateattributes(del, "numeric", {'scalar', 'nonempty', 'integer', 'real', 'finite', 'nonnan', 'nonnegative'}, mfilename, "The maximum delay between sensors and reference point of the array", 2);
+
+ if nargin > 2 && ~isempty(regFac)
+ validateattributes(regFac, "numeric", {'scalar', 'nonempty', 'real', 'finite', 'nonnan', 'nonnegative'}, mfilename, "The regularisation factor", 3);
+ else
+ regFac = 0;
+ end
+
+ if nargin > 3 && ~isempty(sig)
+ validateattributes(sig, "numeric", {'3d', 'nonempty', 'real', 'nonnan', 'finite', 'size', [NaN, NaN, 2]}, mfilename, "Signals to be filtered", 4);
+ else
+ sig = [];
+ end
+
+
+ % ====================================================
+ % Calculate parameters
+ % ====================================================
+ % Calculate right hand side vector
+ rhsVec = zeros(2 * size(sMtx, 2), 1);
+ rhsVec(del) = 1;
+
+
+ % ====================================================
+ % Pre-process arguments
+ % ====================================================
+ % Stack the array manifolds for the two directions into a single matrix
+ sMtx = reshape(sMtx, size(sMtx, 1), []).';
+
+
+ % ====================================================
+ % Calculate filter(s)
+ % ====================================================
+ if regFac ~= 0
+ hMtx = regFac * eye(size(sMtx, 2));
+ else
+ hMtx = 0;
+ end
+
+ hMtx = (sMtx.' * sMtx + hMtx)\(sMtx.') * rhsVec;
+
+ % ====================================================
+ % Return additional output arguments
+ % ====================================================
+ % Return the individual filters
+ if nargout > 1
+ h = reshape(hMtx, [], 2);
+ end
+
+ % Return the filtered signal(s) [output of the DMA]
+ if nargout > 2 && ~isempty(sig)
+ for idx = size(h, 2):-1:1
+ out(:, :, idx) = filter(h(:, idx), 1, sig(:, :, idx));
+ end
+ else
+ out = [];
+ end
+
+ if nargout > 3 && ~isempty(out)
+ outSum = squeeze(sum(out, 3));
+ else
+ outSum = [];
+ end
+end
\ No newline at end of file
--
GitLab
From ce4f60a7dab7252dd69c8e6a7d3afb9930b695c0 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:09:10 +0000
Subject: [PATCH 47/60] Remove arrManTD.m to rename
---
.../MATLAB/Functions/arrManTD.m | 188 ------------------
1 file changed, 188 deletions(-)
delete mode 100644 Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m b/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
deleted file mode 100644
index 1639b2a..0000000
--- a/Signal Processing/Array Processing/MATLAB/Functions/arrManTD.m
+++ /dev/null
@@ -1,188 +0,0 @@
-%% Time-domain array manifold (steering vector)
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 02/11/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Calculate the time-domain equivalent of the array manifold
-% (steering vector) of an arbitrary array.
-% --------------------------------------------------
-% Input
-%
-% mPos [numeric]: The positions of the sensors in Cartesian coordinates.
-% This must be an 3xM matrix where M denotes the number of
-% sensors and for each sensor the coordinates are in
-% [x; y; z] format.
-%
-% dirs [numeric]: These are the directions for which the array manifold
-% will be calculated. This must be an 2xN matrix where N is
-% the number of directions and for each direction the
-% azimuth and elevation/inclination is given. The values
-% must be provided in radians. The directions must be
-% pointing towards the incoming plane waves and not towards
-% the array/origin.
-%
-% filtLen [numeric]: This argument represents the length of the DMA filter
-% if this function is used in the context of time-domain
-% DMA filter design. In general, defines the number of
-% rows the output array will have per direction.
-%
-% delLen [numeric] (Optional): The length of the (fractional) delay
-% filters used to calculate the "array
-% manifold" in samples. If this value is not
-% large enough to accommodate the largest
-% delay between the point of reference of
-% the array (see notes) and the most distant
-% sensor, the length of the filters is set
-% automatically to that number (internal to
-% the "winSincFracDel()" function - see
-% dependencies). [Default: 16].
-%
-% delWin [char/string/numeric] (Optional): The window that will be applied
-% to the (fractional) delay
-% filter. This can be
-% "Rectangular", "Hann", "Hamming"
-% or a positive scalar
-% corresponding to the Kaiser
-% alpha parameter. The string/char
-% values are not case-sensitive.
-% [Default: 3].
-%
-% fs [numeric] (Optional): The sampling rate in Hz. This value is used to
-% calculate the final length of the filters in the
-% returned "array manifold" matrix (the number of
-% columns) and "delLen" in case it is not
-% provided. For the number of columns in the
-% returned "array manifold" see the output
-% arguments. [Default: 8e3].
-%
-% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
-%
-% --------------------------------------------------
-% Output
-%
-% am [numeric]: This is an LxKxMxN array, where L is equal to "filtLen", K
-% is equal to filtLen + delLen + ceil(2 * max(abs(del)))
-% where "del" is the delay between the sensors and the centre
-% of mass of the array. M is the number of sensors and N is
-% the number of directions.
-%
-% maxDel [numeric]: The maximum delay from the centre of mass to the
-% sensors.
-%
-% amMtx [numeric]: This is an LMxKxN array, where LM is equal to L * M.
-% This corresponds to the "vectorised" version of "am"
-% with the "am" of each sensor stacked under each other.
-%
-% --------------------------------------------------
-% Notes
-%
-% - The point of reference for the array is its centre of mass and it is
-% calculated internally with dimensions equal to the mean of each
-% dimension over all the sensor positions.
-%
-% - The implementation is based on "On the Design of Time-Domain
-% Differential Microphone Arrays" by Buchris, Cohen and Benesty.
-%
-% - Dependencies: * winSincFracDel() to calculate the fractional delay
-% filters.
-% --------------------------------------------------
-function [am, maxDel, amMtx] = arrManTD(mPos, dirs, filtLen, delLen, delWin, fs, c)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(3, 7);
- nargoutchk(0, 3);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(mPos, "numeric", {'2d', 'finite', 'nonnan', 'nonempty', 'real', 'nrows', 3}, mfilename, "Position of sensors", 1);
- validateattributes(dirs, "numeric", {'2d', 'finite', 'nonnan', 'nonempty', 'real', 'nrows', 2}, mfilename, "Directions of interest", 2);
- validateattributes(filtLen, "numeric", {'scalar', 'finite', 'nonnan', 'nonempty', 'positive', 'real'}, mfilename, "Number of rows of the 'array manifold'", 3);
-
- % Optional arguments
- if nargin > 3 && ~isempty(delLen)
- validateattributes(delLen, "numeric", {'scalar', 'finite', 'nonnan', 'nonempty', 'positive', 'real'}, mfilename, "The length of the delay filters", 4);
- else
- delLen = 16;
- end
-
- if nargin < 5 || (nargin > 4 && isempty(delWin))
- delWin = 3;
- end
-
- if nargin > 5 && ~isempty(fs)
- validateattributes(fs, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'positive', 'real'}, mfilename, "Sampling rate", 6);
- else
- fs = 8e3;
- end
-
- if nargin > 6 && ~isempty(c)
- validateattributes(c, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'positive', 'real'}, mfilename, "Speed of sound", 7);
- else
- c = 343;
- end
-
-
- % =============================================
- % Calculate delays and filters
- % =============================================
- % Set the centre of mass to the origin
- mPos = mPos - mean(mPos, 2);
-
- % Get Cartesian coordinates of unit vectors pointing in angle directions
- [x, y, z] = sph2cart(dirs(1, :), dirs(2, :), ones(size(dirs(1, :))));
-
- % Calculate each sensor's delay with respect to the centre of mass in samples
- del = -[x; y; z].' * fs * mPos/c;
-
- % Calculate the total length of the filters
- maxDel = max(abs(del), [], "all");
- delLen = ceil(delLen + 2 * maxDel);
-
- % Get the delay filters
- for angIdx = size(del, 1):-1:1
- for mIdx = size(del, 2):-1:1
- h(:, mIdx, angIdx) = winSincFracDel(del(angIdx, mIdx), delLen, delWin);
- end
- end
-
-
- % =============================================
- % Form Sylvester matrices
- % =============================================
- % Turn warnings about the generation of Toeplitz matrices off
- w = warning("off", "MATLAB:toeplitz:DiagonalConflict");
-
- try
- % Generate Sylvester matrices for each mic and angle
- for angIdx = size(h, 3):-1:1
- for mIdx = size(h, 2):-1:1
- am(:, :, mIdx, angIdx) = toeplitz([h(:, mIdx, angIdx); zeros(filtLen - 1, 1)], zeros(filtLen, 1)).';
- end
- end
- catch
- % Make sure to turn the warnings on again if an error occurs
- warning(w.state, w.identifier);
- end
-
- % Reset the warnings state
- warning(w.state, w.identifier);
-
-
- % =============================================
- % Return additional output arguments
- % =============================================
- if nargout > 1
- maxDel = ceil(maxDel);
- end
-
- if nargout > 2
- amMtx = reshape(permute(am, [1, 3, 2, 4]), [], size(am, 2), size(am, 4));
- end
-end
\ No newline at end of file
--
GitLab
From 358b3636181c92fa091f8698571b4fce40a96733 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:09:49 +0000
Subject: [PATCH 48/60] Bring back arrManTd.m after renaming
---
.../MATLAB/Functions/arrManTd.m | 188 ++++++++++++++++++
1 file changed, 188 insertions(+)
create mode 100644 Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m b/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
new file mode 100644
index 0000000..1639b2a
--- /dev/null
+++ b/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
@@ -0,0 +1,188 @@
+%% Time-domain array manifold (steering vector)
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 02/11/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the time-domain equivalent of the array manifold
+% (steering vector) of an arbitrary array.
+% --------------------------------------------------
+% Input
+%
+% mPos [numeric]: The positions of the sensors in Cartesian coordinates.
+% This must be an 3xM matrix where M denotes the number of
+% sensors and for each sensor the coordinates are in
+% [x; y; z] format.
+%
+% dirs [numeric]: These are the directions for which the array manifold
+% will be calculated. This must be an 2xN matrix where N is
+% the number of directions and for each direction the
+% azimuth and elevation/inclination is given. The values
+% must be provided in radians. The directions must be
+% pointing towards the incoming plane waves and not towards
+% the array/origin.
+%
+% filtLen [numeric]: This argument represents the length of the DMA filter
+% if this function is used in the context of time-domain
+% DMA filter design. In general, defines the number of
+% rows the output array will have per direction.
+%
+% delLen [numeric] (Optional): The length of the (fractional) delay
+% filters used to calculate the "array
+% manifold" in samples. If this value is not
+% large enough to accommodate the largest
+% delay between the point of reference of
+% the array (see notes) and the most distant
+% sensor, the length of the filters is set
+% automatically to that number (internal to
+% the "winSincFracDel()" function - see
+% dependencies). [Default: 16].
+%
+% delWin [char/string/numeric] (Optional): The window that will be applied
+% to the (fractional) delay
+% filter. This can be
+% "Rectangular", "Hann", "Hamming"
+% or a positive scalar
+% corresponding to the Kaiser
+% alpha parameter. The string/char
+% values are not case-sensitive.
+% [Default: 3].
+%
+% fs [numeric] (Optional): The sampling rate in Hz. This value is used to
+% calculate the final length of the filters in the
+% returned "array manifold" matrix (the number of
+% columns) and "delLen" in case it is not
+% provided. For the number of columns in the
+% returned "array manifold" see the output
+% arguments. [Default: 8e3].
+%
+% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
+%
+% --------------------------------------------------
+% Output
+%
+% am [numeric]: This is an LxKxMxN array, where L is equal to "filtLen", K
+% is equal to filtLen + delLen + ceil(2 * max(abs(del)))
+% where "del" is the delay between the sensors and the centre
+% of mass of the array. M is the number of sensors and N is
+% the number of directions.
+%
+% maxDel [numeric]: The maximum delay from the centre of mass to the
+% sensors.
+%
+% amMtx [numeric]: This is an LMxKxN array, where LM is equal to L * M.
+% This corresponds to the "vectorised" version of "am"
+% with the "am" of each sensor stacked under each other.
+%
+% --------------------------------------------------
+% Notes
+%
+% - The point of reference for the array is its centre of mass and it is
+% calculated internally with dimensions equal to the mean of each
+% dimension over all the sensor positions.
+%
+% - The implementation is based on "On the Design of Time-Domain
+% Differential Microphone Arrays" by Buchris, Cohen and Benesty.
+%
+% - Dependencies: * winSincFracDel() to calculate the fractional delay
+% filters.
+% --------------------------------------------------
+function [am, maxDel, amMtx] = arrManTD(mPos, dirs, filtLen, delLen, delWin, fs, c)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 7);
+ nargoutchk(0, 3);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(mPos, "numeric", {'2d', 'finite', 'nonnan', 'nonempty', 'real', 'nrows', 3}, mfilename, "Position of sensors", 1);
+ validateattributes(dirs, "numeric", {'2d', 'finite', 'nonnan', 'nonempty', 'real', 'nrows', 2}, mfilename, "Directions of interest", 2);
+ validateattributes(filtLen, "numeric", {'scalar', 'finite', 'nonnan', 'nonempty', 'positive', 'real'}, mfilename, "Number of rows of the 'array manifold'", 3);
+
+ % Optional arguments
+ if nargin > 3 && ~isempty(delLen)
+ validateattributes(delLen, "numeric", {'scalar', 'finite', 'nonnan', 'nonempty', 'positive', 'real'}, mfilename, "The length of the delay filters", 4);
+ else
+ delLen = 16;
+ end
+
+ if nargin < 5 || (nargin > 4 && isempty(delWin))
+ delWin = 3;
+ end
+
+ if nargin > 5 && ~isempty(fs)
+ validateattributes(fs, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'positive', 'real'}, mfilename, "Sampling rate", 6);
+ else
+ fs = 8e3;
+ end
+
+ if nargin > 6 && ~isempty(c)
+ validateattributes(c, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'positive', 'real'}, mfilename, "Speed of sound", 7);
+ else
+ c = 343;
+ end
+
+
+ % =============================================
+ % Calculate delays and filters
+ % =============================================
+ % Set the centre of mass to the origin
+ mPos = mPos - mean(mPos, 2);
+
+ % Get Cartesian coordinates of unit vectors pointing in angle directions
+ [x, y, z] = sph2cart(dirs(1, :), dirs(2, :), ones(size(dirs(1, :))));
+
+ % Calculate each sensor's delay with respect to the centre of mass in samples
+ del = -[x; y; z].' * fs * mPos/c;
+
+ % Calculate the total length of the filters
+ maxDel = max(abs(del), [], "all");
+ delLen = ceil(delLen + 2 * maxDel);
+
+ % Get the delay filters
+ for angIdx = size(del, 1):-1:1
+ for mIdx = size(del, 2):-1:1
+ h(:, mIdx, angIdx) = winSincFracDel(del(angIdx, mIdx), delLen, delWin);
+ end
+ end
+
+
+ % =============================================
+ % Form Sylvester matrices
+ % =============================================
+ % Turn warnings about the generation of Toeplitz matrices off
+ w = warning("off", "MATLAB:toeplitz:DiagonalConflict");
+
+ try
+ % Generate Sylvester matrices for each mic and angle
+ for angIdx = size(h, 3):-1:1
+ for mIdx = size(h, 2):-1:1
+ am(:, :, mIdx, angIdx) = toeplitz([h(:, mIdx, angIdx); zeros(filtLen - 1, 1)], zeros(filtLen, 1)).';
+ end
+ end
+ catch
+ % Make sure to turn the warnings on again if an error occurs
+ warning(w.state, w.identifier);
+ end
+
+ % Reset the warnings state
+ warning(w.state, w.identifier);
+
+
+ % =============================================
+ % Return additional output arguments
+ % =============================================
+ if nargout > 1
+ maxDel = ceil(maxDel);
+ end
+
+ if nargout > 2
+ amMtx = reshape(permute(am, [1, 3, 2, 4]), [], size(am, 2), size(am, 4));
+ end
+end
\ No newline at end of file
--
GitLab
From 6b6fa73fbfffb793560751a95d448cf29d18776f Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:13:28 +0000
Subject: [PATCH 49/60] Update CHANGELOG.md with the addition of
firstOrderDmaTd.m
---
CHANGELOG.md | 3 ++-
1 file changed, 2 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index d629cdb..63ae65c 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -6,7 +6,8 @@
\* Fix error in `invSweep()` throwing an error when the `causIrLen` input argument was not provided.
**Signal Processing - Array Processing**\
-\+ Add function to calculate the time-domain equivalent of the *array manifold* of an arbitrary array.\
+\+ Add function `arrManTd()` to calculate the time domain equivalent of the *array manifold* of an arbitrary array.\
+\+ Add function `firstOrderDmaTd()` to calculate the FIR filter and filtered output(s) for a first order Differential Microphone Array (DMA) in the time domain.\
\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
\* Modify `firstOrderDma()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
--
GitLab
From 6a196311947ad8ff968d74169e25bc434f83bf31 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:40:08 +0000
Subject: [PATCH 50/60] Remove functions from sound fields to rename
---
Sound Fields/MATLAB/Functions/extrpFieldSH.m | 156 ---------------
Sound Fields/MATLAB/Functions/planeWaveSH.m | 160 ----------------
Sound Fields/MATLAB/Functions/ptSrcFieldSH.m | 169 -----------------
Sound Fields/MATLAB/Functions/ptSrcFieldTD.m | 189 -------------------
4 files changed, 674 deletions(-)
delete mode 100644 Sound Fields/MATLAB/Functions/extrpFieldSH.m
delete mode 100644 Sound Fields/MATLAB/Functions/planeWaveSH.m
delete mode 100644 Sound Fields/MATLAB/Functions/ptSrcFieldSH.m
delete mode 100644 Sound Fields/MATLAB/Functions/ptSrcFieldTD.m
diff --git a/Sound Fields/MATLAB/Functions/extrpFieldSH.m b/Sound Fields/MATLAB/Functions/extrpFieldSH.m
deleted file mode 100644
index 20ca9b5..0000000
--- a/Sound Fields/MATLAB/Functions/extrpFieldSH.m
+++ /dev/null
@@ -1,156 +0,0 @@
-%% Extrapolate the sound field in the Spherical Harmonics domain
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 02/06/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Extrapolate sound field pressure field coefficients in the
-% Spherical Harmonics domain.
-% --------------------------------------------------
-% Input
-%
-% shCoeffs [numeric]: The spherical harmonics pressure coefficients at the
-% measurement position. This must be an NxMxK matrix
-% with N being the number of coefficients per position,
-% M the number of positions/measurements and K the
-% number of frequencies of interest.
-%
-% r [numeric]: This is a vector with at least two elements that correspond
-% to the radial distance component of the spherical
-% coordinates of each position the sound field needs to be
-% extrapolated to. The first element corresponds to the radial
-% distance (radius) of the measurement array.
-%
-% k [numeric]: The wavenumbers of the frequencies of interest. This must be
-% a vector with number of elements equal to K, the frequencies
-% of interest.
-%
-% w [numeric] (Optional): Pressure coefficient weights to be applied before
-% the extrapolation process. This can be used to
-% effectively suppress or use only specific modes
-% for the extrapolation process.
-% [Default: ones(size(shCoeffs, 1), 1)].
-%
-% dir [numeric] (Optional): This signifies the "direction" of the sound
-% field. I.e. whether it is an outgoing or
-% ingoing sound field. In practice, an outgoing
-% sound field is one where the radial distance
-% components of the sources are smaller than the
-% corresponding components of the measurement and
-% extrapolated positions. In the opposite case,
-% where the radial distance components are
-% greater than those of the measurement and
-% extrapolated positions, the field is an ingoing
-% field. The possible values of the argument are
-% "In", "Ingoing", "Converging", "Out",
-% "Outgoing", "Diverging" and are NOT
-% case-sensitive. [Default: "In"].
-%
-% --------------------------------------------------
-% Output
-%
-% extrCoeffs [numeric]: This is an NxMxJxK matrix with the extrapolated
-% coefficients, where J is the number of
-% extrapolation positions and equals length(r) - 1.
-% The coefficients may be Inf if the Bessel (or
-% Hankel) functions at the measured radius exhibit a
-% zero (or infinity for Hankel functions).
-%
-% extrCoeffsReduced [numeric]: This argument holds the same values as
-% "extrCoeffs" but with the singular
-% dimensions being reduced. For example, if
-% extrapolation at a single position is needed
-% the array has dimensions NxMxK. If
-% extrapolating at a single frequency, the
-% dimensions are NxMxJ and if extrapolating at
-% a single position and frequency the result
-% is of dimensions NxM. This is equivalent to
-% the command squeeze(extrCoeffs).
-%
-% --------------------------------------------------
-% Notes
-%
-% --------------------------------------------------
-function [extrCoeffs, extrCoeffsReduced] = extrpFieldSH(shCoeffs, r, k, w, dir)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(3, 5);
- nargoutchk(0, 1);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(shCoeffs, "numeric", {'2d', 'nonempty', 'nonnan', 'finite'}, mfilename, "The spherical harmonic pressure coefficients at the measurement positions", 1);
- validateattributes(r, "numeric", {'vector', 'nonempty', 'nonnan', 'finite', 'nonnegative'}, mfilename, "The radial distance components of all positions", 2);
-
- % Check if we have AT LEAST two distances
- if numel(r) < 2
- error("You must provide at least two distances, with the first being the radius of the measurement array");
- end
-
- validateattributes(k, "numeric", {'vector', 'nonempty', 'nonnan', 'finite', 'nonnegative', 'real'}, mfilename, "The wavenumbers of the frequencies of interest", 3);
-
- % Validate optional arguments
- if nargin > 3 && ~isempty(w)
- validateattributes(w, "numeric", {'vector', 'nonempty', 'nonnan', 'numel', size(shCoeffs, 1), 'real', 'finite'}, mfilename, "Coefficient weights", 4);
- else
- w = ones(size(shCoeffs, 1), 1);
- end
-
- if nargin > 4 && ~isempty(dir)
- validateattributes(dir, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "The 'direction' of the sound field", 5);
- validatestring(dir, ["In", "Ingoing", "Converging", "Out", "Outgoing", "Diverging"], mfilename, "The 'direction' of the sound field", 5);
-
- if sum(strcmpi(dir, ["In", "Ingoing", "Converging"])) > 0
- dir = "In";
- end
- else
- dir = "In";
- end
-
- % =============================================
- % Calculate parameters
- % =============================================
- N = sqrt(size(shCoeffs, 1)) - 1; % Highest spherical harmonics order
- nDeg = 1:2:2 * N + 1; % Number of degrees/modes per order
-
-
- % =============================================
- % Extrapolate coefficients
- % =============================================
- % Create some matrices to ease calculations
- if strcmpi(dir, "in")
- Jref = besselj(repmat(0:N, length(k), 1), repmat(r(1) * k(:), 1, N + 1));
- else
- Jref = besselh(repmat(0:N, length(k), 1), repmat(r(1) * k(:), 1, N + 1));
- end
-
- % Go through the frequencies
- for fIdx = size(shCoeffs, 3):-1:1
- % Go through the measurements/positions
- for rIdx = length(r) - 1:-1:1
- tmpIdx = rIdx + 1;
- if strcmpi(dir, "in")
- J = besselj(0:N, repmat(k(fIdx) * r(tmpIdx), 1, N + 1));
- else
- J = besselh(0:N, repmat(k(fIdx) * r(tmpIdx), 1, N + 1));
- end
-
- % Replicate and weight J coefficients
- J = w(:).' .* repelem(J./Jref(fIdx, :), nDeg);
- J = diag(J);
-
- % Calculate extrapolated pressure coefficients
- extrCoeffs(:, :, rIdx, fIdx) = J * shCoeffs(:, :, fIdx);
- end
- end
-
- if nargout > 1
- extrCoeffsReduced = squeeze(extrCoeffs);
- end
-end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/planeWaveSH.m b/Sound Fields/MATLAB/Functions/planeWaveSH.m
deleted file mode 100644
index 3dd6d06..0000000
--- a/Sound Fields/MATLAB/Functions/planeWaveSH.m
+++ /dev/null
@@ -1,160 +0,0 @@
-%% Calculate wavefield generated by plane waves in the spherical harmonics domain
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 29/09/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Calculate the wavefield generated by plane waves at
-% specified directions in the spherical harmonics domain.
-% --------------------------------------------------
-% Input
-%
-% sDir [numeric]: The directions of the sources in spherical coordinates.
-% The variable must be a matrix with dimensions 2xN where N
-% is the number of plane waves and the rows represent the
-% azimuth and inclination respectively. The values must be
-% in radians.
-%
-% rPos [numeric]: The position of the receivers in Cartesian coordinates.
-% The variable must be a matrix with dimensions 3xM where
-% M is the number of receivers and the rows represent the
-% x, y and z coordinates respectively.
-%
-% f [numeric]: The frequencies for which to calculate the wave field. It
-% must be a scalar or vector.
-%
-% N [numeric]: The highest order of the sound field. This must be a
-% real non-negative scalar.
-%
-% Q [numeric] (Optional): The source strength of the plane waves. This
-% variable must be either a scalar corresponding to
-% all sources having the same strength, or a matrix
-% with dimensions NxK, where K is the number of
-% frequencies of interest. The strengths can be
-% complex numbers. [Default: 1].
-%
-% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
-%
-% --------------------------------------------------
-% Output
-%
-% pCoeffs [numeric]: The complex pressure coefficients on a sphere at
-% distance from the origin equal to that of the receiver
-% position distances, as generated by plane waves from
-% directions given in "sDir". The dimensions of the
-% array are ((J + 1)^2)xNxMxF, where J is the highest
-% spherical harmonic order, N the number of plane waves,
-% M the number of receivers and F the number of
-% frequencies.
-%
-% waveFieldDecomp [numeric]: The complex pressures at the position of the
-% receivers as calculate by the Inverse Discrete
-% Spherical Fourier Transform for each plane
-% wave separately. The dimensions of the array
-% are MxNxF, where again M is the number of
-% receivers, N the number of plane waves and F
-% the number of frequencies.
-%
-% waveField [numeric]: The sum of the complex pressures at each receiver
-% position (the "total" sound field at each receiver
-% position). This is an matrix with dimensions MxF,
-% where again M is the number of receivers and F the
-% number of frequencies.
-%
-% --------------------------------------------------
-% Notes
-%
-% --------------------------------------------------
-function [pCoeffs, waveFieldDecomp, waveField] = planeWaveSH(sDir, rPos, f, N, Q, c)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(3, 6);
- nargoutchk(0, 3);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(sDir, "numeric", {'2d', 'rows', 2, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Source positions", 1);
- validateattributes(rPos, "numeric", {'2d', 'rows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Receiver positions", 2);
- validateattributes(f, "numeric", {'vector', 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Frequencies", 3);
- validateattributes(N, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'real', 'nonnegative'}, mfilename, "Highest order of the generated soundfield", 4);
-
- % Validate optional arguments
- if nargin > 4 && ~isempty(Q)
- validateattributes(Q, "numeric", {'2d', 'finite', 'nonempty', 'nonnan'}, mfilename, "Source strength(s)", 5);
-
- % Make sure Q has the correct dimensions
- if ~isscalar(Q) && sum((size(Q) ~= [size(sDir, 2), length(f)])) > 0
- error("Source strength(s) must be either a scalar or its length must match the number of sources.");
- end
- else
- Q = 1;
- end
-
- if nargin > 5 && ~isempty(c)
- validateattributes(c, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Speed of propagation", 6);
- else
- c = 343;
- end
-
-
- % ====================================================
- % Pre-process data
- % ====================================================
- % Make sure the source strength has the correct dimensions
- if isscalar(Q)
- Q = Q * ones(size(sDir, 2), length(f));
- end
-
- % Convert coordinates to spherical coordinate system
- sDir(2, :) = (pi/2) - sDir(2, :);
-
- [rAz, rEl, rR] = cart2sph(rPos(1, :), rPos(2, :), rPos(3, :));
- rEl = (pi/2) - rEl;
-
-
- % =============================================
- % Calculate needed parameters
- % ============================================
- % Wavenumbers
- k = 2 * pi * f/c;
-
- % Spherical harmonics corresponding to the directions of the sources
- sh = sphHarm(sDir, N);
-
-
- % =============================================
- % Calculate variables
- % ============================================
- % Go through the frequencies
- for fIdx = numel(f):-1:1
- % Go through receiver positions
- for rIdx = length(rR):-1:1
- J = 1i.^(0:N) .* besselj(0:N, repelem(k(fIdx) * rR(rIdx), N + 1));
- J = diag(repelem(J, 1:2:2 * N + 1));
-
- pCoeffs(:, :, rIdx, fIdx) = Q(:, fIdx).' .* 4 * pi .* (J * sh');
- end
- end
-
- % Calculate the sound field
- if nargout > 1
- sh = sphHarm([rAz; rEl], N);
-
- for fIdx = numel(k):-1:1
- for rIdx = length(rR):-1:1
- waveFieldDecomp(rIdx, :, fIdx) = idsft(squeeze(pCoeffs(:, :, rIdx, fIdx)), sh(rIdx, :));
- end
- end
- end
-
- % Calculate "total" sound field
- if nargout > 2
- waveField = sum(waveFieldDecomp, 2);
- end
-end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldSH.m b/Sound Fields/MATLAB/Functions/ptSrcFieldSH.m
deleted file mode 100644
index 2b0c561..0000000
--- a/Sound Fields/MATLAB/Functions/ptSrcFieldSH.m
+++ /dev/null
@@ -1,169 +0,0 @@
-%% Calculate wavefield generated by point sources in the spherical harmonics domain
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 29/09/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Calculate the wavefield generated by point sources at
-% specified locations in the spherical harmonics domain.
-% --------------------------------------------------
-% Input
-%
-% sPos [numeric]: The position of the sources in Cartesian coordinates. The
-% variable should be a matrix with dimensions 3xN where N
-% is the number of sources and the rows represent the x, y
-% and z coordinates respectively.
-%
-% rPos [numeric]: The position of the receivers in Cartesian coordinates.
-% The variable should be a matrix with dimensions 3xM where
-% M is the number of receivers and the rows represent the
-% x, y and z coordinates respectively.
-%
-% f [numeric]: The frequencies for which to calculate the wave field. It
-% must be a scalar of vector.
-%
-% N [numeric]: The highest order of the sound field. This must be a
-% real non-negative scalar.
-%
-% Q [numeric] (Optional): The source strength of the plane waves. This
-% variable must be either a scalar corresponding to
-% all sources having the same strength, or a matrix
-% with dimensions NxK, where K is the number of
-% frequencies of interest. The strengths can be
-% complex numbers. [Default: 1].
-%
-% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
-%
-% --------------------------------------------------
-% Output
-%
-% pCoeffs [numeric]: The complex pressure coefficients on a sphere at
-% distance from the origin equal to that of the receiver
-% position distances, as generated by point sources at
-% the source positions. The dimensions of the array are
-% ((J + 1)^2)xNxMxF, where J is the highest spherical
-% harmonic order, N the number of sources, M the number
-% of receivers and F the number of frequencies.
-%
-% waveFieldDecomp [numeric]: The complex pressures at the position of the
-% receivers as calculate by the Inverse Discrete
-% Spherical Fourier Transform for each source
-% separately. The dimensions of the array are
-% MxNxF, where again M is the number of
-% receivers, N the number of sources and F the
-% number of frequencies.
-%
-% waveField [numeric]: The sum of the complex pressures at each receiver
-% position (the "total" sound field at each receiver
-% position). This is an matrix with dimensions MxF,
-% where again M is the number of receivers and F the
-% number of frequencies.
-%
-% --------------------------------------------------
-% Notes
-%
-% --------------------------------------------------
-function [pCoeffs, waveFieldDecomp, waveField] = ptSrcFieldSH(sPos, rPos, f, N, Q, c)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(3, 6);
- nargoutchk(0, 3);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(sPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Source positions", 1);
- validateattributes(rPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Receiver positions", 2);
- validateattributes(f, "numeric", {'vector', 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Frequencies", 3);
- validateattributes(N, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'real', 'nonnegative'}, mfilename, "Highest order of the generated soundfield", 4);
-
- % Validate optional arguments
- if nargin > 4 && ~isempty(Q)
- validateattributes(Q, "numeric", {'2d', 'finite', 'nonempty', 'nonnan'}, mfilename, "Source strength(s)", 5);
-
- % Make sure Q has the correct dimensions
- if ~isscalar(Q) && sum((size(Q) ~= [size(sPos, 1), length(f)])) > 0
- error("Source strength(s) must be either a scalar or its length must match the number of sources.");
- end
- else
- Q = 1;
- end
-
- if nargin > 5 && ~isempty(c)
- validateattributes(c, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Speed of propagation", 6);
- else
- c = 343;
- end
-
-
- % ====================================================
- % Pre-process data
- % ====================================================
- % Make sure the source strength has the correct dimensions
- if isscalar(Q)
- Q = Q * ones(size(sPos, 2), length(f));
- end
-
- % Convert coordinates to spherical coordinate system
- [sAz, sEl, sR] = cart2sph(sPos(1, :), sPos(2, :), sPos(3, :));
- sEl = (pi/2) - sEl;
-
- [rAz, rEl, rR] = cart2sph(rPos(1, :), rPos(2, :), rPos(3, :));
- rEl = (pi/2) - rEl;
-
-
- % =============================================
- % Calculate needed parameters
- % ============================================
- % Wavenumbers
- k = 2 * pi * f/c;
-
- % Spherical harmonics corresponding to the directions of the sources
- sh = sphHarm([sAz; sEl], N);
-
-
- % =============================================
- % Calculate variables
- % ============================================
- % Go through the frequencies
- for fIdx = numel(f):-1:1
- % Go through receiver positions
- for rIdx = length(rR):-1:1
- % Got through the sources
- for sIdx = length(sR):-1:1
- if rR(rIdx) < sR(sIdx)
- Js = besselh(0:N, 2, repelem(k(fIdx) * sR(sIdx), N + 1));
- Jr = besselj(0:N, repelem(k(fIdx) * rR(rIdx), N + 1));
- else
- Js = besselj(0:N, repelem(k(fIdx) * sR(sIdx), N + 1));
- Jr = besselh(0:N, 2, repelem(k(fIdx) * rR(rIdx), N + 1));
- end
- J = repelem(Js, 1:2:2 * N + 1) .* repelem(Jr, 1:2:2 * N + 1);
- J = diag(J);
-
- pCoeffs(:, sIdx, rIdx, fIdx) = -4i * pi * Q(sIdx, fIdx) * k(fIdx) * J * sh(sIdx, :)';
- end
- end
- end
-
- % Calculate the sound field components
- if nargout > 1
- sh = sphHarm([rAz; rEl], N);
-
- for fIdx = numel(k):-1:1
- for rIdx = length(rR):-1:1
- waveFieldDecomp(rIdx, :, fIdx) = idsft(squeeze(pCoeffs(:, :, rIdx, fIdx)), sh(rIdx, :));
- end
- end
- end
-
- % Calculate "total" sound field
- if nargout > 2
- waveField = sum(waveFieldDecomp, 2);
- end
-end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldTD.m b/Sound Fields/MATLAB/Functions/ptSrcFieldTD.m
deleted file mode 100644
index 902e92b..0000000
--- a/Sound Fields/MATLAB/Functions/ptSrcFieldTD.m
+++ /dev/null
@@ -1,189 +0,0 @@
-%% Signal generated by point sources in the time-domain
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 29/09/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Calculate the signals generated by ideal point sources in
-% the time-domain.
-% --------------------------------------------------
-% Input
-%
-% sPos [numeric]: The Cartesian coordinates of the source positions. This
-% must be an 3xN matrix where N is the number of sources.
-%
-% rPos [numeric]: The Cartesian coordinates of the receiver positions. This
-% must be an 3xM matrix where M is the number of receivers.
-%
-% fs [numeric]: The sampling frequency. This is required in order to
-% calculate the time-of-flight in samples and must be a
-% positive real scalar.
-%
-% Q [numeric] (Optional): The source signals. This can be either an IxNxJ
-% array, where I is the length of the source
-% signals in samples and J is the number of
-% trials/sound field realisations, a vector of
-% positive real values holding the source strengths
-% (this is just an amplitude scaling factor), or a
-% real positive integer denoting the common
-% strength of all the source signals. [Default: 1].
-%
-% sigLen [numeric] (Optional): The length of the source signals. This must
-% be real positive integer. The generated
-% signals are zero-mean, (approximately)
-% uniformly distributed and normalised to
-% unity. If "Q" is an array, "sigLen" is used
-% to set the length of the signals. If
-% "sigLen" is smaller than I, the signals are
-% truncated and if it is larger the signals
-% are padded with zeros. Leave empty to have
-% the signals unchanged.
-% [Default: 128 or if Q is matrix size(Q, 1)].
-%
-% nTrials [numeric] (Optional): This is the number of trials/sound field
-% realisations that will be generated. If "Q"
-% is provided, this argument is ignored. It
-% must be a positive real scalar.
-% [Default: 1]
-%
-% c [numeric] (Optional): The speed of sound. [Default: 343].
-%
-% --------------------------------------------------
-% Output
-%
-% rSig [numeric]: The signals at the receiver positions for all
-% trials/sournd field realisations. This is an IxMxJ array.
-%
-% rSigMean [numeric]: The signals at the receiver positions averaged over
-% all trials/sound field realisations. This is an IxM
-% matrix.
-%
-% rSigMtx [numeric]: The signals at each receiver position due to each
-% source for each trial/sound field realisation. This is
-% an IxMxNxJ array.
-%
-% rSigMtxMean [numeric]: The signals at each receiver position due to each
-% source averaged over the trials/sound field
-% realisations. This is an IxMxN array.
-%
-% Q [numeric]: The source signals. If Q is provided, the same variable is
-% returned here. If Q is generated internally, the signals are
-% returned.
-%
-% --------------------------------------------------
-% Notes
-%
-% Dependencies: - twoPtDist(): To calculate the distances between sources
-% and receivers.
-%
-% - It is the responsibility of the user to pick (or provide) adequately
-% long source signals
-%
-% --------------------------------------------------
-function [rSig, rSigMean, rSigMtx, rSigMtxMean, Q] = ptSrcFieldTD(sPos, rPos, fs, Q, sigLen, nTrials, c)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(3, 7);
- nargoutchk(0, 5);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(sPos, "numeric", {'2d', 'real', 'nonnan', 'nonempty', 'finite', 'nrows', 3}, mfilename, "Cartesian coordinates of the source positions", 1);
- validateattributes(rPos, "numeric", {'2d', 'real', 'nonnan', 'nonempty', 'finite', 'nrows', 3}, mfilename, "Cartesian coordinates of the receiver positions", 2);
- validateattributes(fs, "numeric", {'scalar', 'real', 'nonnan', 'nonempty', 'finite', 'positive'}, mfilename, "The sampling frequency", 3);
-
- % Validate optional arguments
- if nargin > 3 && ~isempty(Q)
- if isscalar(Q)
- validateattributes(Q, "numeric", {'scalar', 'real', 'nonnan', 'nonempty', 'finite', 'positive', 'integer'}, mfilename, "Length of source signals in samples", 4);
- elseif isvector(Q)
- validateattributes(Q, "numeric", {'vector', 'real', 'nonnan', 'nonempty', 'finite'}, mfilename, "Source signals", 4);
- else
- validateattributes(Q, "numeric", {'3d', 'real', 'nonnan', 'nonempty', 'finite', 'ncols', size(sPos, 2)}, mfilename, "Source signals", 4);
- end
- else
- Q = 1;
- end
-
- if nargin > 4 && ~isempty(sigLen)
- validateattributes(sigLen, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'real', 'positive', 'integer'}, mfilename, "The length of the source signals", 5)
- else
- if ~isvector(Q)
- sigLen = size(Q, 1);
- elseif isvector(Q) && size(sPos, 2) == 1
- sigLen = length(Q);
- else
- sigLen = 128;
- end
- end
-
-
- if ~ismatrix(Q)
- nTrials = size(Q, 3);
- elseif nargin > 5 && ~isempty(nTrials)
- validateattributes(nTrials, "numeric", {'scalar', 'positive', 'nonnan', 'nonempty', 'finite', 'integer'}, mfilename, "Number of trials/sound field realisations", 6);
- else
- nTrials = 1;
- end
-
- if nargin > 6 && ~isempty(c)
- validateattributes(c, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'real', 'positive'}, mfilename, "The speed of sound", 7);
- else
- c = 343;
- end
-
-
- % ====================================================
- % Pre-process data
- % ====================================================
- % Generate source signals
- if isvector(Q) && length(Q) ~= sigLen
- tmp = rand(sigLen, size(sPos, 2), nTrials);
- tmp = tmp - mean(tmp);
- tmp = tmp./max(abs(tmp));
- Q = Q(:).' .* tmp;
- elseif sigLen ~= size(Q, 1)
- if sigLen > size(Q, 1)
- Q = cat(1, Q, zeros(sigLen - size(Q, 1), size(Q, 2), size(Q, 3)));
- else
- Q = Q(1:sigLen, :, :);
- end
- end
-
- % Calculate source-receiver distances
- dist = twoPtDist(sPos, rPos); % Distances
- del = dist/c; % Delays in seconds
-
- % ====================================================
- % Calculate signals
- % ====================================================
- % Go through the trials/sound field realisations
- for jIdx = nTrials:-1:1
- % Go through the sources (calculate for all receivers in one go)
- for sIdx = size(dist, 1):-1:1
- rSigMtx(:, :, sIdx, jIdx) = (1./dist(sIdx, :)) .* delayseq(Q(:, sIdx, jIdx), del(sIdx, :), fs);
- end
- end
-
- % Sum source signals at each receiver position
- rSig = reshape(sum(rSigMtx, 3), size(rSigMtx, 1), size(rSigMtx, 2), size(rSigMtx, 4));
-
- % ====================================================
- % Calculate output arguments
- % ====================================================
- % Mean receiver signal over all trials/sound field realisations
- if nargout > 1
- rSigMean = squeeze(mean(rSig, 3));
- end
-
- % Mean receiver signal due to each source over all trials/sound field realisations
- if nargout > 3
- rSigMtxMean = squeeze(mean(rSigMtx, 4));
- end
-end
\ No newline at end of file
--
GitLab
From 3764919a7760e076b5335969ac68980f4f8af302 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:44:50 +0000
Subject: [PATCH 51/60] Update names of functions in sound fields folder
---
Sound Fields/MATLAB/Functions/extrpFieldSh.m | 156 +++++++++++++++
Sound Fields/MATLAB/Functions/planeWaveSh.m | 160 ++++++++++++++++
Sound Fields/MATLAB/Functions/ptSrcFieldSh.m | 169 +++++++++++++++++
Sound Fields/MATLAB/Functions/ptSrcFieldTd.m | 189 +++++++++++++++++++
4 files changed, 674 insertions(+)
create mode 100644 Sound Fields/MATLAB/Functions/extrpFieldSh.m
create mode 100644 Sound Fields/MATLAB/Functions/planeWaveSh.m
create mode 100644 Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
create mode 100644 Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
diff --git a/Sound Fields/MATLAB/Functions/extrpFieldSh.m b/Sound Fields/MATLAB/Functions/extrpFieldSh.m
new file mode 100644
index 0000000..20ca9b5
--- /dev/null
+++ b/Sound Fields/MATLAB/Functions/extrpFieldSh.m
@@ -0,0 +1,156 @@
+%% Extrapolate the sound field in the Spherical Harmonics domain
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 02/06/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Extrapolate sound field pressure field coefficients in the
+% Spherical Harmonics domain.
+% --------------------------------------------------
+% Input
+%
+% shCoeffs [numeric]: The spherical harmonics pressure coefficients at the
+% measurement position. This must be an NxMxK matrix
+% with N being the number of coefficients per position,
+% M the number of positions/measurements and K the
+% number of frequencies of interest.
+%
+% r [numeric]: This is a vector with at least two elements that correspond
+% to the radial distance component of the spherical
+% coordinates of each position the sound field needs to be
+% extrapolated to. The first element corresponds to the radial
+% distance (radius) of the measurement array.
+%
+% k [numeric]: The wavenumbers of the frequencies of interest. This must be
+% a vector with number of elements equal to K, the frequencies
+% of interest.
+%
+% w [numeric] (Optional): Pressure coefficient weights to be applied before
+% the extrapolation process. This can be used to
+% effectively suppress or use only specific modes
+% for the extrapolation process.
+% [Default: ones(size(shCoeffs, 1), 1)].
+%
+% dir [numeric] (Optional): This signifies the "direction" of the sound
+% field. I.e. whether it is an outgoing or
+% ingoing sound field. In practice, an outgoing
+% sound field is one where the radial distance
+% components of the sources are smaller than the
+% corresponding components of the measurement and
+% extrapolated positions. In the opposite case,
+% where the radial distance components are
+% greater than those of the measurement and
+% extrapolated positions, the field is an ingoing
+% field. The possible values of the argument are
+% "In", "Ingoing", "Converging", "Out",
+% "Outgoing", "Diverging" and are NOT
+% case-sensitive. [Default: "In"].
+%
+% --------------------------------------------------
+% Output
+%
+% extrCoeffs [numeric]: This is an NxMxJxK matrix with the extrapolated
+% coefficients, where J is the number of
+% extrapolation positions and equals length(r) - 1.
+% The coefficients may be Inf if the Bessel (or
+% Hankel) functions at the measured radius exhibit a
+% zero (or infinity for Hankel functions).
+%
+% extrCoeffsReduced [numeric]: This argument holds the same values as
+% "extrCoeffs" but with the singular
+% dimensions being reduced. For example, if
+% extrapolation at a single position is needed
+% the array has dimensions NxMxK. If
+% extrapolating at a single frequency, the
+% dimensions are NxMxJ and if extrapolating at
+% a single position and frequency the result
+% is of dimensions NxM. This is equivalent to
+% the command squeeze(extrCoeffs).
+%
+% --------------------------------------------------
+% Notes
+%
+% --------------------------------------------------
+function [extrCoeffs, extrCoeffsReduced] = extrpFieldSH(shCoeffs, r, k, w, dir)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 5);
+ nargoutchk(0, 1);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(shCoeffs, "numeric", {'2d', 'nonempty', 'nonnan', 'finite'}, mfilename, "The spherical harmonic pressure coefficients at the measurement positions", 1);
+ validateattributes(r, "numeric", {'vector', 'nonempty', 'nonnan', 'finite', 'nonnegative'}, mfilename, "The radial distance components of all positions", 2);
+
+ % Check if we have AT LEAST two distances
+ if numel(r) < 2
+ error("You must provide at least two distances, with the first being the radius of the measurement array");
+ end
+
+ validateattributes(k, "numeric", {'vector', 'nonempty', 'nonnan', 'finite', 'nonnegative', 'real'}, mfilename, "The wavenumbers of the frequencies of interest", 3);
+
+ % Validate optional arguments
+ if nargin > 3 && ~isempty(w)
+ validateattributes(w, "numeric", {'vector', 'nonempty', 'nonnan', 'numel', size(shCoeffs, 1), 'real', 'finite'}, mfilename, "Coefficient weights", 4);
+ else
+ w = ones(size(shCoeffs, 1), 1);
+ end
+
+ if nargin > 4 && ~isempty(dir)
+ validateattributes(dir, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "The 'direction' of the sound field", 5);
+ validatestring(dir, ["In", "Ingoing", "Converging", "Out", "Outgoing", "Diverging"], mfilename, "The 'direction' of the sound field", 5);
+
+ if sum(strcmpi(dir, ["In", "Ingoing", "Converging"])) > 0
+ dir = "In";
+ end
+ else
+ dir = "In";
+ end
+
+ % =============================================
+ % Calculate parameters
+ % =============================================
+ N = sqrt(size(shCoeffs, 1)) - 1; % Highest spherical harmonics order
+ nDeg = 1:2:2 * N + 1; % Number of degrees/modes per order
+
+
+ % =============================================
+ % Extrapolate coefficients
+ % =============================================
+ % Create some matrices to ease calculations
+ if strcmpi(dir, "in")
+ Jref = besselj(repmat(0:N, length(k), 1), repmat(r(1) * k(:), 1, N + 1));
+ else
+ Jref = besselh(repmat(0:N, length(k), 1), repmat(r(1) * k(:), 1, N + 1));
+ end
+
+ % Go through the frequencies
+ for fIdx = size(shCoeffs, 3):-1:1
+ % Go through the measurements/positions
+ for rIdx = length(r) - 1:-1:1
+ tmpIdx = rIdx + 1;
+ if strcmpi(dir, "in")
+ J = besselj(0:N, repmat(k(fIdx) * r(tmpIdx), 1, N + 1));
+ else
+ J = besselh(0:N, repmat(k(fIdx) * r(tmpIdx), 1, N + 1));
+ end
+
+ % Replicate and weight J coefficients
+ J = w(:).' .* repelem(J./Jref(fIdx, :), nDeg);
+ J = diag(J);
+
+ % Calculate extrapolated pressure coefficients
+ extrCoeffs(:, :, rIdx, fIdx) = J * shCoeffs(:, :, fIdx);
+ end
+ end
+
+ if nargout > 1
+ extrCoeffsReduced = squeeze(extrCoeffs);
+ end
+end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/planeWaveSh.m b/Sound Fields/MATLAB/Functions/planeWaveSh.m
new file mode 100644
index 0000000..3dd6d06
--- /dev/null
+++ b/Sound Fields/MATLAB/Functions/planeWaveSh.m
@@ -0,0 +1,160 @@
+%% Calculate wavefield generated by plane waves in the spherical harmonics domain
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 29/09/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the wavefield generated by plane waves at
+% specified directions in the spherical harmonics domain.
+% --------------------------------------------------
+% Input
+%
+% sDir [numeric]: The directions of the sources in spherical coordinates.
+% The variable must be a matrix with dimensions 2xN where N
+% is the number of plane waves and the rows represent the
+% azimuth and inclination respectively. The values must be
+% in radians.
+%
+% rPos [numeric]: The position of the receivers in Cartesian coordinates.
+% The variable must be a matrix with dimensions 3xM where
+% M is the number of receivers and the rows represent the
+% x, y and z coordinates respectively.
+%
+% f [numeric]: The frequencies for which to calculate the wave field. It
+% must be a scalar or vector.
+%
+% N [numeric]: The highest order of the sound field. This must be a
+% real non-negative scalar.
+%
+% Q [numeric] (Optional): The source strength of the plane waves. This
+% variable must be either a scalar corresponding to
+% all sources having the same strength, or a matrix
+% with dimensions NxK, where K is the number of
+% frequencies of interest. The strengths can be
+% complex numbers. [Default: 1].
+%
+% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
+%
+% --------------------------------------------------
+% Output
+%
+% pCoeffs [numeric]: The complex pressure coefficients on a sphere at
+% distance from the origin equal to that of the receiver
+% position distances, as generated by plane waves from
+% directions given in "sDir". The dimensions of the
+% array are ((J + 1)^2)xNxMxF, where J is the highest
+% spherical harmonic order, N the number of plane waves,
+% M the number of receivers and F the number of
+% frequencies.
+%
+% waveFieldDecomp [numeric]: The complex pressures at the position of the
+% receivers as calculate by the Inverse Discrete
+% Spherical Fourier Transform for each plane
+% wave separately. The dimensions of the array
+% are MxNxF, where again M is the number of
+% receivers, N the number of plane waves and F
+% the number of frequencies.
+%
+% waveField [numeric]: The sum of the complex pressures at each receiver
+% position (the "total" sound field at each receiver
+% position). This is an matrix with dimensions MxF,
+% where again M is the number of receivers and F the
+% number of frequencies.
+%
+% --------------------------------------------------
+% Notes
+%
+% --------------------------------------------------
+function [pCoeffs, waveFieldDecomp, waveField] = planeWaveSH(sDir, rPos, f, N, Q, c)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 6);
+ nargoutchk(0, 3);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(sDir, "numeric", {'2d', 'rows', 2, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Source positions", 1);
+ validateattributes(rPos, "numeric", {'2d', 'rows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Receiver positions", 2);
+ validateattributes(f, "numeric", {'vector', 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Frequencies", 3);
+ validateattributes(N, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'real', 'nonnegative'}, mfilename, "Highest order of the generated soundfield", 4);
+
+ % Validate optional arguments
+ if nargin > 4 && ~isempty(Q)
+ validateattributes(Q, "numeric", {'2d', 'finite', 'nonempty', 'nonnan'}, mfilename, "Source strength(s)", 5);
+
+ % Make sure Q has the correct dimensions
+ if ~isscalar(Q) && sum((size(Q) ~= [size(sDir, 2), length(f)])) > 0
+ error("Source strength(s) must be either a scalar or its length must match the number of sources.");
+ end
+ else
+ Q = 1;
+ end
+
+ if nargin > 5 && ~isempty(c)
+ validateattributes(c, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Speed of propagation", 6);
+ else
+ c = 343;
+ end
+
+
+ % ====================================================
+ % Pre-process data
+ % ====================================================
+ % Make sure the source strength has the correct dimensions
+ if isscalar(Q)
+ Q = Q * ones(size(sDir, 2), length(f));
+ end
+
+ % Convert coordinates to spherical coordinate system
+ sDir(2, :) = (pi/2) - sDir(2, :);
+
+ [rAz, rEl, rR] = cart2sph(rPos(1, :), rPos(2, :), rPos(3, :));
+ rEl = (pi/2) - rEl;
+
+
+ % =============================================
+ % Calculate needed parameters
+ % ============================================
+ % Wavenumbers
+ k = 2 * pi * f/c;
+
+ % Spherical harmonics corresponding to the directions of the sources
+ sh = sphHarm(sDir, N);
+
+
+ % =============================================
+ % Calculate variables
+ % ============================================
+ % Go through the frequencies
+ for fIdx = numel(f):-1:1
+ % Go through receiver positions
+ for rIdx = length(rR):-1:1
+ J = 1i.^(0:N) .* besselj(0:N, repelem(k(fIdx) * rR(rIdx), N + 1));
+ J = diag(repelem(J, 1:2:2 * N + 1));
+
+ pCoeffs(:, :, rIdx, fIdx) = Q(:, fIdx).' .* 4 * pi .* (J * sh');
+ end
+ end
+
+ % Calculate the sound field
+ if nargout > 1
+ sh = sphHarm([rAz; rEl], N);
+
+ for fIdx = numel(k):-1:1
+ for rIdx = length(rR):-1:1
+ waveFieldDecomp(rIdx, :, fIdx) = idsft(squeeze(pCoeffs(:, :, rIdx, fIdx)), sh(rIdx, :));
+ end
+ end
+ end
+
+ % Calculate "total" sound field
+ if nargout > 2
+ waveField = sum(waveFieldDecomp, 2);
+ end
+end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m b/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
new file mode 100644
index 0000000..2b0c561
--- /dev/null
+++ b/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
@@ -0,0 +1,169 @@
+%% Calculate wavefield generated by point sources in the spherical harmonics domain
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 29/09/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the wavefield generated by point sources at
+% specified locations in the spherical harmonics domain.
+% --------------------------------------------------
+% Input
+%
+% sPos [numeric]: The position of the sources in Cartesian coordinates. The
+% variable should be a matrix with dimensions 3xN where N
+% is the number of sources and the rows represent the x, y
+% and z coordinates respectively.
+%
+% rPos [numeric]: The position of the receivers in Cartesian coordinates.
+% The variable should be a matrix with dimensions 3xM where
+% M is the number of receivers and the rows represent the
+% x, y and z coordinates respectively.
+%
+% f [numeric]: The frequencies for which to calculate the wave field. It
+% must be a scalar of vector.
+%
+% N [numeric]: The highest order of the sound field. This must be a
+% real non-negative scalar.
+%
+% Q [numeric] (Optional): The source strength of the plane waves. This
+% variable must be either a scalar corresponding to
+% all sources having the same strength, or a matrix
+% with dimensions NxK, where K is the number of
+% frequencies of interest. The strengths can be
+% complex numbers. [Default: 1].
+%
+% c [numeric] (Optional): The speed of sound in m/s. [Default: 343].
+%
+% --------------------------------------------------
+% Output
+%
+% pCoeffs [numeric]: The complex pressure coefficients on a sphere at
+% distance from the origin equal to that of the receiver
+% position distances, as generated by point sources at
+% the source positions. The dimensions of the array are
+% ((J + 1)^2)xNxMxF, where J is the highest spherical
+% harmonic order, N the number of sources, M the number
+% of receivers and F the number of frequencies.
+%
+% waveFieldDecomp [numeric]: The complex pressures at the position of the
+% receivers as calculate by the Inverse Discrete
+% Spherical Fourier Transform for each source
+% separately. The dimensions of the array are
+% MxNxF, where again M is the number of
+% receivers, N the number of sources and F the
+% number of frequencies.
+%
+% waveField [numeric]: The sum of the complex pressures at each receiver
+% position (the "total" sound field at each receiver
+% position). This is an matrix with dimensions MxF,
+% where again M is the number of receivers and F the
+% number of frequencies.
+%
+% --------------------------------------------------
+% Notes
+%
+% --------------------------------------------------
+function [pCoeffs, waveFieldDecomp, waveField] = ptSrcFieldSH(sPos, rPos, f, N, Q, c)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 6);
+ nargoutchk(0, 3);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(sPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Source positions", 1);
+ validateattributes(rPos, "numeric", {'2d', 'nrows', 3, 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Receiver positions", 2);
+ validateattributes(f, "numeric", {'vector', 'finite', 'nonnan', 'nonempty', 'real'}, mfilename, "Frequencies", 3);
+ validateattributes(N, "numeric", {'scalar', 'finite', 'nonempty', 'nonnan', 'real', 'nonnegative'}, mfilename, "Highest order of the generated soundfield", 4);
+
+ % Validate optional arguments
+ if nargin > 4 && ~isempty(Q)
+ validateattributes(Q, "numeric", {'2d', 'finite', 'nonempty', 'nonnan'}, mfilename, "Source strength(s)", 5);
+
+ % Make sure Q has the correct dimensions
+ if ~isscalar(Q) && sum((size(Q) ~= [size(sPos, 1), length(f)])) > 0
+ error("Source strength(s) must be either a scalar or its length must match the number of sources.");
+ end
+ else
+ Q = 1;
+ end
+
+ if nargin > 5 && ~isempty(c)
+ validateattributes(c, "numeric", {'scalar', 'finite', 'positive', 'nonnan', 'nonempty', 'real'}, mfilename, "Speed of propagation", 6);
+ else
+ c = 343;
+ end
+
+
+ % ====================================================
+ % Pre-process data
+ % ====================================================
+ % Make sure the source strength has the correct dimensions
+ if isscalar(Q)
+ Q = Q * ones(size(sPos, 2), length(f));
+ end
+
+ % Convert coordinates to spherical coordinate system
+ [sAz, sEl, sR] = cart2sph(sPos(1, :), sPos(2, :), sPos(3, :));
+ sEl = (pi/2) - sEl;
+
+ [rAz, rEl, rR] = cart2sph(rPos(1, :), rPos(2, :), rPos(3, :));
+ rEl = (pi/2) - rEl;
+
+
+ % =============================================
+ % Calculate needed parameters
+ % ============================================
+ % Wavenumbers
+ k = 2 * pi * f/c;
+
+ % Spherical harmonics corresponding to the directions of the sources
+ sh = sphHarm([sAz; sEl], N);
+
+
+ % =============================================
+ % Calculate variables
+ % ============================================
+ % Go through the frequencies
+ for fIdx = numel(f):-1:1
+ % Go through receiver positions
+ for rIdx = length(rR):-1:1
+ % Got through the sources
+ for sIdx = length(sR):-1:1
+ if rR(rIdx) < sR(sIdx)
+ Js = besselh(0:N, 2, repelem(k(fIdx) * sR(sIdx), N + 1));
+ Jr = besselj(0:N, repelem(k(fIdx) * rR(rIdx), N + 1));
+ else
+ Js = besselj(0:N, repelem(k(fIdx) * sR(sIdx), N + 1));
+ Jr = besselh(0:N, 2, repelem(k(fIdx) * rR(rIdx), N + 1));
+ end
+ J = repelem(Js, 1:2:2 * N + 1) .* repelem(Jr, 1:2:2 * N + 1);
+ J = diag(J);
+
+ pCoeffs(:, sIdx, rIdx, fIdx) = -4i * pi * Q(sIdx, fIdx) * k(fIdx) * J * sh(sIdx, :)';
+ end
+ end
+ end
+
+ % Calculate the sound field components
+ if nargout > 1
+ sh = sphHarm([rAz; rEl], N);
+
+ for fIdx = numel(k):-1:1
+ for rIdx = length(rR):-1:1
+ waveFieldDecomp(rIdx, :, fIdx) = idsft(squeeze(pCoeffs(:, :, rIdx, fIdx)), sh(rIdx, :));
+ end
+ end
+ end
+
+ % Calculate "total" sound field
+ if nargout > 2
+ waveField = sum(waveFieldDecomp, 2);
+ end
+end
\ No newline at end of file
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m b/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
new file mode 100644
index 0000000..902e92b
--- /dev/null
+++ b/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
@@ -0,0 +1,189 @@
+%% Signal generated by point sources in the time-domain
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 29/09/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate the signals generated by ideal point sources in
+% the time-domain.
+% --------------------------------------------------
+% Input
+%
+% sPos [numeric]: The Cartesian coordinates of the source positions. This
+% must be an 3xN matrix where N is the number of sources.
+%
+% rPos [numeric]: The Cartesian coordinates of the receiver positions. This
+% must be an 3xM matrix where M is the number of receivers.
+%
+% fs [numeric]: The sampling frequency. This is required in order to
+% calculate the time-of-flight in samples and must be a
+% positive real scalar.
+%
+% Q [numeric] (Optional): The source signals. This can be either an IxNxJ
+% array, where I is the length of the source
+% signals in samples and J is the number of
+% trials/sound field realisations, a vector of
+% positive real values holding the source strengths
+% (this is just an amplitude scaling factor), or a
+% real positive integer denoting the common
+% strength of all the source signals. [Default: 1].
+%
+% sigLen [numeric] (Optional): The length of the source signals. This must
+% be real positive integer. The generated
+% signals are zero-mean, (approximately)
+% uniformly distributed and normalised to
+% unity. If "Q" is an array, "sigLen" is used
+% to set the length of the signals. If
+% "sigLen" is smaller than I, the signals are
+% truncated and if it is larger the signals
+% are padded with zeros. Leave empty to have
+% the signals unchanged.
+% [Default: 128 or if Q is matrix size(Q, 1)].
+%
+% nTrials [numeric] (Optional): This is the number of trials/sound field
+% realisations that will be generated. If "Q"
+% is provided, this argument is ignored. It
+% must be a positive real scalar.
+% [Default: 1]
+%
+% c [numeric] (Optional): The speed of sound. [Default: 343].
+%
+% --------------------------------------------------
+% Output
+%
+% rSig [numeric]: The signals at the receiver positions for all
+% trials/sournd field realisations. This is an IxMxJ array.
+%
+% rSigMean [numeric]: The signals at the receiver positions averaged over
+% all trials/sound field realisations. This is an IxM
+% matrix.
+%
+% rSigMtx [numeric]: The signals at each receiver position due to each
+% source for each trial/sound field realisation. This is
+% an IxMxNxJ array.
+%
+% rSigMtxMean [numeric]: The signals at each receiver position due to each
+% source averaged over the trials/sound field
+% realisations. This is an IxMxN array.
+%
+% Q [numeric]: The source signals. If Q is provided, the same variable is
+% returned here. If Q is generated internally, the signals are
+% returned.
+%
+% --------------------------------------------------
+% Notes
+%
+% Dependencies: - twoPtDist(): To calculate the distances between sources
+% and receivers.
+%
+% - It is the responsibility of the user to pick (or provide) adequately
+% long source signals
+%
+% --------------------------------------------------
+function [rSig, rSigMean, rSigMtx, rSigMtxMean, Q] = ptSrcFieldTD(sPos, rPos, fs, Q, sigLen, nTrials, c)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(3, 7);
+ nargoutchk(0, 5);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(sPos, "numeric", {'2d', 'real', 'nonnan', 'nonempty', 'finite', 'nrows', 3}, mfilename, "Cartesian coordinates of the source positions", 1);
+ validateattributes(rPos, "numeric", {'2d', 'real', 'nonnan', 'nonempty', 'finite', 'nrows', 3}, mfilename, "Cartesian coordinates of the receiver positions", 2);
+ validateattributes(fs, "numeric", {'scalar', 'real', 'nonnan', 'nonempty', 'finite', 'positive'}, mfilename, "The sampling frequency", 3);
+
+ % Validate optional arguments
+ if nargin > 3 && ~isempty(Q)
+ if isscalar(Q)
+ validateattributes(Q, "numeric", {'scalar', 'real', 'nonnan', 'nonempty', 'finite', 'positive', 'integer'}, mfilename, "Length of source signals in samples", 4);
+ elseif isvector(Q)
+ validateattributes(Q, "numeric", {'vector', 'real', 'nonnan', 'nonempty', 'finite'}, mfilename, "Source signals", 4);
+ else
+ validateattributes(Q, "numeric", {'3d', 'real', 'nonnan', 'nonempty', 'finite', 'ncols', size(sPos, 2)}, mfilename, "Source signals", 4);
+ end
+ else
+ Q = 1;
+ end
+
+ if nargin > 4 && ~isempty(sigLen)
+ validateattributes(sigLen, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'real', 'positive', 'integer'}, mfilename, "The length of the source signals", 5)
+ else
+ if ~isvector(Q)
+ sigLen = size(Q, 1);
+ elseif isvector(Q) && size(sPos, 2) == 1
+ sigLen = length(Q);
+ else
+ sigLen = 128;
+ end
+ end
+
+
+ if ~ismatrix(Q)
+ nTrials = size(Q, 3);
+ elseif nargin > 5 && ~isempty(nTrials)
+ validateattributes(nTrials, "numeric", {'scalar', 'positive', 'nonnan', 'nonempty', 'finite', 'integer'}, mfilename, "Number of trials/sound field realisations", 6);
+ else
+ nTrials = 1;
+ end
+
+ if nargin > 6 && ~isempty(c)
+ validateattributes(c, "numeric", {'scalar', 'nonempty', 'nonnan', 'finite', 'real', 'positive'}, mfilename, "The speed of sound", 7);
+ else
+ c = 343;
+ end
+
+
+ % ====================================================
+ % Pre-process data
+ % ====================================================
+ % Generate source signals
+ if isvector(Q) && length(Q) ~= sigLen
+ tmp = rand(sigLen, size(sPos, 2), nTrials);
+ tmp = tmp - mean(tmp);
+ tmp = tmp./max(abs(tmp));
+ Q = Q(:).' .* tmp;
+ elseif sigLen ~= size(Q, 1)
+ if sigLen > size(Q, 1)
+ Q = cat(1, Q, zeros(sigLen - size(Q, 1), size(Q, 2), size(Q, 3)));
+ else
+ Q = Q(1:sigLen, :, :);
+ end
+ end
+
+ % Calculate source-receiver distances
+ dist = twoPtDist(sPos, rPos); % Distances
+ del = dist/c; % Delays in seconds
+
+ % ====================================================
+ % Calculate signals
+ % ====================================================
+ % Go through the trials/sound field realisations
+ for jIdx = nTrials:-1:1
+ % Go through the sources (calculate for all receivers in one go)
+ for sIdx = size(dist, 1):-1:1
+ rSigMtx(:, :, sIdx, jIdx) = (1./dist(sIdx, :)) .* delayseq(Q(:, sIdx, jIdx), del(sIdx, :), fs);
+ end
+ end
+
+ % Sum source signals at each receiver position
+ rSig = reshape(sum(rSigMtx, 3), size(rSigMtx, 1), size(rSigMtx, 2), size(rSigMtx, 4));
+
+ % ====================================================
+ % Calculate output arguments
+ % ====================================================
+ % Mean receiver signal over all trials/sound field realisations
+ if nargout > 1
+ rSigMean = squeeze(mean(rSig, 3));
+ end
+
+ % Mean receiver signal due to each source over all trials/sound field realisations
+ if nargout > 3
+ rSigMtxMean = squeeze(mean(rSigMtx, 4));
+ end
+end
\ No newline at end of file
--
GitLab
From 584416c6e400fb3e270ecfafbef72786ad2454d1 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:46:04 +0000
Subject: [PATCH 52/60] Remove virtMicGeo.m to rename
---
.../Geometries/MATLAB/Functions/virtMicGeo.m | 205 ------------------
1 file changed, 205 deletions(-)
delete mode 100644 Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
diff --git a/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m b/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
deleted file mode 100644
index 3b18c17..0000000
--- a/Utilities/Geometries/MATLAB/Functions/virtMicGeo.m
+++ /dev/null
@@ -1,205 +0,0 @@
-%% Calculate virtual microphone positions for specific geometries
-% --------------------------------------------------
-% Author: Achilles Kappis
-% e-mail: axilleaz@protonmail.com
-%
-% Date: 15/12/2024 (DD/MM/YYYY)
-%
-% Copyright: MIT
-% --------------------------------------------------
-% Functionality: Calculate virtual microphone positions for specific
-% geometries.
-% --------------------------------------------------
-% Input
-%
-% gType [char/string]: The type of the source geometry. At the moment the
-% 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)].
-%
-% 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)].
-%
-% trans [numeric] (Optional): This is a real vector with three elements,
-% representing the translation of the geometry
-% along the Cartesian axes.
-% [Default: zeros(3, 1)].
-%
-% rot [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)].
-%
-% rotOrd [char/string] (Optional): The order the rotations will be applied.
-% This can be some permutation of the
-% string "xyz" and is not case-sensistive.
-% The order corresponds to multiplication
-% on the left (i.e. R * v, where "R" is
-% the rotation matrix and "v" a vector to
-% be rotated). [Default: "xyz"].
-%
-% --------------------------------------------------
-% Output
-%
-% vPos [numeric]: The matrix with the source coordinates in the Cartesian
-% system. The matrix has dimensions 3xN, where N is the
-% number of sources and the first dimension corresponds to
-% the x, y and z coodrinates.
-%
-% 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
-% other than "Grid", this is an empty array.
-%
-% --------------------------------------------------
-% Notes
-%
-% - 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.
-%
-% --------------------------------------------------
-function [vPos, vPosMesh] = virtMicGeo(gType, geoDim, nSens, trans, rot, rotOrd)
- % ====================================================
- % Check for number of arguments
- % ====================================================
- narginchk(1, 6);
- nargoutchk(0, 2);
-
- % ====================================================
- % Validate input arguments
- % ====================================================
- % Validate mandatory arguments
- validateattributes(gType, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "Geometry type", 1);
- validatestring(gType, ["Single", "Dual", "Array", "Cube"], mfilename, "Geometry type", 1);
-
- % Validate optional arguments
- if nargin > 1 && ~isempty(geoDim) && ~strcmpi(gType, "Single")
- validateattributes(geoDim, "numeric", {'vector', 'real', 'nonnan', 'finite', 'nonempty', 'nonnegative'}, 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
- 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 > 3 && ~isempty(trans)
- validateattributes(trans, "numeric", {'vector', 'real', 'nonnan', 'finite', 'nonempty', 'numel', 3}, mfilename, "Translation of the geometry", 4);
- else
- trans = zeros(3, 1);
- end
-
- if nargin > 4 && ~isempty(rot) && ~strcmpi(gType, "Single")
- validateattributes(rot, {'numeric'}, {'vector', 'nonempty', 'nonnan', 'finite', 'real', 'numel', 3}, mfilename, "Geometry rotation around the Cartesian axes", 5);
- else
- rot = zeros(3, 1);
- end
-
- if nargin < 6 || (nargin > 5 && isempty(rotOrd))
- rotOrd = "xyz";
- end
-
-
- % ====================================================
- % Calculate virtual microphone positions
- % ====================================================
- switch lower(gType)
- case "single"
- % Create a single point at offset coordinates and return
- vPos = zeros(3, 1);
- case "dual"
- vPos = [-geoDim(1)/2; 0; 0, ...
- geoDim(1)/2; 0; 0];
- case "array"
- vPos = linspace(-geoDim(1)/2, geoDim(1)/2, nSens(1));
- vPos = [vPos; zeros(2, numel(vPos))];
- case "cube"
- % Make sure 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 and translate
- vPos = rotMat3d(rot(1), rot(2), rot(3), rotOrd, true) * vPos + trans(:);
-
- % For "Cube" geometry provide the coordinates in a "mesh format"
- if nargout > 1 && strcmpi(gType, "Cube")
- vPosMesh = reshape(vPos.', size(x, 1), size(x, 2), [], 3);
- elseif nargout > 1
- vPosMesh = [];
- end
-end
\ No newline at end of file
--
GitLab
From 22057ff1716bf0ff7eebcaf98f4e816f1ad18c8e Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:46:54 +0000
Subject: [PATCH 53/60] Change name of virtMicGeo.m to vMicGeo.m
---
.../Geometries/MATLAB/Functions/vMicGeo.m | 205 ++++++++++++++++++
1 file changed, 205 insertions(+)
create mode 100644 Utilities/Geometries/MATLAB/Functions/vMicGeo.m
diff --git a/Utilities/Geometries/MATLAB/Functions/vMicGeo.m b/Utilities/Geometries/MATLAB/Functions/vMicGeo.m
new file mode 100644
index 0000000..3b18c17
--- /dev/null
+++ b/Utilities/Geometries/MATLAB/Functions/vMicGeo.m
@@ -0,0 +1,205 @@
+%% Calculate virtual microphone positions for specific geometries
+% --------------------------------------------------
+% Author: Achilles Kappis
+% e-mail: axilleaz@protonmail.com
+%
+% Date: 15/12/2024 (DD/MM/YYYY)
+%
+% Copyright: MIT
+% --------------------------------------------------
+% Functionality: Calculate virtual microphone positions for specific
+% geometries.
+% --------------------------------------------------
+% Input
+%
+% gType [char/string]: The type of the source geometry. At the moment the
+% 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)].
+%
+% 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)].
+%
+% trans [numeric] (Optional): This is a real vector with three elements,
+% representing the translation of the geometry
+% along the Cartesian axes.
+% [Default: zeros(3, 1)].
+%
+% rot [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)].
+%
+% rotOrd [char/string] (Optional): The order the rotations will be applied.
+% This can be some permutation of the
+% string "xyz" and is not case-sensistive.
+% The order corresponds to multiplication
+% on the left (i.e. R * v, where "R" is
+% the rotation matrix and "v" a vector to
+% be rotated). [Default: "xyz"].
+%
+% --------------------------------------------------
+% Output
+%
+% vPos [numeric]: The matrix with the source coordinates in the Cartesian
+% system. The matrix has dimensions 3xN, where N is the
+% number of sources and the first dimension corresponds to
+% the x, y and z coodrinates.
+%
+% 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
+% other than "Grid", this is an empty array.
+%
+% --------------------------------------------------
+% Notes
+%
+% - 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.
+%
+% --------------------------------------------------
+function [vPos, vPosMesh] = virtMicGeo(gType, geoDim, nSens, trans, rot, rotOrd)
+ % ====================================================
+ % Check for number of arguments
+ % ====================================================
+ narginchk(1, 6);
+ nargoutchk(0, 2);
+
+ % ====================================================
+ % Validate input arguments
+ % ====================================================
+ % Validate mandatory arguments
+ validateattributes(gType, {'char', 'string'}, {'scalartext', 'nonempty'}, mfilename, "Geometry type", 1);
+ validatestring(gType, ["Single", "Dual", "Array", "Cube"], mfilename, "Geometry type", 1);
+
+ % Validate optional arguments
+ if nargin > 1 && ~isempty(geoDim) && ~strcmpi(gType, "Single")
+ validateattributes(geoDim, "numeric", {'vector', 'real', 'nonnan', 'finite', 'nonempty', 'nonnegative'}, 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
+ 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 > 3 && ~isempty(trans)
+ validateattributes(trans, "numeric", {'vector', 'real', 'nonnan', 'finite', 'nonempty', 'numel', 3}, mfilename, "Translation of the geometry", 4);
+ else
+ trans = zeros(3, 1);
+ end
+
+ if nargin > 4 && ~isempty(rot) && ~strcmpi(gType, "Single")
+ validateattributes(rot, {'numeric'}, {'vector', 'nonempty', 'nonnan', 'finite', 'real', 'numel', 3}, mfilename, "Geometry rotation around the Cartesian axes", 5);
+ else
+ rot = zeros(3, 1);
+ end
+
+ if nargin < 6 || (nargin > 5 && isempty(rotOrd))
+ rotOrd = "xyz";
+ end
+
+
+ % ====================================================
+ % Calculate virtual microphone positions
+ % ====================================================
+ switch lower(gType)
+ case "single"
+ % Create a single point at offset coordinates and return
+ vPos = zeros(3, 1);
+ case "dual"
+ vPos = [-geoDim(1)/2; 0; 0, ...
+ geoDim(1)/2; 0; 0];
+ case "array"
+ vPos = linspace(-geoDim(1)/2, geoDim(1)/2, nSens(1));
+ vPos = [vPos; zeros(2, numel(vPos))];
+ case "cube"
+ % Make sure 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 and translate
+ vPos = rotMat3d(rot(1), rot(2), rot(3), rotOrd, true) * vPos + trans(:);
+
+ % For "Cube" geometry provide the coordinates in a "mesh format"
+ if nargout > 1 && strcmpi(gType, "Cube")
+ vPosMesh = reshape(vPos.', size(x, 1), size(x, 2), [], 3);
+ elseif nargout > 1
+ vPosMesh = [];
+ end
+end
\ No newline at end of file
--
GitLab
From afb87528f30fd0dd25ac40928b30c3f093eca304 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sat, 4 Jan 2025 15:47:52 +0000
Subject: [PATCH 54/60] Update CHANGELOG.md with changes in function naming
---
CHANGELOG.md | 8 +++++++-
1 file changed, 7 insertions(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 63ae65c..311bf09 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -19,7 +19,12 @@
**Sound Fields**\
\+ Add `lineSrcField()` to calculate the sound field due to line sources in the frequency domain (tonal).
-\* Parallelize the for-loops in `ptSrcField()` and `lineSrcField()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.
+\* Parallelize the for-loops in `ptSrcField()` and `lineSrcField()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
+\* Rename the following functions (this **breaks** backward compatibility):
+ - `extrpFieldSH()` to `extrpFieldSh()`
+ - `planeWaveSH()` to `planeWaveSh()`
+ - `ptSrcFieldSH()` to `ptSrcFieldSh()`
+ - `ptSrcFieldTD()` to `ptSrcFieldTd()`
**Utilities - Geometries**\
\* Fix error in `rcvGeo()` applying incorrect rotation to the elements in the Uniform Circular Array.\
@@ -31,6 +36,7 @@
\* Change function name from `obsFiltTD()` to `obsFiltTd()`.\
\* Change function name from `obsFiltEstTD()` to `obsFiltEstTd()`.\
\* Update `obsFilt()` to return the condition of the perturbed power spectral density matrix.
+\* Change name of `virtMicGeo()` to `vMicGeo()`. This **breaks** backward compatibility.
-------------------------------------------
--
GitLab
From ee83b113bcfe1196023a5b00e30d2e372999e317 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:33:29 +0000
Subject: [PATCH 55/60] Update the names of the functions whose files where
renamed to match the filename
---
.../Array Processing/MATLAB/Functions/arrManTd.m | 4 ++--
Sound Fields/MATLAB/Functions/extrpFieldSh.m | 4 ++--
Sound Fields/MATLAB/Functions/planeWaveSh.m | 4 ++--
Sound Fields/MATLAB/Functions/ptSrcFieldSh.m | 4 ++--
Sound Fields/MATLAB/Functions/ptSrcFieldTd.m | 4 ++--
Utilities/Geometries/MATLAB/Functions/vMicGeo.m | 4 ++--
.../MATLAB/Functions/obsFiltEstTd.m | 4 ++--
.../Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m | 4 ++--
8 files changed, 16 insertions(+), 16 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m b/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
index 1639b2a..1c92e62 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/arrManTd.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 02/11/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -90,7 +90,7 @@
% - Dependencies: * winSincFracDel() to calculate the fractional delay
% filters.
% --------------------------------------------------
-function [am, maxDel, amMtx] = arrManTD(mPos, dirs, filtLen, delLen, delWin, fs, c)
+function [am, maxDel, amMtx] = arrManTd(mPos, dirs, filtLen, delLen, delWin, fs, c)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Sound Fields/MATLAB/Functions/extrpFieldSh.m b/Sound Fields/MATLAB/Functions/extrpFieldSh.m
index 20ca9b5..124b43b 100644
--- a/Sound Fields/MATLAB/Functions/extrpFieldSh.m
+++ b/Sound Fields/MATLAB/Functions/extrpFieldSh.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 02/06/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -74,7 +74,7 @@
% Notes
%
% --------------------------------------------------
-function [extrCoeffs, extrCoeffsReduced] = extrpFieldSH(shCoeffs, r, k, w, dir)
+function [extrCoeffs, extrCoeffsReduced] = extrpFieldSh(shCoeffs, r, k, w, dir)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Sound Fields/MATLAB/Functions/planeWaveSh.m b/Sound Fields/MATLAB/Functions/planeWaveSh.m
index 3dd6d06..8f48567 100644
--- a/Sound Fields/MATLAB/Functions/planeWaveSh.m
+++ b/Sound Fields/MATLAB/Functions/planeWaveSh.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 29/09/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -68,7 +68,7 @@
% Notes
%
% --------------------------------------------------
-function [pCoeffs, waveFieldDecomp, waveField] = planeWaveSH(sDir, rPos, f, N, Q, c)
+function [pCoeffs, waveFieldDecomp, waveField] = planeWaveSh(sDir, rPos, f, N, Q, c)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m b/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
index 2b0c561..07fec1d 100644
--- a/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
+++ b/Sound Fields/MATLAB/Functions/ptSrcFieldSh.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 29/09/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -66,7 +66,7 @@
% Notes
%
% --------------------------------------------------
-function [pCoeffs, waveFieldDecomp, waveField] = ptSrcFieldSH(sPos, rPos, f, N, Q, c)
+function [pCoeffs, waveFieldDecomp, waveField] = ptSrcFieldSh(sPos, rPos, f, N, Q, c)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m b/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
index 902e92b..cd28066 100644
--- a/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
+++ b/Sound Fields/MATLAB/Functions/ptSrcFieldTd.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 29/09/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -83,7 +83,7 @@
% long source signals
%
% --------------------------------------------------
-function [rSig, rSigMean, rSigMtx, rSigMtxMean, Q] = ptSrcFieldTD(sPos, rPos, fs, Q, sigLen, nTrials, c)
+function [rSig, rSigMean, rSigMtx, rSigMtxMean, Q] = ptSrcFieldTd(sPos, rPos, fs, Q, sigLen, nTrials, c)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Utilities/Geometries/MATLAB/Functions/vMicGeo.m b/Utilities/Geometries/MATLAB/Functions/vMicGeo.m
index 3b18c17..723d296 100644
--- a/Utilities/Geometries/MATLAB/Functions/vMicGeo.m
+++ b/Utilities/Geometries/MATLAB/Functions/vMicGeo.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 15/12/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -88,7 +88,7 @@
% setup is one dimensional.
%
% --------------------------------------------------
-function [vPos, vPosMesh] = virtMicGeo(gType, geoDim, nSens, trans, rot, rotOrd)
+function [vPos, vPosMesh] = vMicGeo(gType, geoDim, nSens, trans, rot, rotOrd)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m
index ee9417f..939ea0b 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltEstTd.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 21/12/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -58,7 +58,7 @@
% Notes
%
% --------------------------------------------------
-function [estPerMic, est, err, estMean, errMean] = obsFiltEstTD(m, O, e)
+function [estPerMic, est, err, estMean, errMean] = obsFiltEstTd(m, O, e)
% ====================================================
% Check for number of arguments
% ====================================================
diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m
index cfeafab..17ca29a 100644
--- a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m
+++ b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTd.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 21/12/2024 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -157,7 +157,7 @@
% for local active sound control" by W. Jung, S. J. Elliott and J. Cheer.
%
% --------------------------------------------------
-function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, condNum, mMtx, Omean, RmeMean, RmmMean, Oopt, RmeMtxMean, RmmMtxMean] = obsFiltTD(e, m, beta, snrVal, filtLen, delay, fs)
+function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, condNum, mMtx, Omean, RmeMean, RmmMean, Oopt, RmeMtxMean, RmmMtxMean] = obsFiltTd(e, m, beta, snrVal, filtLen, delay, fs)
% ====================================================
% Check for number of arguments
% ====================================================
--
GitLab
From e0cb078edb9d95b546b238843d5f2356549f7b9f Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:34:30 +0000
Subject: [PATCH 56/60] Update the name of the firstOrderDmaTd function to
match the filename and add a return argument with the condition of the matrix
to be inverted
---
.../MATLAB/Functions/firstOrderDmaTd.m | 29 +++++++++++++------
1 file changed, 20 insertions(+), 9 deletions(-)
diff --git a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
index 97530c4..f08379a 100644
--- a/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
+++ b/Signal Processing/Array Processing/MATLAB/Functions/firstOrderDmaTd.m
@@ -3,7 +3,7 @@
% Author: Achilles Kappis
% e-mail: axilleaz@protonmail.com
%
-% Date: 04/01/2025 (DD/MM/YYYY)
+% Date: 05/01/2025 (DD/MM/YYYY)
%
% Copyright: MIT
% --------------------------------------------------
@@ -55,6 +55,11 @@
% separately. This is an Lx2 matrix, with L the length of the
% filters.
%
+% cNum [numeric]: The condition number of (sMtx.' * sMtx + regFac * I)
+% where "I" is an identity matrix of appropriate
+% dimensions. This is to provide some indication of the
+% conditioning/sensitivity of the filter.
+%
% out [numeric]: The output of each microphone (after filtering) for each
% signal provided in the input argument "sigs". This is a
% KxPx2 array, with K the length of the signals and P the
@@ -71,12 +76,12 @@
% Differential Microphone Arrays" by Buchris, Cohen and Benesty.
%
% --------------------------------------------------
-function [hMtx, h, out, outSum] = firstOrderDmaTD(sMtx, del, regFac, sig)
+function [hMtx, h, cNum, out, outSum] = firstOrderDmaTd(sMtx, del, regFac, sig)
% ====================================================
% Check for number of arguments
% ====================================================
narginchk(1, 4);
- nargoutchk(0, 4);
+ nargoutchk(0, 5);
% ====================================================
@@ -117,12 +122,12 @@ function [hMtx, h, out, outSum] = firstOrderDmaTD(sMtx, del, regFac, sig)
% Calculate filter(s)
% ====================================================
if regFac ~= 0
- hMtx = regFac * eye(size(sMtx, 2));
+ invQty = regFac * eye(size(sMtx, 2));
else
- hMtx = 0;
+ invQty = 0;
end
-
- hMtx = (sMtx.' * sMtx + hMtx)\(sMtx.') * rhsVec;
+ invQty = (sMtx.' * sMtx + invQty);
+ hMtx = invQty\sMtx.' * rhsVec;
% ====================================================
% Return additional output arguments
@@ -132,8 +137,14 @@ function [hMtx, h, out, outSum] = firstOrderDmaTD(sMtx, del, regFac, sig)
h = reshape(hMtx, [], 2);
end
+ if nargout > 2
+ cNum = cond(invQty);
+ else
+ cNum = [];
+ end
+
% Return the filtered signal(s) [output of the DMA]
- if nargout > 2 && ~isempty(sig)
+ if nargout > 3 && ~isempty(sig)
for idx = size(h, 2):-1:1
out(:, :, idx) = filter(h(:, idx), 1, sig(:, :, idx));
end
@@ -141,7 +152,7 @@ function [hMtx, h, out, outSum] = firstOrderDmaTD(sMtx, del, regFac, sig)
out = [];
end
- if nargout > 3 && ~isempty(out)
+ if nargout > 4 && ~isempty(out)
outSum = squeeze(sum(out, 3));
else
outSum = [];
--
GitLab
From 72cec37f48fc96d6145c05fcb7f55c3b994450b3 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:35:18 +0000
Subject: [PATCH 57/60] Move a couple of changes to their appropriate section
and reformat slightly to fix a couple of errors in CHANGELOG.md
---
CHANGELOG.md | 13 +++++++------
1 file changed, 7 insertions(+), 6 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 311bf09..883aed0 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -11,8 +11,8 @@
\* Changed the name of `firstOrderDMA()` to `firstOrderDma()` (cameltoe). This change **breaks** backward compatibility!\
\* Modify `firstOrderDma()` to return the filters as the first argument, in order to avoid calculating the output if it is not required.\
\* Correct error in `firstOrderDma()` header (documentation) and changed name of local array manifold calculation function to clarify that a local function is called.\
-\* Correct erroneous output calculations in `firstOrderDma()`.
-\* Moved the `input` input argument to the last of the input arguments list and updated the input argument checks. This change **breaks backward compatibility**.
+\* Correct erroneous output calculations in `firstOrderDma()`.\
+\* Moved the `input` input argument to the last of the input arguments list in `firstOrderDma()` and updated the input argument checks. This change **breaks backward compatibility**.
**Optimisation - Memetic Algorithm**\
\* Removed redundant code to handle serial for-loop if parfor is not available.
@@ -28,15 +28,16 @@
**Utilities - Geometries**\
\* Fix error in `rcvGeo()` applying incorrect rotation to the elements in the Uniform Circular Array.\
-\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.
+\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.\
+\* Change name of `virtMicGeo()` to `vMicGeo()`. This **breaks** backward compatibility.
** Virtual Sensing - Remote Microphone Technique**\
\* Parallelise the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
\* Parallelise on of the for-loops in `obsFiltTD()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
-\* Change function name from `obsFiltTD()` to `obsFiltTd()`.\
-\* Change function name from `obsFiltEstTD()` to `obsFiltEstTd()`.\
\* Update `obsFilt()` to return the condition of the perturbed power spectral density matrix.
-\* Change name of `virtMicGeo()` to `vMicGeo()`. This **breaks** backward compatibility.
+\* Rename the following functions (this **breaks** backward compatibility):
+ - `obsFiltTD()` to `obsFiltTd()`
+ - `obsFiltEstTD()` to `obsFiltEstTd()`
-------------------------------------------
--
GitLab
From 4aa4dbd169300ad37353eff788549ef3995ca18f Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:37:16 +0000
Subject: [PATCH 58/60] Fix a small format error in CHANGELOG.md
---
CHANGELOG.md | 2 +-
1 file changed, 1 insertion(+), 1 deletion(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 883aed0..8215022 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -31,7 +31,7 @@
\* Fix error in `virtMicGeo()` not returning the z coordinates of the virtual microphones in the mesh. This fix **breaks backward compatibility**.\
\* Change name of `virtMicGeo()` to `vMicGeo()`. This **breaks** backward compatibility.
-** Virtual Sensing - Remote Microphone Technique**\
+**Virtual Sensing - Remote Microphone Technique**\
\* Parallelise the for-loops in `obsFiltEst()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
\* Parallelise on of the for-loops in `obsFiltTD()`. While execution may be slower for a small number of sources and receiver positions, the overall execution time will not increase significantly due to the functions' generally fast performance.\
\* Update `obsFilt()` to return the condition of the perturbed power spectral density matrix.
--
GitLab
From 9bc2b86acd49a82b729f74e4c8035b9a0277c2b4 Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:39:50 +0000
Subject: [PATCH 59/60] Fix format of version titles in CHANGELOG.md
---
CHANGELOG.md | 30 +++++++++++++++---------------
1 file changed, 15 insertions(+), 15 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 8215022..ab7b55f 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,4 +1,4 @@
-###v0.5.0###
+###v0.5.0###\
**Signal Processing - Generic**\
\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
@@ -42,7 +42,7 @@
-------------------------------------------
-###v0.4.1###
+###v0.4.1###\
** Optimisation - Memetic Algorithm**\
\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
\* Update `selParents()` to use the parallel for-loop (`parfor`) whenever it is available.
@@ -50,7 +50,7 @@
-------------------------------------------
-###v0.4.0###
+###v0.4.0###\
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
@@ -61,7 +61,7 @@
-------------------------------------------
-###v0.3.2###
+###v0.3.2###\
**Signal Processing - Generic**\
\+ Add option to use negative delay values with the `winSincFracDel()` function.
@@ -80,14 +80,14 @@
-------------------------------------------
-###v0.3.1###
+###v0.3.1###\
**Virtual Sensing - Remote Microphone Technique**\
\+ Add the option to add "noise" to the monitoring microphone signals in the time domain estimation of the optimal observation filters in `obsFiltTD()` to replicate the feature of `obsFilt()`. This change **breaks** backwards compatibility.
-------------------------------------------
-###v0.3.0###
+###v0.3.0###\
**General**\
\* Convert the project to accept column position and direction vectors as input. The affected functions are:
- `srcGeo()`
@@ -108,7 +108,7 @@
-------------------------------------------
-###v0.2.8###
+###v0.2.8###\
**Utilities - Generic**\
\+ Add the option to pick the order the rotations will be performed in the `rotMat3d()` function. This results is a **backward incompatible** change.
@@ -118,7 +118,7 @@
-------------------------------------------
-###v0.2.7###
+###v0.2.7###\
**Signal Processing - Measurements**\
\+ Introduction of the *Measurements* "topic" in the Signal Processing part of the codebase.\
\+ Function to estimate the impulse responses from sweep measurements.
@@ -126,14 +126,14 @@
-------------------------------------------
-###v0.2.6###
+###v0.2.6###\
**Signal Processing - Generic**\
\+ Add function to perform fractional octave smoothing of spectra.
-------------------------------------------
-###v0.2.5###
+###v0.2.5###\
**Utilities - Geometries**\
\* Combine translation offsets to a single vector input argument for the `rcvGeo()` MATLAB function.\
\* Combine translation offsets to a single vector input argument for th `virtMicGeo()` MATLAB function.\
@@ -142,7 +142,7 @@
-------------------------------------------
-###v0.2.4###
+###v0.2.4###\
**Virtual Sensing**\
\* Corrected time-domain observation filter calculations.\
\* Removed scaling from the cross-correlation calculations in the time-domain optimal observation filters.\
@@ -157,7 +157,7 @@
-------------------------------------------
-###v0.2.3###
+###v0.2.3###\
**Virtual Sensing**\
\+ Added function to estimate the observation filters in the time-domain.\
\+ Added function to perform estimation with observation filters in the time-domain.
@@ -172,21 +172,21 @@
-------------------------------------------
-###v0.2.2###
+###v0.2.2###\
**Virtual Sensing**\
\* Update the observation filter and estimation with observation filter functions with better input argument checking and changed the name of variables to match those in the literature.
-------------------------------------------
-###v0.2.1###
+###v0.2.1###\
**Virtual Sensing**\
\* Fix a bug where noise was added to the power spectral density matrix for the optimal observation filter calculations in the noiseless case.
-------------------------------------------
-###v0.2.0###
+###v0.2.0###\
**Utilities**\
\+ Added a bisection method for single-valued functions.\
\* Improved rotation matrix calculation function with degrees and radian calculations available.\
--
GitLab
From e8c7fc2c9007634df3c40d5e06015be4f5fb7f5d Mon Sep 17 00:00:00 2001
From: ZaellixA <axilleaz@protonmail.com>
Date: Sun, 5 Jan 2025 18:43:10 +0000
Subject: [PATCH 60/60] Fix the titles of the versions in CHANGELOG.md - this
time it is fixed
---
CHANGELOG.md | 32 ++++++++++++++++----------------
1 file changed, 16 insertions(+), 16 deletions(-)
diff --git a/CHANGELOG.md b/CHANGELOG.md
index ab7b55f..70c69dd 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,4 +1,4 @@
-###v0.5.0###\
+### v0.5.0 ###
**Signal Processing - Generic**\
\* Fix bug returning wrong length filter when negative delays are use in `winSincFracDel()`.
@@ -42,7 +42,7 @@
-------------------------------------------
-###v0.4.1###\
+### v0.4.1 ###
** Optimisation - Memetic Algorithm**\
\* Update `mutPop()` to use the parallel for-loop (`parfor`) whenever it is available.\
\* Update `selParents()` to use the parallel for-loop (`parfor`) whenever it is available.
@@ -50,7 +50,7 @@
-------------------------------------------
-###v0.4.0###\
+### v0.4.0 ###
**Utilities - Generic**\
\* Change the interface of `rotMat3d()`, so that the last argument is a boolean (`logical` in MATLAB), declaring whether the given angles are in degrees (if false, the angles are treated as radians). This is a **backwards incompatible** change.\
\* The following functions are changed internally to apply the changes in `rotMat3d()`:
@@ -61,7 +61,7 @@
-------------------------------------------
-###v0.3.2###\
+### v0.3.2 ###
**Signal Processing - Generic**\
\+ Add option to use negative delay values with the `winSincFracDel()` function.
@@ -80,14 +80,14 @@
-------------------------------------------
-###v0.3.1###\
+### v0.3.1 ###
**Virtual Sensing - Remote Microphone Technique**\
\+ Add the option to add "noise" to the monitoring microphone signals in the time domain estimation of the optimal observation filters in `obsFiltTD()` to replicate the feature of `obsFilt()`. This change **breaks** backwards compatibility.
-------------------------------------------
-###v0.3.0###\
+### v0.3.0 ###
**General**\
\* Convert the project to accept column position and direction vectors as input. The affected functions are:
- `srcGeo()`
@@ -108,7 +108,7 @@
-------------------------------------------
-###v0.2.8###\
+### v0.2.8 ###
**Utilities - Generic**\
\+ Add the option to pick the order the rotations will be performed in the `rotMat3d()` function. This results is a **backward incompatible** change.
@@ -118,7 +118,7 @@
-------------------------------------------
-###v0.2.7###\
+### v0.2.7 ###
**Signal Processing - Measurements**\
\+ Introduction of the *Measurements* "topic" in the Signal Processing part of the codebase.\
\+ Function to estimate the impulse responses from sweep measurements.
@@ -126,14 +126,14 @@
-------------------------------------------
-###v0.2.6###\
+### v0.2.6 ###
**Signal Processing - Generic**\
\+ Add function to perform fractional octave smoothing of spectra.
-------------------------------------------
-###v0.2.5###\
+### v0.2.5 ###
**Utilities - Geometries**\
\* Combine translation offsets to a single vector input argument for the `rcvGeo()` MATLAB function.\
\* Combine translation offsets to a single vector input argument for th `virtMicGeo()` MATLAB function.\
@@ -142,7 +142,7 @@
-------------------------------------------
-###v0.2.4###\
+### v0.2.4 ###
**Virtual Sensing**\
\* Corrected time-domain observation filter calculations.\
\* Removed scaling from the cross-correlation calculations in the time-domain optimal observation filters.\
@@ -157,7 +157,7 @@
-------------------------------------------
-###v0.2.3###\
+### v0.2.3 ###
**Virtual Sensing**\
\+ Added function to estimate the observation filters in the time-domain.\
\+ Added function to perform estimation with observation filters in the time-domain.
@@ -172,21 +172,21 @@
-------------------------------------------
-###v0.2.2###\
+### v0.2.2 ###
**Virtual Sensing**\
\* Update the observation filter and estimation with observation filter functions with better input argument checking and changed the name of variables to match those in the literature.
-------------------------------------------
-###v0.2.1###\
+### v0.2.1 ###
**Virtual Sensing**\
\* Fix a bug where noise was added to the power spectral density matrix for the optimal observation filter calculations in the noiseless case.
-------------------------------------------
-###v0.2.0###\
+### v0.2.0 ###
**Utilities**\
\+ Added a bisection method for single-valued functions.\
\* Improved rotation matrix calculation function with degrees and radian calculations available.\
@@ -237,5 +237,5 @@
-------------------------------------------
-###v0.1.0###
+### v0.1.0 ###
This is the initial version of the project.
--
GitLab