diff --git a/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m b/Virtual Sensing/Remote Microphone Technique/MATLAB/Functions/obsFiltTD.m
index 09a53a42a368f7482eadcaa05b4ac3afa0bad836..117a978d4954f57ebca443d5d753d82c8cd3b763 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: 12/09/2024 (DD/MM/YYYY)
+% Date: 13/09/2024 (DD/MM/YYYY)
 %
 % Copyright: MIT
 % --------------------------------------------------
@@ -37,14 +37,27 @@
 % delay [numeric] (Optional): The delay to be added to the virtual
 %                             microphone signals, effectively implementing
 %                             the delayed Remote Microphone Technique. This
-%                             must be a real non-negative scalar less than
+%                             can be either a scalar, in which case the
+%                             same delay is applied to all virtual
+%                             microphone signals, or a vector with each
+%                             element holding the delay for each virtual
+%                             microphone signal respectively. The delays
+%                             must be real non-negative values less than
 %                             or equal to I, the length of the virtual
-%                             microphone signals. [Default: 0].
+%                             microphone signals. If "fs" is provided,
+%                             these values correspond to seconds otherwise
+%                             they are in samples and must be integer.
+%                             [Default: 0].
+% 
+% fs [numeric] (Optional): The sampling frequency. If this value is given,
+%                          the delay is treated as being in seconds,
+%                          otherwise it is in samples. This argument must
+%                          be a positive real integer.
 % 
 % --------------------------------------------------
 % Output
 % 
-% O [numeric]: The observation filters. This is an firLenxNxMxJ array.
+% O [numeric]: The observation filters. This is an NxfirLenxMxJ array.
 % 
 % Rme [numeric]: The cross-correlation matrices between the monitoring and
 %                virtual microphone signals. This is an filtLenxMxNxJ
@@ -78,10 +91,13 @@
 %                 can be used directly with either "Ovec" or "Oopt"
 %                 arguments to perform the filtering/estimation.
 % 
-% Omean [numeric]: The observation filters averaged over the trials/sound
-%                  field realisations. This is an NxMxfirLen array.
+% Omean [numeric]: The observation filters calculated with the correlation
+%                  matrices RmeMtx and RmmMtx averaged over the
+%                  trials/sound field realisations (these are the
+%                  "RmeMtxMean" and "RmmMtxMean" output arguments). This is
+%                   an NxfirLenxM array.
 % 
-% Rmemean [numeric]: The cross-correlation matrices between the monitoring
+% RmeMean [numeric]: The cross-correlation matrices between the monitoring
 %                    and virtual microphone signals averaged over the
 %                    trials/sound field realisations. This is an
 %                    NxMxfirLen array.
@@ -94,10 +110,11 @@
 % Oopt [numeric]: This is a matrix with the observation filters
 %                 stacked/vectorised, so that they can be applied to a
 %                 stacked/vectorised monitoring microphone measurements
-%                 vector, averaged over the trials/sound field
-%                 realisations. The dimensions of this matrix are
-%                 NxTotFiltLen. These filters are denoted "optimal" since
-%                 they are averaged over many realisations (if J > 1).
+%                 vector, calculated with the correlation matrices RmeMtx
+%                 and RmmMtx averaged over all trials/sound field
+%                 realisations (the "RmeMtxMean" and "RmmMtxMean" output
+%                 arguments). The dimensions of this matrix are
+%                 NxTotFiltLen.
 % 
 % RmeMtxMean [numeric]: This is the matrix with the cross-correlation
 %                       vectors between the monitoring and virtual
@@ -119,11 +136,11 @@
 %   Stephen J. Elliott and Jordan Cheer.
 % 
 % --------------------------------------------------
-function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt, RmeMtxMean, RmmMtxMean] = obsFiltTD(e, m, beta, filtLen, delay)
+function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt, RmeMtxMean, RmmMtxMean] = obsFiltTD(e, m, beta, filtLen, delay, fs)
     % ====================================================
     % Check for number of arguments
     % ====================================================
-    narginchk(2, 5);
+    narginchk(2, 6);
     nargoutchk(0, 13);
 
     % ====================================================
@@ -135,24 +152,50 @@ function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt
 
     % Validate optional arguments
     if nargin > 2 && ~isempty(beta)
-        validateattributes(beta, "numeric", {'scalar', 'nonnegative', 'nonnan', 'nonempty', 'finite'}, mfilename, "Regularisation factor", 3);
+        validateattributes(beta, "numeric", {'scalar', 'nonnegative', 'nonnan', 'nonempty', 'finite', 'real'}, mfilename, "Regularisation factor", 3);
     else
         beta = 0;
     end
 
     if nargin > 3 && ~isempty(filtLen)
-        validateattributes(filtLen, "numeric", {'scalar', 'positive', 'nonnan', 'nonempty', 'finite', '<=', size(m, 1)}, mfilename, "Length of observation filters", 4);
+        validateattributes(filtLen, "numeric", {'scalar', 'positive', 'nonnan', 'real', 'nonempty', 'finite', '<=', size(m, 1)}, mfilename, "Length of observation filters", 4);
     else
         filtLen = size(m, 2);
     end
 
     if nargin > 4 && ~isempty(delay)
-        validateattributes(delay, "numeric", {'scalar', 'nonnegative', 'nonnan', 'finite', 'nonempty', '<=', size(e, 1)}, mfilename, "Delay to be added to virtual microphone signals to implement the delayed Remote Microphone Technique", 5);
+        validateattributes(delay, "numeric", {'vector', 'nonnegative', 'nonnan', 'finite', 'real', 'nonempty', '<=', size(e, 1)}, mfilename, "Delay to be added to virtual microphone signals to implement the delayed Remote Microphone Technique", 5);
+
+        if isscalar(delay)
+            delay = delay * ones(size(e, 2), 1);
+        elseif length(delay) ~= size(e, 2)
+            error("The 'delay' argument must be either a scalar, or a vector with number of elements equal to the number of virtual microphone signals.");
+        end
     else
-        delay = 0;
+        delay = zeros(size(e, 2), 1);
+    end
+
+    if nargin > 5 && ~isempty(fs)
+        validateattributes(fs, "numeric", {'scalar', 'integer', 'positive', 'nonnan', 'nonempty', 'finite', 'real'}, mfilename, "The sampling frequency", 6);
+    else
+        fs = [];
     end
     
 
+    % ====================================================
+    % Pre-process data
+    % ====================================================
+    % Delay the virtual microphone signals
+    if isempty(fs)
+        for jIdx = size(e, 3):-1:1
+            e(:, :, jIdx) = delayseq(e(:, :, jIdx), delay);
+        end
+    else
+        for jIdx = size(e, 3):-1:1
+            e(:, :, jIdx) = delayseq(e(:, :, jIdx), delay, fs);
+        end
+    end
+
     % ====================================================
     % Calculate cross- and auto-correlation matrices
     % ====================================================
@@ -165,7 +208,7 @@ function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt
                 [corr, lags] = xcorr(m(:, mIdx, jIdx), e(:, eIdx, jIdx));
                 lIdx = find(lags == 0);
 
-                Rme(:, mIdx, eIdx, jIdx) = corr(delay + lIdx:-1:lIdx-filtLen + 1);
+                Rme(:, mIdx, eIdx, jIdx) = corr(lIdx:-1:lIdx-filtLen + 1);
             end
 
             % Go through the monitoring microphones to calculate the monitoring microphone correlation matrices
@@ -215,9 +258,17 @@ function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt
         mMtx = reshape(m, prod(size(m, [1, 2])), size(m, 3));
     end
 
-    % Mean observation over trials/sound field realisations
+    % Observation filter calculated with the mean Rme and Rmm over the trials/sound field realisation
     if nargout > 7
-        Omean = mean(O, 4);
+        % Average the RmeMtx and RmmMtx matrices
+        RmeMtxMean = mean(RmeMtx, 3);
+        RmmMtxMean = mean(RmmMtx, 3);
+
+        % Calculate the observation filter
+        Oopt = RmeMtxMean.'/(RmmMtxMean + beta * eye(size(RmmMtxMean)));
+
+        % Reshape
+        Omean = permute(reshape(Oopt.', filtLen, size(m, 2), size(e, 2)), [3, 1, 2]);
     end
     
     % Mean cross-correlations between monitoring and virtual microphones over trials/sound field realisations
@@ -229,19 +280,4 @@ function [O, Rme, Rmm, Ovec, RmeMtx, RmmMtx, mMtx, Omean, RmeMean, RmmMean, Oopt
     if nargout > 9
         RmmMean = mean(Rmm, 5);
     end
-
-    % Average of optimal observation filters over trials/sound field realisations
-    if nargout > 10
-        Oopt = mean(Ovec, 3);
-    end
-
-    % Mean cross-correlation matrix between monitoring and virtual microphones over trials/sound field realisations
-    if nargout > 11
-        RmeMtxMean = mean(RmeMtx, 3);
-    end
-
-    % Mean cross-correlation matrix of monitoring microphones over trials/sound field realisations
-    if nargout > 12
-        RmmMtxMean = mean(RmmMtx, 3);
-    end
 end
\ No newline at end of file