1######################################################################## 2## 3## Copyright (C) 1994-2021 The Octave Project Developers 4## 5## See the file COPYRIGHT.md in the top-level directory of this 6## distribution or <https://octave.org/copyright/>. 7## 8## This file is part of Octave. 9## 10## Octave is free software: you can redistribute it and/or modify it 11## under the terms of the GNU General Public License as published by 12## the Free Software Foundation, either version 3 of the License, or 13## (at your option) any later version. 14## 15## Octave is distributed in the hope that it will be useful, but 16## WITHOUT ANY WARRANTY; without even the implied warranty of 17## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 18## GNU General Public License for more details. 19## 20## You should have received a copy of the GNU General Public License 21## along with Octave; see the file COPYING. If not, see 22## <https://www.gnu.org/licenses/>. 23## 24######################################################################## 25 26## -*- texinfo -*- 27## @deftypefn {} {} fftconv (@var{x}, @var{y}) 28## @deftypefnx {} {} fftconv (@var{x}, @var{y}, @var{n}) 29## Convolve two vectors using the FFT for computation. 30## 31## @code{c = fftconv (@var{x}, @var{y})} returns a vector of length equal to 32## @code{length (@var{x}) + length (@var{y}) - 1}. If @var{x} and @var{y} 33## are the coefficient vectors of two polynomials, the returned value is the 34## coefficient vector of the product polynomial. 35## 36## The computation uses the FFT by calling the function @code{fftfilt}. If 37## the optional argument @var{n} is specified, an N-point FFT is used. 38## @seealso{deconv, conv, conv2} 39## @end deftypefn 40 41function c = fftconv (x, y, n) 42 43 if (nargin < 2 || nargin > 3) 44 print_usage (); 45 endif 46 47 if (! (isvector (x) && isvector (y))) 48 error ("fftconv: both A and B must be vectors"); 49 endif 50 la = length (x); 51 lb = length (y); 52 if ((la == 1) || (lb == 1)) 53 c = x * y; 54 else 55 lc = la + lb - 1; 56 x(lc) = 0; 57 y(lc) = 0; 58 if (nargin == 2) 59 c = fftfilt (x, y); 60 else 61 if (! isscalar (n)) 62 error ("fftconv: N must be a scalar"); 63 endif 64 c = fftfilt (x, y, n); 65 endif 66 endif 67 68endfunction 69 70 71## FIXME: Borrow tests from conv.m. May need tolerance on the assert stmt. 72%!testif HAVE_FFTW 73%! x = ones (3,1); 74%! y = ones (1,3); 75%! b = 2; 76%! c = 3; 77%! assert (fftconv (x, x), [1; 2; 3; 2; 1], 5*eps); 78%! assert (fftconv (y, y), [1, 2, 3, 2, 1], 5*eps); 79%! assert (fftconv (x, y), [1, 2, 3, 2, 1], 5*eps); 80%! assert (fftconv (y, x), [1; 2; 3; 2; 1], 5*eps); 81%! assert (fftconv (c, x), [3; 3; 3], 5*eps); 82%! assert (fftconv (c, y), [3, 3, 3], 5*eps); 83%! assert (fftconv (x, c), [3; 3; 3], 5*eps); 84%! assert (fftconv (y, c), [3, 3, 3], 5*eps); 85%! assert (fftconv (b, c), 6, 5*eps); 86 87%!test 88%! a = 1:10; 89%! b = 1:3; 90%! assert (size (conv (a,b)), [1, numel(a)+numel(b)-1]); 91%! assert (size (conv (b,a)), [1, numel(a)+numel(b)-1]); 92 93%! a = (1:10).'; 94%! b = 1:3; 95%! assert (size (conv (a,b)), [numel(a)+numel(b)-1, 1]); 96%! assert (size (conv (b,a)), [numel(a)+numel(b)-1, 1]); 97 98%!test 99%! a = 1:10; 100%! b = (1:3).'; 101%! assert (size (conv (a,b)), [1, numel(a)+numel(b)-1]); 102%! assert (size (conv (b,a)), [1, numel(a)+numel(b)-1]); 103 104## Test input validation 105%!error fftconv (1) 106%!error fftconv (1,2,3,4) 107%!error fftconv ([1, 2; 3, 4], 3) 108%!error fftconv (2, []) 109%!error fftconv ([1,1], [2,2] , [3, 4]) 110