1 | /* ======================================================================== *\
|
---|
2 | !
|
---|
3 | ! *
|
---|
4 | ! * This file is part of MARS, the MAGIC Analysis and Reconstruction
|
---|
5 | ! * Software. It is distributed to you in the hope that it can be a useful
|
---|
6 | ! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
|
---|
7 | ! * It is distributed WITHOUT ANY WARRANTY.
|
---|
8 | ! *
|
---|
9 | ! * Permission to use, copy, modify and distribute this software and its
|
---|
10 | ! * documentation for any purpose is hereby granted without fee,
|
---|
11 | ! * provided that the above copyright notice appear in all copies and
|
---|
12 | ! * that both that copyright notice and this permission notice appear
|
---|
13 | ! * in supporting documentation. It is provided "as is" without express
|
---|
14 | ! * or implied warranty.
|
---|
15 | ! *
|
---|
16 | !
|
---|
17 | !
|
---|
18 | ! Author(s):
|
---|
19 | !
|
---|
20 | ! Copyright: MAGIC Software Development, 2000-2007
|
---|
21 | !
|
---|
22 | !
|
---|
23 | \* ======================================================================== */
|
---|
24 |
|
---|
25 | /////////////////////////////////////////////////////////////////////////////
|
---|
26 | //
|
---|
27 | // MMcEvt
|
---|
28 | //
|
---|
29 | // This class handles and contains the MonteCarlo information
|
---|
30 | // with which the events have been generated
|
---|
31 | // This information exists for each event.
|
---|
32 | //
|
---|
33 | // Note: All the azimuth Phi angles in this and other MC classes follow
|
---|
34 | // the convention in the Corsika program (see Corsika manual). There, phi
|
---|
35 | // is the azimuth of the momentum vector of particles, and is measured
|
---|
36 | // from the north direction, anticlockwise (i.e, west is phi=90 degrees).
|
---|
37 | // When it refers to the telescope orientation, it is the azimuth of a
|
---|
38 | // vector along the telescope axis, going from the camera to the mirror.
|
---|
39 | // So, fTelescopeTheta=90, fTelescopePhi = 0 means the telescope is
|
---|
40 | // pointing horizontally towards South. For an explanation, see also
|
---|
41 | // TDAS 02-11.
|
---|
42 | //
|
---|
43 | // Version 4:
|
---|
44 | // - Added member fFadcTimeJitter
|
---|
45 | //
|
---|
46 | // Version 5:
|
---|
47 | // - removed fPartId, fEnergy, fImpact, fTelescopeTheta, fTelescopePhi
|
---|
48 | // - derives now from MMcEvtBasic which contains all these values
|
---|
49 | // - moved ParticleId_t to base class MMcEvtBasic
|
---|
50 | //
|
---|
51 | // Version 6:
|
---|
52 | // - added fEventReuse
|
---|
53 | //
|
---|
54 | // Version 7:
|
---|
55 | // - moved fPhi to MMcEvtBasic
|
---|
56 | // - moved fTheta to MmcEvtBasic
|
---|
57 | //
|
---|
58 | /////////////////////////////////////////////////////////////////////////////
|
---|
59 | #include "MMcEvt.hxx"
|
---|
60 |
|
---|
61 | #include "fits.h"
|
---|
62 |
|
---|
63 | #include <TMath.h>
|
---|
64 |
|
---|
65 | #include "MLog.h"
|
---|
66 | #include "MLogManip.h"
|
---|
67 |
|
---|
68 | #include "MString.h"
|
---|
69 |
|
---|
70 | ClassImp(MMcEvt);
|
---|
71 |
|
---|
72 | using namespace std;
|
---|
73 |
|
---|
74 | // --------------------------------------------------------------------------
|
---|
75 | //
|
---|
76 | // Default constructor. Calls Clear()
|
---|
77 | //
|
---|
78 | MMcEvt::MMcEvt()
|
---|
79 | {
|
---|
80 | fName = "MMcEvt";
|
---|
81 | fTitle = "Event info from Monte Carlo";
|
---|
82 |
|
---|
83 | Clear();
|
---|
84 | }
|
---|
85 |
|
---|
86 | // --------------------------------------------------------------------------
|
---|
87 | //
|
---|
88 | // Constructor. Use this to set all data members
|
---|
89 | //
|
---|
90 | // THIS FUNCTION IS FOR THE SIMULATION OLNY.
|
---|
91 | // DON'T USE THIS MEMBERFUNCTION IN THE ANALYSIS.
|
---|
92 | //
|
---|
93 | MMcEvt::MMcEvt(UInt_t fEvtNum, ParticleId_t usPId, Float_t fEner,
|
---|
94 | Float_t fThi0, Float_t fFirTar, Float_t fzFirInt,
|
---|
95 | Float_t fThet, Float_t fPhii, Float_t fCorD,
|
---|
96 | Float_t fCorX, Float_t fCorY, Float_t fImpa,
|
---|
97 | Float_t fTPhii, Float_t fTThet, Float_t fTFirst,
|
---|
98 | Float_t fTLast, Float_t fL_Nmax, Float_t fL_t0,
|
---|
99 | Float_t fL_tmax, Float_t fL_a, Float_t fL_b,
|
---|
100 | Float_t fL_c, Float_t fL_chi2, UInt_t uiPin,
|
---|
101 | UInt_t uiPat, UInt_t uiPre, UInt_t uiPco,
|
---|
102 | UInt_t uiPelS, UInt_t uiPelC, Float_t elec,
|
---|
103 | Float_t muon, Float_t other, Float_t fadc_jitter,
|
---|
104 | Int_t reuse)
|
---|
105 | {
|
---|
106 | fName = "MMcEvt";
|
---|
107 | fTitle = "Event info from Monte Carlo";
|
---|
108 |
|
---|
109 | Fill(fEvtNum, usPId, fEner, fThi0, fFirTar, fzFirInt, fThet,
|
---|
110 | fPhii, fCorD, fCorX, fCorY, fImpa, fTPhii, fTThet, fTFirst,
|
---|
111 | fTLast, fL_Nmax, fL_t0, fL_tmax, fL_a, fL_b, fL_c, fL_chi2,
|
---|
112 | uiPin, uiPat, uiPre, uiPco, uiPelS, uiPelC, elec, muon, other,
|
---|
113 | fadc_jitter, reuse);
|
---|
114 | }
|
---|
115 |
|
---|
116 | // --------------------------------------------------------------------------
|
---|
117 | //
|
---|
118 | // reset all values to values as nonsense as possible
|
---|
119 | //
|
---|
120 | void MMcEvt::Clear(Option_t *opt)
|
---|
121 | {
|
---|
122 | fPartId = kUNDEFINED;
|
---|
123 | fEnergy = -1;
|
---|
124 |
|
---|
125 | fCoreD = 0;
|
---|
126 | fCoreX = 0;
|
---|
127 | fCoreY = 0;
|
---|
128 | fImpact = -1;
|
---|
129 |
|
---|
130 | fPhotIni = 0;
|
---|
131 | fPassPhotAtm = 0;
|
---|
132 | fPassPhotRef = 0;
|
---|
133 | fPassPhotCone = 0;
|
---|
134 | fPhotElfromShower = 0;
|
---|
135 | fPhotElinCamera = 0;
|
---|
136 |
|
---|
137 | fElecCphFraction=0;
|
---|
138 | fMuonCphFraction=0;
|
---|
139 | fOtherCphFraction=0;
|
---|
140 |
|
---|
141 | fEventReuse = 0;
|
---|
142 | }
|
---|
143 |
|
---|
144 | // --------------------------------------------------------------------------
|
---|
145 | //
|
---|
146 | // Use this to set all data members
|
---|
147 | //
|
---|
148 | // THIS FUNCTION IS FOR THE SIMULATION OLNY.
|
---|
149 | // DON'T USE THIS MEMBERFUNCTION IN THE ANALYSIS.
|
---|
150 | //
|
---|
151 | void MMcEvt::Fill( UInt_t fEvtNum, ParticleId_t usPId, Float_t fEner,
|
---|
152 | Float_t fThi0, Float_t fFirTar, Float_t fzFirInt,
|
---|
153 | Float_t fThet, Float_t fPhii, Float_t fCorD,
|
---|
154 | Float_t fCorX, Float_t fCorY, Float_t fImpa,
|
---|
155 | Float_t fTPhii, Float_t fTThet, Float_t fTFirst,
|
---|
156 | Float_t fTLast, Float_t fL_Nmax, Float_t fL_t0,
|
---|
157 | Float_t fL_tmax, Float_t fL_a, Float_t fL_b,
|
---|
158 | Float_t fL_c, Float_t fL_chi2, UInt_t uiPin,
|
---|
159 | UInt_t uiPat, UInt_t uiPre, UInt_t uiPco,
|
---|
160 | UInt_t uiPelS, UInt_t uiPelC, Float_t elec,
|
---|
161 | Float_t muon, Float_t other, Float_t fadc_jitter,
|
---|
162 | Int_t reuse)
|
---|
163 | {
|
---|
164 | //
|
---|
165 | // All datamembers are filled with the correspondin parameters.
|
---|
166 | //
|
---|
167 | // Don't use this memberfunction in analysis
|
---|
168 | //
|
---|
169 | fEvtNumber = fEvtNum;
|
---|
170 | fPartId = usPId ;
|
---|
171 | fEnergy = fEner ;
|
---|
172 | fThick0 = fThi0;
|
---|
173 | fFirstTarget = fFirTar;
|
---|
174 | fZFirstInteraction = fzFirInt;
|
---|
175 |
|
---|
176 | fTheta = fThet ;
|
---|
177 | fPhi = fPhii ;
|
---|
178 |
|
---|
179 | fCoreD = fCorD ;
|
---|
180 | fCoreX = fCorX ;
|
---|
181 | fCoreY = fCorY ;
|
---|
182 | fImpact = fImpa ;
|
---|
183 |
|
---|
184 | fTelescopePhi = fTPhii;
|
---|
185 | fTelescopeTheta = fTThet;
|
---|
186 | fTimeFirst = fTFirst;
|
---|
187 | fTimeLast = fTLast;
|
---|
188 | fLongiNmax = fL_Nmax;
|
---|
189 | fLongit0 = fL_t0;
|
---|
190 | fLongitmax = fL_tmax;
|
---|
191 | fLongia = fL_a;
|
---|
192 | fLongib = fL_b;
|
---|
193 | fLongic = fL_c;
|
---|
194 | fLongichi2 = fL_chi2;
|
---|
195 |
|
---|
196 | fPhotIni = uiPin ;
|
---|
197 | fPassPhotAtm = fPhotIni-uiPat ;
|
---|
198 | fPassPhotRef = fPassPhotAtm-uiPre ;
|
---|
199 | fPassPhotCone = uiPco ;
|
---|
200 | fPhotElfromShower = uiPelS ;
|
---|
201 | fPhotElinCamera = uiPelC ;
|
---|
202 |
|
---|
203 | fElecCphFraction=elec;
|
---|
204 | fMuonCphFraction=muon;
|
---|
205 | fOtherCphFraction=other;
|
---|
206 |
|
---|
207 | fFadcTimeJitter = fadc_jitter;
|
---|
208 |
|
---|
209 | fEventReuse = reuse;
|
---|
210 | }
|
---|
211 |
|
---|
212 | // --------------------------------------------------------------------------
|
---|
213 | //
|
---|
214 | // Print the contents of the container.
|
---|
215 | //
|
---|
216 | // if you specify an option only the requested data members are printed:
|
---|
217 | // allowed options are:
|
---|
218 | // id, energy, impact
|
---|
219 | //
|
---|
220 | void MMcEvt::Print(Option_t *opt) const
|
---|
221 | {
|
---|
222 | MMcEvtBasic::Print(opt);
|
---|
223 |
|
---|
224 | TString str(opt);
|
---|
225 | if (str.IsNull())
|
---|
226 | *fLog << " Photoelectrons: " << fPhotElfromShower << endl;
|
---|
227 | }
|
---|
228 |
|
---|
229 | // --------------------------------------------------------------------------
|
---|
230 | //
|
---|
231 | // Return a proper description of the monte carlo event
|
---|
232 | //
|
---|
233 | TString MMcEvt::GetDescription(const TString &s) const
|
---|
234 | {
|
---|
235 | TString txt("#splitline{");
|
---|
236 |
|
---|
237 | txt += GetParticleName();
|
---|
238 | txt += " ";
|
---|
239 | txt += s;
|
---|
240 | txt += "}{ E=";
|
---|
241 | txt += GetEnergyStr();
|
---|
242 | txt += " r=";
|
---|
243 | txt += TMath::Nint(GetImpact()/100);
|
---|
244 | txt += "m Zd=";
|
---|
245 | txt += MString::Format("%.1f", GetTelescopeTheta()*TMath::RadToDeg());
|
---|
246 | txt += "\\circ ";
|
---|
247 | if (GetPhotElfromShower()>=10000)
|
---|
248 | txt += MString::Format("%dk", TMath::Nint(GetPhotElfromShower()/1000.));
|
---|
249 | else
|
---|
250 | if (GetPhotElfromShower()>=1000)
|
---|
251 | txt += MString::Format("%.1fk", GetPhotElfromShower()/1000.);
|
---|
252 | else
|
---|
253 | txt += GetPhotElfromShower();
|
---|
254 | txt += "PhEl}";
|
---|
255 |
|
---|
256 | return txt;
|
---|
257 | }
|
---|
258 |
|
---|
259 | Bool_t MMcEvt::SetupFits(fits &fin)
|
---|
260 | {
|
---|
261 | if (fName!="MMcEvt")
|
---|
262 | {
|
---|
263 | *fLog << err << "ERROR - SetupFits only supported if name equals MMcEvt." << endl;
|
---|
264 | return kFALSE;
|
---|
265 | }
|
---|
266 |
|
---|
267 | if (!fin.SetRefAddress("MMcEvt.fEvtNumber", fEvtNumber )) return kFALSE;
|
---|
268 | if (!fin.SetRefAddress("MMcEvt.fThick0", fThick0 )) return kFALSE;
|
---|
269 | if (!fin.SetRefAddress("MMcEvt.fFirstTarget", fFirstTarget )) return kFALSE;
|
---|
270 | if (!fin.SetRefAddress("MMcEvt.fZFirstInteraction", fZFirstInteraction)) return kFALSE;
|
---|
271 | if (!fin.SetRefAddress("MMcEvt.fCoreD", fCoreD )) return kFALSE;
|
---|
272 | if (!fin.SetRefAddress("MMcEvt.fCoreX", fCoreX )) return kFALSE;
|
---|
273 | if (!fin.SetRefAddress("MMcEvt.fCoreY", fCoreY )) return kFALSE;
|
---|
274 | if (!fin.SetRefAddress("MMcEvt.fTimeFirst", fTimeFirst )) return kFALSE;
|
---|
275 | if (!fin.SetRefAddress("MMcEvt.fTimeLast", fTimeLast )) return kFALSE;
|
---|
276 | if (!fin.SetRefAddress("MMcEvt.fLongiNmax", fLongiNmax )) return kFALSE;
|
---|
277 | if (!fin.SetRefAddress("MMcEvt.fLongit0", fLongit0 )) return kFALSE;
|
---|
278 | if (!fin.SetRefAddress("MMcEvt.fLongitmax", fLongitmax )) return kFALSE;
|
---|
279 | if (!fin.SetRefAddress("MMcEvt.fLongia", fLongia )) return kFALSE;
|
---|
280 | if (!fin.SetRefAddress("MMcEvt.fLongib", fLongib )) return kFALSE;
|
---|
281 | if (!fin.SetRefAddress("MMcEvt.fLongic", fLongic )) return kFALSE;
|
---|
282 | if (!fin.SetRefAddress("MMcEvt.fLongichi2", fLongichi2 )) return kFALSE;
|
---|
283 | if (!fin.SetRefAddress("MMcEvt.fPhotIni", fPhotIni )) return kFALSE;
|
---|
284 | if (!fin.SetRefAddress("MMcEvt.fPassPhotAtm", fPassPhotAtm )) return kFALSE;
|
---|
285 | if (!fin.SetRefAddress("MMcEvt.fPassPhotRef", fPassPhotRef )) return kFALSE;
|
---|
286 | if (!fin.SetRefAddress("MMcEvt.fPassPhotCone", fPassPhotCone )) return kFALSE;
|
---|
287 | if (!fin.SetRefAddress("MMcEvt.fPhotElfromShower", fPhotElfromShower )) return kFALSE;
|
---|
288 | if (!fin.SetRefAddress("MMcEvt.fPhotElinCamera", fPhotElinCamera )) return kFALSE;
|
---|
289 | if (!fin.SetRefAddress("MMcEvt.fElecCphFraction", fElecCphFraction )) return kFALSE;
|
---|
290 | if (!fin.SetRefAddress("MMcEvt.fMuonCphFraction", fMuonCphFraction )) return kFALSE;
|
---|
291 | if (!fin.SetRefAddress("MMcEvt.fOtherCphFraction", fOtherCphFraction )) return kFALSE;
|
---|
292 | if (!fin.SetRefAddress("MMcEvt.fFadcTimeJitter", fFadcTimeJitter )) return kFALSE;
|
---|
293 | if (!fin.SetRefAddress("MMcEvt.fEventReuse", fEventReuse )) return kFALSE;
|
---|
294 |
|
---|
295 | return MMcEvtBasic::SetupFits(fin);
|
---|
296 | }
|
---|