RandomConeAna/offline/ConeReso.C

0001
0002 #include <iostream>
0003 #include <vector>
0004
0005 #include "sPhenixStyle.h"
0006 #include "sPhenixStyle.C"
0007
0008 int ConeReso()
0009 {
0010
0011     // set sPHENIX style
0012     SetsPhenixStyle();
0013     TH1::SetDefaultSumw2();
0014     TH2::SetDefaultSumw2();
0015
0016
0017     // in/out file names
0018     TString input_file = "/sphenix/user/tmengel/JetPerformancePPG/RandomConeAna/condor/merging/rootfiles/HIJING.root";
0019     TString output_file = "HIJING_ConePlots.root";
0020
0021
0022     // open input file
0023     TFile *f = new TFile(input_file, "READ");
0024     if(!f->IsOpen() || f->IsZombie()){ std::cout << "File " << input_file << " is zombie" << std::endl;  return -1; }
0025
0026     // get tree
0027     TTree *t = (TTree*)f->Get("event_info");
0028     if(!t){ std::cout << "Tree event_info not found in " << input_file << std::endl; return -1; }
0029
0030     // input variables
0031     int m_event_id = 0;
0032     int m_centrality = 0;
0033     float m_impactparam = 0;
0034     float m_zvtx = 0;
0035     double m_mbd_NS = 0;
0036
0037     // rho variables
0038     float m_TowerRho_AREA_rho = 0;
0039     float m_TowerRho_AREA_sigma = 0;
0040     float m_TowerRho_MULT_rho = 0;
0041     float m_TowerRho_MULT_sigma = 0;
0042
0043     // random cone variables
0044     float m_RandomCone_Tower_basic_r04_R = 0;
0045     float m_RandomCone_Tower_basic_r04_eta = 0;
0046     float m_RandomCone_Tower_basic_r04_phi = 0;
0047     float m_RandomCone_Tower_basic_r04_pt = 0;
0048     int m_RandomCone_Tower_basic_r04_nclustered = 0;
0049
0050     float m_RandomCone_Tower_randomize_etaphi_r04_R = 0;
0051     float m_RandomCone_Tower_randomize_etaphi_r04_eta = 0;
0052     float m_RandomCone_Tower_randomize_etaphi_r04_phi = 0;
0053     float m_RandomCone_Tower_randomize_etaphi_r04_pt = 0;
0054     int m_RandomCone_Tower_randomize_etaphi_r04_nclustered = 0;
0055
0056     float m_RandomCone_Tower_Sub1_basic_r04_R = 0;
0057     float m_RandomCone_Tower_Sub1_basic_r04_eta = 0;
0058     float m_RandomCone_Tower_Sub1_basic_r04_phi = 0;
0059     float m_RandomCone_Tower_Sub1_basic_r04_pt = 0;
0060     int m_RandomCone_Tower_Sub1_basic_r04_nclustered = 0;
0061
0062     float m_RandomCone_Tower_Sub1_randomize_etaphi_r04_R = 0;
0063     float m_RandomCone_Tower_Sub1_randomize_etaphi_r04_eta = 0;
0064     float m_RandomCone_Tower_Sub1_randomize_etaphi_r04_phi = 0;
0065     float m_RandomCone_Tower_Sub1_randomize_etaphi_r04_pt = 0;
0066     int m_RandomCone_Tower_Sub1_randomize_etaphi_r04_nclustered = 0;
0067
0068
0069     // set branch addresses
0070     t->SetBranchAddress("event_id", &m_event_id);
0071     t->SetBranchAddress("centrality", &m_centrality);
0072     t->SetBranchAddress("impactparam", &m_impactparam);
0073     t->SetBranchAddress("zvtx", &m_zvtx);
0074     t->SetBranchAddress("mbd_NS", &m_mbd_NS);
0075     t->SetBranchAddress("TowerRho_AREA_rho", &m_TowerRho_AREA_rho);
0076     t->SetBranchAddress("TowerRho_AREA_sigma", &m_TowerRho_AREA_sigma);
0077     t->SetBranchAddress("TowerRho_MULT_rho", &m_TowerRho_MULT_rho);
0078     t->SetBranchAddress("TowerRho_MULT_sigma", &m_TowerRho_MULT_sigma);
0079     t->SetBranchAddress("RandomCone_Tower_basic_r04_R", &m_RandomCone_Tower_basic_r04_R);
0080     t->SetBranchAddress("RandomCone_Tower_basic_r04_eta", &m_RandomCone_Tower_basic_r04_eta);
0081     t->SetBranchAddress("RandomCone_Tower_basic_r04_phi", &m_RandomCone_Tower_basic_r04_phi);
0082     t->SetBranchAddress("RandomCone_Tower_basic_r04_pt", &m_RandomCone_Tower_basic_r04_pt);
0083     t->SetBranchAddress("RandomCone_Tower_basic_r04_nclustered", &m_RandomCone_Tower_basic_r04_nclustered);
0084     t->SetBranchAddress("RandomCone_Tower_randomize_etaphi_r04_R", &m_RandomCone_Tower_randomize_etaphi_r04_R);
0085     t->SetBranchAddress("RandomCone_Tower_randomize_etaphi_r04_eta", &m_RandomCone_Tower_randomize_etaphi_r04_eta);
0086     t->SetBranchAddress("RandomCone_Tower_randomize_etaphi_r04_phi", &m_RandomCone_Tower_randomize_etaphi_r04_phi);
0087     t->SetBranchAddress("RandomCone_Tower_randomize_etaphi_r04_pt", &m_RandomCone_Tower_randomize_etaphi_r04_pt);
0088     t->SetBranchAddress("RandomCone_Tower_randomize_etaphi_r04_nclustered", &m_RandomCone_Tower_randomize_etaphi_r04_nclustered);
0089     t->SetBranchAddress("RandomCone_Tower_Sub1_basic_r04_R", &m_RandomCone_Tower_Sub1_basic_r04_R);
0090     t->SetBranchAddress("RandomCone_Tower_Sub1_basic_r04_eta", &m_RandomCone_Tower_Sub1_basic_r04_eta);
0091     t->SetBranchAddress("RandomCone_Tower_Sub1_basic_r04_phi", &m_RandomCone_Tower_Sub1_basic_r04_phi);
0092     t->SetBranchAddress("RandomCone_Tower_Sub1_basic_r04_pt", &m_RandomCone_Tower_Sub1_basic_r04_pt);
0093     t->SetBranchAddress("RandomCone_Tower_Sub1_basic_r04_nclustered", &m_RandomCone_Tower_Sub1_basic_r04_nclustered);
0094     t->SetBranchAddress("RandomCone_Tower_Sub1_randomize_etaphi_r04_R", &m_RandomCone_Tower_Sub1_randomize_etaphi_r04_R);
0095     t->SetBranchAddress("RandomCone_Tower_Sub1_randomize_etaphi_r04_eta", &m_RandomCone_Tower_Sub1_randomize_etaphi_r04_eta);
0096     t->SetBranchAddress("RandomCone_Tower_Sub1_randomize_etaphi_r04_phi", &m_RandomCone_Tower_Sub1_randomize_etaphi_r04_phi);
0097     t->SetBranchAddress("RandomCone_Tower_Sub1_randomize_etaphi_r04_pt", &m_RandomCone_Tower_Sub1_randomize_etaphi_r04_pt);
0098     t->SetBranchAddress("RandomCone_Tower_Sub1_randomize_etaphi_r04_nclustered", &m_RandomCone_Tower_Sub1_randomize_etaphi_r04_nclustered);
0099
0100     // get number of entries
0101     int nentries = t->GetEntries();
0102
0103     // declare histograms
0104     const int resp_N = 250;
0105     double max_resp = 60;
0106     double min_resp = -60;
0107     double delta_resp = (max_resp - min_resp)/resp_N;
0108     double resp_bins[resp_N+1];
0109     for(int i = 0; i < resp_N+1; i++)
0110     {
0111         resp_bins[i] = min_resp + i*delta_resp;
0112     }
0113
0114     const double cent_bins[] = {0, 10, 20, 30, 40, 50, 60, 80, 100};
0115     const int cent_N = sizeof(cent_bins)/sizeof(double) - 1;
0116
0117     const int rho_N = 250;
0118     double rho_max = 125;
0119     double rho_min = 0;
0120     double delta_rho = (rho_max - rho_min)/rho_N;
0121     double rho_bins[rho_N+1];
0122     for(int i = 0; i < rho_N+1; i++)
0123     {
0124         rho_bins[i] = rho_min + i*delta_rho;
0125     }
0126
0127     const int rho_mult_N = 200;
0128     double rho_mult_max = 2.0;
0129     double rho_mult_min = 0;
0130     double delta_rho_mult = (rho_mult_max - rho_mult_min)/rho_mult_N;
0131     double rho_bins_mult[rho_mult_N+1];
0132     for(int i = 0; i < rho_mult_N+1; i++)
0133     {
0134         rho_bins_mult[i] = rho_mult_min + i*delta_rho_mult;
0135     }
0136
0137     const int nclustered_N = 150;
0138     double nclustered_max = 150;
0139     double nclustered_min = 0;
0140     double delta_nclustered = (nclustered_max - nclustered_min)/nclustered_N;
0141     double nclustered_bins[nclustered_N+1];
0142     for(int i = 0; i < nclustered_N+1; i++)
0143     {
0144         nclustered_bins[i] = nclustered_min + i*delta_nclustered;
0145     }
0146
0147
0148
0149     // cone residual histograms
0150     TH2D * h2_cone_residual_pt_basic_area = new TH2D("h2_cone_residual_pt_basic_area", "h2_cone_residual_pt_basic_area", resp_N, resp_bins, cent_N, cent_bins);
0151     TH2D * h2_cone_residual_pt_basic_mult = new TH2D("h2_cone_residual_pt_basic_mult", "h2_cone_residual_pt_basic_mult", resp_N, resp_bins, cent_N, cent_bins);
0152     TH2D * h2_cone_residual_pt_basic_sub1 = new TH2D("h2_cone_residual_pt_basic_sub1", "h2_cone_residual_pt_basic_sub1", resp_N, resp_bins, cent_N, cent_bins);
0153
0154     TH2D * h2_cone_residual_pt_randomize_etaphi_area = new TH2D("h2_cone_residual_pt_randomize_etaphi_area", "h2_cone_residual_pt_randomize_etaphi_area", resp_N, resp_bins, cent_N, cent_bins);
0155     TH2D * h2_cone_residual_pt_randomize_etaphi_mult = new TH2D("h2_cone_residual_pt_randomize_etaphi_mult", "h2_cone_residual_pt_randomize_etaphi_mult", resp_N, resp_bins, cent_N, cent_bins);
0156     TH2D * h2_cone_residual_pt_randomize_etaphi_sub1 = new TH2D("h2_cone_residual_pt_randomize_etaphi_sub1", "h2_cone_residual_pt_randomize_etaphi_sub1", resp_N, resp_bins, cent_N, cent_bins);
0157
0158     // nclustered histograms
0159     TH2D * h2_nclustered_basic = new TH2D("h2_nclustered_basic", "h2_nclustered_basic", cent_N, cent_bins, nclustered_N, nclustered_bins);
0160     TH2D * h2_nclustered_randomize_etaphi = new TH2D("h2_nclustered_randomize_etaphi", "h2_nclustered_randomize_etaphi", cent_N, cent_bins , nclustered_N, nclustered_bins);
0161
0162     // phi and eta cross check histograms
0163     TH1D * h1_phi_cross_check_basic = new TH1D("h1_phi_cross_check_basic", "h1_phi_cross_check_basic", 100, -TMath::Pi(), TMath::Pi());
0164     TH1D * h1_phi_cross_check_randomize_etaphi = new TH1D("h1_phi_cross_check_randomize_etaphi", "h1_phi_cross_check_randomize_etaphi", 100, -TMath::Pi(), TMath::Pi());
0165     TH1D * h1_eta_cross_check_basic = new TH1D("h1_eta_cross_check_basic", "h1_eta_cross_check_basic", 100, -1.5, 1.5);
0166     TH1D * h1_eta_cross_check_randomize_etaphi = new TH1D("h1_eta_cross_check_randomize_etaphi", "h1_eta_cross_check_randomize_etaphi", 100, -1.5, 1.5);
0167
0168     // rho histograms
0169     TH2D * h2_rho_area_rho = new TH2D("h2_rho_area_rho", "h2_rho_area_rho", rho_N, rho_bins, cent_N, cent_bins);
0170     TH2D * h2_rho_mult_rho = new TH2D("h2_rho_mult_rho", "h2_rho_mult_rho", rho_mult_N, rho_bins_mult, cent_N, cent_bins);
0171     TH2D * h2_rho_area_sigma = new TH2D("h2_rho_area_sigma", "h2_rho_area_sigma", rho_N, rho_bins, cent_N, cent_bins);
0172     TH2D * h2_rho_mult_sigma = new TH2D("h2_rho_mult_sigma", "h2_rho_mult_sigma", rho_mult_N, rho_bins_mult, cent_N, cent_bins);
0173
0174     // mbd, zvtx, impact parameter, and centrality histograms
0175     TH1D * h1_mbd = new TH1D("h1_mbd", "h1_mbd", 100, 0, 1);
0176     TH1D * h1_zvtx = new TH1D("h1_zvtx", "h1_zvtx", 100, -50, 50);
0177     TH1D * h1_impactparam = new TH1D("h1_impactparam", "h1_impactparam", 100, 0, 10);
0178     TH1D * h1_centrality = new TH1D("h1_centrality", "h1_centrality", cent_N, cent_bins);
0179
0180
0181     // loop over entries
0182     for (int i = 0; i < nentries; i++)
0183     {
0184         t->GetEntry(i);
0185
0186         // basic cone area and multiplicity
0187         float pt_basic_area = m_RandomCone_Tower_basic_r04_pt - m_TowerRho_AREA_rho*(TMath::Pi()*m_RandomCone_Tower_basic_r04_R*m_RandomCone_Tower_basic_r04_R);
0188         float pt_basic_mult = m_RandomCone_Tower_basic_r04_pt - m_TowerRho_MULT_rho*m_RandomCone_Tower_basic_r04_nclustered;
0189         float pt_basic_sub1 = m_RandomCone_Tower_Sub1_basic_r04_pt;
0190
0191         float pt_randomize_etaphi_area = m_RandomCone_Tower_randomize_etaphi_r04_pt - m_TowerRho_AREA_rho*(TMath::Pi()*m_RandomCone_Tower_randomize_etaphi_r04_R*m_RandomCone_Tower_randomize_etaphi_r04_R);
0192         float pt_randomize_etaphi_mult = m_RandomCone_Tower_randomize_etaphi_r04_pt - m_TowerRho_MULT_rho*m_RandomCone_Tower_randomize_etaphi_r04_nclustered;
0193         float pt_randomize_etaphi_sub1 = m_RandomCone_Tower_Sub1_randomize_etaphi_r04_pt;
0194
0195         // fill histograms
0196         h2_cone_residual_pt_basic_area->Fill(pt_basic_area, m_centrality);
0197         h2_cone_residual_pt_basic_mult->Fill(pt_basic_mult, m_centrality);
0198         h2_cone_residual_pt_basic_sub1->Fill(pt_basic_sub1, m_centrality);
0199
0200         h2_cone_residual_pt_randomize_etaphi_area->Fill(pt_randomize_etaphi_area, m_centrality);
0201         h2_cone_residual_pt_randomize_etaphi_mult->Fill(pt_randomize_etaphi_mult, m_centrality);
0202         h2_cone_residual_pt_randomize_etaphi_sub1->Fill(pt_randomize_etaphi_sub1, m_centrality);
0203
0204         // nclustered histograms
0205         h2_nclustered_basic->Fill(m_RandomCone_Tower_basic_r04_nclustered, m_centrality);
0206         h2_nclustered_randomize_etaphi->Fill(m_RandomCone_Tower_randomize_etaphi_r04_nclustered, m_centrality);
0207
0208         // phi and eta cross check histograms
0209         h1_phi_cross_check_basic->Fill(m_RandomCone_Tower_basic_r04_phi);
0210         h1_phi_cross_check_randomize_etaphi->Fill(m_RandomCone_Tower_randomize_etaphi_r04_phi);
0211         h1_eta_cross_check_basic->Fill(m_RandomCone_Tower_basic_r04_eta);
0212         h1_eta_cross_check_randomize_etaphi->Fill(m_RandomCone_Tower_randomize_etaphi_r04_eta);
0213
0214         // rho histograms
0215         h2_rho_area_rho->Fill(m_TowerRho_AREA_rho, m_centrality);
0216         h2_rho_mult_rho->Fill(m_TowerRho_MULT_rho, m_centrality);
0217         h2_rho_area_sigma->Fill(m_TowerRho_AREA_sigma, m_centrality);
0218         h2_rho_mult_sigma->Fill(m_TowerRho_MULT_sigma, m_centrality);
0219
0220
0221         // mbd, zvtx, impact parameter, and centrality histograms
0222         h1_mbd->Fill(m_mbd_NS);
0223         h1_zvtx->Fill(m_zvtx);
0224         h1_impactparam->Fill(m_impactparam);
0225         h1_centrality->Fill(m_centrality);
0226
0227
0228     }
0229
0230
0231     // make 1D sigma histograms
0232     TH1D * h1_cone_residual_pt_basic_area_sigma = new TH1D("h1_cone_residual_pt_basic_area_sigma", "h1_cone_residual_pt_basic_area_sigma", cent_N, cent_bins);
0233     TH1D * h1_cone_residual_pt_basic_mult_sigma = new TH1D("h1_cone_residual_pt_basic_mult_sigma", "h1_cone_residual_pt_basic_mult_sigma", cent_N, cent_bins);
0234     TH1D * h1_cone_residual_pt_basic_sub1_sigma = new TH1D("h1_cone_residual_pt_basic_sub1_sigma", "h1_cone_residual_pt_basic_sub1_sigma", cent_N, cent_bins);
0235
0236     TH1D * h1_cone_residual_pt_randomize_etaphi_area_sigma = new TH1D("h1_cone_residual_pt_randomize_etaphi_area_sigma", "h1_cone_residual_pt_randomize_etaphi_area_sigma", cent_N, cent_bins);
0237     TH1D * h1_cone_residual_pt_randomize_etaphi_mult_sigma = new TH1D("h1_cone_residual_pt_randomize_etaphi_mult_sigma", "h1_cone_residual_pt_randomize_etaphi_mult_sigma", cent_N, cent_bins);
0238     TH1D * h1_cone_residual_pt_randomize_etaphi_sub1_sigma = new TH1D("h1_cone_residual_pt_randomize_etaphi_sub1_sigma", "h1_cone_residual_pt_randomize_etaphi_sub1_sigma", cent_N, cent_bins);
0239
0240     // make vectors of 2D histograms for easier access
0241     std::vector<TH2D*> h2_cone_residual_pt_vec = {h2_cone_residual_pt_basic_area,
0242                                                   h2_cone_residual_pt_basic_mult,
0243                                                   h2_cone_residual_pt_basic_sub1,
0244                                                   h2_cone_residual_pt_randomize_etaphi_area,
0245                                                   h2_cone_residual_pt_randomize_etaphi_mult,
0246                                                   h2_cone_residual_pt_randomize_etaphi_sub1};
0247
0248     // 1D vector
0249     std::vector<TH1D*> h1_cone_residual_pt_sigma_vec = {h1_cone_residual_pt_basic_area_sigma,
0250                                                         h1_cone_residual_pt_basic_mult_sigma,
0251                                                         h1_cone_residual_pt_basic_sub1_sigma,
0252                                                         h1_cone_residual_pt_randomize_etaphi_area_sigma,
0253                                                         h1_cone_residual_pt_randomize_etaphi_mult_sigma,
0254                                                         h1_cone_residual_pt_randomize_etaphi_sub1_sigma};
0255
0256     // create output file
0257     TFile *fout = new TFile(output_file, "RECREATE");
0258
0259     // loop over 2D histograms, split by centrality and get RMS
0260     for (unsigned int ihist = 0; ihist < h2_cone_residual_pt_vec.size(); ihist++)
0261     {
0262         // get 2D histogram
0263         TH2D * h2_temp = h2_cone_residual_pt_vec.at(ihist);
0264         h2_temp->GetXaxis()->SetTitle("#delta p_{T}");
0265         h2_temp->GetYaxis()->SetTitle("Centrality");
0266
0267         // get 1D histograms
0268         TH1D * h1_temp = h1_cone_residual_pt_sigma_vec.at(ihist);
0269         h1_temp->GetXaxis()->SetTitle("Centrality");
0270         h1_temp->GetYaxis()->SetTitle("#sigma(#delta p_{T})");
0271
0272         // write 2D histograms
0273         fout->cd();
0274         h2_temp->Write();
0275
0276
0277         for (int icent = 1; icent <= cent_N; icent++)
0278         {
0279             h2_temp->GetYaxis()->SetRange(icent, icent);
0280             std::string hist_name = h2_temp->GetName();
0281             TH1D * h1_temp = (TH1D*)h2_temp->ProjectionX(Form("%s_cent%d", hist_name.c_str(), icent));
0282
0283             h1_temp->GetXaxis()->SetTitle("#delta p_{T}");
0284             h1_temp->GetYaxis()->SetTitle("Counts");
0285
0286             double sigma = h1_temp->GetRMS();
0287             double sigma_err = h1_temp->GetRMSError();
0288
0289             h1_cone_residual_pt_sigma_vec.at(ihist)->SetBinContent(icent, sigma);
0290             h1_cone_residual_pt_sigma_vec.at(ihist)->SetBinError(icent, sigma_err);
0291
0292             // write histograms
0293             fout->cd();
0294             h1_temp->Write();
0295
0296             delete h1_temp;
0297         }
0298
0299         // write 1D sigma histograms
0300         fout->cd();
0301         h1_temp->Write();
0302
0303         delete h2_temp;
0304         delete h1_temp;
0305     }
0306
0307
0308     // loop over centrality bins and get rho histograms
0309     for (int icent = 1; icent <= cent_N; icent++)
0310     {
0311         h2_rho_area_rho->GetYaxis()->SetRange(icent, icent);
0312         h2_rho_mult_rho->GetYaxis()->SetRange(icent, icent);
0313         h2_rho_area_sigma->GetYaxis()->SetRange(icent, icent);
0314         h2_rho_mult_sigma->GetYaxis()->SetRange(icent, icent);
0315
0316         TH1D * h1_rho_area_rho = (TH1D*)h2_rho_area_rho->ProjectionX(Form("%s_cent%d", h2_rho_area_rho->GetName(), icent));
0317         TH1D * h1_rho_mult_rho = (TH1D*)h2_rho_mult_rho->ProjectionX(Form("%s_cent%d", h2_rho_mult_rho->GetName(), icent));
0318         TH1D * h1_rho_area_sigma = (TH1D*)h2_rho_area_sigma->ProjectionX(Form("%s_cent%d", h2_rho_area_sigma->GetName(), icent));
0319         TH1D * h1_rho_mult_sigma = (TH1D*)h2_rho_mult_sigma->ProjectionX(Form("%s_cent%d", h2_rho_mult_sigma->GetName(), icent));
0320
0321         h1_rho_area_rho->GetXaxis()->SetTitle("#rho_{area} (GeV)");
0322         h1_rho_mult_rho->GetXaxis()->SetTitle("#rho_{mult} (GeV)");
0323         h1_rho_area_sigma->GetXaxis()->SetTitle("#rho_{area} (GeV)");
0324         h1_rho_mult_sigma->GetXaxis()->SetTitle("#rho_{mult} (GeV)");
0325
0326         // write histograms
0327         fout->cd();
0328         h1_rho_area_rho->Write();
0329         h1_rho_mult_rho->Write();
0330         h1_rho_area_sigma->Write();
0331         h1_rho_mult_sigma->Write();
0332
0333         delete h1_rho_area_rho;
0334         delete h1_rho_mult_rho;
0335         delete h1_rho_area_sigma;
0336         delete h1_rho_mult_sigma;
0337     }
0338
0339
0340
0341     // project nclustered histograms
0342     TProfile * p_nclustered_basic = h2_nclustered_basic->ProfileX();
0343     TProfile * p_nclustered_randomize_etaphi = h2_nclustered_randomize_etaphi->ProfileX();
0344     TH1D * h1_nclustered_basic = (TH1D*)p_nclustered_basic->ProjectionX(Form("%s_projX", h2_nclustered_basic->GetName()));
0345     TH1D * h1_nclustered_randomize_etaphi = (TH1D*)p_nclustered_randomize_etaphi->ProjectionX(Form("%s_projX", h2_nclustered_randomize_etaphi->GetName()));
0346
0347     // write nclustered histograms
0348     fout->cd();
0349     h2_nclustered_basic->Write();
0350     h2_nclustered_randomize_etaphi->Write();
0351         p_nclustered_basic->Write();
0352     p_nclustered_randomize_etaphi->Write();
0353     h1_nclustered_basic->Write();
0354     h1_nclustered_randomize_etaphi->Write();
0355
0356     // write eta and phi cross check histograms
0357     fout->cd();
0358     h1_phi_cross_check_basic->Write();
0359     h1_phi_cross_check_randomize_etaphi->Write();
0360     h1_eta_cross_check_basic->Write();
0361     h1_eta_cross_check_randomize_etaphi->Write();
0362
0363     // write rho zvtx, impact parameter, and centrality histograms
0364     fout->cd();
0365     h1_mbd->Write();
0366     h1_zvtx->Write();
0367     h1_impactparam->Write();
0368     h1_centrality->Write();
0369
0370
0371
0372     // close output file
0373     fout->Close();
0374
0375     // close input file
0376     f->Close();
0377
0378     return 0;
0379
0380 }