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

#ifndef SPHEREGORSKI_SEEN
#define SPHEREGORSKI_SEEN

#include "sphericalmap.h"
#include "cvector.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 

class SphereGorski : public SphericalMap {
public :


SphereGorski();
SphereGorski(int_4 m);
SphereGorski(char* flnm);
virtual ~SphereGorski();

// ------------ Persistence handling
     enum {classId = 0xF002 };

int_4                  ClassId() const        { return classId; }

virtual void           WriteSelf(POutPersist&) const;
virtual void           ReadSelf(PInPersist&);

// -----------   FITS IO (1D FITS ARRAY)
void ReadFits(char flnm[] );
void WriteFits(char flnm[]);

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

/* Nombre de pixels du decoupage                                         */
virtual  int_4           NbPixels() const;

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

int_4 NbThetaSlices() const;
void  GetThetaSlice(int_4 index, r_4& theta, Vector& phi, Vector& 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; 


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

    // NEST indexing
/* Valeur du contenu du pixel d'indice "NEST" k                                 */
virtual  r_8&       PixValNest(int_4 k);
virtual  r_8 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 );


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


/*    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* ); 


// ------------- méthodes internes ----------------------
          
private :

void          InitNul();
void          Clear();
int           nest2ring(int nside, int ipnest) const;
void          mk_pix2xy(int *pix2x,int *pix2y);
int           ring2nest(int nside, int ipring) const;
void          mk_xy2pix(int *x2pix, int *y2pix);
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 -----------------------
                   private :
class Pix2XY {
  public :
Pix2XY() { 
    pix2x_=new int[1024];
    pix2y_=new int[1024];
    mk_pix2xy(pix2x_,pix2y_); 
  }
~Pix2XY() {
  if (pix2x_)  delete[] pix2x_;
  if (pix2y_)  delete[] pix2y_;
  }
void          pix2ang_nest(int nside, int ipix, double& theta, double& phi) const;

private : 
int *pix2x_;
int *pix2y_;
void          mk_pix2xy(int *pix2x,int *pix2y);

};

int nlmax_;
int nmmax_;
int iseed_;
int  nSide_;
int_4  nPix_;
int *pix2x_;
int *pix2y_;
int *x2pix_;
int *y2pix_;
float fwhm_;
float quadrupole_;
float sym_cut_deg_;
r_8  omeg_;
r_8*    mPix_;
char powFile_[128];
Pix2XY *pix2xy_;
};
#endif