function [w0,s0,d0,targets,sigma] = sardecomp(u0,varargin)
%%
% Usage: [w0,s0,d0,targets,sigma] = sardecomp(u0,Name,Value)
%
% Input(s)/Output(s):
%
%   u0      : (input complex matrix) complex pseudo-raw image
%   w0      : (output complex matrix) speckle component of u0
%   s0      : (output complex matrix) recombination of the extracted targets
%             s0=u0-w0 (linear combination of cardinal sines)
%   d0      : (output complex matrix) recombination of the extracted targets
%             into discrete Diracs
%   targets : (output matrix of double) target(i,:) contains the 7 features
%             of the target #(i-1), we have
%
%                        target(i,:) = [x,y,tx,ty,Ar,Ai,log10nfa],
%
%             where,
%
%             + (x,y)   : integer position of the pixel that contains the target
%                         center (top left corner of the image correspond to
%                         coordinate (0,0))
%
%             + (tx,ty) : relative position of the target center in the pixel
%                         (x,y), the pixel that contains the target center,
%                         which means that the sub-pixellic position of the
%                         target is (x-tx,y-ty)
%
%             + (Ar,Ai) : real and imaginary parts of the complex amplitude
%                         A=Ar+1i*Ai of the target
%
%             + log10nfa: log10 of the NFA-score of the target.
%
%   sigma : (scalar double) return the value of the sigma parameter used in
%           the NFA model (sigma = s when optional input s is set to a
%            positive value)
%
% Optional Name-Value pair arguments:
%
%   'epsilon' eps : (scalar, eps=1.) maximal allowed number of false
%                   detection (or number of false alarms, NFA), the target
%                   extraction process stops when no target with NFA
%                   smaller than eps can be found in the image
%   'K' K         : (scalar, default K=25) locality parameter
%   'Ntransl' N   : (scalar, default N=20) number of tested translations
%   'sigma' s     : (scalar, default s=-1) set a positive value for s
%                   to set sigma = s in the NFA model (otherwise the sigma
%                   parameter of the NFA model is learnt)
%   'nmax' n      : (scalar) if used, the extraction process will stop
%                   after at most n extractions
%   'nfaorder' v  : (scalar logical, default v=false) extract targets
%                   in decreasing NFA order (this operation is time
%                   consuming)
%   'verbose' v   : (scalar logical, default v=false) set v=true to display
%                    progression in the Matlab console
%
% Description: Speckle plus targets decomposition algorithm for pseudo-raw
% SAR images.
%
%   This function is an implementation of the 'SPECKLE PLUS TARGETS
%   DECOMPOSITION' algorithm coming from the research paper [1] (see Algorithm 4
%   therein for a pseudo-code description)
%
%   [1] R. Abergel, L. Denis, S. Ladjal and F. Tupin, ``Sub-pixellic Methods for
%       Sidelobes Suppression and Strong Targets Extraction in Single Look
%       Complex SAR Images'', preprint HAL-IMT, submitted in July 2017.
%
%   This algorithm can be used to decompose a single look complex pseudo raw
%   image (u0) into two components: a speckle-dominated image (w0) and the
%   point-like targets (a set of 3-uple of the type (x,y,A), such as (x,y)
%   denotes the sub-pixellic coordinates of the target center, and A its complex
%   amplitude).
%
%   The algorithm can also produce the images s0 and d0. Image s0 is such as
%
%                                u0 = w0 + s0,
%
%   or equivalently, s0 = u0 -w0. Image s0 corresponds to a linear combination of
%   sampled cardinal sine functions (each cardinal sine represents one extracted
%   target). The image d0 corresponds to the linear combination of discrete
%   Diracs (each discrete Dirac correspond to one extracted target whose
%   sub-pixellic position was projected on the initial sampling grid), so that
%   the image R(u0), defined by
%
%                              R(u0) = w0 + d0,
%
%   is free of sidelobe effects. See more details in [1].

%% Note about the research paper [1]
%
%  The work [1] corresponds to the revision of the manuscript entitled
%  ``Sidelobes Suppression and Strong Targets Extraction in Single Look Complex
%  SAR Images with Several Applications'' that we submitted in July 2017 (we
%  changed the title according to one Reviewer's recommandation).

%% License & copyright
%
% Copyright (C) 2017, Rémy Abergel (Remy.Abergel AT parisdescartes.fr), Loïc
% Denis (Loic.Denis AT univ-st-etienne.fr), Saïd Ladjal (Said.Ladjal AT
% telecom-paristech.fr) and Florence Tupin (Florence.Tupin AT
% telecom-paristech.fr).
%
% GNU GPL License, see <http://www.gnu.org/licenses/>.

%% usage displayer
function display_usage()
  fprintf('\n Usage: [w0,s0,d0,targets,sigma] = sardecomp(u0,Name,Value)\n\n');
  fprintf(' Input(s)/Output(s):\n\n');
  fprintf('   u0      : (input complex matrix) complex pseudo-raw image\n');
  fprintf('   w0      : (output complex matrix) speckle component of u0\n');
  fprintf('   s0      : (output complex matrix) recombination of the extracted targets\n');
  fprintf('             s0=u0-w0 (linear combination of cardinal sines)\n');
  fprintf('   d0      : (output complex matrix) recombination of the extracted targets\n');
  fprintf('             into discrete Diracs\n');
  fprintf('   targets : (output matrix of double) target(i,:) contains the 7 features\n');
  fprintf('             of the target #(i-1), we have\n\n');
  fprintf('                        target(i,:) = [x,y,tx,ty,Ar,Ai,log10nfa],\n');
  fprintf('             where,\n\n');
  fprintf('             + (x,y)   : integer position of the pixel that contains the \n');
  fprintf('                         target center (top left corner of the image \n');
  fprintf('                         correspond to coordinate (0,0)\n\n');
  fprintf('             + (tx,ty) : relative position of the target center in the pixel\n');
  fprintf('                         (x,y), the pixel that contains the target center,\n');
  fprintf('                         which means that the sub-pixellic position of the\n');
  fprintf('                         target is (x-tx,y-ty)\n\n');
  fprintf('             + (Ar,Ai) : real and imaginary parts of the complex amplitude\n');
  fprintf('                         A=Ar+1i*Ai of the target\n\n');
  fprintf('             + log10nfa: log10 of the NFA-score of the target.\n');
  fprintf('   sigma : (scalar double) return the value of the sigma parameter used in\n');
  fprintf('           the NFA model (sigma = s when optional input s is set to a\n');
  fprintf('            positive value)\n\n');
  fprintf(' Optional Name-Value pair arguments:\n\n');
  fprintf('   ''epsilon'' eps : (scalar, eps=1.) maximal allowed number of false\n');
  fprintf('                 detection (or number of false alarms, NFA), the target\n');
  fprintf('                 extraction process stops when no target with NFA\n');
  fprintf('                 smaller than eps can be found in the image\n');
  fprintf('   ''K'' K        : (scalar, default K=25) locality parameter\n');
  fprintf('   ''Ntransl'' N  : (scalar, default N=20) number of tested translations\n');
  fprintf('   ''sigma'' s    : (scalar, default s=-1) set a positive value for s\n');
  fprintf('                    to set sigma = s in the NFA model (otherwise the sigma\n');
  fprintf('                   parameter of the NFA model is learnt)\n');
  fprintf('   ''nmax'' n     : (scalar) if used, the extraction process will stop\n');
  fprintf('                  after at most n extractions\n');
  fprintf('   ''nfaorder'' v : (scalar logical, default v=false) extract targets\n');
  fprintf('                  in decreasing NFA order (this operation is time\n');
  fprintf('                  consuming)\n');
  fprintf('   ''verbose'' v  : (scalar logical, default v=false) set v=true to display\n');
  fprintf('                  progression in the Matlab console\n\n');
  fprintf(' Description: Speckle plus targets decomposition algorithm for pseudo-raw\n');
  fprintf(' SAR images. \n\n');
  fprintf('   This function is an implementation of the ''SPECKLE PLUS TARGETS\n');
  fprintf('   DECOMPOSITION'' algorithm coming from the research paper [1] (see Algorithm 4\n');
  fprintf('   therein for a pseudo-code description)\n\n');
  fprintf('   [1] R. Abergel, L. Denis, S. Ladjal and F. Tupin, ``Sidelobes Suppression and\n');
  fprintf('       Strong Targets Extraction in Single Look Complex SAR Images with Several\n');
  fprintf('       Applications'', preprint HAL-IMT, submitted in July 2017.\n\n');
  fprintf('   This algorithm can be used to decompose a single look complex pseudo raw\n');
  fprintf('   image (u0) into two components: a speckle-dominated image (w0) and the\n');
  fprintf('   point-like targets (a set of 3-uple of the type (x,y,A), such as (x,y)\n');
  fprintf('   denotes the sub-pixellic coordinates of the target center, and A its complex\n');
  fprintf('   amplitude).\n\n');
  fprintf('   The algorithm can also produce the images s0 and d0. Image s0 is such as\n\n');
  fprintf('                                u0 = w0 + s0,\n\n');
  fprintf('   or equivalently, s0 = u0 -w0. Image s0 corresponds to a linear combination of\n');
  fprintf('   sampled cardinal sine functions (each cardinal sine represents one extracted\n');
  fprintf('   target). The image d0 corresponds to the linear combination of discrete\n');
  fprintf('   Diracs (each discrete Dirac correspond to one extracted target whose\n');
  fprintf('   sub-pixellic position was projected on the initial sampling grid), so that\n');
  fprintf('   the image R(u0), defined by\n\n');
  fprintf('                              R(u0) = w0 + d0,\n\n');
  fprintf('   is free of sidelobe effects. See more details in [1].\n\n');
end

if(nargin < 1)
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:batinputnumber','Incorrect number of input(s)');
end

%% parser (consistency checks are done after, to allow precise error messages)
p = inputParser;
p.addRequired('u0');
p.addParameter('epsilon',1);
p.addParameter('Ntransl',20);
p.addParameter('K',25);
p.addParameter('nmax',-1);
p.addParameter('nfaorder',false);
p.addParameter('verbose',false);
p.addParameter('sigma',-1);
parse(p,u0,varargin{:});
epsilon = p.Results.epsilon;
Ntransl = p.Results.Ntransl;
K = p.Results.K;
nmax = p.Results.nmax;
nfaorder = p.Results.nfaorder;
verbose = p.Results.verbose;
sigma = p.Results.sigma;
%% consistency checks
% input u0 (matrix of complex double)
if(~ismatrix(u0) || isscalar(u0) || isreal(u0))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input image u0 must be a (non-scalar) matrix of complex double');
end
% input K
if(~isscalar(K) || ~isreal(K))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''K'' must be a scalar real number');
end
if(K~=floor(K)), fprintf('Warning: input ''K'' will be casted to integer\n'); end
% input Ntransl
if(~isscalar(Ntransl) || ~isreal(Ntransl))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''Ntransl'' must be a scalar real number');
end
if(Ntransl~=floor(Ntransl)), fprintf('Warning: input ''Ntransl'' will be casted to integer\n'); end
% input verbose
if(~isscalar(verbose) || ~islogical(verbose))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''verbose'' must be a logical scalar');
end
% input nfaorder
if(~isscalar(nfaorder) || ~islogical(nfaorder))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''nfaorder'' must be a logical scalar');
end
% input nmax
if(~isscalar(nmax) || ~isreal(nmax))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''nmax'' must be a scalar real number');
end
if(nmax~=floor(nmax)), fprintf('Warning: input ''nmax'' will be casted to integer\n'); end
if(nmax < 0 && nmax ~= -1)
    fprintf('Warning: input ''nmax'' has been set to -1 (all meaningful targets will be extracted) due to incorrect input value\n');
    nmax = -1;
end
% input epsilon (scalar double positive number)
if(~isscalar(epsilon) || ~isreal(epsilon))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''epsilon'' must be a scalar positive real number');
end
if(epsilon <= 0)
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''epsilon'' must be a scalar positive real number');
end
% input sigma (scalar double)
if(~isscalar(sigma) || ~isreal(sigma))
  display_usage();
  error('CPSAR_TOOLBOX:sardecomp:badinputdatatype','input ''sigma'' must be a scalar real number');
end

%% run irregularsampling algorithm (mex interface)
filename_targets = tempname();
[w0,s0,d0,sigma] = sardecomp_mexinterface(double(u0),double(epsilon),int32(K),int32(Ntransl),int32(nmax),int32(nfaorder),int32(verbose),filename_targets,double(sigma));
try
    targets = dlmread(filename_targets);
catch % an error avoids when 'filename_targets' is empty (which occurs when no detection is done)
    targets = [];
end
delete(filename_targets);

end
