1## Copyright (C) 2006 Peter V. Lanspeary <pvl@mecheng.adelaide.edu.au>
2##
3## This program is free software: you can redistribute it and/or modify
4## it under the terms of the GNU General Public License as published by
5## the Free Software Foundation, either version 3 of the License, or
6## (at your option) any later version.
7##
8## This program is distributed in the hope that it will be useful,
9## but WITHOUT ANY WARRANTY; without even the implied warranty of
10## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11## GNU General Public License for more details.
12##
13## You should have received a copy of the GNU General Public License
14## along with this program; see the file COPYING. If not, see
15## <https://www.gnu.org/licenses/>.
16
17## Usage:
18##   [Pxx,freq] = csd(x,y,Nfft,Fs,window,overlap,range,plot_type,detrend)
19##
20##     Estimate cross power spectrum of data "x" and "y" by the Welch (1967)
21##     periodogram/FFT method.  Compatible with Matlab R11 csd and earlier.
22##     See "help pwelch" for description of arguments, hints and references
23##     --- especially hint (7) for Matlab R11 defaults.
24
25function varargout = csd(varargin)
26
27  ## Check fixed argument
28  if ( nargin<2 )
29    error( 'csd: Need at least 2 args. Use help csd.' );
30  endif
31  nvarargin = length(varargin);
32  ## remove any pwelch RESULT args and add 'cross'
33  for iarg=1:nvarargin
34    arg = varargin{iarg};
35    if ( ~isempty(arg) && ischar(arg) && ( strcmp(arg,'power') || ...
36           strcmp(arg,'cross') || strcmp(arg,'trans') || ...
37           strcmp(arg,'coher') || strcmp(arg,'ypower') ))
38      varargin{iarg} = [];
39    endif
40  endfor
41  varargin{nvarargin+1} = 'cross';
42  ##
43  saved_compatib = pwelch('R11-');
44  if ( nargout==0 )
45    pwelch(varargin{:});
46  elseif ( nargout==1 )
47    Pxx = pwelch(varargin{:});
48    varargout{1} = Pxx;
49  elseif ( nargout>=2 )
50    [Pxx,f] = pwelch(varargin{:});
51    varargout{1} = Pxx;
52    varargout{2} = f;
53  endif
54  pwelch(saved_compatib);
55
56endfunction
57