Skip to content

Commit 5d054d8

Browse files
dsekihatdsekihat
committed
PWGEM/Dilepton, Common: update global muon propagation
1 parent d2fee7c commit 5d054d8

File tree

6 files changed

+259
-369
lines changed

6 files changed

+259
-369
lines changed

‎Common/Core/fwdtrackUtilities.h‎

Lines changed: 15 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -18,17 +18,17 @@
1818
#ifndef COMMON_CORE_FWDTRACKUTILITIES_H_
1919
#define COMMON_CORE_FWDTRACKUTILITIES_H_
2020

21+
#include "DetectorsBase/GeometryManager.h"
22+
#include "Field/MagneticField.h"
2123
#include "Framework/AnalysisDataModel.h"
22-
#include <DetectorsBase/GeometryManager.h>
23-
#include <Field/MagneticField.h>
24-
#include <GlobalTracking/MatchGlobalFwd.h>
25-
#include <MCHTracking/TrackExtrap.h>
26-
#include <ReconstructionDataFormats/GlobalFwdTrack.h>
27-
#include <ReconstructionDataFormats/TrackFwd.h>
24+
#include "GlobalTracking/MatchGlobalFwd.h"
25+
#include "MCHTracking/TrackExtrap.h"
26+
#include "ReconstructionDataFormats/GlobalFwdTrack.h"
27+
#include "ReconstructionDataFormats/TrackFwd.h"
2828

29-
#include <Math/MatrixRepresentationsStatic.h>
30-
#include <Math/SMatrix.h>
31-
#include <TGeoGlobalMagField.h>
29+
#include "Math/MatrixRepresentationsStatic.h"
30+
#include "Math/SMatrix.h"
31+
#include "TGeoGlobalMagField.h"
3232

3333
#include <type_traits>
3434
#include <vector>
@@ -51,8 +51,7 @@ using SMatrix5 = ROOT::Math::SVector<double, 5>;
5151
template <typename TFwdTrack, typename TFwdTrackCov>
5252
o2::track::TrackParCovFwd getTrackParCovFwd(TFwdTrack const& track, TFwdTrackCov const& cov)
5353
{
54-
// This function works for both (saMuon, saMuon) and (MFTTrack, MFTTrackCov).
55-
// Don't use covariant matrix of global muons stored in AO2D.root.
54+
// This function works for (glMuon, glMuon), (saMuon, saMuon) and (MFTTrack, MFTTrackCov).
5655

5756
double chi2 = track.chi2();
5857
if constexpr (std::is_same_v<std::decay_t<TFwdTrackCov>, aod::MFTTracksCov::iterator>) {
@@ -128,11 +127,12 @@ o2::dataformats::GlobalFwdTrack propagateTrackParCovFwd(TFwdTrackParCov const& f
128127
// o2::field::MagneticField* field = static_cast<o2::field::MagneticField*>(TGeoGlobalMagField::Instance()->GetField());
129128
// auto Bz = field->getBz(centerMFT); // Get field at centre of MFT in kG.
130129

131-
auto geoMan = o2::base::GeometryManager::meanMaterialBudget(fwdtrack.getX(), fwdtrack.getY(), fwdtrack.getZ(), collision.posX(), collision.posY(), collision.posZ());
132-
auto x2x0 = static_cast<float>(geoMan.meanX2X0);
133-
134130
if (endPoint == propagationPoint::kToVertex) {
135-
fwdtrack.propagateToVtxhelixWithMCS(collision.posZ(), {collision.posX(), collision.posY()}, {collision.covXX(), collision.covYY()}, bzkG, x2x0);
131+
// auto geoMan = o2::base::GeometryManager::meanMaterialBudget(fwdtrack.getX(), fwdtrack.getY(), fwdtrack.getZ(), collision.posX(), collision.posY(), collision.posZ());
132+
// auto x2x0 = static_cast<float>(geoMan.meanX2X0);
133+
// fwdtrack.propagateToVtxhelixWithMCS(collision.posZ(), {collision.posX(), collision.posY()}, {collision.covXX(), collision.covYY()}, bzkG, x2x0);
134+
std::array<double, 3> dcaInfOrig{999.f, 999.f, 999.f};
135+
fwdtrack.propagateToDCAhelix(bzkG, {collision.posX(), collision.posY(), collision.posZ()}, dcaInfOrig);
136136
} else if (endPoint == propagationPoint::kToDCA) {
137137
fwdtrack.propagateToZhelix(collision.posZ(), bzkG);
138138
} else if (endPoint == propagationPoint::kToMatchingPlane) {

‎PWGEM/Dilepton/Core/Dilepton.h‎

Lines changed: 69 additions & 136 deletions
Original file line numberDiff line numberDiff line change
@@ -418,111 +418,84 @@ struct Dilepton {
418418

419419
emh_pos = new TEMH(ndepth);
420420
emh_neg = new TEMH(ndepth);
421-
emh_pair_uls = new MyEMH_pair(ndepth);
422-
emh_pair_lspp = new MyEMH_pair(ndepth);
423-
emh_pair_lsmm = new MyEMH_pair(ndepth);
424421

425-
if (accBins.ConfMllAccBins.value[0] == VARIABLE_WIDTH) {
426-
mll_bin_edges = std::vector<float>(accBins.ConfMllAccBins.value.begin(), accBins.ConfMllAccBins.value.end());
427-
mll_bin_edges.erase(mll_bin_edges.begin());
428-
for (const auto& edge : mll_bin_edges) {
429-
LOGF(info, "VARIABLE_WIDTH: mll_bin_edges = %f", edge);
430-
}
431-
} else {
432-
int nbins = static_cast<int>(accBins.ConfMllAccBins.value[0]);
433-
float xmin = static_cast<float>(accBins.ConfMllAccBins.value[1]);
434-
float xmax = static_cast<float>(accBins.ConfMllAccBins.value[2]);
435-
mll_bin_edges.resize(nbins + 1);
436-
for (int i = 0; i < nbins + 1; i++) {
437-
mll_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
438-
LOGF(info, "FIXED_WIDTH: mll_bin_edges[%d] = %f", i, mll_bin_edges[i]);
439-
}
440-
}
422+
if (cfgAnalysisType == static_cast<int>(o2::aod::pwgem::dilepton::utils::pairutil::DileptonAnalysisType::kPolarization)) { // only for polarization
441423

442-
if (accBins.ConfPtllAccBins.value[0] == VARIABLE_WIDTH) {
443-
ptll_bin_edges = std::vector<float>(accBins.ConfPtllAccBins.value.begin(), accBins.ConfPtllAccBins.value.end());
444-
ptll_bin_edges.erase(ptll_bin_edges.begin());
445-
for (const auto& edge : ptll_bin_edges) {
446-
LOGF(info, "VARIABLE_WIDTH: ptll_bin_edges = %f", edge);
447-
}
448-
} else {
449-
int nbins = static_cast<int>(accBins.ConfPtllAccBins.value[0]);
450-
float xmin = static_cast<float>(accBins.ConfPtllAccBins.value[1]);
451-
float xmax = static_cast<float>(accBins.ConfPtllAccBins.value[2]);
452-
ptll_bin_edges.resize(nbins + 1);
453-
for (int i = 0; i < nbins + 1; i++) {
454-
ptll_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
455-
LOGF(info, "FIXED_WIDTH: ptll_bin_edges[%d] = %f", i, ptll_bin_edges[i]);
456-
}
457-
}
424+
emh_pair_uls = new MyEMH_pair(ndepth);
425+
emh_pair_lspp = new MyEMH_pair(ndepth);
426+
emh_pair_lsmm = new MyEMH_pair(ndepth);
458427

459-
if (accBins.ConfEtallAccBins.value[0] == VARIABLE_WIDTH) {
460-
etall_bin_edges = std::vector<float>(accBins.ConfEtallAccBins.value.begin(), accBins.ConfEtallAccBins.value.end());
461-
etall_bin_edges.erase(etall_bin_edges.begin());
462-
for (const auto& edge : etall_bin_edges) {
463-
LOGF(info, "VARIABLE_WIDTH: etall_bin_edges = %f", edge);
428+
if (accBins.ConfMllAccBins.value[0] == VARIABLE_WIDTH) {
429+
mll_bin_edges = std::vector<float>(accBins.ConfMllAccBins.value.begin(), accBins.ConfMllAccBins.value.end());
430+
mll_bin_edges.erase(mll_bin_edges.begin());
431+
for (const auto& edge : mll_bin_edges) {
432+
LOGF(info, "VARIABLE_WIDTH: mll_bin_edges = %f", edge);
433+
}
434+
} else {
435+
int nbins = static_cast<int>(accBins.ConfMllAccBins.value[0]);
436+
float xmin = static_cast<float>(accBins.ConfMllAccBins.value[1]);
437+
float xmax = static_cast<float>(accBins.ConfMllAccBins.value[2]);
438+
mll_bin_edges.resize(nbins + 1);
439+
for (int i = 0; i < nbins + 1; i++) {
440+
mll_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
441+
LOGF(info, "FIXED_WIDTH: mll_bin_edges[%d] = %f", i, mll_bin_edges[i]);
442+
}
464443
}
465-
} else {
466-
int nbins = static_cast<int>(accBins.ConfEtallAccBins.value[0]);
467-
float xmin = static_cast<float>(accBins.ConfEtallAccBins.value[1]);
468-
float xmax = static_cast<float>(accBins.ConfEtallAccBins.value[2]);
469-
etall_bin_edges.resize(nbins + 1);
470-
for (int i = 0; i < nbins + 1; i++) {
471-
etall_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
472-
LOGF(info, "FIXED_WIDTH: etall_bin_edges[%d] = %f", i, etall_bin_edges[i]);
444+
445+
if (accBins.ConfPtllAccBins.value[0] == VARIABLE_WIDTH) {
446+
ptll_bin_edges = std::vector<float>(accBins.ConfPtllAccBins.value.begin(), accBins.ConfPtllAccBins.value.end());
447+
ptll_bin_edges.erase(ptll_bin_edges.begin());
448+
for (const auto& edge : ptll_bin_edges) {
449+
LOGF(info, "VARIABLE_WIDTH: ptll_bin_edges = %f", edge);
450+
}
451+
} else {
452+
int nbins = static_cast<int>(accBins.ConfPtllAccBins.value[0]);
453+
float xmin = static_cast<float>(accBins.ConfPtllAccBins.value[1]);
454+
float xmax = static_cast<float>(accBins.ConfPtllAccBins.value[2]);
455+
ptll_bin_edges.resize(nbins + 1);
456+
for (int i = 0; i < nbins + 1; i++) {
457+
ptll_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
458+
LOGF(info, "FIXED_WIDTH: ptll_bin_edges[%d] = %f", i, ptll_bin_edges[i]);
459+
}
473460
}
474-
}
475461

476-
if (accBins.ConfPhillAccBins.value[0] == VARIABLE_WIDTH) {
477-
phill_bin_edges = std::vector<float>(accBins.ConfPhillAccBins.value.begin(), accBins.ConfPhillAccBins.value.end());
478-
phill_bin_edges.erase(phill_bin_edges.begin());
479-
for (const auto& edge : phill_bin_edges) {
480-
LOGF(info, "VARIABLE_WIDTH: phill_bin_edges = %f", edge);
462+
if (accBins.ConfEtallAccBins.value[0] == VARIABLE_WIDTH) {
463+
etall_bin_edges = std::vector<float>(accBins.ConfEtallAccBins.value.begin(), accBins.ConfEtallAccBins.value.end());
464+
etall_bin_edges.erase(etall_bin_edges.begin());
465+
for (const auto& edge : etall_bin_edges) {
466+
LOGF(info, "VARIABLE_WIDTH: etall_bin_edges = %f", edge);
467+
}
468+
} else {
469+
int nbins = static_cast<int>(accBins.ConfEtallAccBins.value[0]);
470+
float xmin = static_cast<float>(accBins.ConfEtallAccBins.value[1]);
471+
float xmax = static_cast<float>(accBins.ConfEtallAccBins.value[2]);
472+
etall_bin_edges.resize(nbins + 1);
473+
for (int i = 0; i < nbins + 1; i++) {
474+
etall_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
475+
LOGF(info, "FIXED_WIDTH: etall_bin_edges[%d] = %f", i, etall_bin_edges[i]);
476+
}
481477
}
482-
} else {
483-
int nbins = static_cast<int>(accBins.ConfPhillAccBins.value[0]);
484-
float xmin = static_cast<float>(accBins.ConfPhillAccBins.value[1]);
485-
float xmax = static_cast<float>(accBins.ConfPhillAccBins.value[2]);
486-
phill_bin_edges.resize(nbins + 1);
487-
for (int i = 0; i < nbins + 1; i++) {
488-
phill_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
489-
LOGF(info, "FIXED_WIDTH: phill_bin_edges[%d] = %f", i, phill_bin_edges[i]);
478+
479+
if (accBins.ConfPhillAccBins.value[0] == VARIABLE_WIDTH) {
480+
phill_bin_edges = std::vector<float>(accBins.ConfPhillAccBins.value.begin(), accBins.ConfPhillAccBins.value.end());
481+
phill_bin_edges.erase(phill_bin_edges.begin());
482+
for (const auto& edge : phill_bin_edges) {
483+
LOGF(info, "VARIABLE_WIDTH: phill_bin_edges = %f", edge);
484+
}
485+
} else {
486+
int nbins = static_cast<int>(accBins.ConfPhillAccBins.value[0]);
487+
float xmin = static_cast<float>(accBins.ConfPhillAccBins.value[1]);
488+
float xmax = static_cast<float>(accBins.ConfPhillAccBins.value[2]);
489+
phill_bin_edges.resize(nbins + 1);
490+
for (int i = 0; i < nbins + 1; i++) {
491+
phill_bin_edges[i] = (xmax - xmin) / (nbins)*i + xmin;
492+
LOGF(info, "FIXED_WIDTH: phill_bin_edges[%d] = %f", i, phill_bin_edges[i]);
493+
}
490494
}
491-
}
492495

493-
int nM = mll_bin_edges.size();
494-
int nPt = ptll_bin_edges.size();
495-
int nEta = etall_bin_edges.size();
496-
int nPhi = phill_bin_edges.size();
497-
498-
// emhs_pair_uls.resize(nM);
499-
// emhs_pair_lspp.resize(nM);
500-
// emhs_pair_lsmm.resize(nM);
501-
// for (int im = 0;im<nM;im++) {
502-
// emhs_pair_uls[im].resize(nPt);
503-
// emhs_pair_lspp[im].resize(nPt);
504-
// emhs_pair_lsmm[im].resize(nPt);
505-
// for (int ipt = 0;ipt<nPt;ipt++) {
506-
// emhs_pair_uls[im][ipt].resize(nEta);
507-
// emhs_pair_lspp[im][ipt].resize(nEta);
508-
// emhs_pair_lsmm[im][ipt].resize(nEta);
509-
// for (int ieta = 0;ieta<nEta;ieta++) {
510-
// emhs_pair_uls[im][ipt][ieta].resize(nPhi);
511-
// emhs_pair_lspp[im][ipt][ieta].resize(nPhi);
512-
// emhs_pair_lsmm[im][ipt][ieta].resize(nPhi);
513-
// for (int iphi = 0;iphi<nPhi;iphi++) {
514-
// emhs_pair_uls[im][ipt][ieta][iphi] = new MyEMH_pair(ndepth);
515-
// emhs_pair_lspp[im][ipt][ieta][iphi] = new MyEMH_pair(ndepth);
516-
// emhs_pair_lsmm[im][ipt][ieta][iphi] = new MyEMH_pair(ndepth);
517-
// }
518-
// }
519-
// }
520-
// }
521-
522-
LOGF(info, "nM = %d, nPt = %d, nEta = %d, nPhi = %d", nM, nPt, nEta, nPhi);
523-
524-
std::random_device seed_gen;
525-
engine = std::mt19937(seed_gen());
496+
std::random_device seed_gen;
497+
engine = std::mt19937(seed_gen());
498+
}
526499

527500
DefineEMEventCut();
528501
addhistograms();
@@ -637,30 +610,6 @@ struct Dilepton {
637610
delete emh_pair_lsmm;
638611
emh_pair_lsmm = 0x0;
639612

640-
// for (int im = 0;im<emhs_pair_uls.size();im++) {
641-
// for (int ipt = 0;ipt<emhs_pair_uls[im].size();ipt++) {
642-
// for (int ieta = 0;ieta<emhs_pair_uls[im][ipt].size();ieta++) {
643-
// for (int iphi = 0;iphi<emhs_pair_uls[im][ipt][ieta].size();iphi++) {
644-
// delete emhs_pair_uls[im][ipt][ieta][iphi];
645-
// emhs_pair_uls[im][ipt][ieta][iphi] = 0x0;
646-
647-
// delete emhs_pair_lspp[im][ipt][ieta][iphi];
648-
// emhs_pair_lspp[im][ipt][ieta][iphi] = 0x0;
649-
650-
// delete emhs_pair_lsmm[im][ipt][ieta][iphi];
651-
// emhs_pair_lsmm[im][ipt][ieta][iphi] = 0x0;
652-
// }
653-
// }
654-
// }
655-
// }
656-
657-
// emhs_pair_uls.clear();
658-
// emhs_pair_uls.shrink_to_fit();
659-
// emhs_pair_lspp.clear();
660-
// emhs_pair_lspp.shrink_to_fit();
661-
// emhs_pair_lsmm.clear();
662-
// emhs_pair_lsmm.shrink_to_fit();
663-
664613
used_trackIds_per_col.clear();
665614
used_trackIds_per_col.shrink_to_fit();
666615
map_mixed_eventId_to_globalBC.clear();
@@ -1348,9 +1297,6 @@ struct Dilepton {
13481297
MyEMH_pair* emh_pair_lspp = nullptr;
13491298
MyEMH_pair* emh_pair_lsmm = nullptr;
13501299

1351-
// std::vector<std::vector<std::vector<std::vector<MyEMH_pair*>>>> emhs_pair_uls; // 4D{m, pt, eta, phi}
1352-
// std::vector<std::vector<std::vector<std::vector<MyEMH_pair*>>>> emhs_pair_lspp; // 4D{m, pt, eta, phi}
1353-
// std::vector<std::vector<std::vector<std::vector<MyEMH_pair*>>>> emhs_pair_lsmm; // 4D{m, pt, eta, phi}
13541300
std::map<std::pair<int, int>, uint64_t> map_mixed_eventId_to_globalBC;
13551301

13561302
std::vector<int> used_trackIds_per_col;
@@ -1730,19 +1676,6 @@ struct Dilepton {
17301676
emh_pair_uls->AddCollisionIdAtLast(key_bin, key_df_collision);
17311677
emh_pair_lspp->AddCollisionIdAtLast(key_bin, key_df_collision);
17321678
emh_pair_lsmm->AddCollisionIdAtLast(key_bin, key_df_collision);
1733-
1734-
// for (int im = 0;im<emhs_pair_uls.size();im++) {
1735-
// for (int ipt = 0;ipt<emhs_pair_uls[im].size();ipt++) {
1736-
// for (int ieta = 0;ieta<emhs_pair_uls[im][ipt].size();ieta++) {
1737-
// for (int iphi = 0;iphi<emhs_pair_uls[im][ipt][ieta].size();iphi++) {
1738-
// emhs_pair_uls[im][ipt][ieta][iphi]->AddCollisionIdAtLast(key_bin, key_df_collision);
1739-
// emhs_pair_lspp[im][ipt][ieta][iphi]->AddCollisionIdAtLast(key_bin, key_df_collision);
1740-
// emhs_pair_lsmm[im][ipt][ieta][iphi]->AddCollisionIdAtLast(key_bin, key_df_collision);
1741-
// }
1742-
// }
1743-
// }
1744-
// }
1745-
17461679
} // end of if pair exist
17471680

17481681
} // end of collision loop

‎PWGEM/Dilepton/TableProducer/createEMEventDilepton.cxx‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -181,8 +181,8 @@ struct CreateEMEventDilepton {
181181
q3xft0m = collision.qvecFT0MReVec()[1], q3xft0a = collision.qvecFT0AReVec()[1], q3xft0c = collision.qvecFT0CReVec()[1], q3xfv0a = collision.qvecFV0AReVec()[1], q3xbpos = collision.qvecBPosReVec()[1], q3xbneg = collision.qvecBNegReVec()[1], q3xbtot = collision.qvecBTotReVec()[1];
182182
q3yft0m = collision.qvecFT0MImVec()[1], q3yft0a = collision.qvecFT0AImVec()[1], q3yft0c = collision.qvecFT0CImVec()[1], q3yfv0a = collision.qvecFV0AImVec()[1], q3ybpos = collision.qvecBPosImVec()[1], q3ybneg = collision.qvecBNegImVec()[1], q3ybtot = collision.qvecBTotImVec()[1];
183183
} else if (collision.qvecFT0CReVec().size() >= 1) { // harmonics 2
184-
q2xft0m = collision.qvecFT0MReVec()[0], q2xft0a = collision.qvecFT0AReVec()[0], q2xft0c = collision.qvecFT0CReVec()[0], q2xbpos = collision.qvecBPosReVec()[0], q2xbneg = collision.qvecBNegReVec()[0], q2xbtot = collision.qvecBTotReVec()[0];
185-
q2yft0m = collision.qvecFT0MImVec()[0], q2yft0a = collision.qvecFT0AImVec()[0], q2yft0c = collision.qvecFT0CImVec()[0], q2ybpos = collision.qvecBPosImVec()[0], q2ybneg = collision.qvecBNegImVec()[0], q2ybtot = collision.qvecBTotImVec()[0];
184+
q2xft0m = collision.qvecFT0MReVec()[0], q2xft0a = collision.qvecFT0AReVec()[0], q2xft0c = collision.qvecFT0CReVec()[0], q2xfv0a = collision.qvecFV0AReVec()[0], q2xbpos = collision.qvecBPosReVec()[0], q2xbneg = collision.qvecBNegReVec()[0], q2xbtot = collision.qvecBTotReVec()[0];
185+
q2yft0m = collision.qvecFT0MImVec()[0], q2yft0a = collision.qvecFT0AImVec()[0], q2yft0c = collision.qvecFT0CImVec()[0], q2yfv0a = collision.qvecFV0AImVec()[0], q2ybpos = collision.qvecBPosImVec()[0], q2ybneg = collision.qvecBNegImVec()[0], q2ybtot = collision.qvecBTotImVec()[0];
186186
}
187187
event_qvec(
188188
q2xft0m, q2yft0m, q2xft0a, q2yft0a, q2xft0c, q2yft0c, q2xfv0a, q2yfv0a, q2xbpos, q2ybpos, q2xbneg, q2ybneg, q2xbtot, q2ybtot,

0 commit comments

Comments
 (0)