function [wcoefR,  wcoefI] = FWTC_PO(x,L,cofilt)
% FWT_PO -- Forward Discrete Complex Wavelet Transform
%                             (periodized, orthogonal)
%  Usage
%    [wcoefR,  wcoefI] = FWTC_PO(x,L,cofilt)
%  Inputs
%    x    1-d signal; length(x) = 2^J
%    L    Coarsest Level of V_0;  L << J
%    cofilt  Complex Wavelet ON Filter
%  Outputs
%    wcoefR   The real part of 1-d complex wavelet transform of x.
%    wcoefI   The imaginary part of 1-d complex wavelet transform of x.
%
%  Description
%    1. usually, length(qmf) < 2^(L+1)
%    2. To reconstruct use IWTC_PO
%
% Examples of complex filters...
%  con_filt6 = [ ...
% -0.06629126073624 - 0.08558164961018i  0.11048543456040 - 0.08558164961018i ...
%  0.66291260736239 + 0.17116329922036i  0.66291260736239 + 0.17116329922036i ...
%  0.11048543456040 - 0.08558164961018i -0.06629126073624 - 0.08558164961018i];
%
%con_filt10 = [ ...
%  0.01049245051230 + 0.02059043708702i -0.01712890812780 + 0.00872852869034i ...
 %-0.08063970414533 - 0.11794747353812i  0.15137970843150 - 0.09422365674476i ...
 % 0.64300323451588 + 0.18285216450551i  0.64300323451588 + 0.18285216450551i ...
 % 0.15137970843150 - 0.09422365674476i -0.08063970414533 - 0.11794747353812i ...
 %-0.01712890812780 + 0.00872852869034i  0.01049245051230 + 0.02059043708702i];
%
%con_filt14 = [ ...
 %-0.00195611201186 - 0.00486939130521i  0.00323814597408 - 0.00130081479163i ...
 % 0.01765197050496 + 0.03760151845054i -0.03074642007891 + 0.01404497671030i ...
 %-0.08647851802329 - 0.13085968259941i  0.17151218448326 - 0.09445232865968i ...
 % 0.63388553033831 + 0.17983572219509i  0.63388553033831 + 0.17983572219509i ...
 % 0.17151218448326 - 0.09445232865968i -0.08647851802329 - 0.13085968259941i ...
 %-0.03074642007891 + 0.01404497671030i  0.01765197050496 + 0.03760151845054i ...
 % 0.00323814597408 - 0.00130081479163i -0.00195611201186 - 0.00486939130521i];  
 %
 % Also B.G. Sherlock 1998 wrote a procedure MakeCONFilter(N), N=4m+2, m=1,2,3,...
 % that produces Daubechies Complex ON Filter.
%---------------------------------------------------
  qmfR = real(cofilt);
  qmfI = imag(cofilt);
  [n,J] = dyadlength(x) ;
  wcoefR = zeros(1,n) ;
  wcoefI = zeros(1,n) ;
  betaR = ShapeAsRow(x);  %take samples at finest scale as beta-coeffs
  betaI = 0 .* betaR;
  for j=J-1:-1:L
       alphaRR = DownDyadHi(betaR,qmfR);
       alphaRI = DownDyadHi(betaR,-qmfI);
       alphaIR = DownDyadHi(betaI,qmfR);
       alphaII = DownDyadHi(betaI,-qmfI);
       alphaR = alphaRR - alphaII;
       alphaI = alphaIR + alphaRI;
       wcoefR(dyad(j)) = alphaR;
       wcoefI(dyad(j)) = alphaI;
       betaRR = DownDyadLo(betaR,qmfR) ;
       betaRI = DownDyadLo(betaR,qmfI) ;
       betaIR = DownDyadLo(betaI,qmfR) ;
       betaII = DownDyadLo(betaI,qmfI) ;
       betaR = betaRR - betaII;
       betaI = betaIR + betaRI;
  end
  wcoefR(1:(2^L)) = betaR;
  wcoefI(1:(2^L)) = betaI;
  wcoefR = ShapeLike(wcoefR,x);
  wcoefI = ShapeLike(wcoefI,x);

%
% Copyright (c) 2001  Angelini, C., Katul, G., Lina, J-M., and Vidakovic, B.
%