sparanfilt.m

阵列信号处理的工具箱

function [sigOut, basisMat, doaRanPosOut] = sparanfilt(sigIn, method, model, modelParam, doaRanChoice, taperCoeff, interCorrMat, steVecFlag, extraParam, orthBeamFlag)

%SPARANFILT Space-range linear filtering (fast-time STAP).
%
%--------
%Synopsis:
%  sigOut = sparanfilt(sigIn)
%
%  [sigOut, basisMat, doaRanPosOut] = sparanfilt(sigIn, method,
%    {antenna, waveform}, {propSpeed, interpMethod, storageClass, ...
%    boundaryMethod}, doaRanChoice, taperCoeff, interCorrMat, steVecFlag, ...
%    extraParam, orthBeamFlag)
%
%  [sigOut, basisMat, doaRanPosOut] = sparanfilt(sigIn, keywordList)
%
%Description:
%  Space-range linear filtering (fast time STAP). 
%
%  After the execution of this function, targets in range bin 1 is located
%  in element number 1 in the output signal matrix. This is different from 
%  the function "pulscomp".
%
%Output and Input:
%  sigOut (RxRadarSigT): Output radar signal.
%  basisMat ():
%  doaRanPosOut (RealVectorT.'): The tested DOAs and ranges [rad, m].
%
%  sigIn (RxRadarSigT): Input radar signal.
%
%  keywordList [D](CellArrayT): A cell array with input parameter name
%    (keyword) and value pairs, see the example below. All input parameters
%    of this function following "sigIn" can be used. The use of keywords
%    is to simplify calls to this function with many input parameters.
%
%  method (StringT): Filter method. See below for sections about each
%    method.
%    = 'conv': Non-adaptive pulse compression in the range domain (convolution).
%    = 'aconv': Adaptive pulse compression in the range domain (convolution). 
%      This is pulse compression with adaptive sidelobe cancelation (ASLC).
%    = 'capon': Capon's method.
%    = 'fft': Discrete fourier transform. Cannot be used with this function.
%    = 'ifft': Inverse discrete fourier transform. Cannot be used with this 
%      function.
%    = 'music : MUSIC.
%  antenna [D](AntDefT): Antenna model to use. If this parameters is an
%    empty matrix, the antenna model in the input signal "sigIn" is used
%    instead.
%  waveform [D](WaveformT): Waveform (time properties) model to use. 
%    If this parameters is an empty matrix, the waveform model in the input
%    signal "sigIn" is used instead.
%  propSpeed [D](RealScalarT): Propagation speed of the waves [m/s].
%    Default = speed of light in free space.
%  interpMethod [D](StringT): Interpolation method to use when calculating
%    the range steering matrices. See function "ranstemat" for possible
%    methods.
%  storageClass [D](StringT): Specifies the MATLAB storage class for range
%    steering matrices. See function "ranstemat" for more information.
%  boundaryMethod [D](StringT): How to handle the boundaries of the PRI.
%    See function "ranstemat" for more information.
%  doaRanChoice [D](RealVectorT.'): One or more DOAs and ranges to test[rad,m].
%  taperCoeff [D](CxVectorT): Tapering coefficients to use. They can be 
%    delivered by the function "gettap1".
%  interCorrMat [D](CxMatrixT): Interference antenna signal correlation matrix. 
%    This is used to adaptively suppress the interference.
%  steVecFlag [D](BoolT): Specifies whether the steering vector (= 1, default) 
%    or filter weights (= 0) will be used in the filtering.
%  extraParam [D](): The use and data type depend on filter method. See below 
%    for sections about each method.
%  orthBeamFlag [D](BoolT): 1 = orthogonalize the weight vectors. 
%    0 = Do not do it (default).
%
%--------
%Conventional Convolution.
%Synopsis:
%
%Output and Input:
%
%Description:
%
%Algoritm:
%
%--------
%Conventional Convolution.
%Synopsis:
%
%Output and Input:
%
%Description:
%
%Algoritm:
%
%--------
%Notations:
%  Data type names are shown in parentheses and they start with a capital
%  letter and end with a capital T. Data type definitions can be found in [1]
%  or by "help dbtdata".
%  [D] = This parameter can be omitted and then a default value is used.
%  When the [D]-input parameter is not the last used in the call, it must be
%  given the value [], i.e. an empty matrix.
%  ... = There can be more parameters. They are explained under respective
%  metod or choice.
%
%Examples:
%
%Software Quality:
%  (About what is done to ascertain software quality. What tests are done.)
%  This function is not tested yet.
%
%Known Bugs:
%  *) With methods 'fft' and 'ifft', it is not possible to choose a different 
%  number of output channels than the number of input channels with a
%  correct result.
%  *) This function do not use the function "calcfiltsig".
%
%References:
%  [1]: Bj鰎klund S.: "DBT, A MATLAB Toolbox for Radar Signal Processing.
%    Reference Guide", FOA-D--9x-00xxx-408--SE, To be published.
%
%See Also:
%  sparanstemat, ranfilt, spafilt, pulfilt

%   *  DBT, A Matlab Toolbox for Radar Signal Processing  *
% (c) FOA 1994-2000. See the file dbtright.m for copyright notice.
%
%  Start        : 990301 Svante Bj鰎klund (svabj).
%  Latest change: $Date: 2000/10/16 15:22:03 $ $Author: svabj $.
%  $Revision: 1.6 $
% *****************************************************************************
%Maintenance information:
%  The comment "NFF" stands for "New filter function" and gives tips to
%  necessary changes when using this function as a template for new filter 
%  functions (in other dimensions).
%
%  The functions "pulfilt" and "spapulfilt" contains the latest (and best)
%  code for these filtering function. Look at "pulfilt" and "spapulfilt"
%  for more information.

dbtinfo('sparanfilt: Begin of function.')

% ----------------------------------------------------------------------- %
% Handle input parameters
% ----------------------------------------------------------------------- %

arginNo=1;
if (nargin < arginNo)
  error('DBT-Error: To few input parameters.')
end
arginNo = arginNo +1;

% ****************** Add missing input parameters ******************

arginNo=2;
if (nargin < arginNo)
  method = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  model = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  modelParam = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  doaRanChoice = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  taperCoeff = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  interCorrMat = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  steVecFlag = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  extraParam = [];
end
arginNo = arginNo +1;
if (nargin < arginNo)
  orthBeamFlag = [];
end
arginNo = arginNo +1;


% *************** Evalute the list of keywords  ***************
% If the input parameter "method" is a cell array, it is interpreted as a
% list of input parameter name (keyword) and value pairs.
% The values of these keywords are assigned to the local variables with
% the same name as the keywords. The use of keywords is to make calls to
% this function simpler.

if (iscell(method))
  keywordList = method;
  lenExtra = length(keywordList);
  if (mod(lenExtra,2) ~= 0), error('DBT-Error: Keyword value missing.'),end%if
  for n = 1:2:length(keywordList)
    eval([keywordList{n} '= keywordList{n+1};' ]);
  end%for n
end
arginNo = arginNo +1;


% ****************** Default values ******************

if isempty(method)
  method = 'aconv';
end%if
if isempty(model)
  model = {sigIn.antenna, sigIn.waveform};
end%if

antenna = model{1};
waveform = model{2};
noChannels = antenna.noElem * waveform.noRangeBins;
  % Number of input channels or samples in the dimension we are filtering.
  % NFF: Choose the size of the right dimension.

tempModelParam = {[] [] [] []};
  % NFF: Change to suitable default values.
if isempty(modelParam)
  modelParam = tempModelParam;
else
  tempModelParam(1:length(modelParam)) = modelParam;
  modelParam = tempModelParam;
end%if

propSpeed = modelParam{1};
if isempty(propSpeed)
  propSpeed = speedoflight;
end%if
interpMethod = modelParam{2};
if isempty(interpMethod)
  interpMethod = [];
end%if
storageClass = modelParam{3};
if isempty(storageClass)
  storageClass = [];
end%if
boundaryMethod= modelParam{4};
if isempty(boundaryMethod)
  boundaryMethod = [];
end%if
% NFF: Change to a suitable default values above.

if isempty(doaRanChoice)
  %beamChoice = d2r(-80:2:80);
  beamChoice = linspace(d2r(-90),d2r(90), antenna.noElem+1);
  beamChoice = beamChoice(:,1:(size(beamChoice,2)-1));
  %beamChoice = beamChoice(:,1:4)
  beamChoice = [beamChoice; zeros(1,length(beamChoice))];

  ranBinLen = getranbinlen(waveform, propSpeed);
  noRangeBinToUse = waveform.noRangeBins;
  %noRangeBinToUse = 6;
  rangeChoice = (0.5 : (noRangeBinToUse - 0.5)) * ranBinLen;
    % The center range in meter in each range bin.

  doaRanChoice = {beamChoice,rangeChoice};
    % NFF: Change to a suitable default value for "doaRanChoice".
end%if
if isempty(taperCoeff)
  taperCoeff = [];
end%if
if isempty(interCorrMat)
  interCorrMat = [];
end%if
if isempty(steVecFlag)
  steVecFlag = 1;
end%if
if isempty(extraParam)
  extraParam = [];
end%if
if isempty(orthBeamFlag)
  orthBeamFlag = 0;
end%if


% ************ Add missing input parameters in "xxx" ************


% ****************** Default values in "xxx" ******************


% ************* Pick out some more fields from input parameters. *************

%lambda = waveform.wavelength(1);

if (length(waveform.wavelength) > 1)
  error('DBT-Error: Can only handle one wavelength.')
end%if
lambda = waveform.wavelength;


% ****************** Error check input parameters ******************

chkdtype(sigIn, 'RxRadarSigT')
chkdtype(method, 'StringT')
chkdtype(model, 'CellArrayT')
%chkdtype(modelParam, 'CellArrayT')
chkdtype(taperCoeff, 'CxVectorT')
chkdtype(interCorrMat, 'CxMatrixT', 'EmptyT')
  % NFF: Check any special input parameters.


% ----------------------------------------------------------------------- %
% Preparations before calculations.
% ----------------------------------------------------------------------- %


% ****************** Get signals sizes. ******************
sigSizeIn = sigsize(sigIn);

noTrialsIn  = sigSizeIn(6);
noCpisIn    = sigSizeIn(5);
noExtrasIn  = sigSizeIn(4);
noBeamsIn   = sigSizeIn(3);
noRangesIn  = sigSizeIn(2);
noPulsesIn  = sigSizeIn(1);


% *************** Adjustments for different methods. ***************

if (strcmp(method,'fft') | strcmp(method,'ifft'))
  % Special case when the number of test ranges are not determined by 
  % "doaRanChoice".
  dbterror('FFT or IFFT can not be used with this function.')
    % One reason for the this is that the Fourier tranform uses a unsuitable 
    % basis in the range dimension. The basis should consist of time
    % delayed (flipped back-forth and complex conjugated) versions of
    % the transmitted waveform.
    % A second reason is that the functions "calcfiltinit" and 
    % "calcfiltproc" cannot handle 2D FFT or IFFT.
  if isempty(extraParam)
    extraParam = noBeamsIn * noRangesIn;
      % NFF: Choose the size of the right dimension, i.e. the dimension,
      % in which to filter.
  end%if
  extraParam = extraParam;
  noChanOut = extraParam;
    % Number of ouput channels or samples in the dimension we are filtering.
  steMatIn = [];	% We will not use this parameter.
  interCorrMat = [];	% We will not use this parameter.
else
  dbtinfo('sparanfilt: Before call to sparanstemat.')
  steMatIn = sparanstemat({antenna, waveform}, doaRanChoice, {lambda});
  dbtinfo('sparanfilt: After call to sparanstemat.')
    % NFF: Call the right function with right parameters.
  noChanOut = size(doaRanChoice{1},2) * size(doaRanChoice{1},2);
    % Shall it be doaRanChoice{1} twice or is it an error?
    % Shall it be the number of columns?
    % NFF: Choose the size of the right dimension.
end%if


% ****************** Calculate filtering weights. ******************

sizeSigIn = [noBeamsIn * noRangesIn, noPulsesIn];
  % NFF: Choose the sizes of the right dimensions.
[methodCalc, weightMatrix, basisMat, freqPos] = calcfiltinit(method, ...
   sizeSigIn, steMatIn, taperCoeff, interCorrMat, orthBeamFlag, extraParam);
%basisMat = weightMatrix;


% ****************** Create output variable. ******************

sigOut = sigIn;		% Copy all fields to output signal.
sigOut.signals = zeros(noPulsesIn, noRangesIn, noBeamsIn, noExtrasIn, ...
  noCpisIn, noTrialsIn);
  % NFF: Move "noChanOut" to the right subscript.
  % How to split "noChanOut" into two subscripts in 2D-filtering?


% ----------------------------------------------------------------------- %
% Loops for several range gates and trials.
% ----------------------------------------------------------------------- %

% ****************** Loops. ******************

for trialLoop = 1:noTrialsIn
for cpiLoop = 1:noCpisIn
  for extraLoop = 1:noExtrasIn
  for pulseLoop = 1:noPulsesIn
    % NFF: Loop over suitable dimensions => do not loop
    % over two suitable dimensions.
    %sigMatIn = sigIn.signals(pulseLoop,:,:,extraLoop,cpiLoop, trialLoop);
      % Is this faster? Is this correct?
    sigMatIn = getm3(sigIn.signals, 3, [], pulseLoop,':',':', ...
      extraLoop, cpiLoop, trialLoop);
    size(sigMatIn)
    sigMatIn = sigMatIn(:);
      % NFF: Get the right part of the ND-matrix. Both the
      % "row" input parameters and the ':' input parameters must be correct.
    sigMatOut = calcfiltproc(sigMatIn, methodCalc, weightMatrix, ...
        taperCoeff, orthBeamFlag, extraParam);
disp('Reshaping:')
size(sigMatOut)
noBeamsIn
noRangesIn
    sigMatOut = reshape(sigMatOut, noBeamsIn, noRangesIn);
    sigOut.signals(pulseLoop,':',':',extraLoop,cpiLoop,trialLoop) ...
      = sigMatOut.';
      % NFF: Assign the right part of the ND-matrix.
    infoStr = sprintf('sparanfilt: cpi = %d(%d), pulse = %d(%d)\r',...
      cpiLoop, noCpisIn, pulseLoop, noPulsesIn);
      % NFF: Print the right loop counter.
    dbtinfo(infoStr,1);
  end%for pulseLoop
  end%for extraLoop
end%for cpiLoop
dbtinfo('');	%New line.
end%for trialLoop


% ----------------------------------------------------------------------- %
% Output variables.
% ----------------------------------------------------------------------- %

sigOut = setsparantrans(sigOut, basisMat);
  % NFF: Call the right function.

% For Capon, only the non-adapted steering matrix is stored. See the function
% "calcfiltinit".

doaRanPosOut = doaRanChoice;
beamPosOut  = doaRanChoice{1};
rangePosOut = doaRanChoice{2};
sigOut.doaPos = beamPosOut;
sigOut.ranPos = rangePosOut;
  % NFF: Assign the right field in the output function.

dbtinfo('sparanfilt: End of function.')