1212// / \brief this is a starting point for the Resonances tutorial
1313// / \author junlee kim
1414
15+ #include " PWGLF/DataModel/ReducedHeptaQuarkTables.h"
16+
17+ #include " Common/Core/trackUtilities.h"
18+
19+ #include " CommonConstants/PhysicsConstants.h"
20+ #include " Framework/ASoAHelpers.h"
21+ #include " Framework/AnalysisDataModel.h"
22+ #include " Framework/AnalysisTask.h"
23+ #include " Framework/StepTHn.h"
24+ #include " Framework/runDataProcessing.h"
1525#include < Framework/Configurable.h>
16- # include < TLorentzVector.h >
26+
1727#include < Math/GenVector/Boost.h>
18- #include < Math/Vector4D.h>
1928#include < Math/Vector3D.h>
29+ #include < Math/Vector4D.h>
30+ #include < TLorentzVector.h>
2031#include < TMath.h>
2132#include < TRandom3.h>
33+
2234#include < fairlogger/Logger.h>
35+
36+ #include < algorithm>
2337#include < iostream>
2438#include < iterator>
2539#include < string>
2640#include < vector>
2741
28- #include " Framework/AnalysisTask.h"
29- #include " Framework/ASoAHelpers.h"
30- #include " Framework/runDataProcessing.h"
31- #include " Framework/AnalysisDataModel.h"
32- #include " Framework/StepTHn.h"
33- #include " Common/Core/trackUtilities.h"
34- #include " PWGLF/DataModel/ReducedHeptaQuarkTables.h"
35- #include " CommonConstants/PhysicsConstants.h"
36-
3742using namespace o2 ;
3843using namespace o2 ::framework;
3944using namespace o2 ::framework::expressions;
@@ -47,14 +52,22 @@ struct heptaquark {
4752 Configurable<int > cfgPIDStrategy{" cfgPIDStrategy" , 3 , " PID strategy 1" };
4853 Configurable<float > cfgPIDPrPi{" cfgPIDPrPi" , 3 , " PID selection for proton and pion" };
4954
50- Configurable<float > minPhiMass{ " minPhiMass " , 1.01 , " Minimum phi mass" };
51- Configurable<float > maxPhiMass{ " maxPhiMass " , 1.03 , " Maximum phi mass" };
55+ Configurable<float > cfgMinPhiMass{ " cfgMinPhiMass " , 1.01 , " Minimum phi mass" };
56+ Configurable<float > cfgMaxPhiMass{ " cfgMaxPhiMass " , 1.03 , " Maximum phi mass" };
5257
53- Configurable<float > minLambdaMass{ " minLambdaMass " , 1.1 , " Minimum lambda mass" };
54- Configurable<float > maxLambdaMass{ " maxLambdaMass " , 1.13 , " Maximum lambda mass" };
58+ Configurable<float > cfgMinLambdaMass{ " cfgMinLambdaMass " , 1.1 , " Minimum lambda mass" };
59+ Configurable<float > cfgMaxLambdaMass{ " cfgMaxLambdaMass " , 1.13 , " Maximum lambda mass" };
5560
56- Configurable<float > cutNsigmaTPC{" cutNsigmaTPC" , 2.5 , " nsigma cut TPC" };
57- Configurable<float > cutNsigmaTOF{" cutNsigmaTOF" , 3.0 , " nsigma cut TOF" };
61+ Configurable<float > cfgNsigmaTPC{" cfgNsigmaTPC" , 2.5 , " nsigma cut TPC" };
62+ Configurable<float > cfgNsigmaTOF{" cfgNsigmaTOF" , 3.0 , " nsigma cut TOF" };
63+
64+ Configurable<bool > cfgSelectHQ{" cfgSelectHQ" , true , " switch to select HQ" };
65+
66+ Configurable<float > cfgMinPhiPt{" cfgMinPhiPt" , 0.2 , " Minimum phi pt" };
67+ Configurable<float > cfgMinLambdaPt{" cfgMinLambdaPt" , 0.5 , " Minimum lambda pt" };
68+
69+ Configurable<float > cfgSoftFraction{" cfgSoftFraction" , 0.01 , " Minimum allowed softest fraction" };
70+ Configurable<float > cfgCollinear{" cfgCollinear" , 0.98 , " Maximum allowed collinear selection" };
5871
5972 ConfigurableAxis massAxis{" massAxis" , {600 , 2.8 , 3.4 }, " Invariant mass axis" };
6073 ConfigurableAxis ptAxis{" ptAxis" , {VARIABLE_WIDTH, 0.2 , 0.5 , 1.0 , 1.5 , 2.0 , 2.5 , 3.0 , 4.0 , 5.0 , 6.5 , 8.0 , 10.0 , 100.0 }, " Transverse momentum bins" };
@@ -95,105 +108,115 @@ struct heptaquark {
95108 {
96109 if (PIDStrategy == 0 ) {
97110 if (TOFHit != 1 ) {
98- if (TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
111+ if (TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
99112 return true ;
100113 }
101114 }
102115 if (TOFHit == 1 ) {
103- if (TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
116+ if (TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
104117 return true ;
105118 }
106119 }
107120 }
108121 if (PIDStrategy == 1 ) {
109122 if (ptcand < 0.5 ) {
110- if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
123+ if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
111124 return true ;
112125 }
113- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
126+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
114127 return true ;
115128 }
116129 }
117130 if (ptcand >= 0.5 ) {
118- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
131+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
119132 return true ;
120133 }
121134 }
122135 }
123136 if (PIDStrategy == 2 ) {
124137 if (ptcand < 0.5 ) {
125- if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
138+ if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
126139 return true ;
127140 }
128- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
141+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
129142 return true ;
130143 }
131144 }
132145 if (ptcand >= 0.5 && ptcand < 1.2 ) {
133- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
146+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
134147 return true ;
135148 }
136- if (TOFHit != 1 && nsigmaTPC > -1.5 && nsigmaTPC < cutNsigmaTPC ) {
149+ if (TOFHit != 1 && nsigmaTPC > -1.5 && nsigmaTPC < cfgNsigmaTPC ) {
137150 return true ;
138151 }
139152 }
140153 if (ptcand >= 1.2 ) {
141- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
154+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
142155 return true ;
143156 }
144- if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
157+ if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
145158 return true ;
146159 }
147160 }
148161 }
149162 if (PIDStrategy == 3 ) {
150163 if (ptcand < 0.5 ) {
151- if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
164+ if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
152165 return true ;
153166 }
154- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
167+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
155168 return true ;
156169 }
157170 }
158171 if (ptcand >= 0.5 && ptcand < 1.2 ) {
159- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
172+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
160173 return true ;
161174 }
162175 }
163176 if (ptcand >= 1.2 ) {
164- if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cutNsigmaTOF ) {
177+ if (TOFHit == 1 && TMath::Abs (nsigmaTOF) < cfgNsigmaTOF ) {
165178 return true ;
166179 }
167- if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cutNsigmaTPC ) {
180+ if (TOFHit != 1 && TMath::Abs (nsigmaTPC) < cfgNsigmaTPC ) {
168181 return true ;
169182 }
170183 }
171184 }
172185 return false ;
173186 }
174187
175- template <typename V01, typename V02>
176- ROOT::Math::XYZVector getDCAofV0V0 (V01 const & v01, V02 const & v02)
177- {
178- ROOT::Math::XYZVector v01pos, v02pos, v01mom, v02mom;
179- v01pos.SetXYZ (v01.x (), v01.y (), v01.z ());
180- v02pos.SetXYZ (v02.x (), v02.y (), v02.z ());
181- v01mom.SetXYZ (v01.px (), v01.py (), v01.pz ());
182- v02mom.SetXYZ (v02.px (), v02.py (), v02.pz ());
183-
184- ROOT::Math::XYZVector posdiff = v02pos - v01pos;
185- ROOT::Math::XYZVector cross = v01mom.Cross (v02mom);
186- ROOT::Math::XYZVector dcaVec = (posdiff.Dot (cross) / cross.Mag2 ()) * cross;
187- return dcaVec;
188- }
189-
190- template <typename V01, typename V02>
191- float getCPA (V01 const & v01, V02 const & v02)
188+ template <typename HQ1, typename HQ2, typename HQ3>
189+ int selectHQ (HQ1 const & hq1, HQ2 const & hq2, HQ3 const & hq3)
192190 {
193- ROOT::Math::XYZVector v01mom, v02mom;
194- v01mom.SetXYZ (v01.px () / v01.p (), v01.py () / v01.p (), v01.pz () / v01.p ());
195- v02mom.SetXYZ (v02.px () / v02.p (), v02.py () / v02.p (), v02.pz () / v02.p ());
196- return v01mom.Dot (v02mom);
191+ int selection = 0 ;
192+ if (hq1.Pt () < cfgMinPhiPt || hq2.Pt () < cfgMinPhiPt || hq3.Pt () < cfgMinLambdaPt)
193+ selection += 1 ;
194+
195+ double sumE = hq1.E () + hq2.E () + hq3.E ();
196+ double emin = std::min ({hq1.E (), hq2.E (), hq3.E ()});
197+ double fmin = emin / std::max (1e-9 , sumE);
198+ if (fmin < cfgSoftFraction)
199+ selection += 2 ;
200+
201+ auto ex = hq1 + hq2 + hq3;
202+ TVector3 boost = -ex.BoostVector ();
203+ auto hqphipair_boost = hq1 + hq2;
204+ auto hqlambda_boost = hq3;
205+ hqphipair_boost.Boost (boost);
206+ hqlambda_boost.Boost (boost);
207+ double cosHel = hqlambda_boost.Vect ().Dot (hqphipair_boost.Vect ()) / (hqlambda_boost.Vect ().Mag () * hqphipair_boost.Vect ().Mag ());
208+ if (std::abs (cosHel) > cfgCollinear)
209+ selection += 4 ;
210+ /*
211+ ROOT::Math::XYZVector rPV(col.posX(), col.posY(), col.posZ());
212+ ROOT::Math::XYZVector rSV(hq3.hqx(), hq3.hqy(), hq3.hqz());
213+ ROOT::Math::XYZVector L = rSV - rPV;
214+ ROOT::Math::XYZVector exMom(ex.Px(), ex.Py(), ex.Pz());
215+ double cosPoint = L.Dot(exMom) / (L.R() * pEx.R() + 1e-9);
216+ if (cosPoint < cfgCosPoint)
217+ return 8;
218+ */
219+ return selection;
197220 }
198221
199222 ROOT::Math::PxPyPzMVector DauVec1, DauVec2;
@@ -214,7 +237,7 @@ struct heptaquark {
214237 if (hqtrackd1.hqId () != 333 )
215238 continue ;
216239
217- if (hqtrackd1.hqMass () < minPhiMass || hqtrackd1.hqMass () > maxPhiMass )
240+ if (hqtrackd1.hqMass () < cfgMinPhiMass || hqtrackd1.hqMass () > cfgMaxPhiMass )
218241 continue ;
219242
220243 DauVec1 = ROOT::Math::PxPyPzMVector (hqtrackd1.hqd1Px (), hqtrackd1.hqd1Py (), hqtrackd1.hqd1Pz (), massKa);
@@ -246,7 +269,7 @@ struct heptaquark {
246269 if (hqtrackd2.hqId () != 333 )
247270 continue ;
248271
249- if (hqtrackd2.hqMass () < minPhiMass || hqtrackd2.hqMass () > maxPhiMass )
272+ if (hqtrackd2.hqMass () < cfgMinPhiMass || hqtrackd2.hqMass () > cfgMaxPhiMass )
250273 continue ;
251274
252275 DauVec1 = ROOT::Math::PxPyPzMVector (hqtrackd2.hqd1Px (), hqtrackd2.hqd1Py (), hqtrackd2.hqd1Pz (), massKa);
@@ -278,7 +301,7 @@ struct heptaquark {
278301 if (std::abs (hqtrackd3.hqId ()) != 3122 )
279302 continue ;
280303
281- if (hqtrackd3.hqMass () < minLambdaMass || hqtrackd3.hqMass () > maxLambdaMass )
304+ if (hqtrackd3.hqMass () < cfgMinLambdaMass || hqtrackd3.hqMass () > cfgMaxLambdaMass )
282305 continue ;
283306
284307 int isLambda = static_cast <int >(hqtrackd3.hqId () < 0 );
@@ -321,6 +344,9 @@ struct heptaquark {
321344 HQ12 = HQ1 + HQ2;
322345 HQ13 = HQ1 + HQ3;
323346
347+ if (cfgSelectHQ && selectHQ (HQ1, HQ2, HQ3))
348+ continue ;
349+
324350 histos.fill (HIST (" h_InvMass_same" ), exotic.M (), exotic.Pt (), collision.centrality ());
325351 histos.fill (HIST (" hDalitz" ), HQ12.M2 (), HQ13.M2 (), exotic.M (), exotic.Pt (), isLambda, collision.centrality ());
326352
0 commit comments