123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840 |
- #include "Functions.C"
- void Draw_graphs(TString inFileName="", TString outFileName="./test_graphs.root")
- {
- TFile *fiCorr = new TFile(inFileName.Data(),"read");
- if (inFileName == "" || !fiCorr)
- {
- std::cerr << "No input file was provided!" << std::endl;
- return;
- }
- const std::pair<float, float> bcent_range = {4.06, 8.27};
- const int npid = 6; // h+-, proton, K+, K-, pi+, pi-
- const std::string pidname [] = {"hadrons", "proton", "kaon", "akaon", "pion", "apion"};
- const int ncorrsteps = 4;
- const std::string corrnames [ncorrsteps] = {"PLAIN", "RECENTERED", "TWIST", "RESCALED"};
- const std::string RPcorrname = "PLAIN";
- const int ndatatypes = 3;
- const std::string datatypes [ndatatypes] = {"mc", "reco", "mcreco"};
- // const std::string datatype = "mcreco";
- const int nProj = 2;
- const std::string projv2EP [nProj] = {"cos2", "sin2"};
- const std::string projv1EP [nProj] = {"cos1", "sin1"};
- const std::string projv2SP [nProj] = {"x2", "y2"};
- const std::string projv1SP [nProj] = {"x1", "y1"};
- std::string projPtName = "";
- std::string projYName = "";
- TFile *foGraphs = new TFile(outFileName.Data(),"recreate");
- Qn::DataContainerStatCalculate q2_tpcEP[nProj][nProj], q2_tpcSP[nProj][nProj];
- Qn::DataContainerStatCalculate q1_fhcalEP[nProj][nProj], q1_fhcalSP[nProj][nProj];
- Qn::DataContainerStatCalculate u2_tpcEP1[nProj][nProj][npid], u2_tpcEP2[nProj][nProj][npid];
- Qn::DataContainerStatCalculate u2_tpcSP1[nProj][nProj][npid], u2_tpcSP2[nProj][nProj][npid];
- // Qn::DataContainerStatCalculate u2_rpSP1[nProj][nProj][npid], u2_rpSP2[nProj][nProj][npid];
- Qn::DataContainerStatCalculate u2_fhcalEP1[nProj][nProj][nProj][npid], u2_fhcalEP2[nProj][nProj][nProj][npid];
- Qn::DataContainerStatCalculate u2_fhcalSP1[nProj][nProj][nProj][npid], u2_fhcalSP2[nProj][nProj][nProj][npid];
- Qn::DataContainerStatCalculate u2_rpSP1[nProj][nProj][nProj][npid], u2_rpSP2[nProj][nProj][nProj][npid];
- Qn::DataContainerStatCalculate u1_fhcalEP1[nProj][nProj][npid], u1_fhcalEP2[nProj][nProj][npid];
- Qn::DataContainerStatCalculate u1_fhcalSP1[nProj][nProj][npid], u1_fhcalSP2[nProj][nProj][npid];
- Qn::DataContainerStatCalculate u1_rpSP1[nProj][nProj][npid], u1_rpSP2[nProj][nProj][npid];
- TGraphErrors *graph;
- // Setting up the name of the Qn*Qn correlations
- // Qn*Qn = "dir/nameQn1_correction.nameQn2_correction.component"
- // For example: "QQ/EP/mc_TPC_EP_L_PLAIN.mc_TPC_EP_R_PLAIN.cos2cos2"
- std::string name;
- // Define the correction step that is in the input file
- std::string corrname = "";
- std::string datatype = "";
- for (int idtype=0; idtype<ndatatypes; idtype++)
- {
- for (int icorr=0; icorr<ncorrsteps; icorr++)
- {
- name = datatypes[idtype]+"_TPC_EP_L_"+corrnames[icorr]+"."+datatypes[idtype]+"_TPC_EP_R_"+corrnames[icorr]+".x2x2";
- if (!fiCorr->FindKeyAny(name.c_str())) continue;
- else
- {
- corrname = corrnames[icorr];
- datatype = datatypes[idtype];
- std::cout << "Data type found in the input file : " << datatype << std::endl;
- std::cout << "Correction step found in the input file: " << corrname << std::endl;
- }
- }
- }
- if (corrname == "" || datatype == "")
- {
- std::cerr << "Couldn't find correction or datatype from the file. Exiting..." << std::endl;
- return;
- }
- if (datatype == "mc") projPtName = "_pT";
- if (datatype == "reco") projPtName = "_pT";
- if (datatype == "mcreco") projPtName = "_mc_pT";
- if (projPtName == "")
- {
- std::cerr << "Data type is not correctly set!" << std::endl;
- return;
- }
- if (datatype == "mc") projYName = "_rapidity";
- if (datatype == "reco") projYName = "_rapidity_pdg";
- if (datatype == "mcreco") projYName = "_mc_rapidity";
- if (projYName == "")
- {
- std::cerr << "Data type is not correctly set!" << std::endl;
- return;
- }
- std::string trackBranchName = "";
- if (datatype == "mc") trackBranchName = "McTracks";
- if (datatype == "reco") trackBranchName = "TpcTracks";
- if (datatype == "mcreco") trackBranchName = "TpcTracks";
- if (trackBranchName == "")
- {
- std::cerr << "Data type is not correctly set!" << std::endl;
- return;
- }
- // Getting correlations from the TFile
- // Qn*Qn
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- name = "v2/QQ/EP/"+datatype+"_TPC_EP_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, q2_tpcEP[iproj1][iproj2]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/QQ/SP/"+datatype+"_TPC_EP_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, q2_tpcSP[iproj1][iproj2]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/QQ/EP/"+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, q1_fhcalEP[iproj1][iproj2]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/QQ/SP/"+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, q1_fhcalSP[iproj1][iproj2]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- }
- }
- // un*Qn
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- for (int ipid=0; ipid<npid; ipid++)
- {
- name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u2_tpcEP1[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_TPC_EP_L_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u2_tpcEP2[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u2_tpcSP1[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_TPC_EP_L_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u2_tpcSP2[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
- // if (!GetDCStatCalculate(fiCorr, name, u2_rpSP1[iproj1][iproj2][ipid]))
- // {
- // std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- // }
- // name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
- // if (!GetDCStatCalculate(fiCorr, name, u2_rpSP2[iproj1][iproj2][ipid]))
- // {
- // std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- // }
- name = "v1/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_fhcalEP1[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_fhcalEP2[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_fhcalSP1[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_fhcalSP2[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_rpSP1[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
- if (!GetDCStatCalculate(fiCorr, name, u1_rpSP2[iproj1][iproj2][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- }
- }
- }
- // u2n*Qn*Qn
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- for (int iproj3=0; iproj3<nProj; iproj3++)
- {
- for (int ipid=0; ipid<npid; ipid++)
- {
- name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2EP[iproj1]+projv1EP[iproj2]+projv1EP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_fhcalEP1[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2EP[iproj1]+projv1EP[iproj2]+projv1EP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_fhcalEP2[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_fhcalSP1[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_fhcalSP2[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_rpSP1[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
- if (!GetDCStatCalculate(fiCorr, name, u2_rpSP2[iproj1][iproj2][iproj3][ipid]))
- {
- std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
- }
- }
- }
- }
- }
- // Taking into account that EP2/SP2 is rotated by pi for v1/Q1
- q1_fhcalEP[0][0] = -1. * q1_fhcalEP[0][0];
- q1_fhcalEP[1][1] = -1. * q1_fhcalEP[1][1];
- q1_fhcalSP[0][0] = -1. * q1_fhcalSP[0][0];
- q1_fhcalSP[1][1] = -1. * q1_fhcalSP[1][1];
-
- // Symmetrization for the pT-dependence of v1
- // for (int ipid=0; ipid<npid; ipid++)
- // {
- // u1_fhcalEP2[0][0][ipid] = -1. * u1_fhcalEP2[0][0][ipid];
- // u1_fhcalEP2[1][1][ipid] = -1. * u1_fhcalEP2[1][1][ipid];
- // u1_fhcalSP2[0][0][ipid] = -1. * u1_fhcalSP2[0][0][ipid];
- // u1_fhcalSP2[1][1][ipid] = -1. * u1_fhcalSP2[1][1][ipid];
- // }
- // Calculate resolution: res = Sqrt(<cos(Psi1 - Psi2)>) = Sqrt(<xx> + <yy>)
- Qn::DataContainerStatCalculate res2_tpcEP = Merge(q2_tpcEP[0][0], q2_tpcEP[1][1]);
- Qn::DataContainerStatCalculate res_tpcEP = Sqrt( res2_tpcEP );
- Qn::DataContainerStatCalculate res2_fhcalEP = Merge(q1_fhcalEP[0][0], q1_fhcalEP[1][1]);
- Qn::DataContainerStatCalculate res_fhcalEP = Sqrt( res2_fhcalEP );
- Qn::DataContainerStatCalculate res2_tpcSP = Merge(q2_tpcSP[0][0], q2_tpcSP[1][1]);
- Qn::DataContainerStatCalculate res_tpcSP = Sqrt( res2_tpcSP );
- Qn::DataContainerStatCalculate res2_fhcalSP = Merge(q1_fhcalSP[0][0], q1_fhcalSP[1][1]);
- Qn::DataContainerStatCalculate res_fhcalSP = Sqrt( res2_fhcalSP );
- foGraphs->mkdir("res");
- foGraphs->cd("res");
- name = "res_v2_tpcEP_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(res_tpcEP);
- graph->SetName(Form("%s", name.c_str()));
- name += ";b, fm;Resolution";
- graph->SetTitle(Form("%s", name.c_str()));
- graph->Write();
- name = "res_v1_fhcalEP_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(res_fhcalEP);
- graph->SetName(Form("%s", name.c_str()));
- name += ";b, fm;Resolution";
- graph->SetTitle(Form("%s", name.c_str()));
- graph->Write();
- foGraphs->mkdir("resSP");
- foGraphs->cd("resSP");
- name = "res_v2_tpcSP_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(res_tpcSP);
- graph->SetName(Form("%s", name.c_str()));
- name += ";b, fm;Resolution";
- graph->SetTitle(Form("%s", name.c_str()));
- graph->Write();
- name = "res_v1_fhcalSP_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(res_fhcalSP);
- graph->SetName(Form("%s", name.c_str()));
- name += ";b, fm;Resolution";
- graph->SetTitle(Form("%s", name.c_str()));
- graph->Write();
- Qn::DataContainerStatCalculate vn_pT_obs, vn_pT_div, vn_pT_cent, vn_pT;
- Qn::DataContainerStatCalculate vn_y_obs, vn_y_div, vn_y_cent, vn_y;
- Qn::DataContainerStatCalculate uq_cent, uq_pT, uq_y;
- Qn::DataContainerStatCalculate vn_comp_pT_obs, vn_comp_pT_div, vn_comp_pT_cent, vn_comp_pT;
- Qn::DataContainerStatCalculate vn_comp_y_obs, vn_comp_y_div, vn_comp_y_cent, vn_comp_y;
- // v2 TPC EP
- foGraphs->mkdir("v2tpcEP");
- foGraphs->cd("v2tpcEP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- vn_pT_obs = Merge(Merge(u2_tpcEP1[0][0][ipid], u2_tpcEP1[1][1][ipid]), Merge(u2_tpcEP2[0][0][ipid], u2_tpcEP2[1][1][ipid]));
- vn_pT_div = vn_pT_obs / res_tpcEP;
- vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_pT = vn_pT_cent.Projection({name});
- vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_pT_tpcEP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_y_obs = Merge(Merge(u2_tpcEP1[0][0][ipid], u2_tpcEP1[1][1][ipid]), Merge(u2_tpcEP2[0][0][ipid], u2_tpcEP2[1][1][ipid]));
- vn_y_div = vn_y_obs / res_tpcEP;
- vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_y = vn_y_cent.Projection({name});
- vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_y_tpcEP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_y);
- graph->SetName(name.c_str());
- name += ";y;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- // v2 TPC SP
- foGraphs->mkdir("v2tpcSP");
- foGraphs->cd("v2tpcSP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- vn_pT_obs = Merge(Merge(u2_tpcSP1[0][0][ipid], u2_tpcSP1[1][1][ipid]), Merge(u2_tpcSP2[0][0][ipid], u2_tpcSP2[1][1][ipid]));
- vn_pT_div = vn_pT_obs / res_tpcSP;
- vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_pT = vn_pT_cent.Projection({name});
- vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_pT_tpcSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_y_obs = Merge(Merge(u2_tpcSP1[0][0][ipid], u2_tpcSP1[1][1][ipid]), Merge(u2_tpcSP2[0][0][ipid], u2_tpcSP2[1][1][ipid]));
- vn_y_div = vn_y_obs / res_tpcSP;
- vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_y = vn_y_cent.Projection({name});
- vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_y_tpcSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_y);
- graph->SetName(name.c_str());
- name += ";y;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- // v1 FHCal EP
- foGraphs->mkdir("v1fhcalEP");
- foGraphs->cd("v1fhcalEP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- vn_pT_obs = Merge(Merge(u1_fhcalEP1[0][0][ipid], u1_fhcalEP1[1][1][ipid]), Merge(u1_fhcalEP2[0][0][ipid], u1_fhcalEP2[1][1][ipid]));
- if (corrname == "mc")
- {
- vn_pT_obs = -1.*vn_pT_obs;
- }
- vn_pT_div = vn_pT_obs / res_fhcalEP;
- vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_pT = vn_pT_cent.Projection({name});
- vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_y_obs = Merge(Merge(-1.*u1_fhcalEP1[0][0][ipid], -1.*u1_fhcalEP1[1][1][ipid]), Merge(u1_fhcalEP2[0][0][ipid], u1_fhcalEP2[1][1][ipid]));
- if (corrname == "mc")
- {
- vn_y_obs = -1.*vn_y_obs;
- }
- vn_y_div = vn_y_obs / res_fhcalEP;
- vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_y = vn_y_cent.Projection({name});
- vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- // v1 FHCal SP
- foGraphs->mkdir("v1fhcalSP");
- foGraphs->cd("v1fhcalSP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- vn_pT_obs = Merge(Merge(u1_fhcalSP1[0][0][ipid], u1_fhcalSP1[1][1][ipid]), Merge(u1_fhcalSP2[0][0][ipid], u1_fhcalSP2[1][1][ipid]));
- if (datatype == "mc")
- {
- vn_pT_obs = -1. * vn_pT_obs;
- }
- vn_pT_div = vn_pT_obs / res_fhcalSP;
- vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_pT = vn_pT_cent.Projection({name});
- vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_fhcalSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_y_obs = Merge(Merge(-1.*u1_fhcalSP1[0][0][ipid], -1.*u1_fhcalSP1[1][1][ipid]), Merge(u1_fhcalSP2[0][0][ipid], u1_fhcalSP2[1][1][ipid]));
- if (datatype == "mc")
- {
- vn_y_obs = -1. * vn_y_obs;
- }
- vn_y_div = vn_y_obs / res_fhcalSP;
- vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_y = vn_y_cent.Projection({name});
- vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_fhcalSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- // v1 RP SP
- foGraphs->mkdir("v1rpSP");
- foGraphs->cd("v1rpSP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- vn_pT_obs = Merge(Merge(-1.*u1_rpSP1[0][0][ipid], -1.*u1_rpSP1[1][1][ipid]), Merge(u1_rpSP2[0][0][ipid], u1_rpSP2[1][1][ipid]));
- vn_pT_div = vn_pT_obs;
- vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_pT = vn_pT_cent.Projection({name});
- vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_rpSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_y_obs = Merge(Merge(u1_rpSP1[0][0][ipid], u1_rpSP1[1][1][ipid]), Merge(u1_rpSP2[0][0][ipid], u1_rpSP2[1][1][ipid]));
- vn_y_div = vn_y_obs;
- vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_y = vn_y_cent.Projection({name});
- vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_rpSP_" + pidname[ipid] + "_" + corrname;
- graph = Qn::DataContainerHelper::ToTGraph(vn_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- // QQ correlations
- foGraphs->mkdir("QQ");
- foGraphs->cd("QQ");
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- name = "QQ_tpc_"+ corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(q2_tpcEP[iproj1][iproj2]);
- graph->SetName(name.c_str());
- name += ";b, fm;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = "QQ_tpc_"+ corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(q2_tpcSP[iproj1][iproj2]);
- graph->SetName(name.c_str());
- name += ";b, fm;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = "QQ_fhcal_"+ corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(q1_fhcalEP[iproj1][iproj2]);
- graph->SetName(name.c_str());
- name += ";b, fm;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = "QQ_fhcal_"+ corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(q1_fhcalSP[iproj1][iproj2]);
- graph->SetName(name.c_str());
- name += ";b, fm;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- // uQ/QQ correlations
- foGraphs->mkdir("v1Components_EP");
- foGraphs->cd("v1Components_EP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- vn_comp_pT_obs = Merge(u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
- vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
- vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_comp_pT = vn_comp_pT_cent.Projection({name});
- vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_comp_y_obs = Merge(-1.*u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
- vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
- vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_comp_y = vn_comp_y_cent.Projection({name});
- vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- }
- foGraphs->mkdir("v2Components_EP");
- foGraphs->cd("v2Components_EP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- vn_comp_pT_obs = Merge(u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
- vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
- vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_comp_pT = vn_comp_pT_cent.Projection({name});
- vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_comp_y_obs = Merge(-1.*u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
- vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
- vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_comp_y = vn_comp_y_cent.Projection({name});
- vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- }
- foGraphs->mkdir("v1Components_SP");
- foGraphs->cd("v1Components_SP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- vn_comp_pT_obs = Merge(u1_fhcalSP1[iproj1][iproj2][ipid], u1_fhcalSP2[iproj1][iproj2][ipid]);
- vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalSP[iproj1][iproj2]);
- vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_comp_pT = vn_comp_pT_cent.Projection({name});
- vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_pT_fhcalSP_" + pidname[ipid] + "_" + corrname +"_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_comp_y_obs = Merge(-1.*u1_fhcalSP1[iproj1][iproj2][ipid], u1_fhcalSP2[iproj1][iproj2][ipid]);
- vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalSP[iproj1][iproj2]);
- vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_comp_y = vn_comp_y_cent.Projection({name});
- vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v1_y_fhcalSP_" + pidname[ipid] + "_" + corrname +"_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
- graph->SetName(name.c_str());
- name += ";y;v_{1}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- }
- foGraphs->mkdir("v2Components_SP");
- foGraphs->cd("v2Components_SP");
- for (int ipid=0; ipid<npid; ipid++)
- {
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- vn_comp_pT_obs = Merge(u2_tpcSP1[iproj1][iproj2][ipid], u2_tpcSP2[iproj1][iproj2][ipid]);
- vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q2_tpcSP[iproj1][iproj2]);
- vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- vn_comp_pT = vn_comp_pT_cent.Projection({name});
- vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_pT_tpcSP_" + pidname[ipid] + "_" + corrname +"_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- vn_comp_y_obs = Merge(u2_tpcSP1[iproj1][iproj2][ipid], u2_tpcSP2[iproj1][iproj2][ipid]);
- vn_comp_y_div = vn_comp_y_obs / Sqrt(q2_tpcSP[iproj1][iproj2]);
- vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projYName;
- vn_comp_y = vn_comp_y_cent.Projection({name});
- vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "v2_y_tpcSP_" + pidname[ipid] + "_" + corrname +"_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
- graph->SetName(name.c_str());
- name += ";y;v_{2}";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- }
- // uQ correlations
- foGraphs->mkdir("uQ");
- foGraphs->cd("uQ");
- for (int ipid=0; ipid<npid; ipid++)
- {
- for (int iproj1=0; iproj1<nProj; iproj1++)
- {
- for (int iproj2=0; iproj2<nProj; iproj2++)
- {
- uq_cent = u2_tpcEP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_tpcEP1_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_tpcEP1_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u2_tpcEP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_tpcEP2_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_tpcEP2_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u2_tpcSP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_tpcSP1_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_tpcSP1_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u2_tpcSP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_tpcSP2_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_tpcSP2_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u1_fhcalEP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_fhcalEP1_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_fhcalEP1_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u1_fhcalEP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_fhcalEP2_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_fhcalEP2_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u1_fhcalSP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_fhcalSP1_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_fhcalSP1_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- uq_cent = u1_fhcalSP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
- name = trackBranchName+projPtName;
- uq_pT = uq_cent.Projection({name});
- uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_pT_fhcalSP2_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
- graph->SetName(name.c_str());
- name += ";p_{T}, GeV/c;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- name = trackBranchName+projYName;
- uq_y = uq_cent.Projection({name});
- uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
- name = "uQ_y_fhcalSP2_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
- graph = Qn::DataContainerHelper::ToTGraph(uq_y);
- graph->SetName(name.c_str());
- name += ";y;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
- graph->SetTitle(name.c_str());
- graph->Write();
- }
- }
- }
- foGraphs->Close();
- }
|