| 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): Abelardo Moralejo, 02/2005 <mailto:moralejo@pd.infn.it>
|
|---|
| 19 | !
|
|---|
| 20 | ! Copyright: MAGIC Software Development, 2000-2020
|
|---|
| 21 | !
|
|---|
| 22 | !
|
|---|
| 23 | \* ======================================================================== */
|
|---|
| 24 |
|
|---|
| 25 | /////////////////////////////////////////////////////////////////////////////
|
|---|
| 26 | //
|
|---|
| 27 | // MMcEvtBasic
|
|---|
| 28 | //
|
|---|
| 29 | // This class contains the most basic MonteCarlo information
|
|---|
| 30 | // with which an event has been generated
|
|---|
| 31 | //
|
|---|
| 32 | // Note: The azimuth fTelescopePhi angle in this and other MC classes
|
|---|
| 33 | // follow the convention in the Corsika program (see Corsika manual and
|
|---|
| 34 | // TDAS 02-11). There, phi is the azimuth of the momentum vector of
|
|---|
| 35 | // particles, and is measured from the north direction, anticlockwise
|
|---|
| 36 | // (i.e, west is phi=90 degrees). When it refers to the telescope
|
|---|
| 37 | // orientation, it is the azimuth of a vector along the telescope axis,
|
|---|
| 38 | // going from the camera to the mirror. So, fTelescopeTheta=90,
|
|---|
| 39 | // fTelescopePhi = 0 means the telescope is pointing horizontally towards
|
|---|
| 40 | // South.
|
|---|
| 41 | //
|
|---|
| 42 | //
|
|---|
| 43 | // Version 1:
|
|---|
| 44 | // New container to keep the very basic informations on the
|
|---|
| 45 | // original MC events produced by Corsika
|
|---|
| 46 | //
|
|---|
| 47 | // Version 2:
|
|---|
| 48 | // - added typedef for ParticleId_t from MMcEvt
|
|---|
| 49 | // - replaced MMcEvt::ParticleId_t by ParticleId_t
|
|---|
| 50 | //
|
|---|
| 51 | // Version 3:
|
|---|
| 52 | // - moved fPhi from MMcEvt
|
|---|
| 53 | // - moved fTheta from MMcEvt
|
|---|
| 54 | //
|
|---|
| 55 | // Version 4:
|
|---|
| 56 | // - moved CoreX from MMcEvt
|
|---|
| 57 | // - moved CoreY from MMcEvt
|
|---|
| 58 | //
|
|---|
| 59 | /////////////////////////////////////////////////////////////////////////////
|
|---|
| 60 | #include "MMcEvtBasic.h"
|
|---|
| 61 |
|
|---|
| 62 | #include "fits.h"
|
|---|
| 63 |
|
|---|
| 64 | #include "MString.h"
|
|---|
| 65 |
|
|---|
| 66 | #include "MLog.h"
|
|---|
| 67 | #include "MLogManip.h"
|
|---|
| 68 |
|
|---|
| 69 | ClassImp(MMcEvtBasic);
|
|---|
| 70 |
|
|---|
| 71 | using namespace std;
|
|---|
| 72 |
|
|---|
| 73 | // --------------------------------------------------------------------------
|
|---|
| 74 | //
|
|---|
| 75 | // Default constructor. Calls Clear()
|
|---|
| 76 | //
|
|---|
| 77 | MMcEvtBasic::MMcEvtBasic()
|
|---|
| 78 | {
|
|---|
| 79 | fName = "MMcEvtBasic";
|
|---|
| 80 | fTitle = "Basic event info from Monte Carlo";
|
|---|
| 81 |
|
|---|
| 82 | Clear();
|
|---|
| 83 | }
|
|---|
| 84 |
|
|---|
| 85 | // --------------------------------------------------------------------------
|
|---|
| 86 | //
|
|---|
| 87 | // Constructor. Use this to set all data members
|
|---|
| 88 | //
|
|---|
| 89 | // THIS FUNCTION IS FOR THE SIMULATION OLNY.
|
|---|
| 90 | // DON'T USE THIS MEMBERFUNCTION IN THE ANALYSIS.
|
|---|
| 91 | //
|
|---|
| 92 | MMcEvtBasic::MMcEvtBasic(ParticleId_t usPId, Float_t fEner,
|
|---|
| 93 | Float_t fImpa, Float_t fTPhii, Float_t fTThet)
|
|---|
| 94 | : fPartId(usPId), fEnergy(fEner), fImpact(fImpa),
|
|---|
| 95 | fTelescopePhi(fTPhii), fTelescopeTheta(fTThet), fTheta(0), fPhi(0), fCoreX(0), fCoreY(0)
|
|---|
| 96 | {
|
|---|
| 97 | fName = "MMcEvtBasic";
|
|---|
| 98 | fTitle = "Basic event info from Monte Carlo";
|
|---|
| 99 | }
|
|---|
| 100 |
|
|---|
| 101 | // --------------------------------------------------------------------------
|
|---|
| 102 | //
|
|---|
| 103 | // Copy operator. Copy all data members
|
|---|
| 104 | //
|
|---|
| 105 | void MMcEvtBasic::operator=(const MMcEvtBasic &evt)
|
|---|
| 106 | {
|
|---|
| 107 | fPartId = evt.fPartId;
|
|---|
| 108 | fEnergy = evt.fEnergy;
|
|---|
| 109 | fImpact = evt.fImpact;
|
|---|
| 110 | fTelescopePhi = evt.fTelescopePhi;
|
|---|
| 111 | fTelescopeTheta = evt.fTelescopeTheta;
|
|---|
| 112 | fPhi = evt.fPhi;
|
|---|
| 113 | fTheta = evt.fTheta;
|
|---|
| 114 | fCoreX = evt.fCoreX;
|
|---|
| 115 | fCoreY = evt.fCoreY;
|
|---|
| 116 | }
|
|---|
| 117 |
|
|---|
| 118 | // --------------------------------------------------------------------------
|
|---|
| 119 | //
|
|---|
| 120 | // Reset all values.
|
|---|
| 121 | //
|
|---|
| 122 | void MMcEvtBasic::Clear(Option_t *opt)
|
|---|
| 123 | {
|
|---|
| 124 | fPartId = kUNDEFINED;
|
|---|
| 125 |
|
|---|
| 126 | fEnergy = -1;
|
|---|
| 127 | fImpact = -1;
|
|---|
| 128 |
|
|---|
| 129 | fTelescopePhi = 0;
|
|---|
| 130 | fTelescopeTheta = 0;
|
|---|
| 131 |
|
|---|
| 132 | fTheta = 0;
|
|---|
| 133 | fPhi = 0;
|
|---|
| 134 |
|
|---|
| 135 | fCoreX = 0;
|
|---|
| 136 | fCoreY = 0;
|
|---|
| 137 | }
|
|---|
| 138 |
|
|---|
| 139 | // --------------------------------------------------------------------------
|
|---|
| 140 | //
|
|---|
| 141 | // Fill all data members
|
|---|
| 142 | //
|
|---|
| 143 | void MMcEvtBasic::Fill(ParticleId_t usPId, Float_t fEner,
|
|---|
| 144 | Float_t fImpa, Float_t fTPhii, Float_t fTThet)
|
|---|
| 145 | {
|
|---|
| 146 | fPartId = usPId;
|
|---|
| 147 |
|
|---|
| 148 | fEnergy = fEner;
|
|---|
| 149 | fImpact = fImpa;
|
|---|
| 150 |
|
|---|
| 151 | fTelescopePhi = fTPhii;
|
|---|
| 152 | fTelescopeTheta = fTThet;
|
|---|
| 153 | }
|
|---|
| 154 |
|
|---|
| 155 | TString MMcEvtBasic::GetParticleName(Int_t id)
|
|---|
| 156 | {
|
|---|
| 157 | switch (id)
|
|---|
| 158 | {
|
|---|
| 159 | case kUNDEFINED: return "Undefined";
|
|---|
| 160 | case kGAMMA: return "Gamma";
|
|---|
| 161 | case kPOSITRON: return "Positron";
|
|---|
| 162 | case kELECTRON: return "Electron";
|
|---|
| 163 | case kANTIMUON: return "Anti-Muon";
|
|---|
| 164 | case kMUON: return "Muon";
|
|---|
| 165 | case kPI0: return "Pi-0";
|
|---|
| 166 | case kNEUTRON: return "Neutron";
|
|---|
| 167 | case kPROTON: return "Proton";
|
|---|
| 168 | case kHELIUM: return "Helium";
|
|---|
| 169 | case kOXYGEN: return "Oxygen";
|
|---|
| 170 | case kIRON: return "Iron";
|
|---|
| 171 | case kArtificial: return "Artificial";
|
|---|
| 172 | case kNightSky: return "NightSky";
|
|---|
| 173 | }
|
|---|
| 174 |
|
|---|
| 175 | return MString::Format("Id:%d", id);
|
|---|
| 176 | }
|
|---|
| 177 |
|
|---|
| 178 | TString MMcEvtBasic::GetParticleSymbol(Int_t id)
|
|---|
| 179 | {
|
|---|
| 180 | switch (id)
|
|---|
| 181 | {
|
|---|
| 182 | case kUNDEFINED:return "N/A";
|
|---|
| 183 | case kGAMMA: return "\\gamma";
|
|---|
| 184 | case kPOSITRON: return "e^{+}";
|
|---|
| 185 | case kELECTRON: return "e^{-}";
|
|---|
| 186 | case kANTIMUON: return "\\mu^{+}";
|
|---|
| 187 | case kMUON: return "\\mu^{-}";
|
|---|
| 188 | case kPI0: return "\\pi^{0}";
|
|---|
| 189 | case kNEUTRON: return "n";
|
|---|
| 190 | case kPROTON: return "p";
|
|---|
| 191 | case kHELIUM: return "He";
|
|---|
| 192 | case kOXYGEN: return "O";
|
|---|
| 193 | case kIRON: return "Fe";
|
|---|
| 194 | case kNightSky: return "\\gamma_{NSB}";
|
|---|
| 195 | }
|
|---|
| 196 |
|
|---|
| 197 | return MString::Format("Id:%d", id);
|
|---|
| 198 | }
|
|---|
| 199 |
|
|---|
| 200 | TString MMcEvtBasic::GetEnergyStr(Float_t e)
|
|---|
| 201 | {
|
|---|
| 202 | if (e>=1000)
|
|---|
| 203 | return MString::Format("%.1fTeV", e/1000);
|
|---|
| 204 |
|
|---|
| 205 | if (e>=10)
|
|---|
| 206 | return MString::Format("%dGeV", (Int_t)(e+.5));
|
|---|
| 207 |
|
|---|
| 208 | if (e>=1)
|
|---|
| 209 | return MString::Format("%.1fGeV", e);
|
|---|
| 210 |
|
|---|
| 211 | return MString::Format("%dMeV", (Int_t)(e*1000+.5));
|
|---|
| 212 | }
|
|---|
| 213 |
|
|---|
| 214 |
|
|---|
| 215 | // --------------------------------------------------------------------------
|
|---|
| 216 | //
|
|---|
| 217 | // Print the contents of the container.
|
|---|
| 218 | //
|
|---|
| 219 | // if you specify an option only the requested data members are printed:
|
|---|
| 220 | // allowed options are: id, energy, impact
|
|---|
| 221 | //
|
|---|
| 222 | void MMcEvtBasic::Print(Option_t *opt) const
|
|---|
| 223 | {
|
|---|
| 224 | //
|
|---|
| 225 | // print out the data member on screen
|
|---|
| 226 | //
|
|---|
| 227 | TString str(opt);
|
|---|
| 228 | if (str.IsNull())
|
|---|
| 229 | {
|
|---|
| 230 | *fLog << all << endl;
|
|---|
| 231 | *fLog << "Monte Carlo output:" << endl;
|
|---|
| 232 | if (fPartId>=0)
|
|---|
| 233 | *fLog << " Particle Id: " << GetParticleName() << " [" << fPartId << "]" << endl;
|
|---|
| 234 | *fLog << " Energy: " << fEnergy << "GeV" << endl;
|
|---|
| 235 | *fLog << " Impactparam.: " << fImpact/100 << "m [x=" << fCoreX << ", y=" << fCoreY << "]" << endl;
|
|---|
| 236 | *fLog << endl;
|
|---|
| 237 | return;
|
|---|
| 238 | }
|
|---|
| 239 | if (str.Contains("id", TString::kIgnoreCase))
|
|---|
| 240 | *fLog << "Particle: " << GetParticleName() << endl;
|
|---|
| 241 | if (str.Contains("energy", TString::kIgnoreCase))
|
|---|
| 242 | *fLog << "Energy: " << fEnergy << "GeV" << endl;
|
|---|
| 243 | if (str.Contains("impact", TString::kIgnoreCase))
|
|---|
| 244 | *fLog << "Impact: " << fImpact << "cm [x=" << fCoreX << ", y=" << fCoreY << "]" << endl;
|
|---|
| 245 | }
|
|---|
| 246 |
|
|---|
| 247 | Bool_t MMcEvtBasic::SetupFits(fits &fin)
|
|---|
| 248 | {
|
|---|
| 249 | if (ClassName()==TString("MMcEvtBasic") && fName!="MMcEvtBasic")
|
|---|
| 250 | {
|
|---|
| 251 | *fLog << err << "ERROR - SetupFits only supported if name equals MMcEvtBasic." << endl;
|
|---|
| 252 | return kFALSE;
|
|---|
| 253 | }
|
|---|
| 254 |
|
|---|
| 255 | //if (!fin.SetRefAddress("MMcEvtBasic.fPartId", fPartId )) return kFALSE;
|
|---|
| 256 | if (!fin.SetRefAddress("MMcEvtBasic.fEnergy", fEnergy )) return kFALSE;
|
|---|
| 257 | if (!fin.SetRefAddress("MMcEvtBasic.fImpact", fImpact )) return kFALSE;
|
|---|
| 258 | if (!fin.SetRefAddress("MMcEvtBasic.fTelescopePhi", fTelescopePhi )) return kFALSE;
|
|---|
| 259 | if (!fin.SetRefAddress("MMcEvtBasic.fTelescopeTheta", fTelescopeTheta)) return kFALSE;
|
|---|
| 260 | if (!fin.SetRefAddress("MMcEvtBasic.fTheta", fTheta )) return kFALSE;
|
|---|
| 261 | if (!fin.SetRefAddress("MMcEvtBasic.fPhi", fPhi )) return kFALSE;
|
|---|
| 262 |
|
|---|
| 263 | // Written to MMcEvt for MMcEvtBasic < V4
|
|---|
| 264 | if (fin.HasColumn("MMcEvtBasic.fCoreX"))
|
|---|
| 265 | if (!fin.SetRefAddress("MMcEvtBasic.fCoreX", fCoreX))
|
|---|
| 266 | return kFALSE;
|
|---|
| 267 |
|
|---|
| 268 | if (fin.HasColumn("MMcEvtBasic.fCoreY"))
|
|---|
| 269 | if (!fin.SetRefAddress("MMcEvtBasic.fCoreY", fCoreY))
|
|---|
| 270 | return kFALSE;
|
|---|
| 271 |
|
|---|
| 272 | return kTRUE;
|
|---|
| 273 | }
|
|---|