Changeset 1217 in Sophya for trunk/SophyaLib
- Timestamp:
- Oct 3, 2000, 10:19:18 AM (25 years ago)
- Location:
- trunk/SophyaLib/SkyMap
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/SophyaLib/SkyMap/localmap.cc
r979 r1217 9 9 10 10 11 //***************************************************************************** 12 //++ 13 // Class LocalMap 14 // 15 // include localmap.h 16 // 17 // A local map of a region of the sky, in cartesian coordinates. 18 // It has an origin in (theta0, phi0), mapped to pixel(x0, y0) 19 // (x0, y0 might be outside of this local map) 20 // default value of (x0, y0) is middle of the map, center of 21 // pixel(nx/2, ny/2) 22 // 23 // A local map is a 2 dimensional array, with i as column index and j 24 // as row index. The map is supposed to lie on a plan tangent to the 25 // celestial sphere in a point whose coordinates are (x0,y0) on the local 26 // map and (theta0, phi0) on the sphere. The range of the map is defined 27 // by two values of angles covered respectively by all the pixels in 28 // x direction and all the pixels in y direction (SetSize()). 29 // 30 // A "reference plane" is considered : this plane is tangent to the 31 // celestial sphere in a point with angles theta=Pi/2 and phi=0. This 32 // point is the origine of coordinates is of the reference plane. The 33 // x-axis is the tangent parallel to the equatorial line and oriented 34 // toward the increasing phi's ; the y-axis is parallel to the meridian 35 // line and oriented toward the north pole. 36 // 37 // Internally, a map is first defined within this reference plane and 38 // tranported until the point (theta0, phi0) in such a way that both 39 // axes are kept parallel to meridian and parallel lines of the sphere. 40 // The user can define its own map with axes rotated with respect to 41 // reference axes (this rotation is characterized by angle between 42 // the local parallel line and the wanted x-axis-- see method 43 // SetOrigin(...)) 44 // 45 // 46 47 // 48 //-- 49 //++ 50 // 51 // Links Parents 52 // 53 // PixelMap 54 // 55 //-- 56 //++ 57 // 58 // Links Childs 59 // 60 // 61 //-- 62 //++ 63 // Titre Constructors 64 //-- 65 //++ 11 /*! 12 \class SOPHYA::LocalMap 13 14 A local map has an origin in (theta0, phi0), mapped to pixel(x0, y0) 15 (x0, y0 might be outside of this local map) 16 default value of (x0, y0) is middle of the map, center of pixel(nx/2, ny/2) 17 A local map is a 2 dimensional array, with i as column index and j 18 as row index. The map is supposed to lie on a plan tangent to the 19 celestial sphere in a point whose coordinates are (x0,y0) on the local 20 map and (theta0, phi0) on the sphere. The range of the map is defined 21 by two values of angles covered respectively by all the pixels in 22 x direction and all the pixels in y direction (SetSize()). 23 24 A "reference plane" is considered : this plane is tangent to the 25 celestial sphere in a point with angles theta=Pi/2 and phi=0. This 26 point is the origine of coordinates is of the reference plane. The 27 x-axis is the tangent parallel to the equatorial line and oriented 28 toward the increasing phi's ; the y-axis is parallel to the meridian 29 line and oriented toward the north pole. 30 31 Internally, a map is first defined within this reference plane and 32 tranported until the point (theta0, phi0) in such a way that both 33 axes are kept parallel to meridian and parallel lines of the sphere. 34 The user can define its own map with axes rotated with respect to 35 reference axes (this rotation is characterized by angle between 36 the local parallel line and the wanted x-axis-- see method 37 SetOrigin(...)) 38 */ 66 39 template<class T> 67 40 LocalMap<T>::LocalMap() : pixels_() 68 //69 //70 //--71 41 { 72 42 InitNul(); 73 // pixels_.Reset(); 74 } 75 76 //++ 43 } 44 77 45 template<class T> 78 46 LocalMap<T>::LocalMap(int_4 nx, int_4 ny) : nSzX_(nx), nSzY_(ny) 79 //80 //81 //--82 47 { 83 48 InitNul(); … … 87 52 } 88 53 89 //++90 54 template<class T> 91 55 LocalMap<T>::LocalMap(const LocalMap<T>& lm, bool share) 92 56 : pixels_(lm.pixels_, share) 93 94 //95 // copy constructor96 //--97 57 { 98 58 … … 101 61 } 102 62 103 //++104 63 template<class T> 105 64 LocalMap<T>::LocalMap(const LocalMap<T>& lm) 106 65 : pixels_(lm.pixels_) 107 66 108 //109 // copy constructor110 //--111 67 { 112 68 … … 115 71 } 116 72 117 //++118 // Titre Destructor119 //--120 //++121 73 template<class T> 122 74 LocalMap<T>::~LocalMap() 123 //124 //--125 75 { 126 76 InitNul(); … … 129 79 130 80 131 //++ 132 // Titre Public Methods 133 //-- 134 135 //++ 81 /*! \fn void SOPHYA::LocalMap::ReSize(int_4 nx, int_4 ny) 82 83 Resize storage area for pixels 84 */ 136 85 template<class T> 137 86 void LocalMap<T>::ReSize(int_4 nx, int_4 ny) 138 //139 // Resize storage area for pixels140 //--141 87 { 142 88 InitNul(); … … 188 134 189 135 190 //++ 136 /*! \fn int_4 SOPHYA::LocalMap::NbPixels() const 137 138 \return number of pixels 139 */ 191 140 template<class T> 192 141 int_4 LocalMap<T>::NbPixels() const 193 //194 // Return number of pixels195 //--196 142 { 197 143 return(nPix_); 198 144 } 199 145 200 //++ 146 /*! \fn T& SOPHYA::LocalMap::PixVal(int_4 k) 147 148 \return value of pixel with index k 149 */ 201 150 template<class T> 202 151 T& LocalMap<T>::PixVal(int_4 k) 203 //204 // Return value of pixel with index k205 //--206 152 { 207 153 if((k < 0) || (k >= nPix_)) 208 154 { 209 155 cout << " LocalMap::PIxVal : exceptions a mettre en place" <<endl; 210 // THROW(out_of_range("LocalMap::PIxVal Pixel index out of range"));211 156 THROW(rangeCheckErr); 212 //throw "LocalMap::PIxVal Pixel index out of range";213 157 } 214 158 return(pixels_(k)); 215 159 } 216 160 217 //++ 218 161 /*! \fn T const& SOPHYA::LocalMap::PixVal(int_4 k) const 162 const version of previous method 163 */ 219 164 template<class T> 220 165 T const& LocalMap<T>::PixVal(int_4 k) const 221 //222 // const version of previous method223 //--224 166 { 225 167 if((k < 0) || (k >= nPix_)) 226 168 { 227 169 cout << " LocalMap::PIxVal : exceptions a mettre en place" <<endl; 228 // THROW(out_of_range("LocalMap::PIxVal Pixel index out of range"));229 230 170 throw "LocalMap::PIxVal Pixel index out of range"; 231 171 } … … 233 173 } 234 174 175 /*! \fn bool SOPHYA::LocalMap::ContainsSph(double theta, double phi) const 176 177 \return true if teta,phi in map 178 */ 235 179 template<class T> 236 180 bool LocalMap<T>::ContainsSph(double theta, double phi) const … … 240 184 } 241 185 242 //++ 186 /*! \fn int_4 SOPHYA::LocalMap::PixIndexSph(double theta,double phi) const 187 188 \return index of the pixel with spherical coordinates (theta,phi) 189 */ 243 190 template<class T> 244 191 int_4 LocalMap<T>::PixIndexSph(double theta,double phi) const 245 //246 // Return index of the pixel with spherical coordinates (theta,phi)247 //--248 192 { 249 193 int_4 i,j; … … 284 228 } 285 229 286 //++ 230 /*! \fn void SOPHYA::LocalMap::PixThetaPhi(int_4 k,double& theta,double& phi) const 231 232 \return (theta, phi) coordinates of pixel with index k 233 */ 287 234 template<class T> 288 235 void LocalMap<T>::PixThetaPhi(int_4 k,double& theta,double& phi) const 289 //290 // Return (theta, phi) coordinates of pixel with index k291 //--292 236 { 293 237 if(!(originFlag_) || !(extensFlag_)) … … 311 255 } 312 256 257 /*! \fn T SOPHYA::LocalMap::SetPixels(T v) 258 259 Set all pixels to value v 260 */ 313 261 template <class T> 314 262 T LocalMap<T>::SetPixels(T v) … … 318 266 } 319 267 320 //++ 268 /*! \fn double SOPHYA::LocalMap::PixSolAngle(int_4 k) const 269 270 Pixel Solid angle (steradians) 271 272 All the pixels have not necessarly the same size in (theta, phi) 273 because of the projection scheme which is not yet fixed. 274 */ 321 275 template<class T> 322 276 double LocalMap<T>::PixSolAngle(int_4 k) const 323 //324 // Pixel Solid angle (steradians)325 // All the pixels have not necessarly the same size in (theta, phi)326 // because of the projection scheme which is not yet fixed.327 //--328 277 { 329 278 int_4 i,j; … … 355 304 } 356 305 357 //++ 306 /*! \fn void SOPHYA::LocalMap::SetOrigin(double theta0,double phi0,double angle) 307 308 set the referential of the map (angles in degrees) 309 310 (default x0=siz_x/2, y0=siz_y/2) 311 */ 358 312 template<class T> 359 313 void LocalMap<T>::SetOrigin(double theta0,double phi0,double angle) 360 //361 // set the referential of the map (angles in degrees)362 // (default x0=siz_x/2, y0=siz_y/2)363 //--364 314 { 365 315 theta0_= theta0; … … 374 324 } 375 325 376 //++ 326 /*! \fn void SOPHYA::LocalMap::SetOrigin(double theta0,double phi0,int_4 x0,int_4 y0,double angle) 327 328 set the referential of the map (angles in degrees) 329 */ 377 330 template<class T> 378 331 void LocalMap<T>::SetOrigin(double theta0,double phi0,int_4 x0,int_4 y0,double angle) 379 //380 // set the referential of the map (angles in degrees)381 //--382 332 { 383 333 theta0_= theta0; … … 392 342 } 393 343 394 //++ 344 /*! \fn void SOPHYA::LocalMap::SetSize(double angleX,double angleY) 345 346 angle range of tthe map (angles in degrees) 347 */ 395 348 template<class T> 396 349 void LocalMap<T>::SetSize(double angleX,double angleY) 397 //398 // angle range of tthe map (angles in degrees)399 //--400 350 { 401 351 angleX_= angleX; … … 410 360 } 411 361 412 //++ 362 /*! \fn void SOPHYA::LocalMap::Project(SphericalMap<T>& sphere) const 363 364 Projection to a spherical map 365 */ 413 366 template<class T> 414 367 void LocalMap<T>::Project(SphericalMap<T>& sphere) const 415 //416 // Projection to a spherical map417 //--418 368 { 419 369 for(int m = 0; m < nPix_; m++) … … 440 390 } 441 391 442 //++ 392 393 /*! \fn void SOPHYA::LocalMap::Getij(int_4 k,int_4& i,int_4& j) const 394 395 \return 2 indices corresponding to the pixel number k 396 */ 443 397 template<class T> 444 398 void LocalMap<T>::Getij(int_4 k,int_4& i,int_4& j) const 445 //446 // Return 2 indices corresponding to the pixel number k447 //--448 399 { 449 400 i= (k+1)%nSzX_-1; … … 452 403 } 453 404 454 //++ 405 /*! \fn void SOPHYA::LocalMap::ReferenceToUser(double& theta,double& phi) const 406 407 Transform a pair of coordinates (theta, phi) given in 408 reference coordinates into map coordinates 409 */ 455 410 template<class T> 456 411 void LocalMap<T>::ReferenceToUser(double& theta,double& phi) const 457 //458 // Transform a pair of coordinates (theta, phi) given in459 // reference coordinates into map coordinates460 //--461 412 { 462 413 if(theta > Pi || theta < 0. || phi < 0. || phi >= 2*Pi) … … 491 442 } 492 443 493 //++ 444 /*! \fn void SOPHYA::LocalMap::UserToReference(double& theta,double& phi) const 445 446 Transform a pair of coordinates (theta, phi) given in 447 map coordinates into reference coordinates 448 */ 494 449 template<class T> 495 450 void LocalMap<T>::UserToReference(double& theta,double& phi) const 496 //497 // Transform a pair of coordinates (theta, phi) given in498 // map coordinates into reference coordinates499 //--500 451 { 501 452 if(theta > Pi || theta < 0. || phi < 0. || phi >= 2*Pi) … … 534 485 } 535 486 536 //++ 487 /*! \fn void SOPHYA::LocalMap::PixProjToAngle(double x,double y,double& theta,double& phi) const 488 489 Given coordinates in pixel units in the REFERENCE PLANE, return 490 (theta, phi) in "absolute" referential theta=pi/2 ,phi=0. 491 */ 537 492 template<class T> 538 493 void LocalMap<T>::PixProjToAngle(double x,double y,double& theta,double& phi) const 539 //540 //541 // Given coordinates in pixel units in the REFERENCE PLANE, return542 // (theta, phi) in "absolute" referential theta=pi/2 ,phi=0.543 //--544 494 { 545 495 theta= Pi*0.5-atan(2.*y*tgAngleY_/(double)nSzY_); … … 548 498 } 549 499 550 //++ 500 /*! \fn void SOPHYA::LocalMap::AngleProjToPix(double theta,double phi,double& x,double& y) const 501 502 Given coordinates (theta, phi) in "absolute" referential 503 theta=pi/2 ,phi=0 return pixel indices (i,j) in the REFERENCE PLANE. 504 */ 551 505 template<class T> 552 506 void LocalMap<T>::AngleProjToPix(double theta,double phi,double& x,double& y) const 553 //554 // Given coordinates (theta, phi) in "absolute" referential555 // theta=pi/2 ,phi=0 return pixel indices (i,j) in the REFERENCE PLANE.556 //--557 507 { 558 508 if(phi >= Pi) phi-= DeuxPi; -
trunk/SophyaLib/SkyMap/localmap.h
r1163 r1217 10 10 #include "ppersist.h" 11 11 12 //! A local map of a region of the sky, in cartesian coordinates.13 /*! A local map has an origin in (theta0, phi0), mapped to pixel(x0, y0)14 (x0, y0 might be outside of this local map)15 default value of (x0, y0) is middle of the map, center of pixel(nx/2, ny/2)16 A local map is a 2 dimensional array, with i as column index and j17 as row index. The map is supposed to lie on a plan tangent to the18 celestial sphere in a point whose coordinates are (x0,y0) on the local19 map and (theta0, phi0) on the sphere. The range of the map is defined20 by two values of angles covered respectively by all the pixels in21 x direction and all the pixels in y direction (SetSize()).22 23 A "reference plane" is considered : this plane is tangent to the24 celestial sphere in a point with angles theta=Pi/2 and phi=0. This25 point is the origine of coordinates is of the reference plane. The26 x-axis is the tangent parallel to the equatorial line and oriented27 toward the increasing phi's ; the y-axis is parallel to the meridian28 line and oriented toward the north pole.29 30 Internally, a map is first defined within this reference plane and31 tranported until the point (theta0, phi0) in such a way that both32 axes are kept parallel to meridian and parallel lines of the sphere.33 The user can define its own map with axes rotated with respect to34 reference axes (this rotation is characterized by angle between35 the local parallel line and the wanted x-axis-- see method36 SetOrigin(...))37 */38 //39 // la carte est consideree comme un tableau a deux indices i et j, i etant40 // indice de colonne et j indice de ligne. La carte est supposee resider41 // dans un plan tangent, dont le point de tangence est repere (x0,y0) dans42 // la carte et (theta0, phi0) sur la sphere celeste. L extension de la43 // carte est definie par les valeurs de deux angles couverts respectivement44 // par la totalite des pixels en x de la carte et la totalite des pixels45 // en y. (SetSize()).46 // On considere un "plan de reference" : plan tangent a la sphere celeste47 // aux angles theta=Pi/2 et phi=0. Dans ce plan L origine des coordonnees48 // est le point de tangence. L axe Ox est la tangente parallele a49 // lequateur, dirige vers les phi croissants, l axe Oy est parallele50 // au meridien, dirige vers le pole nord.51 // De maniere interne a la classe une carte est definie dans ce plan de52 // reference et transportee jusqu au point (theta0, phi0) de sorte que les // axes restent paralleles aux meridiens et paralleles. L utilisateur peut53 // definir sa carte selon un repere en rotation par rapport au repere de54 // reference (par l angle entre le parallele et l axe Ox souhaite --55 // methode SetOrigin(...))56 12 57 13 … … 61 17 62 18 namespace SOPHYA { 19 20 21 /* Class LocalMap */ 63 22 64 23 … … 76 35 77 36 inline virtual bool IsTemp(void) const { return pixels_.IsTemp();} 37 78 38 /*! Setting blockdata to temporary (see ndatablock documentation) */ 79 39 inline virtual void SetTemp(bool temp=false) const {pixels_.SetTemp(temp);}; … … 89 49 // ---------- Definition of PixelMap abstract methods ------- 90 50 91 /* return/set the number of pixels */92 /*! Return number of pixels */93 51 virtual int_4 NbPixels() const; // D.Y. int change en int_4 rationalisation Mac 94 52 95 /* return the value of pixel number k */96 /*! Return value of pixel with index k */97 53 virtual T& PixVal(int_4 k); 98 /*! const version of previous method */99 54 virtual T const& PixVal(int_4 k) const; 100 55 101 /* Return true if teta,phi in map */102 56 virtual bool ContainsSph(double theta, double phi) const; 103 /* return the index of pixel at (theta,phi) */104 /*! Return index of the pixel with spherical coordinates (theta,phi) */105 57 virtual int_4 PixIndexSph(double theta,double phi) const; 106 58 107 /* return the spherical coordinates of center of pixel number k */108 /*! Return (theta, phi) coordinates of pixel with index k */109 59 virtual void PixThetaPhi(int_4 k,double& theta,double& phi) const; 110 60 111 /*! Set all pixels to value v */112 61 virtual T SetPixels(T v); 113 62 114 /* return the Pixel Solid angle (steradians) */115 /*! Pixel Solid angle (steradians)116 117 All the pixels have not necessarly the same size in (theta, phi)118 because of the projection scheme which is not yet fixed.119 */120 63 virtual double PixSolAngle(int_4 k) const; 121 64 122 65 // ---------- Specific methods ------------------------------ 123 66 124 /*! Resize storage area for pixels */125 67 void ReSize(int_4 nx, int_4 ny); 126 68 127 69 inline virtual char* TypeOfMap() const {return "LOCAL";}; 128 70 129 /* Origin (with angle between x axis and phi axis, in degrees) x0,y0 the default: middle of map*/130 /*! set the referential of the map (angles in degrees)131 132 (default x0=siz_x/2, y0=siz_y/2)133 */134 71 virtual void SetOrigin(double theta=90.,double phi=0.,double angle=0.); 135 /*! set the referential of the map (angles in degrees) */136 72 virtual void SetOrigin(double theta,double phi,int_4 x0,int_4 y0,double angle=0.); 137 73 138 /* Pixel size (degres) */139 /*! angle range of tthe map (angles in degrees) */140 74 virtual void SetSize(double angleX,double angleY); 141 75 142 /* Check to see if the local mapping is done */76 /*! Check to see if the local mapping is done */ 143 77 inline bool LocalMap_isDone() const {return(originFlag_ && extensFlag_);}; 144 78 145 /*! Projection to a spherical map */146 79 virtual void Project(SphericalMap<T>& sphere) const; 147 80 … … 149 82 -> static method, or separate class */ 150 83 151 /* provides a integer characterizing the pixelization refinement (here : number of pixels) */84 /*! provides a integer characterizing the pixelization refinement (here : number of pixels) */ 152 85 inline virtual int_4 SizeIndex() const {return(nPix_);} 153 86 inline int_4 Size_x() const {return nSzX_;} … … 180 113 181 114 void InitNul(); 182 /*! Return 2 indices corresponding to the pixel number k */183 115 void Getij(int_4 k,int_4& i,int_4& j) const; 184 /*! Transform a pair of coordinates (theta, phi) given in185 reference coordinates into map coordinates186 */187 116 void ReferenceToUser(double& theta,double& phi) const; 188 /*! Transform a pair of coordinates (theta, phi) given in189 map coordinates into reference coordinates190 */191 117 void UserToReference(double& theta,double& phi) const; 192 /*! Given coordinates in pixel units in the REFERENCE PLANE, return193 (theta, phi) in "absolute" referential theta=pi/2 ,phi=0.194 */195 118 void PixProjToAngle(double x,double y,double& theta,double& phi) const; 196 /*! Given coordinates (theta, phi) in "absolute" referential197 theta=pi/2 ,phi=0 return pixel indices (i,j) in the REFERENCE PLANE.198 */199 119 void AngleProjToPix(double theta,double phi,double& x,double& y) const; 200 120 -
trunk/SophyaLib/SkyMap/spherehealpix.cc
r1196 r1217 17 17 using namespace SOPHYA; 18 18 19 //******************************************************************* 20 //++ 21 // Class SphereHEALPix 22 // 23 // include SphereHealpix.h strutil.h 24 // 25 // Pixelisation Gorski 26 // 27 // 28 //| ----------------------------------------------------------------------- 29 //| version 0.8.2 Aug97 TAC Eric Hivon, Kris Gorski 30 //| ----------------------------------------------------------------------- 31 // 32 // the sphere is split in 12 diamond-faces containing nside**2 pixels each 33 // 34 // the numbering of the pixels (in the nested scheme) is similar to 35 // quad-cube 36 // In each face the first pixel is in the lowest corner of the diamond 37 // 38 // the faces are (x,y) coordinate on each face 39 //| . . . . <--- North Pole 40 //| / \ / \ / \ / \ ^ ^ 41 //| . 0 . 1 . 2 . 3 . <--- z = 2/3 \ / 42 //| \ / \ / \ / \ / y \ / x 43 //| 4 . 5 . 6 . 7 . 4 <--- equator \ / 44 //| / \ / \ / \ / \ \/ 45 //| . 8 . 9 .10 .11 . <--- z = -2/3 (0,0) : lowest corner 46 //| \ / \ / \ / \ / 47 //| . . . . <--- South Pole 48 //| 49 // phi:0 2Pi 50 // 51 // in the ring scheme pixels are numbered along the parallels 52 // the first parallel is the one closest to the north pole and so on 53 // on each parallel, pixels are numbered starting from the one closest 54 // to phi = 0 55 // 56 // nside MUST be a power of 2 (<= 8192) 57 //-- 58 //++ 59 // 60 // Links Parents 61 // 62 // SphericalMap 63 //-- 19 /*! 20 \class SOPHYA::SphereHEALPix 21 22 ******************************************************************* 23 Pixelisation Gorski 24 25 \verbatim 26 27 ----------------------------------------------------------------------- 28 version 0.8.2 Aug97 TAC Eric Hivon, Kris Gorski 29 ----------------------------------------------------------------------- 30 31 the sphere is split in 12 diamond-faces containing nside**2 pixels each 32 the numbering of the pixels (in the nested scheme) is similar to 33 quad-cube 34 In each face the first pixel is in the lowest corner of the diamond 35 the faces are (x,y) coordinate on each face 36 37 . . . . <--- North Pole 38 / \ / \ / \ / \ ^ ^ 39 . 0 . 1 . 2 . 3 . <--- z = 2/3 \ / 40 \ / \ / \ / \ / y \ / x 41 4 . 5 . 6 . 7 . 4 <--- equator \ / 42 / \ / \ / \ / \ \/ 43 . 8 . 9 .10 .11 . <--- z = -2/3 (0,0) : lowest corner 44 \ / \ / \ / \ / 45 . . . . <--- South Pole 46 \endverbatim 47 phi:0 2Pi 48 49 in the ring scheme pixels are numbered along the parallels 50 the first parallel is the one closest to the north pole and so on 51 on each parallel, pixels are numbered starting from the one closest 52 to phi = 0 53 54 nside MUST be a power of 2 (<= 8192) 55 56 */ 64 57 65 58 /* --Methode-- */ 66 //++67 // Titre Constructors68 //--69 //++70 59 71 60 template<class T> … … 73 62 sliceLenght_() 74 63 75 //-- 64 76 65 { 77 66 InitNul(); 78 // SetTemp(false); 79 } 80 81 //++ 67 } 68 69 /*! \fn SOPHYA::SphereHEALPix::SphereHEALPix(int_4 m) 70 71 \param <m> : "nside" of the Healpix algorithm 72 73 The total number of pixels will be Npix = 12*nside**2 74 75 nside MUST be a power of 2 (<= 8192) 76 */ 77 82 78 template<class T> 83 79 SphereHEALPix<T>::SphereHEALPix(int_4 m) 84 80 85 // m is the "nside" of the Gorski algorithm86 //87 // The total number of pixels will be Npix = 12*nside**288 //89 // nside MUST be a power of 2 (<= 8192)90 //--91 81 { 92 82 … … 202 192 //-- 203 193 204 //++ 194 /*! \fn void SOPHYA::SphereHEALPix::Resize(int_4 m) 195 \param <m> "nside" of the Gorski algorithm 196 197 The total number of pixels will be Npix = 12*nside**2 198 199 nside MUST be a power of 2 (<= 8192) 200 */ 205 201 template<class T> 206 202 void SphereHEALPix<T>::Resize(int_4 m) 207 203 208 // m is the "nside" of the Gorski algorithm 209 // 210 // The total number of pixels will be Npix = 12*nside**2 211 // 212 // nside MUST be a power of 2 (<= 8192) 213 //-- 204 214 205 { 215 206 if (m<=0 || m> 8192) { … … 270 261 271 262 /* --Methode-- */ 272 //++ 263 /* Nombre de pixels du decoupage */ 264 /*! \fn int_4 SOPHYA::SphereHEALPix::NbPixels() const 265 266 Return number of pixels of the splitting 267 */ 273 268 template<class T> 274 269 int_4 SphereHEALPix<T>::NbPixels() const 275 276 // Retourne le nombre de pixels du découpage277 //--278 270 { 279 271 return(nPix_); 280 272 } 281 273 282 //++ 274 275 /*! \fn uint_4 SOPHYA::SphereHEALPix::NbThetaSlices() const 276 277 \return number of slices in theta direction on the sphere 278 */ 283 279 template<class T> 284 280 uint_4 SphereHEALPix<T>::NbThetaSlices() const 285 286 // Return number of slices in theta direction on the sphere287 //--288 281 { 289 282 uint_4 nbslices = uint_4(4*nSide_-1); … … 296 289 } 297 290 298 //++ 291 /*! \fn void SOPHYA::SphereHEALPix::GetThetaSlice(int_4 index,r_8& theta,TVector<r_8>& phi,TVector<T>& value) const 292 293 For a theta-slice with index 'index', return : 294 295 the corresponding "theta" 296 297 a vector containing the phi's of the pixels of the slice 298 299 a vector containing the corresponding values of pixels 300 */ 299 301 template<class T> 300 302 void SphereHEALPix<T>::GetThetaSlice(int_4 index,r_8& theta,TVector<r_8>& phi,TVector<T>& value) const 301 303 302 // For a theta-slice with index 'index', return :303 //304 // the corresponding "theta"305 //306 // a vector containing the phi's of the pixels of the slice307 //308 // a vector containing the corresponding values of pixels309 //310 //--311 304 { 312 305 313 306 if (index<0 || index >= NbThetaSlices()) 314 307 { 315 // THROW(out_of_range("SphereHEALPix::PIxVal Pixel index out of range"));316 308 cout << " SphereHEALPix::GetThetaSlice : Pixel index out of range" <<endl; 317 309 throw RangeCheckError(" SphereHEALPix::GetThetaSlice : Pixel index out of range"); … … 335 327 theta= TH; 336 328 } 337 //++ 338 //++ 329 /*! \fn void SOPHYA::SphereHEALPix::GetThetaSlice(int_4 sliceIndex,r_8& theta, r_8& phi0, TVector<int_4>& pixelIndices,TVector<T>& value) const 330 331 For a theta-slice with index 'index', return : 332 333 the corresponding "theta" 334 335 the corresponding "phi" for first pixel of the slice 336 337 a vector containing indices of the pixels of the slice 338 339 (equally distributed in phi) 340 341 a vector containing the corresponding values of pixels 342 */ 339 343 340 344 template<class T> 341 345 void SphereHEALPix<T>::GetThetaSlice(int_4 sliceIndex,r_8& theta, r_8& phi0, TVector<int_4>& pixelIndices,TVector<T>& value) const 342 346 343 // For a theta-slice with index 'sliceIndex', return :344 //345 // the corresponding "theta"346 // the corresponding "phi" for first pixel of the slice347 //348 // a vector containing the indices of the pixels of the slice349 // (equally distributed in phi)350 //351 // a vector containing the corresponding values of pixels352 //353 //--354 347 { 355 348 356 349 if (sliceIndex<0 || sliceIndex >= NbThetaSlices()) 357 350 { 358 // THROW(out_of_range("SphereHEALPix::PIxVal Pixel index out of range"));359 351 cout << " SphereHEALPix::GetThetaSlice : Pixel index out of range" <<endl; 360 352 throw RangeCheckError(" SphereHEALPix::GetThetaSlice : Pixel index out of range"); … … 372 364 PixThetaPhi(iring, theta, phi0); 373 365 } 374 //++ 366 375 367 template<class T> 376 368 void SphereHEALPix<T>::SetThetaSlices() 377 378 //--379 369 { 380 370 sliceBeginIndex_.ReSize(4*nSide_-1); … … 399 389 } 400 390 401 /* --Methode-- */ 402 //++ 391 392 /*! \fn T& SOPHYA::SphereHEALPix::PixVal(int_4 k) 393 394 \return value of pixel with "RING" index k 395 */ 403 396 template<class T> 404 397 T& SphereHEALPix<T>::PixVal(int_4 k) 405 398 406 // Return value of pixel with "RING" index k407 //--408 399 { 409 400 if((k < 0) || (k >= nPix_)) … … 414 405 } 415 406 416 /* --Methode-- */ 417 //++ 407 /*! \fn T const& SOPHYA::SphereHEALPix::PixVal(int_4 k) const 408 409 \return value of pixel with "RING" index k 410 */ 418 411 template<class T> 419 412 T const& SphereHEALPix<T>::PixVal(int_4 k) const 420 413 421 // Return value of pixel with "RING" index k422 //--423 414 { 424 415 if((k < 0) || (k >= nPix_)) … … 429 420 } 430 421 431 //++ 422 /*! \fn T& SOPHYA::SphereHEALPix::PixValNest(int_4 k) 423 424 \return value of pixel with "NESTED" index k 425 */ 432 426 template<class T> 433 427 T& SphereHEALPix<T>::PixValNest(int_4 k) 434 428 435 // Return value of pixel with "NESTED" index k436 429 //-- 437 430 { … … 442 435 return pixels_(nest2ring(nSide_,k)); 443 436 } 444 //++ 437 438 /*! \fn T const& SOPHYA::SphereHEALPix::PixValNest(int_4 k) const 439 440 \return value of pixel with "NESTED" index k 441 */ 445 442 446 443 template<class T> 447 444 T const& SphereHEALPix<T>::PixValNest(int_4 k) const 448 445 449 // Return value of pixel with "NESTED" index k450 //--451 446 { 452 447 if((k < 0) || (k >= nPix_)) … … 458 453 } 459 454 455 /*! \fn bool SOPHYA::SphereHEALPix::ContainsSph(double theta, double phi) const 456 457 \return true if teta,phi in map 458 */ 459 template<class T> 460 bool SphereHEALPix<T>::ContainsSph(double theta, double phi) const 461 { 462 return(true); 463 } 464 465 /*! \fn int_4 SOPHYA::SphereHEALPix::PixIndexSph(double theta,double phi) const 466 467 \return "RING" index of the pixel corresponding to direction (theta, phi). 468 */ 469 template<class T> 470 int_4 SphereHEALPix<T>::PixIndexSph(double theta,double phi) const 471 472 { 473 return ang2pix_ring(nSide_,theta,phi); 474 } 475 476 /*! \fn int_4 SOPHYA::SphereHEALPix::PixIndexSphNest(double theta,double phi) const 477 478 \return "NESTED" index of the pixel corresponding to direction (theta, phi). 479 */ 480 template<class T> 481 int_4 SphereHEALPix<T>::PixIndexSphNest(double theta,double phi) const 482 483 { 484 return ang2pix_nest(nSide_,theta,phi); 485 } 486 487 460 488 /* --Methode-- */ 461 //++ 462 template<class T> 463 bool SphereHEALPix<T>::ContainsSph(double /*theta*/, double /*phi*/) const 464 //-- 465 { 466 return(true); 467 } 468 469 /* --Methode-- */ 470 //++ 471 template<class T> 472 int_4 SphereHEALPix<T>::PixIndexSph(double theta,double phi) const 473 474 // Return "RING" index of the pixel corresponding to 475 // direction (theta, phi). 476 //-- 477 { 478 return ang2pix_ring(nSide_,theta,phi); 479 } 480 481 //++ 482 template<class T> 483 int_4 SphereHEALPix<T>::PixIndexSphNest(double theta,double phi) const 484 485 // Return "NESTED" index of the pixel corresponding to 486 // direction (theta, phi). 487 //-- 488 { 489 return ang2pix_nest(nSide_,theta,phi); 490 } 491 492 493 /* --Methode-- */ 494 //++ 489 /*! \fn void SOPHYA::SphereHEALPix::PixThetaPhi(int_4 k,double& theta,double& phi) const 490 \return (theta,phi) coordinates of middle of pixel with "RING" index k 491 */ 495 492 template<class T> 496 493 void SphereHEALPix<T>::PixThetaPhi(int_4 k,double& theta,double& phi) const 497 498 // Return (theta,phi) coordinates of middle of pixel with "RING" index k499 //--500 494 { 501 495 pix2ang_ring(nSide_,k,theta,phi); 502 496 } 503 497 498 /*! \fn T SOPHYA::SphereHEALPix::SetPixels(T v) 499 500 Set all pixels to value v 501 */ 504 502 template <class T> 505 503 T SphereHEALPix<T>::SetPixels(T v) … … 509 507 } 510 508 511 //++ 512 template<class T> 513 double SphereHEALPix<T>::PixSolAngle(int_4 /*dummy*/) const 514 // Pixel Solid angle (steradians) 515 // All the pixels have the same solid angle. The dummy argument is 516 // for compatibility with eventual pixelizations which would not 517 // fulfil this requirement. 518 //-- 519 { 520 return omeg_; 521 } 522 523 //++ 509 510 511 /*! \fn void SOPHYA::SphereHEALPix::PixThetaPhiNest(int_4 k,double& theta,double& phi) const 512 513 \return (theta,phi) coordinates of middle of pixel with "NESTED" index k 514 */ 524 515 template<class T> 525 516 void SphereHEALPix<T>::PixThetaPhiNest(int_4 k,double& theta,double& phi) const 526 517 527 // Return (theta,phi) coordinates of middle of pixel with "NESTED" index k528 //--529 518 { 530 519 pix2ang_nest(nSide_,k,theta,phi); 531 520 } 532 521 533 //++ 522 523 /*! \fn int_4 SOPHYA::SphereHEALPix::NestToRing(int_4 k) const 524 525 translation from NESTED index into RING index 526 */ 534 527 template<class T> 535 528 int_4 SphereHEALPix<T>::NestToRing(int_4 k) const 536 529 537 // translation from NESTED index into RING index538 //539 //--540 530 { 541 531 return nest2ring(nSide_,k); 542 532 } 543 533 544 //++ 534 535 /*! \fn int_4 SOPHYA::SphereHEALPix::RingToNest(int_4 k) const 536 537 translation from RING index into NESTED index 538 */ 545 539 template<class T> 546 540 int_4 SphereHEALPix<T>::RingToNest(int_4 k) const 547 //548 // translation from RING index into NESTED index549 //550 //--551 541 { 552 542 return ring2nest(nSide_,k); … … 571 561 os << endl; 572 562 573 os << endl;574 //const PIXELS_XY& PXY= PIXELS_XY::instance();575 576 //os << endl; os << " contenu des tableaux conversions "<<endl;577 //for(int i=0; i < 5; i++)578 // {579 // os<<PXY.pix2x_(i)<<", "<<PXY.pix2y_(i)<<", "<<PXY.x2pix_(i)<<", "<<PXY.y2pix_(i)<<endl;580 // }581 os << endl;582 563 583 564 } -
trunk/SophyaLib/SkyMap/spherehealpix.h
r1196 r1217 16 16 // ***************** CLASSE SphereHEALPix ***************************** 17 17 18 //! class SphereHEALPix 19 /*! 20 Pixelisation Gorski 21 22 23 ----------------------------------------------------------------------- 24 version 0.8.2 Aug97 TAC Eric Hivon, Kris Gorski 25 ----------------------------------------------------------------------- 26 27 the sphere is split in 12 diamond-faces containing nside**2 pixels each 28 29 the numbering of the pixels (in the nested scheme) is similar to 30 quad-cube 31 In each face the first pixel is in the lowest corner of the diamond 32 33 the faces are (x,y) coordinate on each face 34 \verbatim 35 . . . . <--- North Pole 36 / \ / \ / \ / \ ^ ^ 37 . 0 . 1 . 2 . 3 . <--- z = 2/3 \ / 38 \ / \ / \ / \ / y \ / x 39 4 . 5 . 6 . 7 . 4 <--- equator \ / 40 / \ / \ / \ / \ \/ 41 . 8 . 9 .10 .11 . <--- z = -2/3 (0,0) : lowest corner 42 \ / \ / \ / \ / 43 . . . . <--- South Pole 44 \endverbatim 45 phi:0 2Pi 46 47 in the ring scheme pixels are numbered along the parallels 48 the first parallel is the one closest to the north pole and so on 49 on each parallel, pixels are numbered starting from the one closest 50 to phi = 0 51 52 nside MUST be a power of 2 (<= 8192) 53 54 */ 18 /*! Class SphereHEALPix */ 55 19 56 20 … … 81 45 82 46 SphereHEALPix(); 83 /*!84 m is the "nside" of the Gorski algorithm85 86 The total number of pixels will be Npix = 12*nside**287 88 nside MUST be a power of 2 (<= 8192)89 */90 47 SphereHEALPix(int_4 m); 91 48 SphereHEALPix(const SphereHEALPix<T>& s, bool share); 92 49 SphereHEALPix(const SphereHEALPix<T>& s); 93 //! Destructor94 50 virtual ~SphereHEALPix(); 95 51 96 // Temporaire?97 52 inline virtual bool IsTemp(void) const { 98 53 … … 101 56 return pixels_.IsTemp(); 102 57 } 103 /*! Setting blockdata to temporary (see ndatablock documentation) */ 58 104 59 inline virtual void SetTemp(bool temp=false) const 105 60 { … … 110 65 // ------------------ Definition of PixelMap abstract methods 111 66 112 /* Nombre de pixels du decoupage */113 /*! Return number of pixels of the splitting */114 67 virtual int_4 NbPixels() const; 115 68 116 /* Valeur du contenu du pixel d'indice "RING" k */117 /*! Return value of pixel with "RING" index k */118 69 virtual T& PixVal(int_4 k); 119 70 virtual T const& PixVal(int_4 k) const; 120 71 121 /* Nombre de tranches en theta */122 /*! Return number of slices in theta direction on the sphere */123 72 uint_4 NbThetaSlices() const; 124 /*! For a theta-slice with index 'index', return :125 126 the corresponding "theta"127 128 a vector containing the phi's of the pixels of the slice129 130 a vector containing the corresponding values of pixels131 */132 73 virtual void GetThetaSlice(int_4 index,r_8& theta,TVector<r_8>& phi,TVector<T>& value) const; 133 /*! For a theta-slice with index 'index', return :134 135 the corresponding "theta"136 137 the corresponding "phi" for first pixel of the slice138 139 a vector containing indices of the pixels of the slice140 141 (equally distributed in phi)142 143 a vector containing the corresponding values of pixels144 */145 74 virtual void GetThetaSlice(int_4 sliceIndex,r_8& theta, r_8& phi0, TVector<int_4>& pixelIndices,TVector<T>& value) const ; 146 75 147 /* Return true if teta,phi in map */148 76 virtual bool ContainsSph(double theta, double phi) const; 149 /* Indice "RING" du pixel vers lequel pointe une direction definie par150 ses coordonnees spheriques */151 /*! Return "RING" index of the pixel corresponding to direction (theta, phi).152 */153 77 virtual int_4 PixIndexSph(double theta,double phi) const; 154 78 155 /* Coordonnees spheriques du milieu du pixel d'indice "RING" k */156 79 virtual void PixThetaPhi(int_4 k,double& theta,double& phi) const; 157 80 158 /*! Set all pixels to value v */159 81 virtual T SetPixels(T v); 160 82 161 /* Pixel Solid angle (steradians) */ 162 /*! Pixel Solid angle (steradians) 83 /*! Pixel Solid angle (steradians) 163 84 164 85 All the pixels have the same solid angle. The dummy argument is … … 166 87 fulfil this requirement. 167 88 */ 168 virtual double PixSolAngle(int_4 dummy=0) const; 89 inline virtual double PixSolAngle(int_4 dummy=0) const {return omeg_;} 169 90 170 91 /* Acces to the DataBlock */ … … 174 95 // --------------- Specific methods 175 96 176 /*!177 m is the "nside" of the Gorski algorithm178 179 The total number of pixels will be Npix = 12*nside**2180 181 nside MUST be a power of 2 (<= 8192)182 */183 97 virtual void Resize(int_4 m); 184 98 185 // pour l'instant le tableau est ordonne selon RING, uniquement 99 /*! 100 101 \return type of storage of the map : RING or NESTED 102 103 at the moment, always RING 104 */ 186 105 inline virtual char* TypeOfMap() const {return "RING";}; 187 106 188 107 189 /* Valeur du contenu du pixel d'indice "NEST" k */190 /*! Return value of pixel with "NESTED" index k */191 108 virtual T& PixValNest(int_4 k); 192 /*! Return value of pixel with "NESTED" index k */193 109 virtual T const& PixValNest(int_4 k) const; 194 110 195 /* Indice "NEST" du pixel vers lequel pointe une direction definie par196 ses coordonnees spheriques */197 /*! Return "NESTED" index of the pixel corresponding to direction (theta, phi).198 */199 111 virtual int_4 PixIndexSphNest(double theta,double phi) const; 200 112 201 /* Coordonnees spheriques du milieu du pixel d'indice "NEST" k */202 /*! Return (theta,phi) coordinates of middle of pixel with "NESTED" index k203 */204 113 virtual void PixThetaPhiNest(int_4 k,double& theta,double& phi) const; 205 114 206 /* algorithme de pixelisation */207 115 void Pixelize(int_4); 208 116 209 /* convertit index nested en ring */210 /*! translation from NESTED index into RING index */211 117 int_4 NestToRing(int_4) const; 212 118 213 /* convertit index ring en nested" */214 /*! translation from RING index into NESTED index */215 119 int_4 RingToNest(int_4) const; 216 120 217 121 218 /* retourne la valeur du parametre Gorski*/122 /*! \return value of healpix nside */ 219 123 inline virtual int_4 SizeIndex() const {return(nSide_);} 220 124 221 /* impression */222 125 void print(ostream& os) const; 223 126 … … 256 159 257 160 NDataBlock<T> pixels_; 258 NDataBlock<int_4> sliceBeginIndex_; // Rationalisation Mac.D.Y.161 NDataBlock<int_4> sliceBeginIndex_; // Rationalisation Mac. D.Y. 259 162 NDataBlock<int_4> sliceLenght_; 260 163
Note:
See TracChangeset
for help on using the changeset viewer.