Belle II Software  release-05-02-19
ARICHBackgroundModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Luka Santelj *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 // Own include
11 
12 #include <arich/modules/arichBackground/ARICHBackgroundModule.h>
13 
14 #include <time.h>
15 
16 // Hit classes
17 
18 #include <simulation/dataobjects/BeamBackHit.h>
19 #include <arich/dataobjects/ARICHSimHit.h>
20 #include <arich/dataobjects/ARICHDigit.h>
21 
22 #include <mdst/dataobjects/MCParticle.h>
23 
24 // framework - DataStore
25 #include <framework/datastore/StoreArray.h>
26 #include <framework/datastore/RelationIndex.h>
27 #include <framework/datastore/RelationArray.h>
28 
29 // framework aux
30 #include <framework/logging/Logger.h>
31 
32 using namespace std;
33 using namespace boost;
34 
35 namespace Belle2 {
40  namespace arich {
41  //-----------------------------------------------------------------
42  // Register the Module
43  //-----------------------------------------------------------------
44 
45  REG_MODULE(ARICHBackground)
46 
47 
48  //-----------------------------------------------------------------
49  // Implementation
50  //-----------------------------------------------------------------
51 
52  ARICHBackgroundModule::ARICHBackgroundModule() : Module(), phpos(TVector3()), phmom(TVector3()),
53  phVtx(TVector3()), phMmom(TVector3()), phMvtx(TVector3()), phPvtx(TVector3()), phPmom(TVector3()),
54  phGMvtx(TVector3()), phGMmom(TVector3()), modOrig(TVector3()), source(0), phPDG(0), phMPDG(0), phPPDG(0),
55  phGMPDG(0), type(0), edep(0.0), ttime(0.0), moduleID(0), phnw(0.0), trlen(0.0), en(0.0), ff(NULL),
56  TrHits(NULL)
57  {
58  // Set description()
59  setDescription("ARICHBackground module. Used to extract information relevant for ARICH background from background files");
60 
61  // Add parameters
62  addParam("FileName", m_filename, "output file name", string("mytree.root"));
63  addParam("BkgTag", m_bkgTag, "background type tag (appended to hits in the output tree", 0);
64 
65  }
66 
67  ARICHBackgroundModule::~ARICHBackgroundModule()
68  {
69 
70  }
71 
72  void ARICHBackgroundModule::initialize()
73  {
74  // Print set parameters
75  printModuleParams();
76 
77  ff = new TFile(m_filename.c_str(), "RECREATE");
78 
79  TrHits = new TTree("TrHits", "tree of arich background hits");
80  TrHits->Branch("source", &source, "source/I");
81  TrHits->Branch("phpos", "TVector3", &phpos);
82  TrHits->Branch("phmom", "TVector3", &phmom);
83  TrHits->Branch("phVtx", "TVector3", &phVtx);
84  TrHits->Branch("phMmom", "TVector3", &phMmom);
85  TrHits->Branch("phMvtx", "TVector3", &phMvtx);
86  TrHits->Branch("phPvtx", "TVector3", &phPvtx);
87  TrHits->Branch("phPmom", "TVector3", &phPmom);
88  TrHits->Branch("phGMvtx", "TVector3", &phGMvtx);
89  TrHits->Branch("phGMmom", "TVector3", &phGMmom);
90  TrHits->Branch("modOrig", "TVector3", &modOrig);
91  TrHits->Branch("phPDG", &phPDG, "phPDG/I");
92  TrHits->Branch("phMPDG", &phMPDG, "phMPDG/I");
93  TrHits->Branch("phPPDG", &phPPDG, "phPPDG/I");
94  TrHits->Branch("phGMPDG", &phGMPDG, "phGMPDG/I");
95  TrHits->Branch("type", &type, "type/I");
96  TrHits->Branch("edep", &edep, "edep/D");
97  TrHits->Branch("ttime", &ttime, "ttime/D");
98  TrHits->Branch("moduleID", &moduleID, "moduleID/I");
99  TrHits->Branch("phnw", &phnw, "phnw/D");
100  TrHits->Branch("trlen", &trlen, "trlen/D");
101  TrHits->Branch("en", &en, "en/D");
102 
103  source = m_bkgTag;
104  }
105 
106  void ARICHBackgroundModule::beginRun()
107  {
108  // Print run number
109  B2INFO("ARICHBackground: Processing. ");
110 
111  }
112 
113  void ARICHBackgroundModule::event()
114  {
115 
116  StoreArray<MCParticle> mcParticles;
117 
118  StoreArray<ARICHSimHit> arichSimHits;
119 
120  StoreArray<BeamBackHit> beamBackHits;
121 
122  StoreArray<ARICHDigit> arichDigits;
123 
124  int nHits = beamBackHits.getEntries();
125 
126  for (int iHit = 0; iHit < nHits; ++iHit) {
127 
128  BeamBackHit* arichhit = beamBackHits[iHit];
129  int subdet = arichhit->getSubDet();
130  if (subdet != 4) continue;
131  type = 0;
132  if (arichhit->getIdentifier() == 1) type = 1;
133  edep = arichhit->getEnergyDeposit();
134  ttime = arichhit->getTime();
135  phPDG = arichhit->getPDG();
136  phpos = arichhit->getPosition();
137  phmom = arichhit->getMomentum();
138  moduleID = m_arichgp->getDetectorPlane().pointSlotID(phpos.X(), phpos.Y());
139  double r = m_arichgp->getDetectorPlane().getSlotR(moduleID);
140  double phi = m_arichgp->getDetectorPlane().getSlotPhi(moduleID);
141  modOrig = TVector3(r * cos(phi), r * sin(phi), 0);
142  en = arichhit->getEnergy();
143 
144  if (phPDG == 2112) {
145  phnw = arichhit->getNeutronWeight();
146  trlen = arichhit->getTrackLength();
147  }
148  phVtx = TVector3(0, 0, 0); phMvtx = TVector3(0, 0, 0); phMmom = TVector3(0, 0, 0);
149  phMPDG = -1; phPPDG = -1; phGMPDG = -1;
150  phPvtx = TVector3(0, 0, 0); phPmom = TVector3(0, 0, 0); phGMvtx = TVector3(0, 0, 0); phGMmom = TVector3(0, 0, 0);
151 
152  RelationIndex<MCParticle, BeamBackHit> relBeamBackHitToMCParticle(mcParticles, beamBackHits);
153  if (relBeamBackHitToMCParticle.getFirstElementTo(arichhit)) {
154  const MCParticle* currParticle = relBeamBackHitToMCParticle.getFirstElementTo(arichhit)->from;
155  phVtx = currParticle->getVertex();
156  const MCParticle* mother = currParticle->getMother();
157  int mm = 0;
158  while (mother) {
159  if (mm == 0) {
160  phMPDG = mother->getPDG();
161  phMvtx = mother->getVertex();
162  phMmom = mother->getMomentum();
163  }
164  if (mm == 1) {
165  phGMPDG = mother->getPDG();
166  phGMvtx = mother->getVertex();
167  phGMmom = mother->getMomentum();
168  }
169  const MCParticle* pommother = mother->getMother();
170  if (!pommother) {
171  phPPDG = mother->getPDG();
172  phPvtx = mother->getVertex();
173  phPmom = mother->getMomentum();
174  }
175  mother = pommother;
176  mm++;
177  }
178 
179  }
180  TrHits->Fill();
181  }
182 
183  nHits = arichSimHits.getEntries();
184 
185  for (int iHit = 0; iHit < nHits; ++iHit) {
186  ARICHSimHit* simHit = arichSimHits[iHit];
187  moduleID = simHit->getModuleID();
188  type = 2;
189  phPDG = 0;
190  edep = 0;
191  ttime = simHit->getGlobalTime();
192 
193  phVtx = TVector3(0, 0, 0); phMvtx = TVector3(0, 0, 0); phMmom = TVector3(0, 0, 0);
194  phMPDG = -1; phPPDG = -1; phGMPDG = -1; phmom = TVector3(0, 0, 0);
195  phPvtx = TVector3(0, 0, 0); phPmom = TVector3(0, 0, 0); phGMvtx = TVector3(0, 0, 0); phGMmom = TVector3(0, 0, 0);
196  RelationIndex<MCParticle, ARICHSimHit> relSimHitToMCParticle(mcParticles, arichSimHits);
197  if (relSimHitToMCParticle.getFirstElementTo(simHit)) {
198  const MCParticle* currParticle = relSimHitToMCParticle.getFirstElementTo(simHit)->from;
199  phVtx = currParticle->getVertex();
200  const MCParticle* mother = currParticle->getMother();
201  int mm = 0;
202  while (mother) {
203  if (mm == 0) {
204  phMPDG = mother->getPDG();
205  phMvtx = mother->getVertex();
206  phMmom = mother->getMomentum();
207  }
208  if (mm == 1) {
209  phGMPDG = mother->getPDG();
210  phGMvtx = mother->getVertex();
211  phGMmom = mother->getMomentum();
212  }
213  const MCParticle* pommother = mother->getMother();
214  if (!pommother) {
215  phPPDG = mother->getPDG();
216  phPvtx = mother->getVertex();
217  phPmom = mother->getMomentum();
218  }
219  mother = pommother;
220  mm++;
221  }
222 
223  }
224  TrHits->Fill();
225  }
226  }
227 
228  void ARICHBackgroundModule::endRun()
229  {
230  }
231 
232  void ARICHBackgroundModule::terminate()
233  {
234  // CPU time end
235 
236  // Announce
237  B2INFO("ARICHBackground finished.");
238 
239  //TrHits->Write();
240  ff->Write();
241  ff->Close();
242 
243  }
244 
245  void ARICHBackgroundModule::printModuleParams() const
246  {
247 
248  }
249 
250 
251  } // end arich namespace
253 } // end Belle2 namespace
Belle2::RelationIndex::getFirstElementTo
const Element * getFirstElementTo(const TO &to) const
Return a pointer to the first relation Element of the given object.
Definition: RelationIndex.h:222
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::BeamBackHit::getEnergyDeposit
double getEnergyDeposit() const
Get particle energy deposit in sensitive volume.
Definition: BeamBackHit.h:116
Belle2::BeamBackHit::getMomentum
const TVector3 & getMomentum() const
Get momentum of the particle hit.
Definition: BeamBackHit.h:104
Belle2::RelationIndex
Provides access to fast ( O(log n) ) bi-directional lookups on a specified relation.
Definition: RelationIndex.h:84
Belle2::ARICHSimHit::getModuleID
int getModuleID() const
Get ID number of module that registered hit.
Definition: ARICHSimHit.h:76
Belle2::BeamBackHit::getEnergy
double getEnergy() const
Get energy of the particle.
Definition: BeamBackHit.h:110
Belle2::BeamBackHit::getPDG
int getPDG() const
Get the lund code of the particle that hit the sensitive area.
Definition: BeamBackHit.h:95
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::BeamBackHit::getTrackLength
double getTrackLength() const
the length of the track in the volume
Definition: BeamBackHit.h:119
Belle2::BeamBackHit::getPosition
const TVector3 & getPosition() const
Get global position of the particle hit.
Definition: BeamBackHit.h:101
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::MCParticle::getVertex
TVector3 getVertex() const
Return production vertex position, shorthand for getProductionVertex().
Definition: MCParticle.h:194
Belle2::ARICHSimHit
Class ARICHSimHit - Geant4 simulated hit for ARICH.
Definition: ARICHSimHit.h:39
Belle2::MCParticle::getMother
MCParticle * getMother() const
Returns a pointer to the mother particle.
Definition: MCParticle.h:593
Belle2::BeamBackHit::getTime
double getTime() const
Get the time at which the hit occured.
Definition: BeamBackHit.h:107
Belle2::arich::ARICHBackgroundModule
ARICH digitizer module.
Definition: ARICHBackgroundModule.h:44
Belle2::BeamBackHit::getNeutronWeight
double getNeutronWeight() const
get the effective neutron weigth
Definition: BeamBackHit.h:122
Belle2::BeamBackHit::getIdentifier
int getIdentifier() const
Get the identifier of subdetector component in which hit occured.
Definition: BeamBackHit.h:89
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::StoreArray< MCParticle >
Belle2::BeamBackHit::getSubDet
int getSubDet() const
Det the index of subdetector in which hit occured.
Definition: BeamBackHit.h:92
Belle2::BeamBackHit
Class BeamBackHit - Stores hits from beam backgound simulation.
Definition: BeamBackHit.h:39
Belle2::ARICHSimHit::getGlobalTime
float getGlobalTime() const override
Get global time of hit.
Definition: ARICHSimHit.h:82
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226