Belle II Software  release-05-01-25
B2BIIConvertBeamParamsModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2015-2018 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Ritter *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <b2bii/modules/B2BIIMdstInput/B2BIIConvertBeamParamsModule.h>
12 #include <belle_legacy/benergy/BeamEnergy.h>
13 #include <belle_legacy/ip/IpProfile.h>
14 #include <belle_legacy/panther/panther.h>
15 #include <belle_legacy/tables/belletdf.h>
16 
17 #include <framework/database/Database.h>
18 #include <framework/database/EventDependency.h>
19 #include <framework/database/IntervalOfValidity.h>
20 #include <framework/dbobjects/BeamParameters.h>
21 #include <framework/utilities/FileSystem.h>
22 #include <mdst/dbobjects/BeamSpot.h>
23 #include <mdst/dbobjects/CollisionBoostVector.h>
24 #include <mdst/dbobjects/CollisionInvariantMass.h>
25 
26 #include <list>
27 
28 namespace Belle {
32  int check_beginrun(bool) { return 0; }
33 }
34 
35 namespace {
37  TVector3 CLHEPtoROOT(const HepPoint3D& point)
38  {
39  return TVector3(point.x(), point.y(), point.z());
40  }
42  TMatrixDSym CLHEPtoROOT(const HepSymMatrix& matrix)
43  {
44  TMatrixDSym result(matrix.num_col());
45  for (int i = 0; i < matrix.num_row(); ++i) {
46  for (int j = 0; j < matrix.num_col(); ++j) {
47  result(i, j) = matrix(i + 1, j + 1);
48  }
49  }
50  return result;
51  }
52 }
53 
54 namespace Belle2 {
60  REG_MODULE(B2BIIConvertBeamParams);
61 
63  {
64  addParam("mcFlag", m_mcFlag, "which mc flag to use", m_mcFlag);
65  addParam("missingBenergy", m_missingBenergy, "where to store information about runs with missing database info", m_missingBenergy);
66  addParam("missingIp", m_missingIp, "where to store information about runs with missing IP profile info", m_missingIp);
67  }
68 
70  {
71  BsInit(0);
72  m_event.isRequired();
73  }
74 
76  {
77  // Fudge a runhead record with all the info we need to be able to get the
78  // Beamenergy/IpProfile from the database
79  BsClrTab(BBS_CLEAR);
80  Belle::Belle_runhead_Manager& runMgr = Belle::Belle_runhead_Manager::get_manager();
81  runMgr.add();
82  Belle::Belle_runhead& runhead = runMgr[0];
83  runhead.ExpMC(m_mcFlag);
84  runhead.ExpNo(m_event->getExperiment());
85  runhead.RunNo(m_event->getRun());
86  B2INFO("Obtaining values for exp " << m_event->getExperiment() << ", run " << m_event->getRun());
87  // and then get the values from the database
88  Belle::BeamEnergy::flag_database(0);
89  Belle::IpProfile::begin_run();
90  Belle::Ip_profile_general_Manager& ipMgr = Belle::Ip_profile_general_Manager::get_manager();
91 
92  if (!Belle::BeamEnergy::usable()) {
94  if (!lock.lock()) {
95  B2ERROR("No BeamEnergy for exp " << m_event->getExperiment() << ", run " << m_event->getRun());
96  } else {
97  std::ofstream file(m_missingBenergy.c_str(), std::ios::app);
98  file << m_event->getExperiment() << "," << m_event->getRun() << std::endl;
99  }
100  return;
101  }
102 
103  const double Eher = Belle::BeamEnergy::E_HER();
104  const double Eler = Belle::BeamEnergy::E_LER();
105  const double crossingAngle = Belle::BeamEnergy::Cross_angle();
106  const double angleLer = M_PI; //parallel to negative z axis (different from Belle II!)
107  const double angleHer = crossingAngle; //in positive z and x direction, verified to be consistent with Upsilon(4S) momentum
108 
109  std::vector<double> covariance; //0 entries = no error
110 
111  IntervalOfValidity iov(m_event->getExperiment(), m_event->getRun(), m_event->getExperiment(), m_event->getRun());
112 
113  BeamParameters beamParams;
114  beamParams.setLER(Eler, angleLer, covariance);
115  beamParams.setHER(Eher, angleHer, covariance);
116 
117  CollisionBoostVector collisionBoostVector;
118  CollisionInvariantMass collisionInvM;
119  TLorentzVector cms = beamParams.getLER() + beamParams.getHER();
120  collisionBoostVector.setBoost(cms.BoostVector(), TMatrixTSym<double>(3));
121  //note: maybe we could use Belle::BeamEnergy::E_beam_corr(), Belle::BeamEnergy::E_beam_err()
122  collisionInvM.setMass(cms.M(), 0.0 , 0.0);
123 
124  // Boost vector and invariant mass are not intra-run dependent, store now
125  Database::Instance().storeData("CollisionBoostVector", &collisionBoostVector, iov);
126  Database::Instance().storeData("CollisionInvariantMass", &collisionInvM, iov);
127 
128  // and now we continue with the vertex
129  if (!Belle::IpProfile::usable()) {
131  if (!lock.lock()) {
132  B2ERROR("No IpProfile for exp " << m_event->getExperiment() << ", run " << m_event->getRun());
133  } else {
134  std::ofstream file(m_missingIp.c_str(), std::ios::app);
135  file << m_event->getExperiment() << "," << m_event->getRun() << std::endl;
136  }
137  beamParams.setVertex(TVector3(
138  std::numeric_limits<double>::quiet_NaN(),
139  std::numeric_limits<double>::quiet_NaN(),
140  std::numeric_limits<double>::quiet_NaN()), covariance);
141  Database::Instance().storeData("BeamParameters", &beamParams, iov);
142 
143  BeamSpot beamSpot;
144  beamSpot.setIP(
145  TVector3(std::numeric_limits<double>::quiet_NaN(),
146  std::numeric_limits<double>::quiet_NaN(),
147  std::numeric_limits<double>::quiet_NaN()
148  ), TMatrixTSym<double>(3)
149  );
150  Database::Instance().storeData("BeamSpot", &beamSpot, iov);
151  B2INFO("No IpProfile, created BeamEnergy Payload");
152  return;
153  }
154 
155  int nbins = std::max(ipMgr.count() - 1, 1);
156  B2INFO("exp " << m_event->getExperiment() << ", run " << m_event->getRun() << ": " << nbins << " IpProfile bins");
157 
158  // need to keep the objects alive until we created the payloads
159  std::list<BeamParameters> beamparamsList;
160  std::list<BeamSpot> beamspotList;
161  // and we want them to be intra-run dependent
162  std::unique_ptr<EventDependency> beamparams;
163  std::unique_ptr<EventDependency> beamspots;
164 
165  Belle::Belle_event_Manager& evtMgr = Belle::Belle_event_Manager::get_manager();
166  if (!evtMgr.count()) evtMgr.add();
167  for (int i = 0; i < nbins; ++i) {
168  int evtNr = i * BIN_EVENTS;
169  // fudge belle event record to get the correct run dependent ip profile
170  Belle::Belle_event& evt = evtMgr[0];
171  evt.ExpMC(m_mcFlag);
172  evt.ExpNo(m_event->getExperiment());
173  evt.RunNo(m_event->getRun());
174  evt.EvtNo(evtNr);
175  Belle::IpProfile::set_evtbin_number();
176  B2ASSERT("something wrong: " << Belle::IpProfile::EvtBinNo() << "!=" << i, Belle::IpProfile::EvtBinNo() == i);
177 
178  // and convert
179  HepPoint3D ip = Belle::IpProfile::e_position();
180  HepSymMatrix ipErr = Belle::IpProfile::e_position_err();
181 
182  beamParams.setVertex(CLHEPtoROOT(ip));
183  beamParams.setCovVertex(CLHEPtoROOT(ipErr));
184  beamparamsList.emplace_back(beamParams);
185 
186  BeamSpot beamSpot;
187  beamSpot.setIP(CLHEPtoROOT(ip), CLHEPtoROOT(ipErr));
188  beamspotList.emplace_back(beamSpot);
189 
190  if (!beamparams) {
191  beamparams = std::make_unique<EventDependency>(&beamparamsList.back(), false);
192  beamspots = std::make_unique<EventDependency>(&beamspotList.back(), false);
193  } else {
194  beamparams->add(evtNr, &beamparamsList.back());
195  beamspots->add(evtNr, &beamspotList.back());
196  }
197  }
198  if (beamparamsList.size() < 1) {
199  B2ERROR("Something is wrong with exp " << m_event->getExperiment() << ", run " << m_event->getRun() << ": no bins found");
200  return;
201  }
202  if (beamparamsList.size() == 1) {
203  // just one bin? no need to store event dependency
204  Database::Instance().storeData("BeamParameters", &beamparamsList.back(), iov);
205  Database::Instance().storeData("BeamSpot", &beamspotList.back(), iov);
206  B2INFO("Created event independent payload");
207  } else {
208  // otherwise store full information
209  Database::Instance().storeData("BeamParameters", beamparams.get(), iov);
210  Database::Instance().storeData("BeamSpot", beamspots.get(), iov);
211  B2INFO("Created event dependent payload with " << beamparamsList.size() << " entries");
212  }
213  }
215 } //Belle2 namespace
Belle2::IntervalOfValidity
A class that describes the interval of experiments/runs for which an object in the database is valid.
Definition: IntervalOfValidity.h:35
Belle2::MCInitialParticles::getLER
const TLorentzVector & getLER() const
Get 4vector of the low energy beam.
Definition: MCInitialParticles.h:140
Belle2::B2BIIConvertBeamParamsModule::m_missingBenergy
std::string m_missingBenergy
Where to store information about runs without beam energy information.
Definition: B2BIIConvertBeamParamsModule.h:55
Belle2::CollisionBoostVector::setBoost
void setBoost(const TVector3 &boost, const TMatrixDSym &covariance)
Set the boost vector and its error matrix.
Definition: CollisionBoostVector.h:46
Belle2::FileSystem::Lock::lock
bool lock(int timeout=300, bool ignoreErrors=false)
Try to lock the file.
Definition: FileSystem.cc:190
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::B2BIIConvertBeamParamsModule::initialize
void initialize() override
Initialize phanter banks.
Definition: B2BIIConvertBeamParamsModule.cc:69
Belle2::Database::storeData
bool storeData(const std::string &name, TObject *object, const IntervalOfValidity &iov)
Store an object in the database.
Definition: Database.cc:152
Belle2::B2BIIConvertBeamParamsModule::m_event
StoreObjPtr< EventMetaData > m_event
Event metadata.
Definition: B2BIIConvertBeamParamsModule.h:59
Belle2::MCInitialParticles::getHER
const TLorentzVector & getHER() const
Get 4vector of the high energy beam.
Definition: MCInitialParticles.h:137
Belle2::B2BIIConvertBeamParamsModule::B2BIIConvertBeamParamsModule
B2BIIConvertBeamParamsModule()
Create parameters.
Definition: B2BIIConvertBeamParamsModule.cc:62
Belle2::BeamSpot
This class contains the beam spot position and size modeled as a gaussian distribution in space.
Definition: BeamSpot.h:32
Belle2::BeamParameters::setCovVertex
void setCovVertex(const TMatrixDSym &cov)
Set the covariance matrix of the vertex position.
Definition: BeamParameters.h:80
Belle2::BeamParameters::setVertex
void setVertex(const TVector3 &vertex, const std::vector< double > &cov)
Set the vertex position and error matrix.
Belle2::B2BIIConvertBeamParamsModule::m_mcFlag
int m_mcFlag
mcFlag to use when getting belle database content
Definition: B2BIIConvertBeamParamsModule.h:53
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::BeamSpot::setIP
void setIP(const TVector3 &ipPosition, const TMatrixDSym &covariance)
Set the IP position and its error matrix.
Definition: BeamSpot.h:69
Belle2::BeamParameters::setHER
void setHER(double energy, double angle, const std::vector< double > &cov)
Set the HER FourVector and error matrix from beam energy, angle and covariance entries.
Belle2::CollisionInvariantMass
This class contains the measured average center-of-mass energy, which is equal to the invariant mass ...
Definition: CollisionInvariantMass.h:31
Belle2::B2BIIConvertBeamParamsModule::m_missingIp
std::string m_missingIp
Where to store information about runs without IP profile information.
Definition: B2BIIConvertBeamParamsModule.h:57
Belle2::FileSystem::Lock
Helper class for locking a file.
Definition: FileSystem.h:107
Belle2::Module::addParam
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:562
Belle2::Database::Instance
static Database & Instance()
Instance of a singleton Database.
Definition: Database.cc:54
Belle2::BeamParameters::setLER
void setLER(double energy, double angle, const std::vector< double > &cov)
Set the LER FourVector and error matrix from beam energy, angle and covariance entries.
HepGeom::Point3D< double >
Belle2::CollisionBoostVector
This class contains the measured average boost vector vec(beta) = (beta_x, beta_y,...
Definition: CollisionBoostVector.h:33
Belle2::BeamParameters
This class contains the nominal beam parameters and the parameters used for smearing of the primary v...
Definition: BeamParameters.h:33
Belle2::CollisionInvariantMass::setMass
void setMass(double mass, double error, double spread)
Set the CMS energy and its uncertainty.
Definition: CollisionInvariantMass.h:45
Belle2::B2BIIConvertBeamParamsModule::beginRun
void beginRun() override
Set run info in panther and load IPProfile/Benergy and convert the values to payloads.
Definition: B2BIIConvertBeamParamsModule.cc:75