void sep_macro() { // s t y l e s // ----------- TFile *file_phi = TFile::Open("kpkm_final3.root"); TFile *file_rho = TFile::Open("pippim_final1.root"); TFile *spe = TFile::Open("spe.root", "recreate"); //TFile fRootFile("kpkm_final3","recreate"); double sp, sep_err, xavg, yavg, pavg; TTree *tree = new TTree("tree", "tree"); tree->Branch("sp", &sp, "sp/D"); tree->Branch("sep_err", &sep_err, "sep_err/D"); //TTree tree("stuff","stuff"); //tree.Branch("sp", &sp, "sp/D"); //tree.Branch("sep_err", &sep_err, "sep_err/D"); //tree.Branch("xavg", &xavg, "xavg/D"); //tree.Branch("yavg", &yavg, "yavg/D"); //tree.Branch("pavg", &pavg, "pavg/D"); // gluex_styles auto gluex_style = new TStyle("GlueX","Default GlueX Style"); gluex_style->SetCanvasBorderMode(0); gluex_style->SetPadBorderMode(0); gluex_style->SetPadColor(0); gluex_style->SetCanvasColor(0); gluex_style->SetTitleColor(0); gluex_style->SetStatColor(0); // some default window sizes gluex_style->SetCanvasDefW(800); gluex_style->SetCanvasDefH(800); // let's change the default margins gluex_style->SetPadBottomMargin(0.15); gluex_style->SetPadLeftMargin(0.15); gluex_style->SetPadTopMargin(0.15); gluex_style->SetPadRightMargin(0.08); // axis labels and settings gluex_style->SetLabelSize(0.045,"xyz"); // size of axis value font gluex_style->SetTitleSize(0.05,"xyz"); // size of axis title font gluex_style->SetTitleFont(42,"xyz"); // font option gluex_style->SetTitleColor(kBlack,"xyz"); gluex_style->SetLabelFont(42,"xyz"); gluex_style->SetTitleOffset(1.55,"y"); gluex_style->SetLabelOffset(0.01,"xyz"); // stop collisions of "0"s at the origin // histogram settings gluex_style->SetOptStat(0); // no stats box by default gluex_style->SetOptTitle(0); // no title by default gluex_style->SetHistLineWidth(1); gluex_style->cd(); gROOT->ForceStyle(); // get root files const Int_t nhists = 12; // choose histograms to plot char const *hnames_phi[nhists] = { "Ldiff_KPlus_p1_y1_x1", "Ldiff_KPlus_p1_y1_x2", "Ldiff_KPlus_p1_y1_x3", "Ldiff_KPlus_p1_y2_x1", "Ldiff_KPlus_p1_y2_x2", "Ldiff_KPlus_p1_y2_x3", "Ldiff_KPlus_p1_y3_x1", "Ldiff_KPlus_p1_y3_x2", "Ldiff_KPlus_p1_y3_x3", "Ldiff_KPlus_p1_y4_x1", "Ldiff_KPlus_p1_y4_x2", "Ldiff_KPlus_p1_y4_x3" // Phi(P1) //"Ldiff_KPlus_p2_y1_x1", "Ldiff_KPlus_p2_y1_x2", "Ldiff_KPlus_p2_y1_x3", //"Ldiff_KPlus_p2_y2_x1", "Ldiff_KPlus_p2_y2_x2", "Ldiff_KPlus_p2_y2_x3", //"Ldiff_KPlus_p2_y3_x1", "Ldiff_KPlus_p2_y3_x2", "Ldiff_KPlus_p2_y3_x3", //"Ldiff_KPlus_p2_y4_x1", "Ldiff_KPlus_p2_y4_x2", "Ldiff_KPlus_p2_y4_x3", // Phi(P2) //"Ldiff_KPlus_p3_y1_x1", "Ldiff_KPlus_p3_y1_x2", "Ldiff_KPlus_p3_y1_x3", //"Ldiff_KPlus_p3_y2_x1", "Ldiff_KPlus_p3_y2_x2", "Ldiff_KPlus_p3_y2_x3", //"Ldiff_KPlus_p3_y3_x1", "Ldiff_KPlus_p3_y3_x2", "Ldiff_KPlus_p3_y3_x3", //"Ldiff_KPlus_p3_y4_x1", "Ldiff_KPlus_p3_y4_x2", "Ldiff_KPlus_p3_y4_x3" // Phi(P3) }; char const *hnames_rho[nhists] = { "Cut_PiPlus_p1_y1_x1", "Cut_PiPlus_p1_y1_x2", "Cut_PiPlus_p1_y1_x3", "Cut_PiPlus_p1_y2_x1", "Cut_PiPlus_p1_y2_x2", "Cut_PiPlus_p1_y2_x3", "Cut_PiPlus_p1_y3_x1", "Cut_PiPlus_p1_y3_x2", "Cut_PiPlus_p1_y3_x3", "Cut_PiPlus_p1_y4_x1", "Cut_PiPlus_p1_y4_x2", "Cut_PiPlus_p1_y4_x3" // Rho(P1) //"Cut_PiPlus_p2_y1_x1", "Cut_PiPlus_p2_y1_x2", "Cut_PiPlus_p2_y1_x3", //"Cut_PiPlus_p2_y2_x1", "Cut_PiPlus_p2_y2_x2", "Cut_PiPlus_p2_y2_x3", //"Cut_PiPlus_p2_y3_x1", "Cut_PiPlus_p2_y3_x2", "Cut_PiPlus_p2_y3_x3", //"Cut_PiPlus_p2_y4_x1", "Cut_PiPlus_p2_y4_x2", "Cut_PiPlus_p2_y4_x3", // Rho(P2) //"Cut_PiPlus_p3_y1_x1", "Cut_PiPlus_p3_y1_x2", "Cut_PiPlus_p3_y1_x3", //"Cut_PiPlus_p3_y2_x1", "Cut_PiPlus_p3_y2_x2", "Cut_PiPlus_p3_y2_x3", //"Cut_PiPlus_p3_y3_x1", "Cut_PiPlus_p3_y3_x2", "Cut_PiPlus_p3_y3_x3", //"Cut_PiPlus_p3_y4_x1", "Cut_PiPlus_p3_y4_x2", "Cut_PiPlus_p3_y4_x3" // Rho(P3) }; TH1D* hist_Ldiff_phi_kp[nhists]; TH1F* hist_Ldiff_rho_pip[nhists]; TCanvas* canvas1 = new TCanvas("canv1","histo", 500, 600,500,800); canvas1->Divide(3,4); double scale[nhists]; double y_max[nhists]; double x_avg[nhists]; double y_avg[nhists]; for(int i=0; i< nhists; i++) { hist_Ldiff_phi_kp[i] = (TH1D*) file_phi->Get(hnames_phi[i]); hist_Ldiff_rho_pip[i] = (TH1F*) file_rho->Get(hnames_rho[i]); canvas1->cd(i+1); hist_Ldiff_phi_kp[i]->SetLineColor(2); hist_Ldiff_rho_pip[i]->SetLineColor(4); scale[i] = (hist_Ldiff_rho_pip[i]->GetMaximum()) / (hist_Ldiff_phi_kp[i]->GetMaximum()); y_max[i] = hist_Ldiff_rho_pip[i]->GetMaximum() + 10000; // define Fits TF1 *g1 = new TF1("m1", "gaus", -200, 200); g1->SetLineColor(9); TF1 *g2 = new TF1("m2", "gaus", -200, 200); g2->SetLineColor(46); hist_Ldiff_rho_pip[i]->Draw(); hist_Ldiff_rho_pip[i]->Fit("m1", "R"); hist_Ldiff_phi_kp[i]->Scale(scale[i]); hist_Ldiff_phi_kp[i]->Draw("same"); hist_Ldiff_phi_kp[i]->Fit("m2", "R+"); // plot fit for kaon ldiff after the scaling hist_Ldiff_rho_pip[i]->GetYaxis()->SetRangeUser(0,y_max[i]); // get mean and sigma values of each gaussian double mean_pip = g1->GetParameter(1); double mean_kp = g2->GetParameter(1); double sigma_pip = g1->GetParameter(2); double sigma_kp = g2->GetParameter(2); double dmean_pip = g1->GetParError(1); double dmean_kp = g2->GetParError(1); double dsigma_pip = g1->GetParError(2); double dsigma_kp = g2->GetParError(2); // calculate seperation power: sp = |mean_pi - mean_K| / 0.5*(sigma_pi + sigma_K) double sp = abs(mean_pip - mean_kp) / (0.5*(sigma_pip + sigma_kp)); // calculate seperation power error: double e1,e2,e3,e4; e1 = 2/(sigma_pip + sigma_kp)*dmean_pip; e2 = -2/(sigma_pip + sigma_kp)*dmean_kp; e3 = -2*(mean_pip - mean_kp)/((sigma_pip + sigma_kp)*(sigma_pip + sigma_kp))*dsigma_pip; e4 = -2*(mean_pip - mean_kp)/((sigma_pip + sigma_kp)*(sigma_pip + sigma_kp))*dsigma_kp; double sep_err = sqrt(e1*e1+e2*e2+e3*e3+e4*e4); // find averages pos(x,y) and momentum // xavg = // fill tree tree->Fill(); spe->Write(); //tree.Fill(); // show sep TPaveText* sep = new TPaveText(0.75, 0.65, 0.9, 0.75, "BRNDC"); sep -> AddText(Form("sep power = %.2f #pm %.2f s.d.", sp, sep_err)); sep -> SetTextSize(0.1); sep -> SetTextColor(kBlack); sep -> SetShadowColor(0); sep -> SetFillColor(0); sep -> SetLineColor(0); sep->Draw("same"); } }