|
|
#include "../GlobalFunctions.hh"
//////////////////////////////////////////////////////////////////////
/// OptimizePreSelection()
/// loops over all .root-files created by Preselection() defined in BDTSelection.cc.
/// More kinematic, PID and some more cuts are applied and saved in a new file.
///
THStack *stack(vector<TH1F*> histos){ THStack *tmp = new THStack("IncMC","IncMC"); for(unsigned int h = 0; h < histos.size(); h++){ histos.at(h)->SetFillColor(2+h); tmp->Add(histos.at(h)); } return tmp;
}
int plotIncSample(std::string year = "2012", std::string magnet = "down") {
bool ReferenceChannel = true;
TFile* input = 0; TFile* histos = 0; TTree *tree; // Kst2Kpluspi0Resolved case
if(Kst2Kpluspi0Resolved){ input = TFile::Open(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved%s.root",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") )); histos = TFile::Open(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved_histos%s.root",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") )); } // Kst2Kpluspi0Merged case
if(Kst2Kpluspi0Merged){ //TODO
} // Kst2Kspiplus case
if(Kst2Kspiplus){ //TODO
}
input->cd(); cout << "In input" << endl; tree = (TTree*)input->Get("DecayTree");
cout << "Tree has " << tree->GetEntries() << " entries" << endl; tree->SetBranchStatus("*",1);
Double_t B_plus_M_DTF; //assign variables to TRUEID
Int_t B_plus_TRUEID; Int_t K_star_plus_TRUEID; Int_t K_plus_TRUEID; Int_t K_short_TRUEID; Int_t Ks_pi_plus_TRUEID; Int_t Ks_pi_minus_TRUEID; Int_t pi_plus_TRUEID; Int_t pi_zero_TRUEID; Int_t mu_plus_TRUEID; Int_t mu_minus_TRUEID; Int_t J_psi_TRUEID; Int_t gamma1_TRUEID; Int_t gamma2_TRUEID;
//assign variables to TRUE_MOTHER_ID
Int_t B_plus_MOTHER_ID; Int_t K_star_plus_MOTHER_ID; Int_t K_plus_MOTHER_ID; Int_t K_short_MOTHER_ID; Int_t Ks_pi_plus_MOTHER_ID; Int_t Ks_pi_minus_MOTHER_ID; Int_t pi_plus_MOTHER_ID; Int_t pi_zero_MOTHER_ID; Int_t mu_plus_MOTHER_ID; Int_t mu_minus_MOTHER_ID; Int_t J_psi_MOTHER_ID; Int_t gamma1_MOTHER_ID; Int_t gamma2_MOTHER_ID;
//assign variables to TRUE_GD_MOTHER_ID
Int_t K_plus_GD_MOTHER_ID; Int_t K_plus_GD_GD_MOTHER_ID; Int_t K_short_GD_MOTHER_ID; Int_t Ks_pi_plus_GD_MOTHER_ID; Int_t Ks_pi_minus_GD_MOTHER_ID; Int_t pi_plus_GD_MOTHER_ID; Int_t pi_zero_GD_MOTHER_ID; Int_t mu_plus_GD_MOTHER_ID; Int_t mu_minus_GD_MOTHER_ID; Int_t gamma1_GD_MOTHER_ID; Int_t gamma2_GD_MOTHER_ID;
//assign variables to TRUE_GD_MOTHER_ID
Int_t Ks_pi_plus_GD_GD_MOTHER_ID; Int_t Ks_pi_minus_GD_GD_MOTHER_ID; Int_t gamma1_GD_GD_MOTHER_ID; Int_t gamma2_GD_GD_MOTHER_ID; Int_t pi_zero_merged_MC_GD_GD_MOTHER_ID;
//link variables to branches
tree -> SetBranchAddress( "B_plus_M_DTF" , &B_plus_M_DTF );
tree -> SetBranchAddress( "B_plus_TRUEID" , &B_plus_TRUEID ); tree -> SetBranchAddress( "K_star_plus_TRUEID" , &K_star_plus_TRUEID ); tree -> SetBranchAddress( "mu_plus_TRUEID" , &mu_plus_TRUEID ); tree -> SetBranchAddress( "mu_minus_TRUEID" , &mu_minus_TRUEID );
tree -> SetBranchAddress( "B_plus_MC_MOTHER_ID" , &B_plus_MOTHER_ID ); tree -> SetBranchAddress( "K_star_plus_MC_MOTHER_ID" , &K_star_plus_MOTHER_ID ); tree -> SetBranchAddress( "mu_plus_MC_MOTHER_ID" , &mu_plus_MOTHER_ID ); tree -> SetBranchAddress( "mu_minus_MC_MOTHER_ID" , &mu_minus_MOTHER_ID ); tree -> SetBranchAddress( "mu_plus_MC_GD_MOTHER_ID" , &mu_plus_GD_MOTHER_ID ); tree -> SetBranchAddress( "mu_minus_MC_GD_MOTHER_ID" , &mu_minus_GD_MOTHER_ID );
if(Kst2Kpluspi0Resolved){ tree -> SetBranchAddress( "K_plus_TRUEID" , &K_plus_TRUEID ); tree -> SetBranchAddress( "pi_zero_resolved_TRUEID" , &pi_zero_TRUEID ); tree -> SetBranchAddress( "gamma1_TRUEID" , &gamma1_TRUEID ); tree -> SetBranchAddress( "gamma2_TRUEID" , &gamma2_TRUEID );
tree -> SetBranchAddress( "K_plus_MC_MOTHER_ID" , &K_plus_MOTHER_ID ); tree -> SetBranchAddress( "pi_zero_resolved_MC_MOTHER_ID" , &pi_zero_MOTHER_ID ); tree -> SetBranchAddress( "gamma1_MC_MOTHER_ID" , &gamma1_MOTHER_ID ); tree -> SetBranchAddress( "gamma2_MC_MOTHER_ID" , &gamma2_MOTHER_ID );
tree -> SetBranchAddress( "K_plus_MC_GD_MOTHER_ID" , &K_plus_GD_MOTHER_ID ); tree -> SetBranchAddress( "K_plus_MC_GD_GD_MOTHER_ID" , &K_plus_GD_GD_MOTHER_ID ); tree -> SetBranchAddress( "pi_zero_resolved_MC_GD_MOTHER_ID" , &pi_zero_GD_MOTHER_ID ); tree -> SetBranchAddress( "gamma1_MC_GD_MOTHER_ID" , &gamma1_GD_MOTHER_ID ); tree -> SetBranchAddress( "gamma2_MC_GD_MOTHER_ID" , &gamma2_GD_MOTHER_ID );
tree -> SetBranchAddress( "gamma1_MC_GD_GD_MOTHER_ID" , &gamma1_GD_GD_MOTHER_ID ); tree -> SetBranchAddress( "gamma2_MC_GD_GD_MOTHER_ID" , &gamma2_GD_GD_MOTHER_ID ); }
if(Kst2Kpluspi0Merged){ tree -> SetBranchAddress( "K_plus_TRUEID" , &K_plus_TRUEID ); tree -> SetBranchAddress( "pi_zero_merged_TRUEID" , &pi_zero_TRUEID );
tree -> SetBranchAddress( "K_plus_MC_MOTHER_ID" , &K_plus_MOTHER_ID ); tree -> SetBranchAddress( "pi_zero_merged_MC_MOTHER_ID" , &pi_zero_MOTHER_ID );
tree -> SetBranchAddress( "K_plus_MC_GD_MOTHER_ID" , &K_plus_GD_MOTHER_ID ); tree -> SetBranchAddress( "pi_zero_merged_MC_GD_MOTHER_ID" , &pi_zero_GD_MOTHER_ID );
tree -> SetBranchAddress( "pi_zero_merged_MC_GD_GD_MOTHER_ID" , &pi_zero_merged_MC_GD_GD_MOTHER_ID ); }
if(Kst2Kspiplus){ tree -> SetBranchAddress( "K_short_TRUEID" , &K_short_TRUEID ); tree -> SetBranchAddress( "pi_plus_TRUEID" , &pi_plus_TRUEID ); tree -> SetBranchAddress( "Ks_pi_plus_TRUEID" , &Ks_pi_plus_TRUEID ); tree -> SetBranchAddress( "Ks_pi_minus_TRUEID" , &Ks_pi_minus_TRUEID ); tree -> SetBranchAddress( "K_short_MC_MOTHER_ID" , &K_short_MOTHER_ID ); tree -> SetBranchAddress( "pi_plus_MC_MOTHER_ID" , &pi_plus_MOTHER_ID ); tree -> SetBranchAddress( "Ks_pi_plus_MC_MOTHER_ID" , &Ks_pi_plus_MOTHER_ID ); tree -> SetBranchAddress( "Ks_pi_minus_MC_MOTHER_ID" , &Ks_pi_minus_MOTHER_ID );
tree -> SetBranchAddress( "K_short_MC_GD_MOTHER_ID" , &K_short_GD_MOTHER_ID ); tree -> SetBranchAddress( "pi_plus_MC_GD_MOTHER_ID" , &pi_plus_GD_MOTHER_ID ); tree -> SetBranchAddress( "Ks_pi_plus_MC_GD_MOTHER_ID" , &Ks_pi_plus_GD_MOTHER_ID ); tree -> SetBranchAddress( "Ks_pi_minus_MC_GD_MOTHER_ID" , &Ks_pi_minus_GD_MOTHER_ID );
tree -> SetBranchAddress( "Ks_pi_plus_MC_GD_GD_MOTHER_ID" , &Ks_pi_plus_GD_GD_MOTHER_ID ); tree -> SetBranchAddress( "Ks_pi_minus_MC_GD_GD_MOTHER_ID" , &Ks_pi_minus_GD_GD_MOTHER_ID ); }
tree -> SetBranchAddress( "J_psi_TRUEID" , &J_psi_TRUEID ); tree -> SetBranchAddress( "J_psi_MC_MOTHER_ID" , &J_psi_MOTHER_ID );
//counter
Int_t B_plus_VETO = 0; Int_t K_star_plus_VETO = 0; Int_t K_plus_VETO = 0; Int_t K_short_VETO = 0; Int_t Ks_pi_plus_VETO = 0; Int_t Ks_pi_minus_VETO = 0; Int_t pi_plus_VETO = 0; Int_t pi_zero_VETO = 0; Int_t mu_plus_VETO = 0; Int_t mu_minus_VETO = 0; Int_t J_psi_VETO = 0; Int_t gamma1_VETO = 0; Int_t gamma2_VETO = 0;
Int_t MOTHER_ID_VETO = 0; Int_t GD_MOTHER_ID_VETO = 0; Int_t GD_GD_MOTHER_ID_VETO = 0;
Int_t TRUE_events = 0;
//count the numbers of Kshort daughters, which are pi_zeros
Int_t Kshort_pi_zero_daughters = 0;
//count muons, mis-IDed as pions
Int_t ReconstructedPionsWhichAreMuons = 0;
//======================================
//Get histos of all possible backgrounds
//======================================
//control plots for pi0
Double_t pi_zero_resolved_M = 0.0; if (Kst2Kpluspi0Resolved) tree -> SetBranchAddress( "pi_zero_resolved_M" , &pi_zero_resolved_M );
TH1 * h_pi0_resolved_mass= new TH1D("h_pi0_resolved_mass" , "m_{#pi_{0}} resolved", 100, 75, 175); h_pi0_resolved_mass->GetXaxis()->SetTitle("m_{pi_{0}} resolved [MeV]"); h_pi0_resolved_mass->GetYaxis()->SetTitle("Counts / 1MeV");
TH1 * h_pi0_mismatched_mass= new TH1D("h_pi0_mismatched_mass" , "m_{#pi_{0}} IDed as gamma", 100, 75, 175); h_pi0_mismatched_mass->GetXaxis()->SetTitle("m_{pi_{0}}[MeV] MisIDed"); h_pi0_mismatched_mass->GetYaxis()->SetTitle("Counts / 1MeV");
TH1 * h_B_plus_TRUEID = new TH1D("h_B_plus_TRUEID", "B^{+} TRUEID", 100000, -1, 99999); TH1 * h_K_star_plus_TRUEID = new TH1D("h_K_star_plus_TRUEID", "K^{*+} TRUEID", 100000, -1, 99999); TH1 * h_K_plus_TRUEID = new TH1D("h_K_plus_TRUEID", "K^{+} TRUEID", 100000, -1, 99999); TH1 * h_pi_zero_TRUEID = new TH1D("h_pi_zero_TRUEID", "#pi^{0} TRUEID", 100000, -1, 99999);
TH1 * h_B_plus_TRUEID_rm = new TH1D("h_B_plus_TRUEID_rm", "B^{+} TRUEID", 100000, -1, 99999); TH1 * h_K_star_plus_TRUEID_rm= new TH1D("h_K_star_plus_TRUEID_rm", "K^{*+} TRUEID", 100000, -1, 99999); TH1 * h_K_plus_TRUEID_rm = new TH1D("h_K_plus_TRUEID_rm", "K^{+} TRUEID", 100000, -1, 99999); TH1 * h_pi_zero_TRUEID_rm = new TH1D("h_pi_zero_TRUEID_rm", "#pi^{0} TRUEID", 100000, -1, 99999);
TH1 * h_B_plus_K1 = new TH1D("h_B_plus_K1", "h_B_plus_K1", 100, cut_B_plus_M_low, cut_B_plus_M_high); TH1 * h_B_plus_Kst = new TH1D("h_B_plus_Kst", "h_B_plus_Kst", 100, cut_B_plus_M_low, cut_B_plus_M_high); TH1 * h_B_plus_Kpl = new TH1D("h_B_plus_Kpl", "h_B_plus_Kpl", 100, cut_B_plus_M_low, cut_B_plus_M_high);
vector<std::string> Cuts; string tmp_string = ""; int nBins = 75; TH1 * h_BplusMass_tot= new TH1D("h_BplusMass_tot" , "h_BplusMass_tot", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_tot->GetXaxis()->SetTitle("B^{+} mass [MeV]");
//======================================
// Define trees for all channels
//======================================
TTree *tree_correct = tree->CloneTree(0); tree_correct->SetName("tree_correct"); TH1 * h_BplusMass_correct = new TH1D("h_BplusMass_correct" , "h_BplusMass_correct", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_correct->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_K1_1270 = tree->CloneTree(0); tree_K1_1270->SetName("tree_K1_1270"); TH1 * h_BplusMass_K1_1270 = new TH1D("h_BplusMass_K1_1270" , "h_BplusMass_K1_1270", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_K1_1270->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_K1_1410 = tree->CloneTree(0); tree_K1_1410->SetName("tree_K1_1410"); TH1 * h_BplusMass_K1_1410 = new TH1D("h_BplusMass_K1_1410" , "h_BplusMass_K1_1410", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_K1_1410->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_K1_1400 = tree->CloneTree(0); tree_K1_1400->SetName("tree_K1_1400"); TH1 * h_BplusMass_K1_1400 = new TH1D("h_BplusMass_K1_1400" , "h_BplusMass_K1_1400", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_K1_1400->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_K_plus= tree->CloneTree(0); tree_K_plus->SetName("tree_K_plus"); TH1 * h_BplusMass_K_plus = new TH1D("h_BplusMass_K_plus" , "h_BplusMass_K_plus", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_K_plus->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_K_plus_misID= tree->CloneTree(0); tree_K_plus_misID->SetName("tree_K_plus_misID"); TH1 * h_BplusMass_K_plus_misID = new TH1D("h_BplusMass_K_plus_misID" , "h_BplusMass_K_plus_misID", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_K_plus_misID->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_B_star= tree->CloneTree(0); tree_B_star->SetName("tree_B_star"); TH1 * h_BplusMass_B_star = new TH1D("h_BplusMass_B_star" , "h_BplusMass_B_star", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_B_star->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_Kstar_misID= tree->CloneTree(0); tree_Kstar_misID->SetName("tree_Kstar_misID"); TH1 * h_BplusMass_Kstar_misID = new TH1D("h_BplusMass_Kstar_misID" , "h_BplusMass_Kstar_misID", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_Kstar_misID->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_bquark= tree->CloneTree(0); tree_bquark->SetName("tree_bquark"); TH1 * h_BplusMass_bquark = new TH1D("h_BplusMass_bquark" , "h_BplusMass_bquark", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_bquark->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_combinatorics= tree->CloneTree(0); tree_combinatorics->SetName("tree_combinatorics"); TH1 * h_BplusMass_combinatorics = new TH1D("h_BplusMass_combinatorics" , "h_BplusMass_combinatorics", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_combinatorics->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TTree *tree_others = tree->CloneTree(0); tree_others->SetName("tree_others"); TH1 * h_BplusMass_others = new TH1D("h_BplusMass_others" , "h_BplusMass_others", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_others->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TH1 * h_BplusMass_JpsiVETO = new TH1D("h_BplusMass_JpsiVETO" , "h_BplusMass_JpsiVETO", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_JpsiVETO->GetXaxis()->SetTitle("B^{+} mass [MeV]");
TH1 * h_BplusMass_GDVETO = new TH1D("h_BplusMass_GDVETO" , "h_BplusMass_GDVETO", nBins ,cut_B_plus_M_low, cut_B_plus_M_high); h_BplusMass_GDVETO->GetXaxis()->SetTitle("B^{+} mass [MeV]");
//======================================
//check all particles for correct TRUEID
//======================================
//MC loop
cout << "Tree has " << tree->GetEntries() << " entries" << endl; int nEvents = tree->GetEntries(); std::cout << "Truth matching of " << nEvents << " INC MC events for " << year << magnet << "!" << std::endl; for(int i = 0; i < nEvents; i++){ if(i%2000==0 && i != 0)std::cout << "Loading MC event " << i << " / " << nEvents << std::endl; bool flag_other = true;
//others
tree->GetEntry(i); h_BplusMass_tot->Fill(B_plus_M_DTF); //remove nonsense from B and Jpsi
if(TMath::Abs(B_plus_TRUEID) != TRUEID.B_PLUS && TMath::Abs(B_plus_TRUEID) != 0 && TMath::Abs(B_plus_TRUEID) != 5){ B_plus_VETO++; h_B_plus_TRUEID_rm->Fill(B_plus_TRUEID); continue; } if(ReferenceChannel){ if(TMath::Abs(J_psi_TRUEID) != TRUEID.J_PSI){ J_psi_VETO++; h_BplusMass_JpsiVETO->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_plus_MOTHER_ID) != TRUEID.J_PSI){ MOTHER_ID_VETO++; h_BplusMass_others->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_minus_MOTHER_ID) != TRUEID.J_PSI){ MOTHER_ID_VETO++; h_BplusMass_others->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_plus_GD_MOTHER_ID) != TRUEID.B_PLUS){ GD_MOTHER_ID_VETO++; h_BplusMass_GDVETO->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_minus_GD_MOTHER_ID) != TRUEID.B_PLUS){ GD_MOTHER_ID_VETO++; h_BplusMass_GDVETO->Fill(B_plus_M_DTF); continue; } } else{ if(TMath::Abs(J_psi_TRUEID) == TRUEID.J_PSI){ J_psi_VETO++; h_BplusMass_JpsiVETO->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_plus_MOTHER_ID) != TRUEID.B_PLUS){ MOTHER_ID_VETO++; h_BplusMass_others->Fill(B_plus_M_DTF); continue; } if(TMath::Abs(mu_minus_MOTHER_ID) != TRUEID.B_PLUS){ MOTHER_ID_VETO++; h_BplusMass_others->Fill(B_plus_M_DTF); continue; } }
//correct tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
tree_correct->Fill(); h_BplusMass_correct->Fill(B_plus_M_DTF); continue; } } }
//K1 tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_TRUEID) == 0 || TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS ) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == 0) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
tree_K1_1270->Fill(); h_BplusMass_K1_1270->Fill(B_plus_M_DTF); continue; } } }
//K1 1400 tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS_1400 || TMath::Abs(K_star_plus_TRUEID) == 0|| TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS ) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle ID
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400){//Grandmother IDs
tree_K1_1400->Fill(); h_BplusMass_K1_1410->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } } } //K1 1410 tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS_1410 || TMath::Abs(K_star_plus_TRUEID) == 0 ) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == 0 || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 0 || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410){//Grandmother IDs
tree_K1_1410->Fill(); h_BplusMass_K1_1410->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } }
} //K1 1270 tree 2.0
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
tree_K1_1270->Fill(); h_BplusMass_K1_1270->Fill(B_plus_M_DTF); continue; } } }
//tree_K_plus_misID tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS|| TMath::Abs(K_star_plus_TRUEID) == 0 ) if(TMath::Abs(K_plus_TRUEID) == TRUEID.PI_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
tree_K_plus_misID->Fill(); h_BplusMass_K_plus_misID->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } } } //tree_K_plus_misID tree 2.0 (K+->pi+)
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.PI_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_STAR_PLUS) if (TMath::Abs(K_plus_GD_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
tree_K_plus_misID->Fill(); h_BplusMass_K_plus_misID->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } } }
//Kplus tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS || TMath::Abs(K_star_plus_TRUEID) == 0) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 0){//Grandmother IDs
tree_K_plus->Fill(); h_BplusMass_K_plus->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } } }
//K star misID tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
tree_Kstar_misID->Fill(); h_BplusMass_Kstar_misID->Fill(B_plus_M_DTF); tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF); continue; } } }
//B+ <-> b quark
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 5) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
tree_bquark->Fill(); h_BplusMass_bquark->Fill(B_plus_M_DTF); continue; } }
//B+ <-> b quark 2.0
if (TMath::Abs(K_star_plus_TRUEID) ==5) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 0) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS ){ //Mother IDs
tree_bquark->Fill(); h_BplusMass_bquark->Fill(B_plus_M_DTF); continue; } }
//Bstar tree
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 523) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
tree_B_star->Fill(); h_BplusMass_B_star->Fill(B_plus_M_DTF); continue; } } }
//Bstar tree 2.0
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 523) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
tree_B_star->Fill(); h_BplusMass_B_star->Fill(B_plus_M_DTF); continue; } } }
//Bstar tree 3.0 (K* misID B+)
if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
tree_B_star->Fill(); h_BplusMass_B_star->Fill(B_plus_M_DTF); continue; } } } //Bstar tree 4.0 (K* misID B+)
if (TMath::Abs(K_star_plus_TRUEID) == 0) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 0) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523 ||TMath::Abs(K_plus_GD_MOTHER_ID) == 10521){//Grandmother IDs
tree_B_star->Fill(); h_BplusMass_B_star->Fill(B_plus_M_DTF); continue; } } } //combinatorics tree
if (TMath::Abs(K_star_plus_TRUEID) == 0) if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS ||TMath::Abs(K_plus_TRUEID) == 0 ){ //Paticle IDs
if(TMath::Abs(K_star_plus_MOTHER_ID) == 0) if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS ||TMath::Abs(K_plus_TRUEID) == 0 ){ //Mother IDs
if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
tree_combinatorics->Fill(); h_BplusMass_combinatorics->Fill(B_plus_M_DTF); continue; } } }
//the rest of stuff
tree_others->Fill(); h_BplusMass_others->Fill(B_plus_M_DTF);
} //output truthmatching results
std::cout << "Resulting TruthMatching data: " << TRUE_events << "/" << nEvents << " are found to be correctly reconstructed!" << std::endl; std::cout << "The incorrect reconstructed particles are the following:" << std::endl; std::cout << "(only the first found-incorrect particle is counted)" << std::endl; std::cout << std::endl; //std::cout << (B0 ? "B_zero:\t\t" : "B_plus:\t\t") << B_plus_VETO << std::endl;
//std::cout << (B0 ? "K_star_zero:\t" : ( K1 ? "K_one_plus:\t" : "K_star_plus:\t")) << K_star_plus_VETO << std::endl;
std::cout << "mu_plus:\t" << mu_plus_VETO << std::endl; std::cout << "mu_minus:\t" << mu_minus_VETO << std::endl; if(Kst2Kpluspi0Resolved){ std::cout << "K_plus:\t\t" << K_plus_VETO << std::endl; std::cout << "pi_zero:\t" << pi_zero_VETO << std::endl; std::cout << "gamma1:\t\t" << gamma1_VETO << std::endl; std::cout << "gamma2:\t\t" << gamma2_VETO << std::endl; } if(Kst2Kpluspi0Merged){ std::cout << "K_plus:\t\t" << K_plus_VETO << std::endl; std::cout << "pi_zero:\t" << pi_zero_VETO << std::endl; } if(Kst2Kspiplus){ std::cout << "K_short:\t" << K_short_VETO << std::endl; std::cout << "Ks_pi_plus:\t" << Ks_pi_plus_VETO << std::endl; std::cout << "Ks_pi_minus:\t" << Ks_pi_minus_VETO << std::endl; std::cout << "pi_plus:\t" << pi_plus_VETO << std::endl; } std::cout << "J_psi:\t\t" <<J_psi_VETO << std::endl;
if(Kst2Kspiplus){ std::cout << "Pi0 daughers of the Kshort:\t" << Kshort_pi_zero_daughters << std::endl; std::cout << "Muons mis-IDed as pions:\t" << ReconstructedPionsWhichAreMuons << "/" << Ks_pi_plus_VETO + Ks_pi_minus_VETO + pi_plus_VETO << std::endl; }
std::cout << "Total incorrect family members of all particles: " << MOTHER_ID_VETO + GD_MOTHER_ID_VETO + GD_GD_MOTHER_ID_VETO << std::endl; std::cout << "MotherID:\t" << MOTHER_ID_VETO << std::endl; std::cout << "GDmotherID:\t" << GD_MOTHER_ID_VETO << std::endl; std::cout << "GDGDmotherID:\t" << GD_GD_MOTHER_ID_VETO << std::endl;
//save to file
TFile * outputFile = nullptr;
// Kst2Kpluspi0Resolved case
if(Kst2Kpluspi0Resolved){ outputFile = new TFile(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved_plots.root",path_to_output_KplusPizero.c_str(),year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str()), "RECREATE"); outputFile->cd(); tree_correct->Write(); cout << "Correctly id'd entries " << tree_correct->GetEntries() << endl; tree_K1_1270->Write(); cout << "B+->K1(1270) entries " << tree_K1_1270->GetEntries() << endl; tree_K1_1400->Write(); cout << "B+->K1(1400) entries " << tree_K1_1400->GetEntries() << endl; tree_K1_1410->Write(); cout << "B+->K1(1410) entries " << tree_K1_1410->GetEntries() << endl; tree_K_plus->Write(); cout << "B+->K+pi0 entries " << tree_K_plus->GetEntries() << endl; tree_K_plus_misID->Write(); cout << "B+->K*->pi+pi0 " << tree_K_plus_misID->GetEntries() << endl; tree_B_star->Write(); cout << "B*-> B+->K+ entries " << tree_B_star->GetEntries() << endl; tree_Kstar_misID->Write(); cout << "K* misID entries " << tree_Kstar_misID->GetEntries() << endl; tree_combinatorics->Write(); cout << "b-quark entries " << tree_bquark->GetEntries() << endl; tree_bquark->Write(); cout << "Combinatorics " << tree_combinatorics->GetEntries() << endl; tree_others->Write(); cout << "Other entries " << tree_others->GetEntries() << endl;
THStack *tmp = new THStack("IncMC","IncMC;B^{+} mass [MeV];Entries"); //tmp->GetXaxis()->SetTitle("B^{+} mass [MeV]");
//tmp->GetYaxis()->SetTitle("Entries");
h_BplusMass_others->SetFillColor(kGreen); h_BplusMass_others->SetLineColor(kGreen); tmp->Add(h_BplusMass_others); h_BplusMass_bquark->SetFillColor(kGreen); h_BplusMass_bquark->SetLineColor(kGreen); tmp->Add(h_BplusMass_bquark);
h_BplusMass_K1_1270->SetFillColor(kBlue); h_BplusMass_K1_1270->SetLineColor(kBlue); tmp->Add(h_BplusMass_K1_1270);
// h_BplusMass_K1_1400->SetFillColor(0);
// h_BplusMass_K1_1400->SetLineColor(0);
// tmp->Add(h_BplusMass_K1_1400);
// h_BplusMass_K1_1410->SetFillColor(0);
// h_BplusMass_K1_1410->SetLineColor(0);
// tmp->Add(h_BplusMass_K1_1410);
// h_BplusMass_K_plus->SetFillColor(0);
// h_BplusMass_K_plus->SetLineColor(0);
// tmp->Add(h_BplusMass_K_plus);
// h_BplusMass_K_plus_misID->SetFillColor(0);
// h_BplusMass_K_plus_misID->SetLineColor(0);
// tmp->Add(h_BplusMass_K_plus_misID);
h_BplusMass_B_star->SetFillColor(kYellow); h_BplusMass_B_star->SetLineColor(kYellow); tmp->Add(h_BplusMass_B_star);
// h_BplusMass_Kstar_misID->SetFillColor(0);
// h_BplusMass_Kstar_misID->SetLineColor(0);
// tmp->Add(h_BplusMass_Kstar_misID);
//
// h_BplusMass_combinatorics->SetFillColor(kBlack);
// h_BplusMass_combinatorics->SetLineColor(kBlack);
// tmp->Add(h_BplusMass_combinatorics);
h_BplusMass_JpsiVETO->SetFillColor(kMagenta-10); h_BplusMass_JpsiVETO->SetLineColor(kMagenta-10); tmp->Add(h_BplusMass_JpsiVETO);
h_BplusMass_GDVETO->SetFillColor(kOrange-3); h_BplusMass_GDVETO->SetLineColor(kOrange-3); tmp->Add(h_BplusMass_GDVETO);
h_BplusMass_correct->SetFillColor(kRed); h_BplusMass_correct->SetLineColor(kRed); tmp->Add(h_BplusMass_correct);
tmp->Write();
TCanvas *c1 = new TCanvas("c1","all",900,700); c1->cd(); gStyle->SetOptStat(0); // h_BplusMass_correct->SetFillColor(0);
h_BplusMass_correct->Draw(); // h_BplusMass_others->SetFillColor(0);
h_BplusMass_others->Draw("SAME"); // h_BplusMass_K1_1270->SetFillColor(0);
h_BplusMass_K1_1270->Draw("SAME"); //h_BplusMass_K1_1400->Draw("SAME");
//h_BplusMass_K1_1410->Draw("SAME");
//h_BplusMass_K_plus->Draw("SAME");
//h_BplusMass_K_plus_misID->Draw("SAME");
// h_BplusMass_B_star->SetFillColor(0);
h_BplusMass_B_star->Draw("SAME"); //h_BplusMass_Kstar_misID->Draw("SAME");
//h_BplusMass_combinatorics->Draw("SAME");
// h_BplusMass_bquark->SetFillColor(0);
h_BplusMass_bquark->Draw("SAME"); // h_BplusMass_JpsiVETO->SetFillColor(0);
h_BplusMass_JpsiVETO->Draw("SAME"); // h_BplusMass_GDVETO->SetFillColor(0);
h_BplusMass_GDVETO->Draw("SAME");
c1->Write();
TCanvas *c2 = new TCanvas("c2","al2",900,700); c2->cd(); tmp->Draw(); h_BplusMass_tot->SetLineWidth(3); h_BplusMass_tot->SetLineColor(kBlack); h_BplusMass_tot->Draw("SAMEE"); c2->Write(); c2->SaveAs(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved%s_BKG.eps",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") ),"eps");
} // Kst2Kpluspi0Merged case
else if(Kst2Kpluspi0Merged){ //TODO
}
// Kst2Kspiplus case
else if(Kst2Kspiplus){ //TODO
}
else return 0;
if(!outputFile->IsOpen()){ std::cout << "[ERROR]\t\tFile was not opened succesfully!" << std::endl; return 0; }
outputFile->cd(); //treeTruthMatched->Write("",TObject::kWriteDelete);
if (Kst2Kpluspi0Resolved){ h_pi0_mismatched_mass->SetLineColor(kRed); h_pi0_mismatched_mass->SetLineWidth(3); h_pi0_resolved_mass->SetLineColor(kBlue); h_pi0_resolved_mass->SetLineWidth(3); h_pi0_mismatched_mass->Write("",TObject::kWriteDelete); h_pi0_resolved_mass->Write("",TObject::kWriteDelete);
h_B_plus_TRUEID ->Write("",TObject::kWriteDelete); h_K_star_plus_TRUEID ->Write("",TObject::kWriteDelete); h_K_plus_TRUEID ->Write("",TObject::kWriteDelete); h_pi_zero_TRUEID ->Write("",TObject::kWriteDelete);
h_B_plus_TRUEID_rm ->Write("",TObject::kWriteDelete); h_K_star_plus_TRUEID_rm ->Write("",TObject::kWriteDelete); h_K_plus_TRUEID_rm ->Write("",TObject::kWriteDelete); h_pi_zero_TRUEID_rm ->Write("",TObject::kWriteDelete);
h_B_plus_K1 ->Write("",TObject::kWriteDelete); h_B_plus_Kst ->Write("",TObject::kWriteDelete); h_B_plus_Kpl ->Write("",TObject::kWriteDelete); }
delete outputFile;//->Close();
delete h_pi0_mismatched_mass; delete h_pi0_resolved_mass;
delete tree; std::cout << "[DONE]\t\tFinished Truth-Matching for INC sample: " << year << " " << magnet << std::endl;
return 1;}
int plotIncSampleAll() {
if ( plotIncSample("2011","down") == 0) return 0; if ( plotIncSample("2011","up") == 0) return 0; if ( plotIncSample("2012","down") == 0) return 0; if ( plotIncSample("2012","up") == 0) return 0; if ( plotIncSample("2015","down") == 0) return 0; if ( plotIncSample("2015","up") == 0) return 0; if ( plotIncSample("2016","down") == 0) return 0; if ( plotIncSample("2016","up") == 0) return 0; return 1;
}
|
