pfeiffer/inclusive_detached_dilepton
1
0
Fork 0
Code Issues Pull Requests 1 repo.project_board Releases Wiki Activity

some fitting adjustments

Browse Source
This commit is contained in:
Marius Pfeiffer 2023-11-17 13:18:07 +01:00
parent 5af73a878a
commit aa89e22733
2 changed files with 66 additions and 26 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
inspect_tree.cpp
View File

@ -24,6 +24,30 @@ int inspect_tree() {
auto tree = file->Get<TTree>("BuToKpMuMu23_noPID/DecayTree"); auto tree = file->Get<TTree>("BuToKpMuMu23_noPID/DecayTree");
auto branches = tree->GetListOfBranches(); auto branches = tree->GetListOfBranches();
// /*Lesen*/
// Double_t Bu_M;
// tree->SetBranchAddress("Bu_M", &Bu_M);
// for (size_t i = 0; i < tree->GetEntries(); i++)
// {
// tree->GetEntry(i);
// // Bu_M verwenden
// }
// /*Schreiben*/
// Double_t B_Theta, B_Blubb;
// tree->Branch("B_Theta", &B_Theta, "B_Theta/D");
// tree->Branch("B_Blubb", &B_Blubb, "B_Blubb/D");
// for (size_t i = 0; i < 1000; i++)
// {
// B_Theta = 345435; // ....
// B_Blubb = 348535435435579; // ....
// tree->Fill();
// }
for (size_t i = 0; i < branches->GetSize(); i++) for (size_t i = 0; i < branches->GetSize(); i++)
{ {
TBranch* b = (TBranch*) branches->At(i); TBranch* b = (TBranch*) branches->At(i);

68
train_bdt.cpp
View File

@ -34,7 +34,7 @@
const bool TRAIN_TMVA = false; const bool TRAIN_TMVA = false;
const bool EVALUATE_TMVA = true; const bool EVALUATE_TMVA = true;
const int N_BINS = 150; const int N_BINS = 130;
const double B_PLUS_MASS = 5279.; const double B_PLUS_MASS = 5279.;
const double J_PSI_MASS = 3096.; const double J_PSI_MASS = 3096.;
@ -261,9 +261,9 @@ int train_bdt()
// files to load // files to load
std::vector<std::string> mc_filenames = std::vector<std::string> mc_filenames =
{ {
"/auto/data/pfeiffer/inclusive_detached_dilepton/MC/BuToKpMuMu_rd_btoxll_simulation_12143001_MagDown_v0r0p6316987.root"}; "/auto/data/pfeiffer/inclusive_detached_dilepton/MC/BuToKpMuMu_rd_btoxll_simulation_12143001_MagDown_v0r0p6316365_FULLSTREAM.root"};
TChain *mc_chain = new TChain("BuToKpMuMu23_noPID/DecayTree"); TChain *mc_chain = new TChain("BuToKpMuMu_noPID/DecayTree");
for (unsigned int i = 0; i < mc_filenames.size(); i++) for (unsigned int i = 0; i < mc_filenames.size(); i++)
{ {
mc_chain->Add(mc_filenames.at(i).c_str()); mc_chain->Add(mc_filenames.at(i).c_str());
@ -272,7 +272,9 @@ int train_bdt()
std::vector<TV> vars{ std::vector<TV> vars{
TV::Float("Bplus_PT", "B_PT"), TV::Float("Bplus_PT", "B_PT"),
TV::Float("Bplus_BPVFDCHI2", "B_BPVFDCHI2"), TV::Float("Bplus_BPVFDCHI2", "B_BPVFDCHI2"),
TV::Float("Bplus_BPVDIRA", "B_BPVDIRA"),
TV::Float("Jpsi_BPVIPCHI2", "Jpsi_BPVIPCHI2"), TV::Float("Jpsi_BPVIPCHI2", "Jpsi_BPVIPCHI2"),
TV::Float("Jpsi_BPVDIRA", "Jpsi_BPVDIRA"),
TV::Float("Jpsi_PT", "Jpsi_PT"), TV::Float("Jpsi_PT", "Jpsi_PT"),
TV::Float("Kplus_BPVIPCHI2", "K_BPVIPCHI2"), TV::Float("Kplus_BPVIPCHI2", "K_BPVIPCHI2"),
TV::Float("Kplus_PT", "K_PT"), TV::Float("Kplus_PT", "K_PT"),
@ -378,7 +380,7 @@ int train_bdt()
sig_tree->Fill(); sig_tree->Fill();
added_sig_entries++; added_sig_entries++;
if (added_sig_entries >= added_bck_entries) { if (added_sig_entries >= added_bck_entries * 2) {
break; break;
} }
} }
@ -414,7 +416,7 @@ int train_bdt()
data_loader->AddBackgroundTree(bkg_tree, background_weight); data_loader->AddBackgroundTree(bkg_tree, background_weight);
data_loader->PrepareTrainingAndTestTree("", "", "nTrain_Signal=0:nTrain_Background=0:SplitMode=Random:NormMode=NumEvents:!V"); 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=500: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=600:MinNodeSize=2.5%:CreateMVAPdfs:MaxDepth=3:BoostType=AdaBoost:AdaBoostBeta=0.5:UseBaggedBoost:BaggedSampleFraction=0.5:SeparationType=GiniIndex:nCuts=20");
factory->TrainAllMethods(); factory->TrainAllMethods();
factory->TestAllMethods(); factory->TestAllMethods();
@ -477,7 +479,7 @@ int train_bdt()
double mva_response = reader->EvaluateMVA("BDT"); double mva_response = reader->EvaluateMVA("BDT");
h1_probs->Fill(mva_response); h1_probs->Fill(mva_response);
const double mva_cut_value = 0.01; const double mva_cut_value = 0.09; // -0.02;
if (mva_response > mva_cut_value) if (mva_response > mva_cut_value)
{ {
@ -496,12 +498,12 @@ int train_bdt()
c1->Draw(); c1->Draw();
// TCanvas *c3 = new TCanvas("c3", "Canvas 3", 0, 0, 1000, 600); TCanvas *c3 = new TCanvas("c3", "Canvas 3", 0, 0, 1000, 600);
// auto fitFrame = CreateRooFit(h1_Bplus_M); auto fitFrame = CreateRooFit(h1_Bplus_M);
// fitFrame->Draw(); fitFrame->Draw();
// c3->Draw(); c3->Draw();
} }
else else
{ {
@ -514,31 +516,45 @@ int train_bdt()
RooPlot *CreateRooFit(TH1D *hist) RooPlot *CreateRooFit(TH1D *hist)
{ {
RooRealVar roo_var_mass("roo_var_mass", "B+ Mass Variable", 4700., 6500.); RooRealVar roo_var_mass("roo_var_mass", "B+ Mass Variable", 4700., 6500.);
roo_var_mass.setRange("fitting_range", 5000., 6000.); roo_var_mass.setRange("fitting_range", 4800., 5800.);
roo_var_mass.setRange("plot_range", 4700., 6500.); roo_var_mass.setRange("plot_range", 4700., 6500.);
TString hist_name = "roohist_bplus_M"; TString hist_name = "roohist_bplus_M";
RooDataHist roohist_bplus_M(hist_name, "B Plus Mass Histogram", roo_var_mass, RooFit::Import(*hist)); RooDataHist roohist_bplus_M(hist_name, "B Plus Mass Histogram", roo_var_mass, RooFit::Import(*hist));
// RooRealVar roo_sig_gauss_mean("roo_sig_gauss_mean", "Mass Gauss Mean", 5250., 5100., 5400.); RooRealVar roo_sig_gauss_mean("roo_sig_gauss_mean", "Mass Gauss Mean", 5250., 5100., 5400.);
// RooRealVar roo_sig_gauss_sigma("roo_sig_gauss_sigma", "Mass Gauss Sigma", 60., 0., 150.); RooRealVar roo_sig_gauss_sigma("roo_sig_gauss_sigma", "Mass Gauss Sigma", 30., 20., 40.);
// RooGaussian roo_sig_gauss("roo_sig_gauss", "B+ Mass Signal Gaussian", roo_var_mass, roo_sig_gauss_mean, roo_sig_gauss_sigma); RooGaussian roo_sig_gauss("roo_sig_gauss", "B+ Mass Signal Gaussian", roo_var_mass, roo_sig_gauss_mean, roo_sig_gauss_sigma);
// Crystal Ball for Signal // Crystal Ball for Signal
RooRealVar roo_sig_cb_x0("roo_sig_cry_x0", "Location", B_PLUS_MASS, 5100., 5400.); // RooRealVar roo_sig_cb_x0("roo_sig_cry_x0", "Location", B_PLUS_MASS, 5100., 5400.);
RooRealVar roo_sig_cb_sigmaL("roo_sig_cry_sigmaL", "Sigma L", 30., 0., 60.); // RooRealVar roo_sig_cb_sigmaL("roo_sig_cry_sigmaL", "Sigma L", 30., 0., 60.);
RooRealVar roo_sig_cb_sigmaR("roo_sig_cry_sigmaR", "Sigma R", 30., 0., 60.); // RooRealVar roo_sig_cb_sigmaR("roo_sig_cry_sigmaR", "Sigma R", 30., 0., 60.);
RooRealVar roo_sig_cb_alphaL("roo_sig_cry_alphaL", "Alpha L", 15., 0., 30.); // RooRealVar roo_sig_cb_alphaL("roo_sig_cry_alphaL", "Alpha L", 15., 0., 30.);
RooRealVar roo_sig_cb_nL("roo_sig_cry_nL", "Exponent L", 10., 1., 20.); // RooRealVar roo_sig_cb_nL("roo_sig_cry_nL", "Exponent L", 0., -40., 40.);
RooRealVar roo_sig_cb_alphaR("roo_sig_cry_alphaR", "Alpha R", 15., 0., 30.); // RooRealVar roo_sig_cb_alphaR("roo_sig_cry_alphaR", "Alpha R", 15., 0., 30.);
RooRealVar roo_sig_cb_nR("roo_sig_cry_nR", "Exponent R", 10., 1., 20.); // RooRealVar roo_sig_cb_nR("roo_sig_cry_nR", "Exponent R", 0., -40., 40.);
RooCrystalBall roo_sig_cb("roo_sig_cb", "Signal Crystal Ball", roo_var_mass, roo_sig_cb_x0, roo_sig_cb_sigmaL, roo_sig_cb_sigmaR, roo_sig_cb_alphaL, roo_sig_cb_nL, roo_sig_cb_alphaR, roo_sig_cb_nR); // RooCrystalBall roo_sig_cb("roo_sig_cb", "Signal Crystal Ball", roo_var_mass, roo_sig_cb_x0, roo_sig_cb_sigmaL, roo_sig_cb_sigmaR, roo_sig_cb_alphaL, roo_sig_cb_nL, roo_sig_cb_alphaR, roo_sig_cb_nR);
RooRealVar roo_bkg_exp_c("roo_bkg_exp_c", "Background C", -0.00099021, -0.003, -0.0001); // Double Gauss Signal
// RooRealVar roo_sig_gauss1_mean("roo_sig_gauss1_mean", "Mass Gauss 1 Mean", 5250., 5200., 5300.);
// RooRealVar roo_sig_gauss1_sigma("roo_sig_gauss1_sigma", "Mass Gauss 1 Sigma", 60., 0., 150.);
// RooGaussian roo_sig_gauss1("roo_sig_gauss1", "B+ Mass Signal 1 Gaussian", roo_var_mass, roo_sig_gauss1_mean, roo_sig_gauss1_sigma);
// RooRealVar roo_sig_gauss2_mean("roo_sig_gauss2_mean", "Mass Gauss 2 Mean", 5250., 5200., 5300.);
// RooRealVar roo_sig_gauss2_sigma("roo_sig_gauss2_sigma", "Mass Gauss 2 Sigma", 60., 0., 150.);
// RooGaussian roo_sig_gauss2("roo_sig_gauss2", "B+ Mass Signal 2 Gaussian", roo_var_mass, roo_sig_gauss2_mean, roo_sig_gauss2_sigma);
// RooRealVar roo_sig_double_gauss_frac("roo_sig_double_gauss_frac", "B+ Mass Signal Double Gauss Frac", 0.5);
// RooAddPdf roo_sig_double_gauss("roo_sig_double_gauss", "B+ Mass Signal Double Gauss", roo_sig_gauss1, roo_sig_gauss2, roo_sig_double_gauss_frac);
RooRealVar roo_bkg_exp_c("roo_bkg_exp_c", "Background C", -0.000693147, -0.002, 0.);
RooExponential roo_bkg_exp("roo_bkg_exp", "B+ Mass Background Exp", roo_var_mass, roo_bkg_exp_c); RooExponential roo_bkg_exp("roo_bkg_exp", "B+ Mass Background Exp", roo_var_mass, roo_bkg_exp_c);
RooRealVar roo_var_mass_nsig("roo_var_mass_nsig", "B+ Mass N Signal", 0., hist->GetEntries()); RooRealVar roo_var_mass_nsig("roo_var_mass_nsig", "B+ Mass N Signal", 0., hist->GetEntries());
@ -546,11 +562,11 @@ RooPlot *CreateRooFit(TH1D *hist)
TString pdf_name = "roo_pdf_sig_plus_bkg"; TString pdf_name = "roo_pdf_sig_plus_bkg";
RooAddPdf roo_pdf_sig_plus_bkg(pdf_name, "Sig + Bkg PDF", RooAddPdf roo_pdf_sig_plus_bkg(pdf_name, "Sig + Bkg PDF",
RooArgList(roo_sig_cb, roo_bkg_exp), RooArgList(roo_sig_gauss, roo_bkg_exp),
RooArgList(roo_var_mass_nsig, roo_var_mass_nbkg)); RooArgList(roo_var_mass_nsig, roo_var_mass_nbkg));
RooPlot *roo_frame_mass = roo_var_mass.frame(); RooPlot *roo_frame_mass = roo_var_mass.frame();
roohist_bplus_M.plotOn(roo_frame_mass, RooFit::Binning(N_BINS), RooFit::Name(hist_name), RooFit::MarkerColor(18), RooFit::Range("plot_range")); roohist_bplus_M.plotOn(roo_frame_mass, RooFit::Binning(N_BINS), RooFit::Name(hist_name), RooFit::MarkerColor(15), RooFit::Range("plot_range"));
RooFitResult *fitres = roo_pdf_sig_plus_bkg.fitTo(roohist_bplus_M, RooFit::Save(), RooFit::PrintLevel(1), RooFit::Range("fitting_range")); RooFitResult *fitres = roo_pdf_sig_plus_bkg.fitTo(roohist_bplus_M, RooFit::Save(), RooFit::PrintLevel(1), RooFit::Range("fitting_range"));
@ -561,7 +577,7 @@ RooPlot *CreateRooFit(TH1D *hist)
// roo_sig_cb.plotOn(roo_frame_mass, RooFit::LineColor(kAlpine), RooFit::LineStyle(kDashed), RooFit::Range("fitting_range")); // roo_sig_cb.plotOn(roo_frame_mass, RooFit::LineColor(kAlpine), RooFit::LineStyle(kDashed), RooFit::Range("fitting_range"));
// roo_bkg_exp.plotOn(roo_frame_mass, RooFit::LineColor(kOrange), RooFit::LineStyle(kDashed), RooFit::Range("fitting_range")); // roo_bkg_exp.plotOn(roo_frame_mass, RooFit::LineColor(kOrange), RooFit::LineStyle(kDashed), RooFit::Range("fitting_range"));
roo_sig_cb.paramOn(roo_frame_mass, RooFit::Layout(0.60, 0.99, 0.90)); roo_sig_gauss.paramOn(roo_frame_mass, RooFit::Layout(0.60, 0.99, 0.90));
roo_frame_mass->getAttText()->SetTextSize(0.027); roo_frame_mass->getAttText()->SetTextSize(0.027);
double ymax = roo_frame_mass->GetMaximum(); double ymax = roo_frame_mass->GetMaximum();