Belle II Software  release-05-01-25
SVDCrossTalkCalibrationsCollectorModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: James Webb *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <svd/modules/svdCrossTalkCalibrationsCollector/SVDCrossTalkCalibrationsCollectorModule.h>
12 
13 #include <svd/modules/svdCrossTalkFinder/SVDCrossTalkFinderHelperFunctions.h>
14 
15 #include <vxd/geometry/GeoCache.h>
16 
17 #include <TH1F.h>
18 
19 #include <iostream>
20 
21 
22 using namespace std;
23 using namespace Belle2;
24 
25 
26 REG_MODULE(SVDCrossTalkCalibrationsCollector)
27 
28 SVDCrossTalkCalibrationsCollectorModule::SVDCrossTalkCalibrationsCollectorModule() : CalibrationCollectorModule()
29 {
30  setDescription("This module collects the list of channels exhibiting crossTalk readout.");
31 
33 
34  addParam("SVDShaperDigits", m_svdShaperDigitsName,
35  "SVDShaperDigit collection name", string(""));
36 
37  addParam("uSideOccupancyFactor", m_uSideOccupancyFactor,
38  "Multiple of the average occupancy for high occupancy strip classification", 2);
39 
40  addParam("vSideOccupancyFactor", m_vSideOccupancyFactor,
41  "Multiple of the average occupancy for high occupancy strip classification", 2);
42 
43  addParam("nAPVFactor", m_nAPVFactor,
44  "Number of APV chips with at least one high occupancy strips in event required for cross talk flag", 30);
45 
46 }
47 
49 {
50 
52 
53 
54  m_histogramTree = new TTree("tree", "tree");
55  m_histogramTree->Branch("hist", "TH1F", &m_hist, 32000, 0);
56  m_histogramTree->Branch("layer", &m_layer, "layer/I");
57  m_histogramTree->Branch("ladder", &m_ladder, "ladder/I");
58  m_histogramTree->Branch("sensor", &m_sensor, "sensor/I");
59  m_histogramTree->Branch("side", &m_side, "side/I");
60 
61  registerObject<TTree>("HTreeCrossTalkCalib", m_histogramTree);
62 
63 
64 }
65 
67 {
68 
69  std::string sensorName;
70  TH1F* sensorHist;
71  //Prepare histograms for payload.
73  for (auto& layers : geo.getLayers(VXD::SensorInfoBase::SVD)) {
74  for (auto& ladders : geo.getLadders(layers)) {
75  for (auto& sensors : geo.getSensors(ladders)) {
76  for (int side = 0; side <= 1; side++) {
77  occupancyPDFName(sensors, side, sensorName);
78  if (m_sensorHistograms.count(sensorName) == 0) {
79  if (layers.getLayerNumber() == 3 or side == 1) {
80  sensorHist = new TH1F(sensorName.c_str(), "", 768, 0, 768);
81  } else {
82  sensorHist = new TH1F(sensorName.c_str(), "", 512, 0, 512);
83  }
84  m_sensorHistograms[sensorName] = sensorHist;
85 
86  }
87 
88  }
89  }
90  }
91  }
92 
93 }
94 
96 {
97  //Vectors to hold event info.
98  vector<std::string> highOccChips_uSide ;
99  vector<std::string> highOccChips_vSide ;
100  vector<std::string> highOccChipsStripNum_uSide ;
101  vector<std::string> highOccChipsStripNum_vSide ;
102  vector<std::string> clusterStrips_uSide ;
103  vector<std::string> clusterStrips_vSide ;
104  vector<std::string> clusterChips_uSide ;
105  vector<std::string> clusterChips_vSide ;
106  vector<int> strips_uSide;
107  vector<int> strips_vSide;
108 
109  //loop over ShaperDigits
110  for (auto& svdShaperDigit : m_svdShaperDigits) {
111  //Remove L3 and +fw sensors, not affected by cross-talk
112  if (svdShaperDigit.getSensorID().getLayerNumber() == 3 or svdShaperDigit.getSensorID().getSensorNumber() == 1) {
113  continue;
114  }
115 
116  int side = svdShaperDigit.isUStrip();
117  std::string sensorName;
118  occupancyPDFName(svdShaperDigit.getSensorID(), side, sensorName);
119 
120  double sensorAverage = 0.;
121  calculateAverage(svdShaperDigit.getSensorID(), sensorAverage, side);
122  int stripID = svdShaperDigit.getCellID();
123  std::string sensorStripNum = sensorName + "." + std::to_string(stripID);
124  double stripOccupancy = m_OccupancyCal.getOccupancy(svdShaperDigit.getSensorID(), side, stripID);
125  //Clustering only works assuming digits are ordered.//
126  if (side == 1 && stripOccupancy > (m_uSideOccupancyFactor * sensorAverage)) {
127 
128  int adjacentStrip = 0;
129  if (!strips_uSide.empty()) {
130  adjacentStrip = stripID - strips_uSide.back();
131  }
132  strips_uSide.push_back(stripID);
133  if (!highOccChips_uSide.empty() && sensorName == highOccChips_uSide.back() && adjacentStrip == 1) {
134  clusterChips_uSide.push_back(sensorName);
135  if (clusterStrips_uSide.empty() || clusterStrips_uSide.back() != sensorStripNum) {
136  clusterStrips_uSide.push_back(highOccChipsStripNum_uSide.back());
137  }
138  clusterStrips_uSide.push_back(sensorStripNum);
139  }
140  highOccChipsStripNum_uSide.push_back(sensorStripNum);
141  highOccChips_uSide.push_back(sensorName);
142  }
143 
144  if (side == 0 && stripOccupancy > (m_vSideOccupancyFactor * sensorAverage)) {
145  int adjacentStrip = 0;
146  if (!strips_vSide.empty()) {
147  adjacentStrip = stripID - strips_vSide.back();
148  }
149  strips_vSide.push_back(stripID);
150  if (!highOccChips_vSide.empty() && sensorName == highOccChips_vSide.back() && adjacentStrip == 1) {
151  if (clusterStrips_vSide.empty() || clusterStrips_vSide.back() != sensorStripNum) {
152  clusterStrips_vSide.push_back(highOccChipsStripNum_vSide.back());
153  }
154  clusterStrips_vSide.push_back(sensorStripNum);
155  }
156  highOccChipsStripNum_vSide.push_back(sensorStripNum);
157  highOccChips_vSide.push_back(sensorName);
158  }
159 
160  } //ShaperDigit loop
161 
162  std::sort(clusterChips_uSide.begin(), clusterChips_uSide.end());
163  clusterChips_uSide.erase(unique(clusterChips_uSide.begin(), clusterChips_uSide.end()), clusterChips_uSide.end());
164  int numberOfClusterChips = std::distance(clusterChips_uSide.begin(), std::unique(clusterChips_uSide.begin(),
165  clusterChips_uSide.end()));
166 
167 // Crosstalk events flagged using u-side
168  if (numberOfClusterChips > m_nAPVFactor) {
169  for (auto& svdShaperDigit : m_svdShaperDigits) {
170  std::string sensorID = svdShaperDigit.getSensorID();
171  std::string digitID = sensorID + "." + std::to_string(svdShaperDigit.isUStrip());
172  std::string stripID = digitID + "." + std::to_string(svdShaperDigit.getCellID());
173  if (std::find(clusterStrips_uSide.begin(), clusterStrips_uSide.end(), stripID) != clusterStrips_uSide.end()) {
174  std::string sensorName;
175  occupancyPDFName(svdShaperDigit.getSensorID(), svdShaperDigit.isUStrip(), sensorName);
176  auto xTalkStrip = m_sensorHistograms.at(sensorName);
177  xTalkStrip->Fill(svdShaperDigit.getCellID(), 1);
178  }
179  if (std::find(clusterStrips_vSide.begin(), clusterStrips_vSide.end(), stripID) != clusterStrips_vSide.end()) {
180  std::string sensorName;
181  occupancyPDFName(svdShaperDigit.getSensorID(), svdShaperDigit.isUStrip(), sensorName);
182  auto xTalkStrip = m_sensorHistograms.at(sensorName);
183  xTalkStrip->Fill(svdShaperDigit.getCellID(), 1);
184  }
185 
186 
187  }//shaper digit loop
188  }
189 } //Collector loop
190 
191 
193 {}
194 
196 {
197  std::string sensorName;
199  for (auto& layers : geo.getLayers(VXD::SensorInfoBase::SVD)) {
200  for (auto& ladders : geo.getLadders(layers)) {
201  for (auto& sensors : geo.getSensors(ladders)) {
202  for (int side = 0; side <= 1; side++) {
203 
204  occupancyPDFName(sensors, side, sensorName);
205  auto sensorOnMap = m_sensorHistograms.at(sensorName);
206  m_hist = sensorOnMap;
207  m_layer = sensors.getLayerNumber();
208  m_ladder = sensors.getLadderNumber();
209  m_sensor = sensors.getSensorNumber();
210  m_side = side;
211 
212  getObjectPtr<TTree>("HTreeCrossTalkCalib")->Fill();
213 
214  sensorOnMap->Delete();
215  }
216  }
217  }
218  }
219 
220 }
221 
222 
223 void SVDCrossTalkCalibrationsCollectorModule::calculateAverage(const VxdID& sensorID, double& mean, int side)
224 {
225  double nBins = 0;
226  if (side == 1) nBins = 768; //U-side 768 channels
227  else nBins = 512;
228  double count = 0;
229  for (int i = 0; i < nBins; i++) {
230  count += m_OccupancyCal.getOccupancy(sensorID, side, i);
231  }
232  mean = count / nBins;
233 }
234 
Belle2::CalibrationCollectorModule
Calibration collector module base class.
Definition: CalibrationCollectorModule.h:44
Belle2::SVDCrossTalkCalibrationsCollectorModule::collect
void collect() override final
Event processing.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:95
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_sensor
int m_sensor
Number of sensors.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:72
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_nAPVFactor
int m_nAPVFactor
Parameter to set number of sensors with possible cross-talk clusters required for event flagging.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:88
Belle2::SVDCrossTalkCalibrationsCollectorModule::calculateAverage
void calculateAverage(const VxdID &sensorID, double &mean, int side)
Function to calculate sensor average occupancy.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:223
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition: Module.h:82
Belle2::SVDCrossTalkCalibrationsCollectorModule::prepare
void prepare() override final
Init the module.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:48
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_svdShaperDigitsName
std::string m_svdShaperDigitsName
SVDShaperDigit collection name.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:79
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_svdShaperDigits
StoreArray< SVDShaperDigit > m_svdShaperDigits
The storeArray for svdShaperDigits.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:82
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_uSideOccupancyFactor
int m_uSideOccupancyFactor
Parameter to define high occupancy strips (some multiple above sensor average occupancy)
Definition: SVDCrossTalkCalibrationsCollectorModule.h:84
Belle2::VXD::GeoCache::getLayers
const std::set< Belle2::VxdID > getLayers(SensorInfoBase::SensorType sensortype=SensorInfoBase::VXD)
Return a set of all known Layers.
Definition: GeoCache.cc:177
Belle2::VXD::GeoCache::getSensors
const std::set< Belle2::VxdID > & getSensors(Belle2::VxdID ladder) const
Return a set of all sensor IDs belonging to a given ladder.
Definition: GeoCache.cc:205
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::SVDOccupancyCalibrations::getOccupancy
float getOccupancy(const VxdID &sensorID, const bool &isU, const unsigned short &strip) const
This is the method for getting the occupancy, or the deviation from the average, still to be decided.
Definition: SVDOccupancyCalibrations.h:71
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::SVDCrossTalkCalibrationsCollectorModule::finish
void finish() override
Termination.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:192
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_vSideOccupancyFactor
int m_vSideOccupancyFactor
Parameter to define high occupancy strips (some multiple above sensor average occupancy)
Definition: SVDCrossTalkCalibrationsCollectorModule.h:86
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_layer
int m_layer
Number of layers to define size of TTree.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:70
Belle2::SVDCrossTalkCalibrationsCollectorModule::startRun
void startRun() override final
New run.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:66
Belle2::occupancyPDFName
void occupancyPDFName(const VxdID &sensor, int side, std::string &PDFName)
Function to maintain common naming convention between calibration occupancy file generation and occup...
Definition: SVDCrossTalkFinderHelperFunctions.h:37
Belle2::VXD::SensorInfoBase::SVD
@ SVD
SVD Sensor.
Definition: SensorInfoBase.h:45
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_hist
TH1F * m_hist
Initialisation of crosstalk histogram.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:69
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_OccupancyCal
SVDOccupancyCalibrations m_OccupancyCal
SVDOccupancy calibrations db object.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:92
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_side
int m_side
Number of sides.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:73
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_sensorHistograms
std::map< std::string, TH1F * > m_sensorHistograms
map to store cross-talk strip histograms
Definition: SVDCrossTalkCalibrationsCollectorModule.h:90
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_ladder
int m_ladder
Number of ladders.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:71
Belle2::VXD::GeoCache::getLadders
const std::set< Belle2::VxdID > & getLadders(Belle2::VxdID layer) const
Return a set of all ladder IDs belonging to a given layer.
Definition: GeoCache.cc:194
Belle2::SVDCrossTalkCalibrationsCollectorModule::closeRun
void closeRun() override final
End of run.
Definition: SVDCrossTalkCalibrationsCollectorModule.cc:195
Belle2::SVDCrossTalkCalibrationsCollectorModule::m_histogramTree
TTree * m_histogramTree
Initialisation of TTree object.
Definition: SVDCrossTalkCalibrationsCollectorModule.h:68