Belle II Software development
PhysicsObjectsMiraBelleHadronModule.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 <dqm/modules/PhysicsObjectsMiraBelle/PhysicsObjectsMiraBelleHadronModule.h>
10#include <hlt/softwaretrigger/calculations/utilities.h>
11#include <analysis/dataobjects/ParticleList.h>
12#include <analysis/ContinuumSuppression/FoxWolfram.h>
13#include <analysis/utility/PCmsLabTransform.h>
14#include <framework/datastore/StoreObjPtr.h>
15#include <framework/datastore/StoreArray.h>
16#include <mdst/dataobjects/Track.h>
17#include <mdst/dataobjects/SoftwareTriggerResult.h>
18#include <analysis/ClusterUtility/ClusterUtils.h>
19#include <TDirectory.h>
20#include <map>
21
22using namespace Belle2;
23
24REG_MODULE(PhysicsObjectsMiraBelleHadron);
25
27{
28 setDescription("Monitor Physics Objects Quality");
30
31 addParam("TriggerIdentifier", m_triggerIdentifier,
32 "Trigger identifier string used to select events for the histograms", std::string("software_trigger_cut&skim&accept_hadronb2"));
33 addParam("hadronb2piPListName", m_hadpiPListName, "Name of the pi+ particle list", std::string("pi+:hadb2physMiraBelle"));
34}
35
37{
38 TDirectory* oldDir = gDirectory;
39 oldDir->mkdir("PhysicsObjectsMiraBelleHadron");
40 oldDir->cd("PhysicsObjectsMiraBelleHadron");
41
42 // Mass distributions
43 m_h_nECLClusters = new TH1F("hist_nECLClusters", "hist_nECLClusters", 100, 0, 60);
44 m_h_nECLClusters->SetXTitle("hist_nECLClusters");
45 m_h_visibleEnergyCMSnorm = new TH1F("hist_visibleEnergyCMSnorm", "hist_visibleEnergyCMSnorm", 100, 0, 2);
46 m_h_visibleEnergyCMSnorm->SetXTitle("hist_visibleEnergyCMSnorm");
47 m_h_EsumCMSnorm = new TH1F("hist_EsumCMSnorm", "hist_EsumCMSnorm", 100, 0, 2);
48 m_h_EsumCMSnorm->SetXTitle("hist_EsumCMSnorm");
49 m_h_R2 = new TH1F("hist_R2", "hist_R2", 100, 0, 1);
50 m_h_R2->SetXTitle("hist_R2");
51 m_h_physicsresultsH = new TH1F("hist_physicsresultsH", "hist_physicsresultsH", 10, 0, 10);
52 m_h_physicsresultsH->SetXTitle("hist_physicsresultsH");
53 m_h_physicsresultsH->GetXaxis()->SetBinLabel(2, "Hadronb2");
54 m_h_physicsresultsH->GetXaxis()->SetBinLabel(3, "Hadronb2_tight");
55
56 oldDir->cd();
57}
58
59
61{
62 REG_HISTOGRAM
63
65 result.isOptional();
66}
67
69{
70
71 m_h_nECLClusters->Reset();
73 m_h_EsumCMSnorm->Reset();
74 m_h_R2->Reset();
75 m_h_physicsresultsH->Reset();
76}
77
79{
80
82 if (!result.isValid()) {
83 B2WARNING("SoftwareTriggerResult object not available but needed to select events for the histograms.");
84 return;
85 }
86
87 const std::map<std::string, int>& results = result->getResults();
88 if (results.find(m_triggerIdentifier) == results.end()) {
89 B2WARNING("PhysicsObjectsMiraBelleHadron: Can't find trigger identifier: " << m_triggerIdentifier);
90 return;
91 }
92
93 // apply software trigger
94 const bool accepted = (result->getResult(m_triggerIdentifier) == SoftwareTriggerCutResult::c_accept);
95 if (accepted == false) return;
96 m_h_physicsresultsH->Fill(1);
97 // get pi list
99 std::vector<ROOT::Math::PxPyPzEVector> m_pionHad;
100
101 double EsumPiHad = 0.;
102
103 for (unsigned int i = 0; i < hadpiParticles->getListSize(); i++) {
104 const Particle* parPiHad = hadpiParticles->getParticle(i);
105 ROOT::Math::PxPyPzEVector V4PiHad = PCmsLabTransform::labToCms(parPiHad->get4Vector());
106 m_pionHad.push_back(V4PiHad);
107 EsumPiHad += V4PiHad.E();
108
109 }
110 //nECLClustersLE
111 double neclClusters = -1.;
112 double eneclClusters = 0.;
113 StoreArray<ECLCluster> eclClusters;
114 ClusterUtils Cl;
115 double EsumGamma = 0.;
116
117 if (eclClusters.isValid()) {
118 const unsigned int numberOfECLClusters = std::count_if(eclClusters.begin(), eclClusters.end(),
119 [](const ECLCluster & eclcluster) {
120 return (eclcluster.hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)
121 and eclcluster.getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) > 0.1);
122 });
123 neclClusters = numberOfECLClusters;
124 for (int ncl = 0; ncl < eclClusters.getEntries(); ncl++) {
125 if (eclClusters[ncl]->hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons)
126 && eclClusters[ncl]->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) > 0.1) {
127 eneclClusters += eclClusters[ncl]->getEnergy(ECLCluster::EHypothesisBit::c_nPhotons);
128 if (!eclClusters[ncl]->getRelatedFrom<Track>()) {
129 ROOT::Math::PxPyPzEVector V4Gamma_CMS = PCmsLabTransform::labToCms(Cl.Get4MomentumFromCluster(eclClusters[ncl],
131 EsumGamma += V4Gamma_CMS.E();
132 }
133 }
134 }
135 }
136
137 double visibleEnergyCMSnorm = (EsumPiHad + EsumGamma) / (Belle2::SoftwareTrigger::BeamEnergyCMS() * 2.0);
138 double EsumCMSnorm = eneclClusters / (Belle2::SoftwareTrigger::BeamEnergyCMS() * 2.0);
139 FoxWolfram fw(m_pionHad);
141 double R2 = fw.getR(2);
142
143 m_h_nECLClusters->Fill(neclClusters);
144 m_h_visibleEnergyCMSnorm->Fill(visibleEnergyCMSnorm);
145 m_h_EsumCMSnorm->Fill(EsumCMSnorm);
146 m_h_R2->Fill(R2);
147 bool hadronb_tag = visibleEnergyCMSnorm > 0.4 && EsumCMSnorm > 0.2 && R2 < 0.2;
148 if (hadronb_tag) {
149 m_h_physicsresultsH->Fill(2);
150 }
151
152}
153
154
155
159
163
Class to provide momentum-related information from ECLClusters.
const ROOT::Math::PxPyPzEVector Get4MomentumFromCluster(const ECLCluster *cluster, ECLCluster::EHypothesisBit hypo)
Returns four momentum vector.
ECL cluster data.
Definition ECLCluster.h:27
@ c_nPhotons
CR is split into n photons (N1)
Definition ECLCluster.h:41
Class to calculate the Fox-Wolfram moments up to order 8.
Definition FoxWolfram.h:28
double getR(int i) const
Returns the i-th moment normalized to the 0th-order moment.
Definition FoxWolfram.h:89
void calculateBasicMoments()
Method to perform the calculation of the moments up to order 4, which are the most relevant ones.
Definition FoxWolfram.cc:14
HistoModule()
Constructor.
Definition HistoModule.h:32
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
@ 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:80
static ROOT::Math::PxPyPzMVector labToCms(const ROOT::Math::PxPyPzMVector &vec)
Transforms Lorentz vector into CM System.
Class to store reconstructed particles.
Definition Particle.h:76
ROOT::Math::PxPyPzEVector get4Vector() const
Returns Lorentz vector.
Definition Particle.h:567
void event() override
This method is called for each event.
TH1F * m_h_physicsresultsH
histogram for event results for hadronb2 and hadronb2 tight
void endRun() override
This method is called if the current run ends.
void terminate() override
This method is called at the end of the event processing.
std::string m_triggerIdentifier
Trigger identifier string used to select events for the histograms.
TH1F * m_h_R2
histogram for R2 of hadron events after the hadronb2 selection
std::string m_hadpiPListName
Name of the pi+ particle list.
void beginRun() override
Called when entering a new run.
TH1F * m_h_visibleEnergyCMSnorm
histogram for visible energy of tracks and gammas
TH1F * m_h_EsumCMSnorm
histogram for sum of energy of clusters with E > 0.1
void defineHisto() override
Definition of the histograms.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
bool isValid() const
Check whether the array was registered.
Definition StoreArray.h:288
int getEntries() const
Get the number of objects in the array.
Definition StoreArray.h:216
iterator end()
Return iterator to last entry +1.
Definition StoreArray.h:320
iterator begin()
Return iterator to first entry.
Definition StoreArray.h:318
Type-safe access to single objects in the data store.
Definition StoreObjPtr.h:96
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:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Abstract base class for different kinds of events.