Belle II Software  release-05-01-25
SVDCrossTalkCalibrationsAlgorithm.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/calibration/SVDCrossTalkCalibrationsAlgorithm.h>
12 #include <svd/calibration/SVDCrossTalkStripsCalibrations.h>
13 
14 #include <TH1F.h>
15 
16 #include <iostream>
17 
18 using namespace std;
19 using namespace Belle2;
20 
21 SVDCrossTalkCalibrationsAlgorithm::SVDCrossTalkCalibrationsAlgorithm(const std::string& str) :
22  CalibrationAlgorithm("SVDCrossTalkCalibrationsCollector")
23 {
24  setDescription("Calibration algorithm for SVDCrossTalkCalibrations payload");
25  m_id = str;
26 }
27 
29 {
30  int isCrossTalkCal = 0;
31  auto payload = new Belle2::SVDCrossTalkStripsCalibrations::t_payload(isCrossTalkCal, m_id);
32 
33  auto tree = getObjectPtr<TTree>("HTreeCrossTalkCalib");
34 
35  TH1F* hist = new TH1F("", "", 768, 0, 768);
36 
37  if (!tree) {
38  B2WARNING("No tree");
39  } else if (!tree->GetEntries()) {
40  B2WARNING("Empty tree");
41  }
42 
43  int layer = 0;
44  int ladder = 0;
45  int sensor = 0;
46  int side = 0;
47 
48  tree->SetBranchAddress("hist", &hist);
49  tree->SetBranchAddress("layer", &layer);
50  tree->SetBranchAddress("ladder", &ladder);
51  tree->SetBranchAddress("sensor", &sensor);
52  tree->SetBranchAddress("side", &side);
53 
54 
55 
56  for (int i = 0; i < tree->GetEntries(); i++) {
57  tree->GetEntry(i);
58 
59  int nstrips = 768;
60  if (side == 0 && layer != 3) nstrips = 512;
61 
62  for (int strip = 0; strip < nstrips; strip++) {
63 
64  isCrossTalkCal = hist->GetBinContent(strip + 1);
65 
66  if (layer == 4 && ladder == 1 && sensor == 2 && side == 1
67  && hist->GetEntries() < m_minEntries) { //Check that we have enough events populating the calibration
68  cout << "Not enough Data: " << hist->GetEntries() << " entries found" << endl;
69  return c_NotEnoughData;
70  }
71 
72  payload->set(layer, ladder, sensor, side, strip, isCrossTalkCal);
73  }
74  }
75 
76 
77  saveCalibration(payload, "SVDCrossTalkStripsCalibrations");
78 
79  return c_OK;
80 
81 }
Belle2::SVDCrossTalkCalibrationsAlgorithm::m_id
std::string m_id
Identifier string.
Definition: SVDCrossTalkCalibrationsAlgorithm.h:49
Belle2::CalibrationAlgorithm::saveCalibration
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
Definition: CalibrationAlgorithm.cc:290
Belle2::CalibrationAlgorithm::c_OK
@ c_OK
Finished successfuly =0 in Python.
Definition: CalibrationAlgorithm.h:51
Belle2::CalibrationAlgorithm::setDescription
void setDescription(const std::string &description)
Set algorithm description (in constructor)
Definition: CalibrationAlgorithm.h:331
Belle2::SVDCrossTalkCalibrationsAlgorithm::calibrate
virtual EResult calibrate() override
Run algo on data - pure virtual: needs to be implemented.
Definition: SVDCrossTalkCalibrationsAlgorithm.cc:28
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::SVDCrossTalkCalibrationsAlgorithm::m_minEntries
int m_minEntries
Minimum number of required entries for collector histogram L4.1.2 u-side.
Definition: SVDCrossTalkCalibrationsAlgorithm.h:52
Belle2::CalibrationAlgorithm::EResult
EResult
The result of calibration.
Definition: CalibrationAlgorithm.h:50
Belle2::SVDCrossTalkStripsCalibrations::t_payload
SVDCalibrationsBase< SVDCalibrationsBitmap > t_payload
typedef of the SVDCrossTalkStripsCalibrations payload for all SVD strips
Definition: SVDCrossTalkStripsCalibrations.h:48
Belle2::CalibrationAlgorithm::c_NotEnoughData
@ c_NotEnoughData
Needs more data =2 in Python.
Definition: CalibrationAlgorithm.h:53
Belle2::CalibrationAlgorithm
Base class for calibration algorithms.
Definition: CalibrationAlgorithm.h:47