diff --git a/EventFiltering/PWGEM/HeavyNeutralMesonFilter.cxx b/EventFiltering/PWGEM/HeavyNeutralMesonFilter.cxx index 6909d1561ee..013bee26a88 100644 --- a/EventFiltering/PWGEM/HeavyNeutralMesonFilter.cxx +++ b/EventFiltering/PWGEM/HeavyNeutralMesonFilter.cxx @@ -10,38 +10,52 @@ // or submit itself to any jurisdiction. /// /// \file HeavyNeutralMesonFilter.cxx -/// /// \brief This code loops over collisions to filter events contaning heavy neutral mesons (omega or eta') using EMCal clusters and V0s (PCM) -/// /// \author Nicolas Strangmann (nicolas.strangmann@cern.ch) - Goethe University Frankfurt; Maximilian Korwieser (maximilian.korwieser@cern.ch) - Technical University Munich /// -#include -#include -#include - -#include "Math/GenVector/Boost.h" -#include "Math/Vector4D.h" -#include "TMath.h" -#include "TRandom3.h" - +#include "EventFiltering/filterTables.h" +// #include "PWGEM/PhotonMeson/Utils/HNMUtilities.h" +#include "PWGEM/PhotonMeson/Utils/PairUtilities.h" #include "PWGJE/DataModel/EMCALMatchedCollisions.h" -#include "Common/DataModel/PIDResponse.h" -#include "Common/DataModel/PIDResponseITS.h" -#include "fairlogger/Logger.h" -#include "Framework/Configurable.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisTask.h" -#include "Framework/HistogramRegistry.h" -#include "Framework/runDataProcessing.h" -#include "CommonConstants/MathConstants.h" -#include "Common/Core/TrackSelection.h" -#include "Common/DataModel/TrackSelectionTables.h" +#include "Common/CCDB/TriggerAliases.h" +#include "Common/Core/RecoDecay.h" #include "Common/DataModel/EventSelection.h" #include "Common/DataModel/Multiplicity.h" +#include "Common/DataModel/PIDResponseITS.h" +#include "Common/DataModel/PIDResponseTOF.h" +#include "Common/DataModel/PIDResponseTPC.h" +#include "Common/DataModel/TrackSelectionTables.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include // IWYU pragma: keep +#include +#include +#include + +#include + +#include +#include +#include +#include +#include using namespace o2; using namespace o2::framework; @@ -90,21 +104,22 @@ const std::vector pidCutsName{"TPC min", "TPC max", "TPCTOF max", " const std::vector femtoFilterNames{"PPOmega", "PPEtaPrime", "Omegad", "EtaPrimed", "OmegaP", "EtaPrimeP"}; // configs for tracks -const float pidcutsTable[kNFemtoPartners][kNPIDLimits]{ +// these are need [[maybe_unused]] to silence a warning from clangd, since the compiler will inline them directly to the configs down below and then say: Variable 'X' is not needed and will not be emitted +[[maybe_unused]] const float pidcutsTable[kNFemtoPartners][kNPIDLimits]{ {-4.f, 4.f, 4.f, -99.f, 99.f}, {-4.f, 4.f, 4.f, -6.f, 6.f}, {-4.f, 4.f, 4.f, -99.f, 99.f}}; -const float ptcutsTable[kNFemtoPartners][3]{ +[[maybe_unused]] const float ptcutsTable[kNFemtoPartners][3]{ {0.35f, 6.f, 0.75f}, {0.55f, 2.f, 1.2f}, {0.35f, 6.f, 0.75f}}; -const float nClusterMinTPC[1][kNFemtoPartners]{{80.0f, 80.0f, 80.0f}}; -const float nClusterMinITS[1][kNFemtoPartners]{{4, 4, 4}}; +[[maybe_unused]] const float nClusterMinTPC[1][kNFemtoPartners]{{80.0f, 80.0f, 80.0f}}; +[[maybe_unused]] const float nClusterMinITS[1][kNFemtoPartners]{{4, 4, 4}}; -static const float triggerSwitches[1][kNFemtoTriggers]{{1, 1, 1, 1, 1, 1}}; -const float triggerLimits[1][kNFemtoTriggers]{{1.f, 1.f, 1.f, 1.f, 1.f, 1.f}}; +[[maybe_unused]] static const float triggerSwitches[1][kNFemtoTriggers]{{1, 1, 1, 1, 1, 1}}; +[[maybe_unused]] const float triggerLimits[1][kNFemtoTriggers]{{1.f, 1.f, 1.f, 1.f, 1.f, 1.f}}; } // namespace hnmtrigger struct HeavyNeutralMesonFilter { diff --git a/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h b/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h index 8d2f16d34f6..9dd902396b4 100644 --- a/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h +++ b/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h @@ -9,19 +9,13 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// /// \file Pi0EtaToGammaGamma.h /// \brief This code loops over photons and makes pairs for neutral mesons analyses. -/// /// \author D. Sekihata, daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMA_H_ #define PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMA_H_ -#include "PWGEM/Dilepton/Utils/EMTrack.h" -#include "PWGEM/Dilepton/Utils/EMTrackUtilities.h" -#include "PWGEM/Dilepton/Utils/EventMixingHandler.h" #include "PWGEM/PhotonMeson/Core/DalitzEECut.h" #include "PWGEM/PhotonMeson/Core/EMCPhotonCut.h" #include "PWGEM/PhotonMeson/Core/EMPhotonEventCut.h" @@ -31,200 +25,187 @@ #include "PWGEM/PhotonMeson/Utils/EventHistograms.h" #include "PWGEM/PhotonMeson/Utils/NMHistograms.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" +// Dilepton headers +#include "PWGEM/Dilepton/Utils/EMTrack.h" +#include "PWGEM/Dilepton/Utils/EventMixingHandler.h" -#include "Common/Core/RecoDecay.h" - -#include "CCDB/BasicCCDBManager.h" -#include "DataFormatsParameters/GRPMagField.h" -#include "DataFormatsParameters/GRPObject.h" -#include "DetectorsBase/GeometryManager.h" -#include "EMCALBase/Geometry.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" - -#include "Math/AxisAngle.h" -#include "Math/LorentzRotation.h" -#include "Math/Rotation3D.h" -#include "Math/Vector3D.h" -#include "Math/Vector4D.h" -#include "TString.h" +#include "Common/CCDB/TriggerAliases.h" +#include "Common/DataModel/Centrality.h" +#include "Common/DataModel/EventSelection.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include // IWYU pragma: keep +#include #include -#include +#include +#include #include -#include #include #include #include +#include #include #include -using namespace o2; -using namespace o2::aod; -using namespace o2::framework; -using namespace o2::framework::expressions; -using namespace o2::soa; -using namespace o2::aod::pwgem::photonmeson::photonpair; -using namespace o2::aod::pwgem::photon; -using namespace o2::aod::pwgem::dilepton::utils::emtrackutil; -using namespace o2::aod::pwgem::dilepton::utils; - -using MyCollisions = soa::Join; -using MyCollision = MyCollisions::iterator; - -using MyCollisionsWithJJMC = soa::Join; -using MyCollisionWithJJMC = MyCollisionsWithJJMC::iterator; - -using MyV0Photons = soa::Filtered>; -using MyV0Photon = MyV0Photons::iterator; - -using MyPrimaryElectrons = soa::Filtered>; -using MyPrimaryElectron = MyPrimaryElectrons::iterator; - -using MyEMCClusters = soa::Join; -using MyEMCCluster = MyEMCClusters::iterator; - -using MyPHOSClusters = soa::Join; -using MyPHOSCluster = MyPHOSClusters::iterator; +#include -template +template struct Pi0EtaToGammaGamma { - Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; - Configurable grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"}; - Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; - Configurable skipGRPOquery{"skipGRPOquery", true, "skip grpo query"}; - Configurable d_bz_input{"d_bz_input", -999, "bz field in kG, -999 is automatic"}; - Configurable ndiff_bc_mix{"ndiff_bc_mix", 198, "difference in global BC required in mixed events"}; - - Configurable cfgQvecEstimator{"cfgQvecEstimator", 0, "FT0M:0, FT0A:1, FT0C:2"}; - Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; - Configurable cfgOccupancyEstimator{"cfgOccupancyEstimator", 0, "FT0C:0, Track:1"}; - Configurable cfgCentMin{"cfgCentMin", 0, "min. centrality"}; - Configurable cfgCentMax{"cfgCentMax", 999, "max. centrality"}; - Configurable maxY{"maxY", 0.8, "maximum rapidity for reconstructed particles"}; - Configurable cfgDoMix{"cfgDoMix", true, "flag for event mixing"}; - Configurable ndepth{"ndepth", 10, "depth for event mixing"}; - ConfigurableAxis ConfVtxBins{"ConfVtxBins", {VARIABLE_WIDTH, -10.0f, -8.f, -6.f, -4.f, -2.f, 0.f, 2.f, 4.f, 6.f, 8.f, 10.f}, "Mixing bins - z-vertex"}; - ConfigurableAxis ConfCentBins{"ConfCentBins", {VARIABLE_WIDTH, 0.0f, 5.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.f, 999.f}, "Mixing bins - centrality"}; - ConfigurableAxis ConfEPBins{"ConfEPBins", {VARIABLE_WIDTH, -M_PI / 2, -M_PI / 4, 0.0f, +M_PI / 4, +M_PI / 2}, "Mixing bins - event plane angle"}; - ConfigurableAxis ConfOccupancyBins{"ConfOccupancyBins", {VARIABLE_WIDTH, -1, 1e+10}, "Mixing bins - occupancy"}; + o2::framework::Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; + o2::framework::Configurable grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"}; + o2::framework::Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; + o2::framework::Configurable skipGRPOquery{"skipGRPOquery", true, "skip grpo query"}; + o2::framework::Configurable d_bz_input{"d_bz_input", -999, "bz field in kG, -999 is automatic"}; + o2::framework::Configurable ndiff_bc_mix{"ndiff_bc_mix", 198, "difference in global BC required in mixed events"}; + + o2::framework::Configurable cfgQvecEstimator{"cfgQvecEstimator", 0, "FT0M:0, FT0A:1, FT0C:2"}; + o2::framework::Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; + o2::framework::Configurable cfgOccupancyEstimator{"cfgOccupancyEstimator", 0, "FT0C:0, Track:1"}; + o2::framework::Configurable cfgCentMin{"cfgCentMin", 0, "min. centrality"}; + o2::framework::Configurable cfgCentMax{"cfgCentMax", 999, "max. centrality"}; + o2::framework::Configurable maxY{"maxY", 0.8, "maximum rapidity for reconstructed particles"}; + o2::framework::Configurable cfgDoMix{"cfgDoMix", true, "flag for event mixing"}; + o2::framework::Configurable ndepth{"ndepth", 10, "depth for event mixing"}; + o2::framework::ConfigurableAxis ConfVtxBins{"ConfVtxBins", {o2::framework::VARIABLE_WIDTH, -10.0f, -8.f, -6.f, -4.f, -2.f, 0.f, 2.f, 4.f, 6.f, 8.f, 10.f}, "Mixing bins - z-vertex"}; + o2::framework::ConfigurableAxis ConfCentBins{"ConfCentBins", {o2::framework::VARIABLE_WIDTH, 0.0f, 5.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.f, 999.f}, "Mixing bins - centrality"}; + o2::framework::ConfigurableAxis ConfEPBins{"ConfEPBins", {o2::framework::VARIABLE_WIDTH, -o2::constants::math::PIHalf, -o2::constants::math::PIQuarter, 0.0f, +o2::constants::math::PIQuarter, +o2::constants::math::PIHalf}, "Mixing bins - event plane angle"}; + o2::framework::ConfigurableAxis ConfOccupancyBins{"ConfOccupancyBins", {o2::framework::VARIABLE_WIDTH, -1, 1e+10}, "Mixing bins - occupancy"}; EMPhotonEventCut fEMEventCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "eventcut_group"; - Configurable cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"}; - Configurable cfgZvtxMax{"cfgZvtxMax", +10.f, "max. Zvtx"}; - Configurable cfgRequireSel8{"cfgRequireSel8", true, "require sel8 in event cut"}; - Configurable cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND in event cut"}; - Configurable cfgRequireNoTFB{"cfgRequireNoTFB", false, "require No time frame border in event cut"}; - Configurable cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border in event cut"}; - Configurable cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"}; - Configurable cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV. - Configurable cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"}; - Configurable cfgRequireEMCReadoutInMB{"cfgRequireEMCReadoutInMB", false, "require the EMC to be read out in an MB collision (kTVXinEMC)"}; - Configurable cfgRequireEMCHardwareTriggered{"cfgRequireEMCHardwareTriggered", false, "require the EMC to be hardware triggered (kEMC7 or kDMC7)"}; - Configurable cfgTrackOccupancyMin{"cfgTrackOccupancyMin", -2, "min. occupancy"}; - Configurable cfgTrackOccupancyMax{"cfgTrackOccupancyMax", 1000000000, "max. occupancy"}; - Configurable cfgFT0COccupancyMin{"cfgFT0COccupancyMin", -2, "min. FT0C occupancy"}; - Configurable cfgFT0COccupancyMax{"cfgFT0COccupancyMax", 1000000000, "max. FT0C occupancy"}; - Configurable onlyKeepWeightedEvents{"onlyKeepWeightedEvents", false, "flag to keep only weighted events (for JJ MCs) and remove all MB events (with weight = 1)"}; + o2::framework::Configurable cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"}; + o2::framework::Configurable cfgZvtxMax{"cfgZvtxMax", +10.f, "max. Zvtx"}; + o2::framework::Configurable cfgRequireSel8{"cfgRequireSel8", true, "require sel8 in event cut"}; + o2::framework::Configurable cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND in event cut"}; + o2::framework::Configurable cfgRequireNoTFB{"cfgRequireNoTFB", false, "require No time frame border in event cut"}; + o2::framework::Configurable cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border in event cut"}; + o2::framework::Configurable cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"}; + o2::framework::Configurable cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV. + o2::framework::Configurable cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"}; + o2::framework::Configurable cfgRequireEMCReadoutInMB{"cfgRequireEMCReadoutInMB", false, "require the EMC to be read out in an MB collision (kTVXinEMC)"}; + o2::framework::Configurable cfgRequireEMCHardwareTriggered{"cfgRequireEMCHardwareTriggered", false, "require the EMC to be hardware triggered (kEMC7 or kDMC7)"}; + o2::framework::Configurable cfgTrackOccupancyMin{"cfgTrackOccupancyMin", -2, "min. occupancy"}; + o2::framework::Configurable cfgTrackOccupancyMax{"cfgTrackOccupancyMax", 1000000000, "max. occupancy"}; + o2::framework::Configurable cfgFT0COccupancyMin{"cfgFT0COccupancyMin", -2, "min. FT0C occupancy"}; + o2::framework::Configurable cfgFT0COccupancyMax{"cfgFT0COccupancyMax", 1000000000, "max. FT0C occupancy"}; + o2::framework::Configurable onlyKeepWeightedEvents{"onlyKeepWeightedEvents", false, "flag to keep only weighted events (for JJ MCs) and remove all MB events (with weight = 1)"}; } eventcuts; V0PhotonCut fV0PhotonCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "pcmcut_group"; - Configurable cfg_require_v0_with_itstpc{"cfg_require_v0_with_itstpc", false, "flag to select V0s with ITS-TPC matched tracks"}; - Configurable cfg_require_v0_with_itsonly{"cfg_require_v0_with_itsonly", false, "flag to select V0s with ITSonly tracks"}; - Configurable cfg_require_v0_with_tpconly{"cfg_require_v0_with_tpconly", false, "flag to select V0s with TPConly tracks"}; - Configurable cfg_min_pt_v0{"cfg_min_pt_v0", 0.1, "min pT for v0 photons at PV"}; - Configurable cfg_max_pt_v0{"cfg_max_pt_v0", 1e+10, "max pT for v0 photons at PV"}; - Configurable cfg_min_eta_v0{"cfg_min_eta_v0", -0.8, "min eta for v0 photons at PV"}; - Configurable cfg_max_eta_v0{"cfg_max_eta_v0", 0.8, "max eta for v0 photons at PV"}; - Configurable cfg_min_v0radius{"cfg_min_v0radius", 4.0, "min v0 radius"}; - Configurable cfg_max_v0radius{"cfg_max_v0radius", 90.0, "max v0 radius"}; - Configurable cfg_max_alpha_ap{"cfg_max_alpha_ap", 0.95, "max alpha for AP cut"}; - Configurable cfg_max_qt_ap{"cfg_max_qt_ap", 0.01, "max qT for AP cut"}; - Configurable cfg_min_cospa{"cfg_min_cospa", 0.999, "min V0 CosPA"}; - Configurable cfg_max_pca{"cfg_max_pca", 1.5, "max distance btween 2 legs"}; - Configurable cfg_max_chi2kf{"cfg_max_chi2kf", 1e+10, "max chi2/ndf with KF"}; - Configurable cfg_reject_v0_on_itsib{"cfg_reject_v0_on_itsib", true, "flag to reject V0s on ITSib"}; - Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to V0"}; - - Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; - Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 40, "min ncrossed rows"}; - Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; - Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; - Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; - Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -3.0, "min. TPC n sigma for electron"}; - Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron"}; - Configurable cfg_disable_itsonly_track{"cfg_disable_itsonly_track", false, "flag to disable ITSonly tracks"}; - Configurable cfg_disable_tpconly_track{"cfg_disable_tpconly_track", false, "flag to disable TPConly tracks"}; + o2::framework::Configurable cfg_require_v0_with_itstpc{"cfg_require_v0_with_itstpc", false, "flag to select V0s with ITS-TPC matched tracks"}; + o2::framework::Configurable cfg_require_v0_with_itsonly{"cfg_require_v0_with_itsonly", false, "flag to select V0s with ITSonly tracks"}; + o2::framework::Configurable cfg_require_v0_with_tpconly{"cfg_require_v0_with_tpconly", false, "flag to select V0s with TPConly tracks"}; + o2::framework::Configurable cfg_min_pt_v0{"cfg_min_pt_v0", 0.1, "min pT for v0 photons at PV"}; + o2::framework::Configurable cfg_max_pt_v0{"cfg_max_pt_v0", 1e+10, "max pT for v0 photons at PV"}; + o2::framework::Configurable cfg_min_eta_v0{"cfg_min_eta_v0", -0.8, "min eta for v0 photons at PV"}; + o2::framework::Configurable cfg_max_eta_v0{"cfg_max_eta_v0", 0.8, "max eta for v0 photons at PV"}; + o2::framework::Configurable cfg_min_v0radius{"cfg_min_v0radius", 4.0, "min v0 radius"}; + o2::framework::Configurable cfg_max_v0radius{"cfg_max_v0radius", 90.0, "max v0 radius"}; + o2::framework::Configurable cfg_max_alpha_ap{"cfg_max_alpha_ap", 0.95, "max alpha for AP cut"}; + o2::framework::Configurable cfg_max_qt_ap{"cfg_max_qt_ap", 0.01, "max qT for AP cut"}; + o2::framework::Configurable cfg_min_cospa{"cfg_min_cospa", 0.999, "min V0 CosPA"}; + o2::framework::Configurable cfg_max_pca{"cfg_max_pca", 1.5, "max distance btween 2 legs"}; + o2::framework::Configurable cfg_max_chi2kf{"cfg_max_chi2kf", 1e+10, "max chi2/ndf with KF"}; + o2::framework::Configurable cfg_reject_v0_on_itsib{"cfg_reject_v0_on_itsib", true, "flag to reject V0s on ITSib"}; + o2::framework::Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to V0"}; + + o2::framework::Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; + o2::framework::Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 40, "min ncrossed rows"}; + o2::framework::Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; + o2::framework::Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; + o2::framework::Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; + o2::framework::Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -3.0, "min. TPC n sigma for electron"}; + o2::framework::Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron"}; + o2::framework::Configurable cfg_disable_itsonly_track{"cfg_disable_itsonly_track", false, "flag to disable ITSonly tracks"}; + o2::framework::Configurable cfg_disable_tpconly_track{"cfg_disable_tpconly_track", false, "flag to disable TPConly tracks"}; } pcmcuts; DalitzEECut fDileptonCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "dileptoncut_group"; - Configurable cfg_min_mass{"cfg_min_mass", 0.0, "min mass"}; - Configurable cfg_max_mass{"cfg_max_mass", 0.1, "max mass"}; - Configurable cfg_apply_phiv{"cfg_apply_phiv", true, "flag to apply phiv cut"}; - Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; - Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; - Configurable cfg_phiv_slope{"cfg_phiv_slope", 0.0185, "slope for m vs. phiv"}; - Configurable cfg_phiv_intercept{"cfg_phiv_intercept", -0.0280, "intercept for m vs. phiv"}; - - Configurable cfg_min_pt_track{"cfg_min_pt_track", 0.1, "min pT for single track"}; - Configurable cfg_max_eta_track{"cfg_max_eta_track", 0.8, "max eta for single track"}; - Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; - Configurable cfg_min_ncluster_its{"cfg_min_ncluster_its", 5, "min ncluster its"}; - Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 70, "min ncrossed rows"}; - Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; - Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; - Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.05, "max dca XY for single track in cm"}; - Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.05, "max dca Z for single track in cm"}; - Configurable cfg_max_dca3dsigma_track{"cfg_max_dca3dsigma_track", 1.5, "max DCA 3D in sigma"}; - Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; - Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to electron"}; - Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; - Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; - - Configurable cfg_pid_scheme{"cfg_pid_scheme", static_cast(DalitzEECut::PIDSchemes::kTOFif), "pid scheme [kTOFif : 0, kTPConly : 1]"}; - Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -2.0, "min. TPC n sigma for electron inclusion"}; - Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron inclusion"}; - Configurable cfg_min_TPCNsigmaPi{"cfg_min_TPCNsigmaPi", -0.0, "min. TPC n sigma for pion exclusion"}; - Configurable cfg_max_TPCNsigmaPi{"cfg_max_TPCNsigmaPi", +0.0, "max. TPC n sigma for pion exclusion"}; - Configurable cfg_min_TOFNsigmaEl{"cfg_min_TOFNsigmaEl", -3.0, "min. TOF n sigma for electron inclusion"}; - Configurable cfg_max_TOFNsigmaEl{"cfg_max_TOFNsigmaEl", +3.0, "max. TOF n sigma for electron inclusion"}; + o2::framework::Configurable cfg_min_mass{"cfg_min_mass", 0.0, "min mass"}; + o2::framework::Configurable cfg_max_mass{"cfg_max_mass", 0.1, "max mass"}; + o2::framework::Configurable cfg_apply_phiv{"cfg_apply_phiv", true, "flag to apply phiv cut"}; + o2::framework::Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; + o2::framework::Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; + o2::framework::Configurable cfg_phiv_slope{"cfg_phiv_slope", 0.0185, "slope for m vs. phiv"}; + o2::framework::Configurable cfg_phiv_intercept{"cfg_phiv_intercept", -0.0280, "intercept for m vs. phiv"}; + + o2::framework::Configurable cfg_min_pt_track{"cfg_min_pt_track", 0.1, "min pT for single track"}; + o2::framework::Configurable cfg_max_eta_track{"cfg_max_eta_track", 0.8, "max eta for single track"}; + o2::framework::Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; + o2::framework::Configurable cfg_min_ncluster_its{"cfg_min_ncluster_its", 5, "min ncluster its"}; + o2::framework::Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 70, "min ncrossed rows"}; + o2::framework::Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; + o2::framework::Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; + o2::framework::Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.05, "max dca XY for single track in cm"}; + o2::framework::Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.05, "max dca Z for single track in cm"}; + o2::framework::Configurable cfg_max_dca3dsigma_track{"cfg_max_dca3dsigma_track", 1.5, "max DCA 3D in sigma"}; + o2::framework::Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; + o2::framework::Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to electron"}; + o2::framework::Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; + o2::framework::Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; + + o2::framework::Configurable cfg_pid_scheme{"cfg_pid_scheme", static_cast(DalitzEECut::PIDSchemes::kTOFif), "pid scheme [kTOFif : 0, kTPConly : 1]"}; + o2::framework::Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -2.0, "min. TPC n sigma for electron inclusion"}; + o2::framework::Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron inclusion"}; + o2::framework::Configurable cfg_min_TPCNsigmaPi{"cfg_min_TPCNsigmaPi", -0.0, "min. TPC n sigma for pion exclusion"}; + o2::framework::Configurable cfg_max_TPCNsigmaPi{"cfg_max_TPCNsigmaPi", +0.0, "max. TPC n sigma for pion exclusion"}; + o2::framework::Configurable cfg_min_TOFNsigmaEl{"cfg_min_TOFNsigmaEl", -3.0, "min. TOF n sigma for electron inclusion"}; + o2::framework::Configurable cfg_max_TOFNsigmaEl{"cfg_max_TOFNsigmaEl", +3.0, "max. TOF n sigma for electron inclusion"}; } dileptoncuts; EMCPhotonCut fEMCCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "emccut_group"; - Configurable clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"}; - Configurable minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"}; - Configurable EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"}; - Configurable EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"}; - Configurable EMC_minM02{"EMC_minM02", 0.1, "Minimum M02 for EMCal M02 cut"}; - Configurable EMC_maxM02{"EMC_maxM02", 0.7, "Maximum M02 for EMCal M02 cut"}; - Configurable EMC_minE{"EMC_minE", 0.7, "Minimum cluster energy for EMCal energy cut"}; - Configurable EMC_minNCell{"EMC_minNCell", 1, "Minimum number of cells per cluster for EMCal NCell cut"}; - Configurable> EMC_TM_Eta{"EMC_TM_Eta", {0.01f, 4.07f, -2.5f}, "|eta| <= [0]+(pT+[1])^[2] for EMCal track matching"}; - Configurable> EMC_TM_Phi{"EMC_TM_Phi", {0.015f, 3.65f, -2.f}, "|phi| <= [0]+(pT+[1])^[2] for EMCal track matching"}; - Configurable EMC_Eoverp{"EMC_Eoverp", 1.75, "Minimum cluster energy over track momentum for EMCal track matching"}; - Configurable EMC_UseExoticCut{"EMC_UseExoticCut", true, "FLag to use the EMCal exotic cluster cut"}; - Configurable cfgDistanceToEdge{"cfgDistanceToEdge", 1, "Distance to edge in cells required for rotated cluster to be accepted"}; + o2::framework::Configurable clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"}; + o2::framework::Configurable minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"}; + o2::framework::Configurable EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"}; + o2::framework::Configurable EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"}; + o2::framework::Configurable EMC_minM02{"EMC_minM02", 0.1, "Minimum M02 for EMCal M02 cut"}; + o2::framework::Configurable EMC_maxM02{"EMC_maxM02", 0.7, "Maximum M02 for EMCal M02 cut"}; + o2::framework::Configurable EMC_minE{"EMC_minE", 0.7, "Minimum cluster energy for EMCal energy cut"}; + o2::framework::Configurable EMC_minNCell{"EMC_minNCell", 1, "Minimum number of cells per cluster for EMCal NCell cut"}; + o2::framework::Configurable> EMC_TM_Eta{"EMC_TM_Eta", {0.01f, 4.07f, -2.5f}, "|eta| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + o2::framework::Configurable> EMC_TM_Phi{"EMC_TM_Phi", {0.015f, 3.65f, -2.f}, "|phi| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + o2::framework::Configurable EMC_Eoverp{"EMC_Eoverp", 1.75, "Minimum cluster energy over track momentum for EMCal track matching"}; + o2::framework::Configurable EMC_UseExoticCut{"EMC_UseExoticCut", true, "FLag to use the EMCal exotic cluster cut"}; + o2::framework::Configurable cfgDistanceToEdge{"cfgDistanceToEdge", 1, "Distance to edge in cells required for rotated cluster to be accepted"}; } emccuts; PHOSPhotonCut fPHOSCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "phoscut_group"; - Configurable cfg_min_Ecluster{"cfg_min_Ecluster", 0.3, "Minimum cluster energy for PHOS in GeV"}; + o2::framework::Configurable cfg_min_Ecluster{"cfg_min_Ecluster", 0.3, "Minimum cluster energy for PHOS in GeV"}; } phoscuts; - HistogramRegistry fRegistry{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false}; - static constexpr std::string_view event_types[2] = {"before/", "after/"}; - static constexpr std::string_view event_pair_types[2] = {"same/", "mix/"}; + o2::framework::HistogramRegistry fRegistry{"output", {}, o2::framework::OutputObjHandlingPolicy::AnalysisObject, false, false}; + // static constexpr std::string_view event_types[2] = {"before/", "after/"}; + // static constexpr std::string_view event_pair_types[2] = {"same/", "mix/"}; std::vector zvtx_bin_edges; std::vector cent_bin_edges; @@ -232,12 +213,12 @@ struct Pi0EtaToGammaGamma { std::vector occ_bin_edges; o2::ccdb::CcdbApi ccdbApi; - Service ccdb; + o2::framework::Service ccdb; int mRunNumber; float d_bz; o2::emcal::Geometry* emcalGeom; - void init(InitContext&) + void init(o2::framework::InitContext&) { zvtx_bin_edges = std::vector(ConfVtxBins.value.begin(), ConfVtxBins.value.end()); zvtx_bin_edges.erase(zvtx_bin_edges.begin()); @@ -252,11 +233,11 @@ struct Pi0EtaToGammaGamma { occ_bin_edges = std::vector(ConfOccupancyBins.value.begin(), ConfOccupancyBins.value.end()); occ_bin_edges.erase(occ_bin_edges.begin()); - emh1 = new MyEMH(ndepth); - emh2 = new MyEMH(ndepth); + emh1 = new o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>(ndepth); + emh2 = new o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>(ndepth); o2::aod::pwgem::photonmeson::utils::eventhistogram::addEventHistograms(&fRegistry); - if constexpr (pairtype == PairType::kPCMDalitzEE) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { o2::aod::pwgem::photonmeson::utils::nmhistogram::addNMHistograms(&fRegistry, false, "ee#gamma"); } else { o2::aod::pwgem::photonmeson::utils::nmhistogram::addNMHistograms(&fRegistry, false, "#gamma#gamma"); @@ -267,11 +248,11 @@ struct Pi0EtaToGammaGamma { DefineEMCCut(); DefinePHOSCut(); - if constexpr (pairtype == kEMCEMC) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { fRegistry.addClone("Pair/same/", "Pair/rotation/"); emcalGeom = o2::emcal::Geometry::GetInstanceFromRunNumber(300000); } - fRegistry.add("Pair/mix/hDiffBC", "diff. global BC in mixed event;|BC_{current} - BC_{mixed}|", kTH1D, {{10001, -0.5, 10000.5}}, true); + fRegistry.add("Pair/mix/hDiffBC", "diff. global BC in mixed event;|BC_{current} - BC_{mixed}|", o2::framework::kTH1D, {{10001, -0.5, 10000.5}}, true); mRunNumber = 0; d_bz = 0; @@ -489,7 +470,7 @@ struct Pi0EtaToGammaGamma { if (photons_coll.size() < 3) { return; } - const float rotationAngle = M_PI / 2.0; // rotaion angle 90 degree + const float rotationAngle = o2::constants::math::PIHalf; // rotaion angle 90 degree ROOT::Math::AxisAngle rotationAxis(meson.Vect(), rotationAngle); ROOT::Math::Rotation3D rotationMatrix(rotationAxis); photon1 = rotationMatrix * photon1; @@ -523,7 +504,7 @@ struct Pi0EtaToGammaGamma { // only combine rotated photons with other photons continue; } - if (!(fEMCCut.IsSelected(photon))) { + if (!(fEMCCut.IsSelected::iterator>(photon))) { continue; } @@ -544,18 +525,18 @@ struct Pi0EtaToGammaGamma { return; } - SliceCache cache; - Preslice perCollision_pcm = aod::v0photonkf::emeventId; - Preslice perCollision_emc = aod::emccluster::emeventId; - Preslice perCollision_phos = aod::phoscluster::emeventId; + o2::framework::SliceCache cache; + o2::framework::Preslice>> perCollision_pcm = o2::aod::v0photonkf::emeventId; + o2::framework::Preslice> perCollision_emc = o2::aod::emccluster::emeventId; + o2::framework::Preslice> perCollision_phos = o2::aod::phoscluster::emeventId; - Preslice perCollision_electron = aod::emprimaryelectron::emeventId; - Partition positrons = o2::aod::emprimaryelectron::sign > int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt&& nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; - Partition electrons = o2::aod::emprimaryelectron::sign < int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt && nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; + o2::framework::Preslice>> perCollision_electron = o2::aod::emprimaryelectron::emeventId; + o2::framework::Partition>> positrons = o2::aod::emprimaryelectron::sign > int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt&& nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; + o2::framework::Partition>> electrons = o2::aod::emprimaryelectron::sign < int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt && nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; - using MyEMH = o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, EMTrack>; - MyEMH* emh1 = nullptr; - MyEMH* emh2 = nullptr; + // using MyEMH = o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>; + o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>* emh1 = nullptr; + o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>* emh2 = nullptr; std::vector> used_photonIds; // std::vector> used_dileptonIds; // std::map, uint64_t> map_mixed_eventId_to_globalBC; @@ -570,12 +551,12 @@ struct Pi0EtaToGammaGamma { for (const auto& collision : collisions) { initCCDB(collision); int ndiphoton = 0; - if ((pairtype == PairType::kPHOSPHOS || pairtype == PairType::kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { + if ((pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { continue; } float weight = 1.f; - if constexpr (std::is_same_v, FilteredMyCollisionsWithJJMC>) { + if constexpr (std::is_same_v, o2::soa::Filtered, o2::aod::EMEventsWeight>>>) { weight = collision.weight(); } @@ -639,11 +620,11 @@ struct Pi0EtaToGammaGamma { std::tuple key_bin = std::make_tuple(zbin, centbin, epbin, occbin); std::pair key_df_collision = std::make_pair(ndf, collision.globalIndex()); - if constexpr (pairtype == PairType::kPCMPCM || pairtype == PairType::kPHOSPHOS || pairtype == PairType::kEMCEMC) { // same kinds pairing + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { // same kinds pairing auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto photons2_per_collision = photons2.sliceBy(perCollision2, collision.globalIndex()); - for (const auto& [g1, g2] : combinations(CombinationsStrictlyUpperIndexPolicy(photons1_per_collision, photons2_per_collision))) { + for (const auto& [g1, g2] : o2::soa::combinations(o2::soa::CombinationsStrictlyUpperIndexPolicy(photons1_per_collision, photons2_per_collision))) { if (!cut1.template IsSelected(g1) || !cut2.template IsSelected(g2)) { continue; } @@ -655,7 +636,7 @@ struct Pi0EtaToGammaGamma { continue; } - if (pairtype == PairType::kEMCEMC) { + if (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); if (openingAngle < emccuts.minOpenAngle) { continue; @@ -664,24 +645,24 @@ struct Pi0EtaToGammaGamma { fRegistry.fill(HIST("Pair/same/hs"), v12.M(), v12.Pt(), weight); - if constexpr (pairtype == PairType::kEMCEMC) { - RotationBackground(v12, v1, v2, photons2_per_collision, g1.globalIndex(), g2.globalIndex(), weight); + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { + RotationBackground>(v12, v1, v2, photons2_per_collision, g1.globalIndex(), g2.globalIndex(), weight); } std::pair pair_tmp_id1 = std::make_pair(ndf, g1.globalIndex()); std::pair pair_tmp_id2 = std::make_pair(ndf, g2.globalIndex()); if (std::find(used_photonIds.begin(), used_photonIds.end(), pair_tmp_id1) == used_photonIds.end()) { - emh1->AddTrackToEventPool(key_df_collision, EMTrack(-1, g1.globalIndex(), collision.globalIndex(), g1.globalIndex(), g1.pt(), g1.eta(), g1.phi(), 0)); + emh1->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, g1.globalIndex(), collision.globalIndex(), g1.globalIndex(), g1.pt(), g1.eta(), g1.phi(), 0)); used_photonIds.emplace_back(pair_tmp_id1); } if (std::find(used_photonIds.begin(), used_photonIds.end(), pair_tmp_id2) == used_photonIds.end()) { - emh1->AddTrackToEventPool(key_df_collision, EMTrack(-1, g2.globalIndex(), collision.globalIndex(), g2.globalIndex(), g2.pt(), g2.eta(), g2.phi(), 0)); + emh1->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, g2.globalIndex(), collision.globalIndex(), g2.globalIndex(), g2.pt(), g2.eta(), g2.phi(), 0)); used_photonIds.emplace_back(pair_tmp_id2); } ndiphoton++; } // end of pairing loop - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto positrons_per_collision = positrons->sliceByCached(o2::aod::emprimaryelectron::emeventId, collision.globalIndex(), cache); auto electrons_per_collision = electrons->sliceByCached(o2::aod::emprimaryelectron::emeventId, collision.globalIndex(), cache); @@ -694,7 +675,7 @@ struct Pi0EtaToGammaGamma { auto ele1 = g1.template negTrack_as(); ROOT::Math::PtEtaPhiMVector v_gamma(g1.pt(), g1.eta(), g1.phi(), 0.); - for (const auto& [pos2, ele2] : combinations(CombinationsFullIndexPolicy(positrons_per_collision, electrons_per_collision))) { + for (const auto& [pos2, ele2] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(positrons_per_collision, electrons_per_collision))) { if (pos2.trackId() == ele2.trackId()) { // this is protection against pairing identical 2 tracks. continue; @@ -724,11 +705,11 @@ struct Pi0EtaToGammaGamma { std::pair pair_tmp_id1 = std::make_pair(ndf, g1.globalIndex()); std::tuple tuple_tmp_id2 = std::make_tuple(ndf, collision.globalIndex(), pos2.trackId(), ele2.trackId()); if (std::find(used_photonIds.begin(), used_photonIds.end(), pair_tmp_id1) == used_photonIds.end()) { - emh1->AddTrackToEventPool(key_df_collision, EMTrack(-1, g1.globalIndex(), collision.globalIndex(), -1, g1.pt(), g1.eta(), g1.phi(), 0)); + emh1->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, g1.globalIndex(), collision.globalIndex(), -1, g1.pt(), g1.eta(), g1.phi(), 0)); used_photonIds.emplace_back(pair_tmp_id1); } if (std::find(used_dileptonIds.begin(), used_dileptonIds.end(), tuple_tmp_id2) == used_dileptonIds.end()) { - emh2->AddTrackToEventPool(key_df_collision, EMTrack(-1, -1, collision.globalIndex(), -1, v_ee.Pt(), v_ee.Eta(), v_ee.Phi(), v_ee.M())); + emh2->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, -1, collision.globalIndex(), -1, v_ee.Pt(), v_ee.Eta(), v_ee.Phi(), v_ee.M())); used_dileptonIds.emplace_back(tuple_tmp_id2); } ndiphoton++; @@ -738,7 +719,7 @@ struct Pi0EtaToGammaGamma { auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto photons2_per_collision = photons2.sliceBy(perCollision2, collision.globalIndex()); - for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photons1_per_collision, photons2_per_collision))) { + for (const auto& [g1, g2] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(photons1_per_collision, photons2_per_collision))) { if (!cut1.template IsSelected(g1) || !cut2.template IsSelected(g2)) { continue; } @@ -755,11 +736,11 @@ struct Pi0EtaToGammaGamma { std::pair pair_tmp_id2 = std::make_pair(ndf, g2.globalIndex()); if (std::find(used_photonIds.begin(), used_photonIds.end(), pair_tmp_id1) == used_photonIds.end()) { - emh1->AddTrackToEventPool(key_df_collision, EMTrack(-1, g1.globalIndex(), collision.globalIndex(), -1, g1.pt(), g1.eta(), g1.phi(), 0)); + emh1->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, g1.globalIndex(), collision.globalIndex(), -1, g1.pt(), g1.eta(), g1.phi(), 0)); used_photonIds.emplace_back(pair_tmp_id1); } if (std::find(used_photonIds.begin(), used_photonIds.end(), pair_tmp_id2) == used_photonIds.end()) { - emh2->AddTrackToEventPool(key_df_collision, EMTrack(-1, g2.globalIndex(), collision.globalIndex(), -1, g2.pt(), g2.eta(), g2.phi(), 0)); + emh2->AddTrackToEventPool(key_df_collision, o2::aod::pwgem::dilepton::utils::EMTrack(-1, g2.globalIndex(), collision.globalIndex(), -1, g2.pt(), g2.eta(), g2.phi(), 0)); used_photonIds.emplace_back(pair_tmp_id2); } ndiphoton++; @@ -778,7 +759,7 @@ struct Pi0EtaToGammaGamma { auto collisionIds1_in_mixing_pool = emh1->GetCollisionIdsFromEventPool(key_bin); auto collisionIds2_in_mixing_pool = emh2->GetCollisionIdsFromEventPool(key_bin); - if constexpr (pairtype == PairType::kPCMPCM || pairtype == PairType::kPHOSPHOS || pairtype == PairType::kEMCEMC) { // same kinds pairing + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { // same kinds pairing for (const auto& mix_dfId_collisionId : collisionIds1_in_mixing_pool) { int mix_dfId = mix_dfId_collisionId.first; int64_t mix_collisionId = mix_dfId_collisionId.second; @@ -834,7 +815,7 @@ struct Pi0EtaToGammaGamma { for (const auto& g2 : photons2_from_event_pool) { ROOT::Math::PtEtaPhiMVector v1(g1.pt(), g1.eta(), g1.phi(), 0.); ROOT::Math::PtEtaPhiMVector v2(g2.pt(), g2.eta(), g2.phi(), 0.); - if constexpr (pairtype == PairType::kPCMDalitzEE) { //[photon from event1, dilepton from event2] and [photon from event2, dilepton from event1] + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { //[photon from event1, dilepton from event2] and [photon from event2, dilepton from event1] v2.SetM(g2.mass()); } ROOT::Math::PtEtaPhiMVector v12 = v1 + v2; @@ -867,7 +848,7 @@ struct Pi0EtaToGammaGamma { for (const auto& g2 : photons1_from_event_pool) { ROOT::Math::PtEtaPhiMVector v1(g1.pt(), g1.eta(), g1.phi(), 0.); ROOT::Math::PtEtaPhiMVector v2(g2.pt(), g2.eta(), g2.phi(), 0.); - if constexpr (pairtype == PairType::kPCMDalitzEE) { //[photon from event1, dilepton from event2] and [photon from event2, dilepton from event1] + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { //[photon from event1, dilepton from event2] and [photon from event2, dilepton from event1] v1.SetM(g1.mass()); } ROOT::Math::PtEtaPhiMVector v12 = v1 + v2; @@ -889,32 +870,32 @@ struct Pi0EtaToGammaGamma { } // end of collision loop } - Filter collisionFilter_occupancy_track = eventcuts.cfgTrackOccupancyMin <= o2::aod::evsel::trackOccupancyInTimeRange && o2::aod::evsel::trackOccupancyInTimeRange < eventcuts.cfgTrackOccupancyMax; - Filter collisionFilter_occupancy_ft0c = eventcuts.cfgFT0COccupancyMin <= o2::aod::evsel::ft0cOccupancyInTimeRange && o2::aod::evsel::ft0cOccupancyInTimeRange < eventcuts.cfgFT0COccupancyMax; - Filter collisionFilter_centrality = (cfgCentMin < o2::aod::cent::centFT0M && o2::aod::cent::centFT0M < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0A && o2::aod::cent::centFT0A < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0C && o2::aod::cent::centFT0C < cfgCentMax); - using FilteredMyCollisions = soa::Filtered; + o2::framework::expressions::Filter collisionFilter_occupancy_track = eventcuts.cfgTrackOccupancyMin <= o2::aod::evsel::trackOccupancyInTimeRange && o2::aod::evsel::trackOccupancyInTimeRange < eventcuts.cfgTrackOccupancyMax; + o2::framework::expressions::Filter collisionFilter_occupancy_ft0c = eventcuts.cfgFT0COccupancyMin <= o2::aod::evsel::ft0cOccupancyInTimeRange && o2::aod::evsel::ft0cOccupancyInTimeRange < eventcuts.cfgFT0COccupancyMax; + o2::framework::expressions::Filter collisionFilter_centrality = (cfgCentMin < o2::aod::cent::centFT0M && o2::aod::cent::centFT0M < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0A && o2::aod::cent::centFT0A < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0C && o2::aod::cent::centFT0C < cfgCentMax); + // using FilteredMyCollisions = o2::soa::Filtered>; - Filter prefilter_pcm = ifnode(pcmcuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::v0photonkf::pfbderived == static_cast(0), true); - Filter prefilter_primaryelectron = ifnode(dileptoncuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::emprimaryelectron::pfbderived == static_cast(0), true); + o2::framework::expressions::Filter prefilter_pcm = ifnode(pcmcuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::v0photonkf::pfbderived == static_cast(0), true); + o2::framework::expressions::Filter prefilter_primaryelectron = ifnode(dileptoncuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::emprimaryelectron::pfbderived == static_cast(0), true); int ndf = 0; - void processAnalysis(FilteredMyCollisions const& collisions, Types const&... args) + void processAnalysis(o2::soa::Filtered> const& collisions, Types const&... args) { // LOGF(info, "ndf = %d", ndf); - if constexpr (pairtype == PairType::kPCMPCM) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); runPairing(collisions, v0photons, v0photons, v0legs, v0legs, perCollision_pcm, perCollision_pcm, fV0PhotonCut, fV0PhotonCut); - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); auto emprimaryelectrons = std::get<2>(std::tie(args...)); // LOGF(info, "electrons.size() = %d, positrons.size() = %d", electrons.size(), positrons.size()); runPairing(collisions, v0photons, emprimaryelectrons, v0legs, emprimaryelectrons, perCollision_pcm, perCollision_electron, fV0PhotonCut, fDileptonCut); - } else if constexpr (pairtype == PairType::kEMCEMC) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { auto emcclusters = std::get<0>(std::tie(args...)); runPairing(collisions, emcclusters, emcclusters, nullptr, nullptr, perCollision_emc, perCollision_emc, fEMCCut, fEMCCut); - } else if constexpr (pairtype == PairType::kPHOSPHOS) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS) { auto phosclusters = std::get<0>(std::tie(args...)); runPairing(collisions, phosclusters, phosclusters, nullptr, nullptr, perCollision_phos, perCollision_phos, fPHOSCut, fPHOSCut); } @@ -934,24 +915,24 @@ struct Pi0EtaToGammaGamma { } PROCESS_SWITCH(Pi0EtaToGammaGamma, processAnalysis, "process pair analysis", true); - using FilteredMyCollisionsWithJJMC = soa::Filtered; - void processAnalysisJJMC(FilteredMyCollisionsWithJJMC const& collisions, Types const&... args) + // using FilteredMyCollisionsWithJJMC = o2::soa::Filtered, o2::aod::EMEventsWeight>>; + void processAnalysisJJMC(o2::soa::Filtered, o2::aod::EMEventsWeight>> const& collisions, Types const&... args) { // LOGF(info, "ndf = %d", ndf); - if constexpr (pairtype == PairType::kPCMPCM) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); runPairing(collisions, v0photons, v0photons, v0legs, v0legs, perCollision_pcm, perCollision_pcm, fV0PhotonCut, fV0PhotonCut); - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); auto emprimaryelectrons = std::get<2>(std::tie(args...)); // LOGF(info, "electrons.size() = %d, positrons.size() = %d", electrons.size(), positrons.size()); runPairing(collisions, v0photons, emprimaryelectrons, v0legs, emprimaryelectrons, perCollision_pcm, perCollision_electron, fV0PhotonCut, fDileptonCut); - } else if constexpr (pairtype == PairType::kEMCEMC) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { auto emcclusters = std::get<0>(std::tie(args...)); runPairing(collisions, emcclusters, emcclusters, nullptr, nullptr, perCollision_emc, perCollision_emc, fEMCCut, fEMCCut); - } else if constexpr (pairtype == PairType::kPHOSPHOS) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS) { auto phosclusters = std::get<0>(std::tie(args...)); runPairing(collisions, phosclusters, phosclusters, nullptr, nullptr, perCollision_phos, perCollision_phos, fPHOSCut, fPHOSCut); } @@ -971,7 +952,7 @@ struct Pi0EtaToGammaGamma { } PROCESS_SWITCH(Pi0EtaToGammaGamma, processAnalysisJJMC, "process pair analysis", false); - void processDummy(MyCollisions const&) {} + void processDummy(o2::aod::EMEvents const&) {} PROCESS_SWITCH(Pi0EtaToGammaGamma, processDummy, "Dummy function", false); }; #endif // PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMA_H_ diff --git a/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h b/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h index 16cd3054172..674e663e03c 100644 --- a/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h +++ b/PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h @@ -9,15 +9,13 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch +/// \file Pi0EtaToGammaGammaMC.h +/// \brief This code loops over photons and makes pairs for neutral mesons analyses with MC true info. +/// \author D. Sekihata, daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMAMC_H_ #define PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMAMC_H_ -#include "PWGEM/Dilepton/Utils/MCUtilities.h" #include "PWGEM/PhotonMeson/Core/DalitzEECut.h" #include "PWGEM/PhotonMeson/Core/EMCPhotonCut.h" #include "PWGEM/PhotonMeson/Core/EMPhotonEventCut.h" @@ -28,204 +26,187 @@ #include "PWGEM/PhotonMeson/Utils/MCUtilities.h" #include "PWGEM/PhotonMeson/Utils/NMHistograms.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" +// Dilepton headers +#include "PWGEM/Dilepton/Utils/MCUtilities.h" -#include "Common/Core/RecoDecay.h" - -#include "CCDB/BasicCCDBManager.h" -#include "DataFormatsParameters/GRPMagField.h" -#include "DataFormatsParameters/GRPObject.h" -#include "DetectorsBase/GeometryManager.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" - -#include "Math/Vector4D.h" -#include "TF1.h" -#include "TString.h" - -#include +#include "Common/CCDB/TriggerAliases.h" +#include "Common/DataModel/Centrality.h" +#include "Common/DataModel/EventSelection.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include // IWYU pragma: keep +#include +#include +#include +#include + +#include +#include #include +#include +#include #include -using namespace o2; -using namespace o2::aod; -using namespace o2::framework; -using namespace o2::framework::expressions; -using namespace o2::soa; -using namespace o2::aod::pwgem::photonmeson::photonpair; -using namespace o2::aod::pwgem::photonmeson::utils::mcutil; -using namespace o2::aod::pwgem::dilepton::utils::mcutil; - -using MyCollisions = soa::Join; -using MyCollision = MyCollisions::iterator; - -using MyCollisionsWithJJMC = soa::Join; -using MyCollisionWithJJMC = MyCollisionsWithJJMC::iterator; - -using MyMCCollisions = soa::Join; -using MyMCCollision = MyMCCollisions::iterator; - -using MyV0Photons = soa::Filtered>; -using MyV0Photon = MyV0Photons::iterator; - -using MyEMCClusters = soa::Join; -using MyEMCCluster = MyEMCClusters::iterator; - -using MyPHOSClusters = soa::Join; -using MyPHOSCluster = MyEMCClusters::iterator; - -using MyMCV0Legs = soa::Join; -using MyMCV0Leg = MyMCV0Legs::iterator; - -using MyMCElectrons = soa::Filtered>; -using MyMCElectron = MyMCElectrons::iterator; - -template +template struct Pi0EtaToGammaGammaMC { - Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; - Configurable grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"}; - Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; - Configurable skipGRPOquery{"skipGRPOquery", true, "skip grpo query"}; - Configurable d_bz_input{"d_bz_input", -999, "bz field in kG, -999 is automatic"}; - - Configurable cfgQvecEstimator{"cfgQvecEstimator", 0, "FT0M:0, FT0A:1, FT0C:2"}; - Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; - Configurable cfgCentMin{"cfgCentMin", -1, "min. centrality"}; - Configurable cfgCentMax{"cfgCentMax", 999, "max. centrality"}; - Configurable maxY_rec{"maxY_rec", 0.9, "maximum rapidity for reconstructed particles"}; - Configurable fd_k0s_to_pi0{"fd_k0s_pi0", "1.0", "feed down correction to pi0"}; - Configurable cfgRequireTrueAssociation{"cfgRequireTrueAssociation", false, "flag to require true mc collision association"}; + o2::framework::Configurable ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; + o2::framework::Configurable grpPath{"grpPath", "GLO/GRP/GRP", "Path of the grp file"}; + o2::framework::Configurable grpmagPath{"grpmagPath", "GLO/Config/GRPMagField", "CCDB path of the GRPMagField object"}; + o2::framework::Configurable skipGRPOquery{"skipGRPOquery", true, "skip grpo query"}; + o2::framework::Configurable d_bz_input{"d_bz_input", -999, "bz field in kG, -999 is automatic"}; + + o2::framework::Configurable cfgQvecEstimator{"cfgQvecEstimator", 0, "FT0M:0, FT0A:1, FT0C:2"}; + o2::framework::Configurable cfgCentEstimator{"cfgCentEstimator", 2, "FT0M:0, FT0A:1, FT0C:2"}; + o2::framework::Configurable cfgCentMin{"cfgCentMin", -1, "min. centrality"}; + o2::framework::Configurable cfgCentMax{"cfgCentMax", 999, "max. centrality"}; + o2::framework::Configurable maxY_rec{"maxY_rec", 0.9, "maximum rapidity for reconstructed particles"}; + o2::framework::Configurable fd_k0s_to_pi0{"fd_k0s_pi0", "1.0", "feed down correction to pi0"}; + o2::framework::Configurable cfgRequireTrueAssociation{"cfgRequireTrueAssociation", false, "flag to require true mc collision association"}; EMPhotonEventCut fEMEventCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "eventcut_group"; - Configurable cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"}; - Configurable cfgZvtxMax{"cfgZvtxMax", +10.f, "max. Zvtx"}; - Configurable cfgRequireSel8{"cfgRequireSel8", false, "require sel8 in event cut"}; - Configurable cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND in event cut"}; - Configurable cfgRequireNoTFB{"cfgRequireNoTFB", true, "require No time frame border in event cut"}; - Configurable cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border in event cut"}; - Configurable cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"}; - Configurable cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV. - Configurable cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"}; - Configurable cfgRequireEMCReadoutInMB{"cfgRequireEMCReadoutInMB", false, "require the EMC to be read out in an MB collision (kTVXinEMC)"}; - Configurable cfgRequireEMCHardwareTriggered{"cfgRequireEMCHardwareTriggered", false, "require the EMC to be hardware triggered (kEMC7 or kDMC7)"}; - Configurable cfgTrackOccupancyMin{"cfgTrackOccupancyMin", -2, "min. occupancy"}; - Configurable cfgTrackOccupancyMax{"cfgTrackOccupancyMax", 1000000000, "max. occupancy"}; - Configurable cfgFT0COccupancyMin{"cfgFT0COccupancyMin", -2, "min. FT0C occupancy"}; - Configurable cfgFT0COccupancyMax{"cfgFT0COccupancyMax", 1000000000, "max. FT0C occupancy"}; - Configurable onlyKeepWeightedEvents{"onlyKeepWeightedEvents", false, "flag to keep only weighted events (for JJ MCs) and remove all MB events (with weight = 1)"}; + o2::framework::Configurable cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"}; + o2::framework::Configurable cfgZvtxMax{"cfgZvtxMax", +10.f, "max. Zvtx"}; + o2::framework::Configurable cfgRequireSel8{"cfgRequireSel8", false, "require sel8 in event cut"}; + o2::framework::Configurable cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND in event cut"}; + o2::framework::Configurable cfgRequireNoTFB{"cfgRequireNoTFB", true, "require No time frame border in event cut"}; + o2::framework::Configurable cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border in event cut"}; + o2::framework::Configurable cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"}; + o2::framework::Configurable cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV. + o2::framework::Configurable cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"}; + o2::framework::Configurable cfgRequireEMCReadoutInMB{"cfgRequireEMCReadoutInMB", false, "require the EMC to be read out in an MB collision (kTVXinEMC)"}; + o2::framework::Configurable cfgRequireEMCHardwareTriggered{"cfgRequireEMCHardwareTriggered", false, "require the EMC to be hardware triggered (kEMC7 or kDMC7)"}; + o2::framework::Configurable cfgTrackOccupancyMin{"cfgTrackOccupancyMin", -2, "min. occupancy"}; + o2::framework::Configurable cfgTrackOccupancyMax{"cfgTrackOccupancyMax", 1000000000, "max. occupancy"}; + o2::framework::Configurable cfgFT0COccupancyMin{"cfgFT0COccupancyMin", -2, "min. FT0C occupancy"}; + o2::framework::Configurable cfgFT0COccupancyMax{"cfgFT0COccupancyMax", 1000000000, "max. FT0C occupancy"}; + o2::framework::Configurable onlyKeepWeightedEvents{"onlyKeepWeightedEvents", false, "flag to keep only weighted events (for JJ MCs) and remove all MB events (with weight = 1)"}; } eventcuts; V0PhotonCut fV0PhotonCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "pcmcut_group"; - Configurable cfg_require_v0_with_itstpc{"cfg_require_v0_with_itstpc", false, "flag to select V0s with ITS-TPC matched tracks"}; - Configurable cfg_require_v0_with_itsonly{"cfg_require_v0_with_itsonly", false, "flag to select V0s with ITSonly tracks"}; - Configurable cfg_require_v0_with_tpconly{"cfg_require_v0_with_tpconly", false, "flag to select V0s with TPConly tracks"}; - Configurable cfg_min_pt_v0{"cfg_min_pt_v0", 0.1, "min pT for v0 photons at PV"}; - Configurable cfg_max_pt_v0{"cfg_max_pt_v0", 1e+10, "max pT for v0 photons at PV"}; - Configurable cfg_min_eta_v0{"cfg_min_eta_v0", -0.8, "min eta for v0 photons at PV"}; - Configurable cfg_max_eta_v0{"cfg_max_eta_v0", 0.8, "max eta for v0 photons at PV"}; - Configurable cfg_min_v0radius{"cfg_min_v0radius", 4.0, "min v0 radius"}; - Configurable cfg_max_v0radius{"cfg_max_v0radius", 90.0, "max v0 radius"}; - Configurable cfg_max_alpha_ap{"cfg_max_alpha_ap", 0.95, "max alpha for AP cut"}; - Configurable cfg_max_qt_ap{"cfg_max_qt_ap", 0.01, "max qT for AP cut"}; - Configurable cfg_min_cospa{"cfg_min_cospa", 0.999, "min V0 CosPA"}; - Configurable cfg_max_pca{"cfg_max_pca", 1.5, "max distance btween 2 legs"}; - Configurable cfg_max_chi2kf{"cfg_max_chi2kf", 1e+10, "max chi2/ndf with KF"}; - Configurable cfg_reject_v0_on_itsib{"cfg_reject_v0_on_itsib", true, "flag to reject V0s on ITSib"}; - Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to V0"}; - - Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; - Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 40, "min ncrossed rows"}; - Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; - Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; - Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; - Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -3.0, "min. TPC n sigma for electron"}; - Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron"}; - Configurable cfg_disable_itsonly_track{"cfg_disable_itsonly_track", false, "flag to disable ITSonly tracks"}; - Configurable cfg_disable_tpconly_track{"cfg_disable_tpconly_track", false, "flag to disable TPConly tracks"}; + o2::framework::Configurable cfg_require_v0_with_itstpc{"cfg_require_v0_with_itstpc", false, "flag to select V0s with ITS-TPC matched tracks"}; + o2::framework::Configurable cfg_require_v0_with_itsonly{"cfg_require_v0_with_itsonly", false, "flag to select V0s with ITSonly tracks"}; + o2::framework::Configurable cfg_require_v0_with_tpconly{"cfg_require_v0_with_tpconly", false, "flag to select V0s with TPConly tracks"}; + o2::framework::Configurable cfg_min_pt_v0{"cfg_min_pt_v0", 0.1, "min pT for v0 photons at PV"}; + o2::framework::Configurable cfg_max_pt_v0{"cfg_max_pt_v0", 1e+10, "max pT for v0 photons at PV"}; + o2::framework::Configurable cfg_min_eta_v0{"cfg_min_eta_v0", -0.8, "min eta for v0 photons at PV"}; + o2::framework::Configurable cfg_max_eta_v0{"cfg_max_eta_v0", 0.8, "max eta for v0 photons at PV"}; + o2::framework::Configurable cfg_min_v0radius{"cfg_min_v0radius", 4.0, "min v0 radius"}; + o2::framework::Configurable cfg_max_v0radius{"cfg_max_v0radius", 90.0, "max v0 radius"}; + o2::framework::Configurable cfg_max_alpha_ap{"cfg_max_alpha_ap", 0.95, "max alpha for AP cut"}; + o2::framework::Configurable cfg_max_qt_ap{"cfg_max_qt_ap", 0.01, "max qT for AP cut"}; + o2::framework::Configurable cfg_min_cospa{"cfg_min_cospa", 0.999, "min V0 CosPA"}; + o2::framework::Configurable cfg_max_pca{"cfg_max_pca", 1.5, "max distance btween 2 legs"}; + o2::framework::Configurable cfg_max_chi2kf{"cfg_max_chi2kf", 1e+10, "max chi2/ndf with KF"}; + o2::framework::Configurable cfg_reject_v0_on_itsib{"cfg_reject_v0_on_itsib", true, "flag to reject V0s on ITSib"}; + o2::framework::Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to V0"}; + + o2::framework::Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; + o2::framework::Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 40, "min ncrossed rows"}; + o2::framework::Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; + o2::framework::Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; + o2::framework::Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; + o2::framework::Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -3.0, "min. TPC n sigma for electron"}; + o2::framework::Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron"}; + o2::framework::Configurable cfg_disable_itsonly_track{"cfg_disable_itsonly_track", false, "flag to disable ITSonly tracks"}; + o2::framework::Configurable cfg_disable_tpconly_track{"cfg_disable_tpconly_track", false, "flag to disable TPConly tracks"}; } pcmcuts; DalitzEECut fDileptonCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "dileptoncut_group"; - Configurable cfg_min_mass{"cfg_min_mass", 0.0, "min mass"}; - Configurable cfg_max_mass{"cfg_max_mass", 0.1, "max mass"}; - Configurable cfg_apply_phiv{"cfg_apply_phiv", true, "flag to apply phiv cut"}; - Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; - Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; - Configurable cfg_phiv_slope{"cfg_phiv_slope", 0.0185, "slope for m vs. phiv"}; - Configurable cfg_phiv_intercept{"cfg_phiv_intercept", -0.0280, "intercept for m vs. phiv"}; - - Configurable cfg_min_pt_track{"cfg_min_pt_track", 0.1, "min pT for single track"}; - Configurable cfg_max_eta_track{"cfg_max_eta_track", 0.8, "max eta for single track"}; - Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; - Configurable cfg_min_ncluster_its{"cfg_min_ncluster_its", 5, "min ncluster its"}; - Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 70, "min ncrossed rows"}; - Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; - Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; - Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; - Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.05, "max dca XY for single track in cm"}; - Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.05, "max dca Z for single track in cm"}; - Configurable cfg_max_dca3dsigma_track{"cfg_max_dca3dsigma_track", 1.5, "max DCA 3D in sigma"}; - Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to electron"}; - Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; - Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; - - Configurable cfg_pid_scheme{"cfg_pid_scheme", static_cast(DalitzEECut::PIDSchemes::kTOFif), "pid scheme [kTOFif : 0, kTPConly : 1]"}; - Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -2.0, "min. TPC n sigma for electron inclusion"}; - Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron inclusion"}; - Configurable cfg_min_TPCNsigmaPi{"cfg_min_TPCNsigmaPi", -0.0, "min. TPC n sigma for pion exclusion"}; - Configurable cfg_max_TPCNsigmaPi{"cfg_max_TPCNsigmaPi", +0.0, "max. TPC n sigma for pion exclusion"}; - Configurable cfg_min_TOFNsigmaEl{"cfg_min_TOFNsigmaEl", -3.0, "min. TOF n sigma for electron inclusion"}; - Configurable cfg_max_TOFNsigmaEl{"cfg_max_TOFNsigmaEl", +3.0, "max. TOF n sigma for electron inclusion"}; + o2::framework::Configurable cfg_min_mass{"cfg_min_mass", 0.0, "min mass"}; + o2::framework::Configurable cfg_max_mass{"cfg_max_mass", 0.1, "max mass"}; + o2::framework::Configurable cfg_apply_phiv{"cfg_apply_phiv", true, "flag to apply phiv cut"}; + o2::framework::Configurable cfg_require_itsib_any{"cfg_require_itsib_any", false, "flag to require ITS ib any hits"}; + o2::framework::Configurable cfg_require_itsib_1st{"cfg_require_itsib_1st", true, "flag to require ITS ib 1st hit"}; + o2::framework::Configurable cfg_phiv_slope{"cfg_phiv_slope", 0.0185, "slope for m vs. phiv"}; + o2::framework::Configurable cfg_phiv_intercept{"cfg_phiv_intercept", -0.0280, "intercept for m vs. phiv"}; + + o2::framework::Configurable cfg_min_pt_track{"cfg_min_pt_track", 0.1, "min pT for single track"}; + o2::framework::Configurable cfg_max_eta_track{"cfg_max_eta_track", 0.8, "max eta for single track"}; + o2::framework::Configurable cfg_min_ncluster_tpc{"cfg_min_ncluster_tpc", 0, "min ncluster tpc"}; + o2::framework::Configurable cfg_min_ncluster_its{"cfg_min_ncluster_its", 5, "min ncluster its"}; + o2::framework::Configurable cfg_min_ncrossedrows{"cfg_min_ncrossedrows", 70, "min ncrossed rows"}; + o2::framework::Configurable cfg_max_frac_shared_clusters_tpc{"cfg_max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"}; + o2::framework::Configurable cfg_max_chi2tpc{"cfg_max_chi2tpc", 4.0, "max chi2/NclsTPC"}; + o2::framework::Configurable cfg_max_chi2its{"cfg_max_chi2its", 36.0, "max chi2/NclsITS"}; + o2::framework::Configurable cfg_max_dcaxy{"cfg_max_dcaxy", 0.05, "max dca XY for single track in cm"}; + o2::framework::Configurable cfg_max_dcaz{"cfg_max_dcaz", 0.05, "max dca Z for single track in cm"}; + o2::framework::Configurable cfg_max_dca3dsigma_track{"cfg_max_dca3dsigma_track", 1.5, "max DCA 3D in sigma"}; + o2::framework::Configurable cfg_apply_cuts_from_prefilter_derived{"cfg_apply_cuts_from_prefilter_derived", false, "flag to apply prefilter to electron"}; + o2::framework::Configurable includeITSsa{"includeITSsa", false, "Flag to enable ITSsa tracks"}; + o2::framework::Configurable cfg_max_pt_track_ITSsa{"cfg_max_pt_track_ITSsa", 0.15, "max pt for ITSsa tracks"}; + + o2::framework::Configurable cfg_pid_scheme{"cfg_pid_scheme", static_cast(DalitzEECut::PIDSchemes::kTOFif), "pid scheme [kTOFif : 0, kTPConly : 1]"}; + o2::framework::Configurable cfg_min_TPCNsigmaEl{"cfg_min_TPCNsigmaEl", -2.0, "min. TPC n sigma for electron inclusion"}; + o2::framework::Configurable cfg_max_TPCNsigmaEl{"cfg_max_TPCNsigmaEl", +3.0, "max. TPC n sigma for electron inclusion"}; + o2::framework::Configurable cfg_min_TPCNsigmaPi{"cfg_min_TPCNsigmaPi", -0.0, "min. TPC n sigma for pion exclusion"}; + o2::framework::Configurable cfg_max_TPCNsigmaPi{"cfg_max_TPCNsigmaPi", +0.0, "max. TPC n sigma for pion exclusion"}; + o2::framework::Configurable cfg_min_TOFNsigmaEl{"cfg_min_TOFNsigmaEl", -3.0, "min. TOF n sigma for electron inclusion"}; + o2::framework::Configurable cfg_max_TOFNsigmaEl{"cfg_max_TOFNsigmaEl", +3.0, "max. TOF n sigma for electron inclusion"}; } dileptoncuts; EMCPhotonCut fEMCCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "emccut_group"; - Configurable clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"}; - Configurable minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"}; - Configurable EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"}; - Configurable EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"}; - Configurable EMC_minM02{"EMC_minM02", 0.1, "Minimum M02 for EMCal M02 cut"}; - Configurable EMC_maxM02{"EMC_maxM02", 0.7, "Maximum M02 for EMCal M02 cut"}; - Configurable EMC_minE{"EMC_minE", 0.7, "Minimum cluster energy for EMCal energy cut"}; - Configurable EMC_minNCell{"EMC_minNCell", 1, "Minimum number of cells per cluster for EMCal NCell cut"}; - Configurable> EMC_TM_Eta{"EMC_TM_Eta", {0.01f, 4.07f, -2.5f}, "|eta| <= [0]+(pT+[1])^[2] for EMCal track matching"}; - Configurable> EMC_TM_Phi{"EMC_TM_Phi", {0.015f, 3.65f, -2.f}, "|phi| <= [0]+(pT+[1])^[2] for EMCal track matching"}; - Configurable EMC_Eoverp{"EMC_Eoverp", 1.75, "Minimum cluster energy over track momentum for EMCal track matching"}; - Configurable EMC_UseExoticCut{"EMC_UseExoticCut", true, "FLag to use the EMCal exotic cluster cut"}; + o2::framework::Configurable clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"}; + o2::framework::Configurable minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle"}; + o2::framework::Configurable EMC_minTime{"EMC_minTime", -20., "Minimum cluster time for EMCal time cut"}; + o2::framework::Configurable EMC_maxTime{"EMC_maxTime", +25., "Maximum cluster time for EMCal time cut"}; + o2::framework::Configurable EMC_minM02{"EMC_minM02", 0.1, "Minimum M02 for EMCal M02 cut"}; + o2::framework::Configurable EMC_maxM02{"EMC_maxM02", 0.7, "Maximum M02 for EMCal M02 cut"}; + o2::framework::Configurable EMC_minE{"EMC_minE", 0.7, "Minimum cluster energy for EMCal energy cut"}; + o2::framework::Configurable EMC_minNCell{"EMC_minNCell", 1, "Minimum number of cells per cluster for EMCal NCell cut"}; + o2::framework::Configurable> EMC_TM_Eta{"EMC_TM_Eta", {0.01f, 4.07f, -2.5f}, "|eta| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + o2::framework::Configurable> EMC_TM_Phi{"EMC_TM_Phi", {0.015f, 3.65f, -2.f}, "|phi| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + o2::framework::Configurable EMC_Eoverp{"EMC_Eoverp", 1.75, "Minimum cluster energy over track momentum for EMCal track matching"}; + o2::framework::Configurable EMC_UseExoticCut{"EMC_UseExoticCut", true, "FLag to use the EMCal exotic cluster cut"}; } emccuts; - Configurable maxY_gen{"maxY_gen", 0.9, "maximum rapidity for generated particles"}; // for PCM and dielectron - Configurable maxRgen{"maxRgen", 90.f, "maximum radius for generated particles"}; - Configurable margin_z_mc{"margin_z_mc", 7.0, "margin for z cut in cm for MC"}; + o2::framework::Configurable maxY_gen{"maxY_gen", 0.9, "maximum rapidity for generated particles"}; // for PCM and dielectron + o2::framework::Configurable maxRgen{"maxRgen", 90.f, "maximum radius for generated particles"}; + o2::framework::Configurable margin_z_mc{"margin_z_mc", 7.0, "margin for z cut in cm for MC"}; PHOSPhotonCut fPHOSCut; - struct : ConfigurableGroup { + struct : o2::framework::ConfigurableGroup { std::string prefix = "phoscut_group"; - Configurable cfg_min_Ecluster{"cfg_min_Ecluster", 0.3, "Minimum cluster energy for PHOS in GeV"}; + o2::framework::Configurable cfg_min_Ecluster{"cfg_min_Ecluster", 0.3, "Minimum cluster energy for PHOS in GeV"}; } phoscuts; TF1* f1fd_k0s_to_pi0; - HistogramRegistry fRegistry{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false}; - static constexpr std::string_view event_types[2] = {"before/", "after/"}; - static constexpr std::string_view event_pair_types[2] = {"same/", "mix/"}; - static constexpr std::string_view parnames[2] = {"Pi0/", "Eta/"}; + o2::framework::HistogramRegistry fRegistry{"output", {}, o2::framework::OutputObjHandlingPolicy::AnalysisObject, false, false}; + // static constexpr std::string_view event_types[2] = {"before/", "after/"}; + // static constexpr std::string_view event_pair_types[2] = {"same/", "mix/"}; + static constexpr std::string_view kParnames[2] = {"Pi0/", "Eta/"}; o2::ccdb::CcdbApi ccdbApi; - Service ccdb; + o2::framework::Service ccdb; int mRunNumber; float d_bz; - void init(InitContext&) + void init(o2::framework::InitContext&) { o2::aod::pwgem::photonmeson::utils::eventhistogram::addEventHistograms(&fRegistry); - if constexpr (pairtype == PairType::kPCMDalitzEE) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { o2::aod::pwgem::photonmeson::utils::nmhistogram::addNMHistograms(&fRegistry, true, "ee#gamma"); } else { o2::aod::pwgem::photonmeson::utils::nmhistogram::addNMHistograms(&fRegistry, true, "#gamma#gamma"); @@ -238,7 +219,7 @@ struct Pi0EtaToGammaGammaMC { f1fd_k0s_to_pi0 = new TF1("f1fd_k0s_to_pi0", TString(fd_k0s_to_pi0), 0.f, 100.f); - fRegistry.add("Event/hNrecPerMCCollision", "Nrec per mc collision;N_{rec} collisions per MC collision", kTH1F, {{21, -0.5f, 20.5f}}, false); + fRegistry.add("Event/hNrecPerMCCollision", "Nrec per mc collision;N_{rec} collisions per MC collision", o2::framework::kTH1F, {{21, -0.5f, 20.5f}}, false); mRunNumber = 0; d_bz = 0; @@ -404,14 +385,14 @@ struct Pi0EtaToGammaGammaMC { fPHOSCut.SetEnergyRange(phoscuts.cfg_min_Ecluster, 1e+10); } - SliceCache cache; - Preslice perCollision_pcm = aod::v0photonkf::emeventId; - Preslice perCollision_emc = aod::emccluster::emeventId; - Preslice perCollision_phos = aod::phoscluster::emeventId; + o2::framework::SliceCache cache; + o2::framework::Preslice>> perCollision_pcm = o2::aod::v0photonkf::emeventId; + o2::framework::Preslice> perCollision_emc = o2::aod::emccluster::emeventId; + o2::framework::Preslice> perCollision_phos = o2::aod::phoscluster::emeventId; - Preslice perCollision_electron = aod::emprimaryelectron::emeventId; - Partition positrons = o2::aod::emprimaryelectron::sign > int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt&& nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; - Partition electrons = o2::aod::emprimaryelectron::sign < int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt && nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; + o2::framework::Preslice>> perCollision_electron = o2::aod::emprimaryelectron::emeventId; + o2::framework::Partition>> positrons = o2::aod::emprimaryelectron::sign > int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt&& nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; + o2::framework::Partition>> electrons = o2::aod::emprimaryelectron::sign < int8_t(0) && dileptoncuts.cfg_min_pt_track < o2::aod::track::pt && nabs(o2::aod::track::eta) < dileptoncuts.cfg_max_eta_track; template void runTruePairing(TCollisions const& collisions, @@ -424,12 +405,12 @@ struct Pi0EtaToGammaGammaMC { for (auto& collision : collisions) { initCCDB(collision); - if ((pairtype == PairType::kPHOSPHOS || pairtype == PairType::kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { + if ((pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { continue; } float weight = 1.f; - if constexpr (std::is_same_v, FilteredMyCollisionsWithJJMC>) { + if constexpr (std::is_same_v, o2::soa::Filtered, o2::aod::EMEventsWeight>>>) { weight = collision.weight(); } @@ -451,16 +432,16 @@ struct Pi0EtaToGammaGammaMC { fRegistry.fill(HIST("Event/after/hCollisionCounter"), 12.0, weight); // accepted int photonid1 = -1, photonid2 = -1, pi0id = -1, etaid = -1; - if constexpr (pairtype == PairType::kPCMPCM || pairtype == PairType::kPHOSPHOS || pairtype == PairType::kEMCEMC) { // same kinds pairing + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { // same kinds pairing auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto photons2_per_collision = photons2.sliceBy(perCollision2, collision.globalIndex()); - for (auto& [g1, g2] : combinations(CombinationsStrictlyUpperIndexPolicy(photons1_per_collision, photons2_per_collision))) { + for (auto& [g1, g2] : o2::soa::combinations(o2::soa::CombinationsStrictlyUpperIndexPolicy(photons1_per_collision, photons2_per_collision))) { if (!cut1.template IsSelected(g1) || !cut2.template IsSelected(g2)) { continue; } - if constexpr (pairtype == PairType::kPCMPCM) { // check 2 legs + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { // check 2 legs auto pos1 = g1.template posTrack_as(); auto ele1 = g1.template negTrack_as(); auto pos2 = g2.template posTrack_as(); @@ -471,14 +452,14 @@ struct Pi0EtaToGammaGammaMC { auto pos2mc = pos2.template emmcparticle_as(); auto ele2mc = ele2.template emmcparticle_as(); - photonid1 = FindCommonMotherFrom2Prongs(pos1mc, ele1mc, -11, 11, 22, mcparticles); - photonid2 = FindCommonMotherFrom2Prongs(pos2mc, ele2mc, -11, 11, 22, mcparticles); - } else if constexpr (pairtype == PairType::kEMCEMC) { + photonid1 = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom2Prongs(pos1mc, ele1mc, -11, 11, 22, mcparticles); + photonid2 = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom2Prongs(pos2mc, ele2mc, -11, 11, 22, mcparticles); + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { auto cluster1mcparticle = mcparticles.iteratorAt(g1.emmcparticleId()); auto cluster2mcparticle = mcparticles.iteratorAt(g2.emmcparticleId()); - photonid1 = FindMotherInChain(cluster1mcparticle, mcparticles, std::vector{111, 221}); - photonid2 = FindMotherInChain(cluster2mcparticle, mcparticles, std::vector{111, 221}); + photonid1 = o2::aod::pwgem::photonmeson::utils::mcutil::FindMotherInChain(cluster1mcparticle, mcparticles, std::vector{111, 221}); + photonid2 = o2::aod::pwgem::photonmeson::utils::mcutil::FindMotherInChain(cluster2mcparticle, mcparticles, std::vector{111, 221}); } else { photonid1 = -1; photonid2 = -1; @@ -490,14 +471,14 @@ struct Pi0EtaToGammaGammaMC { auto g1mc = mcparticles.iteratorAt(photonid1); auto g2mc = mcparticles.iteratorAt(photonid2); - if constexpr (pairtype == PairType::kPCMPCM) { - if (!IsConversionPointInAcceptance(g1mc, maxRgen, maxY_gen, margin_z_mc, mcparticles) || !IsConversionPointInAcceptance(g2mc, maxRgen, maxY_gen, margin_z_mc, mcparticles)) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { + if (!o2::aod::pwgem::photonmeson::utils::mcutil::IsConversionPointInAcceptance(g1mc, maxRgen, maxY_gen, margin_z_mc, mcparticles) || !o2::aod::pwgem::photonmeson::utils::mcutil::IsConversionPointInAcceptance(g2mc, maxRgen, maxY_gen, margin_z_mc, mcparticles)) { continue; } } - pi0id = FindCommonMotherFrom2Prongs(g1mc, g2mc, 22, 22, 111, mcparticles); - etaid = FindCommonMotherFrom2Prongs(g1mc, g2mc, 22, 22, 221, mcparticles); + pi0id = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom2Prongs(g1mc, g2mc, 22, 22, 111, mcparticles); + etaid = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom2Prongs(g1mc, g2mc, 22, 22, 221, mcparticles); if (g1mc.globalIndex() != g2mc.globalIndex() && pi0id < 0 && etaid < 0) { // for same gamma no pi0/eta will be found, but we still want to fill the FromSameGamma hist continue; @@ -510,7 +491,7 @@ struct Pi0EtaToGammaGammaMC { continue; } - if (pairtype == PairType::kEMCEMC) { + if (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); if (openingAngle < emccuts.minOpenAngle) { continue; @@ -518,9 +499,9 @@ struct Pi0EtaToGammaGammaMC { } if (g1mc.globalIndex() == g2mc.globalIndex()) { - if (getMotherPDGCode(g1mc, mcparticles) == 111) + if (o2::aod::pwgem::dilepton::utils::mcutil::getMotherPDGCode(g1mc, mcparticles) == 111) fRegistry.fill(HIST("Pair/Pi0/hs_FromSameGamma"), v12.M(), v12.Pt(), weight); - else if (getMotherPDGCode(g1mc, mcparticles) == 221) + else if (o2::aod::pwgem::dilepton::utils::mcutil::getMotherPDGCode(g1mc, mcparticles) == 221) fRegistry.fill(HIST("Pair/Eta/hs_FromSameGamma"), v12.M(), v12.Pt(), weight); continue; } @@ -539,7 +520,7 @@ struct Pi0EtaToGammaGammaMC { o2::aod::pwgem::photonmeson::utils::nmhistogram::fillTruePairInfo(&fRegistry, v12, etamc, mcparticles, mccollisions, f1fd_k0s_to_pi0, weight); } } // end of pairing loop - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto positrons_per_collision = positrons->sliceByCached(o2::aod::emprimaryelectron::emeventId, collision.globalIndex(), cache); auto electrons_per_collision = electrons->sliceByCached(o2::aod::emprimaryelectron::emeventId, collision.globalIndex(), cache); @@ -552,18 +533,18 @@ struct Pi0EtaToGammaGammaMC { auto ele1 = g1.template negTrack_as(); auto pos1mc = pos1.template emmcparticle_as(); auto ele1mc = ele1.template emmcparticle_as(); - photonid1 = FindCommonMotherFrom2Prongs(pos1mc, ele1mc, -11, 11, 22, mcparticles); + photonid1 = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom2Prongs(pos1mc, ele1mc, -11, 11, 22, mcparticles); if (photonid1 < 0) { continue; } auto g1mc = mcparticles.iteratorAt(photonid1); - if (!IsConversionPointInAcceptance(g1mc, maxRgen, maxY_gen, margin_z_mc, mcparticles)) { + if (!o2::aod::pwgem::photonmeson::utils::mcutil::IsConversionPointInAcceptance(g1mc, maxRgen, maxY_gen, margin_z_mc, mcparticles)) { continue; } ROOT::Math::PtEtaPhiMVector v_gamma(g1.pt(), g1.eta(), g1.phi(), 0.f); - for (auto& [pos2, ele2] : combinations(CombinationsFullIndexPolicy(positrons_per_collision, electrons_per_collision))) { // ULS - if (pos2.trackId() == ele2.trackId()) { // this is protection against pairing identical 2 tracks. + for (auto& [pos2, ele2] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(positrons_per_collision, electrons_per_collision))) { // ULS + if (pos2.trackId() == ele2.trackId()) { // this is protection against pairing identical 2 tracks. continue; } if (pos1.trackId() == pos2.trackId() || ele1.trackId() == ele2.trackId()) { @@ -580,8 +561,8 @@ struct Pi0EtaToGammaGammaMC { auto pos2mc = mcparticles.iteratorAt(pos2.emmcparticleId()); auto ele2mc = mcparticles.iteratorAt(ele2.emmcparticleId()); - pi0id = FindCommonMotherFrom3Prongs(g1mc, pos2mc, ele2mc, 22, -11, 11, 111, mcparticles); - etaid = FindCommonMotherFrom3Prongs(g1mc, pos2mc, ele2mc, 22, -11, 11, 221, mcparticles); + pi0id = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom3Prongs(g1mc, pos2mc, ele2mc, 22, -11, 11, 111, mcparticles); + etaid = o2::aod::pwgem::dilepton::utils::mcutil::FindCommonMotherFrom3Prongs(g1mc, pos2mc, ele2mc, 22, -11, 11, 221, mcparticles); if (pi0id < 0 && etaid < 0) { continue; } @@ -610,7 +591,7 @@ struct Pi0EtaToGammaGammaMC { auto photons1_per_collision = photons1.sliceBy(perCollision1, collision.globalIndex()); auto photons2_per_collision = photons2.sliceBy(perCollision2, collision.globalIndex()); - for (auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photons1_per_collision, photons2_per_collision))) { + for (auto& [g1, g2] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(photons1_per_collision, photons2_per_collision))) { if (!cut1.template IsSelected(g1) || !cut2.template IsSelected(g2)) { continue; } @@ -636,8 +617,8 @@ struct Pi0EtaToGammaGammaMC { void fillBinnedData(TBinnedData const& binned_data, const float weight = 1.f) { int xbin = 0, ybin = 0, zbin = 0; - auto hPtY = fRegistry.get(HIST("Generated/") + HIST(parnames[par_id]) + HIST("hPtY")); // 2D - auto hPt = fRegistry.get(HIST("Generated/") + HIST(parnames[par_id]) + HIST("hPt")); // 1D + auto hPtY = fRegistry.get(HIST("Generated/") + HIST(kParnames[par_id]) + HIST("hPtY")); // 2D + auto hPt = fRegistry.get(HIST("Generated/") + HIST(kParnames[par_id]) + HIST("hPt")); // 1D for (int ibin = 0; ibin < hPtY->GetNcells(); ibin++) { int nentry = binned_data[ibin]; @@ -655,8 +636,8 @@ struct Pi0EtaToGammaGammaMC { } } - PresliceUnsorted perMcCollision = aod::emmcparticle::emmceventId; - PresliceUnsorted rec_perMcCollision = aod::emmceventlabel::emmceventId; + o2::framework::PresliceUnsorted perMcCollision = o2::aod::emmcparticle::emmceventId; + o2::framework::PresliceUnsorted> rec_perMcCollision = o2::aod::emmceventlabel::emmceventId; template void runGenInfo(TCollisions const& collisions, TMCCollisions const& mccollisions, TMCParticles const& /*mcparticles*/) @@ -671,12 +652,12 @@ struct Pi0EtaToGammaGammaMC { } for (auto& collision : collisions) { - if ((pairtype == kPHOSPHOS || pairtype == kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { + if ((pairtype == o2::aod::pwgem::photonmeson::photonpair::kPHOSPHOS || pairtype == o2::aod::pwgem::photonmeson::photonpair::kPCMPHOS) && !collision.alias_bit(triggerAliases::kTVXinPHOS)) { continue; // I don't know why this is necessary in simulation. } float weight = 1.f; - if constexpr (std::is_same_v, FilteredMyCollisionsWithJJMC>) { + if constexpr (std::is_same_v, o2::soa::Filtered, o2::aod::EMEventsWeight>>>) { weight = collision.weight(); } @@ -701,45 +682,45 @@ struct Pi0EtaToGammaGammaMC { } // end of collision loop } - Filter collisionFilter_occupancy_track = eventcuts.cfgTrackOccupancyMin <= o2::aod::evsel::trackOccupancyInTimeRange && o2::aod::evsel::trackOccupancyInTimeRange < eventcuts.cfgTrackOccupancyMax; - Filter collisionFilter_occupancy_ft0c = eventcuts.cfgFT0COccupancyMin <= o2::aod::evsel::ft0cOccupancyInTimeRange && o2::aod::evsel::ft0cOccupancyInTimeRange < eventcuts.cfgFT0COccupancyMax; - Filter collisionFilter_centrality = (cfgCentMin < o2::aod::cent::centFT0M && o2::aod::cent::centFT0M < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0A && o2::aod::cent::centFT0A < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0C && o2::aod::cent::centFT0C < cfgCentMax); - using FilteredMyCollisions = soa::Filtered; + o2::framework::expressions::Filter collisionFilter_occupancy_track = eventcuts.cfgTrackOccupancyMin <= o2::aod::evsel::trackOccupancyInTimeRange && o2::aod::evsel::trackOccupancyInTimeRange < eventcuts.cfgTrackOccupancyMax; + o2::framework::expressions::Filter collisionFilter_occupancy_ft0c = eventcuts.cfgFT0COccupancyMin <= o2::aod::evsel::ft0cOccupancyInTimeRange && o2::aod::evsel::ft0cOccupancyInTimeRange < eventcuts.cfgFT0COccupancyMax; + o2::framework::expressions::Filter collisionFilter_centrality = (cfgCentMin < o2::aod::cent::centFT0M && o2::aod::cent::centFT0M < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0A && o2::aod::cent::centFT0A < cfgCentMax) || (cfgCentMin < o2::aod::cent::centFT0C && o2::aod::cent::centFT0C < cfgCentMax); + // using FilteredMyCollisions = o2::soa::Filtered>; - Filter prefilter_pcm = ifnode(pcmcuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::v0photonkf::pfbderived == static_cast(0), true); - Filter prefilter_primaryelectron = ifnode(dileptoncuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::emprimaryelectron::pfbderived == static_cast(0), true); + o2::framework::expressions::Filter prefilter_pcm = ifnode(pcmcuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::v0photonkf::pfbderived == static_cast(0), true); + o2::framework::expressions::Filter prefilter_primaryelectron = ifnode(dileptoncuts.cfg_apply_cuts_from_prefilter_derived.node(), o2::aod::emprimaryelectron::pfbderived == static_cast(0), true); - void processAnalysis(FilteredMyCollisions const& collisions, MyMCCollisions const& mccollisions, aod::EMMCParticles const& mcparticles, Types const&... args) + void processAnalysis(o2::soa::Filtered> const& collisions, o2::soa::Join const& mccollisions, o2::aod::EMMCParticles const& mcparticles, Types const&... args) { - if constexpr (pairtype == PairType::kPCMPCM) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); runTruePairing(collisions, v0photons, v0photons, v0legs, v0legs, perCollision_pcm, perCollision_pcm, fV0PhotonCut, fV0PhotonCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); auto emprimaryelectrons = std::get<2>(std::tie(args...)); // LOGF(info, "electrons.size() = %d, positrons.size() = %d", electrons.size(), positrons.size()); runTruePairing(collisions, v0photons, emprimaryelectrons, v0legs, emprimaryelectrons, perCollision_pcm, perCollision_electron, fV0PhotonCut, fDileptonCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); - } else if constexpr (pairtype == PairType::kEMCEMC) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { auto emcclusters = std::get<0>(std::tie(args...)); runTruePairing(collisions, emcclusters, emcclusters, nullptr, nullptr, perCollision_emc, perCollision_emc, fEMCCut, fEMCCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); } - // else if constexpr (pairtype == PairType::kPHOSPHOS) { + // else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS) { // auto phosclusters = std::get<0>(std::tie(args...)); // runPairing(collisions, phosclusters, phosclusters, nullptr, nullptr, perCollision_phos, perCollision_phos, fPHOSCut, fPHOSCut, nullptr, nullptr); // } - // else if constexpr (pairtype == PairType::kPCMEMC) { + // else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMEMC) { // auto v0photons = std::get<0>(std::tie(args...)); // auto v0legs = std::get<1>(std::tie(args...)); // auto emcclusters = std::get<2>(std::tie(args...)); // auto emcmatchedtracks = std::get<3>(std::tie(args...)); // runPairing(collisions, v0photons, emcclusters, v0legs, nullptr, perCollision_pcm, perCollision_emc, fV0PhotonCut, fEMCCut, emcmatchedtracks, nullptr); - // } else if constexpr (pairtype == PairType::kPCMPHOS) { + // } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPHOS) { // auto v0photons = std::get<0>(std::tie(args...)); // auto v0legs = std::get<1>(std::tie(args...)); // auto phosclusters = std::get<2>(std::tie(args...)); @@ -748,38 +729,38 @@ struct Pi0EtaToGammaGammaMC { } PROCESS_SWITCH(Pi0EtaToGammaGammaMC, processAnalysis, "process pair analysis", true); - using FilteredMyCollisionsWithJJMC = soa::Filtered; - void processAnalysisJJMC(FilteredMyCollisionsWithJJMC const& collisions, MyMCCollisions const& mccollisions, aod::EMMCParticles const& mcparticles, Types const&... args) + // using FilteredMyCollisionsWithJJMC = o2::soa::Filtered, aod::EMEventsWeight>>; + void processAnalysisJJMC(o2::soa::Filtered, o2::aod::EMEventsWeight>> const& collisions, o2::soa::Join const& mccollisions, o2::aod::EMMCParticles const& mcparticles, Types const&... args) { - if constexpr (pairtype == PairType::kPCMPCM) { + if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPCM) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); runTruePairing(collisions, v0photons, v0photons, v0legs, v0legs, perCollision_pcm, perCollision_pcm, fV0PhotonCut, fV0PhotonCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); - } else if constexpr (pairtype == PairType::kPCMDalitzEE) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMDalitzEE) { auto v0photons = std::get<0>(std::tie(args...)); auto v0legs = std::get<1>(std::tie(args...)); auto emprimaryelectrons = std::get<2>(std::tie(args...)); // LOGF(info, "electrons.size() = %d, positrons.size() = %d", electrons.size(), positrons.size()); runTruePairing(collisions, v0photons, emprimaryelectrons, v0legs, emprimaryelectrons, perCollision_pcm, perCollision_electron, fV0PhotonCut, fDileptonCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); - } else if constexpr (pairtype == PairType::kEMCEMC) { + } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kEMCEMC) { auto emcclusters = std::get<0>(std::tie(args...)); runTruePairing(collisions, emcclusters, emcclusters, nullptr, nullptr, perCollision_emc, perCollision_emc, fEMCCut, fEMCCut, mccollisions, mcparticles); runGenInfo(collisions, mccollisions, mcparticles); } - // else if constexpr (pairtype == PairType::kPHOSPHOS) { + // else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPHOSPHOS) { // auto phosclusters = std::get<0>(std::tie(args...)); // runPairing(collisions, phosclusters, phosclusters, nullptr, nullptr, perCollision_phos, perCollision_phos, fPHOSCut, fPHOSCut, nullptr, nullptr); // } - // else if constexpr (pairtype == PairType::kPCMEMC) { + // else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMEMC) { // auto v0photons = std::get<0>(std::tie(args...)); // auto v0legs = std::get<1>(std::tie(args...)); // auto emcclusters = std::get<2>(std::tie(args...)); // auto emcmatchedtracks = std::get<3>(std::tie(args...)); // runPairing(collisions, v0photons, emcclusters, v0legs, nullptr, perCollision_pcm, perCollision_emc, fV0PhotonCut, fEMCCut, emcmatchedtracks, nullptr); - // } else if constexpr (pairtype == PairType::kPCMPHOS) { + // } else if constexpr (pairtype == o2::aod::pwgem::photonmeson::photonpair::PairType::kPCMPHOS) { // auto v0photons = std::get<0>(std::tie(args...)); // auto v0legs = std::get<1>(std::tie(args...)); // auto phosclusters = std::get<2>(std::tie(args...)); @@ -788,7 +769,7 @@ struct Pi0EtaToGammaGammaMC { } PROCESS_SWITCH(Pi0EtaToGammaGammaMC, processAnalysisJJMC, "process pair analysis", false); - void processDummy(MyCollisions const&) {} + void processDummy(o2::aod::EMEvents const&) {} PROCESS_SWITCH(Pi0EtaToGammaGammaMC, processDummy, "Dummy function", false); }; #endif // PWGEM_PHOTONMESON_CORE_PI0ETATOGAMMAGAMMAMC_H_ diff --git a/PWGEM/PhotonMeson/DataModel/gammaTables.h b/PWGEM/PhotonMeson/DataModel/gammaTables.h index 5f3392e1704..b10f56a2299 100644 --- a/PWGEM/PhotonMeson/DataModel/gammaTables.h +++ b/PWGEM/PhotonMeson/DataModel/gammaTables.h @@ -584,6 +584,12 @@ DECLARE_SOA_TABLE(SkimEMCClusters, "AOD", "SKIMEMCCLUSTER", //! table of skimmed emccluster::DeltaEta, emccluster::TrackP, emccluster::TrackPt, emccluster::Pt); using SkimEMCCluster = SkimEMCClusters::iterator; +DECLARE_SOA_TABLE_VERSIONED(SkimEMCClusters_001, "AOD", "SKIMEMCCLUSTER", 1, //! table of skimmed EMCal clusters + o2::soa::Index<>, skimmedcluster::CollisionId, emccluster::Definition, skimmedcluster::E, skimmedcluster::Eta, skimmedcluster::Phi, + skimmedcluster::M02, skimmedcluster::NCells, skimmedcluster::Time, emccluster::IsExotic, emccluster::DeltaPhi, + emccluster::DeltaEta, emccluster::TrackP, emccluster::TrackPt, emccluster::Pt); +using SkimEMCCluster = SkimEMCClusters::iterator; + DECLARE_SOA_TABLE(EMCEMEventIds, "AOD", "EMCEMEVENTID", emccluster::EMEventId); // To be joined with SkimEMCClusters table at analysis level. // iterators using EMCEMEventId = EMCEMEventIds::iterator; diff --git a/PWGEM/PhotonMeson/TableProducer/photonconversionbuilder.cxx b/PWGEM/PhotonMeson/TableProducer/photonconversionbuilder.cxx index 770cd6ad1c3..f84f260bc46 100644 --- a/PWGEM/PhotonMeson/TableProducer/photonconversionbuilder.cxx +++ b/PWGEM/PhotonMeson/TableProducer/photonconversionbuilder.cxx @@ -9,11 +9,10 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// \file photonconversionbuilder.cxx -// \brief this task produces photon data table with KFParticle. -// -// \author Daiki Sekihata , Tokyo + +/// \file photonconversionbuilder.cxx +/// \brief this task produces photon data table with KFParticle. +/// \author Daiki Sekihata , Tokyo #include "PWGEM/Dilepton/Utils/PairUtilities.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" @@ -22,31 +21,42 @@ #include "Common/Core/RecoDecay.h" #include "Common/Core/TPCVDriftManager.h" -#include "Common/Core/TableHelper.h" -#include "Common/Core/TrackSelection.h" #include "Common/Core/trackUtilities.h" -#include "Common/DataModel/TrackSelectionTables.h" +#include "Common/DataModel/EventSelection.h" +#include "Common/DataModel/PIDResponseTPC.h" #include "Tools/KFparticle/KFUtilities.h" -#include "CCDB/BasicCCDBManager.h" -#include "DataFormatsParameters/GRPMagField.h" -#include "DataFormatsParameters/GRPObject.h" -#include "DetectorsBase/GeometryManager.h" -#include "DetectorsBase/Propagator.h" -#include "Framework/ASoAHelpers.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/AnalysisTask.h" -#include "Framework/RunningWorkflowInfo.h" -#include "Framework/runDataProcessing.h" -#include "ReconstructionDataFormats/Track.h" - -#include "Math/Vector4D.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include // IWYU pragma: keep +#include +#include + +#include +#include +#include #include #include #include +#include #include -#include #include #include #include @@ -69,6 +79,12 @@ using MyCollisionsMC = soa::Join; using MyTracksIU = soa::Join; using MyTracksIUMC = soa::Join; +enum MatCorrType { + None = 0, + TGeo = 1, + LUT = 2 +}; + struct PhotonConversionBuilder { Produces v0photonskf; Produces v0legs; @@ -152,7 +168,7 @@ struct PhotonConversionBuilder { {"V0/hConversionPointXY", "conversion point in XY;X (cm);Y (cm)", {HistType::kTH2F, {{400, -100.0f, 100.0f}, {400, -100.f, 100.f}}}}, {"V0/hConversionPointRZ", "conversion point in RZ;Z (cm);R_{xy} (cm)", {HistType::kTH2F, {{200, -100.0f, 100.0f}, {200, 0.f, 100.f}}}}, {"V0/hPt", "pT of V0 at PV;p_{T,#gamma} (GeV/c)", {HistType::kTH1F, {{1000, 0.0f, 10.0f}}}}, - {"V0/hEtaPhi", "#eta vs. #varphi of V0 at PV;#varphi (rad.);#eta", {HistType::kTH2F, {{72, 0.0f, 2 * M_PI}, {200, -1, +1}}}}, + {"V0/hEtaPhi", "#eta vs. #varphi of V0 at PV;#varphi (rad.);#eta", {HistType::kTH2F, {{72, 0.0f, o2::constants::math::TwoPI}, {200, -1, +1}}}}, {"V0/hCosPA", "cosine of pointing angle;cosine of pointing angle", {HistType::kTH1F, {{100, 0.99f, 1.f}}}}, {"V0/hCosPA_Rxy", "cosine of pointing angle;r_{xy} (cm);cosine of pointing angle", {HistType::kTH2F, {{200, 0, 100}, {100, 0.99f, 1.f}}}}, {"V0/hCosPAXY_Rxy", "cosine of pointing angle;r_{xy} (cm);cosine of pointing angle", {HistType::kTH2F, {{200, 0, 100}, {100, 0.99f, 1.f}}}}, @@ -170,7 +186,7 @@ struct PhotonConversionBuilder { {"V0/hPCA_diffX", "PCA vs. trackiu X - R_{xy};distance btween 2 legs (cm);min trackiu X - R_{xy} (cm)", {HistType::kTH2F, {{500, 0.0f, 5.f}, {100, -50.0, 50.0f}}}}, {"V0/hPhiV", "#phi_{V}; #phi_{V} (rad.)", {HistType::kTH1F, {{500, 0.0f, 2 * M_PI}}}}, {"V0Leg/hPt", "pT of leg at SV;p_{T,e} (GeV/c)", {HistType::kTH1F, {{1000, 0.0f, 10.0f}}}}, - {"V0Leg/hEtaPhi", "#eta vs. #varphi of leg at SV;#varphi (rad.);#eta", {HistType::kTH2F, {{72, 0.0f, 2 * M_PI}, {200, -1, +1}}}}, + {"V0Leg/hEtaPhi", "#eta vs. #varphi of leg at SV;#varphi (rad.);#eta", {HistType::kTH2F, {{72, 0.0f, o2::constants::math::TwoPI}, {200, -1, +1}}}}, {"V0Leg/hRelDeltaPt", "pT resolution;p_{T} (GeV/c);#Deltap_{T}/p_{T}", {HistType::kTH2F, {{1000, 0.f, 10.f}, {100, 0, 1}}}}, {"V0Leg/hDCAxyz", "DCA xy vs. z to PV;DCA_{xy} (cm);DCA_{z} (cm)", {HistType::kTH2F, {{200, -50.f, 50.f}, {200, -50.f, +50.f}}}}, {"V0Leg/hdEdx_Pin", "TPC dE/dx vs. p_{in};p_{in} (GeV/c);TPC dE/dx", {HistType::kTH2F, {{1000, 0.f, 10.f}, {200, 0.f, 200.f}}}}, @@ -190,21 +206,21 @@ struct PhotonConversionBuilder { ccdb->setLocalObjectValidityChecking(); ccdb->setFatalWhenNull(false); - if (useMatCorrType == 1) { + if (useMatCorrType == MatCorrType::TGeo) { LOGF(info, "TGeo correction requested, loading geometry"); if (!o2::base::GeometryManager::isGeometryLoaded()) { ccdb->get(geoPath); } } - if (useMatCorrType == 2) { + if (useMatCorrType == MatCorrType::LUT) { LOGF(info, "LUT correction requested, loading LUT"); lut = o2::base::MatLayerCylSet::rectifyPtrFromFile(ccdb->get(lutPath)); } - if (useMatCorrType == 1) { + if (useMatCorrType == MatCorrType::TGeo) { matCorr = o2::base::Propagator::MatCorrType::USEMatCorrTGeo; } - if (useMatCorrType == 2) { + if (useMatCorrType == MatCorrType::LUT) { matCorr = o2::base::Propagator::MatCorrType::USEMatCorrLUT; } } @@ -467,8 +483,8 @@ struct PhotonConversionBuilder { KFPTrack kfp_track_pos = createKFPTrackFromTrackParCov(pTrack, pos.sign(), pos.tpcNClsFound(), pos.tpcChi2NCl()); KFPTrack kfp_track_ele = createKFPTrackFromTrackParCov(nTrack, ele.sign(), ele.tpcNClsFound(), ele.tpcChi2NCl()); - KFParticle kfp_pos(kfp_track_pos, -11); - KFParticle kfp_ele(kfp_track_ele, 11); + KFParticle kfp_pos(kfp_track_pos, kPositron); + KFParticle kfp_ele(kfp_track_ele, kElectron); const KFParticle* GammaDaughters[2] = {&kfp_pos, &kfp_ele}; KFParticle gammaKF; @@ -553,7 +569,7 @@ struct PhotonConversionBuilder { gammaKF_PV.SetProductionVertex(KFPV); float v0pt = RecoDecay::sqrtSumOfSquares(gammaKF_PV.GetPx(), gammaKF_PV.GetPy()); float v0eta = RecoDecay::eta(std::array{gammaKF_PV.GetPx(), gammaKF_PV.GetPy(), gammaKF_PV.GetPz()}); - float v0phi = RecoDecay::phi(gammaKF_PV.GetPx(), gammaKF_PV.GetPy()) > 0.f ? RecoDecay::phi(gammaKF_PV.GetPx(), gammaKF_PV.GetPy()) : RecoDecay::phi(gammaKF_PV.GetPx(), gammaKF_PV.GetPy()) + TMath::TwoPi(); + float v0phi = RecoDecay::constrainAngle(RecoDecay::phi(gammaKF_PV.GetPx(), gammaKF_PV.GetPy())); // KFParticle gammaKF_DecayVtx2 = gammaKF; // gammaKF_DecayVtx2.SetProductionVertex(KFPV); @@ -657,18 +673,18 @@ struct PhotonConversionBuilder { registry.fill(HIST("V0/hCosPAXY_Rxy"), rxy, cospaXY_kf); registry.fill(HIST("V0/hCosPARZ_Rxy"), rxy, cospaRZ_kf); - for (auto& leg : {kfp_pos_DecayVtx, kfp_ele_DecayVtx}) { + for (const auto& leg : {kfp_pos_DecayVtx, kfp_ele_DecayVtx}) { float legpt = RecoDecay::sqrtSumOfSquares(leg.GetPx(), leg.GetPy()); float legeta = RecoDecay::eta(std::array{leg.GetPx(), leg.GetPy(), leg.GetPz()}); - float legphi = RecoDecay::phi(leg.GetPx(), leg.GetPy()) > 0.f ? RecoDecay::phi(leg.GetPx(), leg.GetPy()) : RecoDecay::phi(leg.GetPx(), leg.GetPy()) + TMath::TwoPi(); + float legphi = RecoDecay::constrainAngle(RecoDecay::phi(leg.GetPx(), leg.GetPy())); registry.fill(HIST("V0Leg/hPt"), legpt); registry.fill(HIST("V0Leg/hEtaPhi"), legphi, legeta); } // end of leg loop - for (auto& leg : {pos, ele}) { + for (const auto& leg : {pos, ele}) { registry.fill(HIST("V0Leg/hdEdx_Pin"), leg.tpcInnerParam(), leg.tpcSignal()); registry.fill(HIST("V0Leg/hTPCNsigmaEl"), leg.tpcInnerParam(), leg.tpcNSigmaEl()); } // end of leg loop - for (auto& leg : {pTrack, nTrack}) { + for (const auto& leg : {pTrack, nTrack}) { registry.fill(HIST("V0Leg/hXZ"), leg.getZ(), leg.getX()); registry.fill(HIST("V0Leg/hRelDeltaPt"), leg.getPt(), leg.getPt() * std::sqrt(leg.getSigma1Pt2())); } // end of leg loop @@ -793,7 +809,7 @@ struct PhotonConversionBuilder { } // end of pca_map loop // LOGF(info, "pca_map.size() = %d", pca_map.size()); - for (auto& fullv0Id : stored_fullv0Ids) { + for (const auto& fullv0Id : stored_fullv0Ids) { auto v0Id = std::get<0>(fullv0Id); // auto collisionId = std::get<1>(fullv0Id); // auto posId = std::get<2>(fullv0Id); @@ -812,7 +828,7 @@ struct PhotonConversionBuilder { fillV0Table(v0, true); } // end of fullv0Id loop - for (auto& collision : collisions) { + for (const auto& collision : collisions) { if constexpr (isMC) { if (!collision.has_mcCollision()) { continue; diff --git a/PWGEM/PhotonMeson/TableProducer/skimmerGammaCalo.cxx b/PWGEM/PhotonMeson/TableProducer/skimmerGammaCalo.cxx index 731a9ed9331..cc253c54746 100644 --- a/PWGEM/PhotonMeson/TableProducer/skimmerGammaCalo.cxx +++ b/PWGEM/PhotonMeson/TableProducer/skimmerGammaCalo.cxx @@ -10,9 +10,9 @@ // or submit itself to any jurisdiction. /// \file skimmerGammaCalo.cxx -/// \brief skim cluster information to write photon cluster table in AO2D.root -/// dependencies: emcal-correction-task +/// \brief skim cluster information to write photon cluster table into derived AO2D.root /// \author marvin.hemmer@cern.ch +/// dependencies: emcal-correction-task #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h" @@ -58,9 +58,11 @@ struct SkimmerGammaCalo { Configurable minM02{"minM02", 0.0, "Minimum M02 for M02 cut"}; Configurable maxM02{"maxM02", 1.0, "Maximum M02 for M02 cut"}; Configurable minE{"minE", 0.5, "Minimum energy for energy cut"}; + Configurable removeExotic{"removeExotic", false, "Flag to enable the removal of exotic clusters."}; Configurable> clusterDefinitions{"clusterDefinitions", {0, 1, 2, 10, 11, 12, 13, 20, 21, 22, 30, 40, 41, 42, 43, 44, 45}, "Cluster definitions to be accepted (e.g. 13 for kV3MostSplitLowSeed)"}; Configurable maxdEta{"maxdEta", 0.1, "Set a maximum difference in eta for tracks and cluster to still count as matched"}; Configurable maxdPhi{"maxdPhi", 0.1, "Set a maximum difference in phi for tracks and cluster to still count as matched"}; + Configurable maxEoverP{"maxEoverP", 1.5, "Set a maximum for cluster E / track p for track matching."}; Configurable needEMCTrigger{"needEMCTrigger", false, "flag to only save events which have kTVXinEMC trigger bit. To reduce PbPb derived data size"}; HistogramRegistry historeg{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false}; @@ -69,25 +71,38 @@ struct SkimmerGammaCalo { { historeg.add("DefinitionIn", "Cluster definitions before cuts;#bf{Cluster definition};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{51, -0.5, 50.5}}); historeg.add("DefinitionOut", "Cluster definitions after cuts;#bf{Cluster definition};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{51, -0.5, 50.5}}); - historeg.add("EIn", "Energy of clusters before cuts", gHistoSpec_clusterE); - historeg.add("EOut", "Energy of clusters after cuts", gHistoSpec_clusterE); - historeg.add("MTEtaPhi", "Eta phi of matched tracks", gHistoSpec_clusterTM_dEtadPhi); + historeg.add("EIn", "Energy of clusters before cuts", gHistoSpecClusterE); + historeg.add("EOut", "Energy of clusters after cuts", gHistoSpecClusterE); + historeg.add("MTEtaPhiBeforeTM", "Eta phi of matched tracks before TM cuts", gHistoSpecClusterTMdEtadPhi); + historeg.add("MTEtaPhiAfterTM", "Eta phi of matched tracks after TM cuts", gHistoSpecClusterTMdEtadPhi); + historeg.add("Eoverp", "E/p for cluster E and track p", gHistoSpecTMEoverP); historeg.add("M02In", "Shape of cluster before cuts;#bf{#it{M}_{02}};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{200, 0, 2}}); historeg.add("M02Out", "Shape of cluster after cuts;#bf{#it{M}_{02}};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{200, 0, 2}}); historeg.add("TimeIn", "Time of cluster before cuts;#bf{#it{t} (ns)};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{200, -100, 100}}); historeg.add("TimeOut", "Time of cluster after cuts;#bf{#it{t} (ns)};#bf{#it{N}_{clusters}}", HistType::kTH1F, {{200, -100, 100}}); - auto hCaloClusterFilter = historeg.add("hCaloClusterFilter", "hCaloClusterFilter", kTH1I, {{6, 0, 6}}); + auto hCaloClusterFilter = historeg.add("hCaloClusterFilter", "hCaloClusterFilter", kTH1I, {{7, 0, 7}}); hCaloClusterFilter->GetXaxis()->SetBinLabel(1, "in"); hCaloClusterFilter->GetXaxis()->SetBinLabel(2, "Definition cut"); hCaloClusterFilter->GetXaxis()->SetBinLabel(3, "E cut"); hCaloClusterFilter->GetXaxis()->SetBinLabel(4, "time cut"); hCaloClusterFilter->GetXaxis()->SetBinLabel(5, "M02 cut"); - hCaloClusterFilter->GetXaxis()->SetBinLabel(6, "out"); + hCaloClusterFilter->GetXaxis()->SetBinLabel(6, "exotic cut"); + hCaloClusterFilter->GetXaxis()->SetBinLabel(7, "out"); + + auto hCaloTrackFilter = historeg.add("hCaloTrackFilter", "hCaloTrackFilter", kTH1I, {{4, 0, 4}}); + hCaloTrackFilter->GetXaxis()->SetBinLabel(1, "in"); + hCaloTrackFilter->GetXaxis()->SetBinLabel(2, "#Delta#eta #Delta#varphi"); + hCaloTrackFilter->GetXaxis()->SetBinLabel(3, "E/p cut"); + hCaloTrackFilter->GetXaxis()->SetBinLabel(4, "out"); + LOG(info) << "| EMCal cluster cuts for skimming:"; LOG(info) << "| Timing cut: " << minTime << " < t < " << maxTime; LOG(info) << "| M02 cut: " << minM02 << " < M02 < " << maxM02; LOG(info) << "| E cut: E > " << minE; + LOG(info) << "| TM - dPhi cut: dPhi < " << maxdPhi; + LOG(info) << "| TM - dEta cut: dEta < " << maxdEta; + LOG(info) << "| TM - E/p cut: E/p < " << maxEoverP; } void processRec(soa::Join::iterator const& collision, aod::EMCALClusters const& emcclusters, aod::EMCALClusterCells const& emcclustercells, aod::EMCALMatchedTracks const& emcmatchedtracks, aod::FullTracks const&) @@ -103,16 +118,15 @@ struct SkimmerGammaCalo { historeg.fill(HIST("hCaloClusterFilter"), 0); historeg.fill(HIST("DefinitionIn"), emccluster.definition()); - historeg.fill(HIST("EIn"), emccluster.energy()); historeg.fill(HIST("M02In"), emccluster.m02()); historeg.fill(HIST("TimeIn"), emccluster.time()); + historeg.fill(HIST("EIn"), emccluster.energy()); // Definition cut if (!(std::find(clusterDefinitions.value.begin(), clusterDefinitions.value.end(), emccluster.definition()) != clusterDefinitions.value.end())) { historeg.fill(HIST("hCaloClusterFilter"), 1); continue; } - historeg.fill(HIST("EIn"), emccluster.energy()); // Energy cut if (emccluster.energy() < minE) { historeg.fill(HIST("hCaloClusterFilter"), 2); @@ -128,12 +142,11 @@ struct SkimmerGammaCalo { historeg.fill(HIST("hCaloClusterFilter"), 4); continue; } - historeg.fill(HIST("hCaloClusterFilter"), 5); - - historeg.fill(HIST("DefinitionOut"), emccluster.definition()); - historeg.fill(HIST("EOut"), emccluster.energy()); - historeg.fill(HIST("M02Out"), emccluster.m02()); - historeg.fill(HIST("TimeOut"), emccluster.time()); + if (removeExotic.value && emccluster.isExotic()) { + historeg.fill(HIST("hCaloClusterFilter"), 5); + continue; + } + historeg.fill(HIST("hCaloClusterFilter"), 6); // Skimmed cell table auto groupedCells = emcclustercells.sliceBy(psCellperCluster, emccluster.globalIndex()); @@ -154,10 +167,19 @@ struct SkimmerGammaCalo { vP.reserve(groupedMTs.size()); vPt.reserve(groupedMTs.size()); for (const auto& emcmatchedtrack : groupedMTs) { + historeg.fill(HIST("hCaloTrackFilter"), 0); + historeg.fill(HIST("MTEtaPhiBeforeTM"), emcmatchedtrack.deltaEta(), emcmatchedtrack.deltaPhi()); if (std::abs(emcmatchedtrack.deltaEta()) >= maxdEta || std::abs(emcmatchedtrack.deltaPhi()) >= maxdPhi) { + historeg.fill(HIST("hCaloTrackFilter"), 1); + continue; + } + historeg.fill(HIST("Eoverp"), emccluster.energy(), emccluster.energy() / emcmatchedtrack.track_as().p()); + if (emccluster.energy() / emcmatchedtrack.track_as().p() > maxEoverP) { + historeg.fill(HIST("hCaloTrackFilter"), 2); continue; } - historeg.fill(HIST("MTEtaPhi"), emccluster.eta() - emcmatchedtrack.track_as().trackEtaEmcal(), emccluster.phi() - emcmatchedtrack.track_as().trackPhiEmcal()); + historeg.fill(HIST("hCaloTrackFilter"), 3); + historeg.fill(HIST("MTEtaPhiAfterTM"), emcmatchedtrack.deltaEta(), emcmatchedtrack.deltaPhi()); vTrackIds.emplace_back(emcmatchedtrack.trackId()); vEta.emplace_back(emcmatchedtrack.deltaEta()); vPhi.emplace_back(emcmatchedtrack.deltaPhi()); @@ -165,6 +187,12 @@ struct SkimmerGammaCalo { vPt.emplace_back(emcmatchedtrack.track_as().pt()); } + historeg.fill(HIST("DefinitionOut"), emccluster.definition()); + historeg.fill(HIST("EOut"), emccluster.energy()); + historeg.fill(HIST("M02Out"), emccluster.m02()); + historeg.fill(HIST("TimeOut"), emccluster.time()); + historeg.fill(HIST("hCaloClusterFilter"), 6); + tableGammaEMCReco(emccluster.collisionId(), emccluster.definition(), emccluster.energy(), emccluster.eta(), emccluster.phi(), emccluster.m02(), emccluster.nCells(), emccluster.time(), emccluster.isExotic(), vPhi, vEta, vP, vPt); } diff --git a/PWGEM/PhotonMeson/Tasks/CMakeLists.txt b/PWGEM/PhotonMeson/Tasks/CMakeLists.txt index 332cafc0795..8dc5bdc0e01 100644 --- a/PWGEM/PhotonMeson/Tasks/CMakeLists.txt +++ b/PWGEM/PhotonMeson/Tasks/CMakeLists.txt @@ -172,3 +172,7 @@ o2physics_add_dpl_workflow(diphoton-hadron-mpc-pcmdalitzee PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::PWGEMPhotonMesonCore COMPONENT_NAME Analysis) +o2physics_add_dpl_workflow(calib-task-emc + SOURCES calibTaskEmc.cxx + PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2::EMCALBase O2::EMCALCalib O2Physics::PWGEMPhotonMesonCore + COMPONENT_NAME Analysis) diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaEMCEMC.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaEMCEMC.cxx index 3b7893f7373..b1b657f8a1f 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaEMCEMC.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaEMCEMC.cxx @@ -9,25 +9,24 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaEMCEMC.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses for EMC-EMC. +/// \author D. Sekihata, daiki.sekihata@cern.ch +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; using MyEMCClusters = soa::Join; -using MyEMCCluster = MyEMCClusters::iterator; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCEMCEMC.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCEMCEMC.cxx index 7cd7ba3c492..982eba314e1 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCEMCEMC.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCEMCEMC.cxx @@ -9,25 +9,24 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaMCEMCEMC.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses in MC for EMC-EMC. +/// \author D. Sekihata, daiki.sekihata@cern.ch +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; using MyEMCClusters = soa::Join; -using MyEMCCluster = MyEMCClusters::iterator; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMDalitzEE.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMDalitzEE.cxx index 2d929f1a82a..e36da8a2556 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMDalitzEE.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMDalitzEE.cxx @@ -9,31 +9,26 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaMCPCMDalitzEE.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses in MC for PCM-Dalitz. +/// \author D. Sekihata, daiki.sekihata@cern.ch +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; -// using MyV0Photons = soa::Join; -// using MyV0Photon = MyV0Photons::iterator; -// -// using MyMCV0Legs = soa::Join; -// using MyMCV0Leg = MyMCV0Legs::iterator; -// -// using MyMCElectrons = soa::Join; -// using MyMCElectron = MyMCElectrons::iterator; +using MyV0Photons = soa::Filtered>; +using MyMCV0Legs = soa::Join; +using MyMCElectrons = soa::Filtered>; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMPCM.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMPCM.cxx index 9b79891739c..62c57593fd0 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMPCM.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaMCPCMPCM.cxx @@ -9,28 +9,25 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaMCPCMPCM.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses in MC for PCM-PCM. +/// \author D. Sekihata, daiki.sekihata@cern.ch +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGammaMC.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; -// using MyV0Photons = soa::Join; -// using MyV0Photon = MyV0Photons::iterator; - -// using MyMCV0Legs = soa::Join; -// using MyMCV0Leg = MyMCV0Legs::iterator; +using MyV0Photons = o2::soa::Filtered>; +using MyMCV0Legs = soa::Join; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMDalitzEE.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMDalitzEE.cxx index 15ccbe72bd4..b2c39f2093a 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMDalitzEE.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMDalitzEE.cxx @@ -9,23 +9,25 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaPCMDalitzEE.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses for PCM-Dalitz. +/// \author D. Sekihata, daiki.sekihata@cern.ch -#include "Common/Core/RecoDecay.h" +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; + +using MyV0Photons = o2::soa::Filtered>; +using MyPrimaryElectrons = o2::soa::Filtered>; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMPCM.cxx b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMPCM.cxx index 24557c828a7..6272af75da6 100644 --- a/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMPCM.cxx +++ b/PWGEM/PhotonMeson/Tasks/Pi0EtaToGammaGammaPCMPCM.cxx @@ -9,25 +9,24 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. // -// ======================== -// -// This code loops over photons and makes pairs for neutral mesons analyses. -// Please write to: daiki.sekihata@cern.ch - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/AnalysisDataModel.h" -#include "Framework/ASoAHelpers.h" +/// \file Pi0EtaToGammaGammaPCMPCM.cxx +/// \brief This code loops over photons and makes pairs for neutral mesons analyses for PCM-PCM. +/// \author D. Sekihata, daiki.sekihata@cern.ch +#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "PWGEM/PhotonMeson/Core/Pi0EtaToGammaGamma.h" + +#include +#include +#include using namespace o2; using namespace o2::aod; +using namespace o2::framework; +using namespace o2::aod::pwgem::photonmeson::photonpair; -// using MyV0Photons = soa::Join; -// using MyV0Photon = MyV0Photons::iterator; +using MyV0Photons = o2::soa::Filtered>; WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) { diff --git a/PWGEM/PhotonMeson/Tasks/calibTaskEmc.cxx b/PWGEM/PhotonMeson/Tasks/calibTaskEmc.cxx new file mode 100644 index 00000000000..4df49aad167 --- /dev/null +++ b/PWGEM/PhotonMeson/Tasks/calibTaskEmc.cxx @@ -0,0 +1,1071 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file calibTaskEmc.cxx +/// \brief Task to produce calibration values for EMCal +/// \author M. Hemmer, marvin.hemmer@cern.ch + +#include "PWGEM/PhotonMeson/Core/EMCPhotonCut.h" +#include "PWGEM/PhotonMeson/Core/EMPhotonEventCut.h" +#include "PWGEM/PhotonMeson/Core/V0PhotonCut.h" +#include "PWGEM/PhotonMeson/DataModel/gammaTables.h" +#include "PWGEM/PhotonMeson/Utils/EventHistograms.h" +#include "PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h" +// +#include "PWGEM/Dilepton/Utils/EMTrack.h" +#include "PWGEM/Dilepton/Utils/EventMixingHandler.h" + +#include "Common/Core/RecoDecay.h" +#include "Common/DataModel/Centrality.h" +#include "Common/DataModel/EventSelection.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include // IWYU pragma: keep +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace o2; +using namespace o2::aod; +using namespace o2::framework; +using namespace o2::framework::expressions; +using namespace o2::soa; +using namespace o2::aod::pwgem::photon; + +enum QvecEstimator { + FT0M = 0, + FT0A = 1, + FT0C, + TPCPos, + TPCNeg, + TPCTot +}; + +enum CentralityEstimator { + None = 0, + CFT0A = 1, + CFT0C, + CFT0M, + NCentralityEstimators +}; + +enum Harmonics { + kNone = 0, + kDirect = 1, + kElliptic = 2, + kTriangluar = 3, + kQuadrangular = 4, + kPentagonal = 5, + kHexagonal = 6, + kHeptagonal = 7, + kOctagonal = 8 +}; + +struct CalibTaskEmc { + static constexpr float MinEnergy = 0.7f; + + // configurable for flow + Configurable centEstimator{"centEstimator", 2, "Centrality estimation (FT0A: 1, FT0C: 2, FT0M: 3)"}; + Configurable ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"}; + Configurable cfgDoRotation{"cfgDoRotation", false, "Flag to enable rotation background method"}; + Configurable cfgDownsampling{"cfgDownsampling", 1, "Calculate rotation background only for every collision"}; + Configurable cfgEMCalMapLevelBackground{"cfgEMCalMapLevelBackground", 4, "Different levels of correction for the background, the smaller number includes the level of the higher number (4: none, 3: only inside EMCal, 2: exclude bad channels, 1: remove edges)"}; + Configurable cfgEMCalMapLevelSameEvent{"cfgEMCalMapLevelSameEvent", 4, "Different levels of correction for the same event, the smaller number includes the level of the higher number (4: none, 3: only inside EMCal, 2: exclude bad channels, 1: remove edges)"}; + Configurable cfgRotAngle{"cfgRotAngle", std::move(const_cast(o2::constants::math::PIHalf)), "Angle used for the rotation method"}; + Configurable cfgDistanceToEdge{"cfgDistanceToEdge", 1, "Distance to edge in cells required for rotated cluster to be accepted"}; + Configurable cfgMaxAsymmetry{"cfgMaxAsymmetry", 0.1f, "Maximum allowed asymmetry for photon pairs used in calibration when using EMC-EMC."}; + + // configurable axis + ConfigurableAxis thnConfigAxisInvMass{"thnConfigAxisInvMass", {400, 0.0, 0.8}, ""}; + ConfigurableAxis thnConfigAxisPt{"thnConfigAxisPt", {100, 0., 20.}, ""}; + ConfigurableAxis thnConfigAxisCent{"thnConfigAxisCent", {20, 0., 100.}, ""}; + ConfigurableAxis thnConfigAxisEnergyCalib{"thnConfigAxisEnergyCalib", {200, 0., 20.}, ""}; + + EMPhotonEventCut fEMEventCut; + struct : ConfigurableGroup { + std::string prefix = "eventcuts"; + Configurable cfgZvtxMax{"cfgZvtxMax", 10.f, "max. Zvtx"}; + Configurable cfgRequireSel8{"cfgRequireSel8", true, "require sel8 in event cut"}; + Configurable cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND in event cut"}; + Configurable cfgRequireNoTFB{"cfgRequireNoTFB", false, "require No time frame border in event cut"}; + Configurable cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border in event cut"}; + Configurable cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"}; + Configurable cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV. + Configurable cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"}; + Configurable cfgRequireEMCReadoutInMB{"cfgRequireEMCReadoutInMB", true, "require the EMC to be read out in an MB collision (kTVXinEMC)"}; + Configurable cfgRequireEMCHardwareTriggered{"cfgRequireEMCHardwareTriggered", false, "require the EMC to be hardware triggered (kEMC7 or kDMC7)"}; + Configurable cfgTrackOccupancyMin{"cfgTrackOccupancyMin", -1, "min. track occupancy"}; + Configurable cfgTrackOccupancyMax{"cfgTrackOccupancyMax", 1000000000, "max. track occupancy"}; + Configurable cfgFT0COccupancyMin{"cfgFT0COccupancyMin", -1, "min. FT0C occupancy"}; + Configurable cfgFT0COccupancyMax{"cfgFT0COccupancyMax", 1000000000, "max. FT0C occupancy"}; + Configurable cfgMinCent{"cfgMinCent", 0, "min. centrality (%)"}; + Configurable cfgMaxCent{"cfgMaxCent", 90, "max. centrality (%)"}; + Configurable onlyKeepWeightedEvents{"onlyKeepWeightedEvents", false, "flag to keep only weighted events (for JJ MCs) and remove all MB events (with weight = 1)"}; + Configurable cfgEnableQA{"cfgEnableQA", false, "flag to turn QA plots on/off"}; + } eventcuts; + + EMCPhotonCut fEMCCut; + struct : ConfigurableGroup { + std::string prefix = "emccuts"; + Configurable clusterDefinition{"clusterDefinition", "kV3Default", "Clusterizer to be selected, e.g. V3Default"}; + Configurable cfgEMCminTime{"cfgEMCminTime", -25., "Minimum cluster time for EMCal time cut"}; + Configurable cfgEMCmaxTime{"cfgEMCmaxTime", +30., "Maximum cluster time for EMCal time cut"}; + Configurable cfgEMCminM02{"cfgEMCminM02", 0.1, "Minimum M02 for EMCal M02 cut"}; + Configurable cfgEMCmaxM02{"cfgEMCmaxM02", 0.7, "Maximum M02 for EMCal M02 cut"}; + Configurable cfgEMCminE{"cfgEMCminE", 0.7, "Minimum cluster energy for EMCal energy cut"}; + Configurable cfgEMCminNCell{"cfgEMCminNCell", 1, "Minimum number of cells per cluster for EMCal NCell cut"}; + Configurable> cfgEMCTMEta{"cfgEMCTMEta", {0.01f, 4.07f, -2.5f}, "|eta| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + Configurable> cfgEMCTMPhi{"cfgEMCTMPhi", {0.015f, 3.65f, -2.f}, "|phi| <= [0]+(pT+[1])^[2] for EMCal track matching"}; + Configurable cfgEMCEoverp{"cfgEMCEoverp", 1.75, "Minimum cluster energy over track momentum for EMCal track matching"}; + Configurable cfgEMCUseExoticCut{"cfgEMCUseExoticCut", true, "FLag to use the EMCal exotic cluster cut"}; + Configurable cfgEMCUseTM{"cfgEMCUseTM", true, "flag to use EMCal track matching cut or not"}; + Configurable cfgEnableQA{"cfgEnableQA", false, "flag to turn QA plots on/off"}; + } emccuts; + + V0PhotonCut fPCMPhotonCut; + struct : ConfigurableGroup { + std::string prefix = "pcmcuts"; + Configurable cfgRequireV0WithITSTPC{"cfgRequireV0WithITSTPC", false, "flag to select V0s with ITS-TPC matched tracks"}; + Configurable cfgRequireV0WithITSonly{"cfgRequireV0WithITSonly", false, "flag to select V0s with ITSonly tracks"}; + Configurable cfgRequireV0WithTPConly{"cfgRequireV0WithTPConly", false, "flag to select V0s with TPConly tracks"}; + Configurable cfgMinPtV0{"cfgMinPtV0", 0.1, "min pT for v0 photons at PV"}; + Configurable cfgMaxPtV0{"cfgMaxPtV0", 1e+10, "max pT for v0 photons at PV"}; + Configurable cfgMinEtaV0{"cfgMinEtaV0", -0.8, "min eta for v0 photons at PV"}; + Configurable cfgMaxEtaV0{"cfgMaxEtaV0", 0.8, "max eta for v0 photons at PV"}; + Configurable cfgMinV0Radius{"cfgMinV0Radius", 4.0, "min v0 radius"}; + Configurable cfgMaxV0Radius{"cfgMaxV0Radius", 90.0, "max v0 radius"}; + Configurable cfgMaxAlphaAP{"cfgMaxAlphaAP", 0.95, "max alpha for AP cut"}; + Configurable cfgMaxQtAP{"cfgMaxQtAP", 0.01, "max qT for AP cut"}; + Configurable cfgMinCosPA{"cfgMinCosPA", 0.999, "min V0 CosPA"}; + Configurable cfgMaxPCA{"cfgMaxPCA", 1.5, "max distance btween 2 legs"}; + Configurable cfgMaxChi2KF{"cfgMaxChi2KF", 1e+10, "max chi2/ndf with KF"}; + Configurable cfgRejectV0OnITSib{"cfgRejectV0OnITSib", true, "flag to reject V0s on ITSib"}; + + Configurable cfgMinNClusterTPC{"cfgMinNClusterTPC", 0, "min ncluster tpc"}; + Configurable cfgMinNCrossedRows{"cfgMinNCrossedRows", 40, "min ncrossed rows"}; + Configurable cfgMaxFracSharedClusterTPC{"cfgMaxFracSharedClusterTPC", 999.f, "max fraction of shared clusters in TPC"}; + Configurable cfgMaxChi2TPC{"cfgMaxChi2TPC", 4.0, "max chi2/NclsTPC"}; + Configurable cfgMaxChi2ITS{"cfgMaxChi2ITS", 36.0, "max chi2/NclsITS"}; + Configurable cfgMinTPCNSigmaEl{"cfgMinTPCNSigmaEl", -3.0, "min. TPC n sigma for electron"}; + Configurable cfgMaxTPCNSigmaEl{"cfgMaxTPCNSigmaEl", +3.0, "max. TPC n sigma for electron"}; + Configurable cfgDisableITSOnly{"cfgDisableITSOnly", false, "flag to disable ITSonly tracks"}; + Configurable cfgDisableTPCOnly{"cfgDisableTPCOnly", false, "flag to disable TPConly tracks"}; + } pcmcuts; + + struct : ConfigurableGroup { + std::string prefix = "mesonConfig"; + Configurable minOpenAngle{"minOpenAngle", 0.0202, "apply min opening angle. Default value one EMCal cell"}; + Configurable enableTanThetadPhi{"enableTanThetadPhi", false, "flag to turn cut opening angle in delta theta delta phi on/off"}; + Configurable minTanThetadPhi{"minTanThetadPhi", 4., "apply min opening angle in delta theta delta phi to cut on late conversion"}; + Configurable cfgEnableQA{"cfgEnableQA", false, "flag to turn QA plots on/off"}; + ConfigurableAxis thConfigAxisTanThetaPhi{"thConfigAxisTanThetaPhi", {180, -90.f, 90.f}, ""}; + } mesonConfig; + + struct : ConfigurableGroup { + std::string prefix = "mixingConfig"; + ConfigurableAxis cfgVtxBins{"cfgVtxBins", {VARIABLE_WIDTH, -10.0f, -9.0f, -8.f, -7.0f, -6.f, -5.0f, -4.f, -3.f, -2.f, -1.f, 0.f, 1.f, 2.f, 3.f, 4.f, 5.f, 6.f, 7.f, 8.f, 9.f, 10.f}, "Mixing bins - z-vertex"}; + ConfigurableAxis cfgCentBins{"cfgCentBins", {VARIABLE_WIDTH, 0.0f, 5.0f, 10.0f, 20.0f, 30.0f, 40.0f, 50.0f, 60.0f, 70.0f, 80.0f, 90.0f, 100.f}, "Mixing bins - centrality"}; + Configurable cfgMixingDepth{"cfgMixingDepth", 2, "Mixing depth"}; + } mixingConfig; + + SliceCache cache; + o2::framework::Service ccdb; + int runNow = 0; + int runBefore = -1; + + // Filter clusterFilter = aod::skimmedcluster::time >= emccuts.cfgEMCminTime && aod::skimmedcluster::time <= emccuts.cfgEMCmaxTime && aod::skimmedcluster::m02 >= emccuts.cfgEMCminM02 && aod::skimmedcluster::m02 <= emccuts.cfgEMCmaxM02 && skimmedcluster::e >= emccuts.cfgEMCminE; + // Filter collisionFilter = (nabs(aod::collision::posZ) <= eventcuts.cfgZvtxMax) && (aod::evsel::ft0cOccupancyInTimeRange <= eventcuts.cfgFT0COccupancyMax) && (aod::evsel::ft0cOccupancyInTimeRange >= eventcuts.cfgFT0COccupancyMin); + // using FilteredEMCalPhotons = soa::Filtered>; + using EMCalPhotons = soa::Join; + using PCMPhotons = soa::Join; + using Colls = soa::Join; + + // for event mixing + using MyEMH = o2::aod::pwgem::dilepton::utils::EventMixingHandler, std::pair, o2::aod::pwgem::dilepton::utils::EMTrack>; + MyEMH* emh1 = nullptr; + MyEMH* emh2 = nullptr; + + Preslice perCollisionEMC = aod::emccluster::emeventId; + Preslice perCollisionPCM = aod::v0photonkf::emeventId; + + using BinningType = ColumnBinningPolicy; + BinningType binningOnPositions{{mixingConfig.cfgVtxBins, mixingConfig.cfgCentBins}, true}; + Pair pairPCMEMC{binningOnPositions, mixingConfig.cfgMixingDepth, -1, &cache}; // indicates that mixingConfig.cfgMixingDepth events should be mixed and under/overflow (-1) to be ignored + + HistogramRegistry registry{"registry", {}, OutputObjHandlingPolicy::AnalysisObject, false, false}; + + o2::emcal::Geometry* emcalGeom; + o2::emcal::BadChannelMap* mBadChannels; + // Constants for eta and phi ranges + double etaMin = -0.75, etaMax = 0.75; + int nBinsEta = 150; // 150 bins for eta + + double phiMin = 1.35, phiMax = 5.75; + int nBinsPhi = 440; // (440 bins = 0.01 step size covering most regions) + + std::vector lookupTable1D; + float epsilon = 1.e-8; + + // static constexpr + static constexpr int64_t NMinPhotonRotBkg = 3; + static constexpr int64_t NMinPhotonRotBkgMixed = 2; + + // Usage when cfgEnableNonLin is enabled + std::unique_ptr fEMCalCorrectionFactor; // ("fEMCalCorrectionFactor","(1 + [0]/x + [1]/x^2) / (1 + [2]/x)", 0.3, 100.); + + // To access the 1D array + inline int getIndex(int iEta, int iPhi) + { + return iEta * nBinsPhi + iPhi; + } + + // Function to access the lookup table + inline int8_t checkEtaPhi1D(double eta, double phi) + { + if (eta < etaMin || eta > etaMax || phi < phiMin || phi > phiMax) { + return 3; // Out of bounds + } + + // Compute indices directly + int iEta = static_cast((eta - etaMin) / ((etaMax - etaMin) / nBinsEta)); + int iPhi = static_cast((phi - phiMin) / ((phiMax - phiMin) / nBinsPhi)); + + return lookupTable1D[getIndex(iEta, iPhi)]; + } + + void defineEMEventCut() + { + fEMEventCut = EMPhotonEventCut("fEMEventCut", "fEMEventCut"); + fEMEventCut.SetRequireSel8(eventcuts.cfgRequireSel8); + fEMEventCut.SetRequireFT0AND(eventcuts.cfgRequireFT0AND); + fEMEventCut.SetZvtxRange(-eventcuts.cfgZvtxMax, +eventcuts.cfgZvtxMax); + fEMEventCut.SetRequireNoTFB(eventcuts.cfgRequireNoTFB); + fEMEventCut.SetRequireNoITSROFB(eventcuts.cfgRequireNoITSROFB); + fEMEventCut.SetRequireNoSameBunchPileup(eventcuts.cfgRequireNoSameBunchPileup); + fEMEventCut.SetRequireVertexITSTPC(eventcuts.cfgRequireVertexITSTPC); + fEMEventCut.SetRequireGoodZvtxFT0vsPV(eventcuts.cfgRequireGoodZvtxFT0vsPV); + fEMEventCut.SetRequireEMCReadoutInMB(eventcuts.cfgRequireEMCReadoutInMB); + fEMEventCut.SetRequireEMCHardwareTriggered(eventcuts.cfgRequireEMCHardwareTriggered); + } + + void defineEMCCut() + { + fEMCCut = EMCPhotonCut("fEMCCut", "fEMCCut"); + const float a = emccuts.cfgEMCTMEta->at(0); + const float b = emccuts.cfgEMCTMEta->at(1); + const float c = emccuts.cfgEMCTMEta->at(2); + + const float d = emccuts.cfgEMCTMPhi->at(0); + const float e = emccuts.cfgEMCTMPhi->at(1); + const float f = emccuts.cfgEMCTMPhi->at(2); + LOGF(info, "EMCal track matching parameters : a = %f, b = %f, c = %f, d = %f, e = %f, f = %f", a, b, c, d, e, f); + fEMCCut.SetTrackMatchingEta([a, b, c](float pT) { return a + std::pow(pT + b, c); }); + fEMCCut.SetTrackMatchingPhi([d, e, f](float pT) { return d + std::pow(pT + e, f); }); + fEMCCut.SetMinEoverP(emccuts.cfgEMCEoverp); + + fEMCCut.SetMinE(emccuts.cfgEMCminE); + fEMCCut.SetMinNCell(emccuts.cfgEMCminNCell); + fEMCCut.SetM02Range(emccuts.cfgEMCminM02, emccuts.cfgEMCmaxM02); + fEMCCut.SetTimeRange(emccuts.cfgEMCminTime, emccuts.cfgEMCmaxTime); + fEMCCut.SetUseExoticCut(emccuts.cfgEMCUseExoticCut); + fEMCCut.SetClusterizer(emccuts.clusterDefinition); + } + + void DefinePCMCut() + { + fPCMPhotonCut = V0PhotonCut("fPCMPhotonCut", "fPCMPhotonCut"); + + // for v0 + fPCMPhotonCut.SetV0PtRange(pcmcuts.cfgMinPtV0, pcmcuts.cfgMaxPtV0); + fPCMPhotonCut.SetV0EtaRange(pcmcuts.cfgMinEtaV0, pcmcuts.cfgMaxEtaV0); + fPCMPhotonCut.SetMinCosPA(pcmcuts.cfgMinCosPA); + fPCMPhotonCut.SetMaxPCA(pcmcuts.cfgMaxPCA); + fPCMPhotonCut.SetMaxChi2KF(pcmcuts.cfgMaxChi2KF); + fPCMPhotonCut.SetRxyRange(pcmcuts.cfgMinV0Radius, pcmcuts.cfgMaxV0Radius); + fPCMPhotonCut.SetAPRange(pcmcuts.cfgMaxAlphaAP, pcmcuts.cfgMaxQtAP); + fPCMPhotonCut.RejectITSib(pcmcuts.cfgRejectV0OnITSib); + + // for track + fPCMPhotonCut.SetMinNClustersTPC(pcmcuts.cfgMinNClusterTPC); + fPCMPhotonCut.SetMinNCrossedRowsTPC(pcmcuts.cfgMinNCrossedRows); + fPCMPhotonCut.SetMinNCrossedRowsOverFindableClustersTPC(0.8); + fPCMPhotonCut.SetMaxFracSharedClustersTPC(pcmcuts.cfgMaxFracSharedClusterTPC); + fPCMPhotonCut.SetChi2PerClusterTPC(0.0, pcmcuts.cfgMaxChi2TPC); + fPCMPhotonCut.SetTPCNsigmaElRange(pcmcuts.cfgMinTPCNSigmaEl, pcmcuts.cfgMaxTPCNSigmaEl); + fPCMPhotonCut.SetChi2PerClusterITS(-1e+10, pcmcuts.cfgMaxChi2ITS); + fPCMPhotonCut.SetNClustersITS(0, 7); + fPCMPhotonCut.SetMeanClusterSizeITSob(0.0, 16.0); + fPCMPhotonCut.SetDisableITSonly(pcmcuts.cfgDisableITSOnly); + fPCMPhotonCut.SetDisableTPConly(pcmcuts.cfgDisableTPCOnly); + fPCMPhotonCut.SetRequireITSTPC(pcmcuts.cfgRequireV0WithITSTPC); + fPCMPhotonCut.SetRequireITSonly(pcmcuts.cfgRequireV0WithITSonly); + fPCMPhotonCut.SetRequireTPConly(pcmcuts.cfgRequireV0WithTPConly); + } + + void init(InitContext&) + { + defineEMEventCut(); + defineEMCCut(); + fEMCCut.SetUseTM(emccuts.cfgEMCUseTM); // disables TM + o2::aod::pwgem::photonmeson::utils::eventhistogram::addEventHistograms(®istry); + + const AxisSpec thnAxisInvMass{thnConfigAxisInvMass, "#it{M}_{#gamma#gamma} (GeV/#it{c}^{2})"}; + const AxisSpec thnAxisPt{thnConfigAxisPt, "#it{p}_{T} (GeV/#it{c})"}; + const AxisSpec thnAxisCent{thnConfigAxisCent, "Centrality (%)"}; + const AxisSpec thAxisTanThetaPhi{mesonConfig.thConfigAxisTanThetaPhi, "atan(#Delta#theta/#Delta#varphi)"}; + const AxisSpec thAxisClusterEnergy{thnConfigAxisPt, "#it{E} (GeV)"}; + const AxisSpec thAxisEnergyCalib{thnConfigAxisEnergyCalib, "#it{E}_{clus} (GeV)"}; + const AxisSpec thAxisAlpha{100, -1., +1, "#alpha"}; + const AxisSpec thAxisEnergy{1000, 0., 100., "#it{E}_{clus} (GeV)"}; + + registry.add("hSparseCalibSE", "THn for Calib same event", HistType::kTHnSparseF, {thnAxisInvMass, thAxisEnergyCalib, thnAxisCent}); + registry.add("hSparseCalibBack", "THn for Calib background", HistType::kTHnSparseF, {thnAxisInvMass, thAxisEnergyCalib, thnAxisCent}); + + auto hClusterCuts = registry.add("hClusterCuts", "hClusterCuts;;Counts", kTH1D, {{6, 0.5, 6.5}}, false); + hClusterCuts->GetXaxis()->SetBinLabel(1, "in"); + hClusterCuts->GetXaxis()->SetBinLabel(2, "opening angle"); + hClusterCuts->GetXaxis()->SetBinLabel(3, "#it{M}_{#gamma#gamma}"); + hClusterCuts->GetXaxis()->SetBinLabel(4, "#it{p}_{T}"); + hClusterCuts->GetXaxis()->SetBinLabel(5, "conversion cut"); + hClusterCuts->GetXaxis()->SetBinLabel(6, "out"); + + auto hClusterCutsMixed = registry.add("hClusterCutsMixed", "hClusterCutsMixed;;Counts", kTH1D, {{6, 0.5, 6.5}}, false); + hClusterCutsMixed->GetXaxis()->SetBinLabel(1, "in"); + hClusterCutsMixed->GetXaxis()->SetBinLabel(2, "opening angle"); + hClusterCutsMixed->GetXaxis()->SetBinLabel(3, "#it{M}_{#gamma#gamma}"); + hClusterCutsMixed->GetXaxis()->SetBinLabel(4, "#it{p}_{T}"); + hClusterCutsMixed->GetXaxis()->SetBinLabel(5, "conversion cut"); + hClusterCutsMixed->GetXaxis()->SetBinLabel(6, "out"); + + if (eventcuts.cfgEnableQA) { + auto hCollisionEMCCheck = registry.add("hCollisionEMCCheck", "collision counter;;Counts", kTH1D, {{7, 0.5, 7.5}}, false); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(1, "all"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(2, "EMC MB Readout"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(3, "has clusters"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(4, "EMC MB Readout & has clusters"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(5, "EMC MB Readout but no clusters"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(6, "No EMC MB Readout but has clusters"); + hCollisionEMCCheck->GetXaxis()->SetBinLabel(7, "No EMC MB Readout and no clusters"); + } + + if (emccuts.cfgEnableQA) { + registry.add("hEClusterBefore", "Histo for cluster energy before cuts", HistType::kTH1D, {thAxisClusterEnergy}); + registry.add("hEClusterAfter", "Histo for cluster energy after cuts", HistType::kTH1D, {thAxisClusterEnergy}); + } + + if (mesonConfig.cfgEnableQA) { + registry.add("hInvMassPt", "Histo for inv pair mass vs pt", HistType::kTH2D, {thnAxisInvMass, thnAxisPt}); + registry.add("hTanThetaPhi", "Histo for identification of conversion cluster", HistType::kTH2D, {thnAxisInvMass, thAxisTanThetaPhi}); + registry.add("hAlphaPt", "Histo of meson asymmetry vs pT", HistType::kTH2D, {thAxisAlpha, thnAxisPt}); + registry.add("hInvMassPtMixed", "Histo for inv pair mass vs pt for mixed event", HistType::kTH2D, {thnAxisInvMass, thnAxisPt}); + registry.add("hTanThetaPhiMixed", "Histo for identification of conversion cluster for mixed event", HistType::kTH2D, {thnAxisInvMass, thAxisTanThetaPhi}); + registry.add("hAlphaPtMixed", "Histo of meson asymmetry vs pT for mixed event", HistType::kTH2D, {thAxisAlpha, thnAxisPt}); + } + + ccdb->setURL(ccdbUrl); + ccdb->setCaching(true); + ccdb->setLocalObjectValidityChecking(); + ccdb->setFatalWhenNull(false); + + LOG(info) << "thnConfigAxisInvMass.value[1] = " << thnConfigAxisInvMass.value[1] << " thnConfigAxisInvMass.value.back() = " << thnConfigAxisInvMass.value.back(); + LOG(info) << "thnConfigAxisPt.value[1] = " << thnConfigAxisPt.value[1] << " thnConfigAxisPt.value.back() = " << thnConfigAxisPt.value.back(); + + fEMCalCorrectionFactor = std::make_unique("fEMCalCorrectionFactor", "(1 + [0]/x + [1]/x^2) / (1 + [2]/x)", 0.3, 100.); + fEMCalCorrectionFactor->SetParameters(-5.33426e-01, 1.40144e-02, -5.24434e-01); + }; // end init + + /// Change radians to degree + /// \param angle in radians + /// \return angle in degree + float getAngleDegree(float angle) + { + return angle * 180.f * std::numbers::inv_pi_v; + } + + /// Get the centrality + /// \param collision is the collision with the centrality information + template + float getCentrality(TCollision const& collision) + { + float cent = -999.; + switch (centEstimator) { + case CentralityEstimator::CFT0M: + cent = collision.centFT0M(); + break; + case CentralityEstimator::CFT0A: + cent = collision.centFT0A(); + break; + case CentralityEstimator::CFT0C: + cent = collision.centFT0C(); + break; + default: + LOG(warning) << "Centrality estimator not valid. Possible values are T0M, T0A, T0C. Fallback to T0C"; + cent = collision.centFT0C(); + break; + } + return cent; + } + + bool isTooCloseToEdge(const int cellID, const int DistanceToBorder = 1) + { + if (DistanceToBorder <= 0) { + return false; + } + if (cellID < 0) { + return true; + } + + int iBadCell = -1; + + // check distance to border in case the cell is okay + auto [iSupMod, iMod, iPhi, iEta] = emcalGeom->GetCellIndex(cellID); + auto [irow, icol] = emcalGeom->GetCellPhiEtaIndexInSModule(iSupMod, iMod, iPhi, iEta); + + // Check rows/phi + int iRowLast = 24; + if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::EMCAL_HALF) { + iRowLast /= 2; // 2/3 sm case + } else if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::EMCAL_THIRD) { + iRowLast /= 3; // 1/3 sm case + } else if (emcalGeom->GetSMType(iSupMod) == o2::emcal::EMCALSMType::DCAL_EXT) { + iRowLast /= 3; // 1/3 sm case + } + + if (irow < DistanceToBorder || (iRowLast - irow) <= DistanceToBorder) { + iBadCell = 1; + } + + if (iBadCell > 0) { + return true; + } + return false; + } + + bool isCellMasked(int cellID) + { + bool masked = false; + if (mBadChannels) { + auto maskStatus = mBadChannels->getChannelStatus(cellID); + masked = (maskStatus != o2::emcal::BadChannelMap::MaskType_t::GOOD_CELL); + } + return masked; + } + + template + void initCCDB(TCollision const& collision) + { + // Load EMCal geometry + emcalGeom = o2::emcal::Geometry::GetInstanceFromRunNumber(collision.runNumber()); + // Load Bad Channel map + mBadChannels = ccdb->getForTimeStamp("EMC/Calib/BadChannelMap", collision.timestamp()); + lookupTable1D = std::vector(nBinsEta * nBinsPhi, -1); + double binWidthEta = (etaMax - etaMin) / nBinsEta; + double binWidthPhi = (phiMax - phiMin) / nBinsPhi; + + for (int iEta = 0; iEta < nBinsEta; ++iEta) { + double etaCenter = etaMin + (iEta + 0.5) * binWidthEta; + for (int iPhi = 0; iPhi < nBinsPhi; ++iPhi) { + double phiCenter = phiMin + (iPhi + 0.5) * binWidthPhi; + try { + // Get the cell ID + int cellID = emcalGeom->GetAbsCellIdFromEtaPhi(etaCenter, phiCenter); + + // Check conditions for the cell + if (isTooCloseToEdge(cellID, 1)) { + lookupTable1D[getIndex(iEta, iPhi)] = 2; // Edge + } else if (isCellMasked(cellID)) { + lookupTable1D[getIndex(iEta, iPhi)] = 1; // Bad + } else { + lookupTable1D[getIndex(iEta, iPhi)] = 0; // Good + } + } catch (o2::emcal::InvalidPositionException& e) { + lookupTable1D[getIndex(iEta, iPhi)] = 3; // Outside geometry + } + } + } + } + + /// \brief Calculate background using rotation background method for calib + template + void rotationBackgroundCalib(const ROOT::Math::PtEtaPhiMVector& meson, ROOT::Math::PtEtaPhiMVector photon1, ROOT::Math::PtEtaPhiMVector photon2, TPhotons const& photons_coll, unsigned int ig1, unsigned int ig2, TCollision const& collision) + { + // if less than 3 clusters are present skip event since we need at least 3 clusters + if (photons_coll.size() < NMinPhotonRotBkg) { + return; + } + float cent = getCentrality(collision); + int iCellIDPhoton1 = 0; + int iCellIDPhoton2 = 0; + + ROOT::Math::AxisAngle rotationAxis(meson.Vect(), cfgRotAngle.value); + ROOT::Math::Rotation3D rotationMatrix(rotationAxis); + photon1 = rotationMatrix * photon1; + photon2 = rotationMatrix * photon2; + + if (checkEtaPhi1D(photon1.Eta(), RecoDecay::constrainAngle(photon1.Phi())) >= cfgEMCalMapLevelBackground.value) { + iCellIDPhoton1 = -1; + } + if (checkEtaPhi1D(photon2.Eta(), RecoDecay::constrainAngle(photon2.Phi())) >= cfgEMCalMapLevelBackground.value) { + iCellIDPhoton2 = -1; + } + + if (iCellIDPhoton1 == -1 && iCellIDPhoton2 == -1) { + return; + } + for (const auto& photon : photons_coll) { + if (photon.globalIndex() == ig1 || photon.globalIndex() == ig2) { + // only combine rotated photons with other photons + continue; + } + if (!(fEMCCut.IsSelected(photon))) { + continue; + } + if (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelBackground.value) { + continue; + } + float energyCorrectionFactor = 1.f; + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(photon.e() > MinEnergy ? photon.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector photon3(energyCorrectionFactor * photon.pt(), photon.eta(), photon.phi(), 0.); + if (iCellIDPhoton1 >= 0) { + if (std::fabs((photon1.E() - photon3.E()) / (photon1.E() + photon3.E()) < cfgMaxAsymmetry)) { // only use symmetric decays + ROOT::Math::PtEtaPhiMVector mother1 = photon1 + photon3; + float openingAngle1 = std::acos(photon1.Vect().Dot(photon3.Vect()) / (photon1.P() * photon3.P())); + + if (openingAngle1 > mesonConfig.minOpenAngle && thnConfigAxisInvMass.value[1] <= mother1.M() && thnConfigAxisInvMass.value.back() >= mother1.M() && thnConfigAxisPt.value[1] <= mother1.Pt() && thnConfigAxisPt.value.back() >= mother1.Pt()) { + if (mesonConfig.enableTanThetadPhi) { + float dTheta = photon1.Theta() - photon3.Theta(); + float dPhi = photon1.Phi() - photon3.Phi(); + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hSparseCalibBack"), mother1.M(), mother1.E() / 2., cent); + } + } else { + registry.fill(HIST("hSparseCalibBack"), mother1.M(), mother1.E() / 2., cent); + } + } + } + } + if (iCellIDPhoton2 >= 0) { + if (std::fabs((photon2.E() - photon3.E()) / (photon2.E() + photon3.E()) < cfgMaxAsymmetry)) { // only use symmetric decays + ROOT::Math::PtEtaPhiMVector mother2 = photon2 + photon3; + float openingAngle2 = std::acos(photon2.Vect().Dot(photon3.Vect()) / (photon2.P() * photon3.P())); + + if (openingAngle2 > mesonConfig.minOpenAngle && thnConfigAxisInvMass.value[1] <= mother2.M() && thnConfigAxisInvMass.value.back() >= mother2.M() && thnConfigAxisPt.value[1] <= mother2.Pt() && thnConfigAxisPt.value.back() >= mother2.Pt()) { + if (mesonConfig.enableTanThetadPhi) { + float dTheta = photon2.Theta() - photon3.Theta(); + float dPhi = photon2.Phi() - photon3.Phi(); + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hSparseCalibBack"), mother2.M(), mother2.E() / 2., cent); + } + } else { + registry.fill(HIST("hSparseCalibBack"), mother2.M(), mother2.E() / 2., cent); + } + } + } + } + } // end of loop over third photon + return; + } + + /// \brief Calculate background using rotation background method for calib + /// \param meson mother particle from photon1 & photon2 which defines rotation axis + /// \param photon1 first photon (EMC) which will be rotated + /// \param photon2 second photon (PCM) which will be rotated + /// \param photonsEMC sub table of EMC photons of current event which will be combined with rotated photon2 + /// \param photonsPCM sub table of PCM photons of current event which will be combined with rotated photon1 + /// \param ig1 index of photon1 + /// \param ig2 index of photon2 + /// \param cent current collisions centrality + template + void rotationBackgroundCalibEMCPCM(const ROOT::Math::PtEtaPhiMVector& meson, ROOT::Math::PtEtaPhiMVector photon1, ROOT::Math::PtEtaPhiMVector photon2, TPhotonEMC const& photonsEMC, TPhotonPCM const& photonsPCM, unsigned int ig1, unsigned int ig2, float cent) + { + // we need at least 2 clusters or 2 pcm photons for rotation + if (photonsEMC.size() < NMinPhotonRotBkgMixed || photonsPCM.size() < NMinPhotonRotBkgMixed) { + return; + } + int iCellIDPhoton1 = 0; + int iCellIDPhoton2 = 0; + + ROOT::Math::AxisAngle rotationAxis(meson.Vect(), cfgRotAngle.value); + ROOT::Math::Rotation3D rotationMatrix(rotationAxis); + photon1 = rotationMatrix * photon1; + photon2 = rotationMatrix * photon2; + + if (checkEtaPhi1D(photon1.Eta(), RecoDecay::constrainAngle(photon1.Phi())) >= cfgEMCalMapLevelBackground.value) { + iCellIDPhoton1 = -1; + } + + if (iCellIDPhoton1 == -1 && iCellIDPhoton2 == -1) { + return; + } + // Combining with EMCal photons from event + if (photonsEMC.size() >= NMinPhotonRotBkgMixed) { + for (const auto& photon : photonsEMC) { + if (photon.globalIndex() == ig1) { + continue; + } + if (!(fEMCCut.IsSelected(photon))) { + continue; + } + if (checkEtaPhi1D(photon.eta(), RecoDecay::constrainAngle(photon.phi())) >= cfgEMCalMapLevelBackground.value) { + continue; + } + float energyCorrectionFactor = 1.f; + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(photon.e() > MinEnergy ? photon.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector photon3(energyCorrectionFactor * photon.pt(), photon.eta(), photon.phi(), 0.); + if (iCellIDPhoton2 >= 0) { + ROOT::Math::PtEtaPhiMVector mother = photon2 + photon3; + float openingAngle2 = std::acos(photon2.Vect().Dot(photon3.Vect()) / (photon2.P() * photon3.P())); + + if (openingAngle2 > mesonConfig.minOpenAngle && thnConfigAxisInvMass.value[1] <= mother.M() && thnConfigAxisInvMass.value.back() >= mother.M() && thnConfigAxisPt.value[1] <= mother.Pt() && thnConfigAxisPt.value.back() >= mother.Pt()) { + if (mesonConfig.enableTanThetadPhi) { + float dTheta = photon2.Theta() - photon3.Theta(); + float dPhi = photon2.Phi() - photon3.Phi(); + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hSparseCalibBack"), mother.M(), mother.E() / 2., cent); + } + } else { + registry.fill(HIST("hSparseCalibBack"), mother.M(), mother.E() / 2., cent); + } + } + } + } // end of loop over third photon + } // check that we have at least 2 clusters + + // Combining with PCM photons from event + if (photonsPCM.size() >= NMinPhotonRotBkgMixed) { + for (const auto& photon : photonsPCM) { + if (photon.globalIndex() == ig2) { + continue; + } + if (!(fPCMPhotonCut.IsSelected(photon))) { + continue; + } + ROOT::Math::PtEtaPhiMVector photon3(photon.pt(), photon.eta(), photon.phi(), 0.); + if (iCellIDPhoton1 >= 0) { + ROOT::Math::PtEtaPhiMVector mother = photon1 + photon3; + float openingAngle2 = std::acos(photon1.Vect().Dot(photon3.Vect()) / (photon1.P() * photon3.P())); + + if (openingAngle2 > mesonConfig.minOpenAngle && thnConfigAxisInvMass.value[1] <= mother.M() && thnConfigAxisInvMass.value.back() >= mother.M() && thnConfigAxisPt.value[1] <= mother.Pt() && thnConfigAxisPt.value.back() >= mother.Pt()) { + if (mesonConfig.enableTanThetadPhi) { + float dTheta = photon1.Theta() - photon3.Theta(); + float dPhi = photon1.Phi() - photon3.Phi(); + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hSparseCalibBack"), mother.M(), mother.E() / 2., cent); + } + } else { + registry.fill(HIST("hSparseCalibBack"), mother.M(), mother.E() / 2., cent); + } + } + } + } // end of loop over third photon + } // check that we have at least 2 clusters + return; + } + + // EMCal calibration same event + void processEMCalCalib(Colls const& collisions, EMCalPhotons const& clusters, PCMPhotons const&) + { + float energyCorrectionFactor = 1.f; + int nColl = 1; + for (const auto& collision : collisions) { + auto photonsPerCollision = clusters.sliceBy(perCollisionEMC, collision.globalIndex()); + o2::aod::pwgem::photonmeson::utils::eventhistogram::fillEventInfo<0>(®istry, collision); + if (!(fEMEventCut.IsSelected(collision))) { + // general event selection + continue; + } + if (!(eventcuts.cfgFT0COccupancyMin <= collision.ft0cOccupancyInTimeRange() && collision.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) { + // occupancy selection + continue; + } + float cent = getCentrality(collision); + if (cent < eventcuts.cfgMinCent || cent > eventcuts.cfgMaxCent) { + // event selection + continue; + } + runNow = collision.runNumber(); + if (runNow != runBefore) { + initCCDB(collision); + runBefore = runNow; + } + o2::aod::pwgem::photonmeson::utils::eventhistogram::fillEventInfo<1>(®istry, collision); + registry.fill(HIST("Event/before/hCollisionCounter"), 12.0); // accepted + registry.fill(HIST("Event/after/hCollisionCounter"), 12.0); // accepted + + if (emccuts.cfgEnableQA) { + for (const auto& photon : photonsPerCollision) { + registry.fill(HIST("hEClusterBefore"), photon.e()); // before cuts + if (!(fEMCCut.IsSelected(photon))) { + continue; + } + registry.fill(HIST("hEClusterAfter"), photon.e()); // accepted after cuts + } + } + for (const auto& [g1, g2] : combinations(CombinationsStrictlyUpperIndexPolicy(photonsPerCollision, photonsPerCollision))) { + if (!(fEMCCut.IsSelected(g1)) || !(fEMCCut.IsSelected(g2))) { + continue; + } + + // Cut edge clusters away, similar to rotation method to ensure same acceptance is used + if (cfgDistanceToEdge.value) { + if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelSameEvent.value) { + continue; + } + if (checkEtaPhi1D(g2.eta(), RecoDecay::constrainAngle(g2.phi())) >= cfgEMCalMapLevelSameEvent.value) { + continue; + } + } + + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.); + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g2.e() > MinEnergy ? g2.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v2(energyCorrectionFactor * g2.pt(), g2.eta(), g2.phi(), 0.); + ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2; + + float dTheta = v1.Theta() - v2.Theta(); + float dPhi = v1.Phi() - v2.Phi(); + float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); + registry.fill(HIST("hClusterCuts"), 1); + if (openingAngle <= mesonConfig.minOpenAngle) { + registry.fill(HIST("hClusterCuts"), 2); + continue; + } + if (cfgDoRotation) { + if (nColl % cfgDownsampling.value == 0) { + rotationBackgroundCalib(vMeson, v1, v2, photonsPerCollision, g1.globalIndex(), g2.globalIndex(), collision); + } + } + if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) { + registry.fill(HIST("hClusterCuts"), 3); + continue; + } + if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) { + registry.fill(HIST("hClusterCuts"), 4); + continue; + } + if (mesonConfig.cfgEnableQA) { + registry.fill(HIST("hInvMassPt"), vMeson.M(), vMeson.Pt()); + registry.fill(HIST("hTanThetaPhi"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi))); + registry.fill(HIST("hAlphaPt"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt()); + } + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hClusterCuts"), 5); + continue; + } + if (std::fabs((v1.E() - v2.E()) / (v1.E() + v2.E()) < cfgMaxAsymmetry)) { // only use symmetric decays + registry.fill(HIST("hClusterCuts"), 6); + registry.fill(HIST("hSparseCalibSE"), vMeson.M(), vMeson.E() / 2., getCentrality(collision)); + } + } + if (cfgDoRotation) { + if (nColl % cfgDownsampling.value == 0) { + nColl = 1; // reset counter + } else { + nColl++; + } + } + } + } + PROCESS_SWITCH(CalibTaskEmc, processEMCalCalib, "Process EMCal calibration same event", true); + + // EMCal calibration + void processEMCalPCMCalib(Colls const& collisions, EMCalPhotons const& clusters, PCMPhotons const& photons, aod::V0Legs const&) + { + float energyCorrectionFactor = 1.f; + int nColl = 1; + for (const auto& collision : collisions) { + o2::aod::pwgem::photonmeson::utils::eventhistogram::fillEventInfo<0>(®istry, collision); + if (!(fEMEventCut.IsSelected(collision))) { + // general event selection + continue; + } + if (!(eventcuts.cfgFT0COccupancyMin <= collision.ft0cOccupancyInTimeRange() && collision.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) { + // occupancy selection + continue; + } + float cent = getCentrality(collision); + if (cent < eventcuts.cfgMinCent || cent > eventcuts.cfgMaxCent) { + // event selection + continue; + } + runNow = collision.runNumber(); + if (runNow != runBefore) { + initCCDB(collision); + runBefore = runNow; + } + o2::aod::pwgem::photonmeson::utils::eventhistogram::fillEventInfo<1>(®istry, collision); + registry.fill(HIST("Event/before/hCollisionCounter"), 12.0); // accepted + registry.fill(HIST("Event/after/hCollisionCounter"), 12.0); // accepted + + auto photonsEMCPerCollision = clusters.sliceBy(perCollisionEMC, collision.globalIndex()); + auto photonsPCMPerCollision = photons.sliceBy(perCollisionPCM, collision.globalIndex()); + + if (emccuts.cfgEnableQA) { + for (const auto& photon : photonsEMCPerCollision) { + registry.fill(HIST("hEClusterBefore"), photon.e()); // before cuts + if (!(fEMCCut.IsSelected(photon))) { + continue; + } + registry.fill(HIST("hEClusterAfter"), photon.e()); // accepted after cuts + } + } + for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photonsEMCPerCollision, photonsPCMPerCollision))) { + if (!(fEMCCut.IsSelected(g1)) || !(fPCMPhotonCut.IsSelected(g2))) { + continue; + } + + // Cut edge clusters away, similar to rotation method to ensure same acceptance is used + if (cfgDistanceToEdge.value) { + if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelSameEvent.value) { + continue; + } + } + + // EMCal photon v1 + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.); + // PCM photon v2s + ROOT::Math::PtEtaPhiMVector v2(g2.pt(), g2.eta(), g2.phi(), 0.); + ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2; + + float dTheta = v1.Theta() - v2.Theta(); + float dPhi = v1.Phi() - v2.Phi(); + float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); + registry.fill(HIST("hClusterCuts"), 1); + if (openingAngle <= mesonConfig.minOpenAngle) { + registry.fill(HIST("hClusterCuts"), 2); + continue; + } + if (cfgDoRotation) { + if (nColl % cfgDownsampling.value == 0) { + rotationBackgroundCalibEMCPCM(vMeson, v1, v2, photonsEMCPerCollision, photonsPCMPerCollision, g1.globalIndex(), g2.globalIndex(), cent); + } + } + if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) { + registry.fill(HIST("hClusterCuts"), 3); + continue; + } + if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) { + registry.fill(HIST("hClusterCuts"), 4); + continue; + } + if (mesonConfig.cfgEnableQA) { + registry.fill(HIST("hInvMassPt"), vMeson.M(), vMeson.Pt()); + registry.fill(HIST("hTanThetaPhi"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi))); + registry.fill(HIST("hAlphaPt"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt()); + } + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hClusterCuts"), 5); + continue; + } + registry.fill(HIST("hSparseCalibSE"), vMeson.M(), vMeson.E() / 2., getCentrality(collision)); + } + if (cfgDoRotation) { + if (nColl % cfgDownsampling.value == 0) { + nColl = 1; // reset counter + } else { + nColl++; + } + } + } + } + PROCESS_SWITCH(CalibTaskEmc, processEMCalPCMCalib, "Process EMCal calibration using PCM-EMC same event", true); + + // EMCal calibration mixed event + void processEMCalCalibMixed(Colls const&, EMCalPhotons const&, PCMPhotons const&) + { + float energyCorrectionFactor = 1.f; + + SameKindPair pair{binningOnPositions, mixingConfig.cfgMixingDepth, -1, &cache}; // indicates that 5 events should be mixed and under/overflow (-1) to be ignored + + for (const auto& [c1, clusters1, c2, clusters2] : pair) { + if (!(fEMEventCut.IsSelected(c1)) || !(fEMEventCut.IsSelected(c2))) { + // general event selection + continue; + } + if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) { + // occupancy selection + continue; + } + if (getCentrality(c1) < eventcuts.cfgMinCent || getCentrality(c1) > eventcuts.cfgMaxCent || getCentrality(c2) < eventcuts.cfgMinCent || getCentrality(c2) > eventcuts.cfgMaxCent) { + // event selection + continue; + } + runNow = c1.runNumber(); + if (runNow != runBefore) { + initCCDB(c1); + runBefore = runNow; + } + for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(clusters1, clusters2))) { + if (!(fEMCCut.IsSelected(g1)) || !(fEMCCut.IsSelected(g2))) { + continue; + } + // Cut edge clusters away, similar to rotation method to ensure same acceptance is used + if (cfgDistanceToEdge.value) { + if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelBackground.value) { + continue; + } + if (checkEtaPhi1D(g2.eta(), RecoDecay::constrainAngle(g2.phi())) >= cfgEMCalMapLevelBackground.value) { + continue; + } + } + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.); + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g2.e() > MinEnergy ? g2.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v2(energyCorrectionFactor * g2.pt(), g2.eta(), g2.phi(), 0.); + ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2; + + float dTheta = v1.Theta() - v2.Theta(); + float dPhi = v1.Phi() - v2.Phi(); + float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); + + registry.fill(HIST("hClusterCutsMixed"), 1); + if (openingAngle <= mesonConfig.minOpenAngle) { + registry.fill(HIST("hClusterCutsMixed"), 2); + continue; + } + if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) { + registry.fill(HIST("hClusterCutsMixed"), 3); + continue; + } + if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) { + registry.fill(HIST("hClusterCutsMixed"), 4); + continue; + } + if (mesonConfig.cfgEnableQA) { + registry.fill(HIST("hInvMassPtMixed"), vMeson.M(), vMeson.Pt()); + registry.fill(HIST("hTanThetaPhiMixed"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi))); + registry.fill(HIST("hAlphaPtMixed"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt()); + } + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hClusterCutsMixed"), 5); + continue; + } + if (std::fabs((v1.E() - v2.E()) / (v1.E() + v2.E()) < cfgMaxAsymmetry)) { // only use symmetric decays + registry.fill(HIST("hClusterCutsMixed"), 6); + registry.fill(HIST("hSparseCalibBack"), vMeson.M(), vMeson.E() / 2., getCentrality(c1)); + } + } + } + } + PROCESS_SWITCH(CalibTaskEmc, processEMCalCalibMixed, "Process EMCal calibration mixed event", false); + + // EMCal calibration + void processEMCalPCMCalibMixed(Colls const&, EMCalPhotons const&, PCMPhotons const&, aod::V0Legs const&) + { + float energyCorrectionFactor = 1.f; + + LOG(info) << "Beginning of processEMCalPCMCalibMixed"; + + for (const auto& [c1, photonEMC, c2, photonPCM] : pairPCMEMC) { + if (!(fEMEventCut.IsSelected(c1)) || !(fEMEventCut.IsSelected(c2))) { + // general event selection + continue; + } + if (!(eventcuts.cfgFT0COccupancyMin <= c1.ft0cOccupancyInTimeRange() && c1.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax) || !(eventcuts.cfgFT0COccupancyMin <= c2.ft0cOccupancyInTimeRange() && c2.ft0cOccupancyInTimeRange() < eventcuts.cfgFT0COccupancyMax)) { + // occupancy selection + continue; + } + if (getCentrality(c1) < eventcuts.cfgMinCent || getCentrality(c1) > eventcuts.cfgMaxCent || getCentrality(c2) < eventcuts.cfgMinCent || getCentrality(c2) > eventcuts.cfgMaxCent) { + // event selection + continue; + } + runNow = c1.runNumber(); + if (runNow != runBefore) { + initCCDB(c1); + runBefore = runNow; + } + for (const auto& [g1, g2] : combinations(CombinationsFullIndexPolicy(photonEMC, photonPCM))) { + if (!(fEMCCut.IsSelected(g1)) || !(fPCMPhotonCut.IsSelected(g2))) { + continue; + } + // Cut edge clusters away, similar to rotation method to ensure same acceptance is used + if (cfgDistanceToEdge.value) { + if (checkEtaPhi1D(g1.eta(), RecoDecay::constrainAngle(g1.phi())) >= cfgEMCalMapLevelBackground.value) { + continue; + } + } + energyCorrectionFactor = fEMCalCorrectionFactor->Eval(g1.e() > MinEnergy ? g1.e() : MinEnergy); + ROOT::Math::PtEtaPhiMVector v1(energyCorrectionFactor * g1.pt(), g1.eta(), g1.phi(), 0.); + ROOT::Math::PtEtaPhiMVector v2(g2.pt(), g2.eta(), g2.phi(), 0.); + ROOT::Math::PtEtaPhiMVector vMeson = v1 + v2; + + float dTheta = v1.Theta() - v2.Theta(); + float dPhi = v1.Phi() - v2.Phi(); + float openingAngle = std::acos(v1.Vect().Dot(v2.Vect()) / (v1.P() * v2.P())); + + registry.fill(HIST("hClusterCutsMixed"), 1); + if (openingAngle <= mesonConfig.minOpenAngle) { + registry.fill(HIST("hClusterCutsMixed"), 2); + continue; + } + if (thnConfigAxisInvMass.value[1] > vMeson.M() || thnConfigAxisInvMass.value.back() < vMeson.M()) { + registry.fill(HIST("hClusterCutsMixed"), 3); + continue; + } + if (thnConfigAxisPt.value[1] > vMeson.Pt() || thnConfigAxisPt.value.back() < vMeson.Pt()) { + registry.fill(HIST("hClusterCutsMixed"), 4); + continue; + } + if (mesonConfig.cfgEnableQA) { + registry.fill(HIST("hInvMassPtMixed"), vMeson.M(), vMeson.Pt()); + registry.fill(HIST("hTanThetaPhiMixed"), vMeson.M(), getAngleDegree(std::atan(dTheta / dPhi))); + registry.fill(HIST("hAlphaPtMixed"), (v1.E() - v2.E()) / (v1.E() + v2.E()), vMeson.Pt()); + } + if (mesonConfig.enableTanThetadPhi && mesonConfig.minTanThetadPhi > std::fabs(getAngleDegree(std::atan(dTheta / dPhi)))) { + registry.fill(HIST("hClusterCutsMixed"), 5); + continue; + } + registry.fill(HIST("hSparseCalibBack"), vMeson.M(), vMeson.E() / 2., getCentrality(c1)); + } + } + } + PROCESS_SWITCH(CalibTaskEmc, processEMCalPCMCalibMixed, "Process EMCal calibration with PCM-EMC mixed event", false); + +}; // End struct CalibTaskEmc + +WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +{ + return WorkflowSpec{adaptAnalysisTask(cfgc)}; +} diff --git a/PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx b/PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx index aa837b02233..4b576285531 100644 --- a/PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx +++ b/PWGEM/PhotonMeson/Tasks/taskPi0FlowEMC.cxx @@ -11,14 +11,12 @@ /// \file taskPi0FlowEMC.cxx /// \brief Analysis task for neutral pion flow with EMCal -/// /// \author M. Hemmer, marvin.hemmer@cern.ch #include "PWGEM/PhotonMeson/Core/EMCPhotonCut.h" #include "PWGEM/PhotonMeson/Core/EMPhotonEventCut.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/EventHistograms.h" -#include "PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h" #include "Common/CCDB/TriggerAliases.h" #include "Common/Core/EventPlaneHelper.h" @@ -50,6 +48,7 @@ #include #include #include // IWYU pragma: keep +#include #include #include #include @@ -209,7 +208,7 @@ struct TaskPi0FlowEMC { int runNow = 0; int runBefore = -1; - Filter clusterFilter = aod::skimmedcluster::time >= emccuts.cfgEMCminTime && aod::skimmedcluster::time <= emccuts.cfgEMCmaxTime && aod::skimmedcluster::m02 >= emccuts.cfgEMCminM02 && aod::skimmedcluster::m02 <= emccuts.cfgEMCmaxM02 && skimmedcluster::e >= emccuts.cfgEMCminE; + Filter clusterFilter = aod::skimmedcluster::time >= emccuts.cfgEMCminTime && aod::skimmedcluster::time <= emccuts.cfgEMCmaxTime && aod::skimmedcluster::m02 >= emccuts.cfgEMCminM02 && aod::skimmedcluster::m02 <= emccuts.cfgEMCmaxM02 && aod::skimmedcluster::e >= emccuts.cfgEMCminE; Filter collisionFilter = (nabs(aod::collision::posZ) <= eventcuts.cfgZvtxMax) && (aod::evsel::ft0cOccupancyInTimeRange <= eventcuts.cfgFT0COccupancyMax) && (aod::evsel::ft0cOccupancyInTimeRange >= eventcuts.cfgFT0COccupancyMin); using FilteredEMCalPhotons = soa::Filtered>; using EMCalPhotons = soa::Join; @@ -330,6 +329,7 @@ struct TaskPi0FlowEMC { const AxisSpec thAxisCN{8, 0.5, 8.5, "#it{c}_{n}"}; const AxisSpec thAxisSN{8, 0.5, 8.5, "#it{s}_{n}"}; const AxisSpec thAxisCPUTime{1000, 0, 10000, "#it{t} (#mus)"}; + const AxisSpec thAxisAzimuth{360, -o2::constants::math::PI, o2::constants::math::PI, "#it{#varphi} (rad)"}; const AxisSpec thnAxisMixingVtx{mixingConfig.cfgVtxBins, "#it{z} (cm)"}; const AxisSpec thnAxisMixingCent{mixingConfig.cfgCentBins, "Centrality (%)"}; @@ -399,6 +399,12 @@ struct TaskPi0FlowEMC { registry.add("epReso/hEpSinCoefficientsTPCneg", "hEpSinCoefficientsTPCneg; centrality; s_{n}", HistType::kTProfile2D, {thnAxisCent, thAxisSN}); registry.add("epReso/hEpCosCoefficientsTPCTots", "hEpCosCoefficientsTPCTots; centrality; c_{n}", HistType::kTProfile2D, {thnAxisCent, thAxisCN}); registry.add("epReso/hEpSinCoefficientsTPCTots", "hEpSinCoefficientsTPCTots; centrality; s_{n}", HistType::kTProfile2D, {thnAxisCent, thAxisSN}); + registry.add("QVector/hQVecMeanRVsPhiFT0a", "hQVecMeanRVsPhiFT0a; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); + registry.add("QVector/hQVecMeanRVsPhiFT0c", "hQVecMeanRVsPhiFT0c; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); + registry.add("QVector/hQVecMeanRVsPhiFT0m", "hQVecMeanRVsPhiFT0m; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); + registry.add("QVector/hQVecMeanRVsPhiTPCpos", "hQVecMeanRVsPhiTPCpos; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); + registry.add("QVector/hQVecMeanRVsPhiTPCneg", "hQVecMeanRVsPhiTPCneg; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); + registry.add("QVector/hQVecMeanRVsPhiTPCTot", "hQVecMeanRVsPhiTPCTot; centrality; #it{#varphi} (rad), <#it{r}> (a.u.)", HistType::kTProfile2D, {thnAxisCent, thAxisAzimuth}); } if (eventcuts.cfgEnableQA) { auto hCollisionEMCCheck = registry.add("hCollisionEMCCheck", "collision counter;;Counts", kTH1D, {{7, 0.5, 7.5}}, false); @@ -1239,6 +1245,13 @@ struct TaskPi0FlowEMC { registry.fill(HIST("hEventPlaneAngleTPCpos"), centrality, epBPoss); registry.fill(HIST("hEventPlaneAngleTPCneg"), centrality, epBNegs); + registry.fill(HIST("QVector/hQVecMeanRVsPhiFT0a"), centrality, std::atan2(yQVecFT0a, xQVecFT0a)); + registry.fill(HIST("QVector/hQVecMeanRVsPhiFT0c"), centrality, std::atan2(yQVecFT0c, xQVecFT0c)); + registry.fill(HIST("QVector/hQVecMeanRVsPhiFT0m"), centrality, std::atan2(yQVecFT0m, xQVecFT0m)); + registry.fill(HIST("QVector/hQVecMeanRVsPhiTPCpos"), centrality, std::atan2(yQVecBPos, xQVecBPos)); + registry.fill(HIST("QVector/hQVecMeanRVsPhiTPCneg"), centrality, std::atan2(yQVecBNeg, xQVecBNeg)); + registry.fill(HIST("QVector/hQVecMeanRVsPhiTPCTot"), centrality, std::atan2(yQVecBTot, xQVecBTot)); + registry.fill(HIST("epReso/hEpResoFT0cFT0a"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFT0c, epFT0a))); registry.fill(HIST("epReso/hEpResoFT0cTPCpos"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFT0c, epBPoss))); registry.fill(HIST("epReso/hEpResoFT0cTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFT0c, epBNegs))); @@ -1251,18 +1264,18 @@ struct TaskPi0FlowEMC { registry.fill(HIST("epReso/hEpResoFT0mTPCtot"), centrality, std::cos(harmonic * getDeltaPsiInRange(epFT0m, epBTots))); registry.fill(HIST("epReso/hEpResoTPCposTPCneg"), centrality, std::cos(harmonic * getDeltaPsiInRange(epBPoss, epBNegs))); for (int n = 1; n <= kOctagonal; n++) { - registry.fill(HIST("epReso/hEpCosCoefficientsFT0c"), centrality, n, std::cos(n * epFT0c)); - registry.fill(HIST("epReso/hEpSinCoefficientsFT0c"), centrality, n, std::sin(n * epFT0c)); - registry.fill(HIST("epReso/hEpCosCoefficientsFT0a"), centrality, n, std::cos(n * epFT0a)); - registry.fill(HIST("epReso/hEpSinCoefficientsFT0a"), centrality, n, std::sin(n * epFT0a)); - registry.fill(HIST("epReso/hEpCosCoefficientsFT0m"), centrality, n, std::cos(n * epFT0m)); - registry.fill(HIST("epReso/hEpSinCoefficientsFT0m"), centrality, n, std::sin(n * epFT0m)); - registry.fill(HIST("epReso/hEpCosCoefficientsTPCpos"), centrality, n, std::cos(n * epBPoss)); - registry.fill(HIST("epReso/hEpSinCoefficientsTPCpos"), centrality, n, std::sin(n * epBPoss)); - registry.fill(HIST("epReso/hEpCosCoefficientsTPCneg"), centrality, n, std::cos(n * epBNegs)); - registry.fill(HIST("epReso/hEpSinCoefficientsTPCneg"), centrality, n, std::sin(n * epBNegs)); - registry.fill(HIST("epReso/hEpCosCoefficientsTPCTots"), centrality, n, std::cos(n * epBTots)); - registry.fill(HIST("epReso/hEpSinCoefficientsTPCTots"), centrality, n, std::sin(n * epBTots)); + registry.fill(HIST("epReso/hEpCosCoefficientsFT0c"), centrality, n, std::cos(n * harmonic * epFT0c)); + registry.fill(HIST("epReso/hEpSinCoefficientsFT0c"), centrality, n, std::sin(n * harmonic * epFT0c)); + registry.fill(HIST("epReso/hEpCosCoefficientsFT0a"), centrality, n, std::cos(n * harmonic * epFT0a)); + registry.fill(HIST("epReso/hEpSinCoefficientsFT0a"), centrality, n, std::sin(n * harmonic * epFT0a)); + registry.fill(HIST("epReso/hEpCosCoefficientsFT0m"), centrality, n, std::cos(n * harmonic * epFT0m)); + registry.fill(HIST("epReso/hEpSinCoefficientsFT0m"), centrality, n, std::sin(n * harmonic * epFT0m)); + registry.fill(HIST("epReso/hEpCosCoefficientsTPCpos"), centrality, n, std::cos(n * harmonic * epBPoss)); + registry.fill(HIST("epReso/hEpSinCoefficientsTPCpos"), centrality, n, std::sin(n * harmonic * epBPoss)); + registry.fill(HIST("epReso/hEpCosCoefficientsTPCneg"), centrality, n, std::cos(n * harmonic * epBNegs)); + registry.fill(HIST("epReso/hEpSinCoefficientsTPCneg"), centrality, n, std::sin(n * harmonic * epBNegs)); + registry.fill(HIST("epReso/hEpCosCoefficientsTPCTots"), centrality, n, std::cos(n * harmonic * epBTots)); + registry.fill(HIST("epReso/hEpSinCoefficientsTPCTots"), centrality, n, std::sin(n * harmonic * epBTots)); } } } @@ -1379,7 +1392,7 @@ struct TaskPi0FlowEMC { } } } - PROCESS_SWITCH(TaskPi0FlowEMC, processEMCalCalib, "Process EMCal calibration", false); + PROCESS_SWITCH(TaskPi0FlowEMC, processEMCalCalib, "Process EMCal calibration mixed event", false); // Pi0 from EMCal void processM02(CollsWithQvecs const& collisions, EMCalPhotons const& clusters) diff --git a/PWGEM/PhotonMeson/Utils/ClusterHistograms.h b/PWGEM/PhotonMeson/Utils/ClusterHistograms.h index b5f0eb6db85..71337021f38 100644 --- a/PWGEM/PhotonMeson/Utils/ClusterHistograms.h +++ b/PWGEM/PhotonMeson/Utils/ClusterHistograms.h @@ -16,6 +16,7 @@ #ifndef PWGEM_PHOTONMESON_UTILS_CLUSTERHISTOGRAMS_H_ #define PWGEM_PHOTONMESON_UTILS_CLUSTERHISTOGRAMS_H_ +#include #include #include @@ -26,16 +27,14 @@ #include -using namespace o2::framework; - namespace o2::aod::pwgem::photonmeson::utils::clusterhistogram { -inline void addClusterHistograms(HistogramRegistry* fRegistry, bool do2DQA) +inline void addClusterHistograms(o2::framework::HistogramRegistry* fRegistry, bool do2DQA) { fRegistry->add("Cluster/before/hE", "E_{cluster};#it{E}_{cluster} (GeV);#it{N}_{cluster}", o2::framework::kTH1F, {{500, 0.0f, 50}}, true); fRegistry->add("Cluster/before/hPt", "Transverse momenta of clusters;#it{p}_{T} (GeV/c);#it{N}_{cluster}", o2::framework::kTH1F, {{500, 0.0f, 50}}, true); fRegistry->add("Cluster/before/hNgamma", "Number of #gamma candidates per collision;#it{N}_{#gamma} per collision;#it{N}_{collisions}", o2::framework::kTH1F, {{51, -0.5f, 50.5f}}, true); - fRegistry->add("Cluster/before/hEtaPhi", "#eta vs #varphi;#eta;#varphi (rad.)", o2::framework::kTH2F, {{280, -0.7f, 0.7f}, {180, 0, 2 * M_PI}}, true); + fRegistry->add("Cluster/before/hEtaPhi", "#eta vs #varphi;#eta;#varphi (rad.)", o2::framework::kTH2F, {{280, -0.7f, 0.7f}, {180, 0, o2::constants::math::TwoPI}}, true); fRegistry->add("Cluster/before/hNTracks", "Number of tracks considered for TM;#it{N}_{tracks};#it{N}_{cluster}", o2::framework::kTH1F, {{20, -0.5f, 19.5}}, true); fRegistry->add("Cluster/before/hTrackdEtadPhi", "d#eta vs. d#varphi of matched tracks;d#eta;d#varphi (rad.)", o2::framework::kTH2F, {{200, -0.2f, 0.2f}, {200, -0.2f, 0.2f}}, true); @@ -67,31 +66,31 @@ inline void addClusterHistograms(HistogramRegistry* fRegistry, bool do2DQA) } template -inline void fillClusterHistograms(HistogramRegistry* fRegistry, TCluster cluster, bool do2DQA, float weight = 1.f) +inline void fillClusterHistograms(o2::framework::HistogramRegistry* fRegistry, TCluster cluster, bool do2DQA, float weight = 1.f) { - static constexpr std::string_view cluster_types[2] = {"before/", "after/"}; - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hE"), cluster.e(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hPt"), cluster.pt(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hEtaPhi"), cluster.eta(), cluster.phi(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hNTracks"), cluster.deltaEta().size(), weight); + static constexpr std::string_view kClusterTypes[2] = {"before/", "after/"}; + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hE"), cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hPt"), cluster.pt(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hEtaPhi"), cluster.eta(), cluster.phi(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hNTracks"), cluster.deltaEta().size(), weight); for (size_t itrack = 0; itrack < cluster.deltaEta().size(); itrack++) { // Fill TrackEtaPhi histogram with delta phi and delta eta of all tracks saved in the vectors in skimmerGammaCalo.cxx - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTrackdEtadPhi"), cluster.deltaEta()[itrack], cluster.deltaPhi()[itrack], weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTrackdEtadPhi"), cluster.deltaEta()[itrack], cluster.deltaPhi()[itrack], weight); } if (do2DQA) { - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hNCell"), cluster.nCells(), cluster.e(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hM02"), cluster.m02(), cluster.e(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTime"), cluster.time(), cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hNCell"), cluster.nCells(), cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hM02"), cluster.m02(), cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTime"), cluster.time(), cluster.e(), weight); for (size_t itrack = 0; itrack < cluster.deltaEta().size(); itrack++) { - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTrackEOverP"), cluster.e() / cluster.trackp()[itrack], cluster.e(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTrackdEta"), cluster.deltaEta()[itrack], cluster.e(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTrackdPhi"), cluster.deltaPhi()[itrack], cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTrackEOverP"), cluster.e() / cluster.trackp()[itrack], cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTrackdEta"), cluster.deltaEta()[itrack], cluster.e(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTrackdPhi"), cluster.deltaPhi()[itrack], cluster.e(), weight); } } else { - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hNCell"), cluster.nCells(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hM02"), cluster.m02(), weight); - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTime"), cluster.time(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hNCell"), cluster.nCells(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hM02"), cluster.m02(), weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTime"), cluster.time(), weight); for (size_t itrack = 0; itrack < cluster.deltaEta().size(); itrack++) { - fRegistry->fill(HIST("Cluster/") + HIST(cluster_types[cls_id]) + HIST("hTrackEOverP"), cluster.e() / cluster.trackp()[itrack], weight); + fRegistry->fill(HIST("Cluster/") + HIST(kClusterTypes[cls_id]) + HIST("hTrackEOverP"), cluster.e() / cluster.trackp()[itrack], weight); } } } diff --git a/PWGEM/PhotonMeson/Utils/EventHistograms.h b/PWGEM/PhotonMeson/Utils/EventHistograms.h index ef5f115e216..5a00d9b6da9 100644 --- a/PWGEM/PhotonMeson/Utils/EventHistograms.h +++ b/PWGEM/PhotonMeson/Utils/EventHistograms.h @@ -9,22 +9,29 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// \header file for histograms +/// \file EventHistograms.h +/// \brief header file for event histograms /// \author daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_EVENTHISTOGRAMS_H_ #define PWGEM_PHOTONMESON_UTILS_EVENTHISTOGRAMS_H_ -#include "Framework/HistogramRegistry.h" +#include "Common/CCDB/EventSelectionParams.h" +#include "Common/CCDB/TriggerAliases.h" -using namespace o2::framework; +#include +#include + +#include + +#include namespace o2::aod::pwgem::photonmeson::utils::eventhistogram { -void addEventHistograms(HistogramRegistry* fRegistry) +inline void addEventHistograms(o2::framework::HistogramRegistry* fRegistry) { // event info - auto hCollisionCounter = fRegistry->add("Event/before/hCollisionCounter", "collision counter;;Number of events", kTH1D, {{12, 0.5, 12.5}}, false); + auto hCollisionCounter = fRegistry->add("Event/before/hCollisionCounter", "collision counter;;Number of events", o2::framework::kTH1D, {{12, 0.5, 12.5}}, false); hCollisionCounter->GetXaxis()->SetBinLabel(1, "all"); hCollisionCounter->GetXaxis()->SetBinLabel(2, "No TF border"); hCollisionCounter->GetXaxis()->SetBinLabel(3, "No ITS ROF border"); @@ -38,86 +45,87 @@ void addEventHistograms(HistogramRegistry* fRegistry) hCollisionCounter->GetXaxis()->SetBinLabel(11, "EMC L0 Triggered"); hCollisionCounter->GetXaxis()->SetBinLabel(12, "accepted"); - const AxisSpec axis_cent_ft0m{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, - 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, - 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, - 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, - 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, - "centrality FT0M (%)"}; + const o2::framework::AxisSpec axis_cent_ft0m{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, + 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, + 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, + 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, + 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, + "centrality FT0M (%)"}; - const AxisSpec axis_cent_ft0a{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, - 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, - 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, - 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, - 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, - "centrality FT0A (%)"}; + const o2::framework::AxisSpec axis_cent_ft0a{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, + 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, + 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, + 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, + 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, + "centrality FT0A (%)"}; - const AxisSpec axis_cent_ft0c{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, - 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, - 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, - 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, - 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, - "centrality FT0C (%)"}; + const o2::framework::AxisSpec axis_cent_ft0c{{0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.06, 0.07, 0.08, 0.09, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, + 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, + 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, + 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, + 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110}, + "centrality FT0C (%)"}; - fRegistry->add("Event/before/hZvtx", "vertex z; Z_{vtx} (cm)", kTH1F, {{100, -50, +50}}, false); - fRegistry->add("Event/before/hMultNTracksPV", "hMultNTracksPV; N_{track} to PV", kTH1F, {{6001, -0.5, 6000.5}}, false); - fRegistry->add("Event/before/hMultNTracksPVeta1", "hMultNTracksPVeta1; N_{track} to PV", kTH1F, {{6001, -0.5, 6000.5}}, false); - fRegistry->add("Event/before/hMultFT0", "hMultFT0;mult. FT0A;mult. FT0C", kTH2F, {{200, 0, 200000}, {60, 0, 60000}}, false); - fRegistry->add("Event/before/hCentFT0A", "hCentFT0A;centrality FT0A (%)", kTH1F, {{axis_cent_ft0a}}, false); - fRegistry->add("Event/before/hCentFT0C", "hCentFT0C;centrality FT0C (%)", kTH1F, {{axis_cent_ft0c}}, false); - fRegistry->add("Event/before/hCentFT0M", "hCentFT0M;centrality FT0M (%)", kTH1F, {{axis_cent_ft0m}}, false); - fRegistry->add("Event/before/hCentFT0CvsMultNTracksPV", "hCentFT0CvsMultNTracksPV;centrality FT0C (%);N_{track} to PV", kTH2F, {{110, 0, 110}, {500, 0, 5000}}, false); - fRegistry->add("Event/before/hMultFT0CvsMultNTracksPV", "hMultFT0CvsMultNTracksPV;mult. FT0C;N_{track} to PV", kTH2F, {{60, 0, 60000}, {500, 0, 5000}}, false); - fRegistry->add("Event/before/hMultFT0CvsOccupancy", "hMultFT0CvsOccupancy;mult. FT0C;N_{track} in time range", kTH2F, {{60, 0, 60000}, {2000, 0, 20000}}, false); + fRegistry->add("Event/before/hZvtx", "vertex z; Z_{vtx} (cm)", o2::framework::kTH1F, {{100, -50, +50}}, false); + fRegistry->add("Event/before/hMultNTracksPV", "hMultNTracksPV; N_{track} to PV", o2::framework::kTH1F, {{6001, -0.5, 6000.5}}, false); + fRegistry->add("Event/before/hMultNTracksPVeta1", "hMultNTracksPVeta1; N_{track} to PV", o2::framework::kTH1F, {{6001, -0.5, 6000.5}}, false); + fRegistry->add("Event/before/hMultFT0", "hMultFT0;mult. FT0A;mult. FT0C", o2::framework::kTH2F, {{200, 0, 200000}, {60, 0, 60000}}, false); + fRegistry->add("Event/before/hCentFT0A", "hCentFT0A;centrality FT0A (%)", o2::framework::kTH1F, {{axis_cent_ft0a}}, false); + fRegistry->add("Event/before/hCentFT0C", "hCentFT0C;centrality FT0C (%)", o2::framework::kTH1F, {{axis_cent_ft0c}}, false); + fRegistry->add("Event/before/hCentFT0M", "hCentFT0M;centrality FT0M (%)", o2::framework::kTH1F, {{axis_cent_ft0m}}, false); + fRegistry->add("Event/before/hCentFT0CvsMultNTracksPV", "hCentFT0CvsMultNTracksPV;centrality FT0C (%);N_{track} to PV", o2::framework::kTH2F, {{110, 0, 110}, {500, 0, 5000}}, false); + fRegistry->add("Event/before/hMultFT0CvsMultNTracksPV", "hMultFT0CvsMultNTracksPV;mult. FT0C;N_{track} to PV", o2::framework::kTH2F, {{60, 0, 60000}, {500, 0, 5000}}, false); + fRegistry->add("Event/before/hMultFT0CvsOccupancy", "hMultFT0CvsOccupancy;mult. FT0C;N_{track} in time range", o2::framework::kTH2F, {{60, 0, 60000}, {2000, 0, 20000}}, false); fRegistry->addClone("Event/before/", "Event/after/"); } template -void fillEventInfo(HistogramRegistry* fRegistry, TCollision const& collision, float weight = 1.) +void fillEventInfo(o2::framework::HistogramRegistry* fRegistry, TCollision const& collision, float weight = 1.) { - static constexpr std::string_view event_types[2] = {"before/", "after/"}; - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 1.0, weight); + const float maxZ = 10.f; + static constexpr std::string_view kEventTypes[2] = {"before/", "after/"}; + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 1.0, weight); if (collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 2.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 2.0, weight); } if (collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 3.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 3.0, weight); } if (collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 4.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 4.0, weight); } if (collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 5.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 5.0, weight); } if (collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 6.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 6.0, weight); } if (collision.selection_bit(o2::aod::evsel::kIsTriggerTVX)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 7.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 7.0, weight); } if (collision.sel8()) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 8.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 8.0, weight); } - if (std::abs(collision.posZ()) < 10.0) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 9.0, weight); + if (std::abs(collision.posZ()) < maxZ) { + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 9.0, weight); } - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hZvtx"), collision.posZ(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hZvtx"), collision.posZ(), weight); if (collision.alias_bit(kTVXinEMC)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 10.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 10.0, weight); } if (collision.alias_bit(kEMC7) || collision.alias_bit(kDMC7)) { - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCollisionCounter"), 11.0, weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCollisionCounter"), 11.0, weight); } - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hMultNTracksPV"), collision.multNTracksPV(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hMultNTracksPVeta1"), collision.multNTracksPVeta1(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hMultFT0"), collision.multFT0A(), collision.multFT0C(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCentFT0A"), collision.centFT0A(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCentFT0C"), collision.centFT0C(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCentFT0M"), collision.centFT0M(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hCentFT0CvsMultNTracksPV"), collision.centFT0C(), collision.multNTracksPV(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hMultFT0CvsMultNTracksPV"), collision.multFT0C(), collision.multNTracksPV(), weight); - fRegistry->fill(HIST("Event/") + HIST(event_types[ev_id]) + HIST("hMultFT0CvsOccupancy"), collision.multFT0C(), collision.trackOccupancyInTimeRange(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hMultNTracksPV"), collision.multNTracksPV(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hMultNTracksPVeta1"), collision.multNTracksPVeta1(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hMultFT0"), collision.multFT0A(), collision.multFT0C(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCentFT0A"), collision.centFT0A(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCentFT0C"), collision.centFT0C(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCentFT0M"), collision.centFT0M(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hCentFT0CvsMultNTracksPV"), collision.centFT0C(), collision.multNTracksPV(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hMultFT0CvsMultNTracksPV"), collision.multFT0C(), collision.multNTracksPV(), weight); + fRegistry->fill(HIST("Event/") + HIST(kEventTypes[ev_id]) + HIST("hMultFT0CvsOccupancy"), collision.multFT0C(), collision.trackOccupancyInTimeRange(), weight); } } // namespace o2::aod::pwgem::photonmeson::utils::eventhistogram diff --git a/PWGEM/PhotonMeson/Utils/HNMUtilities.h b/PWGEM/PhotonMeson/Utils/HNMUtilities.h index 80101807641..e4d71bab8f1 100644 --- a/PWGEM/PhotonMeson/Utils/HNMUtilities.h +++ b/PWGEM/PhotonMeson/Utils/HNMUtilities.h @@ -19,22 +19,21 @@ #ifndef PWGEM_PHOTONMESON_UTILS_HNMUTILITIES_H_ #define PWGEM_PHOTONMESON_UTILS_HNMUTILITIES_H_ -#include -#include -#include -#include "TVector3.h" - -#include "Framework/runDataProcessing.h" -#include "Framework/AnalysisTask.h" -#include "Framework/HistogramRegistry.h" -#include "Common/DataModel/EventSelection.h" -#include "EMCALBase/Geometry.h" -#include "PWGJE/DataModel/EMCALClusters.h" #include "PWGEM/PhotonMeson/DataModel/gammaTables.h" #include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "Common/Core/trackUtilities.h" -#include "Common/DataModel/TrackSelectionTables.h" -#include "EventFiltering/filterTables.h" + +#include +#include + +#include // IWYU pragma: keep +#include + +#include + +#include +#include +#include +#include namespace o2::aod::pwgem::photonmeson::hnmutilities { @@ -124,9 +123,9 @@ struct HeavyNeutralMeson { }; const int nSMEdges = 9; -float smPhiEdges[nSMEdges] = {1.75, 2.1, 2.45, 2.8, 3.14, 4., 4.89, 5.24, 5.58}; +inline float smPhiEdges[nSMEdges] = {1.75, 2.1, 2.45, 2.8, 3.14, 4., 4.89, 5.24, 5.58}; -int getSMNumber(float eta, float phi) +inline int getSMNumber(float eta, float phi) { int smNumber = 0; for (int iPhiInterval = 0; iPhiInterval < nSMEdges; iPhiInterval++) { @@ -175,7 +174,7 @@ void storeGammasInVector(C clusters, std::vector& vPhotons, std::array vPhotons, std::vector& vGGs) +inline void reconstructGGs(std::vector vPhotons, std::vector& vGGs) { vGGs.clear(); // loop over all photon combinations and build meson candidates @@ -195,7 +194,7 @@ void reconstructGGs(std::vector vPhotons, std::vector& v } /// \brief Reconstruct heavy neutral mesons from the given pion, antipion and the GG candidates and add them to the vHNMs vector -void reconstructHeavyNeutralMesons(ROOT::Math::PtEtaPhiMVector const& vecPiPlPiMi, std::vector& vGGs, std::vector& vHNMs) +inline void reconstructHeavyNeutralMesons(ROOT::Math::PtEtaPhiMVector const& vecPiPlPiMi, std::vector& vGGs, std::vector& vHNMs) { for (size_t iGG = 0; iGG < vGGs.size(); iGG++) { HeavyNeutralMeson heavyNeutralMeson(&vGGs.at(iGG), vecPiPlPiMi.E(), vecPiPlPiMi.Px(), vecPiPlPiMi.Py(), vecPiPlPiMi.Pz()); diff --git a/PWGEM/PhotonMeson/Utils/MCUtilities.h b/PWGEM/PhotonMeson/Utils/MCUtilities.h index c72a8a0b06d..8164d846c02 100644 --- a/PWGEM/PhotonMeson/Utils/MCUtilities.h +++ b/PWGEM/PhotonMeson/Utils/MCUtilities.h @@ -9,16 +9,18 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// \commonly used for MC analysis. +/// \file MCUtilities.h +/// \brief commonly used for MC analysis. /// \author daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_MCUTILITIES_H_ #define PWGEM_PHOTONMESON_UTILS_MCUTILITIES_H_ -#include "Framework/AnalysisTask.h" -#include "Framework/HistogramRegistry.h" +#include #include +#include +#include #include //_______________________________________________________________________ @@ -50,7 +52,7 @@ int IsFromWD(TCollision const&, T const& mctrack, TMCs const& mcTracks) if (motherid < mcTracks.size()) { // protect against bad mother indices. why is this needed? auto mp = mcTracks.iteratorAt(motherid); int pdg_mother = mp.pdgCode(); - if (abs(pdg_mother) == 310 || abs(pdg_mother) == 130 || abs(pdg_mother) == 3122) { + if (std::abs(pdg_mother) == kK0Short || std::abs(pdg_mother) == kK0Long || std::abs(pdg_mother) == kLambda0) { // LOGF(info, "mctrack.globalIndex() = %d, mp.globalIndex() = %d , pdg_mother = %d", mctrack.globalIndex(), mp.globalIndex(), pdg_mother); return motherid; } @@ -90,32 +92,36 @@ int IsXFromY(T const& mctrack, TMCs const& mcTracks, const int pdgX, const int p // Go up the decay chain of a mcparticle looking for a mother with the given pdg codes, if found return this mothers daughter // E.g. Find the gamma that was created in a pi0 or eta decay template -int FindMotherInChain(T const& mcparticle, TMCs const& mcparticles, TTargetPDGs const& motherpdgs, const int Depth = 50) +int FindMotherInChain(T const& mcparticle, TMCs const& mcparticles, TTargetPDGs const& motherpdgs, const int Depth = 50) // o2-linter: disable=pdg/explicit-code (false positive) { - if (!mcparticle.has_mothers() || Depth < 1) + if (!mcparticle.has_mothers() || Depth < 1) { return -1; + } int motherid = mcparticle.mothersIds()[0]; auto mother = mcparticles.iteratorAt(motherid); - if (std::find(motherpdgs.begin(), motherpdgs.end(), mother.pdgCode()) != motherpdgs.end()) + if (std::find(motherpdgs.begin(), motherpdgs.end(), mother.pdgCode()) != motherpdgs.end()) { return mcparticle.globalIndex(); // The mother has the required pdg code, so return its daughters global mc particle code. - else + } else { return FindMotherInChain(mother, mcparticles, motherpdgs, Depth - 1); + } } //_______________________________________________________________________ template int IsEleFromPC(T const& mctrack, TMCs const& mcTracks) { // is election from photon conversion? returns index of mother photon - if (abs(mctrack.pdgCode()) != 11) + if (std::abs(mctrack.pdgCode()) != kElectron) { return -1; - if (mctrack.producedByGenerator()) + } + if (mctrack.producedByGenerator()) { return -1; + } if (mctrack.has_mothers()) { int motherid = mctrack.mothersIds()[0]; // first mother auto mp = mcTracks.iteratorAt(motherid); int pdg_mother = mp.pdgCode(); - if (pdg_mother == 22) { + if (pdg_mother == kGamma) { return motherid; } } else { @@ -191,7 +197,7 @@ bool IsInAcceptance(TMCParticle const& mcparticle, TMCParticles const& mcparticl } std::vector pdgs; pdgs.reserve(target_pdgs.size()); - for (auto& daughterId : daughtersIds) { + for (const auto& daughterId : daughtersIds) { if (daughterId < 0) { pdgs.clear(); pdgs.shrink_to_fit(); @@ -226,21 +232,21 @@ bool IsInAcceptance(TMCParticle const& mcparticle, TMCParticles const& mcparticl template bool IsConversionPointInAcceptance(TMCPhoton const& mcphoton, const float max_r_gen, const float max_eta_gen, const float margin_z_mc, TMCParticles const& mcparticles) { - if (std::abs(mcphoton.pdgCode()) != 22) { + if (std::abs(mcphoton.pdgCode()) != kGamma) { return false; } auto daughtersIds = mcphoton.daughtersIds(); - if (daughtersIds.size() != 2) { + if (daughtersIds.size() != 2) { // o2-linter: disable=magic-number (2 is not that magic in this context) return false; } - for (auto& daughterId : daughtersIds) { + for (const auto& daughterId : daughtersIds) { if (daughterId < 0) { return false; } auto daughter = mcparticles.iteratorAt(daughterId); - if (std::abs(daughter.pdgCode()) != 11) { + if (std::abs(daughter.pdgCode()) != kElectron) { return false; } @@ -248,9 +254,9 @@ bool IsConversionPointInAcceptance(TMCPhoton const& mcphoton, const float max_r_ return false; } - float rxy_gen_e = std::sqrt(std::pow(daughter.vx(), 2) + pow(daughter.vy(), 2)); + float rxy_gen_e = std::sqrt(std::pow(daughter.vx(), 2) + std::pow(daughter.vy(), 2)); // LOGF(info, "daughterId = %d , pdg = %d , vx = %f , vy = %f , vz = %f, rxy = %f", daughterId, daughter.pdgCode(), daughter.vx(), daughter.vy(), daughter.vz(), rxy_gen_e); - if (rxy_gen_e > max_r_gen || rxy_gen_e < abs(daughter.vz()) * std::tan(2 * std::atan(std::exp(-max_eta_gen))) - margin_z_mc) { + if (rxy_gen_e > max_r_gen || rxy_gen_e < std::abs(daughter.vz()) * std::tan(2 * std::atan(std::exp(-max_eta_gen))) - margin_z_mc) { return false; } } // end of daughter loop @@ -262,11 +268,13 @@ template bool isGammaGammaDecay(TMCParticle mcParticle, TMCParticles mcParticles) { auto daughtersIds = mcParticle.daughtersIds(); - if (daughtersIds.size() != 2) + if (daughtersIds.size() != 2) { // o2-linter: disable=magic-number (2 is not that magic in this context) return false; - for (auto& daughterId : daughtersIds) { - if (mcParticles.iteratorAt(daughterId).pdgCode() != 22) + } + for (const auto& daughterId : daughtersIds) { + if (mcParticles.iteratorAt(daughterId).pdgCode() != kGamma) { return false; + } } return true; } diff --git a/PWGEM/PhotonMeson/Utils/NMHistograms.h b/PWGEM/PhotonMeson/Utils/NMHistograms.h index 8915c6624b5..b04b8a67878 100644 --- a/PWGEM/PhotonMeson/Utils/NMHistograms.h +++ b/PWGEM/PhotonMeson/Utils/NMHistograms.h @@ -9,58 +9,59 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// \header file for histograms +/// \file NMHistograms.h +/// \brief utility for mass vs pT histograms mainly /// \author daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_NMHISTOGRAMS_H_ #define PWGEM_PHOTONMESON_UTILS_NMHISTOGRAMS_H_ #include "PWGEM/PhotonMeson/Utils/MCUtilities.h" -#include "PWGEM/PhotonMeson/Utils/PairUtilities.h" -#include "Framework/HistogramRegistry.h" +#include +#include -#include "TF1.h" +#include +#include +#include +#include #include -using namespace o2::framework; -using namespace o2::aod::pwgem::photonmeson::utils::mcutil; - namespace o2::aod::pwgem::photonmeson::utils::nmhistogram { -void addNMHistograms(HistogramRegistry* fRegistry, bool isMC, const char* pairname = "#gamma#gamma") +inline void addNMHistograms(o2::framework::HistogramRegistry* fRegistry, bool isMC, const char* pairname = "#gamma#gamma") { // !!Don't change pt,eta,y binning. These binnings have to be consistent with binned data at skimming.!! std::vector ptbins; - for (int i = 0; i < 2; i++) { + for (int i = 0; i < 2; i++) { // o2-linter: disable=magic-number (just numbers for binning) ptbins.emplace_back(0.05 * (i - 0) + 0.0); // from 0 to 0.05 GeV/c, every 0.05 GeV/c } - for (int i = 2; i < 51; i++) { + for (int i = 2; i < 51; i++) { // o2-linter: disable=magic-number (just numbers for binning) ptbins.emplace_back(0.1 * (i - 2) + 0.1); // from 0.1 to 4.9 GeV/c, every 0.1 GeV/c } - for (int i = 51; i < 61; i++) { + for (int i = 51; i < 61; i++) { // o2-linter: disable=magic-number (just numbers for binning) ptbins.emplace_back(0.5 * (i - 51) + 5.0); // from 5 to 9.5 GeV/c, every 0.5 GeV/c } - for (int i = 61; i < 72; i++) { + for (int i = 61; i < 72; i++) { // o2-linter: disable=magic-number (just numbers for binning) ptbins.emplace_back(1.0 * (i - 61) + 10.0); // from 10 to 20 GeV/c, every 1 GeV/c } - const AxisSpec axis_pt{ptbins, Form("p_{T,%s} (GeV/c)", pairname)}; - const AxisSpec axis_mass{400, 0, 0.8, Form("m_{%s} (GeV/c^{2})", pairname)}; + const o2::framework::AxisSpec axis_pt{ptbins, Form("p_{T,%s} (GeV/c)", pairname)}; + const o2::framework::AxisSpec axis_mass{400, 0, 0.8, Form("m_{%s} (GeV/c^{2})", pairname)}; if (isMC) { - fRegistry->add("Pair/Pi0/hs_Primary", "rec. true pi0", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Pi0/hs_FromWD", "rec. true pi0 from weak decay", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Pi0/hs_FromHS", "rec. true pi0 from hadronic shower in material", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Pi0/hs_FromSameGamma", "Two clusters from same gamma that is a pi0 daughter (conversion)", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Eta/hs_FromSameGamma", "Two clusters from same gamma that is a eta daughter (conversion)", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Eta/hs_Primary", "rec. true eta", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Eta/hs_FromWD", "rec. true eta from weak decay", kTHnSparseD, {axis_mass, axis_pt}, true); - fRegistry->add("Pair/Eta/hs_FromHS", "rec. true eta from hadronic shower in material", kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Pi0/hs_Primary", "rec. true pi0", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Pi0/hs_FromWD", "rec. true pi0 from weak decay", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Pi0/hs_FromHS", "rec. true pi0 from hadronic shower in material", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Pi0/hs_FromSameGamma", "Two clusters from same gamma that is a pi0 daughter (conversion)", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Eta/hs_Primary", "rec. true eta", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Eta/hs_FromWD", "rec. true eta from weak decay", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Eta/hs_FromHS", "rec. true eta from hadronic shower in material", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/Eta/hs_FromSameGamma", "Two clusters from same gamma that is a eta daughter (conversion)", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); - const AxisSpec axis_rapidity{{0.0, +0.8, +0.9}, "rapidity |y|"}; - fRegistry->add("Generated/Pi0/hPt", "pT;p_{T} (GeV/c)", kTH1F, {axis_pt}, true); - fRegistry->add("Generated/Pi0/hPtY", "Generated info", kTH2F, {axis_pt, axis_rapidity}, true); + const o2::framework::AxisSpec axis_rapidity{{0.0, +0.8, +0.9}, "rapidity |y|"}; + fRegistry->add("Generated/Pi0/hPt", "pT;p_{T} (GeV/c)", o2::framework::kTH1F, {axis_pt}, true); + fRegistry->add("Generated/Pi0/hPtY", "Generated info", o2::framework::kTH2F, {axis_pt, axis_rapidity}, true); fRegistry->addClone("Generated/Pi0/", "Generated/Eta/"); fRegistry->get(HIST("Generated/Pi0/hPt"))->SetXTitle("p_{T} (GeV/c)"); @@ -70,24 +71,24 @@ void addNMHistograms(HistogramRegistry* fRegistry, bool isMC, const char* pairna fRegistry->get(HIST("Generated/Eta/hPtY"))->SetXTitle("p_{T} (GeV/c)"); fRegistry->get(HIST("Generated/Eta/hPtY"))->SetYTitle("rapidity |y|"); } else { - fRegistry->add("Pair/same/hs", "diphoton", kTHnSparseD, {axis_mass, axis_pt}, true); + fRegistry->add("Pair/same/hs", "diphoton", o2::framework::kTHnSparseD, {axis_mass, axis_pt}, true); fRegistry->addClone("Pair/same/", "Pair/mix/"); } } template -void fillTruePairInfo(HistogramRegistry* fRegistry, TDiphoton const& v12, TMCParitlce const& mcparticle, TMCParticles const& mcparticles, TMCCollisions const&, const TF1* f1fd_k0s_to_pi0 = nullptr, float eventWeight = 1.f) +void fillTruePairInfo(o2::framework::HistogramRegistry* fRegistry, TDiphoton const& v12, TMCParitlce const& mcparticle, TMCParticles const& mcparticles, TMCCollisions const&, const TF1* f1fd_k0s_to_pi0 = nullptr, float eventWeight = 1.f) { int pdg = std::abs(mcparticle.pdgCode()); float weight = eventWeight; - int motherid_strhad = IsFromWD(mcparticle.template emmcevent_as(), mcparticle, mcparticles); + int motherid_strhad = o2::aod::pwgem::photonmeson::utils::mcutil::IsFromWD(mcparticle.template emmcevent_as(), mcparticle, mcparticles); switch (pdg) { - case 111: { + case kPi0: { if (mcparticle.isPhysicalPrimary() || mcparticle.producedByGenerator()) { fRegistry->fill(HIST("Pair/Pi0/hs_Primary"), v12.M(), v12.Pt(), weight); } else if (motherid_strhad > 0) { auto str_had = mcparticles.iteratorAt(motherid_strhad); - if (std::abs(str_had.pdgCode()) == 310 && f1fd_k0s_to_pi0 != nullptr) { + if (std::abs(str_had.pdgCode()) == kK0Short && f1fd_k0s_to_pi0 != nullptr) { weight *= f1fd_k0s_to_pi0->Eval(str_had.pt()); } fRegistry->fill(HIST("Pair/Pi0/hs_FromWD"), v12.M(), v12.Pt(), weight); diff --git a/PWGEM/PhotonMeson/Utils/PCMUtilities.h b/PWGEM/PhotonMeson/Utils/PCMUtilities.h index 095179a24ad..269db4f46fb 100644 --- a/PWGEM/PhotonMeson/Utils/PCMUtilities.h +++ b/PWGEM/PhotonMeson/Utils/PCMUtilities.h @@ -9,17 +9,27 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// \commonly used for PCM analyses. +/// \file PCMUtilities.h +/// \brief helper functions commonly used for PCM analyses. /// \author daiki.sekihata@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_PCMUTILITIES_H_ #define PWGEM_PHOTONMESON_UTILS_PCMUTILITIES_H_ -#include -#include "ReconstructionDataFormats/HelixHelper.h" -#include "DetectorsBase/Propagator.h" -#include "Common/Core/trackUtilities.h" #include "Common/Core/RecoDecay.h" +#include "Common/Core/trackUtilities.h" + +#include +#include +#include +#include + +#include + +#include +#include + +#include //_______________________________________________________________________ inline bool checkAP(const float alpha, const float qt, const float alpha_max = 0.95, const float qt_max = 0.05) @@ -64,8 +74,8 @@ inline void Vtx_recalculationParCov(o2::base::Propagator* prop, const o2::track: // I think this calculation gets the closest point on the track to the conversion point // This alpha is a different alpha than the usual alpha and I think it is the angle between X axis and conversion point - float alphaPos = M_PI + std::atan2(-(xyz[1] - helixPos.yC), -(xyz[0] - helixPos.xC)); - float alphaNeg = M_PI + std::atan2(-(xyz[1] - helixNeg.yC), -(xyz[0] - helixNeg.xC)); + float alphaPos = o2::constants::math::PI + std::atan2(-(xyz[1] - helixPos.yC), -(xyz[0] - helixPos.xC)); + float alphaNeg = o2::constants::math::PI + std::atan2(-(xyz[1] - helixNeg.yC), -(xyz[0] - helixNeg.xC)); float vertexXPos = helixPos.xC + helixPos.rC * std::cos(alphaPos); float vertexYPos = helixPos.yC + helixPos.rC * std::sin(alphaPos); diff --git a/PWGEM/PhotonMeson/Utils/TrackSelection.h b/PWGEM/PhotonMeson/Utils/TrackSelection.h index bdc2f5be533..81803c9cee7 100644 --- a/PWGEM/PhotonMeson/Utils/TrackSelection.h +++ b/PWGEM/PhotonMeson/Utils/TrackSelection.h @@ -9,13 +9,16 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. +/// \file TrackSelection.h /// \brief helper functions for pair track selection /// \author felix.schlepper@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_TRACKSELECTION_H_ #define PWGEM_PHOTONMESON_UTILS_TRACKSELECTION_H_ -#include "TPDGCode.h" +#include + +#include namespace o2::pwgem::photonmeson { diff --git a/PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h b/PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h index f3eac77ca04..0fc655d8858 100644 --- a/PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h +++ b/PWGEM/PhotonMeson/Utils/emcalHistoDefinitions.h @@ -9,7 +9,8 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. -/// \brief commonly used histogram (axis) definitions for emcal in PWGEM +/// \file emcalHistoDefinitions.h +/// \brief commonly used histograms and axes definitions for emcal in PWGEM /// \author marvin.hemmer@cern.ch #include @@ -19,32 +20,32 @@ #ifndef PWGEM_PHOTONMESON_UTILS_EMCALHISTODEFINITIONS_H_ #define PWGEM_PHOTONMESON_UTILS_EMCALHISTODEFINITIONS_H_ -using namespace o2::framework; - -AxisSpec const gAxis_dEta{100, -0.2f, 0.2f, "#Delta#it{#eta}"}; -AxisSpec const gAxis_dPhi{100, -0.2f, 0.2f, "#Delta#it{#varphi} (rad)"}; -AxisSpec const gAxis_dR{100, 0.0f, 0.2f, "#Delta#it{R}"}; -AxisSpec const gAxis_MatchedTrack{10, 0.5f, 10.5f, "matched track"}; -AxisSpec const gAxis_Eta{160, -0.8f, 0.8f, "#it{#eta}"}; -AxisSpec const gAxis_Phi{144, 0.f, 2 * 3.14159, "#it{#varphi} (rad)"}; - -std::vector EClusBins = {0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, - 1, 1.2, 1.4, 1.6, 1.8, 2, 2.2, 2.4, 2.6, 2.8, - 3, 3.2, 3.4, 3.6, 3.8, 4, 4.2, 4.4, 4.6, 4.8, - 5, 5.5, 6, 6.5, 7, 7.5, 8, 8.5, 9, 9.5, 10, 11, - 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, - 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, - 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, - 48, 49, 50, 55, 60, 65, 70, 75, 80, 85, 90, 95, - 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200}; -AxisSpec gAxis_EClus = {EClusBins, "#it{E}_{clus} (GeV)"}; - -HistogramConfigSpec gHistoSpec_clusterTM_dEtadPhi({HistType::kTH2F, {gAxis_dEta, gAxis_dPhi}}); -HistogramConfigSpec gHistoSpec_clusterTM_dEtaMT({HistType::kTH2F, {gAxis_dEta, {10, 0.5, 10.5}}}); -HistogramConfigSpec gHistoSpec_clusterTM_dPhiMT({HistType::kTH2F, {gAxis_dPhi, {10, 0.5, 10.5}}}); -HistogramConfigSpec gHistoSpec_clusterTM_dRMT({HistType::kTH2F, {gAxis_dR, {10, 0.5, 10.5}}}); -HistogramConfigSpec gHistoSpec_EtaPhi({HistType::kTH2F, {gAxis_Eta, gAxis_Phi}}); -HistogramConfigSpec gHistoSpec_clusterE({HistType::kTH1F, {gAxis_EClus}}); -HistogramConfigSpec gHistoSpec_clusterECuts({HistType::kTH2F, {gAxis_EClus, {64, -0.5, 63.5}}}); +o2::framework::AxisSpec const gAxisdEta{200, -0.1f, 0.1f, "#Delta#it{#eta}"}; +o2::framework::AxisSpec const gAxisdPhi{200, -0.1f, 0.1f, "#Delta#it{#varphi} (rad)"}; +o2::framework::AxisSpec const gAxisdR{100, 0.0f, 0.2f, "#Delta#it{R}"}; +o2::framework::AxisSpec const gAxisMatchedTrack{10, 0.5f, 10.5f, "matched track"}; +o2::framework::AxisSpec const gAxisEta{160, -0.8f, 0.8f, "#it{#eta}"}; +o2::framework::AxisSpec const gAxisPhi{144, 0.f, 2 * 3.14159, "#it{#varphi} (rad)"}; +o2::framework::AxisSpec const gAxisEoverP{250, 0.f, 2.5f, "#it{E}_{clus.}/#it{p}_{track} (#it{c})"}; + +const std::vector eClusBins = {0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, + 1, 1.2, 1.4, 1.6, 1.8, 2, 2.2, 2.4, 2.6, 2.8, + 3, 3.2, 3.4, 3.6, 3.8, 4, 4.2, 4.4, 4.6, 4.8, + 5, 5.5, 6, 6.5, 7, 7.5, 8, 8.5, 9, 9.5, 10, 11, + 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, + 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, + 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, + 48, 49, 50, 55, 60, 65, 70, 75, 80, 85, 90, 95, + 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200}; +o2::framework::AxisSpec const gAxisEClus = {eClusBins, "#it{E}_{clus} (GeV)"}; + +o2::framework::HistogramConfigSpec const gHistoSpecClusterTMdEtadPhi({o2::framework::HistType::kTH2F, {gAxisdEta, gAxisdPhi}}); +o2::framework::HistogramConfigSpec const gHistoSpecClusterTMdEtaMT({o2::framework::HistType::kTH2F, {gAxisdEta, {10, 0.5, 10.5}}}); +o2::framework::HistogramConfigSpec const gHistoSpecClusterTMdPhiMT({o2::framework::HistType::kTH2F, {gAxisdPhi, {10, 0.5, 10.5}}}); +o2::framework::HistogramConfigSpec const gHistoSpecClusterTMdRMT({o2::framework::HistType::kTH2F, {gAxisdR, {10, 0.5, 10.5}}}); +o2::framework::HistogramConfigSpec const gHistoSpecEtaPhi({o2::framework::HistType::kTH2F, {gAxisEta, gAxisPhi}}); +o2::framework::HistogramConfigSpec const gHistoSpecClusterE({o2::framework::HistType::kTH1F, {gAxisEClus}}); +o2::framework::HistogramConfigSpec const gHistoSpecClusterECuts({o2::framework::HistType::kTH2F, {gAxisEClus, {64, -0.5, 63.5}}}); +o2::framework::HistogramConfigSpec const gHistoSpecTMEoverP({o2::framework::HistType::kTH2D, {gAxisEClus, gAxisEoverP}}); #endif // PWGEM_PHOTONMESON_UTILS_EMCALHISTODEFINITIONS_H_ diff --git a/PWGEM/PhotonMeson/Utils/gammaConvDefinitions.h b/PWGEM/PhotonMeson/Utils/gammaConvDefinitions.h index a9257d90ee4..2eb11525278 100644 --- a/PWGEM/PhotonMeson/Utils/gammaConvDefinitions.h +++ b/PWGEM/PhotonMeson/Utils/gammaConvDefinitions.h @@ -9,13 +9,15 @@ // granted to it by virtue of its status as an Intergovernmental Organization // or submit itself to any jurisdiction. +/// \file gammaConvDefinitions.h /// \brief commonly used definitions for gammaConv tasks /// \author stephan.friedrich.stiefelmaier@cern.ch #ifndef PWGEM_PHOTONMESON_UTILS_GAMMACONVDEFINITIONS_H_ #define PWGEM_PHOTONMESON_UTILS_GAMMACONVDEFINITIONS_H_ -#include "Framework/AnalysisTask.h" +#include +#include using namespace o2::framework; diff --git a/PWGHF/HFL/TableProducer/electronSelectionWithTpcEmcal.cxx b/PWGHF/HFL/TableProducer/electronSelectionWithTpcEmcal.cxx index 0ffff07c55b..59e22fd8b03 100644 --- a/PWGHF/HFL/TableProducer/electronSelectionWithTpcEmcal.cxx +++ b/PWGHF/HFL/TableProducer/electronSelectionWithTpcEmcal.cxx @@ -53,7 +53,7 @@ using namespace o2::framework; using namespace o2::framework::expressions; using namespace o2::soa; -const int kEta = 221; +const int kEtaLocal = 221; struct HfElectronSelectionWithTpcEmcal { @@ -573,15 +573,14 @@ struct HfElectronSelectionWithTpcEmcal { if (particleMc.has_mothers()) { auto const& mother = particleMc.mothers_first_as(); - if (std::abs(mother.pdgCode()) == kEta || std::abs(mother.pdgCode()) == kPi0 || std::abs(mother.pdgCode()) == kGamma) { - + if (std::abs(mother.pdgCode()) == kEtaLocal || std::abs(mother.pdgCode()) == kPi0 || std::abs(mother.pdgCode()) == kGamma) { registry.fill(HIST("hMcgenAllNonHfeElectron"), particleMc.pt()); auto const& gmother = mother.mothers_first_as(); // cases to consider: eta->e, eta->pi0->e, eta->gamma->e, eta->pi0->gamma->e, pi0->e, pi0->gamma->e //================= eta->e ====================================== - if (std::abs(mother.pdgCode()) == kEta) { + if (std::abs(mother.pdgCode()) == kEtaLocal) { if (mother.isPhysicalPrimary()) { if ((std::abs(gmother.pdgCode()) >= pdgCodeCharmMin && std::abs(gmother.pdgCode()) < pdgCodeCharmMax) || @@ -602,7 +601,7 @@ struct HfElectronSelectionWithTpcEmcal { } isEmbPi0 = true; // pi0 -> e } - if (std::abs(gmother.pdgCode()) == kEta) { + if (std::abs(gmother.pdgCode()) == kEtaLocal) { if (gmother.isPhysicalPrimary() || gmother.has_mothers()) { auto const& ggmother = gmother.mothers_first_as(); if ((std::abs(ggmother.pdgCode()) >= pdgCodeCharmMin && std::abs(ggmother.pdgCode()) < pdgCodeCharmMax) || @@ -617,7 +616,7 @@ struct HfElectronSelectionWithTpcEmcal { /// ==================================== eta->gamma->e and eta->pi0->gamma->e============ if (std::abs(mother.pdgCode()) == kGamma) { - if (std::abs(gmother.pdgCode()) == kEta) { + if (std::abs(gmother.pdgCode()) == kEtaLocal) { if (gmother.isPhysicalPrimary() || gmother.has_mothers()) { auto const& ggmother = gmother.mothers_first_as(); if ((std::abs(ggmother.pdgCode()) >= pdgCodeCharmMin && std::abs(ggmother.pdgCode()) < pdgCodeCharmMax) || @@ -638,7 +637,7 @@ struct HfElectronSelectionWithTpcEmcal { } if (gmother.has_mothers()) { auto const& ggmother = gmother.mothers_first_as(); - if (std::abs(ggmother.pdgCode()) == kEta) { + if (std::abs(ggmother.pdgCode()) == kEtaLocal) { if (ggmother.isPhysicalPrimary() || ggmother.has_mothers()) { auto const& gggmother = ggmother.mothers_first_as(); if ((std::abs(gggmother.pdgCode()) >= pdgCodeCharmMin && std::abs(gggmother.pdgCode()) < pdgCodeCharmMax) || diff --git a/PWGHF/HFL/Tasks/taskSingleElectron.cxx b/PWGHF/HFL/Tasks/taskSingleElectron.cxx index b53043cff80..6c4d1f76833 100644 --- a/PWGHF/HFL/Tasks/taskSingleElectron.cxx +++ b/PWGHF/HFL/Tasks/taskSingleElectron.cxx @@ -279,16 +279,16 @@ struct HfTaskSingleElectron { if (grmotherPdg == kPi0) { return GammaPi0; } - if (grmotherPdg == kEta) { + if (grmotherPdg == PdgCode::kEta) { return GammaEta; } - if (grmotherPdg == kOmega) { + if (grmotherPdg == PdgCode::kOmega) { return GammaOmega; } if (grmotherPdg == kPhi) { return GammaPhi; } - if (grmotherPdg == kEtaPrime) { + if (grmotherPdg == PdgCode::kEtaPrime) { return GammaEtaPrime; } if (grmotherPdg == kRho770_0) { @@ -320,12 +320,12 @@ struct HfTaskSingleElectron { mpt = grmotherPt; return GammaPi0; } - if (grmotherPdg == kEta) { + if (grmotherPdg == PdgCode::kEta) { mpdg = grmotherPdg; mpt = grmotherPt; return GammaEta; } - if (grmotherPdg == kOmega) { + if (grmotherPdg == PdgCode::kOmega) { mpdg = grmotherPdg; mpt = grmotherPt; return GammaOmega; @@ -335,7 +335,7 @@ struct HfTaskSingleElectron { mpt = grmotherPt; return GammaPhi; } - if (grmotherPdg == kEtaPrime) { + if (grmotherPdg == PdgCode::kEtaPrime) { mpdg = grmotherPdg; mpt = grmotherPt; return GammaEtaPrime; @@ -355,10 +355,10 @@ struct HfTaskSingleElectron { if (grmothersIdsVec.empty()) { static const std::map pdgToSource = { {kPi0, Pi0}, - {kEta, Eta}, - {kOmega, Omega}, + {PdgCode::kEta, Eta}, + {PdgCode::kOmega, Omega}, {kPhi, Phi}, - {kEtaPrime, EtaPrime}, + {PdgCode::kEtaPrime, EtaPrime}, {kRho770_0, Rho0}, {kKPlus, Ke3}, {kK0Long, K0l}}; @@ -393,16 +393,16 @@ struct HfTaskSingleElectron { mpdg = motherPdg; return Pi0; } - if (motherPdg == kEta) { + if (motherPdg == PdgCode::kEta) { return Eta; } - if (motherPdg == kOmega) { + if (motherPdg == PdgCode::kOmega) { return Omega; } if (motherPdg == kPhi) { return Phi; } - if (motherPdg == kEtaPrime) { + if (motherPdg == PdgCode::kEtaPrime) { return EtaPrime; } if (motherPdg == kRho770_0) {