2020#include < map>
2121#include < tuple>
2222
23- #include " TString.h"
23+ // #include "TString.h"
2424#include " Math/Vector4D.h"
2525#include " Framework/runDataProcessing.h"
2626#include " Framework/AnalysisTask.h"
2727#include " Framework/ASoAHelpers.h"
28- #include " Common/Core/RecoDecay.h"
28+ // #include "Common/Core/RecoDecay.h"
2929#include " Common/Core/trackUtilities.h"
3030
3131#include " DetectorsBase/Propagator.h"
3939#include " PWGEM/PhotonMeson/Core/V0PhotonCut.h"
4040#include " PWGEM/PhotonMeson/Core/EMPhotonEventCut.h"
4141#include " PWGEM/PhotonMeson/Core/DalitzEECut.h"
42+ #include " PWGEM/Dilepton/Utils/PairUtilities.h"
4243
4344using namespace o2 ;
4445using namespace o2 ::aod;
@@ -61,6 +62,12 @@ struct prefilterPhoton {
6162 Produces<aod::EMPrimaryElectronsPrefilterBitDerived> pfb_ele_derived;
6263
6364 // Configurables
65+ Configurable<std::string> ccdburl{" ccdb-url" , " http://alice-ccdb.cern.ch" , " url of the ccdb repository" };
66+ Configurable<std::string> grpPath{" grpPath" , " GLO/GRP/GRP" , " Path of the grp file" };
67+ Configurable<std::string> grpmagPath{" grpmagPath" , " GLO/Config/GRPMagField" , " CCDB path of the GRPMagField object" };
68+ Configurable<bool > skipGRPOquery{" skipGRPOquery" , true , " skip grpo query" };
69+ Configurable<float > d_bz_input{" d_bz_input" , -999 , " bz field in kG, -999 is automatic" };
70+
6471 Configurable<int > cfgCentEstimator{" cfgCentEstimator" , 2 , " FT0M:0, FT0A:1, FT0C:2" };
6572 Configurable<float > cfgCentMin{" cfgCentMin" , -1 , " min. centrality" };
6673 Configurable<float > cfgCentMax{" cfgCentMax" , 999 .f , " max. centrality" };
@@ -115,7 +122,7 @@ struct prefilterPhoton {
115122
116123 Configurable<float > cfg_min_mee{" cfg_min_mee" , 0.0 , " min mass" };
117124 Configurable<float > cfg_max_mee{" cfg_max_mee" , 0.02 , " max mass" };
118- Configurable<bool > cfg_apply_phiv{" cfg_apply_phiv" , false , " flag to apply phiv cut" };
125+ // Configurable<bool> cfg_apply_phiv{"cfg_apply_phiv", false, "flag to apply phiv cut"};
119126 Configurable<bool > cfg_apply_pf{" cfg_apply_pf" , false , " flag to apply phiv prefilter" };
120127 Configurable<bool > cfg_require_itsib_any{" cfg_require_itsib_any" , false , " flag to require ITS ib any hits" };
121128 Configurable<bool > cfg_require_itsib_1st{" cfg_require_itsib_1st" , true , " flag to require ITS ib 1st hit" };
@@ -139,6 +146,8 @@ struct prefilterPhoton {
139146 Configurable<float > cfg_max_TPCNsigmaEl{" cfg_max_TPCNsigmaEl" , +3.0 , " max. TPC n sigma for electron inclusion" };
140147 Configurable<float > cfg_min_TPCNsigmaPi{" cfg_min_TPCNsigmaPi" , 0.0 , " min. TPC n sigma for pion exclusion" };
141148 Configurable<float > cfg_max_TPCNsigmaPi{" cfg_max_TPCNsigmaPi" , 0.0 , " max. TPC n sigma for pion exclusion" };
149+ Configurable<float > cfg_min_TOFNsigmaEl{" cfg_min_TOFNsigmaEl" , -3.0 , " min. TOF n sigma for electron inclusion" };
150+ Configurable<float > cfg_max_TOFNsigmaEl{" cfg_max_TOFNsigmaEl" , +3.0 , " max. TOF n sigma for electron inclusion" };
142151 } dileptoncuts;
143152
144153 struct : ConfigurableGroup {
@@ -155,24 +164,79 @@ struct prefilterPhoton {
155164
156165 HistogramRegistry fRegistry {" output" , {}, OutputObjHandlingPolicy::AnalysisObject, false , false };
157166
167+ o2::ccdb::CcdbApi ccdbApi;
168+ Service<o2::ccdb::BasicCCDBManager> ccdb;
169+ int mRunNumber ;
170+ float d_bz;
171+
158172 void init (InitContext& /* context*/ )
159173 {
160174 DefineEMEventCut ();
161175 DefinePCMCut ();
162176 addhistograms ();
177+
178+ mRunNumber = 0 ;
179+ d_bz = 0 ;
180+
181+ ccdb->setURL (ccdburl);
182+ ccdb->setCaching (true );
183+ ccdb->setLocalObjectValidityChecking ();
184+ ccdb->setFatalWhenNull (false );
163185 }
164186
165187 ~prefilterPhoton () {}
166188
189+ template <typename TCollision>
190+ void initCCDB (TCollision const & collision)
191+ {
192+ if (mRunNumber == collision.runNumber ()) {
193+ return ;
194+ }
195+
196+ // In case override, don't proceed, please - no CCDB access required
197+ if (d_bz_input > -990 ) {
198+ d_bz = d_bz_input;
199+ o2::parameters::GRPMagField grpmag;
200+ if (fabs (d_bz) > 1e-5 ) {
201+ grpmag.setL3Current (30000 .f / (d_bz / 5 .0f ));
202+ }
203+ mRunNumber = collision.runNumber ();
204+ return ;
205+ }
206+
207+ auto run3grp_timestamp = collision.timestamp ();
208+ o2::parameters::GRPObject* grpo = 0x0 ;
209+ o2::parameters::GRPMagField* grpmag = 0x0 ;
210+ if (!skipGRPOquery)
211+ grpo = ccdb->getForTimeStamp <o2::parameters::GRPObject>(grpPath, run3grp_timestamp);
212+ if (grpo) {
213+ // Fetch magnetic field from ccdb for current collision
214+ d_bz = grpo->getNominalL3Field ();
215+ LOG (info) << " Retrieved GRP for timestamp " << run3grp_timestamp << " with magnetic field of " << d_bz << " kZG" ;
216+ } else {
217+ grpmag = ccdb->getForTimeStamp <o2::parameters::GRPMagField>(grpmagPath, run3grp_timestamp);
218+ if (!grpmag) {
219+ LOG (fatal) << " Got nullptr from CCDB for path " << grpmagPath << " of object GRPMagField and " << grpPath << " of object GRPObject for timestamp " << run3grp_timestamp;
220+ }
221+ // Fetch magnetic field from ccdb for current collision
222+ d_bz = std::lround (5 .f * grpmag->getL3Current () / 30000 .f );
223+ LOG (info) << " Retrieved GRP for timestamp " << run3grp_timestamp << " with magnetic field of " << d_bz << " kZG" ;
224+ }
225+ mRunNumber = collision.runNumber ();
226+ }
227+
167228 void addhistograms ()
168229 {
169230 const AxisSpec axis_mass{200 , 0 , 0.8 , " m_{#gamma#gamma} (GeV/c^{2})" };
170231 const AxisSpec axis_pair_pt{100 , 0 , 10 , " p_{T,#gamma#gamma} (GeV/c)" };
232+ const AxisSpec axis_phiv{180 , 0 , M_PI, " #varphi_{V} (rad.)" };
171233
172234 // for pair
173235 fRegistry .add (" Pair/PCMPCM/before/hMvsPt" , " m_{#gamma#gamma} vs. p_{T,#gamma#gamma}" , kTH2D , {axis_mass, axis_pair_pt}, true );
236+ fRegistry .add (" Pair/PCMDalitzEE/before/hMvsPt" , " m_{ee#gamma} vs. p_{T,ee#gamma}" , kTH2D , {axis_mass, axis_pair_pt}, true );
237+ fRegistry .add (" Pair/PCMDalitzEE/before/hMvsPhiV" , " m_{ee} vs. #varphi_{V}" , kTH2D , {{180 , 0 , M_PI}, {100 , 0 , 0.1 }}, true );
174238 fRegistry .addClone (" Pair/PCMPCM/before/" , " Pair/PCMPCM/after/" );
175- fRegistry .addClone (" Pair/PCMPCM/ " , " Pair/PCMDalitzEE/" );
239+ fRegistry .addClone (" Pair/PCMDalitzEE/before/ " , " Pair/PCMDalitzEE/after /" );
176240 }
177241
178242 void DefineEMEventCut ()
@@ -224,7 +288,7 @@ struct prefilterPhoton {
224288 // for pair
225289 fDileptonCut .SetMeeRange (dileptoncuts.cfg_min_mee , dileptoncuts.cfg_max_mee );
226290 fDileptonCut .SetMaxPhivPairMeeDep ([&](float mll) { return (mll - dileptoncuts.cfg_phiv_intercept ) / dileptoncuts.cfg_phiv_slope ; });
227- fDileptonCut .ApplyPhiV (dileptoncuts. cfg_apply_phiv );
291+ fDileptonCut .ApplyPhiV (false );
228292 fDileptonCut .RequireITSibAny (dileptoncuts.cfg_require_itsib_any );
229293 fDileptonCut .RequireITSib1st (dileptoncuts.cfg_require_itsib_1st );
230294
@@ -244,6 +308,7 @@ struct prefilterPhoton {
244308 fDileptonCut .SetPIDScheme (dileptoncuts.cfg_pid_scheme );
245309 fDileptonCut .SetTPCNsigmaElRange (dileptoncuts.cfg_min_TPCNsigmaEl , dileptoncuts.cfg_max_TPCNsigmaEl );
246310 fDileptonCut .SetTPCNsigmaPiRange (dileptoncuts.cfg_min_TPCNsigmaPi , dileptoncuts.cfg_max_TPCNsigmaPi );
311+ fDileptonCut .SetTOFNsigmaElRange (dileptoncuts.cfg_min_TOFNsigmaEl , dileptoncuts.cfg_max_TOFNsigmaEl );
247312 }
248313
249314 template <PairType pairtype, typename TCollisions, typename TPhotons1, typename TPhotons2, typename TSubInfos1, typename TSubInfos2, typename TPreslice1, typename TPreslice2, typename TCut1, typename TCut2>
@@ -268,6 +333,7 @@ struct prefilterPhoton {
268333
269334 if constexpr (pairtype == PairType::kPCMPCM ) {
270335 for (const auto & collision : collisions) {
336+ initCCDB (collision);
271337 const float centralities[3 ] = {collision.centFT0M (), collision.centFT0A (), collision.centFT0C ()};
272338 bool is_cent_ok = true ;
273339 if (centralities[cfgCentEstimator] < cfgCentMin || cfgCentMax < centralities[cfgCentEstimator]) {
@@ -295,13 +361,14 @@ struct prefilterPhoton {
295361 fRegistry .fill (HIST (" Pair/PCMPCM/before/hMvsPt" ), v12.M (), v12.Pt ());
296362
297363 if (ggcuts.cfg_min_mass < v12.M () && v12.M () < ggcuts.cfg_max_mass ) {
298- map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0ggDefault );
299- map_pfb_v0[g2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0ggDefault );
364+ map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0gg );
365+ map_pfb_v0[g2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0gg );
300366 }
301367 } // end of 2photon pairing loop
302368 } // end of collision loop
303369 } else if constexpr (pairtype == PairType::kPCMDalitzEE ) {
304370 for (const auto & collision : collisions) {
371+ initCCDB (collision);
305372 const float centralities[3 ] = {collision.centFT0M (), collision.centFT0A (), collision.centFT0C ()};
306373 bool is_cent_ok = true ;
307374 if (centralities[cfgCentEstimator] < cfgCentMin || cfgCentMax < centralities[cfgCentEstimator]) {
@@ -337,8 +404,8 @@ struct prefilterPhoton {
337404 fRegistry .fill (HIST (" Pair/PCMPCM/before/hMvsPt" ), v12.M (), v12.Pt ());
338405
339406 if (ggcuts.cfg_min_mass < v12.M () && v12.M () < ggcuts.cfg_max_mass ) {
340- map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0ggDefault );
341- map_pfb_v0[g2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0ggDefault );
407+ map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0gg );
408+ map_pfb_v0[g2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0gg );
342409 }
343410 } // end of 2photon pairing loop
344411
@@ -372,12 +439,33 @@ struct prefilterPhoton {
372439 fRegistry .fill (HIST (" Pair/PCMDalitzEE/before/hMvsPt" ), veeg.M (), veeg.Pt ());
373440
374441 if (eegcuts.cfg_min_mass < veeg.M () && veeg.M () < eegcuts.cfg_max_mass ) {
375- map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0eegDefault );
376- map_pfb_ele[pos2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromPi0eegDefault );
377- map_pfb_ele[ele2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromPi0eegDefault );
442+ map_pfb_v0[g1.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::PhotonPrefilterBitDerived::kPhotonFromPi0eeg );
443+ map_pfb_ele[pos2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromPi0eeg );
444+ map_pfb_ele[ele2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromPi0eeg );
378445 }
379446 } // end of dielectron loop
380447 } // end of g1 loop
448+
449+ for (const auto & [pos2, ele2] : combinations (CombinationsFullIndexPolicy (positrons_per_collision, electrons_per_collision))) {
450+ if (pos2.trackId () == ele2.trackId ()) { // this is protection against pairing identical 2 tracks.
451+ continue ;
452+ }
453+
454+ if (!cut2.template IsSelectedTrack <false >(pos2, collision) || !cut2.template IsSelectedTrack <false >(ele2, collision)) {
455+ continue ;
456+ }
457+
458+ ROOT::Math::PtEtaPhiMVector v_pos (pos2.pt (), pos2.eta (), pos2.phi (), o2::constants::physics::MassElectron);
459+ ROOT::Math::PtEtaPhiMVector v_ele (ele2.pt (), ele2.eta (), ele2.phi (), o2::constants::physics::MassElectron);
460+ ROOT::Math::PtEtaPhiMVector v_ee = v_pos + v_ele;
461+ float phiv = o2::aod::pwgem::dilepton::utils::pairutil::getPhivPair (pos2.px (), pos2.py (), pos2.pz (), ele2.px (), ele2.py (), ele2.pz (), pos2.sign (), ele2.sign (), d_bz);
462+ fRegistry .fill (HIST (" Pair/PCMDalitzEE/before/hMvsPhiV" ), phiv, v_ee.M ());
463+
464+ if (v_ee.M () < phiv * dileptoncuts.cfg_phiv_slope + dileptoncuts.cfg_phiv_intercept ) {
465+ map_pfb_ele[pos2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromFakePC );
466+ map_pfb_ele[ele2.globalIndex ()] |= 1 << static_cast <int >(o2::aod::pwgem::photonmeson::utils::pairutil::ElectronPrefilterBitDerived::kElectronFromFakePC );
467+ }
468+ } // end of dielectron loop to reject photon conversion
381469 } // end of collision loop
382470 }
383471
@@ -425,6 +513,7 @@ struct prefilterPhoton {
425513 } // end of collision loop
426514 } else if constexpr (pairtype == PairType::kPCMDalitzEE ) {
427515 for (auto & collision : collisions) {
516+ initCCDB (collision);
428517 const float centralities[3 ] = {collision.centFT0M (), collision.centFT0A (), collision.centFT0C ()};
429518 if (centralities[cfgCentEstimator] < cfgCentMin || cfgCentMax < centralities[cfgCentEstimator]) {
430519 continue ;
@@ -478,11 +567,13 @@ struct prefilterPhoton {
478567 ROOT::Math::PtEtaPhiMVector v_pos (pos2.pt (), pos2.eta (), pos2.phi (), o2::constants::physics::MassElectron);
479568 ROOT::Math::PtEtaPhiMVector v_ele (ele2.pt (), ele2.eta (), ele2.phi (), o2::constants::physics::MassElectron);
480569 ROOT::Math::PtEtaPhiMVector v_ee = v_pos + v_ele;
570+ float phiv = o2::aod::pwgem::dilepton::utils::pairutil::getPhivPair (pos2.px (), pos2.py (), pos2.pz (), ele2.px (), ele2.py (), ele2.pz (), pos2.sign (), ele2.sign (), d_bz);
481571 if (!(dileptoncuts.cfg_min_mee < v_ee.M () && v_ee.M () < dileptoncuts.cfg_max_mee )) {
482572 continue ;
483573 }
484574 ROOT::Math::PtEtaPhiMVector veeg = v_gamma + v_pos + v_ele;
485575 fRegistry .fill (HIST (" Pair/PCMDalitzEE/after/hMvsPt" ), veeg.M (), veeg.Pt ());
576+ fRegistry .fill (HIST (" Pair/PCMDalitzEE/after/hMvsPhiV" ), phiv, v_ee.M ());
486577 } // end of dielectron loop
487578 } // end of g1 loop
488579 } // end of collision loop
0 commit comments