Belle II Software  release-05-02-19
SVDTimeValidationCollectorModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2017 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Luigi Corona, Giulia Casarosa, Giulio Dujany *
7  * *
8  * *
9  * This software is provided "as is" without any warranty. *
10  **************************************************************************/
11 #include <svd/modules/svdTimeValidationCollector/SVDTimeValidationCollectorModule.h>
12 
13 using namespace std;
14 using namespace Belle2;
15 
16 //-----------------------------------------------------------------
17 // Register the Module
18 //-----------------------------------------------------------------
19 REG_MODULE(SVDTimeValidationCollector)
20 
21 //-----------------------------------------------------------------
22 // Implementation
23 //-----------------------------------------------------------------
24 
26 {
27  //Set module properties
28 
29  setDescription("Collector module used to create the histograms needed for the SVD 6-Sample CoG, 3-Sample CoG and 3-Sample ELS Time calibration");
30  setPropertyFlags(c_ParallelProcessingCertified);
31 
32  addParam("SVDClustersName", m_svdClusters, "Name of the SVDClusters list", m_svdClusters);
33  addParam("SVDClustersOnTracksName", m_svdClustersOnTracks, "Name of the SVDClustersOnTracks list", m_svdClustersOnTracks);
34  addParam("EventT0Name", m_eventTime, "Name of the EventT0 list", m_eventTime);
35  addParam("SVDEventInfoName", m_svdEventInfo, "Name of the SVDEventInfo list", m_svdEventInfo);
36  addParam("RecoTracksName", m_recotrack, "Name of the RecoTracks list", m_recotrack);
37  addParam("TracksName", m_track, "Name of the Tracks list", m_track);
38 }
39 
40 void SVDTimeValidationCollectorModule::prepare()
41 {
42 
43  TH1F hClsTimeOnTracks("clsTimeOnTracks__L@layerL@ladderS@sensor@view",
44  "clsTimeOnTracks in @layer.@ladder.@sensor @view/@side",
45  300, -150, 150);
46  hClsTimeOnTracks.GetXaxis()->SetTitle("clsTime_onTracks (ns)");
47  auto _hClsTimeOnTracks = new SVDHistograms<TH1F>(hClsTimeOnTracks);
48 
49  TH1F hClsTimeAll("clsTimeAll__L@layerL@ladderS@sensor@view",
50  "clsTimeAll in @layer.@ladder.@sensor @view/@side",
51  300, -150, 150);
52  hClsTimeAll.GetXaxis()->SetTitle("clsTime_all (ns)");
53  auto _hClsTimeAll = new SVDHistograms<TH1F>(hClsTimeAll);
54 
55  TH1F hClsDiffTimeOnTracks("clsDiffTimeOnTracks__L@layerL@ladderS@sensor@view",
56  "clsDiffTimeOnTracks in @layer.@ladder.@sensor @view/@side",
57  300, -150, 150);
58  hClsDiffTimeOnTracks.GetXaxis()->SetTitle("clsDiffTime_onTracks (ns)");
59  auto _hClsDiffTimeOnTracks = new SVDHistograms<TH1F>(hClsDiffTimeOnTracks);
60 
61  auto hEventT0 = new TH1F("hEventT0", "EventT0", 300, -150, 150);
62  registerObject<TH1F>("hEventT0", hEventT0);
63 
64  m_svdCls.isRequired(m_svdClusters);
65  m_svdClsOnTrk.isRequired(m_svdClustersOnTracks);
66  m_eventT0.isRequired(m_eventTime);
67  m_svdEI.isRequired(m_svdEventInfo);
68  m_recoTrk.isRequired(m_recotrack);
69  m_trk.isRequired(m_track);
70 
71  VXD::GeoCache& geoCache = VXD::GeoCache::getInstance();
72 
73  for (auto layer : geoCache.getLayers(VXD::SensorInfoBase::SVD)) {
74  for (auto ladder : geoCache.getLadders(layer)) {
75  for (Belle2::VxdID sensor : geoCache.getSensors(ladder)) {
76  for (int view = SVDHistograms<TH1F>::VIndex ; view < SVDHistograms<TH1F>::UIndex + 1; view++) {
77  registerObject<TH1F>(_hClsTimeOnTracks->getHistogram(sensor, view)->GetName(), _hClsTimeOnTracks->getHistogram(sensor, view));
78  registerObject<TH1F>(_hClsTimeAll->getHistogram(sensor, view)->GetName(), _hClsTimeAll->getHistogram(sensor, view));
79  registerObject<TH1F>(_hClsDiffTimeOnTracks->getHistogram(sensor, view)->GetName(), _hClsDiffTimeOnTracks->getHistogram(sensor,
80  view));
81  }
82  }
83  }
84  }
85 }
86 
87 
88 void SVDTimeValidationCollectorModule::collect()
89 {
90  if (m_eventT0->hasEventT0()) {
91  float eventT0 = m_eventT0->getEventT0();
92  getObjectPtr<TH1F>("hEventT0")->Fill(eventT0);
93 
94  // Fill histograms clusters on tracks
95  for (const auto& svdCluster : m_svdClsOnTrk) {
96  // get cluster time
97  float clTime = svdCluster.getClsTime();
98 
99  //get cluster layer, ladder sensor and side
100  auto theVxdID = svdCluster.getSensorID();
101  auto layer_num = theVxdID.getLayerNumber();
102  auto ladder_num = theVxdID.getLadderNumber();
103  auto sensor_num = theVxdID.getSensorNumber();
104  char side = 'V';
105  if (svdCluster.isUCluster())
106  side = 'U';
107 
108  auto hClsTimeOnTracks_name = Form("clsTimeOnTracks__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side);
109 
110  auto hClsDiffTimeOnTracks_name = Form("clsDiffTimeOnTracks__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side);
111 
112  getObjectPtr<TH1F>(hClsTimeOnTracks_name)->Fill(clTime);
113  getObjectPtr<TH1F>(hClsDiffTimeOnTracks_name)->Fill(clTime - eventT0);
114  };
115 
116  // Fill histograms with all clusters
117  for (const auto& svdCluster : m_svdCls) {
118  // get cluster time
119  float clTime = svdCluster.getClsTime();
120 
121  //get cluster layer, ladder sensor and side
122  auto theVxdID = svdCluster.getSensorID();
123  auto layer_num = theVxdID.getLayerNumber();
124  auto ladder_num = theVxdID.getLadderNumber();
125  auto sensor_num = theVxdID.getSensorNumber();
126  char side = 'V';
127  if (svdCluster.isUCluster())
128  side = 'U';
129 
130  auto hClsTimeAll_name = Form("clsTimeAll__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side);
131 
132  getObjectPtr<TH1F>(hClsTimeAll_name)->Fill(clTime);
133  };
134  }
135 }
Belle2::CalibrationCollectorModule
Calibration collector module base class.
Definition: CalibrationCollectorModule.h:44
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::SVDHistograms< TH1F >
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::SVDTimeValidationCollectorModule
Collector module used to create the histograms needed for the SVD CoG-Time calibration.
Definition: SVDTimeValidationCollectorModule.h:51
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41