Draw_graphs.C 43 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840
  1. #include "Functions.C"
  2. void Draw_graphs(TString inFileName="", TString outFileName="./test_graphs.root")
  3. {
  4. TFile *fiCorr = new TFile(inFileName.Data(),"read");
  5. if (inFileName == "" || !fiCorr)
  6. {
  7. std::cerr << "No input file was provided!" << std::endl;
  8. return;
  9. }
  10. const std::pair<float, float> bcent_range = {4.06, 8.27};
  11. const int npid = 6; // h+-, proton, K+, K-, pi+, pi-
  12. const std::string pidname [] = {"hadrons", "proton", "kaon", "akaon", "pion", "apion"};
  13. const int ncorrsteps = 4;
  14. const std::string corrnames [ncorrsteps] = {"PLAIN", "RECENTERED", "TWIST", "RESCALED"};
  15. const std::string RPcorrname = "PLAIN";
  16. const int ndatatypes = 3;
  17. const std::string datatypes [ndatatypes] = {"mc", "reco", "mcreco"};
  18. // const std::string datatype = "mcreco";
  19. const int nProj = 2;
  20. const std::string projv2EP [nProj] = {"cos2", "sin2"};
  21. const std::string projv1EP [nProj] = {"cos1", "sin1"};
  22. const std::string projv2SP [nProj] = {"x2", "y2"};
  23. const std::string projv1SP [nProj] = {"x1", "y1"};
  24. std::string projPtName = "";
  25. std::string projYName = "";
  26. TFile *foGraphs = new TFile(outFileName.Data(),"recreate");
  27. Qn::DataContainerStatCalculate q2_tpcEP[nProj][nProj], q2_tpcSP[nProj][nProj];
  28. Qn::DataContainerStatCalculate q1_fhcalEP[nProj][nProj], q1_fhcalSP[nProj][nProj];
  29. Qn::DataContainerStatCalculate u2_tpcEP1[nProj][nProj][npid], u2_tpcEP2[nProj][nProj][npid];
  30. Qn::DataContainerStatCalculate u2_tpcSP1[nProj][nProj][npid], u2_tpcSP2[nProj][nProj][npid];
  31. // Qn::DataContainerStatCalculate u2_rpSP1[nProj][nProj][npid], u2_rpSP2[nProj][nProj][npid];
  32. Qn::DataContainerStatCalculate u2_fhcalEP1[nProj][nProj][nProj][npid], u2_fhcalEP2[nProj][nProj][nProj][npid];
  33. Qn::DataContainerStatCalculate u2_fhcalSP1[nProj][nProj][nProj][npid], u2_fhcalSP2[nProj][nProj][nProj][npid];
  34. Qn::DataContainerStatCalculate u2_rpSP1[nProj][nProj][nProj][npid], u2_rpSP2[nProj][nProj][nProj][npid];
  35. Qn::DataContainerStatCalculate u1_fhcalEP1[nProj][nProj][npid], u1_fhcalEP2[nProj][nProj][npid];
  36. Qn::DataContainerStatCalculate u1_fhcalSP1[nProj][nProj][npid], u1_fhcalSP2[nProj][nProj][npid];
  37. Qn::DataContainerStatCalculate u1_rpSP1[nProj][nProj][npid], u1_rpSP2[nProj][nProj][npid];
  38. TGraphErrors *graph;
  39. // Setting up the name of the Qn*Qn correlations
  40. // Qn*Qn = "dir/nameQn1_correction.nameQn2_correction.component"
  41. // For example: "QQ/EP/mc_TPC_EP_L_PLAIN.mc_TPC_EP_R_PLAIN.cos2cos2"
  42. std::string name;
  43. // Define the correction step that is in the input file
  44. std::string corrname = "";
  45. std::string datatype = "";
  46. for (int idtype=0; idtype<ndatatypes; idtype++)
  47. {
  48. for (int icorr=0; icorr<ncorrsteps; icorr++)
  49. {
  50. name = datatypes[idtype]+"_TPC_EP_L_"+corrnames[icorr]+"."+datatypes[idtype]+"_TPC_EP_R_"+corrnames[icorr]+".x2x2";
  51. if (!fiCorr->FindKeyAny(name.c_str())) continue;
  52. else
  53. {
  54. corrname = corrnames[icorr];
  55. datatype = datatypes[idtype];
  56. std::cout << "Data type found in the input file : " << datatype << std::endl;
  57. std::cout << "Correction step found in the input file: " << corrname << std::endl;
  58. }
  59. }
  60. }
  61. if (corrname == "" || datatype == "")
  62. {
  63. std::cerr << "Couldn't find correction or datatype from the file. Exiting..." << std::endl;
  64. return;
  65. }
  66. if (datatype == "mc") projPtName = "_pT";
  67. if (datatype == "reco") projPtName = "_pT";
  68. if (datatype == "mcreco") projPtName = "_mc_pT";
  69. if (projPtName == "")
  70. {
  71. std::cerr << "Data type is not correctly set!" << std::endl;
  72. return;
  73. }
  74. if (datatype == "mc") projYName = "_rapidity";
  75. if (datatype == "reco") projYName = "_rapidity_pdg";
  76. if (datatype == "mcreco") projYName = "_mc_rapidity";
  77. if (projYName == "")
  78. {
  79. std::cerr << "Data type is not correctly set!" << std::endl;
  80. return;
  81. }
  82. std::string trackBranchName = "";
  83. if (datatype == "mc") trackBranchName = "McTracks";
  84. if (datatype == "reco") trackBranchName = "TpcTracks";
  85. if (datatype == "mcreco") trackBranchName = "TpcTracks";
  86. if (trackBranchName == "")
  87. {
  88. std::cerr << "Data type is not correctly set!" << std::endl;
  89. return;
  90. }
  91. // Getting correlations from the TFile
  92. // Qn*Qn
  93. for (int iproj1=0; iproj1<nProj; iproj1++)
  94. {
  95. for (int iproj2=0; iproj2<nProj; iproj2++)
  96. {
  97. name = "v2/QQ/EP/"+datatype+"_TPC_EP_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
  98. if (!GetDCStatCalculate(fiCorr, name, q2_tpcEP[iproj1][iproj2]))
  99. {
  100. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  101. }
  102. name = "v2/QQ/SP/"+datatype+"_TPC_EP_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
  103. if (!GetDCStatCalculate(fiCorr, name, q2_tpcSP[iproj1][iproj2]))
  104. {
  105. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  106. }
  107. name = "v1/QQ/EP/"+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
  108. if (!GetDCStatCalculate(fiCorr, name, q1_fhcalEP[iproj1][iproj2]))
  109. {
  110. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  111. }
  112. name = "v1/QQ/SP/"+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
  113. if (!GetDCStatCalculate(fiCorr, name, q1_fhcalSP[iproj1][iproj2]))
  114. {
  115. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  116. }
  117. }
  118. }
  119. // un*Qn
  120. for (int iproj1=0; iproj1<nProj; iproj1++)
  121. {
  122. for (int iproj2=0; iproj2<nProj; iproj2++)
  123. {
  124. for (int ipid=0; ipid<npid; ipid++)
  125. {
  126. name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
  127. if (!GetDCStatCalculate(fiCorr, name, u2_tpcEP1[iproj1][iproj2][ipid]))
  128. {
  129. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  130. }
  131. name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_TPC_EP_L_"+corrname+"."+projv2EP[iproj1]+projv2EP[iproj2];
  132. if (!GetDCStatCalculate(fiCorr, name, u2_tpcEP2[iproj1][iproj2][ipid]))
  133. {
  134. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  135. }
  136. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_TPC_EP_R_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
  137. if (!GetDCStatCalculate(fiCorr, name, u2_tpcSP1[iproj1][iproj2][ipid]))
  138. {
  139. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  140. }
  141. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_TPC_EP_L_"+corrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
  142. if (!GetDCStatCalculate(fiCorr, name, u2_tpcSP2[iproj1][iproj2][ipid]))
  143. {
  144. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  145. }
  146. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
  147. // if (!GetDCStatCalculate(fiCorr, name, u2_rpSP1[iproj1][iproj2][ipid]))
  148. // {
  149. // std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  150. // }
  151. // name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv2SP[iproj2];
  152. // if (!GetDCStatCalculate(fiCorr, name, u2_rpSP2[iproj1][iproj2][ipid]))
  153. // {
  154. // std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  155. // }
  156. name = "v1/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
  157. if (!GetDCStatCalculate(fiCorr, name, u1_fhcalEP1[iproj1][iproj2][ipid]))
  158. {
  159. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  160. }
  161. name = "v1/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+projv1EP[iproj1]+projv1EP[iproj2];
  162. if (!GetDCStatCalculate(fiCorr, name, u1_fhcalEP2[iproj1][iproj2][ipid]))
  163. {
  164. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  165. }
  166. name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
  167. if (!GetDCStatCalculate(fiCorr, name, u1_fhcalSP1[iproj1][iproj2][ipid]))
  168. {
  169. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  170. }
  171. name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
  172. if (!GetDCStatCalculate(fiCorr, name, u1_fhcalSP2[iproj1][iproj2][ipid]))
  173. {
  174. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  175. }
  176. name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
  177. if (!GetDCStatCalculate(fiCorr, name, u1_rpSP1[iproj1][iproj2][ipid]))
  178. {
  179. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  180. }
  181. name = "v1/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+projv1SP[iproj1]+projv1SP[iproj2];
  182. if (!GetDCStatCalculate(fiCorr, name, u1_rpSP2[iproj1][iproj2][ipid]))
  183. {
  184. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  185. }
  186. }
  187. }
  188. }
  189. // u2n*Qn*Qn
  190. for (int iproj1=0; iproj1<nProj; iproj1++)
  191. {
  192. for (int iproj2=0; iproj2<nProj; iproj2++)
  193. {
  194. for (int iproj3=0; iproj3<nProj; iproj3++)
  195. {
  196. for (int ipid=0; ipid<npid; ipid++)
  197. {
  198. name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2EP[iproj1]+projv1EP[iproj2]+projv1EP[iproj3];
  199. if (!GetDCStatCalculate(fiCorr, name, u2_fhcalEP1[iproj1][iproj2][iproj3][ipid]))
  200. {
  201. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  202. }
  203. name = "v2/uQ/EP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2EP[iproj1]+projv1EP[iproj2]+projv1EP[iproj3];
  204. if (!GetDCStatCalculate(fiCorr, name, u2_fhcalEP2[iproj1][iproj2][iproj3][ipid]))
  205. {
  206. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  207. }
  208. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
  209. if (!GetDCStatCalculate(fiCorr, name, u2_fhcalSP1[iproj1][iproj2][iproj3][ipid]))
  210. {
  211. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  212. }
  213. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_FHCal_L_"+corrname+"."+datatype+"_FHCal_R_"+corrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
  214. if (!GetDCStatCalculate(fiCorr, name, u2_fhcalSP2[iproj1][iproj2][iproj3][ipid]))
  215. {
  216. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  217. }
  218. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_L_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
  219. if (!GetDCStatCalculate(fiCorr, name, u2_rpSP1[iproj1][iproj2][iproj3][ipid]))
  220. {
  221. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  222. }
  223. name = "v2/uQ/SP/"+datatype+"_" + pidname[ipid] + "_R_"+corrname+"."+datatype+"_RP_"+RPcorrname+"."+datatype+"_RP_"+RPcorrname+"."+projv2SP[iproj1]+projv1SP[iproj2]+projv1SP[iproj3];
  224. if (!GetDCStatCalculate(fiCorr, name, u2_rpSP2[iproj1][iproj2][iproj3][ipid]))
  225. {
  226. std::cerr << "Error: Cannot get " << name << " from the file!" << std::endl;
  227. }
  228. }
  229. }
  230. }
  231. }
  232. // Taking into account that EP2/SP2 is rotated by pi for v1/Q1
  233. q1_fhcalEP[0][0] = -1. * q1_fhcalEP[0][0];
  234. q1_fhcalEP[1][1] = -1. * q1_fhcalEP[1][1];
  235. q1_fhcalSP[0][0] = -1. * q1_fhcalSP[0][0];
  236. q1_fhcalSP[1][1] = -1. * q1_fhcalSP[1][1];
  237. // Symmetrization for the pT-dependence of v1
  238. // for (int ipid=0; ipid<npid; ipid++)
  239. // {
  240. // u1_fhcalEP2[0][0][ipid] = -1. * u1_fhcalEP2[0][0][ipid];
  241. // u1_fhcalEP2[1][1][ipid] = -1. * u1_fhcalEP2[1][1][ipid];
  242. // u1_fhcalSP2[0][0][ipid] = -1. * u1_fhcalSP2[0][0][ipid];
  243. // u1_fhcalSP2[1][1][ipid] = -1. * u1_fhcalSP2[1][1][ipid];
  244. // }
  245. // Calculate resolution: res = Sqrt(<cos(Psi1 - Psi2)>) = Sqrt(<xx> + <yy>)
  246. Qn::DataContainerStatCalculate res2_tpcEP = Merge(q2_tpcEP[0][0], q2_tpcEP[1][1]);
  247. Qn::DataContainerStatCalculate res_tpcEP = Sqrt( res2_tpcEP );
  248. Qn::DataContainerStatCalculate res2_fhcalEP = Merge(q1_fhcalEP[0][0], q1_fhcalEP[1][1]);
  249. Qn::DataContainerStatCalculate res_fhcalEP = Sqrt( res2_fhcalEP );
  250. Qn::DataContainerStatCalculate res2_tpcSP = Merge(q2_tpcSP[0][0], q2_tpcSP[1][1]);
  251. Qn::DataContainerStatCalculate res_tpcSP = Sqrt( res2_tpcSP );
  252. Qn::DataContainerStatCalculate res2_fhcalSP = Merge(q1_fhcalSP[0][0], q1_fhcalSP[1][1]);
  253. Qn::DataContainerStatCalculate res_fhcalSP = Sqrt( res2_fhcalSP );
  254. foGraphs->mkdir("res");
  255. foGraphs->cd("res");
  256. name = "res_v2_tpcEP_" + corrname;
  257. graph = Qn::DataContainerHelper::ToTGraph(res_tpcEP);
  258. graph->SetName(Form("%s", name.c_str()));
  259. name += ";b, fm;Resolution";
  260. graph->SetTitle(Form("%s", name.c_str()));
  261. graph->Write();
  262. name = "res_v1_fhcalEP_" + corrname;
  263. graph = Qn::DataContainerHelper::ToTGraph(res_fhcalEP);
  264. graph->SetName(Form("%s", name.c_str()));
  265. name += ";b, fm;Resolution";
  266. graph->SetTitle(Form("%s", name.c_str()));
  267. graph->Write();
  268. foGraphs->mkdir("resSP");
  269. foGraphs->cd("resSP");
  270. name = "res_v2_tpcSP_" + corrname;
  271. graph = Qn::DataContainerHelper::ToTGraph(res_tpcSP);
  272. graph->SetName(Form("%s", name.c_str()));
  273. name += ";b, fm;Resolution";
  274. graph->SetTitle(Form("%s", name.c_str()));
  275. graph->Write();
  276. name = "res_v1_fhcalSP_" + corrname;
  277. graph = Qn::DataContainerHelper::ToTGraph(res_fhcalSP);
  278. graph->SetName(Form("%s", name.c_str()));
  279. name += ";b, fm;Resolution";
  280. graph->SetTitle(Form("%s", name.c_str()));
  281. graph->Write();
  282. Qn::DataContainerStatCalculate vn_pT_obs, vn_pT_div, vn_pT_cent, vn_pT;
  283. Qn::DataContainerStatCalculate vn_y_obs, vn_y_div, vn_y_cent, vn_y;
  284. Qn::DataContainerStatCalculate uq_cent, uq_pT, uq_y;
  285. Qn::DataContainerStatCalculate vn_comp_pT_obs, vn_comp_pT_div, vn_comp_pT_cent, vn_comp_pT;
  286. Qn::DataContainerStatCalculate vn_comp_y_obs, vn_comp_y_div, vn_comp_y_cent, vn_comp_y;
  287. // v2 TPC EP
  288. foGraphs->mkdir("v2tpcEP");
  289. foGraphs->cd("v2tpcEP");
  290. for (int ipid=0; ipid<npid; ipid++)
  291. {
  292. 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]));
  293. vn_pT_div = vn_pT_obs / res_tpcEP;
  294. vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  295. name = trackBranchName+projPtName;
  296. vn_pT = vn_pT_cent.Projection({name});
  297. vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  298. name = "v2_pT_tpcEP_" + pidname[ipid] + "_" + corrname;
  299. graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
  300. graph->SetName(name.c_str());
  301. name += ";p_{T}, GeV/c;v_{2}";
  302. graph->SetTitle(name.c_str());
  303. graph->Write();
  304. 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]));
  305. vn_y_div = vn_y_obs / res_tpcEP;
  306. vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  307. name = trackBranchName+projYName;
  308. vn_y = vn_y_cent.Projection({name});
  309. vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  310. name = "v2_y_tpcEP_" + pidname[ipid] + "_" + corrname;
  311. graph = Qn::DataContainerHelper::ToTGraph(vn_y);
  312. graph->SetName(name.c_str());
  313. name += ";y;v_{2}";
  314. graph->SetTitle(name.c_str());
  315. graph->Write();
  316. }
  317. // v2 TPC SP
  318. foGraphs->mkdir("v2tpcSP");
  319. foGraphs->cd("v2tpcSP");
  320. for (int ipid=0; ipid<npid; ipid++)
  321. {
  322. 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]));
  323. vn_pT_div = vn_pT_obs / res_tpcSP;
  324. vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  325. name = trackBranchName+projPtName;
  326. vn_pT = vn_pT_cent.Projection({name});
  327. vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  328. name = "v2_pT_tpcSP_" + pidname[ipid] + "_" + corrname;
  329. graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
  330. graph->SetName(name.c_str());
  331. name += ";p_{T}, GeV/c;v_{2}";
  332. graph->SetTitle(name.c_str());
  333. graph->Write();
  334. 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]));
  335. vn_y_div = vn_y_obs / res_tpcSP;
  336. vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  337. name = trackBranchName+projYName;
  338. vn_y = vn_y_cent.Projection({name});
  339. vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  340. name = "v2_y_tpcSP_" + pidname[ipid] + "_" + corrname;
  341. graph = Qn::DataContainerHelper::ToTGraph(vn_y);
  342. graph->SetName(name.c_str());
  343. name += ";y;v_{2}";
  344. graph->SetTitle(name.c_str());
  345. graph->Write();
  346. }
  347. // v1 FHCal EP
  348. foGraphs->mkdir("v1fhcalEP");
  349. foGraphs->cd("v1fhcalEP");
  350. for (int ipid=0; ipid<npid; ipid++)
  351. {
  352. 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]));
  353. if (corrname == "mc")
  354. {
  355. vn_pT_obs = -1.*vn_pT_obs;
  356. }
  357. vn_pT_div = vn_pT_obs / res_fhcalEP;
  358. vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  359. name = trackBranchName+projPtName;
  360. vn_pT = vn_pT_cent.Projection({name});
  361. vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  362. name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname;
  363. graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
  364. graph->SetName(name.c_str());
  365. name += ";p_{T}, GeV/c;v_{1}";
  366. graph->SetTitle(name.c_str());
  367. graph->Write();
  368. 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]));
  369. if (corrname == "mc")
  370. {
  371. vn_y_obs = -1.*vn_y_obs;
  372. }
  373. vn_y_div = vn_y_obs / res_fhcalEP;
  374. vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  375. name = trackBranchName+projYName;
  376. vn_y = vn_y_cent.Projection({name});
  377. vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  378. name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname;
  379. graph = Qn::DataContainerHelper::ToTGraph(vn_y);
  380. graph->SetName(name.c_str());
  381. name += ";y;v_{1}";
  382. graph->SetTitle(name.c_str());
  383. graph->Write();
  384. }
  385. // v1 FHCal SP
  386. foGraphs->mkdir("v1fhcalSP");
  387. foGraphs->cd("v1fhcalSP");
  388. for (int ipid=0; ipid<npid; ipid++)
  389. {
  390. 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]));
  391. if (datatype == "mc")
  392. {
  393. vn_pT_obs = -1. * vn_pT_obs;
  394. }
  395. vn_pT_div = vn_pT_obs / res_fhcalSP;
  396. vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  397. name = trackBranchName+projPtName;
  398. vn_pT = vn_pT_cent.Projection({name});
  399. vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  400. name = "v1_pT_fhcalSP_" + pidname[ipid] + "_" + corrname;
  401. graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
  402. graph->SetName(name.c_str());
  403. name += ";p_{T}, GeV/c;v_{1}";
  404. graph->SetTitle(name.c_str());
  405. graph->Write();
  406. 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]));
  407. if (datatype == "mc")
  408. {
  409. vn_y_obs = -1. * vn_y_obs;
  410. }
  411. vn_y_div = vn_y_obs / res_fhcalSP;
  412. vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  413. name = trackBranchName+projYName;
  414. vn_y = vn_y_cent.Projection({name});
  415. vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  416. name = "v1_y_fhcalSP_" + pidname[ipid] + "_" + corrname;
  417. graph = Qn::DataContainerHelper::ToTGraph(vn_y);
  418. graph->SetName(name.c_str());
  419. name += ";y;v_{1}";
  420. graph->SetTitle(name.c_str());
  421. graph->Write();
  422. }
  423. // v1 RP SP
  424. foGraphs->mkdir("v1rpSP");
  425. foGraphs->cd("v1rpSP");
  426. for (int ipid=0; ipid<npid; ipid++)
  427. {
  428. 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]));
  429. vn_pT_div = vn_pT_obs;
  430. vn_pT_cent = vn_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  431. name = trackBranchName+projPtName;
  432. vn_pT = vn_pT_cent.Projection({name});
  433. vn_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  434. name = "v1_pT_rpSP_" + pidname[ipid] + "_" + corrname;
  435. graph = Qn::DataContainerHelper::ToTGraph(vn_pT);
  436. graph->SetName(name.c_str());
  437. name += ";p_{T}, GeV/c;v_{1}";
  438. graph->SetTitle(name.c_str());
  439. graph->Write();
  440. 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]));
  441. vn_y_div = vn_y_obs;
  442. vn_y_cent = vn_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  443. name = trackBranchName+projYName;
  444. vn_y = vn_y_cent.Projection({name});
  445. vn_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  446. name = "v1_y_rpSP_" + pidname[ipid] + "_" + corrname;
  447. graph = Qn::DataContainerHelper::ToTGraph(vn_y);
  448. graph->SetName(name.c_str());
  449. name += ";y;v_{1}";
  450. graph->SetTitle(name.c_str());
  451. graph->Write();
  452. }
  453. // QQ correlations
  454. foGraphs->mkdir("QQ");
  455. foGraphs->cd("QQ");
  456. for (int iproj1=0; iproj1<nProj; iproj1++)
  457. {
  458. for (int iproj2=0; iproj2<nProj; iproj2++)
  459. {
  460. name = "QQ_tpc_"+ corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
  461. graph = Qn::DataContainerHelper::ToTGraph(q2_tpcEP[iproj1][iproj2]);
  462. graph->SetName(name.c_str());
  463. name += ";b, fm;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
  464. graph->SetTitle(name.c_str());
  465. graph->Write();
  466. name = "QQ_tpc_"+ corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
  467. graph = Qn::DataContainerHelper::ToTGraph(q2_tpcSP[iproj1][iproj2]);
  468. graph->SetName(name.c_str());
  469. name += ";b, fm;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
  470. graph->SetTitle(name.c_str());
  471. graph->Write();
  472. name = "QQ_fhcal_"+ corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
  473. graph = Qn::DataContainerHelper::ToTGraph(q1_fhcalEP[iproj1][iproj2]);
  474. graph->SetName(name.c_str());
  475. name += ";b, fm;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
  476. graph->SetTitle(name.c_str());
  477. graph->Write();
  478. name = "QQ_fhcal_"+ corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
  479. graph = Qn::DataContainerHelper::ToTGraph(q1_fhcalSP[iproj1][iproj2]);
  480. graph->SetName(name.c_str());
  481. name += ";b, fm;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
  482. graph->SetTitle(name.c_str());
  483. graph->Write();
  484. }
  485. }
  486. // uQ/QQ correlations
  487. foGraphs->mkdir("v1Components_EP");
  488. foGraphs->cd("v1Components_EP");
  489. for (int ipid=0; ipid<npid; ipid++)
  490. {
  491. for (int iproj1=0; iproj1<nProj; iproj1++)
  492. {
  493. for (int iproj2=0; iproj2<nProj; iproj2++)
  494. {
  495. vn_comp_pT_obs = Merge(u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
  496. vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
  497. vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  498. name = trackBranchName+projPtName;
  499. vn_comp_pT = vn_comp_pT_cent.Projection({name});
  500. vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  501. name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
  502. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
  503. graph->SetName(name.c_str());
  504. name += ";p_{T}, GeV/c;v_{1}";
  505. graph->SetTitle(name.c_str());
  506. graph->Write();
  507. vn_comp_y_obs = Merge(-1.*u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
  508. vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
  509. vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  510. name = trackBranchName+projYName;
  511. vn_comp_y = vn_comp_y_cent.Projection({name});
  512. vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  513. name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
  514. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
  515. graph->SetName(name.c_str());
  516. name += ";y;v_{1}";
  517. graph->SetTitle(name.c_str());
  518. graph->Write();
  519. }
  520. }
  521. }
  522. foGraphs->mkdir("v2Components_EP");
  523. foGraphs->cd("v2Components_EP");
  524. for (int ipid=0; ipid<npid; ipid++)
  525. {
  526. for (int iproj1=0; iproj1<nProj; iproj1++)
  527. {
  528. for (int iproj2=0; iproj2<nProj; iproj2++)
  529. {
  530. vn_comp_pT_obs = Merge(u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
  531. vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
  532. vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  533. name = trackBranchName+projPtName;
  534. vn_comp_pT = vn_comp_pT_cent.Projection({name});
  535. vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  536. name = "v1_pT_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
  537. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
  538. graph->SetName(name.c_str());
  539. name += ";p_{T}, GeV/c;v_{1}";
  540. graph->SetTitle(name.c_str());
  541. graph->Write();
  542. vn_comp_y_obs = Merge(-1.*u1_fhcalEP1[iproj1][iproj2][ipid], u1_fhcalEP2[iproj1][iproj2][ipid]);
  543. vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalEP[iproj1][iproj2]);
  544. vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  545. name = trackBranchName+projYName;
  546. vn_comp_y = vn_comp_y_cent.Projection({name});
  547. vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  548. name = "v1_y_fhcalEP_" + pidname[ipid] + "_" + corrname +"_"+projv1EP[iproj1]+projv1EP[iproj2];
  549. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
  550. graph->SetName(name.c_str());
  551. name += ";y;v_{1}";
  552. graph->SetTitle(name.c_str());
  553. graph->Write();
  554. }
  555. }
  556. }
  557. foGraphs->mkdir("v1Components_SP");
  558. foGraphs->cd("v1Components_SP");
  559. for (int ipid=0; ipid<npid; ipid++)
  560. {
  561. for (int iproj1=0; iproj1<nProj; iproj1++)
  562. {
  563. for (int iproj2=0; iproj2<nProj; iproj2++)
  564. {
  565. vn_comp_pT_obs = Merge(u1_fhcalSP1[iproj1][iproj2][ipid], u1_fhcalSP2[iproj1][iproj2][ipid]);
  566. vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q1_fhcalSP[iproj1][iproj2]);
  567. vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  568. name = trackBranchName+projPtName;
  569. vn_comp_pT = vn_comp_pT_cent.Projection({name});
  570. vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  571. name = "v1_pT_fhcalSP_" + pidname[ipid] + "_" + corrname +"_"+projv1SP[iproj1]+projv1SP[iproj2];
  572. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
  573. graph->SetName(name.c_str());
  574. name += ";p_{T}, GeV/c;v_{1}";
  575. graph->SetTitle(name.c_str());
  576. graph->Write();
  577. vn_comp_y_obs = Merge(-1.*u1_fhcalSP1[iproj1][iproj2][ipid], u1_fhcalSP2[iproj1][iproj2][ipid]);
  578. vn_comp_y_div = vn_comp_y_obs / Sqrt(q1_fhcalSP[iproj1][iproj2]);
  579. vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  580. name = trackBranchName+projYName;
  581. vn_comp_y = vn_comp_y_cent.Projection({name});
  582. vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  583. name = "v1_y_fhcalSP_" + pidname[ipid] + "_" + corrname +"_"+projv1SP[iproj1]+projv1SP[iproj2];
  584. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
  585. graph->SetName(name.c_str());
  586. name += ";y;v_{1}";
  587. graph->SetTitle(name.c_str());
  588. graph->Write();
  589. }
  590. }
  591. }
  592. foGraphs->mkdir("v2Components_SP");
  593. foGraphs->cd("v2Components_SP");
  594. for (int ipid=0; ipid<npid; ipid++)
  595. {
  596. for (int iproj1=0; iproj1<nProj; iproj1++)
  597. {
  598. for (int iproj2=0; iproj2<nProj; iproj2++)
  599. {
  600. vn_comp_pT_obs = Merge(u2_tpcSP1[iproj1][iproj2][ipid], u2_tpcSP2[iproj1][iproj2][ipid]);
  601. vn_comp_pT_div = vn_comp_pT_obs / Sqrt(q2_tpcSP[iproj1][iproj2]);
  602. vn_comp_pT_cent = vn_comp_pT_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  603. name = trackBranchName+projPtName;
  604. vn_comp_pT = vn_comp_pT_cent.Projection({name});
  605. vn_comp_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  606. name = "v2_pT_tpcSP_" + pidname[ipid] + "_" + corrname +"_"+projv2SP[iproj1]+projv2SP[iproj2];
  607. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_pT);
  608. graph->SetName(name.c_str());
  609. name += ";p_{T}, GeV/c;v_{2}";
  610. graph->SetTitle(name.c_str());
  611. graph->Write();
  612. vn_comp_y_obs = Merge(u2_tpcSP1[iproj1][iproj2][ipid], u2_tpcSP2[iproj1][iproj2][ipid]);
  613. vn_comp_y_div = vn_comp_y_obs / Sqrt(q2_tpcSP[iproj1][iproj2]);
  614. vn_comp_y_cent = vn_comp_y_div.Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  615. name = trackBranchName+projYName;
  616. vn_comp_y = vn_comp_y_cent.Projection({name});
  617. vn_comp_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  618. name = "v2_y_tpcSP_" + pidname[ipid] + "_" + corrname +"_"+projv2SP[iproj1]+projv2SP[iproj2];
  619. graph = Qn::DataContainerHelper::ToTGraph(vn_comp_y);
  620. graph->SetName(name.c_str());
  621. name += ";y;v_{2}";
  622. graph->SetTitle(name.c_str());
  623. graph->Write();
  624. }
  625. }
  626. }
  627. // uQ correlations
  628. foGraphs->mkdir("uQ");
  629. foGraphs->cd("uQ");
  630. for (int ipid=0; ipid<npid; ipid++)
  631. {
  632. for (int iproj1=0; iproj1<nProj; iproj1++)
  633. {
  634. for (int iproj2=0; iproj2<nProj; iproj2++)
  635. {
  636. uq_cent = u2_tpcEP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  637. name = trackBranchName+projPtName;
  638. uq_pT = uq_cent.Projection({name});
  639. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  640. name = "uQ_pT_tpcEP1_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
  641. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  642. graph->SetName(name.c_str());
  643. name += ";p_{T}, GeV/c;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
  644. graph->SetTitle(name.c_str());
  645. graph->Write();
  646. name = trackBranchName+projYName;
  647. uq_y = uq_cent.Projection({name});
  648. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  649. name = "uQ_y_tpcEP1_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
  650. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  651. graph->SetName(name.c_str());
  652. name += ";y;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
  653. graph->SetTitle(name.c_str());
  654. graph->Write();
  655. uq_cent = u2_tpcEP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  656. name = trackBranchName+projPtName;
  657. uq_pT = uq_cent.Projection({name});
  658. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  659. name = "uQ_pT_tpcEP2_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
  660. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  661. graph->SetName(name.c_str());
  662. name += ";p_{T}, GeV/c;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
  663. graph->SetTitle(name.c_str());
  664. graph->Write();
  665. name = trackBranchName+projYName;
  666. uq_y = uq_cent.Projection({name});
  667. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  668. name = "uQ_y_tpcEP2_" + pidname[ipid] + "_" + corrname + "_"+projv2EP[iproj1]+projv2EP[iproj2];
  669. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  670. graph->SetName(name.c_str());
  671. name += ";y;<"+projv2EP[iproj1]+projv2EP[iproj2]+">";
  672. graph->SetTitle(name.c_str());
  673. graph->Write();
  674. uq_cent = u2_tpcSP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  675. name = trackBranchName+projPtName;
  676. uq_pT = uq_cent.Projection({name});
  677. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  678. name = "uQ_pT_tpcSP1_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
  679. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  680. graph->SetName(name.c_str());
  681. name += ";p_{T}, GeV/c;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
  682. graph->SetTitle(name.c_str());
  683. graph->Write();
  684. name = trackBranchName+projYName;
  685. uq_y = uq_cent.Projection({name});
  686. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  687. name = "uQ_y_tpcSP1_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
  688. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  689. graph->SetName(name.c_str());
  690. name += ";y;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
  691. graph->SetTitle(name.c_str());
  692. graph->Write();
  693. uq_cent = u2_tpcSP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  694. name = trackBranchName+projPtName;
  695. uq_pT = uq_cent.Projection({name});
  696. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  697. name = "uQ_pT_tpcSP2_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
  698. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  699. graph->SetName(name.c_str());
  700. name += ";p_{T}, GeV/c;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
  701. graph->SetTitle(name.c_str());
  702. graph->Write();
  703. name = trackBranchName+projYName;
  704. uq_y = uq_cent.Projection({name});
  705. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  706. name = "uQ_y_tpcSP2_" + pidname[ipid] + "_" + corrname + "_"+projv2SP[iproj1]+projv2SP[iproj2];
  707. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  708. graph->SetName(name.c_str());
  709. name += ";y;<"+projv2SP[iproj1]+projv2SP[iproj2]+">";
  710. graph->SetTitle(name.c_str());
  711. graph->Write();
  712. uq_cent = u1_fhcalEP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  713. name = trackBranchName+projPtName;
  714. uq_pT = uq_cent.Projection({name});
  715. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  716. name = "uQ_pT_fhcalEP1_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
  717. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  718. graph->SetName(name.c_str());
  719. name += ";p_{T}, GeV/c;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
  720. graph->SetTitle(name.c_str());
  721. graph->Write();
  722. name = trackBranchName+projYName;
  723. uq_y = uq_cent.Projection({name});
  724. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  725. name = "uQ_y_fhcalEP1_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
  726. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  727. graph->SetName(name.c_str());
  728. name += ";y;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
  729. graph->SetTitle(name.c_str());
  730. graph->Write();
  731. uq_cent = u1_fhcalEP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  732. name = trackBranchName+projPtName;
  733. uq_pT = uq_cent.Projection({name});
  734. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  735. name = "uQ_pT_fhcalEP2_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
  736. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  737. graph->SetName(name.c_str());
  738. name += ";p_{T}, GeV/c;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
  739. graph->SetTitle(name.c_str());
  740. graph->Write();
  741. name = trackBranchName+projYName;
  742. uq_y = uq_cent.Projection({name});
  743. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  744. name = "uQ_y_fhcalEP2_" + pidname[ipid] + "_" + corrname + "_"+projv1EP[iproj1]+projv1EP[iproj2];
  745. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  746. graph->SetName(name.c_str());
  747. name += ";y;<"+projv1EP[iproj1]+projv1EP[iproj2]+">";
  748. graph->SetTitle(name.c_str());
  749. graph->Write();
  750. uq_cent = u1_fhcalSP1[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  751. name = trackBranchName+projPtName;
  752. uq_pT = uq_cent.Projection({name});
  753. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  754. name = "uQ_pT_fhcalSP1_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
  755. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  756. graph->SetName(name.c_str());
  757. name += ";p_{T}, GeV/c;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
  758. graph->SetTitle(name.c_str());
  759. graph->Write();
  760. name = trackBranchName+projYName;
  761. uq_y = uq_cent.Projection({name});
  762. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  763. name = "uQ_y_fhcalSP1_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
  764. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  765. graph->SetName(name.c_str());
  766. name += ";y;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
  767. graph->SetTitle(name.c_str());
  768. graph->Write();
  769. uq_cent = u1_fhcalSP2[iproj1][iproj2][ipid].Rebin({"McEvent_B", 1, bcent_range.first, bcent_range.second});
  770. name = trackBranchName+projPtName;
  771. uq_pT = uq_cent.Projection({name});
  772. uq_pT.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  773. name = "uQ_pT_fhcalSP2_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
  774. graph = Qn::DataContainerHelper::ToTGraph(uq_pT);
  775. graph->SetName(name.c_str());
  776. name += ";p_{T}, GeV/c;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
  777. graph->SetTitle(name.c_str());
  778. graph->Write();
  779. name = trackBranchName+projYName;
  780. uq_y = uq_cent.Projection({name});
  781. uq_y.SetErrors(Qn::Stat::ErrorType::BOOTSTRAP);
  782. name = "uQ_y_fhcalSP2_" + pidname[ipid] + "_" + corrname + "_"+projv1SP[iproj1]+projv1SP[iproj2];
  783. graph = Qn::DataContainerHelper::ToTGraph(uq_y);
  784. graph->SetName(name.c_str());
  785. name += ";y;<"+projv1SP[iproj1]+projv1SP[iproj2]+">";
  786. graph->SetTitle(name.c_str());
  787. graph->Write();
  788. }
  789. }
  790. }
  791. foGraphs->Close();
  792. }