#ifndef LAMBUILDER_H_SEEN
#define LAMBUILDER_H_SEEN

#include <math.h>


 class LambdaBuilder{
 public:
  LambdaBuilder();
  LambdaBuilder(double theta,int lmax, int mmax);
  LambdaBuilder(const LambdaBuilder& lb);
  void reset();   
  double lamlm() const {return lamlm_;}
  double lamlm(int flip);
  double lamlm(int l, int m);
  double lamlm(int l, int m, int flip);
  int l() const {return l_;}
  int m() const {return m_;}
  int lmax() const {return lmax_;}
  int mmax() const {return mmax_;}
  void step();
 protected:
  double par_lm_;
  int lmax_;
  int mmax_;
  int l_;
  int m_;
  double lamlm_;
  double cth_;
  double sth_;
  double lammm_;
  double lam_0_;
  double lam_1_;
  double a_rec_;
  double b_rec_;
  double fm2_;
  //  static const double bignorm=1e-20*HUGE_VAL;
  static const double bignorm=1e-20*1.d288;
};

class Lambda2Builder : public LambdaBuilder {
 public:
  Lambda2Builder();
  Lambda2Builder(double theta,int lmax, int mmax);
  Lambda2Builder(const Lambda2Builder& lb);
  void reset();
  double lam2lmp() const {return -fact_*(w_-ix_);}
  double lam2lmp(int flip);
  double lam2lmp(int l, int m);
  double lam2lmp(int l, int m, int flip);
  double lam2lmm() const {return -fact_*(w_+ix_);}
  double lam2lmm(int flip);
  double lam2lmm(int l, int m);
  double lam2lmm(int l, int m, int flip);
  void step();
 private:
  double fact_;
  double w_;
  double ix_;
};
#endif

/*! \class LambdaBuilder
 \brief A lambda building class

Objects in this class generate
\f[
            \lambda_l^m=\sqrt{\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}}
	    P_l^m(\cos{\theta})
\f]
where \f$P_l^m\f$ are the associated Legendre polynomials.  The interest
of course is that
\f[
            Y_{lm}(\cos{\theta})=\lambda_l^m(\cos{\theta}) e^{im\phi}.
\f]

Each object has a fixed theta, and maximum l and m to be calculated
(lmax and mmax).

The algorithm solves the \f$\lambda_l^m\f$'s iteratively so there
is an optimal sequance of l and m's as follows...

\code
for (int m = 0; m <= nmmax; m++){
  dum=lb.lamlm(m,m) 
  for (int l = m+1; l<= nlmax; l++){
    dum = lb.lamlm(l,m)
  }
}\endcode
The function LambdaBuilder::step() automatically shifts the l and
m's to the next optimal value.
However, \f$\lambda_l^m\f$ can be calculated for any l and m less than
or equal to lmax and lmax by going through the loop or restarting the
loop.

\f$\lambda_l^m(-\cos{\theta})\f$ can be quickly obtained using
symmetry relations, which can be more efficient that creating a new
LambdaBuilder object.

The algorithm is based on the HEALPIX code.
*/

/*!  \fn LambdaBuilder::LambdaBuilder(double theta,int lmax, int mmax)
  \param theta the desired \f$\theta\f$
  \param lmax the maximum value of l
  \param mmax the maximum value of m
*/
/*!
  \fn LambdaBuilder::LambdaBuilder(const LambdaBuilder& lb)
  \param lb the LambdaBuilder object to be copied
*/
/*! \fn void LambdaBuilder::reset()
  \brief Resets the object to l=0, m=0
*/
/*! \fn double LambdaBuilder::lamlm() const
  \returns \f$\lambda_l^m\f$ for current l and m
*/
/*!
  \fn double LambdaBuilder::lamlm(int flip)
  \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the current l and m.
  Uses symmetry relations to get this value.
*/
/*! \fn double LambdaBuilder::lamlm(int l, int m)
  \param l the desired l
  \param m the desired m
  \returns \f$\lambda_l^m\f$ for l and m
*/
/*!   \fn double LambdaBuilder::lamlm(int l, int m, int flip)
  \param l the desired l
  \param m the desired m
   \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the given l and m.
  Uses symmetry relations to get this value.
*/
/*!   \fn int LambdaBuilder::l() const
  \returns current value of l
*/
/*!   \fn int LambdaBuilder::m() const
  \returns current value of m
*/
/*!   \fn int LambdaBuilder::lmax() const
  \returns the maximum value of l
*/
/*!   \fn int LambdaBuilder::mmax() const
  \returns the maximum value of m
*/
/*!   \fn void LambdaBuilder::step()
  \brief Steps to the next l and m
*/

/*! \class Lambda2Builder
 \brief A lambda building class

Objects in this class generate
\f[
            _\pm\lambda_l^m=2\sqrt{\frac{(l-2)!}{(l+2)!}\frac{(2l+1)}{4\pi}\frac{(l-m)!}{(l+m)!}}\left( G^+_{lm} \mp G^-_{lm}\right)
\f]
where
\f[G^+_{lm}(\cos{\theta})=-\left( \frac{l-m^2}{\sin^2{\theta}}+\frac{1}{2}l\left(l-1\right)\right)P_l^m(\cos{\theta})+\left(l+m\right)\frac{\cos{\theta}}{\sin^2{\theta}}P^m_{l-1}(\cos{\theta})
\f]
and
\f[G^-_{lm}(\cos{\theta})=\frac{m}{\sin^2{\theta}}\left(\left(l-1\right)\cos{\theta}P^m_l(\cos{\theta})-\left(l+m\right)P^m_{l-1}(\cos{\theta})\right)
\f]
and \f$P_l^m\f$ are the associated Legendre polynomials.
The interest in all this is that
\f[_{\pm2}Y_l^m(\cos{\theta})=_\pm\lambda_l^m(\cos{\theta})e^{im\phi}
\f]
where \f$_{\pm2}Y_l^m\f$ are the spin-2 spherical harmonics.
*/

/*!   \fn  Lambda2Builder::Lambda2Builder(double theta,int lmax, int mmax)
  \param theta the desired \f$\theta\f$
  \param lmax the maximum value of l
  \param mmax the maximum value of m
*/
/*!   \fn  Lambda2Builder::Lambda2Builder(const Lambda2Builder& lb)
  \param lb the Lambda2Builder object to be copied
*/
/*!   \fn  double Lambda2Builder::lam2lmp() const
  \returns \f$_+\lambda_l^m\f$ for current l and m
*/
/*!   \fn  double Lambda2Builder::lam2lmp(int flip);
  \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$_+\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the current l and m.
  Uses symmetry relations to get this value.
*/
/*!   \fn  double Lambda2Builder::lam2lmp(int l, int m)
  \param l the desired l
  \param m the desired m
  \returns \f$_+\lambda_l^m\f$ for l and m
*/
/*!   \fn  double Lambda2Builder::lam2lmp(int l, int m, int flip)
  \param l the desired l
  \param m the desired m
   \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$_+\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the given l and m.
  Uses symmetry relations to get this value.
*/
/*!   \fn  double Lambda2Builder::lam2lmm() const
  \returns \f$_-\lambda_l^m\f$ for current l and m
*/
/*!   \fn  double Lambda2Builder::lam2lmm(int flip)
  \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$_-\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the current l and m.
  Uses symmetry relations to get this value.
*/
/*!   \fn  double Lambda2Builder::lam2lmm(int l, int m)
  \param l the desired l
  \param m the desired m
  \returns \f$_-\lambda_l^m\f$ for l and m
*/
/*!   \fn  double Lambda2Builder::lam2lmm(int l, int m, int flip)
  \param l the desired l
  \param m the desired m
   \param flip \f$\theta\f$ to be multiplied by sign(flip)
  \returns \f$_-\lambda_l^m(\cos{(sign(flip)\theta}))\f$ for the given l and m.
  Uses symmetry relations to get this value.
)*/