Belle II Software  release-05-02-19
ARICHReconstructorModule.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, Dino Tahirovic *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 // Own include
11 
12 #include <arich/modules/arichReconstruction/ARICHReconstructorModule.h>
13 #include <time.h>
14 
15 // Hit classes
16 #include <mdst/dataobjects/Track.h>
17 #include <tracking/dataobjects/ExtHit.h>
18 #include <mdst/dataobjects/MCParticle.h>
19 #include <arich/dataobjects/ARICHAeroHit.h>
20 #include <arich/dataobjects/ARICHLikelihood.h>
21 #include <arich/dataobjects/ARICHTrack.h>
22 #include <arich/dataobjects/ARICHHit.h>
23 
24 // framework - DataStore
25 #include <framework/datastore/DataStore.h>
26 #include <framework/datastore/StoreArray.h>
27 
28 // framework aux
29 #include <framework/gearbox/Unit.h>
30 #include <framework/gearbox/Const.h>
31 #include <framework/logging/Logger.h>
32 
33 // ROOT
34 #include <TVector3.h>
35 
36 using namespace std;
37 using namespace boost;
38 
39 namespace Belle2 {
45 //-----------------------------------------------------------------
46 // Register the Module
47 //-----------------------------------------------------------------
48  REG_MODULE(ARICHReconstructor)
49 
50 
51 //-----------------------------------------------------------------
52 // Implementation
53 //-----------------------------------------------------------------
54 
56  m_ana(NULL)
57  {
58  // Set description()
59  setDescription("This module calculates the ARICHLikelihood values for all particle id. hypotheses, for all tracks that enter ARICH in the event.");
60 
61  // Set property flags
62  setPropertyFlags(c_ParallelProcessingCertified);
63 
64  // Add parameters
65  addParam("trackPositionResolution", m_trackPositionResolution,
66  "Resolution of track position on aerogel plane (for additional smearing of MC tracks)", 1.0 * Unit::mm);
67  addParam("trackAngleResolution", m_trackAngleResolution,
68  "Resolution of track direction angle on aerogel plane (for additional smearing of MC tracks)", 2.0 * Unit::mrad);
69  addParam("inputTrackType", m_inputTrackType, "Input tracks switch: tracking (0), from AeroHits - MC info (1)", 0);
70  addParam("storePhotons", m_storePhot, "Set to 1 to store reconstructed photon information (Ch. angle,...)", 0);
71  addParam("useAlignment", m_align, "Use ARICH global position alignment constants", true);
72  addParam("useMirrorAlignment", m_alignMirrors, "Use ARICH mirror alignment constants", true);
73  }
74 
75  ARICHReconstructorModule::~ARICHReconstructorModule()
76  {
77  if (m_ana) delete m_ana;
78  }
79 
80  void ARICHReconstructorModule::initialize()
81  {
82  // Initialize variables
83  if (m_ana) delete m_ana;
84  m_ana = new ARICHReconstruction(m_storePhot);
85  m_ana->setTrackPositionResolution(m_trackPositionResolution);
86  m_ana->setTrackAngleResolution(m_trackAngleResolution);
87  m_ana->useMirrorAlignment(m_alignMirrors);
88 
89  StoreArray<ARICHHit> arichHits;
90  arichHits.isRequired();
91 
92  StoreArray<Track> tracks;
93  StoreArray<ExtHit> extHits;
94  StoreArray<ARICHAeroHit> aeroHits;
95  tracks.isOptional();
96  extHits.isOptional();
97  aeroHits.isOptional();
98 
99  if (!aeroHits.isOptional()) {
100  tracks.isRequired();
101  extHits.isRequired();
102  }
103 
104  StoreArray<MCParticle> mcParticles;
105  mcParticles.isOptional();
106 
107  StoreArray<ARICHLikelihood> likelihoods;
108  likelihoods.registerInDataStore();
109 
110  StoreArray<ARICHTrack> arichTracks;
111  arichTracks.registerInDataStore();
112 
113  if (m_inputTrackType) arichTracks.registerRelationTo(likelihoods);
114  else {
115  arichTracks.registerRelationTo(extHits);
116  tracks.registerRelationTo(likelihoods);
117  }
118  //arichTracks.registerInDataStore(DataStore::c_DontWriteOut);
119  //arichTracks.registerRelationTo(likelihoods, DataStore::c_Event, DataStore::c_DontWriteOut);
120  arichTracks.registerRelationTo(aeroHits);
121  printModuleParams();
122  }
123 
124  void ARICHReconstructorModule::beginRun()
125  {
126  m_ana->initialize();
127  }
128 
129  void ARICHReconstructorModule::event()
130  {
131 
132  // Output - log likelihoods
133  StoreArray<ARICHLikelihood> arichLikelihoods;
134 
135  // input AeroHits
136  StoreArray<ARICHTrack> arichTracks;
137 
138  // Input: ARICHDigits
139  StoreArray<ARICHHit> arichHits;
140 
141  // using track information form tracking system (mdst Track)
142  if (m_inputTrackType == 0) {
143 
144  StoreArray<Track> Tracks;
145 
146  Const::EDetector myDetID = Const::EDetector::ARICH; // arich
147  Const::ChargedStable hypothesis = Const::pion;
148  int pdgCode = abs(hypothesis.getPDGCode());
149 
150  for (int itrk = 0; itrk < Tracks.getEntries(); ++itrk) {
151 
152  const Track* track = Tracks[itrk];
153  const TrackFitResult* fitResult = track->getTrackFitResultWithClosestMass(hypothesis);
154  if (!fitResult) {
155  B2ERROR("No TrackFitResult for " << hypothesis.getPDGCode());
156  continue;
157  }
158 
159  const MCParticle* particle = track->getRelated<MCParticle>();
160 
161  ARICHAeroHit* aeroHit = NULL;
162  if (particle) aeroHit = particle->getRelated<ARICHAeroHit>();
163 
164  RelationVector<ExtHit> extHits = DataStore::getRelationsWithObj<ExtHit>(track);
165  ARICHTrack* arichTrack = NULL;
166 
167  //const ExtHit* arich1stHit = NULL;
168  const ExtHit* arich2ndHit = NULL;
169  const ExtHit* arichWinHit = NULL;
170 
171  for (unsigned i = 0; i < extHits.size(); i++) {
172  const ExtHit* extHit = extHits[i];
173  if (abs(extHit->getPdgCode()) != pdgCode) continue;
174  if (extHit->getDetectorID() != myDetID) continue;
175  if (extHit->getStatus() != EXT_EXIT) continue; // particles registered at the EXIT of the Al plate
176  if (extHit->getMomentum().Z() < 0.0) continue; // track passes in backward
177  if (extHit->getCopyID() == 12345) { continue;}
178  if (extHit->getCopyID() == 6789) { arich2ndHit = extHit; continue;}
179  arichWinHit = extHit;
180  }
181 
182  if (arich2ndHit) {
183  // if aeroHit cannot be found using MCParticle check if it was already related to the extHit (by ARICHRelate module)
184  if (!aeroHit) aeroHit = arich2ndHit->getRelated<ARICHAeroHit>();
185 
186  // make new ARICHTrack
187  arichTrack = arichTracks.appendNew(arich2ndHit);
188  if (arichWinHit) arichTrack->setHapdWindowHit(arichWinHit);
189  }
190 
191  // skip if track has no extHit in ARICH
192  if (!arichTrack) continue;
193  // transform track parameters to ARICH local frame
194  m_ana->transformTrackToLocal(*arichTrack, m_align);
195  // make new ARICHLikelihood
196  ARICHLikelihood* like = arichLikelihoods.appendNew();
197  // calculate and set likelihood values
198  m_ana->likelihood2(*arichTrack, arichHits, *like);
199  // make relations
200  track->addRelationTo(like);
201  arichTrack->addRelationTo(arich2ndHit);
202  if (aeroHit) arichTrack->addRelationTo(aeroHit);
203 
204  } // Tracks loop
205  } // input type if
206 
207 
208  // using track information form MC (stored in ARICHAeroHit)
209  else {
210 
211  StoreArray<ARICHAeroHit> aeroHits;
212  int nTracks = aeroHits.getEntries();
213 
214  // Loop over all ARICHAeroHits
215  for (int iTrack = 0; iTrack < nTracks; ++iTrack) {
216  ARICHAeroHit* aeroHit = aeroHits[iTrack];
217 
218  // make new ARICHTrack
219  ARICHTrack* arichTrack = arichTracks.appendNew(aeroHit);
220  // smearing of track parameters (to mimic tracking system resolutions)
221  m_ana->smearTrack(*arichTrack);
222  // transform track parameters to ARICH local frame
223  m_ana->transformTrackToLocal(*arichTrack, m_align);
224  // make associated ARICHLikelihood
225  ARICHLikelihood* like = arichLikelihoods.appendNew();
226  // calculate and set likelihood values
227  m_ana->likelihood2(*arichTrack, arichHits, *like);
228  // make relation
229  arichTrack->addRelationTo(like);
230  arichTrack->addRelationTo(aeroHit);
231  } // for iTrack
232 
233  }
234  }
235 
236  void ARICHReconstructorModule::endRun()
237  {
238  }
239 
240  void ARICHReconstructorModule::terminate()
241  {
242  }
243 
244  void ARICHReconstructorModule::printModuleParams()
245  {
246  if (m_inputTrackType == 0) { B2DEBUG(100, "ARICHReconstructorModule: track infromation is taken from mdst Tracks.");}
247  else B2DEBUG(100, "ARICHReconstructorModule: track information is taken from MC (ARICHAeroHit).");
248  }
249 
251 } // namespace Belle2
Belle2::StoreArray::appendNew
T * appendNew()
Construct a new T object at the end of the array.
Definition: StoreArray.h:256
Belle2::ARICHLikelihood
This is a class to store ARICH likelihoods in the datastore.
Definition: ARICHLikelihood.h:38
Belle2::StoreArray::registerRelationTo
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
Definition: StoreArray.h:150
Belle2::ExtHit::getPdgCode
int getPdgCode() const
Get PDG code of this extrapolation's hypothesis.
Definition: ExtHit.h:126
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::ExtHit::getDetectorID
Const::EDetector getDetectorID() const
Get detector ID of this extrapolation hit.
Definition: ExtHit.h:130
Belle2::ARICHAeroHit
Datastore class that holds information on track parameters at the entrance in aerogel.
Definition: ARICHAeroHit.h:37
Belle2::RelationsInterface::addRelationTo
void addRelationTo(const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
Add a relation from this object to another object (with caching).
Definition: RelationsObject.h:144
Belle2::Const::ParticleType::getPDGCode
int getPDGCode() const
PDG code.
Definition: Const.h:349
Belle2::RelationsInterface::getRelated
T * getRelated(const std::string &name="", const std::string &namedRelation="") const
Get the object to or from which this object has a relation.
Definition: RelationsObject.h:280
Belle2::Const::EDetector
EDetector
Enum for identifying the detector components (detector and subdetector).
Definition: Const.h:44
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::ExtHit
Store one Ext hit as a ROOT object.
Definition: ExtHit.h:40
Belle2::ARICHReconstruction
Internal ARICH track reconstruction.
Definition: ARICHReconstruction.h:60
Belle2::ExtHit::getStatus
ExtHitStatus getStatus() const
Get state of extrapolation at this hit.
Definition: ExtHit.h:138
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::ARICHTrack
Datastore class that holds position and momentum information of tracks that hit ARICH.
Definition: ARICHTrack.h:42
Belle2::ExtHit::getCopyID
int getCopyID() const
Get detector-element ID of sensitive element within detector.
Definition: ExtHit.h:134
Belle2::Const::ChargedStable
Provides a type-safe way to pass members of the chargedStableSet set.
Definition: Const.h:465
Belle2::ExtHit::getMomentum
TVector3 getMomentum() const
Get momentum at this extrapolation hit.
Definition: ExtHit.h:157
Belle2::Track
Class that bundles various TrackFitResults.
Definition: Track.h:35
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: ECLMatchingPerformanceExpertModule.h:33
Belle2::ARICHReconstructorModule
ARICH subdetector main module.
Definition: ARICHReconstructorModule.h:48
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:226