Skip to content
Snippets Groups Projects
Commit d2280a3e authored by Achilles Kappis's avatar Achilles Kappis
Browse files

Update winSincFracDel.m to fix errors and bugs and provide four window functions

parent dd6cd670
No related branches found
No related tags found
1 merge request!3Update to v0.2.3
...@@ -12,8 +12,8 @@ ...@@ -12,8 +12,8 @@
% -------------------------------------------------- % --------------------------------------------------
% Input % Input
% %
% delay [numeric]: The delay in samples. This must be a non-negative real % del [numeric]: The delay in samples. This must be a non-negative real
% scalar. % scalar.
% %
% len [numeric] (Optional): This is the length of the filter. This must be % len [numeric] (Optional): This is the length of the filter. This must be
% a positive real integer. If this values is less % a positive real integer. If this values is less
...@@ -21,18 +21,16 @@ ...@@ -21,18 +21,16 @@
% equals the "delay" value rounded up. % equals the "delay" value rounded up.
% [Default: ceil(delay)]. % [Default: ceil(delay)].
% %
% win [numeric] (Optional): A window to apply to the FIR filter. This can % winFun [string/char, numeric] (Optional): Windowing functions to be
% be a vector in which case the same window will % applied to the truncated sinc
% be applied to the filter for all signals, or it % function. The available options
% can be a matrix in which case each column will % are "Rectangular", "Hann" and
% contain the window corresponding to the filter % "Hamming" and are not
% for each signal. The number of rows must be % case-sensitive. Additionally, a
% equal to "len" and the number of columns to the % non-negative scalar can be
% number of columns of "sig". If left empty, no % given to be used as the
% window is applied. It has to be noted that the % parameter for a Kaiser window.
% provided window has to be shifted by "delay" % [Default: "Rectangular"].
% samples.
% [Default: []].
% %
% sig [numeric] (Optional): The signal(s) to be delayed in the time-domain. % sig [numeric] (Optional): The signal(s) to be delayed in the time-domain.
% This can be a matrix with each column being a % This can be a matrix with each column being a
...@@ -52,74 +50,75 @@ ...@@ -52,74 +50,75 @@
% -------------------------------------------------- % --------------------------------------------------
% Output % Output
% %
% delFilt [numeric]: The fractional delay windowed-sinc FIR filter(s). If % sincFilt [numeric]: The fractional delay windowed-sinc FIR filter(s).This
% multiple windows are provided then a matrix is % argument is a vector with the FIR filter. This is the
% returned with each column corresponding to the filter % causal version of the filter. If the non-causal
% for each signal. Otherwise, this argument is a vector % version is required, shift to the left by "causDel".
% with the FIR filter.
% %
% intDel [numeric]: The causal latency of the filter. This value % causDel [numeric]: The delay required to make the filter causal. This
% corresponds to half the length of the filter. % value corresponds to half the length of the filter.
% %
% causalFilt [numeric]: The filter after left shifting "intDel" samples to % dSig [numeric]: The delayed signal(s). The length of each signal is equal
% the left. % to that of the original. If the delay is longer than the
% % original signal vector, a zeroed vector is returned (with
% delSig [numeric]: The delayed signal(s). The length of each signal is % possible artefacts from the filtering process). Only
% equal to that of the original. If the delay is % returned if a signal is provided, otherwise an empty
% longer than the original signal vector, a zeroed vector % array is returned.
% is returned (with possible artefacts from the filtering
% process). Only returned if a signal is provided,
% otherwise an empty array is returned.
%
% causalSig [numeric]: The delayed signal(s) shifted to the left by the
% causal delay of the filter.
% %
% -------------------------------------------------- % --------------------------------------------------
% Notes % Notes
% %
% - The function uses time-domain convolution, which can be quite slow for % - The function uses time-domain convolution, which can be quite slow for
% long signals. If the "delSig" is not asked for this calculation is % long signals. If the "dSig" is not required this calculation is
% omitted so, when the signals are not needed (or want to perform % omitted so, when the signals are not needed (or want to perform
% frequency-domain convolution manually), use only one output argument. % frequency-domain convolution manually), use only one output argument.
%
% -------------------------------------------------- % --------------------------------------------------
function [delFilt, intDel, causalFilt, delSig, causalSig] = winSincFracDel(delay, len, win, sig, sigLen) function [sincFilt, causDel, dSig] = winSincFracDel(del, len, winFun, sig, sigLen)
% ==================================================== % ====================================================
% Check for number of arguments % Check for number of arguments
% ==================================================== % ====================================================
narginchk(1, 5); narginchk(1, 5);
nargoutchk(0, 5); nargoutchk(0, 3);
% ==================================================== % ====================================================
% Validate input arguments % Validate input arguments
% ==================================================== % ====================================================
% Validate mandatory arguments % Validate mandatory arguments
validateattributes(delay, "numeric", {'scalar', 'real', 'nonnegative', 'nonnan', 'nonempty', 'finite'}, mfilename, "Delay in samples", 1); validateattributes(del, "numeric", {'scalar', 'real', 'nonnegative', 'nonnan', 'nonempty', 'finite'}, mfilename, "Delay in samples", 1);
% Validate optional arguments % Validate optional arguments
if nargin > 1 && ~isempty(len) if nargin > 1 && ~isempty(len)
validateattributes(len, "numeric", {'scalar', 'integer', 'real', 'finite', 'nonnan', 'nonnegative'}, mfilename, "Filter length", 2); validateattributes(len, "numeric", {'scalar', 'integer', 'real', 'finite', 'nonnan', 'nonnegative'}, mfilename, "Filter length", 2);
if len < delay if len < del
warning("winSincFracDel(): The length of the filter cannot be less than the delay, setting equal to 'delay',"); warning("winSincFracDel(): The length of the filter cannot be less than the delay, setting equal to 'delay',");
len = ceil(delay); len = ceil(del);
end end
else else
len = ceil(delay); len = ceil(del);
end end
if nargin > 2 && ~isempty(win) if nargin > 2 && ~isempty(winFun)
validateattributes(win, "numeric", {'nonnan', '2d', 'real', 'finite', 'nrows', len}, mfilename, "Window(s) to be applied", 3); validateattributes(winFun, {'char', 'string', 'numeric'}, {'nonempty'}, mfilename, "Window to be applied to the truncated sinc function (or the alpha parameter for a Kaiser window)", 3);
if ischar(winFun) || isStringScalar(winFun)
validatestring(winFun, ["Rectangular", "Hann", "Hamming"], mfilename, "Window to be applied to the truncated sinc function", 3);
winFun = string(winFun);
elseif isscalar(winFun)
validateattributes(winFun, {'numeric'}, {'nonempty', 'nonnan', 'finite', 'nonnegative', 'real', 'scalar'}, mfilename, "Alpha parameter for a Kaiser window to be applied to the truncated sinc function", 3);
end
else else
win = ones(len, 1); winFun = "rectangular";
end end
if nargin > 3 && ~isempty(sig) if nargin > 3 && ~isempty(sig)
validateattributes(sig, "numeric", {'2d', 'nonempty', 'real'}, mfilename, "Signal(s) to be delayed", 4); validateattributes(sig, "numeric", {'2d', 'nonempty', 'real'}, mfilename, "Signal(s) to be delayed", 4);
if size(win, 2) == 1 if size(winFun, 2) == 1
win = repmat(win, 1, size(sig, 2)); winFun = repmat(winFun, 1, size(sig, 2));
else else
if size(win, 2) ~= size(sig, 2) if size(winFun, 2) ~= size(sig, 2)
error("Number of windows must match the number of signals provided"); error("Number of windows must match the number of signals provided");
end end
end end
...@@ -135,46 +134,50 @@ function [delFilt, intDel, causalFilt, delSig, causalSig] = winSincFracDel(delay ...@@ -135,46 +134,50 @@ function [delFilt, intDel, causalFilt, delSig, causalSig] = winSincFracDel(delay
end end
% ==================================================== % ====================================================
% Generate FIR filter % Generate FIR filter for fractional delay
% ==================================================== % ====================================================
delFilt = -len + 1:len; delFilt = delFilt.'; % Filter sample indices idx = (-len + 1:len).'; % Filter sample indices
intDel = -floor(delFilt(1)/2); % The causal delay of the filter causDel = -floor(idx(1)/2); % The causal delay of the filter
delFilt = sinc(delFilt - delay); % Calculate the filter sincFilt = sinc(idx - del); % Calculate the filter
% Keep only the central part of the filter to match length asked for
delFilt = delFilt(intDel + 1:intDel + len);
% Replicate filter and apply window function
delFilt = repmat(delFilt, 1, size(sig, 2)) .* win;
% ==================================================== % ====================================================
% Filter the signals % Apply windowing
% ==================================================== % ====================================================
% Return the shifted filter fracDel = abs(rem(del, 1)); % Total filter delay
if nargout > 2
causalFilt = delFilt(intDel:end);
end
if nargout > 3 % Generate window
for idx = size(sig, 2):-1:1 if isstring(winFun)
delSig(:, idx) = conv(sig, delFilt, "full"); switch lower(winFun)
case "hann"
win = sin(pi * (idx - fracDel)/len).^2;
case "hamming"
alpha = 25/46;
win = alpha - (1 - alpha) * cos(2 * pi * (idx - fracDel)/len);
otherwise % Rectangular
win = ones(length(idx), 1);
end end
else
% Kaiser window
win = besseli(0, pi * winFun * sqrt(1 - (2 * (idx - del)/len).^2))/besseli(0, pi * winFun);
end end
% Shift signal to the left by the amount of the "causal" delay of the filter % Shift window function to centre on sinc's peak
if nargout > 4 if ~isscalar(winFun)
causalSig = delSig(intDel:end, :); win = delayseq(win, floor(del) + causDel);
end end
% Return the length of the signal asked for % Apply window function and keep only the causal part up to the length of the filter
if strcmpi(sigLen, "Same") sincFilt = sincFilt(causDel + 1:causDel + len) .* win(causDel + 1:causDel + len);
if nargout > 3
delSig = delSig(1:size(sig, 1), :);
end % ====================================================
if nargout > 4 % Filter the signals
causalSig = causalSig(1:size(sig, 1), :); % ====================================================
if nargout > 2
for idx = size(sig, 2):-1:1
dSig(:, idx) = conv(sig, sincFilt, lower(sigLen));
end end
end end
end end
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment