Added extra convolution function

This function uses the same method as the convolve function from LTFAT

SFS_binaural_synthesis/auralize_ir.m

@@ -128,14 +128,8 @@
 if size(ir,2)>2
  warning('Your impulse response has more than two channels.');
 end
-% Use convolve function from LTFAT if possible
-if exist('convolve','file')
- outsig = convolve(ir,content);
-else
- for ii = 1:size(ir,2)
- outsig(:,ii) = conv(ir(:,ii),content);
- end
-end
+% Convolve the impulse responses with the content signal
+outsig = convolution(ir,content);
 % Scale output
 if(usenorm)
  outsig = norm_signal(outsig);

SFS_binaural_synthesis/compensate_headphone.m

@@ -80,14 +80,8 @@
  'compensation filter.']);
  end
  % Apply filter
- % The following is the original code from Sascha, but it will work only if
- % the length of the IR is sufficient greater than the length of the
- % headphone compensation filter.
- %ir(:,1) = conv(hcomp(:,1),ir(1:end-length(hcomp)+1,1));
- %ir(:,2) = conv(hcomp(:,2),ir(1:end-length(hcomp)+1,2));
- % Therefore we use this one
- tmp1 = conv(hcomp(:,1),ir(:,1));
- tmp2 = conv(hcomp(:,2),ir(:,2));
+ tmp1 = convolution(hcomp(:,1),ir(:,1));
+ tmp2 = convolution(hcomp(:,2),ir(:,2));
  len = length(ir(:,1));
  ir(:,1) = tmp1(1:len);
  ir(:,2) = tmp2(1:len);

SFS_binaural_synthesis/ir_generic.m

@@ -102,10 +102,8 @@
  % === Sum up virtual loudspeakers/HRIRs and add loudspeaker time delay ===
  % Also applying the weights of the secondary sources including integration
  % weights or tapering windows etc.
- ir_generic(:,1) = ir_generic(:,1) + ...
- fix_ir_length(conv(ir(:,1),d(:,ii)),N) .* g .* x0(ii,7);
- ir_generic(:,2) = ir_generic(:,2) + ...
- fix_ir_length(conv(ir(:,2),d(:,ii)),N) .* g .* x0(ii,7);
+ ir_generic = ir_generic + ...
+ fix_ir_length(convolution(ir,d(:,ii)),N) .* g .* x0(ii,7);
 
 end
 warning('on','SFS:irs_intpol');

SFS_general/bandpass.m

@@ -77,6 +77,6 @@
 Hm = [0 0 1 1 0 0];
 b = fir2(N,Hf,Hm);
 % filter signal
-sig = conv(sig,b);
+sig = convolution(sig,b);
 % compensate for delay & truncate result
 sig = sig(N/2:end-(N/2)-1);

SFS_general/convolution.m

@@ -0,0 +1,81 @@
+function z = convolution(x,y)
+%CONVOLUTION convolve the signals x and y
+%
+% Usage: z = convolution(x,y)
+%
+% Input parameters:
+% x - matrix/vector with signals as columns
+% y - matrix/vector with signals as columns, note that only one of
+% the signals can be a matrix
+%
+% Output parameters:
+% z - convolved signal
+%
+% CONVOLUTION(x,y) convolves the signals given with x and y. One of the input
+% signals can be a matrix containing the signals as column vectors, the other
+% one has to be a column vector. The convolution is done in the frequency
+% domain and it is checked if we have only real signals to speed up the
+% calculation. The length of z is length(x)+length(y)-1.
+%
+% see also: fft_real, ifft_real, fft, ifft
+
+%*****************************************************************************
+% Copyright (c) 2010-2013 Quality & Usability Lab, together with *
+% Assessment of IP-based Applications *
+% Deutsche Telekom Laboratories, TU Berlin *
+% Ernst-Reuter-Platz 7, 10587 Berlin, Germany *
+% *
+% Copyright (c) 2013 Institut fuer Nachrichtentechnik *
+% Universitaet Rostock *
+% Richard-Wagner-Strasse 31, 18119 Rostock *
+% *
+% This file is part of the Sound Field Synthesis-Toolbox (SFS). *
+% *
+% The SFS is free software: you can redistribute it and/or modify it under *
+% the terms of the GNU General Public License as published by the Free *
+% Software Foundation, either version 3 of the License, or (at your option) *
+% any later version. *
+% *
+% The SFS is distributed in the hope that it will be useful, but WITHOUT ANY *
+% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS *
+% FOR A PARTICULAR PURPOSE. *
+% See the GNU General Public License for more details. *
+% *
+% You should have received a copy of the GNU General Public License along *
+% with this program. If not, see <http:https://www.gnu.org/licenses/>. *
+% *
+% The SFS is a toolbox for Matlab/Octave to simulate and investigate sound *
+% field synthesis methods like wave field synthesis or higher order *
+% ambisonics. *
+% *
+% http:https://dev.qu.tu-berlin.de/projects/sfs-toolbox [email protected] *
+%*****************************************************************************
+
+
+%% ===== Checking input parameters =======================================
+nargmin = 2;
+nargmax = 2;
+narginchk(nargmin,nargmax);
+isargmatrix(x,y);
+% check if only one of the inputs is a matrix
+if all(size(x)>1) && all(size(y)>1)
+ error('%s: Only one of the inputs can be multi-dimensional.', ...
+ upper(mfilename));
+end
+
+
+%% ===== Computation =====================================================
+% if one of the input signals is a matrix repmat the vector of the other signal
+if all(size(x)>1)
+ y = repmat(y,1,size(x,2));
+elseif all(size(y)>1)
+ x = repmat(x,1,size(y,2));
+end
+% length of output signal
+N = size(x,1)+size(y,1)-1;
+% convolve the signals in frequency domain
+if isreal(x) && isreal(y)
+ z = ifft_real(fft_real(postpad(x,N)).*fft_real(postpad(y,N)),N);
+else
+ z = ifft(fft(postpad(x,N)).*fft(postpad(y,N)));
+end

SFS_general/delayline.m

@@ -91,7 +91,7 @@
  sig = delayline(sig,idt,weight,conf2);
  if(abs(dt-idt)>0)
  h = hgls2(Lls,dt-idt,0.90);
- sig = conv(sig,h);
+ sig = convolution(sig,h);
  sig = sig(Lls/2:end-Lls/2);
  end
 

SFS_general/fft_real.m

@@ -0,0 +1,63 @@
+function y = fft_real(x)
+%FFT_REAL computes an fft for real value signals
+%
+% Usage: y = fft_real(x)
+%
+% Input parameters:
+% x - matrix with signals as columns
+%
+% Output parameters:
+% y - matrix with signals as columns
+%
+% FFT_REAL(x) computes the fft of the real signals in x. The signals have
+% to be the columns of x.
+%
+% see also: ifft_real, convolution
+
+%*****************************************************************************
+% Copyright (c) 2010-2013 Quality & Usability Lab, together with *
+% Assessment of IP-based Applications *
+% Deutsche Telekom Laboratories, TU Berlin *
+% Ernst-Reuter-Platz 7, 10587 Berlin, Germany *
+% *
+% Copyright (c) 2013 Institut fuer Nachrichtentechnik *
+% Universitaet Rostock *
+% Richard-Wagner-Strasse 31, 18119 Rostock *
+% *
+% This file is part of the Sound Field Synthesis-Toolbox (SFS). *
+% *
+% The SFS is free software: you can redistribute it and/or modify it under *
+% the terms of the GNU General Public License as published by the Free *
+% Software Foundation, either version 3 of the License, or (at your option) *
+% any later version. *
+% *
+% The SFS is distributed in the hope that it will be useful, but WITHOUT ANY *
+% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS *
+% FOR A PARTICULAR PURPOSE. *
+% See the GNU General Public License for more details. *
+% *
+% You should have received a copy of the GNU General Public License along *
+% with this program. If not, see <http:https://www.gnu.org/licenses/>. *
+% *
+% The SFS is a toolbox for Matlab/Octave to simulate and investigate sound *
+% field synthesis methods like wave field synthesis or higher order *
+% ambisonics. *
+% *
+% http:https://dev.qu.tu-berlin.de/projects/sfs-toolbox [email protected] *
+%*****************************************************************************
+
+
+%% ===== Checking of input parameters ====================================
+nargmin = 1;
+nargmax = 1;
+narginchk(nargmin,nargmax);
+isargmatrix(x);
+
+
+%% ===== Computation =====================================================
+N = size(x,1);
+N2=floor(N/2)+1;
+
+% Force IFFT along dimension 1
+y = fft(x,N,1);
+y = y(1:N2,:);

SFS_general/ifft_real.m

@@ -0,0 +1,71 @@
+function y = ifft_real(x,N)
+%IFFT_REAL computes an ifft for real value signals
+%
+% Usage: y = ifft_real(x,samples)
+%
+% Input parameters:
+% x - matrix x with signals as columns
+% samples - desired length of the output signals
+%
+% Output parameters:
+% y - matrix y with signals as columns
+%
+% IFFT_REAL(x,samples) computes the ifft of the real signals in x. The signals
+% have to be the columns of x.
+%
+% see also: fft_real, convolution
+
+%*****************************************************************************
+% Copyright (c) 2010-2013 Quality & Usability Lab, together with *
+% Assessment of IP-based Applications *
+% Deutsche Telekom Laboratories, TU Berlin *
+% Ernst-Reuter-Platz 7, 10587 Berlin, Germany *
+% *
+% Copyright (c) 2013 Institut fuer Nachrichtentechnik *
+% Universitaet Rostock *
+% Richard-Wagner-Strasse 31, 18119 Rostock *
+% *
+% This file is part of the Sound Field Synthesis-Toolbox (SFS). *
+% *
+% The SFS is free software: you can redistribute it and/or modify it under *
+% the terms of the GNU General Public License as published by the Free *
+% Software Foundation, either version 3 of the License, or (at your option) *
+% any later version. *
+% *
+% The SFS is distributed in the hope that it will be useful, but WITHOUT ANY *
+% WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS *
+% FOR A PARTICULAR PURPOSE. *
+% See the GNU General Public License for more details. *
+% *
+% You should have received a copy of the GNU General Public License along *
+% with this program. If not, see <http:https://www.gnu.org/licenses/>. *
+% *
+% The SFS is a toolbox for Matlab/Octave to simulate and investigate sound *
+% field synthesis methods like wave field synthesis or higher order *
+% ambisonics. *
+% *
+% http:https://dev.qu.tu-berlin.de/projects/sfs-toolbox [email protected] *
+%*****************************************************************************
+
+
+%% ===== Checking of input parameters ====================================
+nargmin = 2;
+nargmax = 2;
+narginchk(nargmin,nargmax);
+isargmatrix(x);
+isargpositivescalar(N);
+
+
+%% ===== Computation =====================================================
+N2=floor(N/2)+1;
+
+% Force IFFT along dimension 1
+if rem(N,2)==0
+ y= [x; flipud(conj(x(2:end-1,:)))];
+else
+ y = [x; flipud(conj(x(2:end,:)))];
+end;
+
+y = real(ifft(y,N,1));
+
+

SFS_time_domain/wfs_preequalization.m

@@ -75,10 +75,7 @@
  % get FIR filter
  hpre = wfs_fir_prefilter(conf);
  % apply filter
- for ii = 1:size(ir,2)
- ir_tmp(:,ii) = conv(hpre,ir(:,ii));
- end
- ir = ir_tmp;
+ ir = convolution(hpre,ir);
 elseif strcmp('IIR',hpretype)
  % get IIR filter
  hpre = wfs_iir_prefilter(conf);