source: trunk/Mars/mmc/MMcEvt.cxx@ 19956

Last change on this file since 19956 was 19955, checked in by tbretz, 5 years ago
Reset fRunNumber in Clear for deft constructor and allow fRunNumber to be missing in FITS file.
File size: 10.7 KB
Line 
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// Version 8:
59// - added fRunNumber
60//
61/////////////////////////////////////////////////////////////////////////////
62#include "MMcEvt.hxx"
63
64#include "fits.h"
65
66#include <TMath.h>
67
68#include "MLog.h"
69#include "MLogManip.h"
70
71#include "MString.h"
72
73ClassImp(MMcEvt);
74
75using namespace std;
76
77// --------------------------------------------------------------------------
78//
79// Default constructor. Calls Clear()
80//
81MMcEvt::MMcEvt()
82{
83 fName = "MMcEvt";
84 fTitle = "Event info from Monte Carlo";
85
86 Clear();
87}
88
89// --------------------------------------------------------------------------
90//
91// Constructor. Use this to set all data members
92//
93// THIS FUNCTION IS FOR THE SIMULATION OLNY.
94// DON'T USE THIS MEMBERFUNCTION IN THE ANALYSIS.
95//
96MMcEvt::MMcEvt(UInt_t fEvtNum, ParticleId_t usPId, Float_t fEner,
97 Float_t fThi0, Float_t fFirTar, Float_t fzFirInt,
98 Float_t fThet, Float_t fPhii, Float_t fCorD,
99 Float_t fCorX, Float_t fCorY, Float_t fImpa,
100 Float_t fTPhii, Float_t fTThet, Float_t fTFirst,
101 Float_t fTLast, Float_t fL_Nmax, Float_t fL_t0,
102 Float_t fL_tmax, Float_t fL_a, Float_t fL_b,
103 Float_t fL_c, Float_t fL_chi2, UInt_t uiPin,
104 UInt_t uiPat, UInt_t uiPre, UInt_t uiPco,
105 UInt_t uiPelS, UInt_t uiPelC, Float_t elec,
106 Float_t muon, Float_t other, Float_t fadc_jitter,
107 Int_t reuse) : fRunNumber(0)
108{
109 fName = "MMcEvt";
110 fTitle = "Event info from Monte Carlo";
111
112 Fill(fEvtNum, usPId, fEner, fThi0, fFirTar, fzFirInt, fThet,
113 fPhii, fCorD, fCorX, fCorY, fImpa, fTPhii, fTThet, fTFirst,
114 fTLast, fL_Nmax, fL_t0, fL_tmax, fL_a, fL_b, fL_c, fL_chi2,
115 uiPin, uiPat, uiPre, uiPco, uiPelS, uiPelC, elec, muon, other,
116 fadc_jitter, reuse);
117}
118
119// --------------------------------------------------------------------------
120//
121// reset all values to values as nonsense as possible
122//
123void MMcEvt::Clear(Option_t *opt)
124{
125 fRunNumber = 0;
126
127 fPartId = kUNDEFINED;
128 fEnergy = -1;
129
130 fCoreD = 0;
131 fCoreX = 0;
132 fCoreY = 0;
133 fImpact = -1;
134
135 fPhotIni = 0;
136 fPassPhotAtm = 0;
137 fPassPhotRef = 0;
138 fPassPhotCone = 0;
139 fPhotElfromShower = 0;
140 fPhotElinCamera = 0;
141
142 fElecCphFraction=0;
143 fMuonCphFraction=0;
144 fOtherCphFraction=0;
145
146 fEventReuse = 0;
147}
148
149// --------------------------------------------------------------------------
150//
151// Use this to set all data members
152//
153// THIS FUNCTION IS FOR THE SIMULATION OLNY.
154// DON'T USE THIS MEMBERFUNCTION IN THE ANALYSIS.
155//
156void MMcEvt::Fill( UInt_t fEvtNum, ParticleId_t usPId, Float_t fEner,
157 Float_t fThi0, Float_t fFirTar, Float_t fzFirInt,
158 Float_t fThet, Float_t fPhii, Float_t fCorD,
159 Float_t fCorX, Float_t fCorY, Float_t fImpa,
160 Float_t fTPhii, Float_t fTThet, Float_t fTFirst,
161 Float_t fTLast, Float_t fL_Nmax, Float_t fL_t0,
162 Float_t fL_tmax, Float_t fL_a, Float_t fL_b,
163 Float_t fL_c, Float_t fL_chi2, UInt_t uiPin,
164 UInt_t uiPat, UInt_t uiPre, UInt_t uiPco,
165 UInt_t uiPelS, UInt_t uiPelC, Float_t elec,
166 Float_t muon, Float_t other, Float_t fadc_jitter,
167 Int_t reuse)
168{
169 //
170 // All datamembers are filled with the correspondin parameters.
171 //
172 // Don't use this memberfunction in analysis
173 //
174 fEvtNumber = fEvtNum;
175 fPartId = usPId ;
176 fEnergy = fEner ;
177 fThick0 = fThi0;
178 fFirstTarget = fFirTar;
179 fZFirstInteraction = fzFirInt;
180
181 fTheta = fThet ;
182 fPhi = fPhii ;
183
184 fCoreD = fCorD ;
185 fCoreX = fCorX ;
186 fCoreY = fCorY ;
187 fImpact = fImpa ;
188
189 fTelescopePhi = fTPhii;
190 fTelescopeTheta = fTThet;
191 fTimeFirst = fTFirst;
192 fTimeLast = fTLast;
193 fLongiNmax = fL_Nmax;
194 fLongit0 = fL_t0;
195 fLongitmax = fL_tmax;
196 fLongia = fL_a;
197 fLongib = fL_b;
198 fLongic = fL_c;
199 fLongichi2 = fL_chi2;
200
201 fPhotIni = uiPin ;
202 fPassPhotAtm = fPhotIni-uiPat ;
203 fPassPhotRef = fPassPhotAtm-uiPre ;
204 fPassPhotCone = uiPco ;
205 fPhotElfromShower = uiPelS ;
206 fPhotElinCamera = uiPelC ;
207
208 fElecCphFraction=elec;
209 fMuonCphFraction=muon;
210 fOtherCphFraction=other;
211
212 fFadcTimeJitter = fadc_jitter;
213
214 fEventReuse = reuse;
215}
216
217// --------------------------------------------------------------------------
218//
219// Print the contents of the container.
220//
221// if you specify an option only the requested data members are printed:
222// allowed options are:
223// id, energy, impact
224//
225void MMcEvt::Print(Option_t *opt) const
226{
227 MMcEvtBasic::Print(opt);
228
229 TString str(opt);
230 if (str.IsNull())
231 {
232 *fLog << " ID: Run=" << fRunNumber << ", Evt=" << fEvtNumber << endl;
233 *fLog << " Photoelectrons: " << fPhotElfromShower << endl;
234 }
235}
236
237// --------------------------------------------------------------------------
238//
239// Return a proper description of the monte carlo event
240//
241TString MMcEvt::GetDescription(const TString &s) const
242{
243 TString txt("#splitline{");
244
245 txt += GetParticleName();
246 txt += " ";
247 txt += s;
248 txt += "}{ E=";
249 txt += GetEnergyStr();
250 txt += " r=";
251 txt += TMath::Nint(GetImpact()/100);
252 txt += "m Zd=";
253 txt += MString::Format("%.1f", GetTelescopeTheta()*TMath::RadToDeg());
254 txt += "\\circ ";
255 if (GetPhotElfromShower()>=10000)
256 txt += MString::Format("%dk", TMath::Nint(GetPhotElfromShower()/1000.));
257 else
258 if (GetPhotElfromShower()>=1000)
259 txt += MString::Format("%.1fk", GetPhotElfromShower()/1000.);
260 else
261 txt += GetPhotElfromShower();
262 txt += "PhEl}";
263
264 return txt;
265}
266
267Bool_t MMcEvt::SetupFits(fits &fin)
268{
269 if (fName!="MMcEvt")
270 {
271 *fLog << err << "ERROR - SetupFits only supported if name equals MMcEvt." << endl;
272 return kFALSE;
273 }
274
275 if (!fin.SetRefAddress("MMcEvt.fEvtNumber", fEvtNumber )) return kFALSE;
276 if (!fin.SetRefAddress("MMcEvt.fThick0", fThick0 )) return kFALSE;
277 if (!fin.SetRefAddress("MMcEvt.fFirstTarget", fFirstTarget )) return kFALSE;
278 if (!fin.SetRefAddress("MMcEvt.fZFirstInteraction", fZFirstInteraction)) return kFALSE;
279 if (!fin.SetRefAddress("MMcEvt.fCoreD", fCoreD )) return kFALSE;
280 if (!fin.SetRefAddress("MMcEvt.fCoreX", fCoreX )) return kFALSE;
281 if (!fin.SetRefAddress("MMcEvt.fCoreY", fCoreY )) return kFALSE;
282 if (!fin.SetRefAddress("MMcEvt.fTimeFirst", fTimeFirst )) return kFALSE;
283 if (!fin.SetRefAddress("MMcEvt.fTimeLast", fTimeLast )) return kFALSE;
284 if (!fin.SetRefAddress("MMcEvt.fLongiNmax", fLongiNmax )) return kFALSE;
285 if (!fin.SetRefAddress("MMcEvt.fLongit0", fLongit0 )) return kFALSE;
286 if (!fin.SetRefAddress("MMcEvt.fLongitmax", fLongitmax )) return kFALSE;
287 if (!fin.SetRefAddress("MMcEvt.fLongia", fLongia )) return kFALSE;
288 if (!fin.SetRefAddress("MMcEvt.fLongib", fLongib )) return kFALSE;
289 if (!fin.SetRefAddress("MMcEvt.fLongic", fLongic )) return kFALSE;
290 if (!fin.SetRefAddress("MMcEvt.fLongichi2", fLongichi2 )) return kFALSE;
291 if (!fin.SetRefAddress("MMcEvt.fPhotIni", fPhotIni )) return kFALSE;
292 if (!fin.SetRefAddress("MMcEvt.fPassPhotAtm", fPassPhotAtm )) return kFALSE;
293 if (!fin.SetRefAddress("MMcEvt.fPassPhotRef", fPassPhotRef )) return kFALSE;
294 if (!fin.SetRefAddress("MMcEvt.fPassPhotCone", fPassPhotCone )) return kFALSE;
295 if (!fin.SetRefAddress("MMcEvt.fPhotElfromShower", fPhotElfromShower )) return kFALSE;
296 if (!fin.SetRefAddress("MMcEvt.fPhotElinCamera", fPhotElinCamera )) return kFALSE;
297 if (!fin.SetRefAddress("MMcEvt.fElecCphFraction", fElecCphFraction )) return kFALSE;
298 if (!fin.SetRefAddress("MMcEvt.fMuonCphFraction", fMuonCphFraction )) return kFALSE;
299 if (!fin.SetRefAddress("MMcEvt.fOtherCphFraction", fOtherCphFraction )) return kFALSE;
300 if (!fin.SetRefAddress("MMcEvt.fFadcTimeJitter", fFadcTimeJitter )) return kFALSE;
301 if (!fin.SetRefAddress("MMcEvt.fEventReuse", fEventReuse )) return kFALSE;
302
303 // Not in files written from MMcEvt < V8
304 if (fin.HasColumn("MMcEvt.fRunNumber"))
305 if (!fin.SetRefAddress("MMcEvt.fRunNumber", fRunNumber))
306 return kFALSE;
307
308 return MMcEvtBasic::SetupFits(fin);
309}
Note: See TracBrowser for help on using the repository browser.