//------------------------------ PERDECOMP ------------------------------

//       Periodic plus Smooth Image Decomposition
//
//               author: Lionel Moisan
//
//   This program is freely available on the web page
//
//   http://www.mi.parisdescartes.fr/~moisan/p+s
//
//   I hope that you will find it useful.
//   If you use it for a publication, please mention 
//   this web page and the paper it makes reference to.
//   If you modify this program, please indicate in the
//   code that you did so and leave this message.
//   You can also report bugs or suggestions to 
//   lionel.moisan [AT] parisdescartes.fr
//
// This function computes the periodic (p) and smooth (s) components
// of an image (2D array) u
// usage:    p = perdecomp(u)    or    [p,s] = perdecomp(u)
//
// note: this function also works for 1D signals (line or column vectors)
//
// v1.0 (07/2012): initial Scilab version
// v1.1 (10/2012): simplified coordinates, use of meshgrid 
// v1.2 (10/2012): removed slow meshgrid


function [p,s] = perdecomp(u)
  [ny,nx] = size(u); X = 1:nx; Y = 1:ny
  v = zeros(u)
  v(1,X)  = u(1,X)-u(ny,X)
  v(ny,X) = -v(1,X)
  v(Y,1 ) = v(Y,1 )+u(Y,1)-u(Y,nx)
  v(Y,nx) = v(Y,nx)-u(Y,1)+u(Y,nx)
  fx = ones(ny,1)*cos(2.*%pi*(X-1)/nx)
  fy = cos(2.*%pi*(Y'-1)/ny)*ones(1,nx)
  fx(1,1)=0.   // avoid division by 0 in the line below
  s = real(ifft(fft(v)*0.5./(2.-fx-fy)))
  p = u-s
endfunction

