Skip to content

[add] Some gen-matched branches to the new skim format #32

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
janice-cat wants to merge 1 commit into
base: main
Choose a base branch
from
Open

[add] Some gen-matched branches to the new skim format #32

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
8 changes: 8 additions & 0 deletions CommonCode/include/Messenger.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -418,6 +418,11 @@ class DzeroGenTreeMessenger
float Gpt[DZEROGENCOUNTMAX];
float Gy[DZEROGENCOUNTMAX];
int GpdgId[DZEROGENCOUNTMAX];
int Gdau1pdgId[DZEROGENCOUNTMAX];
int Gdau2pdgId[DZEROGENCOUNTMAX];
int Gdau3pdgId[DZEROGENCOUNTMAX];
int Gdau4pdgId[DZEROGENCOUNTMAX];
int GnDa[DZEROGENCOUNTMAX];
int GisSignal[DZEROGENCOUNTMAX];
int GcollisionId[DZEROGENCOUNTMAX];
int GSignalType[DZEROGENCOUNTMAX];
Expand Down Expand Up @@ -762,6 +767,9 @@ class DzeroUPCTreeMessenger
std::vector<bool> *DisSignalCalc;
std::vector<bool> *DisSignalCalcPrompt;
std::vector<bool> *DisSignalCalcFeeddown;
std::vector<bool> *DisSwapCalc;
std::vector<bool> *DisD0CalcLeftSideBand; // D0 decay in the left side band (dominated by KK-state)
std::vector<bool> *DisD0CalcRightSideBand; // D0 decay in the right side band (dominated by pipi-state)
//MC only quantities
int Gsize;
std::vector<float> *Gpt;
Expand Down
26 changes: 26 additions & 0 deletions CommonCode/source/Messenger.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1558,6 +1558,11 @@ bool DzeroGenTreeMessenger::Initialize()
Tree->SetBranchAddress("Gpt", &Gpt);
Tree->SetBranchAddress("Gy", &Gy);
Tree->SetBranchAddress("GpdgId", &GpdgId);
Tree->SetBranchAddress("Gdau1pdgId", &Gdau1pdgId);
Tree->SetBranchAddress("Gdau2pdgId", &Gdau2pdgId);
Tree->SetBranchAddress("Gdau3pdgId", &Gdau3pdgId);
Tree->SetBranchAddress("Gdau4pdgId", &Gdau4pdgId);
Tree->SetBranchAddress("GnDa", &GnDa);
Tree->SetBranchAddress("GisSignal", &GisSignal);
Tree->SetBranchAddress("GcollisionId", &GcollisionId);
Tree->SetBranchAddress("GSignalType", &GSignalType);
Expand Down Expand Up @@ -2753,6 +2758,9 @@ DzeroUPCTreeMessenger::~DzeroUPCTreeMessenger()
delete DisSignalCalc;
delete DisSignalCalcPrompt;
delete DisSignalCalcFeeddown;
delete DisSwapCalc;
delete DisD0CalcLeftSideBand;
delete DisD0CalcRightSideBand;
delete Gpt;
delete Gy;
delete GisSignalCalc;
Expand Down Expand Up @@ -2792,6 +2800,9 @@ bool DzeroUPCTreeMessenger::Initialize(bool Debug)
DisSignalCalc = nullptr;
DisSignalCalcPrompt = nullptr;
DisSignalCalcFeeddown = nullptr;
DisSwapCalc = nullptr;
DisD0CalcLeftSideBand = nullptr;
DisD0CalcRightSideBand = nullptr;
Gpt = nullptr;
Gy = nullptr;
GisSignalCalc = nullptr;
Expand Down Expand Up @@ -2843,6 +2854,9 @@ bool DzeroUPCTreeMessenger::Initialize(bool Debug)
Tree->SetBranchAddress("DisSignalCalc", &DisSignalCalc);
Tree->SetBranchAddress("DisSignalCalcPrompt", &DisSignalCalcPrompt);
Tree->SetBranchAddress("DisSignalCalcFeeddown", &DisSignalCalcFeeddown);
Tree->SetBranchAddress("DisSwapCalc", &DisSwapCalc);
Tree->SetBranchAddress("DisD0CalcLeftSideBand", &DisD0CalcLeftSideBand);
Tree->SetBranchAddress("DisD0CalcRightSideBand", &DisD0CalcRightSideBand);
Tree->SetBranchAddress("Gsize", &Gsize);
Tree->SetBranchAddress("Gpt", &Gpt);
Tree->SetBranchAddress("Gy", &Gy);
Expand Down Expand Up @@ -2896,6 +2910,9 @@ bool DzeroUPCTreeMessenger::SetBranch(TTree *T)
DisSignalCalc = new std::vector<bool>();
DisSignalCalcPrompt = new std::vector<bool>();
DisSignalCalcFeeddown = new std::vector<bool>();
DisSwapCalc = new std::vector<bool>();
DisD0CalcLeftSideBand = new std::vector<bool>();
DisD0CalcRightSideBand = new std::vector<bool>();

Gpt = new std::vector<float>();
Gy = new std::vector<float>();
Expand Down Expand Up @@ -2951,6 +2968,9 @@ bool DzeroUPCTreeMessenger::SetBranch(TTree *T)
Tree->Branch("DisSignalCalc", &DisSignalCalc);
Tree->Branch("DisSignalCalcPrompt", &DisSignalCalcPrompt);
Tree->Branch("DisSignalCalcFeeddown", &DisSignalCalcFeeddown);
Tree->Branch("DisSwapCalc", &DisSwapCalc);
Tree->Branch("DisD0CalcLeftSideBand", &DisD0CalcLeftSideBand);
Tree->Branch("DisD0CalcRightSideBand", &DisD0CalcRightSideBand);

Tree->Branch("Gsize", &Gsize);
Tree->Branch("Gpt", &Gpt);
Expand Down Expand Up @@ -3011,6 +3031,9 @@ void DzeroUPCTreeMessenger::Clear()
DisSignalCalc->clear();
DisSignalCalcPrompt->clear();
DisSignalCalcFeeddown->clear();
DisSwapCalc->clear();
DisD0CalcLeftSideBand->clear();
DisD0CalcRightSideBand->clear();
Gsize = 0;
Gpt->clear();
Gy->clear();
Expand Down Expand Up @@ -3066,6 +3089,9 @@ void DzeroUPCTreeMessenger::CopyNonTrack(DzeroUPCTreeMessenger &M)
if(DisSignalCalc != nullptr && M.DisSignalCalc != nullptr) *DisSignalCalc = *(M.DisSignalCalc);
if(DisSignalCalcPrompt != nullptr && M.DisSignalCalcPrompt != nullptr) *DisSignalCalcPrompt = *(M.DisSignalCalcPrompt);
if(DisSignalCalcFeeddown != nullptr && M.DisSignalCalcFeeddown != nullptr) *DisSignalCalcFeeddown = *(M.DisSignalCalcFeeddown);
if(DisSwapCalc != nullptr && M.DisSwapCalc != nullptr) *DisSwapCalc = *(M.DisSwapCalc);
if(DisD0CalcLeftSideBand != nullptr && M.DisD0CalcLeftSideBand != nullptr) *DisD0CalcLeftSideBand = *(M.DisD0CalcLeftSideBand);
if(DisD0CalcRightSideBand != nullptr && M.DisD0CalcRightSideBand != nullptr) *DisD0CalcRightSideBand = *(M.DisD0CalcRightSideBand);
Gsize = M.Gsize;
if(Gpt != nullptr && M.Gpt != nullptr) *Gpt = *(M.Gpt);
if(Gy != nullptr && M.Gy != nullptr) *Gy = *(M.Gy);
Expand Down
74 changes: 72 additions & 2 deletions SampleGeneration/20241027_ForestReducer_DzeroUPC/ReduceForest.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -24,12 +24,15 @@ using namespace std;

#include "include/DmesonSelection.h"

#define D0_MASS 1.8648

bool logical_or_vectBool(std::vector<bool>* vec) {
return std::any_of(vec->begin(), vec->end(), [](bool b) { return b; });
}

int main(int argc, char *argv[]);
double GetMaxEnergyHF(PFTreeMessenger *M, double etaMin, double etaMax);
bool GisSignalwFSR(int GpdgId, int Gdau1pdgId, int Gdau2pdgId, int Gdau3pdgId, int Gdau4pdgId, int GnDa, int GisSignal);

int main(int argc, char *argv[]) {
string VersionString = "V8";
Expand Down Expand Up @@ -137,7 +140,10 @@ int main(int argc, char *argv[]) {
MDzeroUPC.Gpt->push_back(MDzeroGen.Gpt[iDGen]);
MDzeroUPC.Gy->push_back(MDzeroGen.Gy[iDGen]);
bool isSignalGen =
(MDzeroGen.GisSignal[iDGen] == 1 || MDzeroGen.GisSignal[iDGen] == 2) && MDzeroGen.Gpt[iDGen] > 0.;
MDzeroGen.Gpt[iDGen] > 0. && \
GisSignalwFSR(MDzeroGen.GpdgId[iDGen], MDzeroGen.Gdau1pdgId[iDGen], MDzeroGen.Gdau2pdgId[iDGen],
MDzeroGen.Gdau3pdgId[iDGen], MDzeroGen.Gdau4pdgId[iDGen], MDzeroGen.GnDa[iDGen],
MDzeroGen.GisSignal[iDGen]);
bool isPromptGen = MDzeroGen.GBAncestorpdgId[iDGen] == 0;
bool isFeeddownGen = (MDzeroGen.GBAncestorpdgId[iDGen] >= 500 && MDzeroGen.GBAncestorpdgId[iDGen] < 600) ||
(MDzeroGen.GBAncestorpdgId[iDGen] > -600 && MDzeroGen.GBAncestorpdgId[iDGen] <= -500);
Expand Down Expand Up @@ -247,13 +253,24 @@ int main(int argc, char *argv[]) {
MDzeroUPC.DpassCut->push_back(DmesonSelectionPrelim23(MDzero,iD));
if (IsData == false) {
MDzeroUPC.Dgen->push_back(MDzero.Dgen[iD]);
bool isSignalGenMatched = MDzero.Dgen[iD] == 23333 && MDzero.Dgenpt[iD] > 0.;
bool isSignalGenMatched = MDzero.Dgenpt[iD] > 0. && \
(MDzero.Dgen[iD] == 23333 || MDzero.Dgen[iD] == 41022 || MDzero.Dgen[iD] == 41044);
bool isPromptGenMatched = MDzero.DgenBAncestorpdgId[iD] == 0;
bool isFeeddownGenMatched = (MDzero.DgenBAncestorpdgId[iD] >= 500 && MDzero.DgenBAncestorpdgId[iD] < 600) ||
(MDzero.DgenBAncestorpdgId[iD] > -600 && MDzero.DgenBAncestorpdgId[iD] <= -500);
bool isSwapGenMatched = MDzero.Dgenpt[iD] > 0. && \
(MDzero.Dgen[iD] == 23344 || MDzero.Dgen[iD] == 41122 || MDzero.Dgen[iD] == 41144);
bool isD0DecayLeftSideBand = MDzero.Dgenpt[iD] > 0. && \
(MDzero.Dgen[iD] == 333 && MDzero.Dmass[iD] < D0_MASS);
bool isD0DecayRightSideBand = MDzero.Dgenpt[iD] > 0. && \
(MDzero.Dgen[iD] == 333 && MDzero.Dmass[iD] > D0_MASS);

MDzeroUPC.DisSignalCalc->push_back(isSignalGenMatched);
MDzeroUPC.DisSignalCalcPrompt->push_back(isSignalGenMatched && isPromptGenMatched);
MDzeroUPC.DisSignalCalcFeeddown->push_back(isSignalGenMatched && isFeeddownGenMatched);
MDzeroUPC.DisSwapCalc->push_back(isSwapGenMatched);
MDzeroUPC.DisD0CalcLeftSideBand->push_back(isD0DecayLeftSideBand);
MDzeroUPC.DisD0CalcRightSideBand->push_back(isD0DecayRightSideBand);
}
}
MDzeroUPC.Dsize = countSelDzero;
Expand Down Expand Up @@ -294,3 +311,56 @@ double GetMaxEnergyHF(PFTreeMessenger *M, double etaMin = 3., double etaMax = 5.
}
return EMax;
}

// ============================================================================ //
// Function to Get a Tag for Signal D Decay process with FSR at the Gen Level
// ============================================================================ //
bool GisSignalwFSR(int GpdgId, int Gdau1pdgId, int Gdau2pdgId, int Gdau3pdgId, int Gdau4pdgId, int GnDa, int GisSignal)
{
const int PION_PDGID = 211;
const int KAON_PDGID = 321;
const int DZERO_PDGID = 421;
const int GAMMA_PDGID = 22;
constexpr bool match_all_sign = true;

// 421 -> -321, 211
bool isKplusPiminusDecay = (GpdgId == DZERO_PDGID) &&
(Gdau1pdgId == (-KAON_PDGID) && Gdau2pdgId == PION_PDGID ||
Gdau2pdgId == (-KAON_PDGID) && Gdau1pdgId == PION_PDGID);
// -421 -> 321, -211
bool isKminusPiplusDecay = (GpdgId == (-DZERO_PDGID)) &&
(Gdau1pdgId == KAON_PDGID && Gdau2pdgId == (-PION_PDGID) ||
Gdau2pdgId == KAON_PDGID && Gdau1pdgId == (-PION_PDGID));
if (match_all_sign) {
// 421 -> 321, 211
isKplusPiminusDecay = (TMath::Abs(GpdgId) == DZERO_PDGID) &&
(TMath::Abs(Gdau1pdgId) == KAON_PDGID && TMath::Abs(Gdau2pdgId) == PION_PDGID ||
TMath::Abs(Gdau2pdgId) == KAON_PDGID && TMath::Abs(Gdau1pdgId) == PION_PDGID);
isKminusPiplusDecay = isKplusPiminusDecay;
}
bool isOneGamma = (Gdau3pdgId == GAMMA_PDGID);
bool isTwoGamma = (Gdau4pdgId == GAMMA_PDGID);

switch (GnDa) {
case 2: {
if ((GisSignal==1 || GisSignal==2) !=
(isKplusPiminusDecay || isKminusPiplusDecay)) {
std::cout << "Warning! " << __PRETTY_FUNCTION__ <<
"doesn't match previous definition" << "\n";
}
return isKplusPiminusDecay || isKminusPiplusDecay;
break;
}
case 3: {
return (isKplusPiminusDecay || isKminusPiplusDecay) && isOneGamma;
break;
}
case 4: {
return (isKplusPiminusDecay || isKminusPiplusDecay) && isOneGamma && isTwoGamma;
break;
}
default:
return false;
break;
}
}
12 changes: 9 additions & 3 deletions SkimValidation/20241110_SkimValidation/DumpNewSkim.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -132,13 +132,19 @@ int main(int argc, char *argv[]) {
MDzeroUPC.DsvpvDisErr_2D->at(iD), MDzeroUPC.Ddtheta->at(iD));

if (!IsData) {
fprintf(outfile, "Dgen %d DisSignalCalc %d DisSignalCalcPrompt %d DisSignalCalcFeeddown %d\n",
fprintf(outfile,
"Dgen %d DisSignalCalc %d DisSignalCalcPrompt %d DisSignalCalcFeeddown %d\n"
"DisSwapCalc %d DisD0CalcLeftSideBand %d DisD0CalcRightSideBand %d",
MDzeroUPC.Dgen->at(iD), (int)(MDzeroUPC.DisSignalCalc->at(iD)),
(int)(MDzeroUPC.DisSignalCalcPrompt->at(iD)), (int)(MDzeroUPC.DisSignalCalcFeeddown->at(iD)));
(int)(MDzeroUPC.DisSignalCalcPrompt->at(iD)), (int)(MDzeroUPC.DisSignalCalcFeeddown->at(iD)),
(int)(MDzeroUPC.DisSwapCalc->at(iD)), (int)(MDzeroUPC.DisD0CalcLeftSideBand->at(iD)),
(int)(MDzeroUPC.DisD0CalcRightSideBand->at(iD)));

v.FillRecoDGenMatchedInfo(MDzeroUPC.Dgen->at(iD), (int)(MDzeroUPC.DisSignalCalc->at(iD)),
(int)(MDzeroUPC.DisSignalCalcPrompt->at(iD)),
(int)(MDzeroUPC.DisSignalCalcFeeddown->at(iD)));
(int)(MDzeroUPC.DisSignalCalcFeeddown->at(iD)), (int)(MDzeroUPC.DisSwapCalc->at(iD)),
(int)(MDzeroUPC.DisD0CalcLeftSideBand->at(iD)),
(int)(MDzeroUPC.DisD0CalcRightSideBand->at(iD)));
}
}

Expand Down
73 changes: 68 additions & 5 deletions SkimValidation/20241110_SkimValidation/DumpOldSkim.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,12 +20,16 @@ using namespace std;

#include "ValidateHist.h"

#define D0_MASS 1.8648

// copy from main analysis
void calcRapidityGapsInput(std::vector<float> *pfE, std::vector<float> *pfPt, std::vector<float> *pfEta,
std::vector<int> *pfId, int &nPFHFMinus_, int &nPFHFPlus_);

bool tightsel(int j, DzeroTreeMessenger &MDzero, std::string varySel = "");

bool GisSignalwFSR(int GpdgId, int Gdau1pdgId, int Gdau2pdgId, int Gdau3pdgId, int Gdau4pdgId, int GnDa, int GisSignal);

int findBestVertex(PbPbUPCTrackTreeMessenger &MTrackPbPbUPC) {
int BestVertex = -1;

Expand Down Expand Up @@ -219,23 +223,35 @@ int main(int argc, char *argv[]) {
MDzero.DsvpvDistance_2D[iD], MDzero.DsvpvDisErr_2D[iD], MDzero.Ddtheta[iD]);

if (!IsData) {
bool isSignalGenMatched = MDzero.Dgen[iD] == 23333 && MDzero.Dgenpt[iD] > 0.;
bool isSignalGenMatched = MDzero.Dgenpt[iD] > 0. &&
(MDzero.Dgen[iD] == 23333 || MDzero.Dgen[iD] == 41022 || MDzero.Dgen[iD] == 41044);
bool isPromptGenMatched = isSignalGenMatched && (MDzero.DgenBAncestorpdgId[iD] == 0);
bool isFeeddownGenMatched =
isSignalGenMatched && ((MDzero.DgenBAncestorpdgId[iD] >= 500 && MDzero.DgenBAncestorpdgId[iD] < 600) ||
(MDzero.DgenBAncestorpdgId[iD] > -600 && MDzero.DgenBAncestorpdgId[iD] <= -500));
bool isSwapGenMatched = MDzero.Dgenpt[iD] > 0. &&
(MDzero.Dgen[iD] == 23344 || MDzero.Dgen[iD] == 41122 || MDzero.Dgen[iD] == 41144);
bool isD0DecayLeftSideBand = MDzero.Dgenpt[iD] > 0. && (MDzero.Dgen[iD] == 333 && MDzero.Dmass[iD] < D0_MASS);
bool isD0DecayRightSideBand = MDzero.Dgenpt[iD] > 0. && (MDzero.Dgen[iD] == 333 && MDzero.Dmass[iD] > D0_MASS);

fprintf(outfile, "Dgen %d DisSignalCalc %o DisSignalCalcPrompt %o DisSignalCalcFeeddown %o\n", MDzero.Dgen[iD],
isSignalGenMatched, isPromptGenMatched, isFeeddownGenMatched);
fprintf(outfile,
"Dgen %d DisSignalCalc %o DisSignalCalcPrompt %o DisSignalCalcFeeddown %o\n"
"DisSwapCalc %d DisD0CalcLeftSideBand %d DisD0CalcRightSideBand %d",
MDzero.Dgen[iD], isSignalGenMatched, isPromptGenMatched, isFeeddownGenMatched, isSwapGenMatched,
isD0DecayLeftSideBand, isD0DecayRightSideBand);

v.FillRecoDGenMatchedInfo(MDzero.Dgen[iD], (int)isSignalGenMatched, (int)isPromptGenMatched,
(int)isFeeddownGenMatched);
(int)isFeeddownGenMatched, (int)isSwapGenMatched, (int)isD0DecayLeftSideBand,
(int)isD0DecayRightSideBand);
}
}

if (!IsData) {
for (int iG = 0; iG < MDzeroGen.Gsize; ++iG) {
bool isSignalGen = (MDzeroGen.GisSignal[iG] == 1 || MDzeroGen.GisSignal[iG] == 2) && MDzeroGen.Gpt[iG] > 0.;
bool isSignalGen = MDzeroGen.Gpt[iG] > 0. &&
GisSignalwFSR(MDzeroGen.GpdgId[iG], MDzeroGen.Gdau1pdgId[iG], MDzeroGen.Gdau2pdgId[iG],
MDzeroGen.Gdau3pdgId[iG], MDzeroGen.Gdau4pdgId[iG], MDzeroGen.GnDa[iG],
MDzeroGen.GisSignal[iG]);
bool isPromptGen = isSignalGen && MDzeroGen.GBAncestorpdgId[iG] == 0;
bool isFeeddownGen =
isSignalGen && ((MDzeroGen.GBAncestorpdgId[iG] >= 500 && MDzeroGen.GBAncestorpdgId[iG] < 600) ||
Expand Down Expand Up @@ -444,3 +460,50 @@ bool tightsel(int j, DzeroTreeMessenger &MDzero, std::string varySel) {
// if(ishi_) cut = true; // wtf is this
return cut;
}

bool GisSignalwFSR(int GpdgId, int Gdau1pdgId, int Gdau2pdgId, int Gdau3pdgId, int Gdau4pdgId, int GnDa,
int GisSignal) {
const int PION_PDGID = 211;
const int KAON_PDGID = 321;
const int DZERO_PDGID = 421;
const int GAMMA_PDGID = 22;
constexpr bool match_all_sign = true;

// 421 -> -321, 211
bool isKplusPiminusDecay = (GpdgId == DZERO_PDGID) && (Gdau1pdgId == (-KAON_PDGID) && Gdau2pdgId == PION_PDGID ||
Gdau2pdgId == (-KAON_PDGID) && Gdau1pdgId == PION_PDGID);
// -421 -> 321, -211
bool isKminusPiplusDecay = (GpdgId == (-DZERO_PDGID)) && (Gdau1pdgId == KAON_PDGID && Gdau2pdgId == (-PION_PDGID) ||
Gdau2pdgId == KAON_PDGID && Gdau1pdgId == (-PION_PDGID));
if (match_all_sign) {
// 421 -> 321, 211
isKplusPiminusDecay = (TMath::Abs(GpdgId) == DZERO_PDGID) &&
(TMath::Abs(Gdau1pdgId) == KAON_PDGID && TMath::Abs(Gdau2pdgId) == PION_PDGID ||
TMath::Abs(Gdau2pdgId) == KAON_PDGID && TMath::Abs(Gdau1pdgId) == PION_PDGID);
isKminusPiplusDecay = isKplusPiminusDecay;
}
bool isOneGamma = (Gdau3pdgId == GAMMA_PDGID);
bool isTwoGamma = (Gdau4pdgId == GAMMA_PDGID);

switch (GnDa) {
case 2: {
if ((GisSignal == 1 || GisSignal == 2) != (isKplusPiminusDecay || isKminusPiplusDecay)) {
std::cout << "Warning! " << __PRETTY_FUNCTION__ << "doesn't match previous definition"
<< "\n";
}
return isKplusPiminusDecay || isKminusPiplusDecay;
break;
}
case 3: {
return (isKplusPiminusDecay || isKminusPiplusDecay) && isOneGamma;
break;
}
case 4: {
return (isKplusPiminusDecay || isKminusPiplusDecay) && isOneGamma && isTwoGamma;
break;
}
default:
return false;
break;
}
}
Loading