//   Class for representing spherical harmonic coefficients
//              G. Le Meur      2000
// DAPNIA/SPP (Saclay) / CEA    LAL - IN2P3/CNRS  (Orsay)

#ifndef ALM_SEEN
#define ALM_SEEN

#include "nbrandom.h"
#include "nbmath.h"
#include "triangmtx.h"
#include "tvector.h"

namespace SOPHYA {  

/*! class for the coefficients \f$a_{lm}\f$ of the development of a function defined on a sphere, in spherical harmonics */
template <class T>
class Alm : public TriangularMatrix<complex<T> >
  {
    public :
 
Alm() : TriangularMatrix<complex<T> >()  {;};
    /* instanciate the class from the maximum value of the index \f$l\f$ in the sequence of the \f$a_{lm}\f$'s */
Alm(int nlmax) : TriangularMatrix<complex<T> >(nlmax+1) {;}
Alm(const Alm<T>& a,  bool share=false)  : TriangularMatrix<complex<T> >(a, share)  {;}

Alm(const TVector<T>& clin, const r_8 fwhm) ;
/*! resize with a new lmax */
inline void ReSizeToLmax(int_4 nlmax) {ReSizeRow(nlmax+1);}
inline int_4 Lmax() const {return rowNumber()-1;}
TVector<T> powerSpectrum() const;

inline  Alm<T>& SetComplexT(complex<T> a)
   {
     this->SetT(a);
     return *this;
   }
 inline Alm<T>& operator = (complex<T> a)
 {
   return SetComplexT(a);
 } 

};
/*! class for a vector with an index running from \f$-m_{max}\f$ to \f$+m_{max}\f$ (then the size of the vector will be actually \f$2m_{max}+1)\f$ */
template <class T>
class Bm
  {
    public :
Bm(): nmmax_(0) {;};
    /* instanciate from the maximum absolute value of the index */
Bm(int mmax) : nmmax_((uint_4)mmax)   { bm_.ReSize( (uint_4)(2*mmax+1) );}
 Bm(const Bm<T>& b,  bool share=false) : nmmax_(b.nmmax_), bm_(b.bm_, share) {;}
 /*! resize with a new value of mmax */
inline void ReSizeToMmax(int_4 mmax) 
  {
    nmmax_=(uint_4)mmax;
    bm_.ReSize(2* nmmax_+1);
  }
//inline int_4 Size() const {return 2*nmmax_+1;};
inline T& operator()(int m) {return bm_( adr_i(m));};
inline T const& operator()(int m) const {return *(bm_.Begin()+ adr_i(m));};
/*! return the current value of the maximum absolute value of the index */
inline int_4 Mmax() const 
   {
     return (int_4)nmmax_;
   }
inline Bm<T>& operator = (const Bm<T>& b)
 {
   if (this != &b)
     {
       nmmax_= b.nmmax_;
       bm_=  b.bm_;
     }
   return *this;
 }

  private:
/*! return the address in the array representing the vector */
inline uint_4 adr_i(int i) const 
{
  return(uint_4)(i+nmmax_);
}



 uint_4 nmmax_;
 NDataBlock<T> bm_;

  };

} // namespace SOPHYA

#endif