Skip to content

Commit

Permalink
[PWGEM/Dilepton] update table producers (AliceO2Group#9127)
Browse files Browse the repository at this point in the history
  • Loading branch information
@dsekihat
dsekihat authored Dec 25, 2024
1 parent 31d5088 commit 5880c77
Show file tree
Hide file tree
Showing 4 changed files with 95 additions and 206 deletions.
118 changes: 0 additions & 118 deletions PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -329,128 +329,10 @@ struct AssociateDileptonToEMEvent {
PROCESS_SWITCH(AssociateDileptonToEMEvent, processFwdMuon, "process forward muon indexing", false);
PROCESS_SWITCH(AssociateDileptonToEMEvent, processDummy, "process dummy", true);
};
struct EMEventPropertyTask {
SliceCache cache;
Preslice<aod::Tracks> perCollision = aod::track::collisionId;
using Run3Tracks = soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA>;

Configurable<bool> fillQAHistogram{"fillQAHistogram", false, "flag to fill QA histograms"};

Produces<aod::EMEventsProperty> evprop;
struct : ConfigurableGroup {
std::string prefix = "spherocity_cutgroup";
Configurable<bool> require_isPVContributor{"require_isPVContributor", false, "require tracks to be PV contributors"};
Configurable<int> min_ntrack{"min_ntrack", 3, "min. number of tracks"};
Configurable<float> min_pt{"min_pt", 0.15, "min. pT of track in GeV/c"};
Configurable<float> min_eta{"min_eta", -0.8, "min. eta of track"};
Configurable<float> max_eta{"max_eta", +0.8, "max. eta of track"};
Configurable<float> max_dcaxy{"max_dcaxy", 2.4, "max. DCAxy of track in cm"};
Configurable<float> max_dcaz{"max_dcaz", 3.2, "max. DCAz of track in cm"};
Configurable<int> min_ncluster_tpc{"min_ncluster_tpc", 50, "min. number of TPC clusters"};
Configurable<float> max_chi2tpc{"max_chi2tpc", 4.0, "max. chi2/ncls TPC"};
} spherocity_cuts;

HistogramRegistry fRegistry{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false};
void init(InitContext&)
{
if (fillQAHistogram) {
fRegistry.add("Spherocity/hPt", "pT;p_{T} (GeV/c)", kTH1F, {{200, 0.0f, 10}}, false);
fRegistry.add("Spherocity/hEtaPhi", "#eta vs. #varphi;#varphi (rad.);#eta", kTH2F, {{180, 0, 2 * M_PI}, {40, -1.0f, 1.0f}}, false);
fRegistry.add("Spherocity/hSpherocity_ptweighted", "spherocity;Number of used tracks;spherocity", kTH2F, {{101, -0.5, 100.5}, {100, 0.0f, 1}}, false);
fRegistry.add("Spherocity/hSpherocity_ptunweighted", "spherocity;Number of used tracks;spherocity", kTH2F, {{101, -0.5, 100.5}, {100, 0.0f, 1}}, false);
}
}

template <typename TTracks>
int getSpherocity(TTracks const& tracks, float& spherocity_ptweighted, float& spherocity_ptunweighted)
{
// Reference for spherocity : https://arxiv.org/pdf/1905.07208, https://arxiv.org/abs/2310.20406, https://arxiv.org/abs/1205.3963
int used_ntrack_spherocity = 0;
float sum_pt = 0.f, sum_ntrack = 0.f;

for (auto const& track : tracks) {
if (spherocity_cuts.require_isPVContributor && !track.isPVContributor()) {
continue;
}
if (track.tpcNClsFound() < spherocity_cuts.min_ncluster_tpc) {
continue;
}

sum_pt += track.pt();
sum_ntrack += 1.f;

if (fillQAHistogram) {
fRegistry.fill(HIST("Spherocity/hPt"), track.pt());
fRegistry.fill(HIST("Spherocity/hEtaPhi"), track.phi(), track.eta());
}
used_ntrack_spherocity++;
} // end of track loop per collision

float tempSph = 1.f, tempSph_pt1 = 1.f;
for (int i = 0; i < 360 / 0.1; i++) {
float nx = std::cos(M_PI / 180.f * i * 0.1);
float ny = std::sin(M_PI / 180.f * i * 0.1);
float sum_crossprod = 0.f, sum_crossprod_pt1 = 0.f;
for (auto const& track : tracks) {
if (spherocity_cuts.require_isPVContributor && !track.isPVContributor()) {
continue;
}
if (track.tpcNClsFound() < spherocity_cuts.min_ncluster_tpc) {
continue;
}
float px = track.px();
float py = track.py();
sum_crossprod += abs(px * ny - py * nx);
sum_crossprod_pt1 += abs(std::cos(track.phi()) * ny - std::sin(track.phi()) * nx);
}
float sph = std::pow(sum_crossprod / sum_pt, 2);
float sph_pt1 = std::pow(sum_crossprod_pt1 / sum_ntrack, 2);
if (sph < tempSph) {
tempSph = sph;
}
if (sph_pt1 < tempSph_pt1) {
tempSph_pt1 = sph_pt1;
}
} // end of track loop per collision
spherocity_ptweighted = std::pow(M_PI_2, 2) * tempSph;
spherocity_ptunweighted = std::pow(M_PI_2, 2) * tempSph_pt1;
if (used_ntrack_spherocity < spherocity_cuts.min_ntrack) {
spherocity_ptweighted = -1.f;
spherocity_ptunweighted = -1.f;
}
return used_ntrack_spherocity;
}

Partition<Run3Tracks> tracks_for_spherocity = spherocity_cuts.min_pt < aod::track::pt && spherocity_cuts.min_eta < o2::aod::track::eta && o2::aod::track::eta < spherocity_cuts.max_eta && nabs(o2::aod::track::dcaXY) < spherocity_cuts.max_dcaxy && nabs(o2::aod::track::dcaZ) < spherocity_cuts.max_dcaz && ncheckbit(aod::track::v001::detectorMap, (uint8_t)o2::aod::track::ITS) == true && ncheckbit(aod::track::v001::detectorMap, (uint8_t)o2::aod::track::TPC) == true && o2::aod::track::tpcChi2NCl < spherocity_cuts.max_chi2tpc; // ITS-TPC matched tracks

void processProp(aod::EMEvents const& collisions, Run3Tracks const&)
{
for (auto& collision : collisions) {
auto tracks_for_spherocity_per_collision = tracks_for_spherocity->sliceByCached(o2::aod::track::collisionId, collision.collisionId(), cache);
float spherocity_ptweighted = -1.f, spherocity_ptunweighted = -1.f;
int ntrack = getSpherocity(tracks_for_spherocity_per_collision, spherocity_ptweighted, spherocity_ptunweighted);
if (fillQAHistogram) {
fRegistry.fill(HIST("Spherocity/hSpherocity_ptweighted"), ntrack, spherocity_ptweighted);
fRegistry.fill(HIST("Spherocity/hSpherocity_ptunweighted"), ntrack, spherocity_ptunweighted);
}
evprop(spherocity_ptweighted, spherocity_ptunweighted, ntrack);
} // end of collision loop
}
PROCESS_SWITCH(EMEventPropertyTask, processProp, "process event property", true);

void processDummy(aod::EMEvents const& collisions)
{
for (int i = 0; i < collisions.size(); i++) {
evprop(-1.f, -1.f, 0);
} // end of collision loop
}
PROCESS_SWITCH(EMEventPropertyTask, processDummy, "process dummy", false);
};
WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
{
return WorkflowSpec{
adaptAnalysisTask<CreateEMEventDilepton>(cfgc, TaskName{"create-emevent-dilepton"}),
adaptAnalysisTask<AssociateDileptonToEMEvent>(cfgc, TaskName{"associate-dilepton-to-emevent"}),
adaptAnalysisTask<EMEventPropertyTask>(cfgc, TaskName{"emevent-property"}),
};
}
37 changes: 12 additions & 25 deletions PWGEM/Dilepton/TableProducer/filterDielectronEvent.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -74,10 +74,10 @@ struct filterDielectronEvent {
Configurable<int> min_ncluster_tpc{"min_ncluster_tpc", 0, "min ncluster tpc"};
Configurable<int> mincrossedrows{"mincrossedrows", 80, "min. crossed rows"};
Configurable<float> min_tpc_cr_findable_ratio{"min_tpc_cr_findable_ratio", 0.8, "min. TPC Ncr/Nf ratio"};
Configurable<float> max_mean_its_cluster_size{"max_mean_its_cluster_size", 16.f, "max. <ITS cluster size> x cos(lambda)"};
Configurable<float> max_p_for_its_cluster_size{"max_p_for_its_cluster_size", 0, "its cluster size cut is applied below this p"};
Configurable<float> max_pin_for_pion_rejection{"max_pin_for_pion_rejection", -1, "pion rejection is applied below this pin"};
Configurable<int> minitsncls{"minitsncls", 4, "min. number of ITS clusters"};
Configurable<float> max_pin_for_pion_rejection{"max_pin_for_pion_rejection", 1e+10, "pion rejection is applied below this pin"};
Configurable<float> max_frac_shared_clusters_tpc{"max_frac_shared_clusters_tpc", 999.f, "max fraction of shared clusters in TPC"};
Configurable<int> min_ncluster_its{"min_ncluster_its", 4, "min ncluster its"};
Configurable<int> min_ncluster_itsib{"min_ncluster_itsib", 1, "min ncluster itsib"};
Configurable<float> maxchi2tpc{"maxchi2tpc", 5.0, "max. chi2/NclsTPC"};
Configurable<float> maxchi2its{"maxchi2its", 6.0, "max. chi2/NclsITS"};
Configurable<float> minpt{"minpt", 0.15, "min pt for track"};
Expand All @@ -102,8 +102,6 @@ struct filterDielectronEvent {

HistogramRegistry fRegistry{"output", {}, OutputObjHandlingPolicy::AnalysisObject, false, false};

std::pair<int8_t, std::set<uint8_t>> itsRequirement = {1, {0, 1, 2}}; // any hits on 3 ITS ib layers.

int mRunNumber;
float d_bz;
Service<o2::ccdb::BasicCCDBManager> ccdb;
Expand Down Expand Up @@ -220,25 +218,10 @@ struct filterDielectronEvent {
if (!track.hasITS() || !track.hasTPC()) {
return false;
}
if (track.itsNCls() < minitsncls) {
if (track.itsNCls() < min_ncluster_its) {
return false;
}

auto hits = std::count_if(itsRequirement.second.begin(), itsRequirement.second.end(), [&](auto&& requiredLayer) { return track.itsClusterMap() & (1 << requiredLayer); });
if (hits < itsRequirement.first) {
return false;
}

uint32_t itsClusterSizes = track.itsClusterSizes();
int total_cluster_size = 0, nl = 0;
for (unsigned int layer = 0; layer < 7; layer++) {
int cluster_size_per_layer = (itsClusterSizes >> (layer * 4)) & 0xf;
if (cluster_size_per_layer > 0) {
nl++;
}
total_cluster_size += cluster_size_per_layer;
}
if (static_cast<float>(total_cluster_size) / static_cast<float>(nl) * std::cos(std::atan(track.tgl())) > max_mean_its_cluster_size && track.p() < max_p_for_its_cluster_size) {
if (track.itsNClsInnerBarrel() < min_ncluster_itsib) {
return false;
}

Expand All @@ -254,6 +237,10 @@ struct filterDielectronEvent {
return false;
}

if (track.tpcFractionSharedCls() > max_frac_shared_clusters_tpc) {
return false;
}

gpu::gpustd::array<float, 2> dcaInfo;
auto track_par_cov_recalc = getTrackParCov(track);
track_par_cov_recalc.setPID(o2::track::PID::Electron);
Expand Down Expand Up @@ -298,7 +285,7 @@ struct filterDielectronEvent {
if (track.tpcNSigmaEl() < minTPCNsigmaEl || maxTPCNsigmaEl < track.tpcNSigmaEl()) {
return false;
}
if (minTPCNsigmaPi < track.tpcNSigmaPi() && track.tpcNSigmaPi() < maxTPCNsigmaPi) {
if (minTPCNsigmaPi < track.tpcNSigmaPi() && track.tpcNSigmaPi() < maxTPCNsigmaPi && track.tpcInnerParam() < max_pin_for_pion_rejection) {
return false;
}
if (minTPCNsigmaKa < track.tpcNSigmaKa() && track.tpcNSigmaKa() < maxTPCNsigmaKa) {
Expand All @@ -316,7 +303,7 @@ struct filterDielectronEvent {
template <typename TTrack>
bool isElectron_TOFreq(TTrack const& track)
{
if (minTPCNsigmaPi < track.tpcNSigmaPi() && track.tpcNSigmaPi() < maxTPCNsigmaPi) {
if (minTPCNsigmaPi < track.tpcNSigmaPi() && track.tpcNSigmaPi() < maxTPCNsigmaPi && track.tpcInnerParam() < max_pin_for_pion_rejection) {
return false;
}
return minTPCNsigmaEl < track.tpcNSigmaEl() && track.tpcNSigmaEl() < maxTPCNsigmaEl && fabs(track.tofNSigmaEl()) < maxTOFNsigmaEl;
Expand Down
Loading

0 comments on commit 5880c77

Please sign in to comment.