Belle II Software development
KLMStripEfficiencyAlgorithm.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/* Own header. */
10#include <klm/calibration/KLMStripEfficiencyAlgorithm.h>
11
12/* KLM headers. */
13#include <klm/dataobjects/KLMChannelIndex.h>
14
15/* ROOT headers. */
16#include <TFile.h>
17#include <TH1F.h>
18
19using namespace Belle2;
20
24
26{
27 m_AchievedPrecision = results.m_AchievedPrecision;
28 m_MatchedDigits = results.m_MatchedDigits;
29 m_ExtHits = results.m_ExtHits;
31 m_Efficiency = new float[nPlanes];
32 m_ExtHitsPlane = new int[nPlanes];
33 std::memcpy(m_Efficiency, results.m_Efficiency, nPlanes * sizeof(float));
34 std::memcpy(m_ExtHitsPlane, results.m_ExtHitsPlane, nPlanes * sizeof(int));
35}
36
38{
39 if (m_Efficiency != nullptr)
40 delete m_Efficiency;
41 if (m_ExtHitsPlane != nullptr)
42 delete m_ExtHitsPlane;
43}
44
46 CalibrationAlgorithm("KLMStripEfficiencyCollector"),
47 m_ElementNumbers(&(KLMElementNumbers::Instance())),
49 m_StripEfficiency(nullptr)
50{
51 int nPlanes = m_PlaneArrayIndex->getNElements();
52 m_Results.m_Efficiency = new float[nPlanes];
53 m_Results.m_ExtHitsPlane = new int[nPlanes];
54}
55
61
63{
64 /* Get collected results. */
65 int nPlanes = m_PlaneArrayIndex->getNElements();
66 TH1F* efficiencyHistogram =
67 new TH1F("plane_efficiency", "KLM plane efficiency",
68 nPlanes, -0.5, double(nPlanes) - 0.5);
69 std::shared_ptr<TH1F> matchedDigitsInPlane;
70 matchedDigitsInPlane = getObjectPtr<TH1F>("matchedDigitsInPlane");
71 m_Results.m_MatchedDigits = matchedDigitsInPlane->Integral();
72 std::shared_ptr<TH1F> allExtHitsInPlane;
73 allExtHitsInPlane = getObjectPtr<TH1F>("allExtHitsInPlane");
74 m_Results.m_ExtHits = allExtHitsInPlane->Integral();
75 matchedDigitsInPlane.get()->Sumw2();
76 allExtHitsInPlane.get()->Sumw2();
77 efficiencyHistogram->Divide(matchedDigitsInPlane.get(),
78 allExtHitsInPlane.get(), 1, 1, "B");
79 for (int i = 0; i < nPlanes; ++i) {
80 m_Results.m_Efficiency[i] = efficiencyHistogram->GetBinContent(i + 1);
81 m_Results.m_ExtHitsPlane[i] = allExtHitsInPlane->GetBinContent(i + 1);
82 }
83 /* Check whether the amount of data is sufficient. */
84 bool notEnoughData = false;
85 m_Results.m_AchievedPrecision = 0;
87 for (KLMChannelIndex& klmPlane : klmPlanes) {
88 KLMPlaneNumber plane = klmPlane.getKLMPlaneNumber();
89 KLMPlaneNumber planeIndex = m_PlaneArrayIndex->getIndex(plane);
90 int extHits = allExtHitsInPlane->GetBinContent(planeIndex + 1);
91 float efficiencyError = efficiencyHistogram->GetBinError(planeIndex + 1);
92 if (efficiencyError > m_Results.m_AchievedPrecision)
93 m_Results.m_AchievedPrecision = efficiencyError;
94 /*
95 * No hits is not considered as "not enough data", because this can
96 * happen in case KLM is excluded.
97 */
98 switch (m_CalibrationStage) {
100 if (extHits != 0 && extHits < m_MinimalExtHits)
101 notEnoughData = true;
102 break;
104 if (efficiencyError > m_RequestedPrecision)
105 notEnoughData = true;
106 break;
107 }
108 }
109 /*
110 * Fill the payload. A new object is created, because saveCalibration()
111 * stores a pointer to KLMStripEfficiency, and it is necessary to save
112 * the payloads to commit them at the end of calibration.
113 */
115 (!notEnoughData || m_ForcedCalibration)) {
117 KLMChannelIndex klmChannels;
118 for (KLMChannelIndex& klmChannel : klmChannels) {
119 int subdetector = klmChannel.getSubdetector();
120 int section = klmChannel.getSection();
121 int sector = klmChannel.getSector();
122 int layer = klmChannel.getLayer();
123 int plane = klmChannel.getPlane();
124 int strip = klmChannel.getStrip();
125 KLMPlaneNumber planeKLM = 0;
126 if (subdetector == KLMElementNumbers::c_BKLM) {
127 planeKLM = m_ElementNumbers->planeNumberBKLM(
128 section, sector, layer, plane);
129 } else {
130 planeKLM = m_ElementNumbers->planeNumberEKLM(
131 section, sector, layer, plane);
132 }
133 KLMPlaneNumber planeIndex = m_PlaneArrayIndex->getIndex(planeKLM);
134 float efficiency = efficiencyHistogram->GetBinContent(planeIndex + 1);
135 float efficiencyError = efficiencyHistogram->GetBinError(planeIndex + 1);
136 /* Fill the efficiency for this strip. */
137 if (subdetector == KLMElementNumbers::c_BKLM) {
138 m_StripEfficiency->setBarrelEfficiency(
139 section, sector, layer, plane, strip, efficiency, efficiencyError);
140 } else {
141 m_StripEfficiency->setEndcapEfficiency(
142 section, sector, layer, plane, strip, efficiency, efficiencyError);
143 }
144 }
145 saveCalibration(m_StripEfficiency, "KLMStripEfficiency");
146 }
147 /* Write histograms to output file. */
148 TFile* outputFile = new TFile(m_OutputFileName.c_str(), "recreate");
149 outputFile->cd();
150 matchedDigitsInPlane.get()->Write();
151 allExtHitsInPlane.get()->Write();
152 efficiencyHistogram->Write();
153 delete efficiencyHistogram;
154 delete outputFile;
155 /* Set output status. */
156 if (notEnoughData && !m_ForcedCalibration)
159}
160
162 const float* efficiency) const
163{
164 const int nPlanes = KLMPlaneArrayIndex::Instance().getNElements();
165 int newPlanes = 0;
166 for (int i = 0; i < nPlanes; ++i) {
167 if (m_Efficiency[i] > 0 && efficiency[i] == 0)
168 newPlanes++;
169 }
170 return newPlanes;
171}
172
174 const int* extHitsPlane) const
175{
176 const int nPlanes = KLMPlaneArrayIndex::Instance().getNElements();
177 int newPlanes = 0;
178 for (int i = 0; i < nPlanes; ++i) {
179 if (m_ExtHitsPlane[i] > 0 && extHitsPlane[i] == 0)
180 newPlanes++;
181 }
182 return newPlanes;
183}
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
EResult
The result of calibration.
@ c_OK
Finished successfully =0 in Python.
@ c_NotEnoughData
Needs more data =2 in Python.
CalibrationAlgorithm(const std::string &collectorModuleName)
Constructor - sets the prefix for collected objects (won't be accesses until execute(....
KLM channel index.
uint16_t getNElements() const
Get number of elements.
KLM plane array index.
static const KLMPlaneArrayIndex & Instance()
Instantiation.
int newMeasuredPlanes(const float *efficiency) const
Get number of new measured planes.
float m_AchievedPrecision
Achieved precision of efficiency measurement.
int newExtHitsPlanes(const int *extHitsPlane) const
Get number of new measured planes with ExtHits.
@ c_MeasurablePlaneCheck
Check of set of planes with determined efficiency.
enum CalibrationStage m_CalibrationStage
Calibration stage.
bool m_ForcedCalibration
Whether the calibration is forced.
const KLMElementNumbers * m_ElementNumbers
Element numbers.
const KLMPlaneArrayIndex * m_PlaneArrayIndex
Plane array index.
KLMStripEfficiency * m_StripEfficiency
Efficiency data object.
CalibrationAlgorithm::EResult calibrate() override
Calibration.
int m_MinimalExtHits
Minimal number of ExtHits per plane.
float m_RequestedPrecision
Requested precision of efficiency measurement.
DBObject used to store the efficiencies of KLM strips.
std::shared_ptr< T > getObjectPtr(const std::string &name, const std::vector< Calibration::ExpRun > &requestedRuns)
Get calibration data object by name and list of runs, the Merge function will be called to generate t...
uint16_t KLMPlaneNumber
Plane number.
Abstract base class for different kinds of events.