source: Sophya/trunk/SophyaLib/Samba/spheregorski.h@ 470

#ifndef SPHEREGORSKI_SEEN
#define SPHEREGORSKI_SEEN

#include "sphericalmap.h"
#include "tvector.h"
#include "ndatablock.h"

#include "anydataobj.h"
#include "ppersist.h"

// les attributs de classe pix2x_ et pix2y_ sont relatifs a la traduction des 
// indices  RING en indices NESTED (ou l inverse). Ce sont des tableaux 
// d'entiers initialises et remplis par le constructeur. Dans la version 
// FORTRAN de healpix, il s'agissait de tableaux mis en COMMON. Ils etaient 
// initialises au premier appel d'une conversion d indice. Je ne peux pas 
// garder cette procedure car, par exemple, la fonction PixValNest() const 
// n'autorisera pas la constitution de ces tableaux (a cause du const). Ainsi, 
// des la creation d un objet SphereGorski ces tableaux sont calcules, ce qui 
// est, certes, inutile si on ne se sert pas des indices NESTED. Ca ne me 
// rejouit pas, mais c est le seul moyen que j ai trouve pour temir compte de 
// toutes les demandes et des options prises.
// 
//                                                G. Le Meur 

// ***************** CLASSE SphereGorski *****************************

template<class T>
class SphereGorski : public SphericalMap<T>,  public AnyDataObj 
{

public :

SphereGorski();
SphereGorski(int_4 m);
SphereGorski(const SphereGorski<T>& s);
virtual ~SphereGorski();

// ------------------ Definition of PixelMap abstract methods

/* Nombre de pixels du decoupage */
virtual int_4 NbPixels() const;
inline void setNbPixels(int_4 n) {nPix_= n;}

/* Valeur du contenu du pixel d'indice "RING" k  */
virtual T& PixVal(int_4 k);
virtual T const& PixVal(int_4 k) const;

/* Nombre de tranches en theta */ 
int_4 NbThetaSlices() const;
void GetThetaSlice(int_4 index,r_4& theta,TVector<float>& phi,TVector<T>& value) const;

/* Indice "RING" du pixel vers lequel pointe une direction definie par 
ses  coordonnees spheriques */                                    
virtual int_4 PixIndexSph(float theta, float phi) const;

/* Coordonnees spheriques du milieu du pixel d'indice "RING" k   */
virtual void PixThetaPhi(int_4 k, float& teta, float& phi) const;

/* Pixel Solid angle  (steradians) */
virtual r_8 PixSolAngle(int_4 dummy) const;
inline void setPixSolAngle(r_8 x) {omeg_= x;}

// --------------- Specific methods 

virtual void Resize(int_4 m);

inline virtual char* TypeOfMap() const {return "RING";};

/* Valeur du contenu du pixel d'indice "NEST" k                                 */
virtual T& PixValNest(int_4 k);
virtual T const& PixValNest(int_4 k) const;

/* Indice "NEST" du pixel vers lequel pointe une direction definie par 
ses  coordonnees spheriques */                                    
virtual int_4 PixIndexSphNest(float theta, float phi) const;

/* Coordonnees spheriques du milieu du pixel d'indice "NEST" k       */
virtual void PixThetaPhiNest(int_4 k, float& teta, float& phi) const;

/* algorithme de pixelisation */
void Pixelize(int_4); 

/* convertit index nested en ring  */
int_4 NestToRing(int_4 ) const;

/* convertit index ring en nested" */
int_4 RingToNest(int_4 ) const;

/*    analyse en harmoniques spheriques des valeurs des pixels de la 
   sphere : appel du module anafast (Gorski-Hivon) */
//void anharm(int, float, float*);
 
/*synthese  des valeurs des pixels de la sphere par l'intermediaire 
  des coefficients en harmoniques spheriques reconstitues apartir d'un 
  spectre en puissance : appel du module synfast (Gorski-Hivon) */
//void synharm(int, int, float, float*); 

/* retourne/fixe la valeur du parametre Gorski */
inline virtual int_4 SizeIndex() const {return(nSide_);}
inline void setSizeIndex(int_4 n) {nSide_= n;}

/* retourne les pointeurs /remplit les tableaux */
inline const NDataBlock<T>* getDataBlock() const { return (&pixels_); }
inline void setDataBlock(T* data, int_4 m) { pixels_.FillFrom(m,data); }

/* retourne/fixe les parametres des modules anafast et synfast */
void getParafast(int_4& nlmax,int_4& nmmax,int_4& iseed,float& fwhm,float& quadr,float& cut) const;
void setParafast(int_4 nlmax,int_4 nmmax,int_4 iseed,float fwhm,float quadr,float cut,char* filename);

/* retourne/fixe le nom du fichier qui contient le spectre de puissance */
void powfile(char filename[]) const;

/* impression */ 
void print(ostream& os) const;
         
private :

// ------------- méthodes internes ----------------------
void InitNul();

int  nest2ring(int nside, int ipnest) const;
int  ring2nest(int nside, int ipring) const;

int  ang2pix_ring(int nside, double theta, double phi) const;
int  ang2pix_nest(int nside, double theta, double phi) const;
void pix2ang_ring(int nside, int ipix, double& theta,  double& phi) const;
void pix2ang_nest(int nside, int ipix, double& theta, double& phi) const;

// ------------- variables internes -----------------------
int_4 nSide_;
int_4 nPix_;
r_8 omeg_;

NDataBlock<T> pixels_;

int nlmax_;
int nmmax_;
int iseed_;
float fwhm_;
float quadrupole_;
float sym_cut_deg_;
char powFile_[128];
};

//
// ------------- Classe pour la gestion de persistance --
//
template <class T>
class FIO_SphereGorski : public PPersist  
{
public:

FIO_SphereGorski();
FIO_SphereGorski(string const & filename); 
FIO_SphereGorski(const SphereGorski<T>& obj);
FIO_SphereGorski(SphereGorski<T>* obj);
virtual ~FIO_SphereGorski();
virtual AnyDataObj* DataObj();
inline operator SphereGorski<T>() { return(*dobj); }
inline SphereGorski<T> getObj() { return(*dobj); }

protected :

virtual void ReadSelf(PInPersist&);           
virtual void WriteSelf(POutPersist&) const;  
SphereGorski<T>* dobj;
bool ownobj;
};

//
// ------------- Classe PIXELS_XY -----------------------
//
class PIXELS_XY
{

public :

static PIXELS_XY& instance();

NDataBlock<int> pix2x_;
NDataBlock<int> pix2y_;
NDataBlock<int> x2pix_;
NDataBlock<int> y2pix_;

private :

PIXELS_XY();
void mk_pix2xy();
void mk_xy2pix();
};
#endif