Belle II Software  release-05-01-25
SVDMaxStripTTreeModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa *
7  * *
8  * The module is used to create a TTree to study the number of strips *
9  * fired per event per APV chip. *
10  * *
11  * This software is provided "as is" without any warranty. *
12  **************************************************************************/
13 
14 #include <svd/modules/svdPerformance/SVDMaxStripTTreeModule.h>
15 #include <hlt/softwaretrigger/core/FinalTriggerDecisionCalculator.h>
16 
17 #include <vxd/dataobjects/VxdID.h>
18 #include <vxd/geometry/GeoCache.h>
19 
20 using namespace Belle2;
21 using namespace std;
22 using namespace SoftwareTrigger;
23 
24 REG_MODULE(SVDMaxStripTTree)
25 
27 {
28  setDescription("The module is used to create a TTree to study the number of strips per APV per event.");
29  addParam("outputFileName", m_rootFileName, "Name of output root file.", std::string("SVDMaxStripTTree.root"));
30  addParam("ShaperDigits", m_shapersStoreArrayName, "StoreArray name of the input ShaperDigits.",
31  m_shapersStoreArrayName);
32  addParam("skipHLTRejectedEvents", m_skipRejectedEvents, "If TRUE skip events rejected by HLT", bool(true));
33 }
34 
36 {
37  m_shapers.isRequired(m_shapersStoreArrayName);
38 
39  m_rootFilePtr = new TFile(m_rootFileName.c_str(), "RECREATE");
40 
41  //Tree for SVD u and v strips
42  m_t = new TTree("tree", "Tree for SVD u/v-strips");
43  m_t->Branch("evt", &m_event, "evt/i");
44  m_t->Branch("svdLayer", &m_svdLayer, "svdLayer/i");
45  m_t->Branch("svdLadder", &m_svdLadder, "svdLadder/i");
46  m_t->Branch("svdSensor", &m_svdSensor, "svdSensor/i");
47  m_t->Branch("svdSide", &m_svdSide, "svdSide/i");
48  m_t->Branch("svdChip", &m_svdChip, "svdChip/i");
49  m_t->Branch("svdHits", &m_svdHits, "svdHits/i");
50 
51  m_event = 0;
52 }
53 
55 {
56 
57  TH1F hHits("nHits_L@layerL@ladderS@sensor@view@apv",
58  "Number of Hits per Event in @layer.@ladder.@sensor chip @apv on the @view/@side side",
59  2, -0.5 , 1.5);
60 
61  m_hHits = new SVDAPVHistograms<TH1F>(hHits);
62 
63 }
64 
65 
67 {
68 
69  if (m_skipRejectedEvents && (m_resultStoreObjectPointer.isValid())) {
70  const bool eventAccepted = FinalTriggerDecisionCalculator::getFinalTriggerDecision(*m_resultStoreObjectPointer);
71  if (!eventAccepted) return;
72  }
73 
74  //count the number of strips per APV in the event
75  for (const auto& shaper : m_shapers)
76  m_hHits->fill(shaper.getSensorID(), shaper.isUStrip(), shaper.getCellID() / 128, 0);
77 
78  //loop on geometry and fill the tree for this event
80 
81  for (auto layer : geoCache.getLayers(VXD::SensorInfoBase::SVD))
82  for (auto ladder : geoCache.getLadders(layer))
83  for (Belle2::VxdID sensor : geoCache.getSensors(ladder))
84  for (int view = SVDAPVHistograms<TH1F>::VIndex ; view < SVDAPVHistograms<TH1F>::UIndex + 1; view++)
85  for (int apv = 0; apv < 6; apv ++) {
86  m_svdLayer = sensor.getLayerNumber();
87  m_svdLadder = sensor.getLadderNumber();
88  m_svdSensor = sensor.getSensorNumber();
89  m_svdSide = view;
90  m_svdChip = apv;
91  m_svdHits = (m_hHits->getHistogram(sensor, view, apv))->GetEntries();
92  m_t->Fill();
93 
94  //reset the histogram used as counters
95  (m_hHits->getHistogram(sensor, view, apv))->Reset();
96 
97  }
98 
99  m_event++;
100 }
101 
102 
104 {
105 
106  if (m_rootFilePtr != NULL) {
107  m_rootFilePtr->cd();
108  m_t->Write();
109  m_rootFilePtr->Close();
110  }
111 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::SVDMaxStripTTreeModule::initialize
void initialize() override
Register input and output data.
Definition: SVDMaxStripTTreeModule.cc:35
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::SVDMaxStripTTreeModule::event
void event() override
Compute the variables and fill the tree.
Definition: SVDMaxStripTTreeModule.cc:66
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::SVDMaxStripTTreeModule::terminate
void terminate() override
Write the TTrees to the file.
Definition: SVDMaxStripTTreeModule.cc:103
Belle2::SVDMaxStripTTreeModule
The module is used to create a TTree to study the number of strips fired per event per APV chip.
Definition: SVDMaxStripTTreeModule.h:41
Belle2::VXD::SensorInfoBase::SVD
@ SVD
SVD Sensor.
Definition: SensorInfoBase.h:45
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::SVDAPVHistograms< TH1F >
Belle2::SVDMaxStripTTreeModule::beginRun
void beginRun() override
Define APVHistogram.
Definition: SVDMaxStripTTreeModule.cc:54