Skip to content

Please consider the following formatting changes to #14737 #1

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
alibuild wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Please consider the following formatting changes to #14737 #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
67 changes: 35 additions & 32 deletions PWGDQ/Core/VarManager.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1104,9 +1104,10 @@ class VarManager : public TObject
return false;
}

// Flag to set PV recalculation via KF
static void SetPVrecalculationKF(const bool pvRecalKF) {
fgPVrecalKF = pvRecalKF;
// Flag to set PV recalculation via KF
static void SetPVrecalculationKF(const bool pvRecalKF)
{
fgPVrecalKF = pvRecalKF;
}

// Setup the collision system
Expand Down Expand Up @@ -1311,7 +1312,7 @@ class VarManager : public TObject
static void FillQuadMC(T1 const& t1, T2 const& t2, T2 const& t3, float* values = nullptr);
template <int pairType, uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
static void FillPairVertexing(C const& collision, T const& t1, T const& t2, bool propToSV = false, float* values = nullptr);
template <int pairType, uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
template <int pairType, uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
static void FillPairVertexingRecomputePV(C const& /*collision*/, T const& t1, T const& t2, o2::dataformats::VertexBase pvRefitted, float* values = nullptr);
template <uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
static void FillTripletVertexing(C const& collision, T const& t1, T const& t2, T const& t3, PairCandidateType tripletType, float* values = nullptr);
Expand Down Expand Up @@ -4202,7 +4203,7 @@ void VarManager::FillPairVertexing(C const& collision, T const& t1, T const& t2,

if constexpr ((pairType == kDecayToEE || pairType == kDecayToKPi) && trackHasCov) {
secondaryVertex = fgFitterTwoProngBarrel.getPCACandidate();
// printf("secVtx (first) %f %f %f \n",secondaryVertex[0],secondaryVertex[1],secondaryVertex[2]);
// printf("secVtx (first) %f %f %f \n",secondaryVertex[0],secondaryVertex[1],secondaryVertex[2]);
covMatrixPCA = fgFitterTwoProngBarrel.calcPCACovMatrixFlat();
auto chi2PCA = fgFitterTwoProngBarrel.getChi2AtPCACandidate();
auto trackParVar0 = fgFitterTwoProngBarrel.getTrack(0);
Expand All @@ -4211,7 +4212,8 @@ void VarManager::FillPairVertexing(C const& collision, T const& t1, T const& t2,
v1 = {trackParVar0.getPt(), trackParVar0.getEta(), trackParVar0.getPhi(), m1};
v2 = {trackParVar1.getPt(), trackParVar1.getEta(), trackParVar1.getPhi(), m2};
v12 = v1 + v2;
if(fgPVrecalKF) primaryVertexNew = RecalculatePrimaryVertex(t1, t2, collision);
if (fgPVrecalKF)
primaryVertexNew = RecalculatePrimaryVertex(t1, t2, collision);

} else if constexpr (pairType == kDecayToMuMu && muonHasCov) {
// Get pca candidate from forward DCA fitter
Expand Down Expand Up @@ -4269,14 +4271,15 @@ void VarManager::FillPairVertexing(C const& collision, T const& t1, T const& t2,
values[kVertexingLxyProjected] = values[kVertexingLxyProjected] / TMath::Sqrt((v12.Px() * v12.Px()) + (v12.Py() * v12.Py()));
values[kVertexingLxyzProjected] = ((secondaryVertex[0] - collision.posX()) * v12.Px()) + ((secondaryVertex[1] - collision.posY()) * v12.Py()) + ((secondaryVertex[2] - collision.posZ()) * v12.Pz());
values[kVertexingLxyzProjected] = values[kVertexingLxyzProjected] / TMath::Sqrt((v12.Px() * v12.Px()) + (v12.Py() * v12.Py()) + (v12.Pz() * v12.Pz()));
if(fgPVrecalKF){
values[kVertexingLxyProjectedRecalculatePV] = (secondaryVertex[0] - primaryVertexNew.getX()) * v12.Px() + (secondaryVertex[1] - primaryVertexNew.getY()) * v12.Py();
values[kVertexingLxyProjectedRecalculatePV] = values[kVertexingLxyProjectedRecalculatePV] / v12.Pt();
if (fgPVrecalKF) {
values[kVertexingLxyProjectedRecalculatePV] = (secondaryVertex[0] - primaryVertexNew.getX()) * v12.Px() + (secondaryVertex[1] - primaryVertexNew.getY()) * v12.Py();
values[kVertexingLxyProjectedRecalculatePV] = values[kVertexingLxyProjectedRecalculatePV] / v12.Pt();
}
values[kVertexingTauxyProjected] = values[kVertexingLxyProjected] * v12.M() / (v12.Pt());
values[kVertexingTauxyProjectedPoleJPsiMass] = values[kVertexingLxyProjected] * o2::constants::physics::MassJPsi / (v12.Pt());
values[kVertexingTauxyProjectedNs] = values[kVertexingTauxyProjected] / o2::constants::physics::LightSpeedCm2NS;
if(fgPVrecalKF) values[kVertexingTauxyProjectedPoleJPsiMassRecalculatePV] = values[kVertexingLxyProjectedRecalculatePV] * o2::constants::physics::MassJPsi / (v12.Pt());
if (fgPVrecalKF)
values[kVertexingTauxyProjectedPoleJPsiMassRecalculatePV] = values[kVertexingLxyProjectedRecalculatePV] * o2::constants::physics::MassJPsi / (v12.Pt());
values[kVertexingTauzProjected] = values[kVertexingLzProjected] * v12.M() / TMath::Abs(v12.Pz());
values[kVertexingTauxyzProjected] = values[kVertexingLxyzProjected] * v12.M() / (v12.P());
}
Expand Down Expand Up @@ -4495,37 +4498,38 @@ void VarManager::FillPairVertexing(C const& collision, T const& t1, T const& t2,
}

template <int pairType, uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
void VarManager::FillPairVertexingRecomputePV(C const& /*collision*/, T const& t1, T const& t2, o2::dataformats::VertexBase pvRefitted, float* values )
{
// recompute decay lenght variables using updated primary vertex
void VarManager::FillPairVertexingRecomputePV(C const& /*collision*/, T const& t1, T const& t2, o2::dataformats::VertexBase pvRefitted, float* values)
{
// recompute decay lenght variables using updated primary vertex

// check at compile time that the event and cov matrix have the cov matrix
constexpr bool eventHasVtxCov = ((collFillMap & Collision) > 0 || (collFillMap & ReducedEventVtxCov) > 0);
constexpr bool trackHasCov = ((fillMap & TrackCov) > 0 || (fillMap & ReducedTrackBarrelCov) > 0);
constexpr bool muonHasCov = ((fillMap & MuonCov) > 0 || (fillMap & ReducedMuonCov) > 0);
// check at compile time that the event and cov matrix have the cov matrix
constexpr bool eventHasVtxCov = ((collFillMap & Collision) > 0 || (collFillMap & ReducedEventVtxCov) > 0);
constexpr bool trackHasCov = ((fillMap & TrackCov) > 0 || (fillMap & ReducedTrackBarrelCov) > 0);
constexpr bool muonHasCov = ((fillMap & MuonCov) > 0 || (fillMap & ReducedMuonCov) > 0);

if (!values) {
if (!values) {
values = fgValues;
}
}

float m1 = o2::constants::physics::MassElectron;
float m2 = o2::constants::physics::MassElectron;
if constexpr (pairType == kDecayToKPi) {
float m1 = o2::constants::physics::MassElectron;
float m2 = o2::constants::physics::MassElectron;
if constexpr (pairType == kDecayToKPi) {
m1 = o2::constants::physics::MassKaonCharged;
m2 = o2::constants::physics::MassPionCharged;
}
if constexpr (pairType == kDecayToMuMu && muonHasCov) {
}
if constexpr (pairType == kDecayToMuMu && muonHasCov) {
m1 = o2::constants::physics::MassMuon;
m2 = o2::constants::physics::MassMuon;
}
ROOT::Math::PtEtaPhiMVector v1(t1.pt(), t1.eta(), t1.phi(), m1);
ROOT::Math::PtEtaPhiMVector v2(t2.pt(), t2.eta(), t2.phi(), m2);
ROOT::Math::PtEtaPhiMVector v12 = v1 + v2;
}
ROOT::Math::PtEtaPhiMVector v1(t1.pt(), t1.eta(), t1.phi(), m1);
ROOT::Math::PtEtaPhiMVector v2(t2.pt(), t2.eta(), t2.phi(), m2);
ROOT::Math::PtEtaPhiMVector v12 = v1 + v2;

if(fgFitterTwoProngBarrel.getNCandidates() == 0) return;
Vec3D secondaryVertex;
if (fgFitterTwoProngBarrel.getNCandidates() == 0)
return;
Vec3D secondaryVertex;

if (!fgUsedKF) { // to be updated when seconday vertex is computed with KF
if (!fgUsedKF) { // to be updated when seconday vertex is computed with KF
if constexpr (eventHasVtxCov) {

if constexpr ((pairType == kDecayToEE || pairType == kDecayToKPi) && trackHasCov) {
Expand Down Expand Up @@ -4554,7 +4558,6 @@ void VarManager::FillPairVertexingRecomputePV(C const& /*collision*/, T const& t
values[kVertexingTauxyProjectedPoleJPsiMassRecalculatePV] = values[kVertexingLxyProjectedRecalculatePV] * o2::constants::physics::MassJPsi / (v12.Pt());
}
}

}

template <uint32_t collFillMap, uint32_t fillMap, typename C, typename T>
Expand Down
122 changes: 62 additions & 60 deletions PWGDQ/Tasks/dqEfficiency_withAssoc_direct.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,9 +23,6 @@
#include "PWGDQ/Core/MixingLibrary.h"
#include "PWGDQ/Core/VarManager.h"
#include "PWGDQ/DataModel/ReducedInfoTables.h"
#include <DetectorsVertexing/PVertexer.h>
#include <ReconstructionDataFormats/Track.h>
#include <Common/Core/trackUtilities.h>

#include "Common/Core/PID/PIDTOFParamService.h"
#include "Common/Core/TableHelper.h"
Expand All @@ -34,6 +31,7 @@
#include "Common/DataModel/McCollisionExtra.h"
#include "Common/DataModel/Multiplicity.h"
#include "Common/DataModel/TrackSelectionTables.h"
#include <Common/Core/trackUtilities.h>

#include "CCDB/BasicCCDBManager.h"
#include "DataFormatsParameters/GRPMagField.h"
Expand All @@ -46,6 +44,8 @@
#include "Framework/AnalysisHelpers.h"
#include "Framework/AnalysisTask.h"
#include "Framework/runDataProcessing.h"
#include <DetectorsVertexing/PVertexer.h>
#include <ReconstructionDataFormats/Track.h>

#include "TGeoGlobalMagField.h"
#include <TH1F.h>
Expand Down Expand Up @@ -1163,7 +1163,6 @@ struct AnalysisSameEventPairing {
std::vector<o2::track::TrackParCov> pvContribTrackPars;
std::vector<bool> vec_useTrk_PVrefit;


// keep histogram class names in maps, so we don't have to buld their names in the pair loops
std::map<int, std::vector<TString>> fTrackHistNames;
std::map<int, std::vector<TString>> fBarrelHistNamesMCmatched;
Expand Down Expand Up @@ -1509,17 +1508,16 @@ struct AnalysisSameEventPairing {
cout << "AnalysisSameEventPairing::init() completed" << endl;
}


void initParamsFromCCDB(uint64_t timestamp, bool withTwoProngFitter = true)
{
cout << "AnalysisSameEventPairing::initParamsFromCCDB() called for timestamp " << timestamp << endl;
if (fConfigOptions.useRemoteField.value) {
o2::parameters::GRPMagField* grpmag = fCCDB->getForTimeStamp<o2::parameters::GRPMagField>(fConfigCCDB.grpMagPath, timestamp);
o2::base::MatLayerCylSet* lut = o2::base::MatLayerCylSet::rectifyPtrFromFile(fCCDB->get<o2::base::MatLayerCylSet>(fConfigCCDB.lutPath));
o2::base::MatLayerCylSet* lut = o2::base::MatLayerCylSet::rectifyPtrFromFile(fCCDB->get<o2::base::MatLayerCylSet>(fConfigCCDB.lutPath));
float magField = 0.0;
if (grpmag != nullptr) {
magField = grpmag->getNominalL3Field();
o2::base::Propagator::initFieldFromGRP(grpmag);
o2::base::Propagator::initFieldFromGRP(grpmag);
o2::base::Propagator::Instance()->setMatLUT(lut);
} else {
LOGF(fatal, "GRP object is not available in CCDB at timestamp=%llu", timestamp);
Expand Down Expand Up @@ -1549,53 +1547,56 @@ struct AnalysisSameEventPairing {
cout << "AnalysisSameEventPairing::initParamsFromCCDB() completed" << endl;
}

template <typename Events, typename TTracks, typename Tracks>
bool refitPVWithPVertexer(Events const& collision, TTracks const& tracks, Tracks const& t1, Tracks const& t2, o2::dataformats::VertexBase& pvRefitted)
{
// --- build PV contributor list ---
pvContribGlobIDs.clear();
pvContribTrackPars.clear();
// int nMyPVContrib = 0; int nMyPVContribOrig = 0;
for (auto const& trk : tracks) {
// check if it is PV contributor
if (!trk.isPVContributor()) continue;
// check if it contributes to the vtx of this collision
if (trk.collisionId() != collision.globalIndex()) continue;
// nMyPVContribOrig++;
// --- remove t1 and t2 if they are PV contributors ---
if (trk.globalIndex() == t1.globalIndex() || trk.globalIndex() == t2.globalIndex()) continue;
// add tracks and parameters to the list
pvContribGlobIDs.push_back(trk.globalIndex());
pvContribTrackPars.push_back(getTrackParCov(trk));
// nMyPVContrib++;
}
template <typename Events, typename TTracks, typename Tracks>
bool refitPVWithPVertexer(Events const& collision, TTracks const& tracks, Tracks const& t1, Tracks const& t2, o2::dataformats::VertexBase& pvRefitted)
{
// --- build PV contributor list ---
pvContribGlobIDs.clear();
pvContribTrackPars.clear();
// int nMyPVContrib = 0; int nMyPVContribOrig = 0;
for (auto const& trk : tracks) {
// check if it is PV contributor
if (!trk.isPVContributor())
continue;
// check if it contributes to the vtx of this collision
if (trk.collisionId() != collision.globalIndex())
continue;
// nMyPVContribOrig++;
// --- remove t1 and t2 if they are PV contributors ---
if (trk.globalIndex() == t1.globalIndex() || trk.globalIndex() == t2.globalIndex())
continue;
// add tracks and parameters to the list
pvContribGlobIDs.push_back(trk.globalIndex());
pvContribTrackPars.push_back(getTrackParCov(trk));
// nMyPVContrib++;
}

// cout << "contributors from collision: " << collision.numContrib() << " - from refitting: before -> " << nMyPVContribOrig << " after -> " << nMyPVContrib << endl;
vec_useTrk_PVrefit.assign(pvContribGlobIDs.size(), true);
// --- build VertexBase from event collision ---
o2::dataformats::VertexBase Pvtx;
Pvtx.setX(collision.posX());
Pvtx.setY(collision.posY());
Pvtx.setZ(collision.posZ());
Pvtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(),
collision.covXZ(), collision.covYZ(), collision.covZZ());

// --- configure vertexer ---
o2::vertexing::PVertexer vertexer;
if (fConfigOptions.removeDiamondConstrPV) {
o2::conf::ConfigurableParam::updateFromString("pvertexer.useMeanVertexConstraint=false");
}
vertexer.init();
// cout << "contributors from collision: " << collision.numContrib() << " - from refitting: before -> " << nMyPVContribOrig << " after -> " << nMyPVContrib << endl;
vec_useTrk_PVrefit.assign(pvContribGlobIDs.size(), true);
// --- build VertexBase from event collision ---
o2::dataformats::VertexBase Pvtx;
Pvtx.setX(collision.posX());
Pvtx.setY(collision.posY());
Pvtx.setZ(collision.posZ());
Pvtx.setCov(collision.covXX(), collision.covXY(), collision.covYY(),
collision.covXZ(), collision.covYZ(), collision.covZZ());

bool PVrefit_doable = vertexer.prepareVertexRefit(pvContribTrackPars, Pvtx);
if (!PVrefit_doable) return false;
// --- configure vertexer ---
o2::vertexing::PVertexer vertexer;
if (fConfigOptions.removeDiamondConstrPV) {
o2::conf::ConfigurableParam::updateFromString("pvertexer.useMeanVertexConstraint=false");
}
vertexer.init();

// --- do the refit ---
pvRefitted = vertexer.refitVertex(vec_useTrk_PVrefit, Pvtx);
bool PVrefit_doable = vertexer.prepareVertexRefit(pvContribTrackPars, Pvtx);
if (!PVrefit_doable)
return false;

return true;
}
// --- do the refit ---
pvRefitted = vertexer.refitVertex(vec_useTrk_PVrefit, Pvtx);

return true;
}

// Template function to run same event pairing (barrel-barrel, muon-muon, barrel-muon)
template <bool TTwoProngFitter, int TPairType, uint32_t TEventFillMap, uint32_t TTrackFillMap, typename TEvents, typename TTracks>
Expand Down Expand Up @@ -1710,17 +1711,18 @@ bool refitPVWithPVertexer(Events const& collision, TTracks const& tracks, Track
}
if constexpr (TTwoProngFitter) {
VarManager::FillPairVertexing<TPairType, TEventFillMap, TTrackFillMap>(event, t1, t2, fConfigOptions.propToPCA);
if(fConfigOptions.recomputePV){
VarManager::SetPVrecalculationKF(false);
VarManager::ResetValues(VarManager::kVertexingLxyProjectedRecalculatePV, VarManager::kVertexingLxyProjectedRecalculatePV+1);
VarManager::ResetValues(VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV, VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV+1);
// cout << "primary vertex (before): x -> " << event.posX() << " y -> " << event.posY() << " z -> " << event.posZ() << endl;
o2::dataformats::VertexBase pvRefit;
bool ok = refitPVWithPVertexer(event, tracks, t1, t2, pvRefit);
if(ok) VarManager::FillPairVertexingRecomputePV<TPairType, TEventFillMap, TTrackFillMap>(event, t1, t2, pvRefit);
// cout << "primary vertex (after): ok -> " << ok << " x -> " << pvRefit.getX() << " y -> " << pvRefit.getY() << " z -> " << pvRefit.getZ() << endl;
if (fConfigOptions.recomputePV) {
VarManager::SetPVrecalculationKF(false);
VarManager::ResetValues(VarManager::kVertexingLxyProjectedRecalculatePV, VarManager::kVertexingLxyProjectedRecalculatePV + 1);
VarManager::ResetValues(VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV, VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV + 1);
// cout << "primary vertex (before): x -> " << event.posX() << " y -> " << event.posY() << " z -> " << event.posZ() << endl;
o2::dataformats::VertexBase pvRefit;
bool ok = refitPVWithPVertexer(event, tracks, t1, t2, pvRefit);
if (ok)
VarManager::FillPairVertexingRecomputePV<TPairType, TEventFillMap, TTrackFillMap>(event, t1, t2, pvRefit);
// cout << "primary vertex (after): ok -> " << ok << " x -> " << pvRefit.getX() << " y -> " << pvRefit.getY() << " z -> " << pvRefit.getZ() << endl;
}
}
}
if constexpr (eventHasQvector) {
VarManager::FillPairVn<TPairType>(t1, t2);
}
Expand Down Expand Up @@ -1858,7 +1860,7 @@ bool refitPVWithPVertexer(Events const& collision, TTracks const& tracks, Track
fHistMan->FillHistClass(histNames[icut][0].Data(), VarManager::fgValues); // reconstructed, unmatched
for (unsigned int isig = 0; isig < fRecMCSignals.size(); isig++) { // loop over MC signals
if (mcDecision & (static_cast<uint32_t>(1) << isig)) {
PromptNonPromptSepTable(VarManager::fgValues[VarManager::kMass], VarManager::fgValues[VarManager::kPt], VarManager::fgValues[VarManager::kEta], VarManager::fgValues[VarManager::kRap], VarManager::fgValues[VarManager::kPhi], VarManager::fgValues[VarManager::kVertexingTauxyProjected], VarManager::fgValues[VarManager::kVertexingTauxyProjectedPoleJPsiMass], VarManager::fgValues[VarManager::kVertexingTauzProjected], VarManager::fgValues[VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV], isAmbiInBunch, isAmbiOutOfBunch, isCorrect_pair, VarManager::fgValues[VarManager::kMultFT0A], VarManager::fgValues[VarManager::kMultFT0C], VarManager::fgValues[VarManager::kCentFT0M], VarManager::fgValues[VarManager::kVtxNcontribReal]);
PromptNonPromptSepTable(VarManager::fgValues[VarManager::kMass], VarManager::fgValues[VarManager::kPt], VarManager::fgValues[VarManager::kEta], VarManager::fgValues[VarManager::kRap], VarManager::fgValues[VarManager::kPhi], VarManager::fgValues[VarManager::kVertexingTauxyProjected], VarManager::fgValues[VarManager::kVertexingTauxyProjectedPoleJPsiMass], VarManager::fgValues[VarManager::kVertexingTauzProjected], VarManager::fgValues[VarManager::kVertexingTauxyProjectedPoleJPsiMassRecalculatePV], isAmbiInBunch, isAmbiOutOfBunch, isCorrect_pair, VarManager::fgValues[VarManager::kMultFT0A], VarManager::fgValues[VarManager::kMultFT0C], VarManager::fgValues[VarManager::kCentFT0M], VarManager::fgValues[VarManager::kVtxNcontribReal]);
fHistMan->FillHistClass(histNamesMC[icut * fRecMCSignals.size() + isig][0].Data(), VarManager::fgValues); // matched signal
/*if (fConfigOptions.fConfigMiniTree) {
if constexpr (TPairType == VarManager::kDecayToMuMu) {
Expand Down
Loading