// prepare final histograms for human analysis
// of reconstructed info
// 01/05/2014 A.Karavdina

void AnaAna(TString file_sig_name= "/panda/scrut14/macro/scrut/fastsim_hc/hc_rho0pi0_fast.root_taskana.root", TString path="/panda/scrut14/macro/scrut/fastsim_hc/"){
  // cout<<"path="<<path<<endl;
  TString fout_pdf_o = path+"/plotRes.pdf(";
  TString fout_pdf_c = path+"/plotRes.pdf)";
  TString fout_root = path+"/plotRes.root";
  TFile *finSIG = new TFile(file_sig_name.Data(),"READ");
  TTree *tree_sig_raw = (TTree*)finSIG->Get("ntp3");
  TTree *tree_sig_4Cfit = (TTree*)finSIG->Get("ntp4");
  TTree *tree_sig_rho = (TTree*)finSIG->Get("ntp2");
  TTree *tree_sig_pi0 = (TTree*)finSIG->Get("ntp1");
  TTree *tree_sig_eta = (TTree*)finSIG->Get("ntp5");
  TTree *tree_sig_pipl = (TTree*)finSIG->Get("ntp6");
  TTree *tree_sig_pimn = (TTree*)finSIG->Get("ntp7");
  TTree *tree_sig_mc = (TTree*)finSIG->Get("nmc");
  TH1 *hchi2 = (TH1D*)finSIG->Get("hchi2_4C");

  TH1 *hhcTh_mc = new TH1D("hhcTh_mc","Efficiency (MC); #theta, ^#circ",2e2,0,2);
  tree_sig_mc->Project("hhcTh_mc","(tht*TMath::RadToDeg())","pdg==10443");
  TH1 *hhcTh_raw = new TH1D("hhcTh_raw","Efficiency (REC, raw); #theta, ^#circ",2e2,0,2);
  tree_sig_raw->Project("hhcTh_raw","(hctht*TMath::RadToDeg())");
  //  hhcTh_raw->Divide(hhcTh_mc);
  TH1 *hhcPh_mc = new TH1D("hhcPh_mc","Efficiency (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hhcPh_mc","(phi*TMath::RadToDeg())","pdg==10443");
  TH1 *hhcPh_raw = new TH1D("hhcPh_raw","Efficiency (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_raw->Project("hhcPh_raw","(hcphi*TMath::RadToDeg())");
  //  hhcPh_raw->Divide(hhcPh_mc);

  TH1 *hrhoTh_mc = new TH1D("hrhoTh_mc","Efficiency #pi^{+}#pi^{-} (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hrhoTh_mc","(tht*TMath::RadToDeg())","pdg==223 || pdg==113");
  TH1 *hrhoTh_raw = new TH1D("hrhoTh_raw","Efficiency #pi^{+}#pi^{-} (REC, raw); #theta, ^#circ",2e2,0,180);
  TH1 *hrhoTh_raw2 = new TH1D("hrhoTh_raw2","Efficiency #pi^{+}#pi^{-} (REC, raw); #theta, ^#circ",2e2,0,180);
  // TH1 *hrhoTh_fit = new TH1D("hrhoTh_fit","Efficiency #pi^{+}#pi^{-} (REC, 4C fit); #theta, ^#circ",2e2,0,180);
  tree_sig_rho->Project("hrhoTh_raw","(rhotht*TMath::RadToDeg())");
  tree_sig_rho->Project("hrhoTh_raw2","(rhotht*TMath::RadToDeg())");
  //  tree_sig_4Cfit->Project("hrhoTh_fit","(hcfitd0tht*TMath::RadToDeg())");
  hrhoTh_raw2->Divide(hrhoTh_mc);
  //hrhoTh_fit->Divide(hrhoTh_mc);

  TH1 *hrhoPh_mc = new TH1D("hrhoPh_mc","Efficiency #pi^{+}#pi^{-} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hrhoPh_mc","(phi*TMath::RadToDeg())","pdg==223 || pdg==113");
  TH1 *hrhoPh_raw = new TH1D("hrhoPh_raw","Efficiency #pi^{+}#pi^{-} (REC, raw);#phi, ^#circ",2e2,-180,180);
  TH1 *hrhoPh_raw2 = new TH1D("hrhoPh_raw2","Efficiency #pi^{+}#pi^{-} (REC, raw);#phi, ^#circ",2e2,-180,180);
  //  TH1 *hrhoPh_fit = new TH1D("hrhoPh_fit","Efficiency #pi^{+}#pi^{-} (REC, 4C fit);#phi, ^#circ",2e2,-180,180);
  //  tree_sig_4Cfit->Project("hrhoPh_fit","(hcfitd0phi*TMath::RadToDeg())");
  tree_sig_rho->Project("hrhoPh_raw","(rhophi*TMath::RadToDeg())");
  tree_sig_rho->Project("hrhoPh_raw2","(rhophi*TMath::RadToDeg())");
  hrhoPh_raw2->Divide(hrhoPh_mc);
  //  hrhoPh_fit->Divide(hrhoPh_mc);

  TH2D *hrhoThPh_mc = new TH2D("hrhoThPh_mc","Efficiency #pi^{+}#pi^{-} (MC);#phi, ^#circ",2e2,-180,180,2e2,0,180);
  tree_sig_mc->Project("hrhoThPh_mc","(tht*TMath::RadToDeg()):(phi*TMath::RadToDeg())","pdg==223 || pdg==113");
 TH2D *hrhoThPh_raw = new TH2D("hrhoThPh_raw","Efficiency #pi^{+}#pi^{-} (MC);#phi, ^#circ",2e2,-180,180,2e2,0,180);
 tree_sig_rho->Project("hrhoThPh_raw","(rhotht*TMath::RadToDeg()):(rhophi*TMath::RadToDeg())");
 TH2D *hrhoThPh_raw2 = new TH2D("hrhoThPh_raw2","Efficiency #pi^{+}#pi^{-} (MC);#phi, ^#circ",2e2,-180,180,2e2,0,180);
 tree_sig_rho->Project("hrhoThPh_raw2","(rhotht*TMath::RadToDeg()):(rhophi*TMath::RadToDeg())");
 hrhoThPh_raw2->Divide(hrhoThPh_mc);

 TH1 *hetaTh_mc = new TH1D("hetaTh_mc","Efficiency #pi^{+}#pi^{-}#pi^{0} (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hetaTh_mc","(tht*TMath::RadToDeg())","pdg==221");
  TH1 *hetaTh_raw = new TH1D("hetaTh_raw","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_eta->Project("hetaTh_raw","(etatht*TMath::RadToDeg())");
  TH1 *hetaTh_raw2 = new TH1D("hetaTh_raw2","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_eta->Project("hetaTh_raw2","(etatht*TMath::RadToDeg())");
  //  TH1 *hetaTh_fit = new TH1D("hetaTh_fit","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, 4C fit); #theta, ^#circ",2e2,0,180);
  //  tree_sig_4Cfit->Project("hetaTh_fit","(hcfitd1tht*TMath::RadToDeg())");
  hetaTh_raw2->Divide(hetaTh_mc);
  //  hetaTh_fit->Divide(hetaTh_mc);

  TH1 *hetaPh_mc = new TH1D("hetaPh_mc","Efficiency #pi^{+}#pi^{-}#pi^{0} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hetaPh_mc","(phi*TMath::RadToDeg())","pdg==221");
  TH1 *hetaPh_raw = new TH1D("hetaPh_raw","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_eta->Project("hetaPh_raw","(etaphi*TMath::RadToDeg())");
  TH1 *hetaPh_raw2 = new TH1D("hetaPh_raw2","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_eta->Project("hetaPh_raw2","(etaphi*TMath::RadToDeg())");

  //TH1 *hetaPh_fit = new TH1D("hetaPh_fit","Efficiency #pi^{+}#pi^{-}#pi^{0} (REC, 4C fit);#phi, ^#circ",2e2,-180,180);
  // tree_sig_4Cfit->Project("hetaPh_fit","(hcfitd1phi*TMath::RadToDeg())");
  // tree_sig_eta->Project("hetaPh_raw","(etaphi*TMath::RadToDeg())");
  hetaPh_raw2->Divide(hetaPh_mc);
 // hetaPh_fit->Divide(hetaPh_mc);
 TH2 *hetaThPh_mc = new TH2D("hetaThPh_mc","Efficiency #pi^{+}#pi^{-}#pi^{0} (MC);#phi, #circ; #theta, #circ",2e2,-180,180,2e2,0,180);
 TH2 *hetaThPh_raw = new TH2D("hetaThPh_raw","Efficiency #pi^{+}#pi^{-}#pi^{0} (MC);#phi, #circ; #theta, #circ",2e2,-180,180,2e2,0,180);
 TH2 *hetaThPh_raw2 = new TH2D("hetaThPh_raw2","Efficiency #pi^{+}#pi^{-}#pi^{0} (MC);#phi, #circ; #theta, #circ",2e2,-180,180,2e2,0,180);
 tree_sig_mc->Project("hetaThPh_mc","(tht*TMath::RadToDeg()):(phi*TMath::RadToDeg())","pdg==221");
 tree_sig_eta->Project("hetaThPh_raw","(etatht*TMath::RadToDeg()):(etaphi*TMath::RadToDeg())");
 tree_sig_eta->Project("hetaThPh_raw2","(etatht*TMath::RadToDeg()):(etaphi*TMath::RadToDeg())");
 hetaThPh_raw2->Divide(hetaThPh_mc);

  TH1 *hgammaTh_mc = new TH1D("hgammaTh_mc","Efficiency #gamma (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hgammaTh_mc","(tht*TMath::RadToDeg())","pdg==22");
  TH1 *hgammaTh_raw = new TH1D("hgammaTh_raw","Efficiency #gamma (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pi0->Project("hgammaTh_raw","(pi0g1tht*TMath::RadToDeg())");
  TH1 *hgammaTh_raw2 = new TH1D("hgammaTh_raw2","Efficiency #gamma (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pi0->Project("hgammaTh_raw2","(pi0g2tht*TMath::RadToDeg())");
  hgammaTh_raw->Add(hgammaTh_raw2);
 TH1 *hgammaPh_mc = new TH1D("hgammaPh_mc","Efficiency #pi^{0} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hgammaPh_mc","(phi*TMath::RadToDeg())","pdg==111");
  TH1 *hgammaPh_raw = new TH1D("hgammaPh_raw","Efficiency #pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pi0->Project("hgammaPh_raw","(pi0g1phi*TMath::RadToDeg())");
  TH1 *hgammaPh_raw2 = new TH1D("hgammaPh_raw2","Efficiency #pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pi0->Project("hgammaPh_raw2","(pi0g2phi*TMath::RadToDeg())");
  hgammaPh_raw->Add(hgammaPh_raw2);

  TH1 *hpi0Th_mc = new TH1D("hpi0Th_mc","Efficiency #pi^{0} (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hpi0Th_mc","(tht*TMath::RadToDeg())","pdg==111");
  TH1 *hpi0Th_raw = new TH1D("hpi0Th_raw","Efficiency #pi^{0} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pi0->Project("hpi0Th_raw","(pi0tht*TMath::RadToDeg())");
  TH1 *hpi0Th_raw2 = new TH1D("hpi0Th_raw2","Efficiency #pi^{0} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pi0->Project("hpi0Th_raw2","(pi0tht*TMath::RadToDeg())");
  hpi0Th_raw2->Divide(hpi0Th_mc);
  //  TH1 *hpi0Th_fit = new TH1D("hpi0Th_fit","Efficiency #pi^{0} (REC, 4C fit);#theta, ^#circ",2e2,0,180);
  //  tree_sig_4Cfit->Project("hpi0Th_fit","(hcfitd0tht*TMath::RadToDeg())");
  //  hpi0Th_fit->Divide(hpi0Th_mc);

  TH1 *hpi0Ph_mc = new TH1D("hpi0Ph_mc","Efficiency #pi^{0} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hpi0Ph_mc","(phi*TMath::RadToDeg())","pdg==111");
  TH1 *hpi0Ph_raw = new TH1D("hpi0Ph_raw","Efficiency #pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pi0->Project("hpi0Ph_raw","(pi0phi*TMath::RadToDeg())");
  TH1 *hpi0Ph_raw2 = new TH1D("hpi0Ph_raw2","Efficiency #pi^{0} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pi0->Project("hpi0Ph_raw2","(pi0phi*TMath::RadToDeg())");
  hpi0Ph_raw2->Divide(hpi0Ph_mc);
  //  TH1 *hpi0Ph_fit = new TH1D("hpi0Ph_fit","Efficiency #pi^{0} (REC, 4C fit);#phi, ^#circ",2e2,-180,180);
  //  tree_sig_4Cfit->Project("hpi0Ph_fit","(hcfitd0phi*TMath::RadToDeg())");
  //  hpi0Ph_fit->Divide(hpi0Ph_mc);

  TH2D *hpi0ThPh_mc = new TH2D("hpi0ThPh_mc","Efficiency #pi^{0} (MC);#phi, #circ;#theta,#circ",2e2,-180,180,2e2,0,180);
  tree_sig_mc->Project("hpi0ThPh_mc","(tht*TMath::RadToDeg()):(phi*TMath::RadToDeg())","pdg==111");
  TH2D *hpi0ThPh_raw = new TH2D("hpi0ThPh_raw","Efficiency #pi^{0} (MC);#phi, #circ;#theta,#circ",2e2,-180,180,2e2,0,180);
  tree_sig_pi0->Project("hpi0ThPh_raw","(pi0tht*TMath::RadToDeg()):(pi0phi*TMath::RadToDeg())");
  TH2D *hpi0ThPh_raw2 = new TH2D("hpi0ThPh_raw2","Efficiency #pi^{0} (MC);#phi, #circ;#theta,#circ",2e2,-180,180,2e2,0,180);
  tree_sig_pi0->Project("hpi0ThPh_raw2","(pi0tht*TMath::RadToDeg()):(pi0phi*TMath::RadToDeg())");
  hpi0ThPh_raw2->Divide(hpi0ThPh_mc);

 TH1 *hpiplPh_mc = new TH1D("hpiplPh_mc","Efficiency #pi^{+} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hpiplPh_mc","(phi*TMath::RadToDeg())","pdg==211");
  TH1 *hpiplPh_raw = new TH1D("hpiplPh_raw","Efficiency #pi^{+} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pipl->Project("hpiplPh_raw","(piplphi*TMath::RadToDeg())");
  TH1 *hpiplPh_raw2 = new TH1D("hpiplPh_raw2","Efficiency #pi^{+} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pipl->Project("hpiplPh_raw2","(piplphi*TMath::RadToDeg())");
  hpiplPh_raw2->Divide(hpiplPh_mc);

 TH1 *hpimnPh_mc = new TH1D("hpimnPh_mc","Efficiency #pi^{-} (MC);#phi, ^#circ",2e2,-180,180);
  tree_sig_mc->Project("hpimnPh_mc","(phi*TMath::RadToDeg())","pdg==-211");
  TH1 *hpimnPh_raw = new TH1D("hpimnPh_raw","Efficiency #pi^{-} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pimn->Project("hpimnPh_raw","(pimnphi*TMath::RadToDeg())");
  TH1 *hpimnPh_raw2 = new TH1D("hpimnPh_raw2","Efficiency #pi^{-} (REC, raw);#phi, ^#circ",2e2,-180,180);
  tree_sig_pimn->Project("hpimnPh_raw2","(pimnphi*TMath::RadToDeg())");
  hpimnPh_raw2->Divide(hpimnPh_mc);


 TH1 *hpiplTh_mc = new TH1D("hpiplTh_mc","Efficiency #pi^{+} (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hpiplTh_mc","(tht*TMath::RadToDeg())","pdg==211");
  TH1 *hpiplTh_raw = new TH1D("hpiplTh_raw","Efficiency #pi^{+} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pipl->Project("hpiplTh_raw","(pipltht*TMath::RadToDeg())");
  TH1 *hpiplTh_raw2 = new TH1D("hpiplTh_raw2","Efficiency #pi^{+} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pipl->Project("hpiplTh_raw2","(pipltht*TMath::RadToDeg())");
  hpiplTh_raw2->Divide(hpiplTh_mc);

  TH1 *hpimnTh_mc = new TH1D("hpimnTh_mc","Efficiency #pi^{-} (MC); #theta, ^#circ",2e2,0,180);
  tree_sig_mc->Project("hpimnTh_mc","(tht*TMath::RadToDeg())","pdg==-211");
  TH1 *hpimnTh_raw = new TH1D("hpimnTh_raw","Efficiency #pi^{-} (REC, raw); #theta, ^#circ",2e2,0,180);
  TH1 *hpimnTh_raw2= new TH1D("hpimnTh_raw2","Efficiency #pi^{-} (REC, raw); #theta, ^#circ",2e2,0,180);
  tree_sig_pimn->Project("hpimnTh_raw","(pimntht*TMath::RadToDeg())");
  tree_sig_pimn->Project("hpimnTh_raw2","(pimntht*TMath::RadToDeg())");
  hpimnTh_raw2->Divide(hpimnTh_mc);

  TH2 *hpiplThPh_mc = new TH2D("hpiplThPh_mc","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  TH2 *hpimnThPh_mc = new TH2D("hpimnThPh_mc","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  tree_sig_mc->Project("hpiplThPh_mc","(tht*TMath::RadToDeg()):(phi*TMath::RadToDeg())","pdg==211");
  tree_sig_mc->Project("hpimnThPh_mc","(tht*TMath::RadToDeg()):(phi*TMath::RadToDeg())","pdg==-211");
  TH2 *hpiplThPh_raw = new TH2D("hpiplThPh_raw","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  TH2 *hpimnThPh_raw = new TH2D("hpimnThPh_raw","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  TH2 *hpiplThPh_raw2 = new TH2D("hpiplThPh_raw2","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  TH2 *hpimnThPh_raw2 = new TH2D("hpimnThPh_raw2","Efficiency #pi^{+} (MC);#phi, #circ;#theta, #circ",2e2,-180,180,2e2,0,180);
  tree_sig_pipl->Project("hpiplThPh_raw","(pipltht*TMath::RadToDeg()):(piplphi*TMath::RadToDeg())");
  tree_sig_pipl->Project("hpiplThPh_raw2","(pipltht*TMath::RadToDeg()):(piplphi*TMath::RadToDeg())");
  tree_sig_pimn->Project("hpimnThPh_raw","(pimntht*TMath::RadToDeg()):(pimnphi*TMath::RadToDeg())");
  tree_sig_pimn->Project("hpimnThPh_raw2","(pimntht*TMath::RadToDeg()):(pimnphi*TMath::RadToDeg())");
  hpiplThPh_raw2->Divide(hpiplThPh_mc);
  hpimnThPh_raw2->Divide(hpimnThPh_mc);

  TH1 *hhcmraw = new TH1D("hhcmraw","Signal ;m, GeV",3e2,3.,4.);
  tree_sig_raw->Project("hhcmraw","hcm");
  TH2 *hDalitz = new TH2D("hDalitz","; m^{2}_{#pi^{+}#pi^{-}#pi^{0}}, (GeV/c^{2})^{2};m^{2}_{#pi^{+}#pi^{-}}, (GeV/c^{2})^{2}",3e2,0.,1.,3e2,0.,5.0);
  tree_sig_raw->Project("hDalitz","hcd0m*hcd0m:hcd1m*hcd1m");
  TH1 *hhcmfit = new TH1D("hhcmfit","Signal ;m, GeV",1e4,3.,4.);
  tree_sig_4Cfit->Project("hhcmfit","hcfitm");

  TH1 *hrhom = new TH1D("hrhom","Signal ; m, GeV",3e2,0.,5.0);
  tree_sig_rho->Project("hrhom","rhom");
  TH2 *hPeyrou_rho = new TH2D("hPeyrou_rho","Peyrou Diagram: #rho (SIG); p_{z}, GeV/c ; p_{#perp} , GeV/c",3e2,-1,7,3e2,0,2);
  tree_sig_rho->Project("hPeyrou_rho","rhopt:rhopz");
  TH1 *hpi0m = new TH1D("hpi0m","Signal ; m, GeV",3e2,0,0.3);
  tree_sig_pi0->Project("hpi0m","pi0m");
  TH2 *hPeyrou_pi0 = new TH2D("hPeyrou_pi0","Peyrou Diagram: #pi^{0} (SIG); p_{z}, GeV/c; p_{#perp} , GeV/c",3e2,-1,7,3e2,0,2);
  tree_sig_pi0->Project("hPeyrou_pi0","pi0pt:pi0pz");
  

  TH2 *hPeyrou_hc = new TH2D("hPeyrou_hc","Peyrou Diagram: h_{c} (SIG); p_{z}, GeV/c; p_{#perp} , GeV/c",3e2,5,7,3e2,0,0.2);
  tree_sig_raw->Project("hPeyrou_hc","hcpt:hcpz");
  //  TH2 *hPeyrou_hc_fit = new TH2D("hPeyrou_hc_fit","Peyrou Diagram: h_{c} (SIG, 4C fit); p_{z}, GeV/c; p_{#perp}, GeV/c",3e2,5,7,3e2,0,4e-11);
  //  tree_sig_4Cfit->Project("hPeyrou_hc_fit","hcfitpt:hcfitpz");

  TH1 *hetam = new TH1D("hetam","Signal ; m, GeV",3e2,0,5.0);
  tree_sig_eta->Project("hetam","etam");
  TH2 *hPeyrou_eta = new TH2D("hPeyrou_eta","Peyrou Diagram: #eta (SIG); p_{z}, GeV/c; p_{#perp}, GeV/c",3e2,-1,7,3e2,0,2);
  tree_sig_eta->Project("hPeyrou_eta","etapt:etapz");
 
  TCanvas c1;
  c1.SetLogy(1);
  c1.Divide(3,2);
  c1.cd(1);
  hpi0m->Draw();
  c1.cd(2);
  hetam->Draw();
  c1.cd(3);
  hrhom->Draw();
  c1.cd(4);
  hhcmraw->Draw();
  c1.cd(5);
  hhcmfit->Draw();
  c1.cd(6);
  hchi2->Draw();
  c1.Print(fout_pdf_o);

  TCanvas c2;
  c2.Divide(2,2);
  c2.cd(1);
  //hPeyrou_pi0->Draw("colz");
  hPeyrou_rho->Draw("colz");
  c2.cd(2);
  hPeyrou_eta->Draw("colz");
  c2.cd(3);
  hPeyrou_hc->Draw("colz");
  c2.cd(4);
  hDalitz->Draw("colz");
  //  hPeyrou_hc_fit->Draw("colz");
  c2.Print(fout_pdf_o);

  TCanvas c3;
  c3.Divide(2,3);
  c3.cd(1);
  hpi0Th_raw2->SetMinimum(0);
  hpi0Th_raw2->SetMaximum(1.2);
  hpi0Th_raw2->Draw();
  c3.cd(2);
  hpi0Ph_raw2->SetMinimum(0);
  hpi0Ph_raw2->SetMaximum(1.2);
  hpi0Ph_raw2->Draw();
  c3.cd(3);
  hpiplTh_raw2->SetMinimum(0);
  hpiplTh_raw2->SetMaximum(1.2);
  hpiplTh_raw2->Draw();
  c3.cd(4);
  hpiplPh_raw2->SetMinimum(0);
  hpiplPh_raw2->SetMaximum(1.2);
  hpiplPh_raw2->Draw();
  c3.cd(5);
  hpimnTh_raw2->SetMinimum(0);
  hpimnTh_raw2->SetMaximum(1.2);
  hpimnTh_raw2->Draw();
  c3.cd(6);
  hpimnPh_raw2->SetMinimum(0);
  hpimnPh_raw2->SetMaximum(1.2);
  hpimnPh_raw2->Draw();
  c3.Print(fout_pdf_o);

  TCanvas c5;
  c5.Divide(3,2);
  c5.cd(1);
  hpi0ThPh_raw2->Draw("colz");
  c5.cd(2);
  hpiplThPh_raw2->Draw("colz");
  c5.cd(3);
  hpimnThPh_raw2->Draw("colz");
  c5.cd(4);
  hrhoThPh_raw2->Draw("colz");
  c5.cd(5);
  hetaThPh_raw2->Draw("colz");
  c5.Print(fout_pdf_o);

 TCanvas c4;
 c4.Divide(2,2);
 c4.cd(1);
 hetaTh_raw2->SetMinimum(0);
 hetaTh_raw2->SetMaximum(1.2);
  hetaTh_raw2->Draw();
  c4.cd(2);
  hetaPh_raw2->SetMinimum(0);
  hetaPh_raw2->SetMaximum(1.2);
  hetaPh_raw2->Draw();
  c4.cd(3);
  hrhoTh_raw2->SetMinimum(0);
  hrhoTh_raw2->SetMaximum(1.2);
  hrhoTh_raw2->Draw();
  c4.cd(4);
  hrhoPh_raw2->SetMinimum(0);
  hrhoPh_raw2->SetMaximum(1.2);
  hrhoPh_raw2->Draw();
  c4.Print(fout_pdf_c);


  // Compare reconstructed info with true one
  //  TH2 *hPeyrou_rho_tr = new TH2D("hPeyrou_rho_tr","Peyrou Diagram: #rho (SIG); p_{z}, GeV/c ; p_{#perp} , GeV/c",3e2,-1,7,3e2,0,2);
  TH2 *hPeyrou_rho_tr = new TH2D("hPeyrou_rho_tr","Peyrou Diagram: #rho (SIG); p_{z}, GeV/c ; p_{#perp} , GeV/c",3e2,-1,7,3e2,0,2);
  tree_sig_rho->Project("hPeyrou_rho_tr","rhopt:rhopz","trrhom>0");
  //  TH2 *hPeyrou_eta_tr = new TH2D("hPeyrou_eta_tr","Peyrou Diagram: #eta (SIG); p_{z}, GeV/c; p_{#perp}, GeV/c",3e2,-1,7,3e2,0,2);
  TH2 *hPeyrou_eta_tr = new TH2D("hPeyrou_eta_tr","Peyrou Diagram: #eta (SIG); p_{z}, GeV/c; p_{#perp}, GeV/c",3e2,-1,7,3e2,0,2);
  tree_sig_eta->Project("hPeyrou_eta_tr","etapt:etapz","tretam>0");

  // TCanvas c4;
  // c4.Divide(2,2);
  // c4.cd(1);
  // //hPeyrou_pi0->Draw("colz");
  // hPeyrou_rho_tr->Draw("colz");
  // c4.cd(2);
  // hPeyrou_eta_tr->Draw("colz");
  // c4.Print(fout_pdf_c);


  //  TString out="plotRes.root";
  TFile *fout = new TFile(fout_root,"RECREATE");
  hpi0ThPh_raw->Write();
  hpi0ThPh_mc->Write();
  hrhoThPh_raw->Write();
  hrhoThPh_mc->Write();
  hpiplThPh_mc->Write();
  hpiplThPh_raw->Write();
  hpimnThPh_mc->Write();
  hpimnThPh_raw->Write();
  hetaThPh_raw->Write();
  hetaThPh_mc->Write();
   hpi0m->Write();
  hetam->Write();
  hrhom->Write();
  hhcmraw->Write();
  hhcmfit->Write();
  hPeyrou_pi0->Write();
  hPeyrou_rho->Write();
  hPeyrou_eta->Write();
  hPeyrou_hc->Write();
  //  hPeyrou_hc_fit->Write();
  hDalitz->Write();
  hhcTh_mc->Write();
  hhcTh_raw->Write();
  hhcPh_mc->Write();
  hhcPh_raw->Write();
  hrhoTh_mc->Write();
  hrhoPh_mc->Write();
  hrhoTh_raw->Write();
  hrhoPh_raw->Write();
  hetaTh_mc->Write();
  hetaPh_mc->Write();
  hetaTh_raw->Write();
  hetaPh_raw->Write();
  hpi0Th_mc->Write();
  hpi0Ph_mc->Write();
  hpi0Th_raw->Write();
  hpi0Ph_raw->Write();

  // hrhoTh_fit->Write();
  // hrhoPh_fit->Write(); 
  // hetaTh_fit->Write();
  // hetaPh_fit->Write();
  // hpi0Th_fit->Write();
  // hpi0Ph_fit->Write();
  hPeyrou_rho_tr->Write();
  hPeyrou_eta_tr->Write();
  hpiplPh_raw->Write();
  hpiplPh_mc->Write();
  hpimnPh_raw->Write();
  hpimnPh_mc->Write();
  hpiplTh_raw->Write();
  hpiplTh_mc->Write();
  hpimnTh_raw->Write();
  hpimnTh_mc->Write();
  hgammaTh_raw->Write();
  hgammaPh_raw->Write();
  hgammaTh_mc->Write();
  hgammaPh_mc->Write();
  // c1.Write();
  // c2.Write();
  // c3.Write();
  // c4.Write();

  hpi0Th_raw2->Write();
  hpi0Ph_raw2->Write();
  hpiplTh_raw2->Write();
  hpiplPh_raw2->Write();
  hpimnTh_raw2->Write();
  hpimnPh_raw2->Write();
  hetaTh_raw2->Write();
  hetaPh_raw2->Write();
  hrhoTh_raw2->Write();
  hrhoPh_raw2->Write();
  fout->Write();
  fout->Close();
}