[PWGEM/Dilepton] add a simple task to check PV in MC

PWGEM/Dilepton/Tasks/CMakeLists.txt

@@ -10,22 +10,22 @@
 # or submit itself to any jurisdiction.
 
 
 o2physics_add_dpl_workflow(efficiency-ee
  SOURCES emEfficiencyEE.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2::DetectorsBase O2Physics::PWGDQCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(lmee-lf-cocktail
  SOURCES lmeeLFCocktail.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(lmee-hf-cocktail
  SOURCES lmeeHFCocktail.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(mc-templates
  SOURCES MCtemplates.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2::DetectorsBase O2Physics::PWGDQCore
  COMPONENT_NAME Analysis)
@@ -55,32 +55,32 @@
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2::DetectorsBase O2Physics::PWGDQCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(create-ttp
  SOURCES createTTP.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(event-qc
  SOURCES eventQC.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore O2Physics::EventFilteringUtils
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(single-electron-qc
  SOURCES singleElectronQC.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2Physics::AnalysisCore O2Physics::PWGEMDileptonCore O2Physics::EventFilteringUtils
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(single-electron-qc-mc
  SOURCES singleElectronQCMC.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2Physics::AnalysisCore O2Physics::PWGEMDileptonCore
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(single-muon-qc
  SOURCES singleMuonQC.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2Physics::AnalysisCore O2Physics::PWGEMDileptonCore O2Physics::EventFilteringUtils
  COMPONENT_NAME Analysis)
 
 o2physics_add_dpl_workflow(single-muon-qc-mc
  SOURCES singleMuonQCMC.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsBase O2Physics::AnalysisCore O2Physics::PWGEMDileptonCore
  COMPONENT_NAME Analysis)
@@ -170,6 +170,11 @@
  PUBLIC_LINK_LIBRARIES O2::Framework O2::DetectorsVertexing O2::DCAFitter O2Physics::AnalysisCore
  COMPONENT_NAME Analysis)
 
+o2physics_add_dpl_workflow(test-pv
+ SOURCES testPV.cxx
+ PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore
+ COMPONENT_NAME Analysis)
+
 o2physics_add_dpl_workflow(evaluate-acceptance
  SOURCES evaluateAcceptance.cxx
  PUBLIC_LINK_LIBRARIES O2::Framework O2Physics::AnalysisCore

PWGEM/Dilepton/Tasks/testPV.cxx

@@ -0,0 +1,205 @@
+// 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 testPV.cxx
+/// \brief a task to check PV quality in MC
+/// \author daiki.sekihata@cern.ch
+
+#include "Common/CCDB/EventSelectionParams.h"
+#include "Common/DataModel/EventSelection.h"
+
+#include <Framework/AnalysisDataModel.h>
+#include <Framework/AnalysisHelpers.h>
+#include <Framework/AnalysisTask.h>
+#include <Framework/Configurable.h>
+#include <Framework/HistogramRegistry.h>
+#include <Framework/HistogramSpec.h>
+#include <Framework/InitContext.h>
+#include <Framework/runDataProcessing.h>
+
+#include <string>
+
+using namespace o2;
+using namespace o2::soa;
+using namespace o2::framework;
+using namespace o2::framework::expressions;
+using namespace o2::constants::physics;
+
+struct testPV {
+
+ // Configurables
+ // Configurable<std::string> ccdburl{"ccdb-url", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
+
+ struct : ConfigurableGroup {
+ std::string prefix = "eventCut";
+ Configurable<int> cfgEventGeneratorId{"cfgEventGeneratorId", -1, "event generator index. e.g. select gap/signal events"};
+ Configurable<float> cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"};
+ Configurable<float> cfgZvtxMax{"cfgZvtxMax", 10.f, "max. Zvtx"};
+ Configurable<bool> cfgRequireFT0AND{"cfgRequireFT0AND", true, "require FT0AND"};
+ Configurable<bool> cfgRequireNoTFB{"cfgRequireNoTFB", true, "require No time frame border"};
+ Configurable<bool> cfgRequireNoITSROFB{"cfgRequireNoITSROFB", false, "require no ITS readout frame border"};
+ Configurable<bool> cfgRequireNoSameBunchPileup{"cfgRequireNoSameBunchPileup", false, "require no same bunch pileup in event cut"};
+ Configurable<bool> cfgRequireGoodZvtxFT0vsPV{"cfgRequireGoodZvtxFT0vsPV", false, "require good Zvtx between FT0 vs. PV in event cut"};
+ Configurable<bool> cfgRequireVertexITSTPC{"cfgRequireVertexITSTPC", false, "require Vertex ITSTPC in event cut"}; // ITS-TPC matched track contributes PV.
+ Configurable<bool> cfgRequireVertexTOFmatched{"cfgRequireVertexTOFmatched", false, "require Vertex TOFmatched in event cut"}; // ITS-TPC-TOF matched track contributes PV.
+ // for RCT
+ o2::framework::Configurable<bool> cfgRequireGoodRCT{"cfgRequireGoodRCT", true, "require good detector flag in run condtion table"};
+ o2::framework::Configurable<std::string> cfgRCTLabel{"cfgRCTLabel", "CBT", "select 1 [CBT, CBT_hadronPID, CBT_muon_glo] see O2Physics/Common/CCDB/RCTSelectionFlags.h"};
+ o2::framework::Configurable<bool> cfgCheckZDC{"cfgCheckZDC", false, "set ZDC flag for AA"};
+ o2::framework::Configurable<bool> cfgTreatLimitedAcceptanceAsBad{"cfgTreatLimitedAcceptanceAsBad", false, "reject all events where the detectors relevant for the specified Runlist are flagged as LimitedAcceptance"};
+ } eventCut;
+
+ HistogramRegistry fRegistry{"fRegistry"};
+
+ void init(o2::framework::InitContext&)
+ {
+ // ccdb->setURL(ccdburl);
+ // ccdb->setCaching(true);
+ // ccdb->setLocalObjectValidityChecking();
+ // ccdb->setFatalWhenNull(false);
+
+ mRunNumber = 0;
+
+ addHistograms();
+ }
+
+ int mRunNumber{0};
+ // Service<o2::ccdb::BasicCCDBManager> ccdb;
+
+ template <typename TBC>
+ void initCCDB(TBC const& bc)
+ {
+ if (mRunNumber == bc.runNumber()) {
+ return;
+ }
+
+ mRunNumber = bc.runNumber();
+ }
+
+ void addHistograms()
+ {
+ auto hCollisionCounter = fRegistry.add<TH1>("hCollisionCounter", "collision counter", kTH1D, {{2, 0.5f, 2.5f}}, false);
+ hCollisionCounter->GetXaxis()->SetBinLabel(1, "all");
+ hCollisionCounter->GetXaxis()->SetBinLabel(2, "accepted");
+
+ fRegistry.add("Vertex/hZvtx", "vertex z; Z_{vtx} (cm)", kTH1F, {{100, -50, +50}}, false);
+ fRegistry.add("Vertex/hNContrib", "Number of PV contributors;N_{contrib}", kTH1F, {{101, -0.5, 100.5}}, false);
+ fRegistry.add("Vertex/hChi2", "vertex chi2;#chi^{2}/N_{contrib}", kTH1F, {{100, 0, 10}}, false);
+ fRegistry.add("Vertex/hChi2vsNContrib", "vertex #chi^{2}/N_{contrib} vs. N_{contrib};N_{contrib};#chi^{2}/N_{contrib}", kTH2F, {{101, -0.5, 100.5}, {100, 0, 10}}, false);
+
+ fRegistry.add("Vertex/hSigmaX", "vertex #sigma_{X} vs. N_{contrib};N_{contrib};#sigma_{X} (#mum)", kTH2F, {{101, -0.5, 100.5}, {1000, 0, 100}}, false);
+ fRegistry.add("Vertex/hSigmaY", "vertex #sigma_{Y} vs. N_{contrib};N_{contrib};#sigma_{Y} (#mum)", kTH2F, {{101, -0.5, 100.5}, {1000, 0, 100}}, false);
+ fRegistry.add("Vertex/hSigmaZ", "vertex #sigma_{Z} vs. N_{contrib};N_{contrib};#sigma_{Z} (#mum)", kTH2F, {{101, -0.5, 100.5}, {1000, 0, 100}}, false);
+
+ fRegistry.add("Vertex/hDeltaX", "vertex #DeltaX vs. N_{contrib};N_{contrib};#DeltaX = (X_{rec} #minus X_{gen})/#sigma_{X}", kTH2F, {{101, -0.5, 100.5}, {200, -10, 10}}, false);
+ fRegistry.add("Vertex/hDeltaY", "vertex #DeltaY vs. N_{contrib};N_{contrib};#DeltaY = (Y_{rec} #minus Y_{gen})/#sigma_{Y}", kTH2F, {{101, -0.5, 100.5}, {200, -10, 10}}, false);
+ fRegistry.add("Vertex/hDeltaZ", "vertex #DeltaZ vs. N_{contrib};N_{contrib};#DeltaZ = (Z_{rec} #minus Z_{gen})/#sigma_{Z}", kTH2F, {{101, -0.5, 100.5}, {200, -10, 10}}, false);
+
+ fRegistry.add("Vertex/hCollisionTime", "vertex time;N_{contrib};collision time (ns)", kTH2F, {{101, -0.5, 100.5}, {500, -25, 25}}, false);
+ fRegistry.add("Vertex/hCollisionTimeRes", "vertex time resolution;N_{contrib};collision time resolution (ns)", kTH2F, {{101, -0.5, 100.5}, {250, 0, 25}}, false);
+ }
+
+ template <typename TCollision>
+ bool isSelectedCollision(TCollision const& collision)
+ {
+ if (!(eventCut.cfgZvtxMin < collision.posZ() && collision.posZ() < eventCut.cfgZvtxMax)) {
+ return false;
+ }
+ if (eventCut.cfgRequireFT0AND && !collision.selection_bit(o2::aod::evsel::kIsTriggerTVX)) {
+ return false;
+ }
+ if (eventCut.cfgRequireNoTFB && !collision.selection_bit(o2::aod::evsel::kNoTimeFrameBorder)) {
+ return false;
+ }
+ if (eventCut.cfgRequireNoITSROFB && !collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder)) {
+ return false;
+ }
+ if (eventCut.cfgRequireNoSameBunchPileup && !collision.selection_bit(o2::aod::evsel::kNoSameBunchPileup)) {
+ return false;
+ }
+ if (eventCut.cfgRequireGoodZvtxFT0vsPV && !collision.selection_bit(o2::aod::evsel::kIsGoodZvtxFT0vsPV)) {
+ return false;
+ }
+ if (eventCut.cfgRequireVertexITSTPC && !collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC)) {
+ return false;
+ }
+ if (eventCut.cfgRequireVertexTOFmatched && !collision.selection_bit(o2::aod::evsel::kIsVertexTOFmatched)) {
+ return false;
+ }
+ return true;
+ }
+
+ template <typename TCollision, typename TMCCollision>
+ void fillVertexHistograms(TCollision const& collision, TMCCollision const& mcCollision)
+ {
+ fRegistry.fill(HIST("Vertex/hZvtx"), collision.posZ());
+ fRegistry.fill(HIST("Vertex/hNContrib"), collision.numContrib());
+ fRegistry.fill(HIST("Vertex/hChi2"), collision.chi2() / collision.numContrib());
+ fRegistry.fill(HIST("Vertex/hChi2vsNContrib"), collision.numContrib(), collision.chi2() / collision.numContrib());
+
+ fRegistry.fill(HIST("Vertex/hSigmaX"), collision.numContrib(), std::sqrt(collision.covXX()) * 1e+4); // convert cm to um
+ fRegistry.fill(HIST("Vertex/hSigmaY"), collision.numContrib(), std::sqrt(collision.covYY()) * 1e+4); // convert cm to um
+ fRegistry.fill(HIST("Vertex/hSigmaZ"), collision.numContrib(), std::sqrt(collision.covZZ()) * 1e+4); // convert cm to um
+
+ fRegistry.fill(HIST("Vertex/hDeltaX"), collision.numContrib(), (collision.posX() - mcCollision.posX()) / std::sqrt(collision.covXX()));
+ fRegistry.fill(HIST("Vertex/hDeltaY"), collision.numContrib(), (collision.posY() - mcCollision.posY()) / std::sqrt(collision.covYY()));
+ fRegistry.fill(HIST("Vertex/hDeltaZ"), collision.numContrib(), (collision.posZ() - mcCollision.posZ()) / std::sqrt(collision.covZZ()));
+
+ fRegistry.fill(HIST("Vertex/hCollisionTime"), collision.numContrib(), collision.collisionTime());
+ fRegistry.fill(HIST("Vertex/hCollisionTimeRes"), collision.numContrib(), collision.collisionTimeRes());
+ }
+
+ template <typename TBCs, typename TCollisions, typename TMCCollisions, typename TMCParticles>
+ void run(TBCs const&, TCollisions const& collisions, TMCCollisions const&, TMCParticles const&)
+ {
+ for (const auto& collision : collisions) {
+ if (!collision.has_mcCollision()) {
+ continue;
+ }
+ auto bc = collision.template bc_as<TBCs>();
+ initCCDB(bc);
+
+ auto mcCollision = collision.template mcCollision_as<aod::McCollisions>();
+ if (eventCut.cfgEventGeneratorId > -1 && mcCollision.getSubGeneratorId() != eventCut.cfgEventGeneratorId) {
+ continue;
+ }
+
+ fRegistry.fill(HIST("hCollisionCounter"), 1);
+
+ if (!isSelectedCollision(collision)) {
+ continue;
+ }
+ fRegistry.fill(HIST("hCollisionCounter"), 2);
+
+ fillVertexHistograms(collision, mcCollision);
+
+ } // end of collision loop
+ }
+
+ using MyCollisions = soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels>;
+ using MyBCs = soa::Join<aod::BCsWithTimestamps, aod::BcSels>;
+
+ Filter collisionFilter_evsel = eventCut.cfgZvtxMin < o2::aod::collision::posZ && o2::aod::collision::posZ < eventCut.cfgZvtxMax;
+ using FilteredMyCollisions = soa::Filtered<MyCollisions>;
+
+ void processMC(FilteredMyCollisions const& collisions, MyBCs const& bcs, aod::McCollisions const& mcCollisions, aod::McParticles const& mcParticles)
+ {
+ run(bcs, collisions, mcCollisions, mcParticles);
+ }
+ PROCESS_SWITCH(testPV, processMC, "processMC", true);
+
+ void processDummy(MyCollisions const&) {}
+ PROCESS_SWITCH(testPV, processDummy, "processDummy", false);
+};
+WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
+{
+ return WorkflowSpec{adaptAnalysisTask<testPV>(cfgc, TaskName{"test-pv"})};
+}