[PWGUD] fix efficiency and cuts (#16001)

Author: miedema-11

PWGUD/Tasks/flowCorrelationsUpc.cxx

@@ -153,28 +153,20 @@ struct FlowCorrelationsUpc {
   O2_DEFINE_CONFIGURABLE(cfgCutMerging, float, 0.02, "Merging cut on track merge")
   O2_DEFINE_CONFIGURABLE(cfgRadiusLow, float, 0.8, "Low radius for merging cut")
   O2_DEFINE_CONFIGURABLE(cfgRadiusHigh, float, 2.5, "High radius for merging cut")
-  O2_DEFINE_CONFIGURABLE(cfgIsGoodItsLayers, bool, false, "whether choose itslayers")
   O2_DEFINE_CONFIGURABLE(cfgDcaxy, bool, true, "choose dcaxy")
   O2_DEFINE_CONFIGURABLE(cfgDcaz, bool, false, "choose dcaz")
   O2_DEFINE_CONFIGURABLE(cfgDcazCut, float, 10.0, "dcaz cut")
-  O2_DEFINE_CONFIGURABLE(cfgItsClusterSize, unsigned int, 5, "ITS cluster size")
   O2_DEFINE_CONFIGURABLE(cfgMaxTPCChi2NCl, int, 4, "tpcchi2")
-  O2_DEFINE_CONFIGURABLE(cfgEvSelOccupancy, bool, true, "Occupancy cut")
+  O2_DEFINE_CONFIGURABLE(cfgCutOccupancy, bool, true, "Occupancy cut")
   O2_DEFINE_CONFIGURABLE(cfgCutOccupancyHigh, int, 1000, "High cut on TPC occupancy")
   O2_DEFINE_CONFIGURABLE(cfgCutOccupancyLow, int, 0, "Low cut on TPC occupancy")
   O2_DEFINE_CONFIGURABLE(cfgCutTPCCrossedRows, float, 70.0f, "minimum number of crossed TPC Rows")
   O2_DEFINE_CONFIGURABLE(cfgCutTPCclu, float, 50.0f, "minimum number of found TPC clusters")
   O2_DEFINE_CONFIGURABLE(cfgCutITSclu, float, 5.0f, "minimum number of ITS clusters")
-  O2_DEFINE_CONFIGURABLE(cfgCutTPCChi2NCl, int, 4, "max chi2 per TPC clusters")
   O2_DEFINE_CONFIGURABLE(cfgGlobalTrack, bool, true, "require TPC+ITS track")
   O2_DEFINE_CONFIGURABLE(cfgUseNchCorrected, bool, true, "use corrected Nch for X axis")
-  O2_DEFINE_CONFIGURABLE(cfgOutputNUAWeights, bool, false, "Fill and output NUA weights")
-  O2_DEFINE_CONFIGURABLE(cfgOutputNUAWeightsRefPt, bool, false, "NUA weights are filled in ref pt bins")
-  O2_DEFINE_CONFIGURABLE(cfgOutputNUAWeightsRunbyRun, bool, false, "NUA weights are filled run-by-run")
   O2_DEFINE_CONFIGURABLE(cfgEfficiency, std::string, "", "CCDB path to efficiency object")
-  O2_DEFINE_CONFIGURABLE(cfgAcceptance, std::string, "", "CCDB path to acceptance object")
   O2_DEFINE_CONFIGURABLE(cfgUseEventWeights, bool, false, "Use event weights for mixed event")
-  O2_DEFINE_CONFIGURABLE(cfgLocalEfficiency, bool, false, "Use local efficiency object")
 
   ConfigurableAxis axisVertex{"axisVertex", {10, -10, 10}, "vertex axis for histograms"};
   ConfigurableAxis axisEta{"axisEta", {40, -1., 1.}, "eta axis for histograms"};
@@ -203,25 +195,22 @@ struct FlowCorrelationsUpc {
   ConfigurableAxis axisNch{"axisNch", {300, 0, 300}, "N_{ch}"};
 
   // Corrections
-  TH1D* mEfficiency = nullptr;
-  GFWWeights* mAcceptance = nullptr;
+  TH3D* mEfficiency = nullptr;
   bool correctionsLoaded = false;
 
   // make the filters and cuts.
   Filter trackFilter = (aod::udtrack::isPVContributor == true);
-  Filter collisionFilter = ((aod::udcollision::gapSide == (uint8_t)1 || aod::udcollision::gapSide == (uint8_t)0) && (cfgIfVertex == false || aod::collision::posZ < cfgZVtxCut) && (aod::udcollision::occupancyInTime > 0 && aod::udcollision::occupancyInTime < cfgCutOccupancyHigh) && (aod::flowcorrupc::truegapside == 1 || aod::flowcorrupc::truegapside == 0));
+  Filter collisionFilter = ((aod::udcollision::gapSide == (uint8_t)1 || aod::udcollision::gapSide == (uint8_t)0) && (cfgIfVertex == false || aod::collision::posZ < cfgZVtxCut) && (!cfgCutOccupancy || (aod::udcollision::occupancyInTime > 0 && aod::udcollision::occupancyInTime < cfgCutOccupancyHigh)) && (aod::flowcorrupc::truegapside == 1 || aod::flowcorrupc::truegapside == 0));
 
   // Connect to ccdb
   Service<ccdb::BasicCCDBManager> ccdb;
   Configurable<std::string> ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
 
   OutputObj<GFWWeights> fWeights{GFWWeights("weights")};
-  OutputObj<GFWWeights> fWeightsMc{GFWWeights("weightsMC")};
 
   TAxis* fPtAxis = nullptr;
   int lastRunNumber = -1;
-  std::vector<int> runNumbers;
-  std::map<int, std::shared_ptr<TH3>> th3sPerRun; // map of TH3 histograms for all runs
+  std::vector<int> runNumbers; // map of TH3 histograms for all runs
   std::vector<float> efficiencyCache;
 
   using UdTracks = soa::Filtered<soa::Join<aod::UDTracks, aod::UDTracksExtra, aod::UDTracksPID>>;
@@ -248,9 +237,6 @@ struct FlowCorrelationsUpc {
     registry.add("pT", "pT", {HistType::kTH1D, {axisPtTrigger}});
     registry.add("Nch", "N_{ch}", {HistType::kTH1D, {axisMultiplicity}});
     registry.add("zVtx", "zVtx", {HistType::kTH1D, {axisVertex}});
-    registry.add("Nch_vs_zVtx", "Nch vs zVtx", {HistType::kTH2D, {axisVertex, axisMultiplicity}});
-    registry.add("Nch_same", "Nch same event", {HistType::kTH1D, {axisMultiplicity}});
-    registry.add("Nch_mixed", "Nch mixed event", {HistType::kTH1D, {axisMultiplicity}});
     registry.add("EtaCorrected", "Eta corrected", {HistType::kTH1D, {axisEta}});
     registry.add("pTCorrected", "pT corrected", {HistType::kTH1D, {axisPtTrigger}});
 
@@ -273,13 +259,6 @@ struct FlowCorrelationsUpc {
     double* ptBins = &(axis.binEdges)[0];
     fPtAxis = new TAxis(nPtBins, ptBins);
 
-    if (cfgOutputNUAWeights) {
-      fWeights->setPtBins(nPtBins, ptBins);
-      fWeights->init(true, false);
-      fWeightsMc->setPtBins(nPtBins, ptBins);
-      fWeightsMc->init(true, false);
-    }
-
     std::vector<AxisSpec> corrAxis = {{axisSample, "Sample"},
                                       {axisVertex, "z-vtx (cm)"},
                                       {axisIndependent, "Independent (N_{ch} corrected)"},
@@ -333,6 +312,9 @@ struct FlowCorrelationsUpc {
     if (track.pt() < cfgPtCutMin || track.pt() > cfgPtCutMax) {
       return false;
     }
+    if (cfgGlobalTrack && !(track.hasITS() && track.hasTPC())) {
+      return false;
+    }
     // registry.fill(HIST("hTrackCount"), 1.5);
     if (cfgDcaz && !(std::fabs(track.dcaZ()) < cfgDcazCut)) {
       return false;
@@ -361,28 +343,13 @@ struct FlowCorrelationsUpc {
     return true;
   }
 
-  void createOutputObjectsForRun(int runNumber)
-  {
-    const AxisSpec axisPhi{60, 0.0, constants::math::TwoPI, "#varphi"};
-    std::shared_ptr<TH3> histPhiEtaVtxz = registry.add<TH3>(Form("%d/hPhiEtaVtxz", runNumber), ";#varphi;#eta;v_{z}", {HistType::kTH3D, {axisPhi, {64, -1.8, 1.8}, {40, -10, 10}}});
-    th3sPerRun.insert(std::make_pair(runNumber, histPhiEtaVtxz));
-  }
-
   void loadCorrections(uint64_t timestamp)
   {
     if (correctionsLoaded) {
       return;
     }
-    if (cfgAcceptance.value.empty() == false) {
-      mAcceptance = ccdb->getForTimeStamp<GFWWeights>(cfgAcceptance, timestamp);
-      if (mAcceptance) {
-        LOGF(info, "Loaded acceptance weights from %s (%p)", cfgAcceptance.value.c_str(), (void*)mAcceptance);
-      } else {
-        LOGF(warning, "Could not load acceptance weights from %s (%p)", cfgAcceptance.value.c_str(), (void*)mAcceptance);
-      }
-    }
     if (cfgEfficiency.value.empty() == false) {
-      mEfficiency = ccdb->getForTimeStamp<TH1D>(cfgEfficiency, timestamp);
+      mEfficiency = ccdb->getForTimeStamp<TH3D>(cfgEfficiency, timestamp);
       if (mEfficiency == nullptr) {
         LOGF(fatal, "Could not load efficiency histogram for trigger particles from %s", cfgEfficiency.value.c_str());
       }
@@ -408,7 +375,7 @@ struct FlowCorrelationsUpc {
     return true;
   }
 
-  bool setCurrentParticleWeights(float& weight_nue, float& weight_nua, float phi, float eta, float pt, float vtxz)
+  bool setCurrentParticleWeights(float& weight_nue, float& weight_nua, float pt)
   {
     float eff = 1.;
     if (mEfficiency)
@@ -418,17 +385,13 @@ struct FlowCorrelationsUpc {
     if (eff == 0)
       return false;
     weight_nue = 1. / eff;
-
-    if (mAcceptance)
-      weight_nua = mAcceptance->getNUA(phi, eta, vtxz);
-    else
-      weight_nua = 1;
+    weight_nua = 1.; // Set to 1 as NUA weight is not being used
     return true;
   }
 
   // fill multiple histograms
   template <typename TCollision, typename TTracks>
-  void fillYield(TCollision collision, TTracks tracks, int runNumber, float vtxz) // function to fill the yield and etaphi histograms.
+  void fillYield(TCollision collision, TTracks tracks, float vtxz) // function to fill the yield and etaphi histograms.
   {
     registry.fill(HIST("Nch"), tracks.size());
     registry.fill(HIST("zVtx"), collision.posZ());
@@ -449,20 +412,6 @@ struct FlowCorrelationsUpc {
       registry.fill(HIST("pT"), pt);
       registry.fill(HIST("EtaCorrected"), eta, weff);
       registry.fill(HIST("pTCorrected"), pt, weff);
-
-      if (cfgOutputNUAWeights) {
-        if (cfgOutputNUAWeightsRefPt) {
-          if (pt > cfgPtCutMin && pt < cfgPtCutMax) {
-            fWeights->fill(phi, eta, vtxz, pt, 0., 0);
-            if (cfgOutputNUAWeightsRunbyRun)
-              th3sPerRun[runNumber]->Fill(phi, eta, vtxz);
-          }
-        } else {
-          fWeights->fill(phi, eta, vtxz, pt, 0., 0);
-          if (cfgOutputNUAWeightsRunbyRun)
-            th3sPerRun[runNumber]->Fill(phi, eta, vtxz);
-        }
-      }
     }
   }
 
@@ -497,7 +446,7 @@ struct FlowCorrelationsUpc {
 
       // 计算track1的权重
       float weff1 = 1., wacc1 = 1.;
-      if (!setCurrentParticleWeights(weff1, wacc1, phi1, eta1, pt1, vtxz))
+      if (!setCurrentParticleWeights(weff1, wacc1, pt1))
         continue;
 
       if (system == SameEvent) {
@@ -520,24 +469,16 @@ struct FlowCorrelationsUpc {
         double phi2 = RecoDecay::phi(momentum);
         double eta2 = RecoDecay::eta(momentum);
 
-        // 计算track2的权重
         float weff2 = 1., wacc2 = 1.;
         if (mEfficiency) {
           weff2 = efficiencyCache[track2.filteredIndex()];
         } else {
           getEfficiencyCorrection(weff2, eta2, pt2, vtxz);
         }
 
-        if (mAcceptance) {
-          wacc2 = mAcceptance->getNUA(phi2, eta2, vtxz);
-        } else {
-          wacc2 = 1;
-        }
-
         float deltaPhi = RecoDecay::constrainAngle(phi1 - phi2, -PIHalf);
         float deltaEta = eta1 - eta2;
 
-        // 计算组合权重
         float weight = eventWeight * weff1 * weff2 * wacc1 * wacc2;
 
         // Merging cut
@@ -585,19 +526,6 @@ struct FlowCorrelationsUpc {
 
     loadCorrections(runDuration.first);
 
-    if (cfgOutputNUAWeightsRunbyRun && currentRunNumber != lastRunNumber) {
-      lastRunNumber = currentRunNumber;
-      if (std::find(runNumbers.begin(), runNumbers.end(), currentRunNumber) == runNumbers.end()) {
-        createOutputObjectsForRun(currentRunNumber);
-        runNumbers.push_back(currentRunNumber);
-      }
-
-      if (th3sPerRun.find(currentRunNumber) == th3sPerRun.end()) {
-        LOGF(fatal, "RunNumber %d not found in th3sPerRun", currentRunNumber);
-        return;
-      }
-    }
-
     registry.fill(HIST("eventcont"), 3.5);
 
     //-----------independent---------------
@@ -628,7 +556,7 @@ struct FlowCorrelationsUpc {
       independent = nTracksCorrected;
     }
 
-    fillYield(collision, tracks, currentRunNumber, vtxz);
+    fillYield(collision, tracks, vtxz);
 
     fillCorrelations<CorrelationContainer::kCFStepReconstructed>(
       tracks, tracks, collision.posZ(), SameEvent,