Belle II Software  release-06-02-00
VXDDQMExpressRecoModule.cc
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 
9 #include "vxd/modules/vxdDQM/VXDDQMExpressRecoModule.h"
10 
11 #include <framework/datastore/StoreArray.h>
12 #include <framework/datastore/RelationArray.h>
13 
14 #include <svd/dataobjects/SVDShaperDigit.h>
15 #include <svd/dataobjects/SVDCluster.h>
16 #include <svd/geometry/SensorInfo.h>
17 
18 #include <pxd/dataobjects/PXDDigit.h>
19 #include <pxd/dataobjects/PXDCluster.h>
20 #include <pxd/geometry/SensorInfo.h>
21 
22 #include <vxd/geometry/SensorInfoBase.h>
23 #include <vxd/geometry/GeoTools.h>
24 
25 #include <TMath.h>
26 
27 #include <boost/format.hpp>
28 
29 #include "TDirectory.h"
30 
31 using namespace std;
32 using boost::format;
33 using namespace Belle2;
34 
35 //-----------------------------------------------------------------
36 // Register the Module
37 //-----------------------------------------------------------------
38 REG_MODULE(VXDDQMExpressReco)
39 
40 
41 //-----------------------------------------------------------------
42 // Implementation
43 //-----------------------------------------------------------------
44 
45 VXDDQMExpressRecoModule::VXDDQMExpressRecoModule() : HistoModule()
46 {
47  //Set module properties
48  setDescription("VXD DQM module for Express Reco "
49  "Recommended Number of events for monitor is 40 kEvents or more to fill all histograms "
50  );
51 
52  setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
53  addParam("UseDigits", m_UseDigits,
54  "flag <0,1> for using digits only, no cluster information will be required, default = 0 ", m_UseDigits);
55  addParam("CorrelationGranulation", m_CorrelationGranulation,
56  "Set granulation of histogram plots, default is 1 degree, min = 0.02 degree, max = 1 degree ", m_CorrelationGranulation);
57  addParam("SwapPXD", m_SwapPXD, "flag <0,1> very special case for swap of phi-theta coordinates, default = 0 ", m_SwapPXD);
58  addParam("CutCorrelationSigPXD", m_CutCorrelationSigPXD,
59  "Cut threshold of PXD signal for accepting to correlations, default = 0 ADU ", m_CutCorrelationSigPXD);
60  addParam("CutCorrelationSigUSVD", m_CutCorrelationSigUSVD,
61  "Cut threshold of SVD signal for accepting to correlations in u, default = 0 ADU ", m_CutCorrelationSigUSVD);
62  addParam("CutCorrelationSigVSVD", m_CutCorrelationSigVSVD,
63  "Cut threshold of SVD signal for accepting to correlations in v, default = 0 ADU ", m_CutCorrelationSigVSVD);
64  addParam("CutCorrelationTimeSVD", m_CutCorrelationTimeSVD,
65  "Cut threshold of SVD time window for accepting to correlations, default = 70 ns ", m_CutCorrelationTimeSVD);
66 
67  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
68  std::string("VXDExpReco"));
69 }
70 
71 
72 VXDDQMExpressRecoModule::~VXDDQMExpressRecoModule()
73 {
74 }
75 
76 //------------------------------------------------------------------
77 // Function to define histograms
78 //-----------------------------------------------------------------
79 
80 void VXDDQMExpressRecoModule::defineHisto()
81 {
82  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
83  if (gTools->getNumberOfLayers() == 0) {
84  B2WARNING("Missing geometry for VXD, check steering file.");
85  return;
86  }
87  if (gTools->getNumberOfPXDLayers() == 0) {
88  B2WARNING("Missing geometry for PXD.");
89  }
90  if (gTools->getNumberOfSVDLayers() == 0) {
91  B2WARNING("Missing geometry for SVD.");
92  }
93 
94  // Create a separate histogram directories and cd into it.
95  TDirectory* oldDir = gDirectory;
96  if (m_histogramDirectoryName != "") {
97  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
98  oldDir->cd(m_histogramDirectoryName.c_str());
99  }
100 
101  if (m_CorrelationGranulation > 5.0) m_CorrelationGranulation = 5.0; // set maximum of gramularity to 1 degree
102  if (m_CorrelationGranulation < 0.02) m_CorrelationGranulation = 0.02; // set minimum of gramularity to 0.02 degree
103 
104  // if binning go over h_MaxBins it decrease preset of range
105  int h_MaxBins = 2000; //maximal size of histogram binning:
106 
107  VXD::GeoCache& geo = VXD::GeoCache::getInstance();
108 
109  // basic constants presets:
110  int nVXDLayers = gTools->getNumberOfLayers();
111 
112  // Create basic histograms:
113  m_correlationsSP = (TH2F**) new TH2F*[nVXDLayers * nVXDLayers];
114  m_correlationsSP1DPhi = (TH1F**) new TH1F*[nVXDLayers * nVXDLayers];
115  m_correlationsSP1DTheta = (TH1F**) new TH1F*[nVXDLayers * nVXDLayers];
116 
117 
118  for (int i = 0; i < nVXDLayers; i++) {
119  for (int j = 0; j < nVXDLayers; j++) {
120  m_correlationsSP[nVXDLayers * j + i] = NULL;
121  m_correlationsSP1DPhi[nVXDLayers * j + i] = NULL;
122  m_correlationsSP1DTheta[nVXDLayers * j + i] = NULL;
123  }
124  }
125 
126  string Diru = str(format("Phi"));
127  string Dirv = str(format("Theta"));
128  string Unit = str(format("degree"));
129  string AxisTitPhi = str(format("%1% position [%2%]") % Diru.c_str() % Unit.c_str()); // 0..360 degree, from x+ to y+...
130  string AxisTitTheta = str(format("%1% position [%2%]") % Dirv.c_str() % Unit.c_str()); // 0..180 degree, z+ to z-
131  if (m_UseDigits) {
132  AxisTitPhi = str(format("From digits: %1%") % AxisTitPhi);
133  AxisTitTheta = str(format("From digits: %1%") % AxisTitTheta);
134  }
135  for (VxdID layer1 : geo.getLayers()) {
136  int i = layer1.getLayerNumber() - gTools->getFirstLayer();
137  float uSize1s, uSize1e;
138  float vSize1s, vSize1e;
139  int nStripsU1, nStripsV1;
140  uSize1s = -180.0; // degree
141  uSize1e = 180.0; // degree
142  vSize1s = 0.0; // degree
143  vSize1e = 180.0; // degree
144  nStripsU1 = (uSize1e - uSize1s) / m_CorrelationGranulation;
145  if (nStripsU1 > h_MaxBins) nStripsU1 = h_MaxBins;
146  nStripsV1 = nStripsU1;
147  for (VxdID layer2 : geo.getLayers()) {
148  int j = layer2.getLayerNumber() - gTools->getFirstLayer();
149  float uSize2s, uSize2e;
150  float vSize2s, vSize2e;
151  int nStripsU2, nStripsV2;
152  uSize2s = -180.0; // degree
153  uSize2e = 180.0; // degree
154  vSize2s = 0.0; // degree
155  vSize2e = 180.0; // degree
156  nStripsU2 = nStripsU1;
157  nStripsV2 = nStripsV1;
158 
159  if (i == j) { // hit maps
160  string nameSP = str(format("DQMER_VXD_Hitmap_L%1%") % layer2.getLayerNumber());
161  string titleSP = str(format("DQM ER VXD Hitmap, layer %1%") % layer2.getLayerNumber());
162  m_correlationsSP[nVXDLayers * j + i] = new TH2F(nameSP.c_str(), titleSP.c_str(),
163  nStripsU2, uSize2s, uSize2e,
164  nStripsV2, vSize2s, vSize2e);
165  m_correlationsSP[nVXDLayers * j + i]->GetXaxis()->SetTitle(AxisTitPhi.c_str());
166  m_correlationsSP[nVXDLayers * j + i]->GetYaxis()->SetTitle(AxisTitTheta.c_str());
167  m_correlationsSP[nVXDLayers * j + i]->GetZaxis()->SetTitle("hits");
168 
169  nameSP = str(format("DQMER_VXD_Hitmap_%1%_L%2%") % Diru.c_str() % layer2.getLayerNumber());
170  titleSP = str(format("DQM ER VXD Hitmap in %1%, layer %2%") % Diru.c_str() % layer2.getLayerNumber());
171  m_correlationsSP1DPhi[nVXDLayers * j + i] = new TH1F(nameSP.c_str(), titleSP.c_str(),
172  nStripsU2, uSize2s, uSize2e);
173  m_correlationsSP1DPhi[nVXDLayers * j + i]->GetXaxis()->SetTitle(AxisTitPhi.c_str());
174  m_correlationsSP1DPhi[nVXDLayers * j + i]->GetYaxis()->SetTitle("hits");
175 
176  nameSP = str(format("DQMER_VXD_Hitmap_%1%_L%2%") % Dirv.c_str() % layer2.getLayerNumber());
177  titleSP = str(format("DQM ER VXD Hitmap in %1%, layer %2%") % Dirv.c_str() % layer2.getLayerNumber());
178  m_correlationsSP1DTheta[nVXDLayers * j + i] = new TH1F(nameSP.c_str(), titleSP.c_str(),
179  nStripsV2, vSize2s, vSize2e);
180  m_correlationsSP1DTheta[nVXDLayers * j + i]->GetXaxis()->SetTitle(AxisTitTheta.c_str());
181  m_correlationsSP1DTheta[nVXDLayers * j + i]->GetYaxis()->SetTitle("hits");
182 
183  } else if (i < j) { // correlations for Phi
184  string nameSP = str(format("DQMER_VXD_Correlations_%1%_L%2%_L%3%") % Diru.c_str() % layer1.getLayerNumber() %
185  layer2.getLayerNumber());
186  string titleSP = str(format("DQM ER VXD Correlations in %1%, layers %2% %3%") % Diru.c_str() % layer1.getLayerNumber() %
187  layer2.getLayerNumber());
188  m_correlationsSP[nVXDLayers * j + i] = new TH2F(nameSP.c_str(), titleSP.c_str(),
189  nStripsU1, uSize1s, uSize1e,
190  nStripsU2, uSize2s, uSize2e);
191  string axisxtitle = str(format("%1%, layer %2%") % AxisTitPhi.c_str() % layer1.getLayerNumber());
192  string axisytitle = str(format("%1%, layer %2%") % AxisTitPhi.c_str() % layer2.getLayerNumber());
193  m_correlationsSP[nVXDLayers * j + i]->GetXaxis()->SetTitle(axisxtitle.c_str());
194  m_correlationsSP[nVXDLayers * j + i]->GetYaxis()->SetTitle(axisytitle.c_str());
195  m_correlationsSP[nVXDLayers * j + i]->GetZaxis()->SetTitle("hits");
196 
197  nameSP = str(format("DQMER_VXD_1D_Correlations_%1%_L%2%_L%3%") % Diru.c_str() % layer1.getLayerNumber() % layer2.getLayerNumber());
198  titleSP = str(format("DQM ER VXD 1D Correlations in %1%, layers %2% %3%") % Diru.c_str() % layer1.getLayerNumber() %
199  layer2.getLayerNumber());
200  m_correlationsSP1DPhi[nVXDLayers * j + i] = new TH1F(nameSP.c_str(), titleSP.c_str(),
201  nStripsU1, uSize1s, uSize1e);
202  axisxtitle = str(format("%1%, layer %2% - %3%") % AxisTitPhi.c_str() % layer1.getLayerNumber() % layer2.getLayerNumber());
203  m_correlationsSP1DPhi[nVXDLayers * j + i]->GetXaxis()->SetTitle(axisxtitle.c_str());
204  m_correlationsSP1DPhi[nVXDLayers * j + i]->GetYaxis()->SetTitle("hits");
205  } else { // correlations for Theta
206  string nameSP = str(format("DQMER_VXD_Correlations_%1%_L%2%_L%3%") % Dirv.c_str() % layer1.getLayerNumber() %
207  layer2.getLayerNumber());
208  string titleSP = str(format("DQM ER VXD Correlations in %1%, layers %2% %3%") % Dirv.c_str() % layer1.getLayerNumber() %
209  layer2.getLayerNumber());
210  m_correlationsSP[nVXDLayers * j + i] = new TH2F(nameSP.c_str(), titleSP.c_str(),
211  nStripsV1, vSize1s, vSize1e,
212  nStripsV2, vSize2s, vSize2e);
213  string axisxtitle = str(format("%1%, layer %2%") % AxisTitTheta.c_str() % layer1.getLayerNumber());
214  string axisytitle = str(format("%1%, layer %2%") % AxisTitTheta.c_str() % layer2.getLayerNumber());
215  m_correlationsSP[nVXDLayers * j + i]->GetXaxis()->SetTitle(axisxtitle.c_str());
216  m_correlationsSP[nVXDLayers * j + i]->GetYaxis()->SetTitle(axisytitle.c_str());
217  m_correlationsSP[nVXDLayers * j + i]->GetZaxis()->SetTitle("hits");
218 
219  nameSP = str(format("DQMER_VXD_1D_Correlations_%1%_L%2%_L%3%") % Dirv.c_str() % layer1.getLayerNumber() % layer2.getLayerNumber());
220  titleSP = str(format("DQM ER VXD 1D Correlations in %1%, layers %2% %3%") % Dirv.c_str() % layer1.getLayerNumber() %
221  layer2.getLayerNumber());
222  m_correlationsSP1DTheta[nVXDLayers * j + i] = new TH1F(nameSP.c_str(), titleSP.c_str(),
223  nStripsV1, -vSize1e, vSize1e);
224  axisxtitle = str(format("%1%, layer %2% - %3%") % AxisTitTheta.c_str() % layer1.getLayerNumber() % layer2.getLayerNumber());
225  m_correlationsSP1DTheta[nVXDLayers * j + i]->GetXaxis()->SetTitle(axisxtitle.c_str());
226  m_correlationsSP1DTheta[nVXDLayers * j + i]->GetYaxis()->SetTitle("hits");
227  }
228  }
229  }
230 
231  oldDir->cd();
232 }
233 
234 
235 void VXDDQMExpressRecoModule::initialize()
236 {
237  // Register histograms (calls back defineHisto)
238  REG_HISTOGRAM
239 
240  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
241  if (gTools->getNumberOfLayers() != 0) {
242  //Register collections
243  StoreArray<PXDDigit> storePXDDigits(m_storePXDDigitsName);
244  StoreArray<SVDShaperDigit> storeSVDShaperDigits(m_storeSVDShaperDigitsName);
245  StoreArray<PXDCluster> storePXDClusters(m_storePXDClustersName);
246  StoreArray<SVDCluster> storeSVDClusters(m_storeSVDClustersName);
247  RelationArray relPXDClusterDigits(storePXDClusters, storePXDDigits);
248  RelationArray relSVDClusterDigits(storeSVDClusters, storeSVDShaperDigits);
249  m_storePXDClustersName = storePXDClusters.getName();
250  m_relPXDClusterDigitName = relPXDClusterDigits.getName();
251  m_storeSVDClustersName = storeSVDClusters.getName();
252  m_relSVDClusterDigitName = relSVDClusterDigits.getName();
253 
254  storePXDDigits.isRequired();
255  storeSVDShaperDigits.isRequired();
256 
257  //Store names to speed up creation later
258  m_storePXDDigitsName = storePXDDigits.getName();
259  m_storeSVDShaperDigitsName = storeSVDShaperDigits.getName();
260  }
261 }
262 
263 void VXDDQMExpressRecoModule::beginRun()
264 {
265  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
266  if (gTools->getNumberOfLayers() == 0) return;
267 
268  // Just to make sure, reset all the histograms.
269  int nVXDLayers = gTools->getNumberOfLayers();
270  for (int i = 0; i < nVXDLayers; i++) {
271  for (int j = 0; j < nVXDLayers; j++) {
272  if (m_correlationsSP[nVXDLayers * j + i] != NULL) m_correlationsSP[nVXDLayers * j + i]->Reset();
273  if (m_correlationsSP1DPhi[nVXDLayers * j + i] != NULL) m_correlationsSP1DPhi[nVXDLayers * j + i]->Reset();
274  if (m_correlationsSP1DTheta[nVXDLayers * j + i] != NULL) m_correlationsSP1DTheta[nVXDLayers * j + i]->Reset();
275  }
276  }
277 }
278 
279 
280 void VXDDQMExpressRecoModule::event()
281 {
282  auto gTools = VXD::GeoCache::getInstance().getGeoTools();
283  if (gTools->getNumberOfLayers() == 0) return;
284 
285  const StoreArray<PXDDigit> storePXDDigits(m_storePXDDigitsName);
286  const StoreArray<SVDShaperDigit> storeSVDShaperDigits(m_storeSVDShaperDigitsName);
287 
288  const StoreArray<SVDCluster> storeSVDClusters(m_storeSVDClustersName);
289  const StoreArray<PXDCluster> storePXDClusters(m_storePXDClustersName);
290 
291  const RelationArray relPXDClusterDigits(storePXDClusters, storePXDDigits, m_relPXDClusterDigitName);
292  const RelationArray relSVDClusterDigits(storeSVDClusters, storeSVDShaperDigits, m_relSVDClusterDigitName);
293 
294  // If there are no digits, leave
295  if (!storePXDDigits && !storeSVDShaperDigits) return;
296 
297  // basic constants presets:
298  int nVXDLayers = gTools->getNumberOfLayers();
299  int firstVXDLayer = gTools->getFirstLayer();
300  int firstPXDLayer = gTools->getFirstPXDLayer();
301  int lastPXDLayer = gTools->getLastPXDLayer();
302  int firstSVDLayer = gTools->getFirstSVDLayer();
303  int lastSVDLayer = gTools->getLastSVDLayer();
304 
305  int MaxHits = 0;
306  if (m_UseDigits) {
307  MaxHits = storeSVDShaperDigits.getEntries() + storePXDDigits.getEntries();
308  } else {
309  MaxHits = storeSVDClusters.getEntries() + storePXDClusters.getEntries();
310  }
311 
312  // If there are no hits, leave
313  if (MaxHits == 0) return;
314 
315  for (int i1 = 0; i1 < MaxHits; i1++) {
316  // preparing of first value for correlation plots with postfix "1":
317  float fTime1 = 0.0;
318  float fPosSPU1 = 0.0;
319  float fPosSPV1 = 0.0;
320  int iIsPXD1 = 0;
321  int index1 = 0;
322  int iIsU1 = 0;
323  int iIsV1 = 0;
324  int iLayer1 = 0;
325  if (m_UseDigits) { // DIGITS:
326  if (i1 < storePXDDigits.getEntries()) { // PXD digits/clusters:
327  const PXDDigit& digitPXD1 = *storePXDDigits[i1];
328  iLayer1 = digitPXD1.getSensorID().getLayerNumber();
329  if ((iLayer1 < firstPXDLayer) || (iLayer1 > lastPXDLayer)) continue;
330  index1 = iLayer1 - firstVXDLayer;
331  float fCharge1 = digitPXD1.getCharge();
332  if (fCharge1 < m_CutCorrelationSigPXD) continue;
333  VxdID sensorID1 = digitPXD1.getSensorID();
334  auto info = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(sensorID1));
335  TVector3 rLocal1(info.getUCellPosition(digitPXD1.getUCellID()), info.getVCellPosition(digitPXD1.getVCellID()), 0);
336  TVector3 ral1 = info.pointToGlobal(rLocal1);
337  iIsPXD1 = 1;
338  iIsU1 = 1;
339  iIsV1 = 1;
340  fPosSPU1 = ral1.Phi() / TMath::Pi() * 180;
341  fPosSPV1 = ral1.Theta() / TMath::Pi() * 180;
342  if (m_SwapPXD) {
343  fPosSPV1 = ral1.Phi() / TMath::Pi() * 180;
344  fPosSPU1 = ral1.Theta() / TMath::Pi() * 180;
345  }
346  } else { // SVD digits/clusters:
347  const SVDShaperDigit& digitSVD1 = *storeSVDShaperDigits[i1 - storePXDDigits.getEntries()];
348  iLayer1 = digitSVD1.getSensorID().getLayerNumber();
349  if ((iLayer1 < firstSVDLayer) || (iLayer1 > lastSVDLayer)) continue;
350  index1 = iLayer1 - firstVXDLayer;
351  fTime1 = digitSVD1.getFADCTime();
352  VxdID sensorID1 = digitSVD1.getSensorID();
353  auto info = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(sensorID1));
354  SVDShaperDigit::APVFloatSamples samples = digitSVD1.getSamples();
355  if (digitSVD1.isUStrip()) {
356  int iCont = 0;
357  for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
358  float fCharge1 = samples[i];
359  if (fCharge1 > m_CutCorrelationSigUSVD) iCont = 1;
360  }
361  if (iCont == 0) continue;
362  float possi = info.getUCellPosition(digitSVD1.getCellID());
363  TVector3 rLocal1(possi, 0 , 0);
364  TVector3 ral1 = info.pointToGlobal(rLocal1);
365  iIsU1 = 1;
366  fPosSPU1 = ral1.Phi() / TMath::Pi() * 180;
367  } else {
368  int iCont = 0;
369  for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
370  float fCharge1 = samples[i];
371  if (fCharge1 > m_CutCorrelationSigVSVD) iCont = 1;
372  }
373  if (iCont == 0) continue;
374 
375  //float possi = digitSVD1.getCellPosition(); // is not work anymore
376  float possi = info.getVCellPosition(digitSVD1.getCellID());
377 
378  TVector3 rLocal1(0, possi, 0);
379  TVector3 ral1 = info.pointToGlobal(rLocal1);
380  iIsV1 = 1;
381  fPosSPV1 = ral1.Theta() / TMath::Pi() * 180;
382  }
383  }
384  } else { // CLUSTERS:
385  if (i1 < storePXDClusters.getEntries()) { // PXD digits/clusters:
386  const PXDCluster& clusterPXD1 = *storePXDClusters[i1];
387  iLayer1 = clusterPXD1.getSensorID().getLayerNumber();
388  if ((iLayer1 < firstPXDLayer) || (iLayer1 > lastPXDLayer)) continue;
389  index1 = iLayer1 - firstVXDLayer;
390  float fCharge1 = clusterPXD1.getCharge();
391  if (fCharge1 < m_CutCorrelationSigPXD) continue;
392  VxdID sensorID1 = clusterPXD1.getSensorID();
393  auto info = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(sensorID1));
394  TVector3 rLocal1(clusterPXD1.getU(), clusterPXD1.getV(), 0);
395  TVector3 ral1 = info.pointToGlobal(rLocal1);
396  iIsPXD1 = 1;
397  iIsU1 = 1;
398  iIsV1 = 1;
399  fPosSPU1 = ral1.Phi() / TMath::Pi() * 180;
400  fPosSPV1 = ral1.Theta() / TMath::Pi() * 180;
401  if (m_SwapPXD) {
402  fPosSPV1 = ral1.Phi() / TMath::Pi() * 180;
403  fPosSPU1 = ral1.Theta() / TMath::Pi() * 180;
404  }
405  } else { // SVD digits/clusters:
406  const SVDCluster& clusterSVD1 = *storeSVDClusters[i1 - storePXDClusters.getEntries()];
407  iLayer1 = clusterSVD1.getSensorID().getLayerNumber();
408  if ((iLayer1 < firstSVDLayer) || (iLayer1 > lastSVDLayer)) continue;
409  index1 = iLayer1 - firstVXDLayer;
410  float fCharge1 = clusterSVD1.getCharge();
411  fTime1 = clusterSVD1.getClsTime();
412  VxdID sensorID1 = clusterSVD1.getSensorID();
413  auto info = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(sensorID1));
414  if (clusterSVD1.isUCluster()) {
415  if (fCharge1 < m_CutCorrelationSigUSVD * 200) continue; // in electrons
416  TVector3 rLocal1(clusterSVD1.getPosition(), 0 , 0);
417  TVector3 ral1 = info.pointToGlobal(rLocal1);
418  iIsU1 = 1;
419  fPosSPU1 = ral1.Phi() / TMath::Pi() * 180;
420  } else {
421  if (fCharge1 < m_CutCorrelationSigVSVD * 200) continue; // in electrons
422  TVector3 rLocal1(0, clusterSVD1.getPosition(), 0);
423  TVector3 ral1 = info.pointToGlobal(rLocal1);
424  iIsV1 = 1;
425  fPosSPV1 = ral1.Theta() / TMath::Pi() * 180;
426  }
427  }
428  }
429  // hit maps for PXD:
430  if ((iIsU1 == 1) && (iIsV1 == 1)) {
431  if (m_correlationsSP[nVXDLayers * index1 + index1] != NULL)
432  m_correlationsSP[nVXDLayers * index1 + index1]->Fill(fPosSPU1, fPosSPV1);
433  if (m_correlationsSP1DPhi[nVXDLayers * index1 + index1] != NULL)
434  m_correlationsSP1DPhi[nVXDLayers * index1 + index1]->Fill(fPosSPU1);
435  if (m_correlationsSP1DTheta[nVXDLayers * index1 + index1] != NULL)
436  m_correlationsSP1DTheta[nVXDLayers * index1 + index1]->Fill(fPosSPV1);
437  }
438  for (int i2 = 0; i2 < MaxHits; i2++) {
439  // preparing of second value for correlation plots with postfix "2":
440  float fTime2 = 0.0;
441  float fPosSPU2 = 0.0;
442  float fPosSPV2 = 0.0;
443  int iIsPXD2 = 0;
444  int index2 = 0;
445  int iIsU2 = 0;
446  int iIsV2 = 0;
447  int iLayer2 = 0;
448  if (m_UseDigits) { // DIGITS:
449  if (i2 < storePXDDigits.getEntries()) { // PXD digits/clusters:
450  const PXDDigit& digitPXD2 = *storePXDDigits[i2];
451  iLayer2 = digitPXD2.getSensorID().getLayerNumber();
452  if ((iLayer2 < firstPXDLayer) || (iLayer2 > lastPXDLayer)) continue;
453  index2 = iLayer2 - firstVXDLayer;
454  float fCharge2 = digitPXD2.getCharge();
455  if (fCharge2 < m_CutCorrelationSigPXD) continue;
456  VxdID sensorID2 = digitPXD2.getSensorID();
457  auto info = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(sensorID2));
458  TVector3 rLocal2(info.getUCellPosition(digitPXD2.getUCellID()), info.getVCellPosition(digitPXD2.getVCellID()), 0);
459  TVector3 ral2 = info.pointToGlobal(rLocal2);
460  iIsPXD2 = 1;
461  iIsU2 = 1;
462  iIsV2 = 1;
463  fPosSPU2 = ral2.Phi() / TMath::Pi() * 180;
464  fPosSPV2 = ral2.Theta() / TMath::Pi() * 180;
465  if (m_SwapPXD) {
466  fPosSPV2 = ral2.Phi() / TMath::Pi() * 180;
467  fPosSPU2 = ral2.Theta() / TMath::Pi() * 180;
468  }
469  } else { // SVD digits/clusters:
470  const SVDShaperDigit& digitSVD2 = *storeSVDShaperDigits[i2 - storePXDDigits.getEntries()];
471  iLayer2 = digitSVD2.getSensorID().getLayerNumber();
472  if ((iLayer2 < firstSVDLayer) || (iLayer2 > lastSVDLayer)) continue;
473  index2 = iLayer2 - firstVXDLayer;
474  fTime2 = digitSVD2.getFADCTime();
475  VxdID sensorID2 = digitSVD2.getSensorID();
476  auto info = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(sensorID2));
477 
478  SVDShaperDigit::APVFloatSamples samples = digitSVD2.getSamples();
479  if (digitSVD2.isUStrip()) {
480  int iCont = 0;
481  for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
482  float fCharge2 = samples[i];
483  if (fCharge2 > m_CutCorrelationSigUSVD) iCont = 1;
484  }
485  if (iCont == 0) continue;
486  float possi = info.getUCellPosition(digitSVD2.getCellID());
487  TVector3 rLocal2(possi, 0 , 0);
488  TVector3 ral2 = info.pointToGlobal(rLocal2);
489  iIsU2 = 1;
490  fPosSPU2 = ral2.Phi() / TMath::Pi() * 180;
491  } else {
492  int iCont = 0;
493  for (size_t i = 0; i < SVDShaperDigit::c_nAPVSamples; ++i) {
494  float fCharge2 = samples[i];
495  if (fCharge2 > m_CutCorrelationSigVSVD) iCont = 1;
496  }
497  if (iCont == 0) continue;
498  float possi = info.getVCellPosition(digitSVD2.getCellID());
499  TVector3 rLocal2(0, possi, 0);
500  TVector3 ral2 = info.pointToGlobal(rLocal2);
501  iIsV2 = 1;
502  fPosSPV2 = ral2.Theta() / TMath::Pi() * 180;
503  }
504  }
505  } else { // CLUSTERS:
506  if (i2 < storePXDClusters.getEntries()) { // PXD digits/clusters:
507  const PXDCluster& clusterPXD2 = *storePXDClusters[i2];
508  iLayer2 = clusterPXD2.getSensorID().getLayerNumber();
509  if ((iLayer2 < firstPXDLayer) || (iLayer2 > lastPXDLayer)) continue;
510  index2 = iLayer2 - firstVXDLayer;
511  float fCharge2 = clusterPXD2.getCharge();
512  if (fCharge2 < m_CutCorrelationSigPXD) continue;
513  VxdID sensorID2 = clusterPXD2.getSensorID();
514  auto info = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::get(sensorID2));
515  TVector3 rLocal2(clusterPXD2.getU(), clusterPXD2.getV(), 0);
516  TVector3 ral2 = info.pointToGlobal(rLocal2);
517  iIsPXD2 = 1;
518  iIsU2 = 1;
519  iIsV2 = 1;
520  fPosSPU2 = ral2.Phi() / TMath::Pi() * 180;
521  fPosSPV2 = ral2.Theta() / TMath::Pi() * 180;
522  if (m_SwapPXD) {
523  fPosSPV2 = ral2.Phi() / TMath::Pi() * 180;
524  fPosSPU2 = ral2.Theta() / TMath::Pi() * 180;
525  }
526  } else { // SVD digits/clusters:
527  const SVDCluster& clusterSVD2 = *storeSVDClusters[i2 - storePXDClusters.getEntries()];
528  iLayer2 = clusterSVD2.getSensorID().getLayerNumber();
529  if ((iLayer2 < firstSVDLayer) || (iLayer2 > lastSVDLayer)) continue;
530  index2 = iLayer2 - firstVXDLayer;
531  float fCharge2 = clusterSVD2.getCharge();
532  fTime2 = clusterSVD2.getClsTime();
533  VxdID sensorID2 = clusterSVD2.getSensorID();
534  auto info = dynamic_cast<const SVD::SensorInfo&>(VXD::GeoCache::get(sensorID2));
535  if (clusterSVD2.isUCluster()) {
536  if (fCharge2 < m_CutCorrelationSigUSVD * 200) continue; // in electrons
537  TVector3 rLocal2(clusterSVD2.getPosition(), 0 , 0);
538  TVector3 ral2 = info.pointToGlobal(rLocal2);
539  iIsU2 = 1;
540  fPosSPU2 = ral2.Phi() / TMath::Pi() * 180;
541  } else {
542  if (fCharge2 < m_CutCorrelationSigVSVD * 200) continue; // in electrons
543  TVector3 rLocal2(0, clusterSVD2.getPosition(), 0);
544  TVector3 ral2 = info.pointToGlobal(rLocal2);
545  iIsV2 = 1;
546  fPosSPV2 = ral2.Theta() / TMath::Pi() * 180;
547  }
548  }
549  }
550  if ((iIsPXD1 == 0) && (iIsPXD2 == 0))
551  if ((fabs(fTime1 - fTime2)) > m_CutCorrelationTimeSVD) continue;
552  // ready to fill correlation histograms and hit maps:
553  if ((index1 == index2) && (iIsU1 == 1) && (iIsV2 == 1) && (iIsPXD1 == 0) && (iIsPXD2 == 0)) {
554  // hit maps for SVD:
555  if (m_correlationsSP[nVXDLayers * index2 + index1] != NULL)
556  m_correlationsSP[nVXDLayers * index2 + index1]->Fill(fPosSPU1, fPosSPV2);
557  if (m_correlationsSP1DPhi[nVXDLayers * index2 + index1] != NULL)
558  m_correlationsSP1DPhi[nVXDLayers * index2 + index1]->Fill(fPosSPU1);
559  if (m_correlationsSP1DTheta[nVXDLayers * index2 + index1] != NULL)
560  m_correlationsSP1DTheta[nVXDLayers * index2 + index1]->Fill(fPosSPV2);
561  } else if ((index1 < index2) && (iIsU1 == iIsU2) && (iIsU1 == 1)) {
562  // correlations for u
563  if (m_correlationsSP[nVXDLayers * index2 + index1] != NULL)
564  m_correlationsSP[nVXDLayers * index2 + index1]->Fill(fPosSPU1, fPosSPU2);
565  if (m_correlationsSP1DPhi[nVXDLayers * index2 + index1] != NULL)
566  m_correlationsSP1DPhi[nVXDLayers * index2 + index1]->Fill(fPosSPU2 - fPosSPU1);
567  } else if ((index1 > index2) && (iIsV1 == iIsV2) && (iIsV1 == 1)) {
568  // correlations for v
569  if (m_correlationsSP[nVXDLayers * index2 + index1] != NULL)
570  m_correlationsSP[nVXDLayers * index2 + index1]->Fill(fPosSPV2, fPosSPV1);
571  if (m_correlationsSP1DTheta[nVXDLayers * index2 + index1] != NULL)
572  m_correlationsSP1DTheta[nVXDLayers * index2 + index1]->Fill(fPosSPV2 - fPosSPV1);
573  }
574  }
575  }
576 
577 }
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
Belle2::PXDCluster
The PXD Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition: PXDCluster.h:30
Belle2::PXDCluster::getV
float getV() const
Get v coordinate of hit position.
Definition: PXDCluster.h:136
Belle2::PXDCluster::getCharge
unsigned short getCharge() const
Get collected charge.
Definition: PXDCluster.h:156
Belle2::PXDCluster::getSensorID
VxdID getSensorID() const
Get the sensor ID.
Definition: PXDCluster.h:126
Belle2::PXDCluster::getU
float getU() const
Get u coordinate of hit position.
Definition: PXDCluster.h:131
Belle2::PXDDigit
The PXD digit class.
Definition: PXDDigit.h:27
Belle2::PXDDigit::getVCellID
unsigned short getVCellID() const
Get cell ID in v.
Definition: PXDDigit.h:74
Belle2::PXDDigit::getUCellID
unsigned short getUCellID() const
Get cell ID in u.
Definition: PXDDigit.h:69
Belle2::PXDDigit::getCharge
unsigned short getCharge() const
Get collected charge.
Definition: PXDDigit.h:79
Belle2::PXDDigit::getSensorID
VxdID getSensorID() const
Get the sensor ID.
Definition: PXDDigit.h:64
Belle2::PXD::SensorInfo
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition: SensorInfo.h:23
Belle2::RelationArray
Low-level class to create/modify relations between StoreArrays.
Definition: RelationArray.h:62
Belle2::SVDCluster
The SVD Cluster class This class stores all information about reconstructed SVD clusters.
Definition: SVDCluster.h:28
Belle2::SVDCluster::getClsTime
float getClsTime() const
Get average of waveform maximum times of cluster strip signals.
Definition: SVDCluster.h:133
Belle2::SVDCluster::getCharge
float getCharge() const
Get collected charge.
Definition: SVDCluster.h:143
Belle2::SVDCluster::getSensorID
VxdID getSensorID() const
Get the sensor ID.
Definition: SVDCluster.h:101
Belle2::SVDCluster::isUCluster
bool isUCluster() const
Get the direction of strips.
Definition: SVDCluster.h:109
Belle2::SVDCluster::getPosition
float getPosition(double v=0) const
Get the coordinate of reconstructed hit.
Definition: SVDCluster.h:116
Belle2::SVDShaperDigit
The SVD ShaperDigit class.
Definition: SVDShaperDigit.h:32
Belle2::SVDShaperDigit::APVFloatSamples
std::array< APVFloatSampleType, c_nAPVSamples > APVFloatSamples
array of APVFloatSampleType objects
Definition: SVDShaperDigit.h:49
Belle2::SVDShaperDigit::getSensorID
VxdID getSensorID() const
Get the sensor ID.
Definition: SVDShaperDigit.h:98
Belle2::SVDShaperDigit::getSamples
APVFloatSamples getSamples() const
Get array of samples.
Definition: SVDShaperDigit.h:119
Belle2::SVDShaperDigit::getCellID
short int getCellID() const
Get strip ID.
Definition: SVDShaperDigit.h:114
Belle2::SVDShaperDigit::isUStrip
bool isUStrip() const
Get strip direction.
Definition: SVDShaperDigit.h:109
Belle2::SVDShaperDigit::getFADCTime
float getFADCTime() const
Get digit FADCTime estimate.
Definition: SVDShaperDigit.h:155
Belle2::SVD::SensorInfo
Specific implementation of SensorInfo for SVD Sensors which provides additional sensor specific infor...
Definition: SensorInfo.h:25
Belle2::StoreAccessorBase::isRequired
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.
Definition: StoreAccessorBase.h:78
Belle2::StoreAccessorBase::getName
const std::string & getName() const
Return name under which the object is saved in the DataStore.
Definition: StoreAccessorBase.h:128
Belle2::StoreArray
Accessor to arrays stored in the data store.
Definition: StoreArray.h:113
Belle2::StoreArray::getEntries
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:216
Belle2::Unit
The Unit class.
Definition: Unit.h:40
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:39
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:175
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:33
Belle2::VxdID::getLayerNumber
baseType getLayerNumber() const
Get the layer id.
Definition: VxdID.h:96
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:650
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:17