9#include <analysis/modules/ContinuumSuppressionBuilder/ContinuumSuppressionBuilderModule.h>
11#include <analysis/ContinuumSuppression/Thrust.h>
12#include <analysis/ContinuumSuppression/KsfwMoments.h>
13#include <analysis/ContinuumSuppression/FoxWolfram.h>
14#include <analysis/ContinuumSuppression/CleoCones.h>
16#include <analysis/dataobjects/Particle.h>
17#include <analysis/dataobjects/RestOfEvent.h>
19#include <analysis/utility/CLHEPToROOT.h>
20#include <analysis/utility/PCmsLabTransform.h>
21#include <analysis/utility/ROOTToCLHEP.h>
23#include <analysis/VertexFitting/KFit/MakeMotherKFit.h>
24#include <analysis/VertexFitting/KFit/VertexFitKFit.h>
26#include <framework/dataobjects/Helix.h>
28#include <framework/geometry/BFieldManager.h>
30#include <mdst/dataobjects/TrackFitResult.h>
33#include <Math/Vector3D.h>
49 setDescription(
"Creates for each Particle in the given ParticleLists a ContinuumSuppression dataobject and makes basf2 relation between them.");
55 addParam(
"performIPProfileFit",
m_ipProfileFit,
"switch to turn on vertex fit of tracks with IP profile constraint",
false);
69 B2ERROR(
"The ROE mask for the continuum suppression must not be called " <<
m_ROEMask);
86 for (
unsigned i = 0; i <
m_plist->getListSize(); i++) {
97 particle->addRelationTo(qqVars);
99 std::vector<ROOT::Math::PxPyPzEVector> p_cms_sigB, p_cms_roe, p_cms_all;
101 std::vector<std::pair<ROOT::Math::PxPyPzEVector, int>> p_cms_q_sigA;
102 std::vector<std::pair<ROOT::Math::PxPyPzEVector, int>> p_cms_q_sigB;
103 std::vector<std::pair<ROOT::Math::PxPyPzEVector, int>> p_cms_q_roe;
105 std::vector<float> ksfwFS0;
106 std::vector<float> ksfwFS1;
108 std::vector<float> cleoConesAll;
109 std::vector<float> cleoConesROE;
113 ROOT::Math::XYZVector thrustB;
114 ROOT::Math::XYZVector thrustO;
127 ROOT::Math::PxPyPzEVector p_cms_missA(0, 0, 0, 2 * BeamEnergy);
128 ROOT::Math::PxPyPzEVector p_cms_missB(0, 0, 0, 2 * BeamEnergy);
132 std::vector<Belle2::Particle*> signalDaughters = particle->getDaughters();
135 ROOT::Math::PxPyPzEVector p_cms = T.rotateLabToCms() * sigFS0->get4Vector();
137 p_cms_q_sigA.emplace_back(p_cms, sigFS0->getCharge());
139 p_cms_missA -= p_cms;
144 std::vector<const Belle2::Particle*> signalFSParticles = particle->getFinalStateDaughters();
147 ROOT::Math::PxPyPzEVector p_cms = T.rotateLabToCms() * sigFS1->get4Vector();
149 p_cms_all.push_back(p_cms);
150 p_cms_sigB.push_back(p_cms);
152 p_cms_q_sigB.emplace_back(p_cms, sigFS1->getCharge());
154 p_cms_missB -= p_cms;
167 for (
const Particle* chargedROEParticle : chargedROEParticles) {
169 ROOT::Math::PxPyPzEVector p =
ipProfileFit(chargedROEParticle);
171 ROOT::Math::PxPyPzEVector p_cms = T.rotateLabToCms() * p;
173 p_cms_all.push_back(p_cms);
174 p_cms_roe.push_back(p_cms);
176 p_cms_q_roe.emplace_back(p_cms, chargedROEParticle->getCharge());
178 p_cms_missA -= p_cms;
179 p_cms_missB -= p_cms;
188 for (
const Particle* photon : roePhotons) {
192 ROOT::Math::PxPyPzEVector p_cms = T.rotateLabToCms() * photon->get4Vector();
193 p_cms_all.push_back(p_cms);
194 p_cms_roe.push_back(p_cms);
196 p_cms_q_roe.emplace_back(p_cms, photon->getCharge());
198 p_cms_missA -= p_cms;
199 p_cms_missB -= p_cms;
208 thrustBm = thrustB.R();
209 thrustOm = thrustO.R();
210 cosTBTO = fabs(thrustB.Unit().Dot(thrustO.Unit()));
211 cosTBz = fabs(cos(thrustB.Theta()));
214 CleoCones cc(p_cms_all, p_cms_roe, thrustB,
true,
true);
224 ROOT::Math::PxPyPzEVector p_cms_B = T.rotateLabToCms() * particle->get4Vector();
225 double Hso0_max(2 * (2 * BeamEnergy - p_cms_B.E()));
235 ksfwFS0.push_back(KsfwM.
mm2());
236 ksfwFS0.push_back(KsfwM.
et());
237 ksfwFS0.push_back(KsfwM.
Hso(0, 0));
238 ksfwFS0.push_back(KsfwM.
Hso(0, 1));
239 ksfwFS0.push_back(KsfwM.
Hso(0, 2));
240 ksfwFS0.push_back(KsfwM.
Hso(0, 3));
241 ksfwFS0.push_back(KsfwM.
Hso(0, 4));
242 ksfwFS0.push_back(KsfwM.
Hso(1, 0));
243 ksfwFS0.push_back(KsfwM.
Hso(1, 2));
244 ksfwFS0.push_back(KsfwM.
Hso(1, 4));
245 ksfwFS0.push_back(KsfwM.
Hso(2, 0));
246 ksfwFS0.push_back(KsfwM.
Hso(2, 2));
247 ksfwFS0.push_back(KsfwM.
Hso(2, 4));
248 ksfwFS0.push_back(KsfwM.
Hoo(0));
249 ksfwFS0.push_back(KsfwM.
Hoo(1));
250 ksfwFS0.push_back(KsfwM.
Hoo(2));
251 ksfwFS0.push_back(KsfwM.
Hoo(3));
252 ksfwFS0.push_back(KsfwM.
Hoo(4));
255 ksfwFS1.push_back(KsfwM.
mm2());
256 ksfwFS1.push_back(KsfwM.
et());
257 ksfwFS1.push_back(KsfwM.
Hso(0, 0));
258 ksfwFS1.push_back(KsfwM.
Hso(0, 1));
259 ksfwFS1.push_back(KsfwM.
Hso(0, 2));
260 ksfwFS1.push_back(KsfwM.
Hso(0, 3));
261 ksfwFS1.push_back(KsfwM.
Hso(0, 4));
262 ksfwFS1.push_back(KsfwM.
Hso(1, 0));
263 ksfwFS1.push_back(KsfwM.
Hso(1, 2));
264 ksfwFS1.push_back(KsfwM.
Hso(1, 4));
265 ksfwFS1.push_back(KsfwM.
Hso(2, 0));
266 ksfwFS1.push_back(KsfwM.
Hso(2, 2));
267 ksfwFS1.push_back(KsfwM.
Hso(2, 4));
268 ksfwFS1.push_back(KsfwM.
Hoo(0));
269 ksfwFS1.push_back(KsfwM.
Hoo(1));
270 ksfwFS1.push_back(KsfwM.
Hoo(2));
271 ksfwFS1.push_back(KsfwM.
Hoo(3));
272 ksfwFS1.push_back(KsfwM.
Hoo(4));
300 const TrackFitResult* trackFitResult = particle->getTrackFitResult();
301 if (!trackFitResult)
return particle->get4Vector();
308 if (std::abs(helix.getD0()) > 5.0 || std::abs(helix.getZ0()) > 10.0)
return particle->get4Vector();
311 ROOT::Math::XYZVector p_helix = helix.getMomentum(
m_Bfield);
312 double shiftedParticleMass = particle->getPDGMass();
313 double shiftedParticleEnergy = std::sqrt(p_helix.Mag2() + shiftedParticleMass * shiftedParticleMass);
314 ROOT::Math::PxPyPzEVector momentumOfShiftedParticle(p_helix.X(), p_helix.y(), p_helix.Z(), shiftedParticleEnergy);
315 const Particle shiftedParticle(momentumOfShiftedParticle, particle->getPDGCode());
static ROOT::Math::XYZVector getFieldInTesla(const ROOT::Math::XYZVector &pos)
return the magnetic field at a given position in Tesla.
Class to calculate the Cleo clone variables.
std::vector< float > cleo_cone_with_roe()
Returns calculated Cleo Cones constructed from only ROE tracks.
std::vector< float > cleo_cone_with_all()
Returns calculated Cleo Cones constructed from all tracks.
TMatrixDSym m_beamSpotCov
Beam spot covariance matrix.
ContinuumSuppressionBuilderModule()
constructor
std::string m_ROEMask
ROE mask.
ROOT::Math::PxPyPzEVector ipProfileFit(const Particle *particle)
get 4Vector for CS calculation with or without IP profile fit
std::string m_particleListName
Name of the ParticleList.
virtual void initialize() override
initialize the module (setup the data store)
virtual void event() override
process event
bool m_ipProfileFit
Switch to turn on / off vertex fit of tracks with IP profile constraint.
ROOT::Math::XYZVector m_BeamSpotCenter
Beam spot position.
DBObjPtr< BeamSpot > m_beamSpotDB
Beam spot database object.
void addContinuumSuppression(const Particle *particle)
calculate continuum suppression quantities
StoreArray< ContinuumSuppression > m_csarray
StoreArray of ContinuumSuppression.
double m_Bfield
magnetic field from data base
StoreObjPtr< ParticleList > m_plist
input particle list
This is a class for collecting variables used in continuum suppression.
void addR2(float R2)
Add reduced Fox-Wolfram moment R2.
void addThrustBm(float thrustBm)
Add magnitude of B thrust axis.
void addCleoConesALL(const std::vector< float > &cleoConesALL)
Add vector of Cleo Cones constructed of all final state particles.
void addCleoConesROE(const std::vector< float > &cleoConesROE)
Add vector of Cleo Cones constructed of only ROE particles.
void addThrustB(const ROOT::Math::XYZVector &thrustB)
Add ROE thrust axis.
void addCosTBz(float cosTBz)
Add cosine of the angle between the thrust axis of the B and the z-axis.
void addThrustOm(float thrustOm)
Add magnitude of ROE thrust axis.
void addKsfwFS0(const std::vector< float > &ksfwFS0)
Add vector of KSFW moments, Et, and mm2 for final state = 0.
void addThrustO(const ROOT::Math::XYZVector &thrustO)
Add ROE thrust axis.
void addKsfwFS1(const std::vector< float > &ksfwFS1)
Add vector of KSFW moments, Et, and mm2 for final state = 1.
void addCosTBTO(float cosTBTO)
Add cosine of the angle between the thrust axis of the B and the thrust axis of the ROE.
@ c_nPhotons
CR is split into n photons (N1)
Class to calculate the Fox-Wolfram moments up to order 8.
double getR(int i) const
Returns the i-th moment normalized to the 0th-order moment.
void calculateBasicMoments()
Method to perform the calculation of the moments up to order 4, which are the most relevant ones.
Moment-calculation of the k_sfw improved Super-Fox-Wolfram moments.
int usefinal(int uf)
Sets the flag that specifiies we are using the finalstate for signal.
double et(int uf=-1) const
Returns calculated transverse energy.
double Hso(int i, int j, int uf=-1) const
Returns calculated KSFW Moments.
double mm2(int uf=-1) const
Returns calculated missing mass squared.
double Hoo(int i, int uf=-1) const
Returns calculated KSFW Moments.
void setDescription(const std::string &description)
Sets the description of the module.
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Class to store reconstructed particles.
This is a general purpose class for collecting reconstructed MDST data objects that are not used in r...
std::vector< const Particle * > getChargedParticles(const std::string &maskName=c_defaultMaskName, unsigned int pdg=0, bool unpackComposite=true) const
Get charged particles from ROE mask.
static constexpr const char * c_defaultMaskName
Default mask name.
std::vector< const Particle * > getPhotons(const std::string &maskName=c_defaultMaskName, bool unpackComposite=true) const
Get photons from ROE mask.
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.
Accessor to arrays stored in the data store.
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
static ROOT::Math::XYZVector calculateThrust(const std::vector< ROOT::Math::PxPyPzEVector > &momenta)
calculates the thrust axis
Values of the result of a track fit with a given particle hypothesis.
Helix getHelix() const
Conversion to framework Helix (without covariance).
const CLHEP::HepSymMatrix getTrackError(const int id) const
Get an error matrix of the track.
const CLHEP::HepLorentzVector getTrackMomentum(const int id) const
Get a Lorentz vector of the track.
const HepPoint3D getTrackPosition(const int id) const
Get a position of the track.
enum KFitError::ECode setMagneticField(const double mf)
Change a magnetic field from the default value KFitConst::kDefaultMagneticField.
const KFitTrack getTrack(const int id) const
Get a specified track object.
enum KFitError::ECode addParticle(const Particle *particle)
Add a particle to the fitter.
ECode
ECode is a error code enumerate.
double getCharge(void) const
Get a charge of the track.
MakeMotherKFit is a class to build mother particle from kinematically fitted daughters.
enum KFitError::ECode setVertex(const HepPoint3D &v)
Set a vertex position of the mother particle.
enum KFitError::ECode addTrack(const KFitTrack &kp)
Add a track to the make-mother object.
enum KFitError::ECode doMake(void)
Perform a reconstruction of mother particle.
enum KFitError::ECode setVertexError(const CLHEP::HepSymMatrix &e)
Set a vertex error matrix of the mother particle.
enum KFitError::ECode setTrackVertexError(const CLHEP::HepMatrix &e)
Set a vertex error matrix of the child particle in the addTrack'ed order.
const CLHEP::HepLorentzVector getMotherMomentum(void) const
Get a Lorentz vector of the mother particle.
enum KFitError::ECode setMagneticField(const double mf)
Change a magnetic field from the default value KFitConst::kDefaultMagneticField.
VertexFitKFit is a derived class from KFitBase to perform vertex-constraint kinematical fit.
const CLHEP::HepSymMatrix getVertexError(void) const
Get a fitted vertex error matrix.
enum KFitError::ECode doFit(void)
Perform a vertex-constraint fit.
enum KFitError::ECode setIpProfile(const HepPoint3D &ip, const CLHEP::HepSymMatrix &ipe)
Set an IP-ellipsoid shape for the vertex constraint fit.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
const CLHEP::HepMatrix getTrackVertexError(const int id) const
Get a vertex error matrix of the track.
void addParam(const std::string &name, T ¶mVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Abstract base class for different kinds of events.