// $Id: WaveDWT.cc,v 0.1 2001/08/10 19:10:00 klimenko Exp $

#define WAVEDWT_CC

#include <strstream>
#include <stdexcept>
#include "WaveDWT.hh"

namespace datacondAPIwat {

// constructors

template<class DataType_t>
WaveDWT<DataType_t>::
WaveDWT(int mH=1,int mL=1,int tree=0, enum BORDER border=B_CYCLE) :
   Wavelet(mH, mL, tree, border), pWWS(NULL), nWWS(0)
{ }

template<class DataType_t>
WaveDWT<DataType_t>::WaveDWT(const Wavelet &w) : 
   Wavelet(w), pWWS(NULL), nWWS(0)
{ }

template<class DataType_t>
WaveDWT<DataType_t>::WaveDWT(const WaveDWT<DataType_t> &w) : 
   Wavelet(w), pWWS(NULL), nWWS(0)
{ }

// destructor
template<class DataType_t>
WaveDWT<DataType_t>::~WaveDWT()
{ 
   release(); 
}

template<class DataType_t>
WaveDWT<DataType_t>* WaveDWT<DataType_t>::Clone() const
{
  return new WaveDWT<DataType_t>(*this);
}

template<class DataType_t>
void WaveDWT<DataType_t>::decompose(int level, int layer)
{ }	

// - Virtual procedure for forward wavelet transform. 
// - Real code for decompose will appear in a derivative 
//   classe, since it depends of the wavelet type.
// - Only ldeep steps of decomposition will be done.
// - By default ldeep=1, which means do one step.
template<class DataType_t>
void WaveDWT<DataType_t>::t2w(int ldeep)
{
   int maxLevel = getMaxLevel();

   int levs = m_Level;
   int levf = m_Level+ldeep;
   if((ldeep == -1) || (levf > maxLevel)) levf = maxLevel;
   
   int layf;
   for(int level=levs; level<levf; level++)
   {
      layf = (m_TreeType==1) ? 1<<level : 1;

      for(int layer=0; layer<layf; layer++)
	 decompose(level,layer);
      
      m_Level=level+1;
      
   }
   
   m_Level=levf;
}

template<class DataType_t>
void WaveDWT<DataType_t>::reconstruct(int level, int layer)
{ }	

// - Virtual procedure for inverse wavelet transform. 
// - Real code for reconstruct will appear in a derivative 
//   classe, since it depends of the type of wavelet.
// - Only ldeep steps of reconstruction will be done.
// - By default ldeep=1, which means do one step.
template<class DataType_t>
void WaveDWT<DataType_t>::w2t(int ldeep)
{
   int levs = m_Level;
   int levf = m_Level-ldeep;
   if((ldeep == -1) || (levf < 0)) levf = 0;
   
   int layf;
   for(int level=levs-1; level>=levf; level--)
   {
      layf = (m_TreeType==1) ? 1<<level : 1;

      for(int layer=0; layer<layf; layer++)
	 reconstruct(level,layer);
      
      m_Level=level;
      
   }
   m_Level=levf;
}


// access function
// convert layer index or frequency index to slice
template<class DataType_t>
slice WaveDWT<DataType_t>::getSlice(const int index=0)
{
   int level = m_Level; 
   int layer = abs(index);

   int maxlayer = (BinaryTree()) ? (1<<level)-1 : level;

   if(layer > maxlayer){ 
      layer = maxlayer;     

      std::ostrstream oss;
      oss << "Wavelet::getSlice(): "
	  << "argument 2 is set to " << layer << endl;
  
      std::invalid_argument exception(oss.str());
      oss.freeze(0);
      throw exception;
   }

   if(BinaryTree()){
      if(index>0) layer = convertF2L(level,layer);
   } 
   else {
      
      if(layer) {               // detail layers
	 level -= layer-1; 
	 layer = 1;
      }
      else{                     // approximation layer
	 layer = 0;
      }
   }

   return getSlice(level,layer);
}

// convert (level,layer) to slice
template<class DataType_t>
slice WaveDWT<DataType_t>::getSlice(const int level, const int layer)
{
   if(!allocate()){
      std::invalid_argument("WaveDWT::getSlice(): data is not allocated");
      return slice(0,1,1);
   }      

   size_t m = nWWS>>level;                   // number of elements
   size_t s = 1<<level;                      // slice step
   size_t i = getOffset(level,layer);        // first element

   if(i>=m){
      std::invalid_argument("WaveDWT::getSlice(): invalide arguments");
      return slice(0,1,1);
   }      

   return slice(i,m,s);
}


// Calculate maximal level of wavelet decomposition,
// which depends on the length of wavelet filters.
template<class DataType_t>
int WaveDWT<DataType_t>::getMaxLevel()
{
   int maxLevel = 0;
   if(!allocate()) return 0;
   int n = nWWS;;
   
   for(; (n>m_H) && (n>m_L); n=n>>1) maxLevel++;

   return maxLevel;
}

// allocate input data
template<class DataType_t>
bool WaveDWT<DataType_t>::
allocate(unsigned int n, DataType_t *p)
{
   bool allocate = false;
   if(pWWS == NULL && n>0){
      allocate = true;
      pWWS = p;
      nWWS = n;
   }
   return allocate;
}

// check allocation status
template<class DataType_t>
bool WaveDWT<DataType_t>::allocate()
{
   return (pWWS == NULL || nWWS == 0) ? false : true;
}

// release input data
template<class DataType_t>
void WaveDWT<DataType_t>::release()
{
   pWWS = NULL;
   nWWS = 0;
}

// instantiations

#define CLASS_INSTANTIATION(class_) template class WaveDWT< class_ >;

CLASS_INSTANTIATION(float)
CLASS_INSTANTIATION(double)
//CLASS_INSTANTIATION(std::complex<float>)
//CLASS_INSTANTIATION(std::complex<double>)

#undef CLASS_INSTANTIATION

template WaveDWT<float>::
WaveDWT(const WaveDWT<float> &);
template WaveDWT<double>::
WaveDWT(const WaveDWT<double> &);

template bool WaveDWT<float>::
allocate(unsigned int, float *);
template bool WaveDWT<double>::
allocate(unsigned int, double *);

} // namespace datacondAPIwat