Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 12 additions & 8 deletions PWGEM/PhotonMeson/DataModel/bcWiseTables.h
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.

Check failure on line 1 in PWGEM/PhotonMeson/DataModel/bcWiseTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/file-cpp]

Use lowerCamelCase or UpperCamelCase for names of C++ files. See the O2 naming conventions for details.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
Expand All @@ -19,10 +19,10 @@
#ifndef PWGEM_PHOTONMESON_DATAMODEL_BCWISETABLES_H_
#define PWGEM_PHOTONMESON_DATAMODEL_BCWISETABLES_H_

#include <limits>

#include "Framework/AnalysisDataModel.h"

#include <limits>

namespace o2::aod
{

Expand Down Expand Up @@ -72,14 +72,14 @@
DECLARE_SOA_COLUMN(StoredMu, storedMu, uint16_t); //! probability of TVX collision per BC (x1000)

DECLARE_SOA_DYNAMIC_COLUMN(Centrality, centrality, [](uint8_t storedcentrality) -> float { return std::nextafter(storedcentrality / emdownscaling::downscalingFactors[emdownscaling::kFT0MCent], std::numeric_limits<float>::infinity()); }); //! Centrality (0-100)
DECLARE_SOA_DYNAMIC_COLUMN(FT0MAmplitude, ft0Amplitude, [](uint16_t storedFT0MAmplitude) -> float { return std::nextafter(storedFT0MAmplitude / emdownscaling::downscalingFactors[emdownscaling::kFT0Amp], std::numeric_limits<float>::infinity()); }); //! FT0M amplitude

Check failure on line 75 in PWGEM/PhotonMeson/DataModel/bcWiseTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column]

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.
DECLARE_SOA_DYNAMIC_COLUMN(Mu, mu, [](uint16_t storedMu) -> float { return std::nextafter(storedMu / emdownscaling::downscalingFactors[emdownscaling::kMu], std::numeric_limits<float>::infinity()); }); //! probability of TVX collision per BC
} // namespace bcwisebc
DECLARE_SOA_TABLE(BCWiseBCs, "AOD", "BCWISEBC", //! table of bc wise centrality estimation and event selection input
o2::soa::Index<>, bcwisebc::HasFT0, bcwisebc::HasTVX, bcwisebc::HaskTVXinEMC, bcwisebc::HasEMCCell, bcwisebc::HasNoTFROFBorder, bcwisebc::StoredCentrality,
bcwisebc::StoredFT0MAmplitude, bcwisebc::StoredMu, bcwisebc::Centrality<bcwisebc::StoredCentrality>, bcwisebc::FT0MAmplitude<bcwisebc::StoredFT0MAmplitude>, bcwisebc::Mu<bcwisebc::StoredMu>);

DECLARE_SOA_INDEX_COLUMN(BCWiseBC, bcWiseBC); //! bunch crossing ID used as index

Check failure on line 82 in PWGEM/PhotonMeson/DataModel/bcWiseTables.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column]

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.

namespace bcwisecollision
{
Expand Down Expand Up @@ -121,30 +121,34 @@
bcwisecluster::Definition<bcwisecluster::StoredDefinition>, bcwisecluster::E<bcwisecluster::StoredE>, bcwisecluster::Eta<bcwisecluster::StoredEta>, bcwisecluster::Phi<bcwisecluster::StoredPhi>, bcwisecluster::NCells<bcwisecluster::StoredNCells>, bcwisecluster::M02<bcwisecluster::StoredM02>, bcwisecluster::Time<bcwisecluster::StoredTime>, bcwisecluster::IsExotic<bcwisecluster::StoredIsExotic>,
bcwisecluster::Pt<bcwisecluster::StoredE, bcwisecluster::StoredEta>);

namespace bcwisemcpi0s
namespace bcwisemcmesons
{
DECLARE_SOA_COLUMN(StoredPt, storedPt, uint16_t); //! Transverse momentum of generated pi0 (1 MeV -> Maximum pi0 pT of ~65 GeV)
DECLARE_SOA_COLUMN(IsAccepted, isAccepted, bool); //! Both decay photons are within the EMCal acceptance
DECLARE_SOA_COLUMN(IsPrimary, isPrimary, bool); //! mcParticle.isPhysicalPrimary() || mcParticle.producedByGenerator()
DECLARE_SOA_COLUMN(IsFromWD, isFromWD, bool); //! Pi0 from a weak decay according to pwgem::photonmeson::utils::mcutil::IsFromWD

DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, [](uint16_t storedpt) -> float { return std::nextafter(storedpt / emdownscaling::downscalingFactors[emdownscaling::kpT], std::numeric_limits<float>::infinity()); }); //! pT of pi0 (GeV)
} // namespace bcwisemcpi0s
} // namespace bcwisemcmesons

DECLARE_SOA_TABLE(BCWiseMCPi0s, "AOD", "BCWISEMCPI0", //! table of pi0s on MC level
o2::soa::Index<>, BCWiseBCId, bcwisemcpi0s::StoredPt, bcwisemcpi0s::IsAccepted, bcwisemcpi0s::IsPrimary, bcwisemcpi0s::IsFromWD,
bcwisemcpi0s::Pt<bcwisemcpi0s::StoredPt>);
o2::soa::Index<>, BCWiseBCId, bcwisemcmesons::StoredPt, bcwisemcmesons::IsAccepted, bcwisemcmesons::IsPrimary, bcwisemcmesons::IsFromWD,
bcwisemcmesons::Pt<bcwisemcmesons::StoredPt>);
DECLARE_SOA_TABLE(BCWiseMCEtas, "AOD", "BCWISEMCETA", //! table of eta mesons on MC level
o2::soa::Index<>, BCWiseBCId, bcwisemcmesons::StoredPt, bcwisemcmesons::IsAccepted, bcwisemcmesons::IsPrimary, bcwisemcmesons::IsFromWD,
bcwisemcmesons::Pt<bcwisemcmesons::StoredPt>);

namespace bcwisemccluster
{
DECLARE_SOA_COLUMN(Pi0ID, pi0ID, int32_t); //! Index of the mother pi0 (-1 if not from pi0)
DECLARE_SOA_COLUMN(MesonID, mesonID, int32_t); //! Index of the mother mesom (-1 if not from a pi0 or eta)
DECLARE_SOA_COLUMN(IsEta, isEta, bool); //! Boolean flag to indicate if the cluster is from an eta meson, otherwise it is from a pi0
DECLARE_SOA_COLUMN(StoredTrueE, storedTrueE, uint16_t); //! energy of cluster inducing particle (1 MeV -> Maximum cluster energy of ~65 GeV)

DECLARE_SOA_DYNAMIC_COLUMN(TrueE, trueE, [](uint16_t storedTrueE) -> float { return std::nextafter(storedTrueE / emdownscaling::downscalingFactors[emdownscaling::kEnergy], std::numeric_limits<float>::infinity()); }); //! energy of cluster inducing particle (GeV)
} // namespace bcwisemccluster

DECLARE_SOA_TABLE(BCWiseMCClusters, "AOD", "BCWISEMCCLS", //! table of MC information for clusters -> To be joined with the cluster table
o2::soa::Index<>, BCWiseBCId, bcwisemccluster::Pi0ID, bcwisemccluster::StoredTrueE,
o2::soa::Index<>, BCWiseBCId, bcwisemccluster::MesonID, bcwisemccluster::IsEta, bcwisemccluster::StoredTrueE,
bcwisemccluster::TrueE<bcwisemccluster::StoredTrueE>);

} // namespace o2::aod
Expand Down
82 changes: 51 additions & 31 deletions PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -16,25 +16,25 @@
/// \author Nicolas Strangmann (nicolas.strangmann@cern.ch) - Goethe University Frankfurt
///

#include <limits>
#include <vector>
#include <map>
#include <string>
#include "PWGEM/PhotonMeson/DataModel/bcWiseTables.h"
#include "PWGEM/PhotonMeson/Utils/MCUtilities.h"
#include "PWGJE/DataModel/EMCALClusters.h"

#include "Framework/runDataProcessing.h"
#include "Framework/AnalysisTask.h"
#include "Common/CCDB/ctpRateFetcher.h"
#include "Common/DataModel/Centrality.h"
#include "Common/DataModel/EventSelection.h"

#include "CCDB/BasicCCDBManager.h"
#include "DataFormatsParameters/GRPLHCIFData.h"
#include "DetectorsBase/GeometryManager.h"
#include "EMCALBase/Geometry.h"
#include "Framework/AnalysisTask.h"
#include "Framework/runDataProcessing.h"

#include "Common/DataModel/EventSelection.h"
#include "Common/DataModel/Centrality.h"
#include "Common/CCDB/ctpRateFetcher.h"
#include "CCDB/BasicCCDBManager.h"
#include "DataFormatsParameters/GRPLHCIFData.h"
#include "PWGJE/DataModel/EMCALClusters.h"
#include "PWGEM/PhotonMeson/Utils/MCUtilities.h"
#include "PWGEM/PhotonMeson/DataModel/bcWiseTables.h"
#include <limits>
#include <map>
#include <string>
#include <vector>

using namespace o2;
using namespace o2::aod::emdownscaling;
Expand All @@ -56,6 +56,7 @@
Produces<aod::BCWiseClusters> clusterTable;
Produces<aod::BCWiseCollisions> collisionTable;
Produces<aod::BCWiseMCPi0s> mcpi0Table;
Produces<aod::BCWiseMCEtas> mcetaTable;
Produces<aod::BCWiseMCClusters> mcclusterTable;

PresliceUnsorted<MyCollisions> perFoundBC = aod::evsel::foundBCId;
Expand All @@ -70,7 +71,7 @@
Configurable<float> cfgMinTime{"cfgMinTime", -25, "Minimum time of selected clusters (ns)"};
Configurable<float> cfgMaxTime{"cfgMaxTime", 25, "Maximum time of selected clusters (ns)"};
Configurable<float> cfgRapidityCut{"cfgRapidityCut", 0.8f, "Maximum absolute rapidity of counted generated particles"};
Configurable<float> cfgMinPtGenPi0{"cfgMinPtGenPi0", 0., "Minimum pT for stored generated pi0s (reduce disk space of derived data)"};
Configurable<float> cfgMinPtGen{"cfgMinPtGen", 0., "Minimum pT for stored generated mesons (reduce disk space of derived data)"};

Configurable<bool> cfgRequirekTVXinEMC{"cfgRequirekTVXinEMC", false, "Only store kTVXinEMC triggered BCs"};
Configurable<bool> cfgRequireGoodRCTQuality{"cfgRequireGoodRCTQuality", false, "Only store BCs with good quality of T0 and EMC in RCT"};
Expand All @@ -90,6 +91,7 @@
HistogramRegistry mHistManager{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false};

std::map<int32_t, int32_t> fMapPi0Index; // Map to connect the MC index of the pi0 to the one saved in the derived table
std::map<int32_t, int32_t> fMapEtaIndex; // Map to connect the MC index of the eta to the one saved in the derived table

void init(o2::framework::InitContext&)
{
Expand All @@ -108,7 +110,7 @@
LOG(info) << "| Shape cut: " << cfgMinM02 << " < M02 < " << cfgMaxM02;
LOG(info) << "| Energy cut: " << cfgMinClusterEnergy << " < E < " << cfgMaxClusterEnergy;
LOG(info) << "| Rapidity cut: |y| < " << cfgRapidityCut;
LOG(info) << "| Min gen pi0 pt: pT > " << cfgMinPtGenPi0;
LOG(info) << "| Min gen pt: pT > " << cfgMinPtGen;

o2::emcal::Geometry::GetInstanceFromRunNumber(300000);
if (cfgRequireGoodRCTQuality)
Expand Down Expand Up @@ -154,31 +156,42 @@
{
for (const auto& cluster : clusters) {
float clusterInducerEnergy = 0.;
int32_t pi0MCIndex = -1;
int32_t mesonMCIndex = -1;
if (cluster.amplitudeA().size() > 0) {
int clusterInducerId = cluster.mcParticleIds()[0];
auto clusterInducer = mcParticles.iteratorAt(clusterInducerId);
clusterInducerEnergy = clusterInducer.e();
int daughterId = aod::pwgem::photonmeson::utils::mcutil::FindMotherInChain(clusterInducer, mcParticles, std::vector<int>{111});
int daughterId = aod::pwgem::photonmeson::utils::mcutil::FindMotherInChain(clusterInducer, mcParticles, std::vector<int>{111, 221});
if (daughterId > 0) {
pi0MCIndex = mcParticles.iteratorAt(daughterId).mothersIds()[0];
if (mcParticles.iteratorAt(pi0MCIndex).pt() < cfgMinPtGenPi0)
pi0MCIndex = -1;
mesonMCIndex = mcParticles.iteratorAt(daughterId).mothersIds()[0];
if (mcParticles.iteratorAt(mesonMCIndex).pt() < cfgMinPtGen)
mesonMCIndex = -1;
}
}
if (pi0MCIndex >= 0) {
if (fMapPi0Index.find(pi0MCIndex) != fMapPi0Index.end()) // Some pi0s might not be found (not gg decay or too large y)
pi0MCIndex = fMapPi0Index[pi0MCIndex]; // If pi0 was stored in table, change index from the MC index to the pi0 index from this task
else // If pi0 was not stored, treat photon as if not from pi0
pi0MCIndex = -1;
bool isEta = false;
if (mesonMCIndex >= 0) {
if (mcParticles.iteratorAt(mesonMCIndex).pdgCode() == 111) {

Check failure on line 173 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number]

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.

Check failure on line 173 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
if (fMapPi0Index.find(mesonMCIndex) != fMapPi0Index.end()) // Some pi0s might not be found (not gg decay or too large y)
mesonMCIndex = fMapPi0Index[mesonMCIndex]; // If pi0 was stored in table, change index from the MC index to the pi0 index from this task
else // If pi0 was not stored, treat photon as if not from pi0
mesonMCIndex = -1;
} else if (mcParticles.iteratorAt(mesonMCIndex).pdgCode() == 221) {

Check failure on line 178 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
isEta = true;
if (fMapEtaIndex.find(mesonMCIndex) != fMapEtaIndex.end()) // Some etas might not be found (not gg decay or too large y)
mesonMCIndex = fMapEtaIndex[mesonMCIndex]; // If eta was stored in table, change index from the MC index to the eta index from this task
else // If eta was not stored, treat photon as if not from eta
mesonMCIndex = -1;
} else {
mesonMCIndex = -1; // Not a pi0 or eta
}
}
mcclusterTable(bcID, pi0MCIndex, convertForStorage<uint16_t>(clusterInducerEnergy, kEnergy));
mcclusterTable(bcID, mesonMCIndex, isEta, convertForStorage<uint16_t>(clusterInducerEnergy, kEnergy));
}
}

bool isBCSelected(const auto& bc)
{
if (cfgRequirekTVXinEMC && !bc.selection_bit(aod::evsel::kIsTriggerTVX))
if (cfgRequirekTVXinEMC && !bc.alias_bit(kTVXinEMC))
return false;
if (cfgRequireGoodRCTQuality && !isFT0EMCGoodRCTChecker(bc))
return false;
Expand Down Expand Up @@ -253,7 +266,7 @@
if (daughtersIds.size() != 2)
return false;
for (const auto& daughterId : daughtersIds) {
if (mcParticles.iteratorAt(daughterId).pdgCode() != 22)

Check failure on line 269 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
return false;
}
return true;
Expand All @@ -266,7 +279,7 @@
if (daughtersIds.size() != 2)
return false;
for (const auto& daughterId : daughtersIds) {
if (mcParticles.iteratorAt(daughterId).pdgCode() != 22)

Check failure on line 282 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
return false;
int iCellID = -1;
try {
Expand Down Expand Up @@ -319,12 +332,18 @@
auto mcParticlesInColl = mcParticles.sliceBy(perMcCollision, mcCollision.globalIndex());
mHistManager.fill(HIST("CentralityVsGenMultiplicity"), bc.centFT0M(), mcParticlesInColl.size());
for (const auto& mcParticle : mcParticlesInColl) {
if (mcParticle.pdgCode() != 111 || std::abs(mcParticle.y()) > cfgRapidityCut || !isGammaGammaDecay(mcParticle, mcParticles) || mcParticle.pt() < cfgMinPtGenPi0)
if (std::abs(mcParticle.y()) > cfgRapidityCut || !isGammaGammaDecay(mcParticle, mcParticles) || mcParticle.pt() < cfgMinPtGen)
continue;
bool isPrimary = mcParticle.isPhysicalPrimary() || mcParticle.producedByGenerator();
bool isFromWD = (aod::pwgem::photonmeson::utils::mcutil::IsFromWD(mcCollision, mcParticle, mcParticles)) > 0;
mcpi0Table(bc.globalIndex(), convertForStorage<uint16_t>(mcParticle.pt(), kpT), isAccepted(mcParticle, mcParticles), isPrimary, isFromWD);
fMapPi0Index[mcParticle.globalIndex()] = static_cast<int32_t>(mcpi0Table.lastIndex());

if (mcParticle.pdgCode() == 111) {

Check failure on line 340 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
mcpi0Table(bc.globalIndex(), convertForStorage<uint16_t>(mcParticle.pt(), kpT), isAccepted(mcParticle, mcParticles), isPrimary, isFromWD);
fMapPi0Index[mcParticle.globalIndex()] = static_cast<int32_t>(mcpi0Table.lastIndex());
} else if (mcParticle.pdgCode() == 221) {

Check failure on line 343 in PWGEM/PhotonMeson/TableProducer/bcWiseClusterSkimmer.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[pdg/explicit-code]

Avoid hard-coded PDG codes. Use named values from PDG_t or o2::constants::physics::Pdg instead.
mcetaTable(bc.globalIndex(), convertForStorage<uint16_t>(mcParticle.pt(), kpT), isAccepted(mcParticle, mcParticles), isPrimary, isFromWD);
fMapEtaIndex[mcParticle.globalIndex()] = static_cast<int32_t>(mcetaTable.lastIndex());
}
}
}

Expand All @@ -338,6 +357,7 @@
processClusterMCInfo(clustersInBC, bc.globalIndex(), mcParticles);
}
fMapPi0Index.clear();
fMapEtaIndex.clear();
}
}
PROCESS_SWITCH(bcWiseClusterSkimmer, processMC, "Run skimming for MC", false);
Expand Down
Loading
Loading