diff --git a/PWGCF/EbyEFluctuations/Tasks/netchargeFluctuations.cxx b/PWGCF/EbyEFluctuations/Tasks/netchargeFluctuations.cxx index cf7c6a46628..be1644c1387 100644 --- a/PWGCF/EbyEFluctuations/Tasks/netchargeFluctuations.cxx +++ b/PWGCF/EbyEFluctuations/Tasks/netchargeFluctuations.cxx @@ -15,35 +15,37 @@ /// For RUN-3 /// /// \author Nida Malik -#include // Include for std::vector +#include "PWGCF/Core/CorrelationContainer.h" +#include "PWGCF/Core/PairCuts.h" -#include "Framework/AnalysisTask.h" -#include "Framework/runDataProcessing.h" +#include "Common/CCDB/EventSelectionParams.h" +#include "Common/CCDB/TriggerAliases.h" +#include "Common/Core/TrackSelection.h" +#include "Common/Core/trackUtilities.h" +#include "Common/DataModel/Centrality.h" #include "Common/DataModel/EventSelection.h" +#include "Common/DataModel/FT0Corrected.h" #include "Common/DataModel/Multiplicity.h" #include "Common/DataModel/PIDResponse.h" -#include "Common/Core/trackUtilities.h" -#include "Common/CCDB/EventSelectionParams.h" -#include "Common/Core/TrackSelection.h" #include "Common/DataModel/TrackSelectionTables.h" -#include "Common/DataModel/Centrality.h" + +#include "CCDB/BasicCCDBManager.h" #include "CommonConstants/MathConstants.h" -#include "Common/DataModel/FT0Corrected.h" -#include "Framework/AnalysisDataModel.h" +#include "CommonConstants/PhysicsConstants.h" #include "Framework/ASoAHelpers.h" -#include "Framework/RunningWorkflowInfo.h" -#include "PWGCF/Core/CorrelationContainer.h" -#include "PWGCF/Core/PairCuts.h" -#include "Common/CCDB/TriggerAliases.h" +#include "Framework/AnalysisDataModel.h" +#include "Framework/AnalysisTask.h" #include "Framework/HistogramRegistry.h" #include "Framework/O2DatabasePDGPlugin.h" -#include "CommonConstants/PhysicsConstants.h" -#include "Common/CCDB/TriggerAliases.h" -#include "CCDB/BasicCCDBManager.h" +#include "Framework/RunningWorkflowInfo.h" +#include "Framework/runDataProcessing.h" + #include "TProfile.h" #include "TProfile2D.h" #include "TRandom3.h" +#include // Include for std::vector + using namespace o2; using namespace o2::framework; using namespace o2::framework::expressions; @@ -103,8 +105,6 @@ struct NetchargeFluctuations { Configurable cfgNSubsample{"cfgNSubsample", 30, "Number of subsamples for Error"}; Configurable deltaEta{"deltaEta", 8, "Delta eta bin count"}; Configurable threshold{"threshold", 1e-6, "Delta eta bin count"}; - - // Event selections Configurable cSel8Trig{"cSel8Trig", true, "Sel8 (T0A + T0C) Selection Run3"}; // sel8 @@ -116,13 +116,13 @@ struct NetchargeFluctuations { Configurable cPileupReject{"cPileupReject", false, "Pileup rejection"}; // pileup Configurable cZVtxTimeDiff{"cZVtxTimeDiff", false, "z-vtx time diff selection"}; // pileup Configurable cfgUseGoodItsLayerAllCut{"cfgUseGoodItsLayerAllCut", false, "Good ITS Layers All"}; // pileup - Configurable cDcaXy{"cDcaXy", false,"Dca XY cut"}; - Configurable cDcaZ{"cDcaZ", false,"Dca Z cut"}; - Configurable cTpcCr{"cTpcCr", false,"tpc crossrows"}; - Configurable cItsChi{"cItsChi", false,"ITS chi"}; - Configurable cTpcChi{"cTpcChi", false,"TPC chi"}; - Configurable cCorrectedReco{"cCorrectedReco", true,"Corrected Reco"}; - Configurable cCorrectedData{"cCorrectedData", false,"Corrected Data"}; + Configurable cDcaXy{"cDcaXy", false, "Dca XY cut"}; + Configurable cDcaZ{"cDcaZ", false, "Dca Z cut"}; + Configurable cTpcCr{"cTpcCr", false, "tpc crossrows"}; + Configurable cItsChi{"cItsChi", false, "ITS chi"}; + Configurable cTpcChi{"cTpcChi", false, "TPC chi"}; + Configurable cCorrectedReco{"cCorrectedReco", true, "Corrected Reco"}; + Configurable cCorrectedData{"cCorrectedData", false, "Corrected Data"}; // CCDB efficiency histograms TH2D* efficiency = nullptr; @@ -147,9 +147,8 @@ struct NetchargeFluctuations { const AxisSpec nch1Axis = {1500, 0, 1500, "Nch"}; const AxisSpec nchpAxis = {50000, 0, 50000, "Nch"}; - - std::vector centBining = {0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90,100}; - AxisSpec cent1Axis = {centBining, "Multiplicity percentile from FT0M (%)"}; + std::vector centBining = {0, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100}; + AxisSpec cent1Axis = {centBining, "Multiplicity percentile from FT0M (%)"}; auto noSubsample = static_cast(cfgNSubsample); float maxSubsample = 1.0 * noSubsample; @@ -186,7 +185,7 @@ struct NetchargeFluctuations { histogramRegistry.add("gen/hEta_cent", "", kTH2F, {centAxis, etaAxis}); histogramRegistry.add("gen/hSign", "", kTH1F, {signAxis}); histogramRegistry.add("gen/hPt", "", kTH1F, {ptAxis}); - histogramRegistry.add("gen/hPt_cent", "", kTH2F, {centAxis,ptAxis}); + histogramRegistry.add("gen/hPt_cent", "", kTH2F, {centAxis, ptAxis}); histogramRegistry.add("gen/nch", "", kTH1F, {nchAxis}); histogramRegistry.add("mult_dist/nch", "", kTH1D, {nchAxis}); @@ -241,30 +240,30 @@ struct NetchargeFluctuations { histogramRegistry.add("subsample/neg_sq", "", kTProfile2D, {cent1Axis, subsampleAxis}); histogramRegistry.add("subsample/posneg", "", kTProfile2D, {cent1Axis, subsampleAxis}); - - if (cfgLoadEff) { + if (cfgLoadEff) { ccdb->setURL(cfgUrlCCDB.value); ccdb->setCaching(true); ccdb->setLocalObjectValidityChecking(); - //ccdb->setCreatedNotAfter(ccdbNoLaterThan.value); - //LOGF(info, "Getting object %s", ccdbPath.value.data()); - + // ccdb->setCreatedNotAfter(ccdbNoLaterThan.value); + // LOGF(info, "Getting object %s", ccdbPath.value.data()); + TList* list = ccdb->getForTimeStamp(cfgPathCCDB.value, -1); efficiency = reinterpret_cast(list->FindObject("efficiency_Run3")); if (!efficiency) { - LOGF(info, "FATAL!! Could not find required histograms in CCDB");} - } - + LOGF(info, "FATAL!! Could not find required histograms in CCDB"); + } + } } - -template + + template bool selCollision(C const& coll, float& cent, float& mult) { - if (std::abs(coll.posZ()) > vertexZcut) return false; + if (std::abs(coll.posZ()) > vertexZcut) + return false; - if (run == kRun3) { + if (run == kRun3) { if (cSel8Trig && !coll.sel8()) { return false; } // require min bias trigger @@ -281,362 +280,393 @@ template mult = coll.multFV0M(); // multiplicity for run2 } - if (cNoItsROBorder && !coll.selection_bit(aod::evsel::kNoITSROFrameBorder)) return false; - if (cTFBorder && !coll.selection_bit(aod::evsel::kNoTimeFrameBorder)) return false; - if (cPileupReject && !coll.selection_bit(aod::evsel::kNoSameBunchPileup)) return false; - if (cZVtxTimeDiff && !coll.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) return false; - if (cItsTpcVtx && !coll.selection_bit(aod::evsel::kIsVertexITSTPC)) return false; - if (cfgUseGoodItsLayerAllCut && !(collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll))) return false; - + if (cNoItsROBorder && !coll.selection_bit(aod::evsel::kNoITSROFrameBorder)) + return false; + if (cTFBorder && !coll.selection_bit(aod::evsel::kNoTimeFrameBorder)) + return false; + if (cPileupReject && !coll.selection_bit(aod::evsel::kNoSameBunchPileup)) + return false; + if (cZVtxTimeDiff && !coll.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV)) + return false; + if (cItsTpcVtx && !coll.selection_bit(aod::evsel::kIsVertexITSTPC)) + return false; + if (cfgUseGoodItsLayerAllCut && !(collision.selection_bit(o2::aod::evsel::kIsGoodITSLayersAll))) + return false; + return true; } -template -void fillBeforeQA(T const& track) -{ - histogramRegistry.fill(HIST("data/hTPCchi2perCluster_before"), track.tpcChi2NCl()); - histogramRegistry.fill(HIST("data/hITSchi2perCluster_before"), track.itsChi2NCl()); - histogramRegistry.fill(HIST("data/hTPCCrossedrows_before"), track.tpcNClsCrossedRows()); - histogramRegistry.fill(HIST("data/hDcaXY_before"), track.dcaXY()); - histogramRegistry.fill(HIST("data/hDcaZ_before"), track.dcaZ()); - histogramRegistry.fill(HIST("data/hPtDcaXY_before"), track.pt(), track.dcaXY()); - histogramRegistry.fill(HIST("data/hPtDcaZ_before"), track.pt(), track.dcaZ()); -} - -template -void fillAfterQA(T const& track) -{ - histogramRegistry.fill(HIST("data/hDcaXY_after"), track.dcaXY()); - histogramRegistry.fill(HIST("data/hDcaZ_after"), track.dcaZ()); - histogramRegistry.fill(HIST("data/hPt"), track.pt()); - histogramRegistry.fill(HIST("data/hEta"), track.eta()); - histogramRegistry.fill(HIST("data/hPt_eta"), track.pt(), track.eta()); - histogramRegistry.fill(HIST("data/hPtDcaXY_after"), track.pt(), track.dcaXY()); - histogramRegistry.fill(HIST("data/hPtDcaZ_after"), track.pt(), track.dcaZ()); - histogramRegistry.fill(HIST("data/hTPCCrossedrows_after"), track.tpcNClsCrossedRows()); - histogramRegistry.fill(HIST("data/hTPCchi2perCluster_after"), track.tpcChi2NCl()); - histogramRegistry.fill(HIST("data/hITSchi2perCluster_after"), track.itsChi2NCl()); + template + void fillBeforeQA(T const& track) + { + histogramRegistry.fill(HIST("data/hTPCchi2perCluster_before"), track.tpcChi2NCl()); + histogramRegistry.fill(HIST("data/hITSchi2perCluster_before"), track.itsChi2NCl()); + histogramRegistry.fill(HIST("data/hTPCCrossedrows_before"), track.tpcNClsCrossedRows()); + histogramRegistry.fill(HIST("data/hDcaXY_before"), track.dcaXY()); + histogramRegistry.fill(HIST("data/hDcaZ_before"), track.dcaZ()); + histogramRegistry.fill(HIST("data/hPtDcaXY_before"), track.pt(), track.dcaXY()); + histogramRegistry.fill(HIST("data/hPtDcaZ_before"), track.pt(), track.dcaZ()); } + template + void fillAfterQA(T const& track) + { + histogramRegistry.fill(HIST("data/hDcaXY_after"), track.dcaXY()); + histogramRegistry.fill(HIST("data/hDcaZ_after"), track.dcaZ()); + histogramRegistry.fill(HIST("data/hPt"), track.pt()); + histogramRegistry.fill(HIST("data/hEta"), track.eta()); + histogramRegistry.fill(HIST("data/hPt_eta"), track.pt(), track.eta()); + histogramRegistry.fill(HIST("data/hPtDcaXY_after"), track.pt(), track.dcaXY()); + histogramRegistry.fill(HIST("data/hPtDcaZ_after"), track.pt(), track.dcaZ()); + histogramRegistry.fill(HIST("data/hTPCCrossedrows_after"), track.tpcNClsCrossedRows()); + histogramRegistry.fill(HIST("data/hTPCchi2perCluster_after"), track.tpcChi2NCl()); + histogramRegistry.fill(HIST("data/hITSchi2perCluster_after"), track.itsChi2NCl()); + } -template + template bool selTrack(T const& track) { - if (!track.isGlobalTrack()) return false ; - if (std::fabs(track.eta()) >= etaCut) return false; - if (track.pt() <= ptMinCut || track.pt() >= ptMaxCut) return false; - if(track.sign() == 0) return false; - if (cDcaXy && std::fabs(track.dcaXY()) > dcaXYCut) return false; - if (cDcaZ && std::fabs(track.dcaZ()) > dcaZCut) return false ; - if (cTpcCr && track.tpcNClsCrossedRows() < tpcCrossCut) return false; - if (cItsChi && track.itsChi2NCl() > itsChiCut) return false; - if (cTpcChi && track.tpcChi2NCl() > tpcChiCut) return false; - -return true; - } + if (!track.isGlobalTrack()) + return false; + if (std::fabs(track.eta()) >= etaCut) + return false; + if (track.pt() <= ptMinCut || track.pt() >= ptMaxCut) + return false; + if (track.sign() == 0) + return false; + if (cDcaXy && std::fabs(track.dcaXY()) > dcaXYCut) + return false; + if (cDcaZ && std::fabs(track.dcaZ()) > dcaZCut) + return false; + if (cTpcCr && track.tpcNClsCrossedRows() < tpcCrossCut) + return false; + if (cItsChi && track.itsChi2NCl() > itsChiCut) + return false; + if (cTpcChi && track.tpcChi2NCl() > tpcChiCut) + return false; -double getEfficiency(double pt, double eta, TH2D* hEff) -{ - if (!hEff) { - LOGF(error, "Efficiency histogram is null — check CCDB loading."); return 1e-6; + return true; } - int binX = hEff->GetXaxis()->FindBin(pt); - int binY = hEff->GetYaxis()->FindBin(eta); - if (binX < 1 || binX > hEff->GetNbinsX() || binY < 1 || binY > hEff->GetNbinsY()) { - LOGF(warn, "pt or eta out of histogram bounds: pt = %f, eta = %f", pt, eta); return 1e-6; + + double getEfficiency(double pt, double eta, TH2D* hEff) + { + if (!hEff) { + LOGF(error, "Efficiency histogram is null — check CCDB loading."); + return 1e-6; + } + int binX = hEff->GetXaxis()->FindBin(pt); + int binY = hEff->GetYaxis()->FindBin(eta); + if (binX < 1 || binX > hEff->GetNbinsX() || binY < 1 || binY > hEff->GetNbinsY()) { + LOGF(warn, "pt or eta out of histogram bounds: pt = %f, eta = %f", pt, eta); + return 1e-6; + } + double eff = hEff->GetBinContent(binX, binY); + return eff; } - double eff = hEff->GetBinContent(binX, binY); - return eff; -} -void fillHistograms( float nch, float cent, float fpos, float fneg, float posneg, float termp, float termn) -{ - histogramRegistry.fill(HIST("mult_dist/nch"), nch); - histogramRegistry.fill(HIST("mult_dist/nch_cent"), cent, nch); - histogramRegistry.fill(HIST("mult_dist/nch_pos"), fpos); - histogramRegistry.fill(HIST("mult_dist/nch_pos_cent"), cent, fpos); - histogramRegistry.fill(HIST("mult_dist/nch_neg"), fneg); - histogramRegistry.fill(HIST("mult_dist/nch_neg_cent"), cent, fneg); - histogramRegistry.fill(HIST("mult_dist/nch_negpos"), posneg); - histogramRegistry.fill(HIST("mult_dist/nch_negpos_cent"), cent, posneg); - - histogramRegistry.fill(HIST("cent/pos"), cent, fpos); - histogramRegistry.fill(HIST("cent/neg"), cent, fneg); - histogramRegistry.fill(HIST("cent/termp"), cent, termp); - histogramRegistry.fill(HIST("cent/termn"), cent, termn); - histogramRegistry.fill(HIST("cent/pos_sq"), cent, fpos * fpos); - histogramRegistry.fill(HIST("cent/neg_sq"), cent, fneg * fneg); - histogramRegistry.fill(HIST("cent/posneg"), cent, posneg); - - float lRandom = fRndm->Rndm(); + void fillHistograms(float nch, float cent, float fpos, float fneg, float posneg, float termp, float termn) + { + histogramRegistry.fill(HIST("mult_dist/nch"), nch); + histogramRegistry.fill(HIST("mult_dist/nch_cent"), cent, nch); + histogramRegistry.fill(HIST("mult_dist/nch_pos"), fpos); + histogramRegistry.fill(HIST("mult_dist/nch_pos_cent"), cent, fpos); + histogramRegistry.fill(HIST("mult_dist/nch_neg"), fneg); + histogramRegistry.fill(HIST("mult_dist/nch_neg_cent"), cent, fneg); + histogramRegistry.fill(HIST("mult_dist/nch_negpos"), posneg); + histogramRegistry.fill(HIST("mult_dist/nch_negpos_cent"), cent, posneg); + + histogramRegistry.fill(HIST("cent/pos"), cent, fpos); + histogramRegistry.fill(HIST("cent/neg"), cent, fneg); + histogramRegistry.fill(HIST("cent/termp"), cent, termp); + histogramRegistry.fill(HIST("cent/termn"), cent, termn); + histogramRegistry.fill(HIST("cent/pos_sq"), cent, fpos * fpos); + histogramRegistry.fill(HIST("cent/neg_sq"), cent, fneg * fneg); + histogramRegistry.fill(HIST("cent/posneg"), cent, posneg); + + float lRandom = fRndm->Rndm(); int sampleIndex = static_cast(cfgNSubsample * lRandom); - histogramRegistry.fill(HIST("subsample/pos"), cent,sampleIndex, fpos); - histogramRegistry.fill(HIST("subsample/neg"), cent,sampleIndex, fneg); - histogramRegistry.fill(HIST("subsample/termp"), cent,sampleIndex, termp); - histogramRegistry.fill(HIST("subsample/termn"), cent,sampleIndex, termn); - histogramRegistry.fill(HIST("subsample/pos_sq"), cent,sampleIndex, fpos *fpos); - histogramRegistry.fill(HIST("subsample/neg_sq"), cent,sampleIndex, fneg*fneg); - histogramRegistry.fill(HIST("subsample/posneg"), cent,sampleIndex, posneg); - -} + histogramRegistry.fill(HIST("subsample/pos"), cent, sampleIndex, fpos); + histogramRegistry.fill(HIST("subsample/neg"), cent, sampleIndex, fneg); + histogramRegistry.fill(HIST("subsample/termp"), cent, sampleIndex, termp); + histogramRegistry.fill(HIST("subsample/termn"), cent, sampleIndex, termn); + histogramRegistry.fill(HIST("subsample/pos_sq"), cent, sampleIndex, fpos * fpos); + histogramRegistry.fill(HIST("subsample/neg_sq"), cent, sampleIndex, fneg * fneg); + histogramRegistry.fill(HIST("subsample/posneg"), cent, sampleIndex, posneg); + } -template + template void calculationData(C const& coll, T const& tracks) { - float cent = -1, mult = -1; - histogramRegistry.fill(HIST("data/hVtxZ_before"), coll.posZ()); + float cent = -1, mult = -1; + histogramRegistry.fill(HIST("data/hVtxZ_before"), coll.posZ()); if (!selCollision(coll, cent, mult)) { - return; + return; } histogramRegistry.fill(HIST("data/hVtxZ_after"), coll.posZ()); histogramRegistry.fill(HIST("data/hCentrality"), cent); histogramRegistry.fill(HIST("data/hMultiplicity"), mult); - int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; int nch =0, nchCor = 0; + int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; + int nch = 0, nchCor = 0; double posWeight = 0, negWeight = 0; for (const auto& track : tracks) { - fillBeforeQA(track); - if (!selTrack(track)) continue; - nch +=1; - fillAfterQA(track); - histogramRegistry.fill(HIST("data/hEta_cent"), cent, track.eta()); - histogramRegistry.fill(HIST("data/hPt_cent"),cent, track.pt()); - - double eff = getEfficiency(track.pt(), track.eta(), efficiency); - if (eff < threshold) continue; - double weight = 1.0 / eff; - - histogramRegistry.fill(HIST("cor/hPt_cor"), track.pt(), weight); - histogramRegistry.fill(HIST("cor/hEta_cor"), track.eta(), weight); - - nchCor += weight; - if (track.sign() == 1) { - fpos += 1; - posWeight +=weight; - } - else if (track.sign() == -1) { - fneg += 1; - negWeight += weight; - } - - }//track + fillBeforeQA(track); + if (!selTrack(track)) + continue; + nch += 1; + fillAfterQA(track); + histogramRegistry.fill(HIST("data/hEta_cent"), cent, track.eta()); + histogramRegistry.fill(HIST("data/hPt_cent"), cent, track.pt()); + + double eff = getEfficiency(track.pt(), track.eta(), efficiency); + if (eff < threshold) + continue; + double weight = 1.0 / eff; + + histogramRegistry.fill(HIST("cor/hPt_cor"), track.pt(), weight); + histogramRegistry.fill(HIST("cor/hEta_cor"), track.eta(), weight); + + nchCor += weight; + if (track.sign() == 1) { + fpos += 1; + posWeight += weight; + } else if (track.sign() == -1) { + fneg += 1; + negWeight += weight; + } + + } // track termp = fpos * (fpos - 1); termn = fneg * (fneg - 1); - posneg = fpos * fneg; + posneg = fpos * fneg; histogramRegistry.fill(HIST("cor/nch_vs_nchCor"), nch, nchCor); histogramRegistry.fill(HIST("cor/nchCor"), nchCor); histogramRegistry.fill(HIST("cor/cent_nchCor"), cent, nchCor); - histogramRegistry.fill(HIST("cor/fpos_cent"), cent, posWeight); + histogramRegistry.fill(HIST("cor/fpos_cent"), cent, posWeight); histogramRegistry.fill(HIST("cor/fneg_cent"), cent, negWeight); - fillHistograms(nch, cent, fpos, fneg, posneg, termp, termn); - -template - void calculationMc(C const& coll, T const& inputTracks, M const& mcCollisions, P const& mcParticles) - { - (void) mcCollisions; - if (!coll.has_mcCollision()) { return;} - histogramRegistry.fill(HIST("gen/hVtxZ_before"), coll.mcCollision().posZ()); - float cent = -1, mult = -1; - histogramRegistry.fill(HIST("data/hVtxZ_before"), coll.posZ()); - if (!selCollision(coll, cent, mult)) {return;} - histogramRegistry.fill(HIST("data/hVtxZ_after"), coll.posZ()); - histogramRegistry.fill(HIST("data/hCentrality"), cent); - histogramRegistry.fill(HIST("data/hMultiplicity"), mult); - - int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; int nch =0, nchCor = 0; - double posRecWeight = 0, negRecWeight = 0; - - for (const auto& track : inputTracks) { - fillBeforeQA(track); - if (!selTrack(track)) continue; - nch +=1; - fillAfterQA(track); - histogramRegistry.fill(HIST("data/hEta_cent"), cent, track.eta()); - histogramRegistry.fill(HIST("data/hPt_cent"),cent, track.pt()); - - double eff = getEfficiency(track.pt(), track.eta(), efficiency); - if (eff < threshold) continue; - double weight = 1.0 / eff; - histogramRegistry.fill(HIST("cor/hPt_cor"), track.pt(), weight); - histogramRegistry.fill(HIST("cor/hEta_cor"), track.eta(), weight); - nchCor += weight; - - if (track.sign() == 1) { - fpos += 1; - posRecWeight += weight; - } - else if (track.sign() == -1) { - fneg += 1; - negRecWeight += weight; - } - }//track - termp = fpos * (fpos - 1); - - termn = fneg * (fneg - 1); - - posneg = fpos * fneg; - histogramRegistry.fill(HIST("cor/nch_vs_nchCor"), nch, nchCor); - histogramRegistry.fill(HIST("cor/nchCor"), nchCor); - histogramRegistry.fill(HIST("cor/cent_nchCor"), cent, nchCor); - histogramRegistry.fill(HIST("cor/fpos_cent"), cent, posRecWeight); - histogramRegistry.fill(HIST("cor/fneg_cent"), cent,negRecWeight); - fillHistograms(nch, cent, fpos, fneg, posneg, termp, termn); - - int posGen = 0, negGen = 0, posNegGen = 0, termNGen = 0, termPGen = 0, nchGen = 0; - - const auto& mccolgen = coll.template mcCollision_as(); - if (std::abs(mccolgen.posZ()) > vertexZcut) return; - const auto& mcpartgen = mcParticles.sliceByCached(aod::mcparticle::mcCollisionId, mccolgen.globalIndex(), cache); - histogramRegistry.fill(HIST("gen/hVtxZ_after"), mccolgen.posZ()); - for (const auto& mcpart : mcpartgen) { - if (!mcpart.isPhysicalPrimary()) continue; - int pid = mcpart.pdgCode(); - auto sign = 0; - auto* pd = pdgService->GetParticle(pid); - if (pd != nullptr) { sign = pd->Charge() / 3.; } - if (sign == 0) continue; - if (std::abs(pid) != kElectron && std::abs(pid) != kMuonMinus && std::abs(pid) != kPiPlus && std::abs(pid) != kKPlus && std::abs(pid) != kProton) continue; - if (std::fabs(mcpart.eta()) > etaCut) continue; - if ((mcpart.pt() <= ptMinCut) || (mcpart.pt() >= ptMaxCut)) continue; - histogramRegistry.fill(HIST("gen/hPt"), mcpart.pt()); - histogramRegistry.fill(HIST("gen/hPt_cent"),cent, mcpart.pt()); - histogramRegistry.fill(HIST("gen/hEta"), mcpart.eta()); - histogramRegistry.fill(HIST("gen/hEta_cent"),cent, mcpart.eta()); - histogramRegistry.fill(HIST("gen/hSign"), sign); - histogramRegistry.fill(HIST("gen/hPt_eta"), mcpart.pt(), mcpart.eta()); - nchGen +=1; - if (sign == 1) { + + template + void calculationMc(C const& coll, T const& inputTracks, M const& mcCollisions, P const& mcParticles) + { + (void)mcCollisions; + if (!coll.has_mcCollision()) { + return; + } + histogramRegistry.fill(HIST("gen/hVtxZ_before"), coll.mcCollision().posZ()); + float cent = -1, mult = -1; + histogramRegistry.fill(HIST("data/hVtxZ_before"), coll.posZ()); + if (!selCollision(coll, cent, mult)) { + return; + } + histogramRegistry.fill(HIST("data/hVtxZ_after"), coll.posZ()); + histogramRegistry.fill(HIST("data/hCentrality"), cent); + histogramRegistry.fill(HIST("data/hMultiplicity"), mult); + + int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; + int nch = 0, nchCor = 0; + double posRecWeight = 0, negRecWeight = 0; + + for (const auto& track : inputTracks) { + fillBeforeQA(track); + if (!selTrack(track)) + continue; + nch += 1; + fillAfterQA(track); + histogramRegistry.fill(HIST("data/hEta_cent"), cent, track.eta()); + histogramRegistry.fill(HIST("data/hPt_cent"), cent, track.pt()); + + double eff = getEfficiency(track.pt(), track.eta(), efficiency); + if (eff < threshold) + continue; + double weight = 1.0 / eff; + histogramRegistry.fill(HIST("cor/hPt_cor"), track.pt(), weight); + histogramRegistry.fill(HIST("cor/hEta_cor"), track.eta(), weight); + nchCor += weight; + + if (track.sign() == 1) { + fpos += 1; + posRecWeight += weight; + } else if (track.sign() == -1) { + fneg += 1; + negRecWeight += weight; + } + } // track + termp = fpos * (fpos - 1); + + termn = fneg * (fneg - 1); + + posneg = fpos * fneg; + histogramRegistry.fill(HIST("cor/nch_vs_nchCor"), nch, nchCor); + histogramRegistry.fill(HIST("cor/nchCor"), nchCor); + histogramRegistry.fill(HIST("cor/cent_nchCor"), cent, nchCor); + histogramRegistry.fill(HIST("cor/fpos_cent"), cent, posRecWeight); + histogramRegistry.fill(HIST("cor/fneg_cent"), cent, negRecWeight); + + fillHistograms(nch, cent, fpos, fneg, posneg, termp, termn); + + int posGen = 0, negGen = 0, posNegGen = 0, termNGen = 0, termPGen = 0, nchGen = 0; + + const auto& mccolgen = coll.template mcCollision_as(); + if (std::abs(mccolgen.posZ()) > vertexZcut) + return; + const auto& mcpartgen = mcParticles.sliceByCached(aod::mcparticle::mcCollisionId, mccolgen.globalIndex(), cache); + histogramRegistry.fill(HIST("gen/hVtxZ_after"), mccolgen.posZ()); + for (const auto& mcpart : mcpartgen) { + if (!mcpart.isPhysicalPrimary()) + continue; + int pid = mcpart.pdgCode(); + auto sign = 0; + auto* pd = pdgService->GetParticle(pid); + if (pd != nullptr) { + sign = pd->Charge() / 3.; + } + if (sign == 0) + continue; + if (std::abs(pid) != kElectron && std::abs(pid) != kMuonMinus && std::abs(pid) != kPiPlus && std::abs(pid) != kKPlus && std::abs(pid) != kProton) + continue; + if (std::fabs(mcpart.eta()) > etaCut) + continue; + if ((mcpart.pt() <= ptMinCut) || (mcpart.pt() >= ptMaxCut)) + continue; + histogramRegistry.fill(HIST("gen/hPt"), mcpart.pt()); + histogramRegistry.fill(HIST("gen/hPt_cent"), cent, mcpart.pt()); + histogramRegistry.fill(HIST("gen/hEta"), mcpart.eta()); + histogramRegistry.fill(HIST("gen/hEta_cent"), cent, mcpart.eta()); + histogramRegistry.fill(HIST("gen/hSign"), sign); + histogramRegistry.fill(HIST("gen/hPt_eta"), mcpart.pt(), mcpart.eta()); + nchGen += 1; + if (sign == 1) { posGen += 1; - } - if (sign == -1) { + } + if (sign == -1) { negGen += 1; - } - }//particle + } + } // particle termPGen = posGen * (posGen - 1); termNGen = negGen * (negGen - 1); posNegGen = posGen * negGen; - histogramRegistry.fill(HIST("cent/gen_pos"), cent, posGen); - histogramRegistry.fill(HIST("cent/gen_neg"), cent, negGen); - histogramRegistry.fill(HIST("cent/gen_termp"), cent, termPGen); - histogramRegistry.fill(HIST("cent/gen_termn"), cent, termNGen); - histogramRegistry.fill(HIST("cent/gen_pos_sq"), cent, posGen*posGen); - histogramRegistry.fill(HIST("cent/gen_neg_sq"), cent, negGen*negGen); - histogramRegistry.fill(HIST("cent/gen_posneg"), cent, posNegGen); - histogramRegistry.fill(HIST("cent/gen_nch"), cent, nchGen); - histogramRegistry.fill(HIST("gen/nch"), nchGen); - - }// void - - -template -void calculationDeltaEta(C const& coll, T const& tracks, float deta1, float deta2) -{ - float cent = -1, mult = -1; - if (!selCollision(coll, cent, mult)) return; - if (!(cent >= centMin && cent < centMax)) return; - histogramRegistry.fill(HIST("delta_eta/cent"), cent); - - int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; - for (const auto& track : tracks) { - if (!selTrack(track)) continue; - double eta = track.eta(); - if (eta < deta1 || eta > deta2) continue; - - histogramRegistry.fill(HIST("delta_eta/track_eta"), eta); - - if (track.sign() == 1) fpos++; - else if (track.sign() == -1) fneg++; - } - termp = fpos * (fpos - 1); - termn = fneg * (fneg - 1); - posneg = fpos * fneg; - - float deltaEtaWidth = deta2 - deta1 + 1e-5f; - - histogramRegistry.fill(HIST("delta_eta/pos"), deltaEtaWidth, fpos); - histogramRegistry.fill(HIST("delta_eta/neg"), deltaEtaWidth, fneg); - histogramRegistry.fill(HIST("delta_eta/termp"), deltaEtaWidth, termp); - histogramRegistry.fill(HIST("delta_eta/termn"), deltaEtaWidth, termn); - histogramRegistry.fill(HIST("delta_eta/pos_sq"), deltaEtaWidth, fpos * fpos); - histogramRegistry.fill(HIST("delta_eta/neg_sq"), deltaEtaWidth, fneg * fneg); - histogramRegistry.fill(HIST("delta_eta/posneg"), deltaEtaWidth, posneg); -} + histogramRegistry.fill(HIST("cent/gen_pos"), cent, posGen); + histogramRegistry.fill(HIST("cent/gen_neg"), cent, negGen); + histogramRegistry.fill(HIST("cent/gen_termp"), cent, termPGen); + histogramRegistry.fill(HIST("cent/gen_termn"), cent, termNGen); + histogramRegistry.fill(HIST("cent/gen_pos_sq"), cent, posGen * posGen); + histogramRegistry.fill(HIST("cent/gen_neg_sq"), cent, negGen * negGen); + histogramRegistry.fill(HIST("cent/gen_posneg"), cent, posNegGen); + histogramRegistry.fill(HIST("cent/gen_nch"), cent, nchGen); + histogramRegistry.fill(HIST("gen/nch"), nchGen); + + } // void + + template + void calculationDeltaEta(C const& coll, T const& tracks, float deta1, float deta2) + { + float cent = -1, mult = -1; + if (!selCollision(coll, cent, mult)) + return; + if (!(cent >= centMin && cent < centMax)) + return; + histogramRegistry.fill(HIST("delta_eta/cent"), cent); + + int fpos = 0, fneg = 0, posneg = 0, termn = 0, termp = 0; + for (const auto& track : tracks) { + if (!selTrack(track)) + continue; + double eta = track.eta(); + if (eta < deta1 || eta > deta2) + continue; + + histogramRegistry.fill(HIST("delta_eta/track_eta"), eta); + + if (track.sign() == 1) + fpos++; + else if (track.sign() == -1) + fneg++; + } + termp = fpos * (fpos - 1); + termn = fneg * (fneg - 1); + posneg = fpos * fneg; + + float deltaEtaWidth = deta2 - deta1 + 1e-5f; + + histogramRegistry.fill(HIST("delta_eta/pos"), deltaEtaWidth, fpos); + histogramRegistry.fill(HIST("delta_eta/neg"), deltaEtaWidth, fneg); + histogramRegistry.fill(HIST("delta_eta/termp"), deltaEtaWidth, termp); + histogramRegistry.fill(HIST("delta_eta/termn"), deltaEtaWidth, termn); + histogramRegistry.fill(HIST("delta_eta/pos_sq"), deltaEtaWidth, fpos * fpos); + histogramRegistry.fill(HIST("delta_eta/neg_sq"), deltaEtaWidth, fneg * fneg); + histogramRegistry.fill(HIST("delta_eta/posneg"), deltaEtaWidth, posneg); + } SliceCache cache; Preslice mcTrack = o2::aod::mcparticle::mcCollisionId; -// process function for Data Run3 + // process function for Data Run3 void processDataRun3(aod::MyCollisionRun3 const& coll, aod::MyTracks const& tracks) { calculationData(coll, tracks); for (int ii = 0; ii < deltaEta; ii++) { - float etaMin = -0.1f * (ii + 1); - float etaMax = 0.1f * (ii + 1); + float etaMin = -0.1f * (ii + 1); + float etaMax = 0.1f * (ii + 1); - calculationDeltaEta(coll, tracks, etaMin, etaMax); + calculationDeltaEta(coll, tracks, etaMin, etaMax); } } PROCESS_SWITCH(NetchargeFluctuations, processDataRun3, "Process for Run3 DATA", false); -// process function for Data Run2 + // process function for Data Run2 void processDataRun2(aod::MyCollisionRun2 const& coll, aod::MyTracks const& tracks) { calculationData(coll, tracks); for (int ii = 0; ii < deltaEta; ii++) { - float etaMin = -0.1f * (ii + 1); // -0.1, -0.2, ..., -0.8 - float etaMax = 0.1f * (ii + 1); // +0.1, +0.2, ..., +0.8 + float etaMin = -0.1f * (ii + 1); // -0.1, -0.2, ..., -0.8 + float etaMax = 0.1f * (ii + 1); // +0.1, +0.2, ..., +0.8 - calculationDeltaEta(coll, tracks, etaMin, etaMax); + calculationDeltaEta(coll, tracks, etaMin, etaMax); } } PROCESS_SWITCH(NetchargeFluctuations, processDataRun2, "Process for Run2 DATA", false); -// process function for MC Run3 + // process function for MC Run3 void processMcRun3(aod::MyMCCollisionRun3 const& coll, aod::MyMCTracks const& inputTracks, aod::McCollisions const& mcCollisions, aod::McParticles const& mcParticles) { calculationMc(coll, inputTracks, mcCollisions, mcParticles); - + for (int ii = 0; ii < deltaEta; ii++) { - float etaMin = -0.1f * (ii + 1); - float etaMax = 0.1f * (ii + 1); + float etaMin = -0.1f * (ii + 1); + float etaMax = 0.1f * (ii + 1); - calculationDeltaEta(coll, inputTracks, etaMin, etaMax); + calculationDeltaEta(coll, inputTracks, etaMin, etaMax); } - } PROCESS_SWITCH(NetchargeFluctuations, processMcRun3, "Process reconstructed", false); -// process function for MC Run2 + // process function for MC Run2 void processMcRun2(aod::MyMCCollisionRun2 const& coll, aod::MyMCTracks const& inputTracks, aod::McCollisions const& mcCollisions, aod::McParticles const& mcParticles) { calculationMc(coll, inputTracks, mcCollisions, mcParticles); - for (int ii = 0; ii < deltaEta; ii++) { - float etaMin = -0.1f * (ii + 1); // -0.1, -0.2, ..., -0.8 - float etaMax = 0.1f * (ii + 1); // +0.1, +0.2, ..., +0.8 + for (int ii = 0; ii < deltaEta; ii++) { + float etaMin = -0.1f * (ii + 1); // -0.1, -0.2, ..., -0.8 + float etaMax = 0.1f * (ii + 1); // +0.1, +0.2, ..., +0.8 - calculationDeltaEta(coll, inputTracks, etaMin, etaMax); + calculationDeltaEta(coll, inputTracks, etaMin, etaMax); } } PROCESS_SWITCH(NetchargeFluctuations, processMcRun2, "Process reconstructed", true); - -}; - -// struct -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) -{ - return WorkflowSpec{ - {adaptAnalysisTask(cfgc)} }; -} - + // struct + WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) + { + return WorkflowSpec{ + {adaptAnalysisTask(cfgc)}}; +}