Line data Source code
1 : //----------------------------------------------------------------------
2 : // This file is part of PARSEC (http://physik.rwth-aachen.de/parsec)
3 : // a parametrized simulation engine for cosmic rays.
4 : //
5 : // Copyright (C) 2011 Martin Erdmann, Peter Schiffer, Tobias Winchen
6 : // RWTH Aachen University, Germany
7 : // Contact: winchen@physik.rwth-aachen.de
8 : //
9 : // This program is free software: you can redistribute it and/or
10 : // modify it under the terms of the GNU General Public License as
11 : // published by the Free Software Foundation, either version 3 of
12 : // the License, or (at your option) any later version.
13 : //
14 : // This program is distributed in the hope that it will be useful,
15 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : // GNU General Public License for more details.
18 : //
19 : // You should have received a copy of the GNU General Public License
20 : // along with this program. If not, see <http://www.gnu.org/licenses/>.
21 : //----------------------------------------------------------------------
22 :
23 :
24 : #ifndef PIXELIZATION_HH
25 : #define PIXELIZATION_HH
26 :
27 : #include "healpix_base/healpix_base.h"
28 : #include <cmath>
29 : #include <stdint.h>
30 :
31 : namespace crpropa
32 : {
33 : /**
34 : * \addtogroup MagneticLenses
35 : * @{
36 : */
37 :
38 : /// Helpers to makes work with Healpix smooth
39 : const uint8_t _nOrder_max = 13;
40 : const uint32_t _nPix[] =
41 : {
42 : 48,
43 : 192,
44 : 768,
45 : 3072,
46 : 12288,
47 : 49152,
48 : 196608,
49 : 786432,
50 : 3145728,
51 : 12582912,
52 : 50331648,
53 : 201326592,
54 : 805306368
55 : };
56 :
57 : /// Every communication with healpix is done through this class to avoid
58 : /// bugs with missmatching coordinates (and make python hooks easier)
59 : class Pixelization
60 : {
61 : public:
62 0 : Pixelization()
63 0 : {
64 0 : _healpix = new healpix::T_Healpix_Base<int>(6, healpix::RING);
65 0 : }
66 :
67 : /// Constructor creating Pixelization with healpix order 6 (about
68 : /// 50000 pixels)
69 10 : Pixelization(uint8_t order)
70 10 : {
71 10 : _healpix = new healpix::T_Healpix_Base<int>(order, healpix::RING);
72 10 : }
73 :
74 : ~Pixelization()
75 : {
76 10 : delete _healpix;
77 10 : }
78 :
79 : /// Returns the number of the pixel which includes the direction (phi,theta)
80 : /// phi in [-pi, pi], theta in [-pi/2, pi/2]
81 : uint32_t direction2Pix(double longitude, double latitude) const;
82 :
83 : /// Returns the number of pixels of the pixelization
84 : uint32_t nPix() const
85 : {
86 86026 : return _healpix->Npix();
87 : }
88 :
89 : /// Returns the number of pixels given by healpix order
90 : static uint32_t nPix(uint8_t order);
91 :
92 : /// Returns the number of pixels of the pixelization
93 : int getNumberOfPixels()
94 : {
95 10317737 : return _healpix->Npix();
96 : }
97 :
98 : /// Returns the order, a given pixel number corresponds to. 0 if no
99 : /// match!
100 : static uint8_t pix2Order(uint32_t pix);
101 :
102 : /// Gives the center of pixel i in longitude [rad] and latitude [rad]
103 : void pix2Direction(uint32_t i, double &longitude, double &latitude) const;
104 :
105 : /// Calculate the angle [rad] between the vectors pointing to pixels i and j
106 : double angularDistance(uint32_t i, uint32_t j) const;
107 :
108 : /// Returns the maximum possible pixelization order
109 : uint8_t getMaxOrder() const
110 : {
111 : return _nOrder_max;
112 : }
113 :
114 : /// Returns healpix order
115 : uint8_t getOrder() const
116 : {
117 0 : return _healpix->Order();
118 : }
119 :
120 : void getRandomDirectionInPixel(uint32_t pixel, double &longitude, double &latitude);
121 :
122 0 : void getPixelsInCone(double longitude, double latitude,double radius, std::vector<int>& listpix)
123 : {
124 : healpix::vec3 v;
125 0 : spherCo2Vec(longitude, latitude, v);
126 : healpix::pointing p;
127 0 : p.from_vec3(v);
128 0 : _healpix->query_disc(p, radius, listpix);
129 0 : }
130 :
131 : private:
132 : void spherCo2Vec(double phi, double theta, healpix::vec3 &V) const;
133 : void vec2SphereCo(double &phi , double &theta, const healpix::vec3 &V) const;
134 : healpix::T_Healpix_Base<int> *_healpix;
135 : static healpix::T_Healpix_Base<healpix::int64> _healpix_nest;
136 : };
137 :
138 :
139 : /** @}*/
140 : } // namespace
141 : #endif // PIXELIZATION_HH
|
