Belle II Software prerelease-10-00-00a
PhysicsObjectsMiraBelleEcmsBBModule.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#include <dqm/modules/PhysicsObjectsMiraBelle/PhysicsObjectsMiraBelleEcmsBBModule.h>
9
10#include <analysis/variables/ContinuumSuppressionVariables.h>
11#include <analysis/utility/PCmsLabTransform.h>
12#include <framework/particledb/EvtGenDatabasePDG.h>
13#include <mdst/dataobjects/SoftwareTriggerResult.h>
14#include <framework/datastore/StoreObjPtr.h>
15#include <analysis/dataobjects/ParticleList.h>
16
17#include <TDirectory.h>
18#include <TH1D.h>
19
20#include <cmath>
21
22using namespace Belle2;
23
24
25//-----------------------------------------------------------------
26// Register the Module
27//-----------------------------------------------------------------
28REG_MODULE(PhysicsObjectsMiraBelleEcmsBB);
29
30//-----------------------------------------------------------------
31// Implementation
32//-----------------------------------------------------------------
33
35{
36 //Set module properties
37 setDescription("Monitor of the CMS collision energy based on hadronic B decays");
39
40 addParam("TriggerIdentifier", m_triggerIdentifier,
41 "Trigger identifier string used to select events for the histograms", std::string("software_trigger_cut&skim&accept_btocharm"));
42 addParam("BmListName", m_BmListName, "Name of the B- particle list", std::string("B-:combined"));
43 addParam("B0ListName", m_B0ListName, "Name of the B0 particle list", std::string("B0:combined"));
44}
45
46
51
53{
54 TDirectory* oldDir = gDirectory;
55 oldDir->mkdir("PhysicsObjectsMiraBelleEcmsBB");
56 oldDir->cd("PhysicsObjectsMiraBelleEcmsBB");
57
58 const double cMBp = EvtGenDatabasePDG::Instance()->GetParticle("B+")->Mass();
59
60 // In the original ML-based analysis there are 40 bins for visualisation
61 // Since binned fit is used here, 80 bins are used instead.
62 // For simplicity, the same histogram lower bound (B+ mass) is used
63 m_hB0 = new TH1D("hB0", "", 80, cMBp, 5.37);
64 m_hBp = new TH1D("hBp", "", 80, cMBp, 5.37);
65 oldDir->cd();
66}
67
69{
70 if (m_hBp) m_hBp->Reset();
71 if (m_hB0) m_hB0->Reset();
72}
73
74
75
76
77
79{
81 if (!result.isValid()) {
82 B2WARNING("SoftwareTriggerResult object not available but needed to select events for the histograms.");
83 return;
84 }
85
86 const bool accepted = (result->getResult(m_triggerIdentifier) == SoftwareTriggerCutResult::c_accept);
87 if (accepted == false) return;
88
91
92 //put all the B candidates into the vector
93 std::vector<const Particle*> Bparts;
94
95 if (B0.isValid()) {
96 for (unsigned i = 0; i < B0->getListSize(); ++i)
97 if (B0->getParticle(i))
98 Bparts.push_back(B0->getParticle(i));
99 }
100 if (Bm.isValid()) {
101 for (unsigned i = 0; i < Bm->getListSize(); ++i)
102 if (Bm->getParticle(i))
103 Bparts.push_back(Bm->getParticle(i));
104 }
105
106
107 if (Bparts.size() == 0) return;
108
109 for (unsigned i = 0; i < Bparts.size(); ++i) {
110 const Particle* Bpart = Bparts[i];
111
112 const Particle* D = nullptr;
113 double dmDstar = std::numeric_limits<double>::quiet_NaN();
114
115 static const int c_PdgD0 = std::abs(EvtGenDatabasePDG::Instance()->GetParticle("D0")->PdgCode());
116 static const int c_PdgDplus = std::abs(EvtGenDatabasePDG::Instance()->GetParticle("D+")->PdgCode());
117
118 static const int c_PdgDstar0 = std::abs(EvtGenDatabasePDG::Instance()->GetParticle("D*0")->PdgCode());
119 static const int c_PdgDstarPlus = std::abs(EvtGenDatabasePDG::Instance()->GetParticle("D*+")->PdgCode());
120
121 //if D0 or D+ meson
122 if (std::abs(Bpart->getDaughter(0)->getPDGCode()) == c_PdgD0 || std::abs(Bpart->getDaughter(0)->getPDGCode()) == c_PdgDplus) {
123 D = Bpart->getDaughter(0);
124 } else if (std::abs(Bpart->getDaughter(0)->getPDGCode()) == c_PdgDstar0
125 || std::abs(Bpart->getDaughter(0)->getPDGCode()) == c_PdgDstarPlus) {
126 const Particle* Dstar = Bpart->getDaughter(0);
127 D = Dstar->getDaughter(0);
128 dmDstar = Dstar->getMass() - D->getMass();
129 } else {
130 B2INFO("No D meson found");
131 }
132 double mD = D ? D->getMass() : std::numeric_limits<double>::quiet_NaN();
133
134 //Convert mBC and deltaE to the Y4S reference
135 double pBcms = PCmsLabTransform::labToCms(Bpart->get4Vector()).P();
136 double mInv = Bpart->getMass();
137 double pdg = Bpart->getPDGCode();
138 double R2 = Variable::R2(Bpart);
139
140 //get mass of B+- or B0
141 double mB = EvtGenDatabasePDG::Instance()->GetParticle(std::abs(pdg))->Mass();
142
143
144 // Filling the histograms
145 if (c_mDmin < mD && mD < c_mDmax)
146 if (std::abs(mInv - mB) < c_mBwindow)
147 if (R2 < c_R2max)
148 if (std::isnan(dmDstar) || (c_dmDstarMin < dmDstar && dmDstar < c_dmDstarMax)) {
149 double eBC = sqrt(pBcms * pBcms + mB * mB); // beam constrained energy
150 if (std::abs(pdg) == 511) {
151 m_hB0->Fill(eBC);
152 } else {
153 m_hBp->Fill(eBC);
154 }
155 }
156 }
157}
static EvtGenDatabasePDG * Instance()
Instance method that loads the EvtGen table.
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
int getPDGCode(void) const
Returns PDG code.
Definition Particle.h:465
ROOT::Math::PxPyPzEVector get4Vector() const
Returns Lorentz vector.
Definition Particle.h:567
const Particle * getDaughter(unsigned i) const
Returns a pointer to the i-th daughter particle.
Definition Particle.cc:662
double getMass() const
Returns invariant mass (= nominal for FS particles)
Definition Particle.h:527
void initialize() override final
Register the histograms.
std::string m_B0ListName
List name for neutral B candidates.
static constexpr double c_mDmax
Maximal value of the D meson inv mass.
void defineHisto() override final
Initialize the histograms.
std::string m_triggerIdentifier
Trigger identifier string used to select events for the histograms.
void event() override final
Event processor Filling of the histograms.
static constexpr double c_mDmin
Minimal value of the D meson inv mass.
static constexpr double c_dmDstarMin
Minimal value of the m(D*)-(mD)
std::string m_BmListName
List name for charged B candidates.
static constexpr double c_dmDstarMax
Maximal value of the m(D*)-(mD)
static constexpr double c_mBwindow
Maximal deviation of B meson inv mass from PDG value.
static constexpr double c_R2max
Maximal allowed R2 value (to suppress continuum)
Type-safe access to single objects in the data store.
Definition StoreObjPtr.h:96
bool isValid() const
Check whether the object was created.
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
double sqrt(double a)
sqrt for double
Definition beamHelpers.h:28
Abstract base class for different kinds of events.