diff --git a/B0ToHpHmMuMu_results.txt b/B0ToHpHmMuMu_results.txt index 1dabbae..2c3d791 100644 --- a/B0ToHpHmMuMu_results.txt +++ b/B0ToHpHmMuMu_results.txt @@ -1,11 +1,11 @@ -#### B0ToHpHmMuMu @ Tue Apr 2 10:44:44 2024 #### -J/Psi Mode: 341 #pm 23 / 185 #pm 23 - Sig/Bkg: 1.85 #pm 0.26 -Psi(2S) Mode: 23 #pm 6 / 43 #pm 11 - Sig/Bkg: 0.54 #pm 0.19 -Mode Yield Ratio: 0.068 #pm 0.018 +#### B0ToHpHmMuMu @ Fri Apr 5 14:10:32 2024 #### +J/Psi Mode: 252 #pm 19 / 59 #pm 14 + Sig/Bkg: 4.25 #pm 1.02 +Psi(2S) Mode: 15 #pm 4 / 3 #pm 4 + Sig/Bkg: 4.59 #pm 4.45 +Mode Yield Ratio: 0.060 #pm 0.017 Rel Br Frac MuMu: 7.701 #pm 0.804 -Rel Br Frac: 0.522 #pm 0.149 +Rel Br Frac: 0.463 #pm 0.136 Fitted Parameters: Simulation @@ -19,14 +19,14 @@ n_{R} = 5.67 \pm 0.73 Fitted Parameters: J/PSI -N_{Sig} = 341 \pm 23 -N_{Bkg} = 185 \pm 23 -N_{Bkg,2#sigma} = 97 \pm 6 -N_{Sig}/N_{Bkg} = 1.85 \pm 0.26 -N_{Sig}/N_{Bkg,2#sigma} = 3.53 \pm 0.32 -#lambda = -0.00217 \pm 0.00066 -#mu = 5270.03 \pm 1.92 -#sigma_{LR} = 26.26 \pm 1.98 +N_{Sig} = 252 \pm 19 +N_{Bkg} = 59 \pm 14 +N_{Bkg,2#sigma} = 79 \pm 6 +N_{Sig}/N_{Bkg} = 4.25 \pm 1.02 +N_{Sig}/N_{Bkg,2#sigma} = 3.19 \pm 0.33 +#lambda = -0.00400 \pm 0.00395 +#mu = 5267.82 \pm 2.42 +#sigma_{LR} = 29.94 \pm 2.37 #alpha_{L} = 1.72 (c) n_{L} = 3.34 (c) #alpha_{R} = 1.76 (c) @@ -34,14 +34,14 @@ n_{R} = 5.67 (c) Fitted Parameters: PSI(2S) -N_{Sig} = 23 \pm 6 -N_{Bkg} = 43 \pm 10 -N_{Bkg,2#sigma} = 12 \pm 2 -N_{Sig}/N_{Bkg} = 0.54 \pm 0.19 -N_{Sig}/N_{Bkg,2#sigma} = 1.93 \pm 0.60 -#lambda = -0.00202 \pm 0.00120 -#mu = 5267.35 \pm 7.06 -#sigma_{LR} = 26.26 (c) +N_{Sig} = 15 \pm 4 +N_{Bkg} = 3 \pm 3 +N_{Bkg,2#sigma} = 2 \pm 1 +N_{Sig}/N_{Bkg} = 4.59 \pm 4.44 +N_{Sig}/N_{Bkg,2#sigma} = 6.17 \pm 2.34 +#lambda = -0.00000 \pm 0.00320 +#mu = 5275.52 \pm 8.60 +#sigma_{LR} = 29.94 (c) #alpha_{L} = 1.72 (c) n_{L} = 3.34 (c) #alpha_{R} = 1.76 (c) @@ -50,11 +50,11 @@ n_{R} = 5.67 (c) # J/psi \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$341 \pm 23$ & $185 \pm 23$ +$252 \pm 19$ & $59 \pm 14$ \end{tabular} # psi(2S) \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$23 \pm 6$ & $43 \pm 10$ +$15 \pm 4$ & $3 \pm 3$ \end{tabular} diff --git a/BuToHpMuMu_results.txt b/BuToHpMuMu_results.txt index 368349e..1f6c90e 100644 --- a/BuToHpMuMu_results.txt +++ b/BuToHpMuMu_results.txt @@ -1,18 +1,18 @@ -#### BuToHpMuMu @ Tue Apr 2 10:44:14 2024 #### -J/Psi Mode: 897 #pm 34 / 402 #pm 30 - Sig/Bkg: 2.23 #pm 0.19 -Psi(2S) Mode: 69 #pm 12 / 233 #pm 26 - Sig/Bkg: 0.30 #pm 0.06 -Mode Yield Ratio: 0.077 #pm 0.013 +#### BuToHpMuMu @ Wed Apr 10 11:06:16 2024 #### +J/Psi Mode: 733 #pm 30 / 174 #pm 20 + Sig/Bkg: 4.20 #pm 0.49 +Psi(2S) Mode: 64 #pm 9 / 43 #pm 11 + Sig/Bkg: 1.50 #pm 0.42 +Mode Yield Ratio: 0.087 #pm 0.013 Rel Br Frac MuMu: 7.701 #pm 0.804 -Rel Br Frac: 0.594 #pm 0.116 +Rel Br Frac: 0.673 #pm 0.121 Params from Sim: aL: 1.76 nL: 4.02 aR: 1.93 nR: 4.33 -S: 23.48 +S: 23.72 #mu = 5278.16 \pm 0.05 #sigma_{LR} = 16.43 \pm 0.05 @@ -23,14 +23,14 @@ n_{R} = 4.33 \pm 0.19 Fitted Parameters: J/PSI -N_{Sig} = 897 \pm 34 -N_{Bkg} = 402 \pm 30 -N_{Bkg,2#sigma} = 246 \pm 9 -N_{Sig}/N_{Bkg} = 2.23 \pm 0.19 -N_{Sig}/N_{Bkg,2#sigma} = 3.64 \pm 0.19 -#lambda = -0.00350 \pm 0.00041 -#mu = 5271.03 \pm 0.96 -#sigma_{LR} = 23.48 \pm 0.90 +N_{Sig} = 733 \pm 29 +N_{Bkg} = 174 \pm 19 +N_{Bkg,2#sigma} = 179 \pm 7 +N_{Sig}/N_{Bkg} = 4.20 \pm 0.49 +N_{Sig}/N_{Bkg,2#sigma} = 4.08 \pm 0.23 +#lambda = -0.00400 \pm 0.00018 +#mu = 5271.39 \pm 1.03 +#sigma_{LR} = 23.72 \pm 0.94 #alpha_{L} = 1.76 (c) n_{L} = 4.02 (c) #alpha_{R} = 1.93 (c) @@ -38,14 +38,14 @@ n_{R} = 4.33 (c) Fitted Parameters: PSI(2S) -N_{Sig} = 69 \pm 11 -N_{Bkg} = 233 \pm 26 -N_{Bkg,2#sigma} = 44 \pm 4 -N_{Sig}/N_{Bkg} = 0.30 \pm 0.06 -N_{Sig}/N_{Bkg,2#sigma} = 1.56 \pm 0.29 -#lambda = -0.00142 \pm 0.00056 -#mu = 5263.32 \pm 4.14 -#sigma_{LR} = 23.48 (c) +N_{Sig} = 64 \pm 9 +N_{Bkg} = 43 \pm 10 +N_{Bkg,2#sigma} = 18 \pm 2 +N_{Sig}/N_{Bkg} = 1.50 \pm 0.42 +N_{Sig}/N_{Bkg,2#sigma} = 3.56 \pm 0.67 +#lambda = -0.00217 \pm 0.00121 +#mu = 5264.03 \pm 3.67 +#sigma_{LR} = 23.72 (c) #alpha_{L} = 1.76 (c) n_{L} = 4.02 (c) #alpha_{R} = 1.93 (c) @@ -54,11 +54,11 @@ n_{R} = 4.33 (c) # J/psi \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$897 \pm 34$ & $402 \pm 30$ +$733 \pm 29$ & $174 \pm 19$ \end{tabular} # psi(2S) \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$69 \pm 11$ & $233 \pm 26$ +$64 \pm 9$ & $43 \pm 10$ \end{tabular} diff --git a/BuToKpMuMu_results.txt b/BuToKpMuMu_results.txt index 1d3f300..f67f124 100644 --- a/BuToKpMuMu_results.txt +++ b/BuToKpMuMu_results.txt @@ -1,7 +1,7 @@ -#### BuToKpMuMu @ Tue Apr 2 10:45:35 2024 #### -J/Psi Mode: 815 #pm 35 / 743 #pm 45 - Sig/Bkg: 1.10 #pm 0.09 -Psi(2S) Mode: 69 #pm 12 / 290 #pm 30 +#### BuToKpMuMu @ Wed Apr 10 10:51:08 2024 #### +J/Psi Mode: 813 #pm 35 / 734 #pm 44 + Sig/Bkg: 1.11 #pm 0.09 +Psi(2S) Mode: 69 #pm 12 / 283 #pm 30 Sig/Bkg: 0.24 #pm 0.05 Mode Yield Ratio: 0.085 #pm 0.015 Rel Br Frac MuMu: 7.701 #pm 0.804 @@ -12,7 +12,7 @@ aL: 1.77 nL: 3.86 aR: 1.96 nR: 4.43 -S: 24.12 +S: 24.08 #mu = 5278.10 \pm 0.05 #sigma_{LR} = 16.35 \pm 0.05 @@ -23,14 +23,14 @@ n_{R} = 4.43 \pm 0.20 Fitted Parameters: J/PSI -N_{Sig} = 815 \pm 34 -N_{Bkg} = 743 \pm 44 -N_{Bkg,2#sigma} = 259 \pm 9 -N_{Sig}/N_{Bkg} = 1.10 \pm 0.08 -N_{Sig}/N_{Bkg,2#sigma} = 3.15 \pm 0.18 -#lambda = -0.00207 \pm 0.00031 -#mu = 5271.47 \pm 1.08 -#sigma_{LR} = 24.12 \pm 1.02 +N_{Sig} = 813 \pm 34 +N_{Bkg} = 734 \pm 44 +N_{Bkg,2#sigma} = 257 \pm 9 +N_{Sig}/N_{Bkg} = 1.11 \pm 0.08 +N_{Sig}/N_{Bkg,2#sigma} = 3.16 \pm 0.18 +#lambda = -0.00209 \pm 0.00031 +#mu = 5271.35 \pm 1.08 +#sigma_{LR} = 24.08 \pm 1.02 #alpha_{L} = 1.77 (c) n_{L} = 3.86 (c) #alpha_{R} = 1.96 (c) @@ -38,14 +38,14 @@ n_{R} = 4.43 (c) Fitted Parameters: PSI(2S) -N_{Sig} = 69 \pm 12 -N_{Bkg} = 290 \pm 30 -N_{Bkg,2#sigma} = 53 \pm 5 +N_{Sig} = 69 \pm 11 +N_{Bkg} = 283 \pm 29 +N_{Bkg,2#sigma} = 52 \pm 5 N_{Sig}/N_{Bkg} = 0.24 \pm 0.05 -N_{Sig}/N_{Bkg,2#sigma} = 1.32 \pm 0.25 -#lambda = -0.00130 \pm 0.00051 -#mu = 5268.10 \pm 4.25 -#sigma_{LR} = 24.12 (c) +N_{Sig}/N_{Bkg,2#sigma} = 1.33 \pm 0.25 +#lambda = -0.00136 \pm 0.00051 +#mu = 5268.31 \pm 4.17 +#sigma_{LR} = 24.08 (c) #alpha_{L} = 1.77 (c) n_{L} = 3.86 (c) #alpha_{R} = 1.96 (c) @@ -54,11 +54,11 @@ n_{R} = 4.43 (c) # J/psi \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$815 \pm 34$ & $743 \pm 44$ +$813 \pm 34$ & $734 \pm 44$ \end{tabular} # psi(2S) \begin{tabular}{c|c} $N_{Sig}$ & $N_{Bkg}$\\\hline -$69 \pm 12$ & $290 \pm 30$ +$69 \pm 11$ & $283 \pm 29$ \end{tabular} diff --git a/bdt_classification.h b/bdt_classification.h index d5b9c3d..6906c58 100644 --- a/bdt_classification.h +++ b/bdt_classification.h @@ -229,7 +229,7 @@ void TrainBDT(std::vector vars, const char* unique_id, TTree *sig_tree, TT data_loader->AddBackgroundTree(bkg_tree, background_weight); data_loader->PrepareTrainingAndTestTree("", "", "nTrain_Signal=0:nTrain_Background=0:SplitMode=Random:NormMode=NumEvents:V"); - factory->BookMethod(data_loader, TMVA::Types::kBDT, "BDT", "!H:!V:NTrees=600:MinNodeSize=2.5%:CreateMVAPdfs:MaxDepth=3:BoostType=AdaBoost:AdaBoostBeta=0.5:UseBaggedBoost:BaggedSampleFraction=0.5:SeparationType=GiniIndex:nCuts=20"); + factory->BookMethod(data_loader, TMVA::Types::kBDT, "BDT", "!H:!V:NTrees=400:MinNodeSize=2.5%:CreateMVAPdfs:MaxDepth=3:BoostType=AdaBoost:AdaBoostBeta=0.5:UseBaggedBoost:BaggedSampleFraction=0.5:SeparationType=GiniIndex:nCuts=20"); factory->TrainAllMethods(); factory->TestAllMethods(); diff --git a/new_analysis_b02hphmmumu.cpp b/new_analysis_b02hphmmumu.cpp index efd7798..55ab354 100644 --- a/new_analysis_b02hphmmumu.cpp +++ b/new_analysis_b02hphmmumu.cpp @@ -51,6 +51,7 @@ int new_analysis_b02hphmmumu() const char *sim_tree_name = "B0ToHpHmMuMu_noPID_mapped"; const char *end_state_mass_literal = "m(#pi^{+}#pi^{-}_{(#rightarrow K^{-})}#mu^{+}#mu^{-} & #pi^{+}_{(#rightarrow K^{+})}#pi^{-}#mu^{+}#mu^{-})"; const bool retrain_bdt = false; + const bool skip_fit = false; TChain *data_chain = new TChain(TString::Format("%s/DecayTree", data_tree_name)); data_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/data_samples/Collision23_Beam6800GeV-VeloClosed-MagDown-Excl-UT_RealData_Sprucing23r1_90000000_RD.root"); @@ -108,8 +109,10 @@ int new_analysis_b02hphmmumu() TV::Float("B0_PT", "B0_PT"), TV::Float("B0_BPVFDCHI2", "B0_BPVFDCHI2"), TV::Float("B0_BPVDIRA", "B0_BPVDIRA"), + TV::Double("B0_CHI2", "B0_CHI2"), TV::Float("Jpsi_BPVIPCHI2", "Jpsi_BPVIPCHI2"), TV::Float("Jpsi_PT", "Jpsi_PT"), + TV::Double("Jpsi_CHI2", "Jpsi_CHI2"), TV::Float("Kst0_BPVIPCHI2", "Kst0_BPVIPCHI2"), TV::Float("Kst0_PT", "Kst0_PT"), TV::Float("Kplus_BPVIPCHI2", "Kplus_BPVIPCHI2"), @@ -117,10 +120,10 @@ int new_analysis_b02hphmmumu() TV::Float("piminus_BPVIPCHI2", "piminus_BPVIPCHI2"), TV::Float("piminus_PT", "piminus_PT"), kplus_pid_k_var, - TV::Float("muminus_BPVIPCHI2", "muminus_BPVIPCHI2"), - // TV::Float("muminus_PT", "muminus_PT"), - TV::Float("muplus_BPVIPCHI2", "muplus_BPVIPCHI2"), - // TV::Float("muplus_PT", "muplus_PT"), + //TV::Float("muminus_BPVIPCHI2", "muminus_BPVIPCHI2"), + TV::Float("muminus_PT", "muminus_PT"), + //TV::Float("muplus_BPVIPCHI2", "muplus_BPVIPCHI2"), + TV::Float("muplus_PT", "muplus_PT"), }; TTree *sig_tree = new TTree("TreeS", "tree containing signal data"); @@ -147,7 +150,8 @@ int new_analysis_b02hphmmumu() if (std::all_of(vars.begin(), vars.end(), [](TV *v) { return v->IsDataFinite(); })) { - if (reconstructed_B_Mass > 5500. && ((TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) || (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.))) + if (reconstructed_B_Mass > 5500. && ((TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) || (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.)) + && kplus_pid_k_var->GetDataDouble() > -3 && (kplus_pid_k_var->GetDataDouble() - piminus_PID_K) > 0) { bkg_tree->Fill(); bkg_events++; @@ -194,6 +198,11 @@ int new_analysis_b02hphmmumu() std::cout << "# Finished BDT retrain." << std::endl; } + if (skip_fit) { + std::cout << "# Skipping evaluation of data." << std::endl; + return 0; + } + std::cout << "# Starting evaluation of data." << std::endl; Float_t *train_vars = new Float_t[vars.size()]; diff --git a/new_analysis_bu2hpmumu.cpp b/new_analysis_bu2hpmumu.cpp index 48b29cd..0ae6784 100644 --- a/new_analysis_bu2hpmumu.cpp +++ b/new_analysis_bu2hpmumu.cpp @@ -50,7 +50,8 @@ int new_analysis_bu2hpmumu() const char *data_tree_name = "SpruceRD_BuToHpMuMu"; const char *sim_tree_name = "BuToHpMuMu_noPID_mapped"; const char *end_state_mass_literal = "m(#pi^{+}_{(#rightarrow K^{+})}#mu^{+}#mu^{-})"; - const bool retrain_bdt = false; + const bool retrain_bdt = true; + const bool skip_fit = false; TChain *data_chain = new TChain(TString::Format("%s/DecayTree", data_tree_name)); data_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/data_samples/Collision23_Beam6800GeV-VeloClosed-MagDown-Excl-UT_RealData_Sprucing23r1_90000000_RD.root"); @@ -60,9 +61,9 @@ int new_analysis_bu2hpmumu() FourVect *l24v_data = FourVect::Init(data_chain, "muplus"); FourVect *hp4v_data = FourVect::Init(data_chain, "Kplus"); - // Double_t Kplus_PID_K; + Double_t Kplus_PID_K; - // data_chain->SetBranchAddress("Kplus_PID_K", &Kplus_PID_K); + data_chain->SetBranchAddress("Kplus_PID_K", &Kplus_PID_K); TChain *sim_chain = new TChain(TString::Format("%s/DecayTree", sim_tree_name)); sim_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/MC/BuToHpMuMu_mapped_mc.root"); @@ -100,17 +101,19 @@ int new_analysis_bu2hpmumu() auto hlt1_decision_histos = CreateHlt1DecisionHistos(analysis_name); std::map exclusive_hits{}; - TV* kplus_pid_k_var = TV::Double("Kplus_PID_K", "Kplus_PID_K"); + // TV* kplus_pid_k_var = TV::Double("Kplus_PID_K", "Kplus_PID_K"); std::vector vars{ TV::Float("B_PT", "B_PT"), TV::Float("B_BPVFDCHI2", "B_BPVFDCHI2"), TV::Float("B_BPVDIRA", "B_BPVDIRA"), + TV::Double("B_CHI2", "B_CHI2"), TV::Float("Jpsi_BPVIPCHI2", "Jpsi_BPVIPCHI2"), TV::Float("Jpsi_PT", "Jpsi_PT"), + TV::Double("Jpsi_CHI2", "Jpsi_CHI2"), TV::Float("Kplus_BPVIPCHI2", "Kplus_BPVIPCHI2"), TV::Float("Kplus_PT", "Kplus_PT"), - kplus_pid_k_var, + // kplus_pid_k_var, TV::Float("muminus_BPVIPCHI2", "muminus_BPVIPCHI2"), // TV::Float("muminus_PT", "muminus_PT"), TV::Float("muplus_BPVIPCHI2", "muplus_BPVIPCHI2"), @@ -142,7 +145,7 @@ int new_analysis_bu2hpmumu() if (std::all_of(vars.begin(), vars.end(), [](TV *v) { return v->IsDataFinite(); })) { - if (reconstructed_B_Mass > 5500. && ((TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) || (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.))) + if (reconstructed_B_Mass > 5500. && ((TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) || (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.)) && Kplus_PID_K > -3) { bkg_tree->Fill(); bkg_events++; @@ -186,6 +189,12 @@ int new_analysis_bu2hpmumu() std::cout << "# Finished BDT retrain." << std::endl; } + if (skip_fit) + { + std::cout << "# Skipping evaluation of data." << std::endl; + return 0; + } + std::cout << "# Starting evaluation of data." << std::endl; Float_t *train_vars = new Float_t[vars.size()]; @@ -222,23 +231,26 @@ int new_analysis_bu2hpmumu() FillHlt1DecisionHistos(hlt1_decision_histos, reconstructed_B_Mass); } - double mva_response = reader->EvaluateMVA("BDT"); - h1_bdt_probs->Fill(mva_response); + if (Kplus_PID_K > -3 && ((TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) || (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.))) + { + double mva_response = reader->EvaluateMVA("BDT"); + h1_bdt_probs->Fill(mva_response); - h1_B_Mass_unf->Fill(reconstructed_B_Mass); + h1_B_Mass_unf->Fill(reconstructed_B_Mass); - if (mva_response > mva_cut_value && kplus_pid_k_var->GetDataDouble() > -3) - { - h1_B_Mass_bdtf->Fill(reconstructed_B_Mass); - if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) - { - B_Mass_jpsi_var = reconstructed_B_Mass; - tree_B_Mass_jpsi->Fill(); - } - else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.) + if (mva_response > mva_cut_value) { - B_Mass_psi2s_var = reconstructed_B_Mass; - tree_B_Mass_psi2s->Fill(); + h1_B_Mass_bdtf->Fill(reconstructed_B_Mass); + if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) + { + B_Mass_jpsi_var = reconstructed_B_Mass; + tree_B_Mass_jpsi->Fill(); + } + else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.) + { + B_Mass_psi2s_var = reconstructed_B_Mass; + tree_B_Mass_psi2s->Fill(); + } } } @@ -290,36 +302,41 @@ int new_analysis_bu2hpmumu() res_file << "#### " << analysis_name << " @ " << std::put_time(&tm, "%c") << " ####" << std::endl; res_file << "J/Psi Mode: " << ErrToStr(roofit_hist_jpsi_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_jpsi_fitsum.background_yield, 0) << std::endl; res_file << " Sig/Bkg: " << ErrToStr(jpsi_sigobkg, 2) << std::endl; - res_file << "Psi(2S) Mode: " << ErrToStr(roofit_hist_psi2s_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_psi2s_fitsum.background_yield, 0) << std::endl; + res_file << "Psi(2S) Mode: " << ErrToStr(roofit_hist_psi2s_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_psi2s_fitsum.background_yield, 0) << std::endl; res_file << " Sig/Bkg: " << ErrToStr(psi2s_sigobkg, 2) << std::endl; res_file << "Mode Yield Ratio: " << ErrToStr(signal_ratio, 3) << std::endl; res_file << "Rel Br Frac MuMu: " << ErrToStr(mumu_br_frac, 3) << std::endl; auto br_frac = MultWithErr(signal_ratio.first, signal_ratio.second, mumu_br_frac.first, mumu_br_frac.second); - res_file << "Rel Br Frac: " << ErrToStr(br_frac, 3) << std::endl << std::endl; + res_file << "Rel Br Frac: " << ErrToStr(br_frac, 3) << std::endl + << std::endl; - res_file << "Params from Sim:" << std::endl << roofit_hist_sim.shape_parameters.ToString() << std::endl; + res_file << "Params from Sim:" << std::endl + << roofit_hist_sim.shape_parameters.ToString() << std::endl; for (const auto &par : roofit_hist_sim.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - res_file << std::endl << "Fitted Parameters: J/PSI" << std::endl << std::endl; + res_file << std::endl + << "Fitted Parameters: J/PSI" << std::endl + << std::endl; for (const auto &par : roofit_hist_jpsi_fitsum.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - res_file << std::endl << "Fitted Parameters: PSI(2S)" << std::endl << std::endl; + res_file << std::endl + << "Fitted Parameters: PSI(2S)" << std::endl + << std::endl; for (const auto &par : roofit_hist_psi2s_fitsum.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - auto print_table = [&res_file](std::string name, std::pair sig, std::pair bkg) { res_file << std::endl; diff --git a/new_analysis_bu2kpmumu.cpp b/new_analysis_bu2kpmumu.cpp index b034d6d..f6888a6 100644 --- a/new_analysis_bu2kpmumu.cpp +++ b/new_analysis_bu2kpmumu.cpp @@ -58,6 +58,10 @@ int new_analysis_bu2kpmumu() FourVect *l24v_data = FourVect::Init(data_chain, "muminus"); FourVect *hp4v_data = FourVect::Init(data_chain, "Kplus"); + Double_t Kplus_PID_K; + + data_chain->SetBranchAddress("Kplus_PID_K", &Kplus_PID_K); + TChain *sim_chain = new TChain(TString::Format("%s/DecayTree", sim_tree_name)); sim_chain->Add("/auto/data/pfeiffer/inclusive_detached_dilepton/MC/rd_btoxll_simulation_turbo_v0r0p6657752_BuToKpMuMu_12143001_magdown.root"); @@ -134,15 +138,18 @@ int new_analysis_bu2kpmumu() h1_B_Mass_unf->Fill(reconstructed_B_Mass); - if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) + if (Kplus_PID_K > -3) { - B_Mass_jpsi_var = reconstructed_B_Mass; - tree_B_Mass_jpsi->Fill(); - } - else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.) - { - B_Mass_psi2s_var = reconstructed_B_Mass; - tree_B_Mass_psi2s->Fill(); + if (TMath::Abs(dimuon.M() - JPSI_MASS) < 100.) + { + B_Mass_jpsi_var = reconstructed_B_Mass; + tree_B_Mass_jpsi->Fill(); + } + else if (TMath::Abs(dimuon.M() - PSI2S_MASS) < 100.) + { + B_Mass_psi2s_var = reconstructed_B_Mass; + tree_B_Mass_psi2s->Fill(); + } } PrintProgress(TString::Format("%s BDT Evaluation", analysis_name), data_entries, 10000, i); @@ -161,7 +168,7 @@ int new_analysis_bu2kpmumu() auto roofit_hist_jpsi_fitsum = CreateRooDataSetAndFitCB(tree_B_Mass_jpsi, B_Mass_jpsi_var_name, end_state_mass_literal, true, true, roofit_hist_sim.shape_parameters); roofit_hist_sim.shape_parameters.sigma_lr = roofit_hist_jpsi_fitsum.shape_parameters.sigma_lr; - + auto roofit_hist_psi2s_fitsum = CreateRooDataSetAndFitCB(tree_B_Mass_psi2s, B_Mass_psi2s_var_name, end_state_mass_literal, true, true, roofit_hist_sim.shape_parameters, true); DrawInDefaultCanvas(roofit_hist_jpsi_fitsum, analysis_name); @@ -185,36 +192,41 @@ int new_analysis_bu2kpmumu() res_file << "#### " << analysis_name << " @ " << std::put_time(&tm, "%c") << " ####" << std::endl; res_file << "J/Psi Mode: " << ErrToStr(roofit_hist_jpsi_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_jpsi_fitsum.background_yield, 0) << std::endl; res_file << " Sig/Bkg: " << ErrToStr(jpsi_sigobkg, 2) << std::endl; - res_file << "Psi(2S) Mode: " << ErrToStr(roofit_hist_psi2s_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_psi2s_fitsum.background_yield, 0) << std::endl; + res_file << "Psi(2S) Mode: " << ErrToStr(roofit_hist_psi2s_fitsum.signal_yield, 0) << " / " << ErrToStr(roofit_hist_psi2s_fitsum.background_yield, 0) << std::endl; res_file << " Sig/Bkg: " << ErrToStr(psi2s_sigobkg, 2) << std::endl; res_file << "Mode Yield Ratio: " << ErrToStr(signal_ratio, 3) << std::endl; res_file << "Rel Br Frac MuMu: " << ErrToStr(mumu_br_frac, 3) << std::endl; auto br_frac = MultWithErr(signal_ratio.first, signal_ratio.second, mumu_br_frac.first, mumu_br_frac.second); - res_file << "Rel Br Frac: " << ErrToStr(br_frac, 3) << std::endl << std::endl; + res_file << "Rel Br Frac: " << ErrToStr(br_frac, 3) << std::endl + << std::endl; - res_file << "Params from Sim:" << std::endl << roofit_hist_sim.shape_parameters.ToString() << std::endl; + res_file << "Params from Sim:" << std::endl + << roofit_hist_sim.shape_parameters.ToString() << std::endl; for (const auto &par : roofit_hist_sim.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - res_file << std::endl << "Fitted Parameters: J/PSI" << std::endl << std::endl; + res_file << std::endl + << "Fitted Parameters: J/PSI" << std::endl + << std::endl; for (const auto &par : roofit_hist_jpsi_fitsum.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - res_file << std::endl << "Fitted Parameters: PSI(2S)" << std::endl << std::endl; + res_file << std::endl + << "Fitted Parameters: PSI(2S)" << std::endl + << std::endl; for (const auto &par : roofit_hist_psi2s_fitsum.fitted_params) { res_file << par.ToString(true).c_str() << std::endl; } - auto print_table = [&res_file](std::string name, std::pair sig, std::pair bkg) { res_file << std::endl;