Skip to content

Commit 224eb87

Browse files
adding weights for pT shape Study
1 parent 2c61921 commit 224eb87

File tree

2 files changed

+98
-62
lines changed

2 files changed

+98
-62
lines changed

PWGHF/D2H/Tasks/taskDstarToD0Pi.cxx

Lines changed: 97 additions & 61 deletions
Original file line numberDiff line numberDiff line change
@@ -61,15 +61,16 @@ struct HfTaskDstarToD0Pi {
6161
Configurable<bool> isCentStudy{"isCentStudy", true, "Flag to select centrality study"};
6262
Configurable<bool> qaEnabled{"qaEnabled", true, "Flag to enable QA histograms"};
6363
Configurable<bool> studyD0ToPiKPi0{"studyD0ToPiKPi0", false, "Flag to study D*->D0(piKpi0)pi channel"};
64+
Configurable<bool> ptShapeStudy{"ptShapeStudy", true, "Flag to enable pT shape study"};
6465

6566
// CCDB configuration
66-
Configurable<bool> useWeight{"useWeight", true, "Flag to use weights from CCDB"};
67-
Configurable<int> nWeights{"nWeights", 6, "Number of weights to be used from CCDB"};
67+
// Configurable<bool> useWeight{"useWeight", true, "Flag to use weights from CCDB"};
68+
// Configurable<int> nWeights{"nWeights", 6, "Number of weights to be used from CCDB"};
6869
Configurable<std::string> ccdbUrl{"ccdbUrl", "http://alice-ccdb.cern.ch", "url of the ccdb repository"};
69-
Configurable<std::string> ccdbPathForWeight{"ccdbPathForWeight", "Users/d/desharma/weights", "CCDB path for PVContrib weights"};
70-
Configurable<int64_t> timestampCCDB{"timestampCCDB", -1, "CCDB timestamp for the weights"};
71-
Configurable<std::string> weightFileName{"weightFileName", "Weights.root", "Name of the weight file to be used for PVContrib"};
72-
Configurable<std::vector<double>> centRangesForWeights{"centRangesForWeights", {0.0, 5.0, 10.0, 30.0, 50.0, 70.0, 100.0}, "Centrality ranges for weights. Size of ranges should be equal to nWeights + 1"};
70+
Configurable<std::string> ccdbPathForWeight{"ccdbPathForWeight", "", "CCDB path for pt shape weights"};
71+
Configurable<int64_t> timestampCCDB{"timestampCCDB", -1, "CCDB timestamp for pt shape weights"};
72+
Configurable<std::string> weightFileName{"weightFileName", "Weights.root", "Name of the weight file to be used for pt shape study"};
73+
// Configurable<std::vector<double>> centRangesForWeights{"centRangesForWeights", {0.0, 5.0, 10.0, 30.0, 50.0, 70.0, 100.0}, "Centrality ranges for weights. Size of ranges should be equal to nWeights + 1"};
7374

7475
Configurable<double> yCandDstarRecoMax{"yCandDstarRecoMax", 0.8, "max. candidate Dstar rapidity"};
7576
Configurable<double> yCandDstarGenMax{"yCandDstarGenMax", 0.5, "max. rapidity of Generator level Particle"};
@@ -80,7 +81,19 @@ struct HfTaskDstarToD0Pi {
8081

8182
o2::ccdb::CcdbApi ccdbApi;
8283
SliceCache cache;
83-
std::vector<TH2F*> hWeights;
84+
std::vector<TH1D*> hWeights;
85+
int const nWeights = 2; // prompt and non-prompt weights
86+
std::vector<std::string> const weightHistNames = {"promptWeightVsPt", "nonPromptWeightVsPt"};
87+
enum weightType {
88+
Prompt = 0,
89+
NonPrompt = 1
90+
};
91+
92+
// follwing commented axes can be removed if not needed
93+
// std::vector<AxisSpec> axesPtVsCentVsBDTVsPvContribVsPtB;
94+
// std::vector<AxisSpec> axesPtVsCentVsPvContribVsPtB;
95+
// std::vector<AxisSpec> axesPtVsBDTVsPtB;
96+
8497
std::vector<AxisSpec> axesPtVsCentVsBDTVsPvContrib;
8598
std::vector<AxisSpec> axesPtVsCentVsPvContrib;
8699
std::vector<AxisSpec> axesPtVsBDT;
@@ -114,6 +127,8 @@ struct HfTaskDstarToD0Pi {
114127
ConfigurableAxis binningBkgBDTScore{"binningBkgBDTScore", {100, 0.0f, 1.0f}, "Bins for background BDT Score"};
115128
ConfigurableAxis binningSigBDTScore{"binningSigBDTScore", {100, 0.0f, 1.0f}, "Bins for Signal (Prompts + Non Prompt) BDT Score"};
116129
ConfigurableAxis binningPvContrib{"binningPvContrib", {100, 0.0f, 300.0f}, "Bins for PVContrib"};
130+
// can be removed if not needed
131+
// ConfigurableAxis binningPtFine{"binningPtFine", {100, 0.0f, 100.0f}, "Bins for pT shape study"};
117132

118133
HistogramRegistry registry{"registry", {}};
119134

@@ -136,10 +151,25 @@ struct HfTaskDstarToD0Pi {
136151
AxisSpec axisPvContrib = {binningPvContrib, "PV Contribution"};
137152
AxisSpec const axisPt = {vecPtBins, "#it{p}_{T} (GeV/#it{c})"};
138153

154+
// can be removed if not needed
155+
// AxisSpec axisPtFine = {vecPtBins, "#it{p}_{T} (GeV/#it{c})"};
156+
// if (ptShapeStudy){
157+
// axisPtFine = {binningPtFine, "pT (GeV/c)"};
158+
// }
159+
// // std::cout<< "fine pt bins for pt shape study: binningPtFine"<< binningPtFine->data() << std::endl;
160+
// // std::cout<< "fine pt bins for pt shape study: axisPtFine"<< axisPtFine.binEdges[2] << std::endl;
161+
139162
axesPtVsCentVsBDTVsPvContrib = {axisPt, axisCentrality, axisBDTScoreBackground, axisBDTScorePrompt, axisBDTScoreNonPrompt, axisPvContrib};
140163
axesPtVsCentVsPvContrib = {axisPt, axisCentrality, axisPvContrib};
141164
axesPtVsBDT = {axisPt, axisBDTScoreBackground, axisBDTScorePrompt, axisBDTScoreNonPrompt};
142165

166+
// can be removed if
167+
// if (ptShapeStudy){
168+
// axesPtVsCentVsBDTVsPvContribVsPtB = {axisPt, axisCentrality, axisBDTScoreBackground, axisBDTScorePrompt, axisBDTScoreNonPrompt, axisPvContrib, axisPtFine};
169+
// axesPtVsCentVsPvContribVsPtB = {axisPt, axisCentrality, axisPvContrib, axisPtFine};
170+
// axesPtVsBDTVsPtB = {axisPt, axisBDTScoreBackground, axisBDTScorePrompt, axisBDTScoreNonPrompt, axisPtFine};
171+
// }
172+
143173
if (qaEnabled) {
144174
// only QA
145175
registry.add("QA/hPtDstar", "Dstar Candidates; Dstar candidate #it{p}_{T} (GeV/#it{c}); entries", {HistType::kTH1F, {{360, 0., 36.}}});
@@ -289,8 +319,8 @@ struct HfTaskDstarToD0Pi {
289319
}
290320
}
291321

292-
// if weights to be applied
293-
if (useWeight) {
322+
// if weights to be applied for pt shape study
323+
if (ptShapeStudy) {
294324
ccdbApi.init(ccdbUrl);
295325
std::map<std::string, std::string> const metadata;
296326
// Retrieve the file from CCDB
@@ -300,28 +330,30 @@ struct HfTaskDstarToD0Pi {
300330
return;
301331
}
302332

303-
if (isCentStudy) {
304-
// Open the ROOT file
305-
TFile* weightFile = TFile::Open(weightFileName.value.c_str(), "READ");
306-
if ((weightFile != nullptr) && !weightFile->IsZombie()) {
307-
// Ensure hWeights is properly sized
308-
hWeights.resize(nWeights);
309-
for (int ithWeight = 0; ithWeight < nWeights; ++ithWeight) {
310-
std::string const histName = "hMult" + std::to_string(ithWeight + 1) + "_Weight";
311-
hWeights[ithWeight] = dynamic_cast<TH2F*>(weightFile->Get(histName.c_str()));
312-
if (hWeights[ithWeight] == nullptr) {
313-
LOGF(fatal, "Histogram %s not found in weight file!", histName.c_str());
314-
return;
315-
}
316-
hWeights[ithWeight]->SetDirectory(nullptr);
317-
hWeights[ithWeight]->SetName(("hWeight" + std::to_string(ithWeight + 1)).c_str());
318-
}
319-
weightFile->Close();
320-
delete weightFile;
321-
} else {
322-
LOGF(fatal, "Failed to open weight file from CCDB: %s", weightFileName.value.c_str());
333+
// Open the ROOT file to intialise weight hists
334+
TFile* weightFile = TFile::Open(weightFileName.value.c_str(), "READ");
335+
if ((weightFile != nullptr) && !weightFile->IsZombie()) {
336+
// Ensure hWeights is properly sized
337+
hWeights.resize(nWeights); // prompt and non-prompt
338+
339+
hWeights[weightType::Prompt] = dynamic_cast<TH1D*>(weightFile->Get(weightHistNames[weightType::Prompt].c_str()));
340+
hWeights[weightType::NonPrompt] = dynamic_cast<TH1D*>(weightFile->Get(weightHistNames[weightType::NonPrompt].c_str()));
341+
if (hWeights[weightType::Prompt] == nullptr) {
342+
LOGF(fatal, "Histogram %s not found in weight file!", weightHistNames[weightType::Prompt].c_str());
323343
return;
324344
}
345+
if (hWeights[weightType::NonPrompt] == nullptr) {
346+
LOGF(fatal, "Histogram %s not found in weight file!", weightHistNames[weightType::NonPrompt].c_str());
347+
return;
348+
}
349+
hWeights[weightType::Prompt]->SetDirectory(nullptr);
350+
hWeights[weightType::NonPrompt]->SetDirectory(nullptr);
351+
352+
weightFile->Close();
353+
delete weightFile;
354+
} else {
355+
LOGF(fatal, "Failed to open weight file from CCDB: %s", weightFileName.value.c_str());
356+
return;
325357
}
326358
}
327359
}
@@ -473,9 +505,10 @@ struct HfTaskDstarToD0Pi {
473505
auto prong0 = candDstarMcRec.template prong0_as<aod::TracksWMc>();
474506
auto indexMother = RecoDecay::getMother(rowsMcPartilces, prong0.template mcParticle_as<CandDstarMcGen>(), o2::constants::physics::Pdg::kDStar, true, &signDstar, 2);
475507
auto particleMother = rowsMcPartilces.rawIteratorAt(indexMother); // What is difference between rawIterator() or iteratorAt() methods?
508+
auto ptMother = particleMother.pt();
476509
if (qaEnabled) {
477-
registry.fill(HIST("QA/hPtSkimDstarGenSig"), particleMother.pt()); // generator level pt
478-
registry.fill(HIST("QA/hPtVsCentSkimDstarGenSig"), particleMother.pt(), centrality);
510+
registry.fill(HIST("QA/hPtSkimDstarGenSig"), ptMother); // generator level pt
511+
registry.fill(HIST("QA/hPtVsCentSkimDstarGenSig"), ptMother, centrality);
479512
registry.fill(HIST("QA/hPtVsYSkimDstarRecSig"), ptDstarRecSig, yDstarRecSig); // Skimed at level of trackIndexSkimCreator
480513
if (candDstarMcRec.isRecoTopol()) { // if Topological selection are passed
481514
registry.fill(HIST("QA/hPtVsYRecoTopolDstarRecSig"), ptDstarRecSig, yDstarRecSig);
@@ -487,18 +520,17 @@ struct HfTaskDstarToD0Pi {
487520

488521
if (candDstarMcRec.isSelDstarToD0Pi()) { // if all selection passed
489522
float weightValue = 1.0;
490-
if (useWeight && (hWeights.empty() || hWeights[0] == nullptr)) {
523+
std::vector<float> ptShapeWeightValues(nWeights, 1.0); // Assuming two weights: one for prompt and one for non-prompt
524+
if (ptShapeStudy && (hWeights.empty() || hWeights[0] == nullptr)) {
491525
LOGF(fatal, "Weight histograms are not initialized or empty. Check CCDB path or weight file.");
492526
return;
493527
}
494-
if (useWeight && isCentStudy) {
495-
for (int ithWeight = 0; ithWeight < nWeights; ++ithWeight) {
496-
if (centrality > centRangesForWeights.value[ithWeight] && centrality <= centRangesForWeights.value[ithWeight + 1]) {
497-
weightValue = hWeights[ithWeight]->GetBinContent(hWeights[ithWeight]->FindBin(nPVContributors));
498-
break;
499-
}
500-
}
528+
529+
if (ptShapeStudy) {
530+
ptShapeWeightValues[weightType::Prompt] = hWeights[weightType::Prompt]->GetBinContent(hWeights[weightType::Prompt]->FindBin(ptDstarRecSig));
531+
ptShapeWeightValues[weightType::NonPrompt] = hWeights[weightType::NonPrompt]->GetBinContent(hWeights[weightType::NonPrompt]->FindBin(ptMother));
501532
}
533+
502534
if (qaEnabled) {
503535
registry.fill(HIST("QA/hPtFullRecoDstarRecSig"), ptDstarRecSig);
504536
}
@@ -508,33 +540,33 @@ struct HfTaskDstarToD0Pi {
508540
auto bdtScore = candDstarMcRec.mlProbDstarToD0Pi();
509541
registry.fill(HIST("Efficiency/hPtVsCentVsBDTScoreVsPvContribRecSig"), ptDstarRecSig, centrality, bdtScore[0], bdtScore[1], bdtScore[2], nPVContributors, weightValue);
510542
if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::Prompt) {
511-
registry.fill(HIST("Efficiency/hPtPromptVsCentVsBDTScorePvContribRecSig"), ptDstarRecSig, centrality, bdtScore[0], bdtScore[1], bdtScore[2], nPVContributors, weightValue);
543+
registry.fill(HIST("Efficiency/hPtPromptVsCentVsBDTScorePvContribRecSig"), ptDstarRecSig, centrality, bdtScore[0], bdtScore[1], bdtScore[2], nPVContributors, ptShapeWeightValues[weightType::Prompt]);
512544
} else if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::NonPrompt) {
513-
registry.fill(HIST("Efficiency/hPtNonPrompRectVsCentVsBDTScorePvContribRecSig"), ptDstarRecSig, centrality, bdtScore[0], bdtScore[1], bdtScore[2], nPVContributors, weightValue);
545+
registry.fill(HIST("Efficiency/hPtNonPrompRectVsCentVsBDTScorePvContribRecSig"), ptDstarRecSig, centrality, bdtScore[0], bdtScore[1], bdtScore[2], nPVContributors, ptShapeWeightValues[weightType::NonPrompt]);
514546
}
515547
} else {
516548
auto bdtScore = candDstarMcRec.mlProbDstarToD0Pi();
517549
registry.fill(HIST("Efficiency/hPtVsBDTScoreRecSig"), ptDstarRecSig, bdtScore[0], bdtScore[1], bdtScore[2], weightValue);
518550
if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::Prompt) {
519-
registry.fill(HIST("Efficiency/hPtPromptVsBDTScoreRecSig"), ptDstarRecSig, bdtScore[0], bdtScore[1], bdtScore[2], weightValue);
551+
registry.fill(HIST("Efficiency/hPtPromptVsBDTScoreRecSig"), ptDstarRecSig, bdtScore[0], bdtScore[1], bdtScore[2], ptShapeWeightValues[weightType::Prompt]);
520552
} else if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::NonPrompt) {
521-
registry.fill(HIST("Efficiency/hPtNonPromptVsBDTScoreRecSig"), ptDstarRecSig, bdtScore[0], bdtScore[1], bdtScore[2], weightValue);
553+
registry.fill(HIST("Efficiency/hPtNonPromptVsBDTScoreRecSig"), ptDstarRecSig, bdtScore[0], bdtScore[1], bdtScore[2], ptShapeWeightValues[weightType::NonPrompt]);
522554
}
523555
}
524556
} else { // All efficiency histograms at reconstruction level w/o ml
525557
if (isCentStudy) {
526558
registry.fill(HIST("Efficiency/hPtVsCentVsPvContribRecSig"), ptDstarRecSig, centrality, nPVContributors, weightValue);
527559
if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::Prompt) {
528-
registry.fill(HIST("Efficiency/hPtPromptVsCentVsPvContribRecSig"), ptDstarRecSig, centrality, nPVContributors, weightValue);
560+
registry.fill(HIST("Efficiency/hPtPromptVsCentVsPvContribRecSig"), ptDstarRecSig, centrality, nPVContributors, ptShapeWeightValues[weightType::Prompt]);
529561
} else if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::NonPrompt) {
530-
registry.fill(HIST("Efficiency/hPtNonPromptVsCentVsPvContribRecSig"), ptDstarRecSig, centrality, nPVContributors, weightValue);
562+
registry.fill(HIST("Efficiency/hPtNonPromptVsCentVsPvContribRecSig"), ptDstarRecSig, centrality, nPVContributors, ptShapeWeightValues[weightType::NonPrompt]);
531563
}
532564
} else {
533565
registry.fill(HIST("Efficiency/hPtVsPvContribRecSig"), ptDstarRecSig, nPVContributors, weightValue);
534566
if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::Prompt) {
535-
registry.fill(HIST("Efficiency/hPtPromptVsPvContribRecSig"), ptDstarRecSig, nPVContributors, weightValue);
567+
registry.fill(HIST("Efficiency/hPtPromptVsPvContribRecSig"), ptDstarRecSig, nPVContributors, ptShapeWeightValues[weightType::Prompt]);
536568
} else if (candDstarMcRec.originMcRec() == RecoDecay::OriginType::NonPrompt) {
537-
registry.fill(HIST("Efficiency/hPtNonPromptVsPvContribRecSig"), ptDstarRecSig, nPVContributors, weightValue);
569+
registry.fill(HIST("Efficiency/hPtNonPromptVsPvContribRecSig"), ptDstarRecSig, nPVContributors, ptShapeWeightValues[weightType::NonPrompt]);
538570
}
539571
}
540572
}
@@ -631,6 +663,12 @@ struct HfTaskDstarToD0Pi {
631663
for (auto const& mcParticle : rowsMcPartilces) {
632664
if (std::abs(mcParticle.flagMcMatchGen()) == hf_decay::hf_cand_dstar::DecayChannelMain::DstarToPiKPi) { // MC Matching is successful at Generator Level
633665
auto ptGen = mcParticle.pt();
666+
667+
// mother information
668+
auto idxBhadMotherPart = mcParticle.idxBhadMotherPart();
669+
auto bMother = rowsMcPartilces.rawIteratorAt(idxBhadMotherPart);
670+
auto ptBMother = bMother.pt();
671+
634672
auto yGen = RecoDecay::y(mcParticle.pVector(), o2::constants::physics::MassDStar);
635673
if (yCandDstarGenMax >= 0. && std::abs(yGen) > yCandDstarGenMax) {
636674
continue;
@@ -668,17 +706,15 @@ struct HfTaskDstarToD0Pi {
668706
}
669707

670708
float weightValue = 1.0;
671-
if (useWeight && (hWeights.empty() || hWeights[0] == nullptr)) {
709+
std::vector<float> ptShapeWeightValues(nWeights, 1.0); // Assuming two weights: one for prompt and one for non-prompt
710+
if (ptShapeStudy && (hWeights.empty() || hWeights[0] == nullptr)) {
672711
LOGF(fatal, "Weight histograms are not initialized or empty. Check CCDB path or weight file.");
673712
return;
674713
}
675-
if (useWeight && isCentStudy) {
676-
for (int ithWeight = 0; ithWeight < nWeights; ++ithWeight) {
677-
if (centFT0MGen > centRangesForWeights.value[ithWeight] && centFT0MGen <= centRangesForWeights.value[ithWeight + 1]) {
678-
weightValue = hWeights[ithWeight]->GetBinContent(hWeights[ithWeight]->FindBin(centFT0MGen, pvContributors));
679-
break;
680-
}
681-
}
714+
715+
if (ptShapeStudy) {
716+
ptShapeWeightValues[weightType::Prompt] = hWeights[weightType::Prompt]->GetBinContent(hWeights[weightType::Prompt]->FindBin(ptGen));
717+
ptShapeWeightValues[weightType::NonPrompt] = hWeights[weightType::NonPrompt]->GetBinContent(hWeights[weightType::NonPrompt]->FindBin(ptBMother));
682718
}
683719

684720
registry.fill(HIST("Efficiency/hPtVsYDstarGen"), ptGen, yGen, weightValue);
@@ -693,20 +729,20 @@ struct HfTaskDstarToD0Pi {
693729

694730
// Prompt
695731
if (mcParticle.originMcGen() == RecoDecay::OriginType::Prompt) {
696-
registry.fill(HIST("Efficiency/hPtVsYPromptDstarGen"), ptGen, yGen, weightValue);
732+
registry.fill(HIST("Efficiency/hPtVsYPromptDstarGen"), ptGen, yGen, ptShapeWeightValues[weightType::Prompt]);
697733
if (isCentStudy) {
698-
registry.fill(HIST("Efficiency/hPtPromptVsCentVsPvContribGen"), ptGen, centFT0MGen, pvContributors, weightValue);
734+
registry.fill(HIST("Efficiency/hPtPromptVsCentVsPvContribGen"), ptGen, centFT0MGen, pvContributors, ptShapeWeightValues[weightType::Prompt]);
699735
} else {
700-
registry.fill(HIST("Efficiency/hPtPromptVsGen"), ptGen, weightValue);
736+
registry.fill(HIST("Efficiency/hPtPromptVsGen"), ptGen, ptShapeWeightValues[weightType::Prompt]);
701737
}
702738
// Non-Prompt
703739
} else if (mcParticle.originMcGen() == RecoDecay::OriginType::NonPrompt) {
704-
registry.fill(HIST("Efficiency/hPtVsYNonPromptDstarGen"), ptGen, yGen, weightValue);
740+
registry.fill(HIST("Efficiency/hPtVsYNonPromptDstarGen"), ptGen, yGen, ptShapeWeightValues[weightType::NonPrompt]);
705741
if (isCentStudy) {
706-
registry.fill(HIST("Efficiency/hPtNonPromptVsCentVsPvContribGen"), ptGen, centFT0MGen, pvContributors, weightValue);
742+
registry.fill(HIST("Efficiency/hPtNonPromptVsCentVsPvContribGen"), ptGen, centFT0MGen, pvContributors, ptShapeWeightValues[weightType::NonPrompt]);
707743

708744
} else {
709-
registry.fill(HIST("Efficiency/hPtNonPromptVsGen"), ptGen, weightValue);
745+
registry.fill(HIST("Efficiency/hPtNonPromptVsGen"), ptGen, ptShapeWeightValues[weightType::NonPrompt]);
710746
}
711747
}
712748
}

PWGHF/TableProducer/candidateCreatorDstar.cxx

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -720,7 +720,7 @@ struct HfCandidateCreatorDstarExpressions {
720720

721721
// check wether the particle is non-promt (from a B0 hadron)
722722
if (flagDstar != 0) {
723-
auto particleDstar = mcParticles.iteratorAt(indexRecDstar);
723+
auto particleDstar = mcParticles.rawIteratorAt(indexRecDstar);
724724
originDstar = RecoDecay::getCharmHadronOrigin(mcParticles, particleDstar, false, &idxBhadMothers);
725725
}
726726
if (originDstar == RecoDecay::OriginType::NonPrompt) {

0 commit comments

Comments
 (0)