octave_packages/m/signal/stft.m

   1 ## Copyright (C) 1995-2012 Andreas Weingessel
   2 ##
   3 ## This file is part of Octave.
   4 ##
   5 ## Octave is free software; you can redistribute it and/or modify it
   6 ## under the terms of the GNU General Public License as published by
   7 ## the Free Software Foundation; either version 3 of the License, or (at
   8 ## your option) any later version.
   9 ##
  10 ## Octave is distributed in the hope that it will be useful, but
  11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 ## General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with Octave; see the file COPYING.  If not, see
  17 ## <http://www.gnu.org/licenses/>.
  18
  19 ## -*- texinfo -*-
  20 ## @deftypefn {Function File} {[@var{y}, @var{c}] =} stft (@var{x}, @var{win_size}, @var{inc}, @var{num_coef}, @var{win_type})
  21 ## Compute the short-time Fourier transform of the vector @var{x} with
  22 ## @var{num_coef} coefficients by applying a window of @var{win_size} data
  23 ## points and an increment of @var{inc} points.
  24 ##
  25 ## Before computing the Fourier transform, one of the following windows
  26 ## is applied:
  27 ##
  28 ## @table @asis
  29 ## @item @nospell{hanning}
  30 ## win_type = 1
  31 ##
  32 ## @item @nospell{hamming}
  33 ## win_type = 2
  34 ##
  35 ## @item rectangle
  36 ## win_type = 3
  37 ## @end table
  38 ##
  39 ## The window names can be passed as strings or by the @var{win_type} number.
  40 ##
  41 ## If not all arguments are specified, the following defaults are used:
  42 ## @var{win_size} = 80, @var{inc} = 24, @var{num_coef} = 64, and
  43 ## @var{win_type} = 1.
  44 ##
  45 ## @code{@var{y} = stft (@var{x}, @dots{})} returns the absolute values
  46 ## of the Fourier coefficients according to the @var{num_coef} positive
  47 ## frequencies.
  48 ##
  49 ## @code{[@var{y}, @var{c}] = stft (@code{x}, @dots{})} returns the
  50 ## entire STFT-matrix @var{y} and a 3-element vector @var{c} containing
  51 ## the window size, increment, and window type, which is needed by the
  52 ## synthesis function.
  53 ## @end deftypefn
  54
  55 ## Author: AW <Andreas.Weingessel@ci.tuwien.ac.at>
  56 ## Description: Short-Time Fourier Transform
  57
  58 function [y, c] = stft(x, win_size, inc, num_coef, win_type)
  59
  60   ## Default values of unspecified arguments.
  61   if (nargin < 5)
  62     win_type = 1;
  63     if (nargin < 4)
  64       num_coef = 64;
  65       if (nargin < 3)
  66         inc = 24;
  67         if (nargin < 2)
  68           win_size = 80;
  69         endif
  70       endif
  71     endif
  72   elseif (nargin == 5)
  73     if (ischar (win_type))
  74       if (strcmp (win_type, "hanning"))
  75         win_type = 1;
  76       elseif (strcmp (win_type, "hamming"))
  77         win_type = 2;
  78       elseif (strcmp (win_type, "rectangle"))
  79         win_type = 3;
  80       else
  81         error ("stft: unknown window type `%s'", win_type);
  82       endif
  83     endif
  84   else
  85     print_usage ();
  86   endif
  87
  88   ## Check whether X is a vector.
  89   [nr, nc] = size (x);
  90   if (nc != 1)
  91     if (nr == 1)
  92       x = x';
  93       nr = nc;
  94     else
  95       error ("stft: X must be a vector");
  96     endif
  97   endif
  98
  99   ncoef = 2 * num_coef;
 100   if (win_size > ncoef)
 101     win_size = ncoef;
 102     printf ("stft: window size adjusted to %f\n", win_size);
 103   endif
 104   num_win = fix ((nr - win_size) / inc);
 105
 106   ## compute the window coefficients
 107   if (win_type == 3)
 108     ## Rectangular window.
 109     win_coef = ones (win_size, 1);
 110   elseif (win_type == 2)
 111     ## Hamming window.
 112     win_coef = hamming (win_size);
 113   else
 114     ## Hanning window.
 115     win_coef = hanning (win_size);
 116   endif
 117
 118   ## Create a matrix Z whose columns contain the windowed time-slices.
 119   z = zeros (ncoef, num_win + 1);
 120   start = 1;
 121   for i = 0:num_win
 122     z(1:win_size, i+1) = x(start:start+win_size-1) .* win_coef;
 123     start = start + inc;
 124   endfor
 125
 126   y = fft (z);
 127
 128   if (nargout == 1)
 129     y = abs (y(1:num_coef, :));
 130   else
 131     c = [win_size, inc, win_type];
 132   endif
 133
 134 endfunction