From e127a71bdda9a0b70b3e91162862f20ef63689b2 Mon Sep 17 00:00:00 2001 From: Chiara Zampolli Date: Tue, 6 Aug 2024 15:53:35 +0200 Subject: [PATCH 1/4] Adding plots vs TPC occupancy --- Detectors/GLOQC/CMakeLists.txt | 4 +- Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h | 23 +++- Detectors/GLOQC/src/MatchITSTPCQC.cxx | 126 ++++++++++++++++-- GPU/GPUTracking/CMakeLists.txt | 1 + 4 files changed, 138 insertions(+), 16 deletions(-) diff --git a/Detectors/GLOQC/CMakeLists.txt b/Detectors/GLOQC/CMakeLists.txt index cdb307a6515c8..9d2db85460b69 100644 --- a/Detectors/GLOQC/CMakeLists.txt +++ b/Detectors/GLOQC/CMakeLists.txt @@ -15,7 +15,9 @@ o2_add_library(GLOQC SOURCES src/MatchITSTPCQC.cxx src/ITSTPCMatchingQCParams.cxx - PUBLIC_LINK_LIBRARIES O2::DetectorsVertexing) + PUBLIC_LINK_LIBRARIES O2::DetectorsVertexing + PRIVATE_LINK_LIBRARIES O2::GPUO2Interface + O2::GPUTracking) o2_target_root_dictionary(GLOQC HEADERS include/GLOQC/MatchITSTPCQC.h diff --git a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h index 035aada3a765f..83da8e42f3110 100644 --- a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h +++ b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h @@ -19,6 +19,7 @@ #include #include #include +#include #include #include #include "DataFormatsGlobalTracking/RecoContainer.h" @@ -28,6 +29,11 @@ #include "Steer/MCKinematicsReader.h" #include "ReconstructionDataFormats/PID.h" #include "DCAFitter/DCAFitterN.h" +#include "GPUO2InterfaceConfiguration.h" +// #include "GPUSettingsO2.h" +#include "GPUParam.h" +#include "GPUParam.inc" + #include #include #include @@ -63,7 +69,7 @@ class MatchITSTPCQC void setDataRequest(const std::shared_ptr& dr) { mDataRequest = dr; } void finalize(); void reset(); - bool processV0(int iv, o2::globaltracking::RecoContainer& recoData); + bool processV0(int iv, o2::globaltracking::RecoContainer& recoData, std::vector& mTBinClOcc, float pvTime); bool refitV0(const o2::dataformats::V0Index& id, o2::dataformats::V0& v0, o2::globaltracking::RecoContainer& recoData); TH1D* getHistoPtNum(matchType m) const { return mPtNum[m]; } @@ -130,7 +136,7 @@ class MatchITSTPCQC TH1D* getHisto1OverPtPhysPrimDen(matchType m) const { return m1OverPtPhysPrimDen[m]; } TEfficiency* getFractionITSTPCmatchPhysPrim1OverPt(matchType m) const { return mFractionITSTPCmatchPhysPrim1OverPt[m]; } - TH2F* getHistoK0MassVsPt() const { return mK0MassVsPt; } + TH3F* getHistoK0MassVsPtVsOcc() const { return mK0MassVsPtVsOcc; } void getHistos(TObjArray& objar); @@ -228,7 +234,7 @@ class MatchITSTPCQC publisher->startPublishing(mDCArVsPtDen); publisher->startPublishing(mFractionITSTPCmatchDCArVsPt); if (mDoK0QC) { - publisher->startPublishing(mK0MassVsPt); + publisher->startPublishing(mK0MassVsPtVsOcc); } } @@ -408,12 +414,21 @@ class MatchITSTPCQC // for V0s o2::vertexing::DCAFitterN<2> mFitterV0; - TH2F* mK0MassVsPt = nullptr; + TH3F* mK0MassVsPtVsOcc = nullptr; bool mDoK0QC = false; // whether to fill the K0 QC plot(s) float mCutK0Mass = 0.05; // cut on the difference between the K0 mass and the PDG mass bool mRefit = false; // whether to refit or not float mMaxEtaK0 = 0.8; // cut on the K0 eta long int mTimestamp = -1; // timestamp used to load the SVertexParam object: if differnt from -1, we don't load (it means we already did it) + std::unique_ptr mConfig; + std::unique_ptr mConfParam; + // std::unique_ptr mParam; + std::shared_ptr mParam = nullptr; + int mNHBPerTF = 0; + int mNTPCOccBinLength = 0; ///< TPC occ. histo bin length in TBs + float mNTPCOccBinLengthInv; + std::vector mTBinClOcc; ///< TPC occupancy histo: i-th entry is the integrated occupancy for ~1 orbit starting from the TB = i*mNTPCOccBinLength + gsl::span mTPCRefitterOccMap; ///< externally set TPC clusters occupancy map ClassDefNV(MatchITSTPCQC, 3); }; diff --git a/Detectors/GLOQC/src/MatchITSTPCQC.cxx b/Detectors/GLOQC/src/MatchITSTPCQC.cxx index a0ae896105d2d..e95410275a138 100644 --- a/Detectors/GLOQC/src/MatchITSTPCQC.cxx +++ b/Detectors/GLOQC/src/MatchITSTPCQC.cxx @@ -26,6 +26,11 @@ #include "DetectorsVertexing/SVertexerParams.h" #include "Framework/InputRecord.h" #include "Framework/TimingInfo.h" +#include "GPUO2InterfaceUtils.h" +#include "CommonConstants/LHCConstants.h" +#include "DataFormatsTPC/Constants.h" + +#include "GPUO2InterfaceRefit.h" using namespace o2::gloqc; using namespace o2::mcutils; @@ -42,6 +47,7 @@ MatchITSTPCQC::~MatchITSTPCQC() void MatchITSTPCQC::deleteHistograms() { + LOG(info) << "Deleting histos..."; for (int i = 0; i < matchType::SIZE; ++i) { // Pt delete mPtNum[i]; @@ -124,7 +130,7 @@ void MatchITSTPCQC::deleteHistograms() delete mFractionITSTPCmatchDCArVsPt; // K0 - delete mK0MassVsPt; + delete mK0MassVsPtVsOcc; } //__________________________________________________________ @@ -205,7 +211,7 @@ void MatchITSTPCQC::reset() // K0 if (mDoK0QC) { - mK0MassVsPt->Reset(); + mK0MassVsPtVsOcc->Reset(); } } @@ -375,9 +381,39 @@ bool MatchITSTPCQC::init() } } + // log binning for pT for K0s + const Int_t nbinsPtK0 = 10; + const Double_t xminPtK0 = 0.01; + const Double_t xmaxPtK0 = 20; + Double_t* xbinsPtK0 = new Double_t[nbinsPtK0 + 1]; + Double_t xlogminPtK0 = TMath::Log10(xminPtK0); + Double_t xlogmaxPtK0 = TMath::Log10(xmaxPtK0); + Double_t dlogxPtK0 = (xlogmaxPtK0 - xlogminPtK0) / nbinsPtK0; + for (int i = 0; i <= nbinsPtK0; i++) { + Double_t xlogPtK0 = xlogminPtK0 + i * dlogxPtK0; + xbinsPtK0[i] = TMath::Exp(TMath::Log(10) * xlogPtK0); + } + // the other bins + const Int_t nbinsMassK0 = 100; + Double_t* ybinsMassK0 = new Double_t[nbinsMassK0 + 1]; + Double_t yminMassK0 = 0.4; + Double_t ymaxMassK0 = 0.6; + Double_t dyMassK0 = (ymaxMassK0 - yminMassK0) / nbinsMassK0; + for (int i = 0; i <= nbinsMassK0; i++) { + ybinsMassK0[i] = yminMassK0 + i * dyMassK0; + } + const Int_t nbinsMultK0 = 5; + Double_t* zbinsMultK0 = new Double_t[nbinsMultK0 + 1]; + Double_t zminMultK0 = 100000.; + Double_t zmaxMultK0 = 1000000.; + Double_t dzMultK0 = (zmaxMultK0 - zminMultK0) / nbinsMultK0; + for (int i = 0; i <= nbinsMultK0; i++) { + zbinsMultK0[i] = zminMultK0 + i * dzMultK0; + } + if (mDoK0QC) { // V0s - mK0MassVsPt = new TH2F("mK0MassVsPt", "K0 invariant mass vs Pt; Pt [GeV/c]; K0s mass [GeV/c^2]", 100, 0.f, 20.f, 100, 0.3, 0.7); + mK0MassVsPtVsOcc = new TH3F("mK0MassVsPt", "K0 invariant mass vs Pt; Pt [GeV/c]; K0s mass [GeV/c^2]", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0, zbinsMultK0); } LOG(info) << "Printing configuration cuts"; @@ -404,6 +440,7 @@ void MatchITSTPCQC::initDataRequest() if (mDoK0QC) { mDataRequest->requestPrimaryVertices(mUseMC); mDataRequest->requestSecondaryVertices(mUseMC); + mDataRequest->requestTPCClusters(false); } } @@ -900,13 +937,74 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) const auto pvertices = mRecoCont.getPrimaryVertices(); LOG(info) << "Found " << pvertices.size() << " primary vertices"; + // getting occupancy estimator + mNHBPerTF = o2::base::GRPGeomHelper::instance().getGRPECS()->getNHBFPerTF(); + if (!mParam) { + // for occupancy estimator + mParam = o2::gpu::GPUO2InterfaceUtils::getFullParamShared(0.f, mNHBPerTF); + } + size_t occupancyMapSizeBytes = o2::gpu::GPUO2InterfaceRefit::fillOccupancyMapGetSize(mNHBPerTF, mParam.get()); + LOG(debug) << "occupancyMapSizeBytes = " << occupancyMapSizeBytes; + mTPCRefitterOccMap = mRecoCont.occupancyMapTPC; + o2::gpu::GPUO2InterfaceUtils::paramUseExternalOccupancyMap(mParam.get(), mNHBPerTF, mTPCRefitterOccMap.data(), occupancyMapSizeBytes); + + std::vector mTBinClOcc; ///< TPC occupancy histo: i-th entry is the integrated occupancy for ~1 orbit starting from the TB = i * mNTPCOccBinLength + mTBinClOcc.clear(); + int mNTPCOccBinLength = mParam->rec.tpc.occupancyMapTimeBins; + mNTPCOccBinLengthInv = 1. / mNTPCOccBinLength; + if (mNTPCOccBinLength > 1 && mTPCRefitterOccMap.size()) { + int nTPCBinsInTF = mNHBPerTF * o2::constants::lhc::LHCMaxBunches / 8; // number of TPC time bins in 1 TF, considering that 1 TPC time bin is 8 bunches + int ninteg = 0; + int nTPCOccBinsInTF = nTPCBinsInTF * mNTPCOccBinLengthInv; // how many occupancy bins in 1 TF; mNTPCOccBinLengthInv is the inverse of the length of an occupancy bin + int sumBins = std::max(1, int(o2::constants::lhc::LHCMaxBunches / 8 * mNTPCOccBinLengthInv)); // we will integrate occupancy at max for this number of bins: the max between 1 and the number of occupancy bins in 1 orbit + LOG(debug) << "number of TPC TB in 1 TF = nTPCBinsInTF = " << nTPCBinsInTF << " ; number of occupancy bins in 1 TF = nTPCOccBinsInTF = " << nTPCOccBinsInTF; + LOG(debug) << "bins to integrate = sumBins = " << sumBins; + mTBinClOcc.resize(nTPCOccBinsInTF); + std::vector mltHistTB(nTPCOccBinsInTF); + float sm = 0., tb = 0.5 * mNTPCOccBinLength; + bool foundNotZero = false; + for (int i = 0; i < nTPCOccBinsInTF; i++) { // for every occupancy bin in the TF + mltHistTB[i] = mParam->GetUnscaledMult(tb); + if (mParam->GetUnscaledMult(tb) != 0) { + LOG(debug) << "i = " << i << " tb = " << tb << " mltHistTB[" << i << "] = " << mltHistTB[i]; + foundNotZero = true; + } + tb += mNTPCOccBinLength; + } + if (!foundNotZero) { + LOG(debug) << "No mult bin was found different from 0!"; + } + foundNotZero = false; + // now we fill the occupancy map; we integrate the sumBins after the current one, but when we are at the last 27 bins of the TF, where we integrate what we have left till the end of the TF; for practical reasons, we start from the end, adding all the time, and then also removing the last bin, when we have enough, so that we always add together sumBins bins (except, as said, for the last part of the TF) + for (int i = nTPCOccBinsInTF; i--;) { + if (mltHistTB[i] != 0) { + foundNotZero = true; + } + LOG(debug) << "i = " << i << " sm before = " << sm; + sm += mltHistTB[i]; + LOG(debug) << "i = " << i << " sm after = " << sm; + if (i + sumBins < nTPCOccBinsInTF) { + LOG(debug) << "i = " << i << " sumBins = " << sumBins << " nTPCOccBinsInTF = " << nTPCOccBinsInTF << " we have to decrease sm by = " << mltHistTB[i + sumBins]; + sm -= mltHistTB[i + sumBins]; + LOG(debug) << "i = " << i << " sm after 2 = " << sm; + } + mTBinClOcc[i] = sm; + LOG(debug) << "i = " << i << " mTBinClOcc[" << i << "] = " << mTBinClOcc[i]; + } + if (!foundNotZero) { + LOG(debug) << "No mult bin was found different from 0! sm = " << sm; + } + } else { + mTBinClOcc.resize(1); + } + auto v0IDs = mRecoCont.getV0sIdx(); auto nv0 = v0IDs.size(); if (nv0 > mRecoCont.getV0s().size()) { mRefit = true; } - LOG(info) << "Found " << mRecoCont.getV0s().size() << " V0s in reco container"; - LOG(info) << "Found " << nv0 << " V0s ids"; + LOG(debug) << "Found " << mRecoCont.getV0s().size() << " V0s in reco container"; + LOG(debug) << "Found " << nv0 << " V0s ids"; // associating sec vtxs to prim vtx std::map> pv2sv; static int tfID = 0; @@ -922,18 +1020,20 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) if (pvID < 0 || vv.size() == 0) { continue; } + const auto& pv = mRecoCont.getPrimaryVertex(pvID); + float pvTime = pv.getTimeStamp().getTimeStamp(); // in \mus for (int iv0 : vv) { - nV0sOk += processV0(iv0, mRecoCont) ? 1 : 0; + nV0sOk += processV0(iv0, mRecoCont, mTBinClOcc, pvTime) ? 1 : 0; } } - LOG(info) << "Processed " << nV0sOk << " V0s"; + LOG(debug) << "Processed " << nV0sOk << " V0s"; } evCount++; } //__________________________________________________________ -bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoData) +bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoData, std::vector& mTBinClOcc, float pvTime) { o2::dataformats::V0 v0; auto v0s = recoData.getV0s(); @@ -948,11 +1048,15 @@ bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoDat if (mMaxEtaK0 < std::abs(v0sel.getEta())) { return false; } - LOG(info) << "Find K0 with mass " << std::sqrt(v0sel.calcMass2AsK0()); if (mCutK0Mass > 0 && std::abs(std::sqrt(v0sel.calcMass2AsK0()) - 0.497) > mCutK0Mass) { return false; } - mK0MassVsPt->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0())); + // get the corresponding PV + int tb = pvTime / (8 * o2::constants::lhc::LHCBunchSpacingMUS) * mNTPCOccBinLengthInv; // V0 time in TPC time bins + LOG(debug) << "pvTime = " << pvTime << " tb = " << tb; + float mltTPC = tb < 0 ? mTBinClOcc[0] : (tb >= mTBinClOcc.size() ? mTBinClOcc.back() : mTBinClOcc[tb]); + LOG(debug) << "Filling plot with pt = " << v0sel.getPt() << " mass = " << std::sqrt(v0sel.calcMass2AsK0()) << " mult TPC = " << mltTPC; + mK0MassVsPtVsOcc->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); return true; } @@ -1228,5 +1332,5 @@ void MatchITSTPCQC::getHistos(TObjArray& objar) objar.Add(mFractionITSTPCmatchDCArVsPt); // V0 - objar.Add(mK0MassVsPt); + objar.Add(mK0MassVsPtVsOcc); } diff --git a/GPU/GPUTracking/CMakeLists.txt b/GPU/GPUTracking/CMakeLists.txt index a4ca1f126f013..2fd512e8f7cfb 100644 --- a/GPU/GPUTracking/CMakeLists.txt +++ b/GPU/GPUTracking/CMakeLists.txt @@ -166,6 +166,7 @@ set(HDRS_INSTALL TRDTracking/GPUTRDTrackletLabels.h TRDTracking/GPUTRDTrackPoint.h TRDTracking/GPUTRDTrackPoint.h + DataTypes/GPUTPCGMPolynomialField.h ) # Sources for O2 and for Standalone if requested in config file From b28079d619a009759937728f2d9d4c82a2cc7bc6 Mon Sep 17 00:00:00 2001 From: shahoian Date: Thu, 3 Oct 2024 15:42:03 +0200 Subject: [PATCH 2/4] Unconditionally load PVs if requested --- .../Detectors/GlobalTracking/src/RecoContainer.cxx | 12 ++++-------- 1 file changed, 4 insertions(+), 8 deletions(-) diff --git a/DataFormats/Detectors/GlobalTracking/src/RecoContainer.cxx b/DataFormats/Detectors/GlobalTracking/src/RecoContainer.cxx index e634395f6d723..60c18b966abed 100644 --- a/DataFormats/Detectors/GlobalTracking/src/RecoContainer.cxx +++ b/DataFormats/Detectors/GlobalTracking/src/RecoContainer.cxx @@ -833,13 +833,11 @@ void RecoContainer::addSVertices(ProcessingContext& pc, bool) //____________________________________________________________ void RecoContainer::addPVertices(ProcessingContext& pc, bool mc) { - if (!pvtxPool.isLoaded(PVTX)) { // in case was loaded via addPVerticesTMP - pvtxPool.registerContainer(pc.inputs().get>("pvtx"), PVTX); - } + pvtxPool.registerContainer(pc.inputs().get>("pvtx"), PVTX); pvtxPool.registerContainer(pc.inputs().get>("pvtx_trmtc"), PVTX_TRMTC); pvtxPool.registerContainer(pc.inputs().get>("pvtx_tref"), PVTX_TRMTCREFS); - if (mc && !pvtxPool.isLoaded(PVTX_MCTR)) { // in case was loaded via addPVerticesTMP + if (mc) { // in case was loaded via addPVerticesTMP pvtxPool.registerContainer(pc.inputs().get>("pvtx_mc"), PVTX_MCTR); } } @@ -856,13 +854,11 @@ void RecoContainer::addStrangeTracks(ProcessingContext& pc, bool mc) //____________________________________________________________ void RecoContainer::addPVerticesTMP(ProcessingContext& pc, bool mc) { - if (!pvtxPool.isLoaded(PVTX)) { // in case was loaded via addPVertices - pvtxPool.registerContainer(pc.inputs().get>("pvtx"), PVTX); - } + pvtxPool.registerContainer(pc.inputs().get>("pvtx"), PVTX); pvtxPool.registerContainer(pc.inputs().get>("pvtx_cont"), PVTX_CONTID); pvtxPool.registerContainer(pc.inputs().get>("pvtx_contref"), PVTX_CONTIDREFS); - if (mc && !pvtxPool.isLoaded(PVTX_MCTR)) { // in case was loaded via addPVertices + if (mc) { // in case was loaded via addPVertices pvtxPool.registerContainer(pc.inputs().get>("pvtx_mc"), PVTX_MCTR); } } From a46c8c402fc4f2877d16ebe32d01ff91531241f5 Mon Sep 17 00:00:00 2001 From: Chiara Zampolli Date: Thu, 3 Oct 2024 21:48:20 +0200 Subject: [PATCH 3/4] Different histos in case of pp or PbPb (waiting for a feature in QC) and sampling --- Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h | 14 +++- Detectors/GLOQC/src/MatchITSTPCQC.cxx | 80 ++++++++++++++----- 2 files changed, 70 insertions(+), 24 deletions(-) diff --git a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h index 83da8e42f3110..1aaf34bc244fb 100644 --- a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h +++ b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h @@ -136,7 +136,8 @@ class MatchITSTPCQC TH1D* getHisto1OverPtPhysPrimDen(matchType m) const { return m1OverPtPhysPrimDen[m]; } TEfficiency* getFractionITSTPCmatchPhysPrim1OverPt(matchType m) const { return mFractionITSTPCmatchPhysPrim1OverPt[m]; } - TH3F* getHistoK0MassVsPtVsOcc() const { return mK0MassVsPtVsOcc; } + TH3F* getHistoK0MassVsPtVsOccpp() const { return mK0MassVsPtVsOccpp; } + TH3F* getHistoK0MassVsPtVsOccPbPb() const { return mK0MassVsPtVsOccPbPb; } void getHistos(TObjArray& objar); @@ -234,7 +235,8 @@ class MatchITSTPCQC publisher->startPublishing(mDCArVsPtDen); publisher->startPublishing(mFractionITSTPCmatchDCArVsPt); if (mDoK0QC) { - publisher->startPublishing(mK0MassVsPtVsOcc); + publisher->startPublishing(mK0MassVsPtVsOccpp); + publisher->startPublishing(mK0MassVsPtVsOccPbPb); } } @@ -247,6 +249,8 @@ class MatchITSTPCQC void setBz(float bz) { mBz = bz; } void setDoK0QC(bool v) { mDoK0QC = v; } bool getDoK0QC() const { return mDoK0QC; } + void setK0Scaling(float v) { mK0Scaling = v; } + float getK0Scaling() const { return mK0Scaling; } // ITS track void setMinPtITSCut(float v) { mPtITSCut = v; }; @@ -414,7 +418,8 @@ class MatchITSTPCQC // for V0s o2::vertexing::DCAFitterN<2> mFitterV0; - TH3F* mK0MassVsPtVsOcc = nullptr; + TH3F* mK0MassVsPtVsOccpp = nullptr; + TH3F* mK0MassVsPtVsOccPbPb = nullptr; bool mDoK0QC = false; // whether to fill the K0 QC plot(s) float mCutK0Mass = 0.05; // cut on the difference between the K0 mass and the PDG mass bool mRefit = false; // whether to refit or not @@ -429,6 +434,9 @@ class MatchITSTPCQC float mNTPCOccBinLengthInv; std::vector mTBinClOcc; ///< TPC occupancy histo: i-th entry is the integrated occupancy for ~1 orbit starting from the TB = i*mNTPCOccBinLength gsl::span mTPCRefitterOccMap; ///< externally set TPC clusters occupancy map + bool mIsHI = false; + float mK0Scaling = 1.f; // permill that we want to keep of K0S + uint64_t mNK0 = 0; // number of found V0s ClassDefNV(MatchITSTPCQC, 3); }; diff --git a/Detectors/GLOQC/src/MatchITSTPCQC.cxx b/Detectors/GLOQC/src/MatchITSTPCQC.cxx index e95410275a138..81b5743a58893 100644 --- a/Detectors/GLOQC/src/MatchITSTPCQC.cxx +++ b/Detectors/GLOQC/src/MatchITSTPCQC.cxx @@ -47,7 +47,7 @@ MatchITSTPCQC::~MatchITSTPCQC() void MatchITSTPCQC::deleteHistograms() { - LOG(info) << "Deleting histos..."; + LOG(debug) << "Deleting histos..."; for (int i = 0; i < matchType::SIZE; ++i) { // Pt delete mPtNum[i]; @@ -130,7 +130,8 @@ void MatchITSTPCQC::deleteHistograms() delete mFractionITSTPCmatchDCArVsPt; // K0 - delete mK0MassVsPtVsOcc; + delete mK0MassVsPtVsOccpp; + delete mK0MassVsPtVsOccPbPb; } //__________________________________________________________ @@ -211,7 +212,8 @@ void MatchITSTPCQC::reset() // K0 if (mDoK0QC) { - mK0MassVsPtVsOcc->Reset(); + mK0MassVsPtVsOccpp->Reset(); + mK0MassVsPtVsOccPbPb->Reset(); } } @@ -403,17 +405,26 @@ bool MatchITSTPCQC::init() ybinsMassK0[i] = yminMassK0 + i * dyMassK0; } const Int_t nbinsMultK0 = 5; - Double_t* zbinsMultK0 = new Double_t[nbinsMultK0 + 1]; - Double_t zminMultK0 = 100000.; - Double_t zmaxMultK0 = 1000000.; - Double_t dzMultK0 = (zmaxMultK0 - zminMultK0) / nbinsMultK0; + Double_t* zbinsMultK0pp = new Double_t[nbinsMultK0 + 1]; + Double_t* zbinsMultK0PbPb = new Double_t[nbinsMultK0 + 1]; + Double_t zminMultK0pp = 100000.; + Double_t zmaxMultK0pp = 1000000.; + Double_t zminMultK0PbPb = 1.e6; + Double_t zmaxMultK0PbPb = 6.e6; + Double_t dzMultK0pp = (zmaxMultK0pp - zminMultK0pp) / nbinsMultK0; for (int i = 0; i <= nbinsMultK0; i++) { - zbinsMultK0[i] = zminMultK0 + i * dzMultK0; + zbinsMultK0pp[i] = zminMultK0pp + i * dzMultK0pp; + } + Double_t dzMultK0PbPb = (zmaxMultK0PbPb - zminMultK0PbPb) / nbinsMultK0; + for (int i = 0; i <= nbinsMultK0; i++) { + zbinsMultK0PbPb[i] = zminMultK0PbPb + i * dzMultK0PbPb; } if (mDoK0QC) { // V0s - mK0MassVsPtVsOcc = new TH3F("mK0MassVsPt", "K0 invariant mass vs Pt; Pt [GeV/c]; K0s mass [GeV/c^2]", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0, zbinsMultK0); + mK0MassVsPtVsOccpp = new TH3F("mK0MassVsPtVsOccpp", "K0 invariant mass vs Pt; Pt [GeV/c]; K0s mass [GeV/c^2]", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0, zbinsMultK0pp); + + mK0MassVsPtVsOccPbPb = new TH3F("mK0MassVsPtVsOccPbPb", "K0 invariant mass vs Pt; Pt [GeV/c]; K0s mass [GeV/c^2]", nbinsPtK0, xbinsPtK0, nbinsMassK0, ybinsMassK0, nbinsMultK0, zbinsMultK0PbPb); } LOG(info) << "Printing configuration cuts"; @@ -457,6 +468,14 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) // we have not yet initialized the SVertexer params; let's do it ctx.inputs().get("SVParam"); mTimestamp = ctx.services().get().creation; + auto grplhcif = o2::base::GRPGeomHelper::instance().getGRPLHCIF(); + if (grplhcif->getBeamZ(0) != 1 || grplhcif->getBeamZ(1) != 1) { + LOG(info) << "We are in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1); + mIsHI = true; + } + else { + LOG(info) << "We are not in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1); + } } static int evCount = 0; @@ -569,7 +588,7 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) } } } - LOG(info) << "number of entries in map for nominator (without duplicates) = " << mMapLabels.size(); + LOG(debug) << "number of entries in map for nominator (without duplicates) = " << mMapLabels.size(); // now we use only the tracks in the map to fill the histograms (--> tracks have passed the // track selection and there are no duplicated tracks wrt the same MC label) for (int i = 0; i < matchType::SIZE; ++i) { @@ -823,8 +842,8 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) } } } - LOG(info) << "number of entries in map for denominator of TPC tracks (without duplicates) = " << mMapRefLabels[matchType::TPC].size() + mMapLabels[matchType::TPC].size(); - LOG(info) << "number of entries in map for denominator of ITS tracks (without duplicates) = " << mMapRefLabels[matchType::ITS].size() + mMapLabels[matchType::ITS].size(); + LOG(debug) << "number of entries in map for denominator of TPC tracks (without duplicates) = " << mMapRefLabels[matchType::TPC].size() + mMapLabels[matchType::TPC].size(); + LOG(debug) << "number of entries in map for denominator of ITS tracks (without duplicates) = " << mMapRefLabels[matchType::ITS].size() + mMapLabels[matchType::ITS].size(); // now we use only the tracks in the map to fill the histograms (--> tracks have passed the // track selection and there are no duplicated tracks wrt the same MC label) for (auto const& el : mMapRefLabels[matchType::TPC]) { @@ -932,10 +951,10 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) } } - if (mDoK0QC) { + if (mDoK0QC && mRecoCont.getPrimaryVertices().size() > 0) { // now doing K0S const auto pvertices = mRecoCont.getPrimaryVertices(); - LOG(info) << "Found " << pvertices.size() << " primary vertices"; + LOG(debug) << "****** Number of PVs = " << pvertices.size(); // getting occupancy estimator mNHBPerTF = o2::base::GRPGeomHelper::instance().getGRPECS()->getNHBFPerTF(); @@ -951,6 +970,7 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) std::vector mTBinClOcc; ///< TPC occupancy histo: i-th entry is the integrated occupancy for ~1 orbit starting from the TB = i * mNTPCOccBinLength mTBinClOcc.clear(); int mNTPCOccBinLength = mParam->rec.tpc.occupancyMapTimeBins; + LOG(debug) << "mNTPCOccBinLength = " << mNTPCOccBinLength; mNTPCOccBinLengthInv = 1. / mNTPCOccBinLength; if (mNTPCOccBinLength > 1 && mTPCRefitterOccMap.size()) { int nTPCBinsInTF = mNHBPerTF * o2::constants::lhc::LHCMaxBunches / 8; // number of TPC time bins in 1 TF, considering that 1 TPC time bin is 8 bunches @@ -1014,6 +1034,7 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) } int nV0sOk = 0; // processing every sec vtx for each prim vtx + int myCount = 0; for (auto it : pv2sv) { int pvID = it.first; auto& vv = it.second; @@ -1025,6 +1046,7 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) for (int iv0 : vv) { nV0sOk += processV0(iv0, mRecoCont, mTBinClOcc, pvTime) ? 1 : 0; } + ++myCount; } LOG(debug) << "Processed " << nV0sOk << " V0s"; @@ -1041,6 +1063,14 @@ bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoDat static int tfID = 0; const auto& v0id = v0IDs[iv]; + ++mNK0; + if (mNK0 % int(1/mK0Scaling) == 0) { + LOG(debug) << "Checking " << mNK0 << "th V0: refitting it, since we keep " << mK0Scaling * 100 << "% of all V0s"; + } + else { + LOG(debug) << "Checking " << mNK0 << "th K0: NOT refitting it, but skipping it, since we keep " << mK0Scaling * 100 << "% of all V0s"; + return false; + } if (mRefit && !refitV0(v0id, v0, recoData)) { return false; } @@ -1055,14 +1085,21 @@ bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoDat int tb = pvTime / (8 * o2::constants::lhc::LHCBunchSpacingMUS) * mNTPCOccBinLengthInv; // V0 time in TPC time bins LOG(debug) << "pvTime = " << pvTime << " tb = " << tb; float mltTPC = tb < 0 ? mTBinClOcc[0] : (tb >= mTBinClOcc.size() ? mTBinClOcc.back() : mTBinClOcc[tb]); - LOG(debug) << "Filling plot with pt = " << v0sel.getPt() << " mass = " << std::sqrt(v0sel.calcMass2AsK0()) << " mult TPC = " << mltTPC; - mK0MassVsPtVsOcc->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); + ++mNK0; + LOG(debug) << "Filling K0 histogram with pt = " << v0sel.getPt() << " mass = " << std::sqrt(v0sel.calcMass2AsK0()) << " mult TPC = " << mltTPC; + if (!mIsHI) { + mK0MassVsPtVsOccpp->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); + } + else { + mK0MassVsPtVsOccPbPb->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); + } return true; } //__________________________________________________________ bool MatchITSTPCQC::refitV0(const o2::dataformats::V0Index& id, o2::dataformats::V0& v0, o2::globaltracking::RecoContainer& recoData) { + LOG(debug) << "Refitting V0"; if (!recoData.isTrackSourceLoaded(id.getProngID(0).getSource()) || !recoData.isTrackSourceLoaded(id.getProngID(1).getSource())) { return false; } @@ -1198,10 +1235,10 @@ void MatchITSTPCQC::setEfficiency(TEfficiency* eff, TH1* hnum, TH1* hden, bool i // we need to force to replace the total histogram, otherwise it will compare it to the previous passed one, and it might get an error of inconsistency in the bin contents if constexpr (false) { // checking - bool bad{false}; - LOG(info) << "Setting efficiency " << eff->GetName() << " from " << hnum->GetName() << " and " << hden->GetName(); - LOG(info) << "Num " << hnum->GetName() << " " << hnum->GetNbinsX() << " " << hnum->GetNbinsY() << " with " << hnum->GetEntries() << " entries"; - LOG(info) << "Den " << hden->GetName() << " " << hden->GetNbinsX() << " " << hden->GetNbinsY() << " with " << hden->GetEntries() << " entries"; + bool bad{false}; + LOG(debug) << "Setting efficiency " << eff->GetName() << " from " << hnum->GetName() << " and " << hden->GetName(); + LOG(debug) << "Num " << hnum->GetName() << " " << hnum->GetNbinsX() << " " << hnum->GetNbinsY() << " with " << hnum->GetEntries() << " entries"; + LOG(debug) << "Den " << hden->GetName() << " " << hden->GetNbinsX() << " " << hden->GetNbinsY() << " with " << hden->GetEntries() << " entries"; if (hnum->GetDimension() != hden->GetDimension()) { LOGP(warning, "Histograms have different dimensions (num={} to den={})", hnum->GetDimension(), hden->GetDimension()); bad = true; @@ -1332,5 +1369,6 @@ void MatchITSTPCQC::getHistos(TObjArray& objar) objar.Add(mFractionITSTPCmatchDCArVsPt); // V0 - objar.Add(mK0MassVsPtVsOcc); + objar.Add(mK0MassVsPtVsOccpp); + objar.Add(mK0MassVsPtVsOccPbPb); } From 3d06bbb3e2bdc7056a39378742f00ba09707fcf2 Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Thu, 3 Oct 2024 19:49:01 +0000 Subject: [PATCH 4/4] Please consider the following formatting changes --- Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h | 4 ++-- Detectors/GLOQC/src/MatchITSTPCQC.cxx | 13 +++++-------- 2 files changed, 7 insertions(+), 10 deletions(-) diff --git a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h index 1aaf34bc244fb..00862dd3c734d 100644 --- a/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h +++ b/Detectors/GLOQC/include/GLOQC/MatchITSTPCQC.h @@ -435,8 +435,8 @@ class MatchITSTPCQC std::vector mTBinClOcc; ///< TPC occupancy histo: i-th entry is the integrated occupancy for ~1 orbit starting from the TB = i*mNTPCOccBinLength gsl::span mTPCRefitterOccMap; ///< externally set TPC clusters occupancy map bool mIsHI = false; - float mK0Scaling = 1.f; // permill that we want to keep of K0S - uint64_t mNK0 = 0; // number of found V0s + float mK0Scaling = 1.f; // permill that we want to keep of K0S + uint64_t mNK0 = 0; // number of found V0s ClassDefNV(MatchITSTPCQC, 3); }; diff --git a/Detectors/GLOQC/src/MatchITSTPCQC.cxx b/Detectors/GLOQC/src/MatchITSTPCQC.cxx index 81b5743a58893..d010d7dab96b5 100644 --- a/Detectors/GLOQC/src/MatchITSTPCQC.cxx +++ b/Detectors/GLOQC/src/MatchITSTPCQC.cxx @@ -472,8 +472,7 @@ void MatchITSTPCQC::run(o2::framework::ProcessingContext& ctx) if (grplhcif->getBeamZ(0) != 1 || grplhcif->getBeamZ(1) != 1) { LOG(info) << "We are in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1); mIsHI = true; - } - else { + } else { LOG(info) << "We are not in Heavy Ion: Z for beam 0 = " << grplhcif->getBeamZ(0) << " ; Z for beam 1 = " << grplhcif->getBeamZ(1); } } @@ -1064,10 +1063,9 @@ bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoDat const auto& v0id = v0IDs[iv]; ++mNK0; - if (mNK0 % int(1/mK0Scaling) == 0) { + if (mNK0 % int(1 / mK0Scaling) == 0) { LOG(debug) << "Checking " << mNK0 << "th V0: refitting it, since we keep " << mK0Scaling * 100 << "% of all V0s"; - } - else { + } else { LOG(debug) << "Checking " << mNK0 << "th K0: NOT refitting it, but skipping it, since we keep " << mK0Scaling * 100 << "% of all V0s"; return false; } @@ -1089,8 +1087,7 @@ bool MatchITSTPCQC::processV0(int iv, o2::globaltracking::RecoContainer& recoDat LOG(debug) << "Filling K0 histogram with pt = " << v0sel.getPt() << " mass = " << std::sqrt(v0sel.calcMass2AsK0()) << " mult TPC = " << mltTPC; if (!mIsHI) { mK0MassVsPtVsOccpp->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); - } - else { + } else { mK0MassVsPtVsOccPbPb->Fill(v0sel.getPt(), std::sqrt(v0sel.calcMass2AsK0()), mltTPC); } return true; @@ -1235,7 +1232,7 @@ void MatchITSTPCQC::setEfficiency(TEfficiency* eff, TH1* hnum, TH1* hden, bool i // we need to force to replace the total histogram, otherwise it will compare it to the previous passed one, and it might get an error of inconsistency in the bin contents if constexpr (false) { // checking - bool bad{false}; + bool bad{false}; LOG(debug) << "Setting efficiency " << eff->GetName() << " from " << hnum->GetName() << " and " << hden->GetName(); LOG(debug) << "Num " << hnum->GetName() << " " << hnum->GetNbinsX() << " " << hnum->GetNbinsY() << " with " << hnum->GetEntries() << " entries"; LOG(debug) << "Den " << hden->GetName() << " " << hden->GetNbinsX() << " " << hden->GetNbinsY() << " with " << hden->GetEntries() << " entries";