1 | /////////////////////////////////////////////////////////////////
|
---|
2 | //
|
---|
3 | // MCEventHeader_2
|
---|
4 | //
|
---|
5 | // Created: Sat Oct 19 12:11:43 CEST 2002
|
---|
6 | // Author: Oscar Blanch Bigas
|
---|
7 | // Purpose: Base class for reflector 0.6 EventHeader
|
---|
8 | // Notes:
|
---|
9 | //
|
---|
10 | /////////////////////////////////////////////////////////////////
|
---|
11 |
|
---|
12 | // @T \newpage
|
---|
13 |
|
---|
14 | // @section Source code of {\tt MCEventHeader_2.hxx}
|
---|
15 |
|
---|
16 | /* @text
|
---|
17 | This section shows the include file {\tt MCEventHeader_2.hxx}
|
---|
18 | @endtext */
|
---|
19 |
|
---|
20 | #ifndef MCEventHeader_2_Class
|
---|
21 | #define MCEventHeader_2_Class
|
---|
22 |
|
---|
23 | // @subsection Include files
|
---|
24 |
|
---|
25 | // @code
|
---|
26 | #include "TROOT.h"
|
---|
27 | #include "TObject.h"
|
---|
28 |
|
---|
29 | #include <iostream>
|
---|
30 | #include <iomanip>
|
---|
31 | #include <fstream>
|
---|
32 | #include <stdlib.h>
|
---|
33 | #include <math.h>
|
---|
34 |
|
---|
35 | #include "COREventHeader.hxx"
|
---|
36 |
|
---|
37 | #define SIZE_OF_MCEVENTHEADER_2 181 /* floats */
|
---|
38 |
|
---|
39 | // @endcode
|
---|
40 |
|
---|
41 | // @subsection Class {\em MCEventHeader_2}: Definition
|
---|
42 |
|
---|
43 | // @code
|
---|
44 | class MCEventHeader_2 {
|
---|
45 |
|
---|
46 | protected:
|
---|
47 | Float_t EvtNumber;
|
---|
48 | Float_t PrimaryID;
|
---|
49 | Float_t Etotal;
|
---|
50 | Float_t Thick0;
|
---|
51 | Float_t FirstTarget;
|
---|
52 | Float_t zFirstInt;
|
---|
53 | Float_t p[3];
|
---|
54 | Float_t Theta;
|
---|
55 | Float_t Phi;
|
---|
56 |
|
---|
57 | Float_t NumRndSeq;
|
---|
58 | Float_t RndData[10][3];
|
---|
59 |
|
---|
60 | Float_t RunNumber;
|
---|
61 | Float_t DateRun;
|
---|
62 | Float_t VersionPGM;
|
---|
63 |
|
---|
64 | Float_t NumObsLev;
|
---|
65 | Float_t HeightLev[10];
|
---|
66 |
|
---|
67 | Float_t SlopeSpec;
|
---|
68 | Float_t ELowLim;
|
---|
69 | Float_t EUppLim;
|
---|
70 |
|
---|
71 | Float_t Ecutoffh;
|
---|
72 | Float_t Ecutoffm;
|
---|
73 | Float_t Ecutoffe;
|
---|
74 | Float_t Ecutoffg;
|
---|
75 |
|
---|
76 | Float_t NFLAIN;
|
---|
77 | Float_t NFLDIF;
|
---|
78 | Float_t NFLPI0;
|
---|
79 | Float_t NFLPIF;
|
---|
80 | Float_t NFLCHE;
|
---|
81 | Float_t NFRAGM;
|
---|
82 |
|
---|
83 | Float_t Bx;
|
---|
84 | Float_t By;
|
---|
85 |
|
---|
86 | Float_t EGS4yn;
|
---|
87 | Float_t NKGyn;
|
---|
88 | Float_t GHEISHAyn;
|
---|
89 | Float_t VENUSyn;
|
---|
90 | Float_t CERENKOVyn;
|
---|
91 | Float_t NEUTRINOyn;
|
---|
92 | Float_t HORIZONTyn;
|
---|
93 | Float_t COMPUTER;
|
---|
94 |
|
---|
95 | Float_t ThetaMin;
|
---|
96 | Float_t ThetaMax;
|
---|
97 | Float_t PhiMin;
|
---|
98 | Float_t PhiMax;
|
---|
99 |
|
---|
100 | Float_t CBunchSize;
|
---|
101 | Float_t CDetInX,CDetInY;
|
---|
102 | Float_t CSpacInX,CSpacInY;
|
---|
103 | Float_t CLenInX,CLenInY;
|
---|
104 | Float_t COutput;
|
---|
105 |
|
---|
106 | Float_t AngleNorthX;
|
---|
107 | Float_t MuonInfo;
|
---|
108 |
|
---|
109 | Float_t StepLength;
|
---|
110 | Float_t CWaveLower;
|
---|
111 | Float_t CWaveUpper;
|
---|
112 | Float_t Multipl;
|
---|
113 | Float_t CorePos[2][20];
|
---|
114 | Float_t SIBYLL[2];
|
---|
115 | Float_t QGSJET[2];
|
---|
116 | Float_t DPMJET[2];
|
---|
117 | Float_t VENUS_cross;
|
---|
118 | Float_t mu_mult_scat;
|
---|
119 | Float_t NKG_range;
|
---|
120 | Float_t EFRCTHN[2];
|
---|
121 | Float_t WMAX[2];
|
---|
122 | Float_t rthin_rmax;
|
---|
123 |
|
---|
124 | Float_t viewcone_angles[2];
|
---|
125 | /* (degrees) Inner and outer angles in
|
---|
126 | * Corsika's VIEWCONE option. This is
|
---|
127 | * only possible with Corsika>6 versions. In
|
---|
128 | * that case, PhiMin=PhiMax and
|
---|
129 | * ThetaMin=ThetaMax (also in this header)
|
---|
130 | * indicate the axis of this cone.
|
---|
131 | */
|
---|
132 |
|
---|
133 | /* ^^^ Up to here, the info from the CORSIKA event header. */
|
---|
134 |
|
---|
135 | /* Telescope orientation: */
|
---|
136 | Float_t telescopePhi; /* rad */
|
---|
137 | Float_t telescopeTheta; /* rad */
|
---|
138 |
|
---|
139 | /* Time of first and last photon: */
|
---|
140 | Float_t TimeFirst;
|
---|
141 | Float_t TimeLast;
|
---|
142 |
|
---|
143 | /* 6 parameters and chi2 of the NKG fit to the longitudinal
|
---|
144 | * particle distribution (see CORSIKA manual for explanation):
|
---|
145 | */
|
---|
146 | Float_t longi_Nmax;
|
---|
147 | Float_t longi_t0;
|
---|
148 | Float_t longi_tmax;
|
---|
149 | Float_t longi_a;
|
---|
150 | Float_t longi_b;
|
---|
151 | Float_t longi_c;
|
---|
152 | Float_t longi_chi2;
|
---|
153 |
|
---|
154 | Float_t CORSIKAPhs; // Original photons written by Corsika
|
---|
155 | Float_t AtmAbsPhs; // Photons absorved by the atmosphere
|
---|
156 | Float_t MirrAbsPhs; // Photons absorved by the mirrors
|
---|
157 | Float_t OutOfMirrPhs; // Photons outside the mirror
|
---|
158 | Float_t BlackSpotPhs; // Photons lost in the "black spot"
|
---|
159 | Float_t OutOfChamPhs; // Photons outside the chamber
|
---|
160 | Float_t CPhotons; // Photons reaching the chamber
|
---|
161 |
|
---|
162 | /* Now follow the fraction of photons reaching the camera produced by
|
---|
163 | * electrons, muons and other particles respectively:
|
---|
164 | */
|
---|
165 | Float_t elec_cph_fraction;
|
---|
166 | Float_t muon_cph_fraction;
|
---|
167 | Float_t other_cph_fraction;
|
---|
168 |
|
---|
169 | Float_t dummy[7]; /* not used */
|
---|
170 |
|
---|
171 | public:
|
---|
172 | MCEventHeader_2(void) {} // default constructor
|
---|
173 |
|
---|
174 | virtual ~MCEventHeader_2(void) {} // default destructor
|
---|
175 |
|
---|
176 | // reads EventHeader from binary input stream
|
---|
177 | Int_t read ( ifstream &is ) {
|
---|
178 | is.read ( (char *)this, mysize() );
|
---|
179 | return is.gcount();
|
---|
180 | }
|
---|
181 |
|
---|
182 | Int_t read (FILE* f) {
|
---|
183 | return fread((char*)&EvtNumber, (SIZE_OF_MCEVENTHEADER_2)*sizeof(Float_t), 1, f );
|
---|
184 | }
|
---|
185 |
|
---|
186 | // writes EventHeader to binary output stream
|
---|
187 | Int_t write ( ofstream &os ) {
|
---|
188 | os.write ( (char *)this, mysize() );
|
---|
189 | return 0;
|
---|
190 | }
|
---|
191 |
|
---|
192 | // get information about the EventHeader
|
---|
193 |
|
---|
194 | // get the event number
|
---|
195 | inline Float_t get_evt_number (void) {return EvtNumber;}
|
---|
196 |
|
---|
197 | // get the primary type (GEANT code)
|
---|
198 | inline Float_t get_primary ( void ) { return ( PrimaryID ); }
|
---|
199 |
|
---|
200 | // get the total primary energy
|
---|
201 | inline Float_t get_energy ( void ) { return ( Etotal ); }
|
---|
202 |
|
---|
203 | // get Thick0
|
---|
204 | inline Float_t get_thick0 (void) { return Thick0;}
|
---|
205 |
|
---|
206 | // get First Target
|
---|
207 | inline Float_t get_first_target (void) {return FirstTarget;}
|
---|
208 |
|
---|
209 | // get z first interaction
|
---|
210 | inline Float_t get_z_first_int (void) {return zFirstInt;}
|
---|
211 |
|
---|
212 | // get the initial zenith angle
|
---|
213 | inline Float_t get_theta ( void ) { return ( Theta ); }
|
---|
214 |
|
---|
215 | // get the initial phi angle
|
---|
216 | inline Float_t get_phi ( void ) { return ( Phi ); }
|
---|
217 |
|
---|
218 | // get the telescope zenith angle
|
---|
219 | inline Float_t get_theta_CT ( void ) { return ( telescopeTheta ); }
|
---|
220 |
|
---|
221 | // get the telescope phi angle
|
---|
222 | inline Float_t get_phi_CT ( void ) { return ( telescopePhi ); }
|
---|
223 |
|
---|
224 | // get the spectrum slope
|
---|
225 | inline Float_t get_slope ( void ) { return ( SlopeSpec ); }
|
---|
226 |
|
---|
227 | // get height of first interaction (in cm)
|
---|
228 | inline Float_t get_height ( void ) { return ( zFirstInt ); }
|
---|
229 |
|
---|
230 | // get the energy range of this run
|
---|
231 | inline void get_energy_range ( Float_t *elow, Float_t *eup ) {
|
---|
232 | *elow = ELowLim;
|
---|
233 | *eup = EUppLim;
|
---|
234 | }
|
---|
235 |
|
---|
236 | // get the theta range of this run
|
---|
237 | inline void get_theta_range ( Float_t *thlow, Float_t *thup ) {
|
---|
238 | *thlow = ThetaMin;
|
---|
239 | *thup = ThetaMax;
|
---|
240 | }
|
---|
241 |
|
---|
242 | inline void get_viewcone ( Float_t *rmin, Float_t *rmax){
|
---|
243 | *rmin=viewcone_angles[0];
|
---|
244 | *rmax=viewcone_angles[1];
|
---|
245 | }
|
---|
246 |
|
---|
247 | // get the energy range of this run
|
---|
248 | inline void get_phi_range ( Float_t *plow, Float_t *pup ) {
|
---|
249 | *plow = PhiMin;
|
---|
250 | *pup = PhiMax;
|
---|
251 | }
|
---|
252 |
|
---|
253 | // get the core position
|
---|
254 | inline Float_t get_core ( Float_t *x, Float_t *y, Int_t ncore = 0 ) {
|
---|
255 | *x = CorePos[0][ncore];
|
---|
256 | *y = CorePos[1][ncore];
|
---|
257 | return ( (Float_t) sqrt ( (*x)*(*x) + (*y)*(*y) ) );
|
---|
258 | }
|
---|
259 |
|
---|
260 | // get the core position
|
---|
261 | inline Float_t get_core ( Int_t ncore = 0 ) {
|
---|
262 | return ( (Float_t) sqrt ( ((CorePos[0][ncore])*(CorePos[0][ncore]))+
|
---|
263 | ((CorePos[1][ncore])*(CorePos[1][ncore])) ) );
|
---|
264 | }
|
---|
265 |
|
---|
266 | // get the core position in X
|
---|
267 | inline Float_t get_coreX ( Int_t ncore = 0 ) {
|
---|
268 | return ( (Float_t) CorePos[0][ncore] );
|
---|
269 | }
|
---|
270 |
|
---|
271 | // get the core position in Y
|
---|
272 | inline Float_t get_coreY ( Int_t ncore = 0 ) {
|
---|
273 | return ( (Float_t) CorePos[1][ncore] );
|
---|
274 | }
|
---|
275 |
|
---|
276 | //get photons at different stages
|
---|
277 | inline Float_t get_CORSIKA (){ return CORSIKAPhs;}
|
---|
278 |
|
---|
279 | inline Float_t get_AtmAbs (){ return AtmAbsPhs;}
|
---|
280 |
|
---|
281 | inline Float_t get_MirrAbs (){ return MirrAbsPhs;}
|
---|
282 |
|
---|
283 | inline Float_t get_OutOfMirr (){ return OutOfMirrPhs;}
|
---|
284 |
|
---|
285 | inline Float_t get_BlackSpot (){ return BlackSpotPhs;}
|
---|
286 |
|
---|
287 | inline Float_t get_OutOfCham (){ return OutOfChamPhs;}
|
---|
288 |
|
---|
289 | inline Float_t get_CPhotons (){ return CPhotons;}
|
---|
290 |
|
---|
291 | inline Float_t get_NumObsLev (){return NumObsLev;}
|
---|
292 |
|
---|
293 | inline Float_t get_HeightLev (Int_t i) {return HeightLev[i];}
|
---|
294 |
|
---|
295 | inline Float_t get_VersionPGM () {return VersionPGM;}
|
---|
296 |
|
---|
297 | inline Float_t get_RunNumber () {return RunNumber;}
|
---|
298 |
|
---|
299 | inline Float_t get_DateRun () {return DateRun;}
|
---|
300 |
|
---|
301 | inline Float_t get_NumRndSeq () {return NumRndSeq;}
|
---|
302 |
|
---|
303 | inline Float_t get_CWaveLower() {return CWaveLower;}
|
---|
304 |
|
---|
305 | inline Float_t get_CWaveUpper() {return CWaveUpper;}
|
---|
306 |
|
---|
307 | inline Float_t get_ElecFraction() {return elec_cph_fraction;}
|
---|
308 |
|
---|
309 | inline Float_t get_MuonFraction() {return muon_cph_fraction;}
|
---|
310 |
|
---|
311 | inline Float_t get_OtherFraction() {return other_cph_fraction;}
|
---|
312 |
|
---|
313 | // transport from COREventHeader to MCEventHeader_2
|
---|
314 | void transport ( COREventHeader *e );
|
---|
315 |
|
---|
316 | // write extreme times
|
---|
317 | inline void put_times ( Float_t t1, Float_t t2 ) {
|
---|
318 | TimeFirst = t1;
|
---|
319 | TimeLast = t2;
|
---|
320 | }
|
---|
321 |
|
---|
322 | // get extreme times
|
---|
323 | inline Float_t get_times ( Float_t *t1, Float_t *t2 ) {
|
---|
324 | *t1 = TimeFirst;
|
---|
325 | *t2 = TimeLast;
|
---|
326 | return ( TimeLast - TimeFirst );
|
---|
327 | }
|
---|
328 |
|
---|
329 | inline Float_t get_NKGfit(Float_t *Nmax, Float_t *t0, Float_t *tmax,
|
---|
330 | Float_t *a, Float_t *b, Float_t *c){
|
---|
331 | *Nmax=longi_Nmax; *t0=longi_t0; *tmax=longi_tmax;
|
---|
332 | *a=longi_a; *b=longi_b; *c=longi_c;
|
---|
333 |
|
---|
334 | return longi_chi2;
|
---|
335 |
|
---|
336 | }
|
---|
337 |
|
---|
338 | inline int mysize(void)
|
---|
339 | { return ( sizeof( float ) * SIZE_OF_MCEVENTHEADER_2 ); }
|
---|
340 |
|
---|
341 | // get deviations of the CT from the original shower direction
|
---|
342 | inline Float_t get_deviations ( Float_t *t1, Float_t *t2 ) {
|
---|
343 |
|
---|
344 | float ct1,st1,cp1,sp1;
|
---|
345 | float ct2,st2,cp2,sp2;
|
---|
346 | float x1,y1,z1;
|
---|
347 | float x2,y2,z2;
|
---|
348 |
|
---|
349 | ct1 = cos(Theta);
|
---|
350 | st1 = sin(Theta);
|
---|
351 | cp1 = cos(Phi);
|
---|
352 | sp1 = sin(Phi);
|
---|
353 |
|
---|
354 | ct2 = cos(telescopeTheta);
|
---|
355 | st2 = sin(telescopeTheta);
|
---|
356 | cp2 = cos(telescopePhi);
|
---|
357 | sp2 = sin(telescopePhi);
|
---|
358 |
|
---|
359 | x1 = st1*cp1; y1 = st1*sp1; z1 = ct1;
|
---|
360 | x2 = st2*cp2; y2 = st2*sp2; z2 = ct2;
|
---|
361 |
|
---|
362 | *t1 = (telescopeTheta-Theta);
|
---|
363 | *t2 = (telescopePhi-Phi);
|
---|
364 |
|
---|
365 | return ( acos(x1*x2 + y1*y2 + z1*z2) );
|
---|
366 | }
|
---|
367 |
|
---|
368 | inline void print ( void ) {
|
---|
369 | float *ptr = (float *)this;
|
---|
370 | for(int i=0; i<SIZE_OF_MCEVENTHEADER_2; ++i,++ptr)
|
---|
371 | cerr << i << ':' << *ptr << '\n';
|
---|
372 | }
|
---|
373 |
|
---|
374 | };
|
---|
375 | // @endcode
|
---|
376 |
|
---|
377 | #endif // not defined MCEventHeader_2_Class
|
---|