[u/mrichter/AliRoot.git] / PWGPP / analysisQA / processJETriggerQA.C

/***************************************************
processJETriggerQA:
To procees JE Triggere QA wagon's output

****************************************************/


const Float_t ptmin =  0. ; //lower cutoff of jet pt spectrum

void processJETriggerQA(TString strFileIn = "AnalysisResults.root", 
			TString suftype="eps",
			Float_t jetR         = 0.2, 
			Float_t minTrkPT     = 0.15, 
			Float_t minClusterET = 0.3, 
			Int_t run            = 0, 
			const char* outfile   = "JETriggerQA_outfile.root"
			){

  gStyle->SetOptStat(0);
  gStyle->SetOptTitle(0);

  TString prefix = "fig_je_TriggerQA_";

  TFile * f1 = TFile::Open(strFileIn.Data());

  //Load histogram list
  TList *histList = 0x0;
  histList =  (TList*)f1->Get(Form("TriggerQA_Jet_AKTFullR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_Jet_AKTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TC/TriggerQA_Jet_AKTFullR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_Jet_AKTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TC", 
       TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000),
       TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000),   
       TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000),  
       TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000)));

  //Load a second histogram list in order to incorporate results from the Rho Task as well
  TList *histListRho = 0x0;
  histListRho = (TList*)f1->Get(Form("Rho_Jet_KTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TPC_histos", 
       TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000)));

  //---------------------------------------------------------------------------------------------------
  //       jet histograms
  //---------------------------------------------------------------------------------------------------

  const Int_t kJetType = 2;
  TString suffix [kJetType] =  {"Charged","Full"};
  TH3F *h3PtEtaPhiJet[kJetType];
  TH1F *hPtJet[kJetType];
  TH2F *hEtaPhiJet[kJetType];
  TH2F *hRhoCent;

  TH1F *hNEventSel = histList->FindObject("fhNEvents");
  Float_t nEvents  = hNEventSel->GetBinContent(2); 

  for(Int_t itype = 0; itype < kJetType; itype++){

     h3PtEtaPhiJet[itype] = (TH3F*)  histList->FindObject(Form("fh3PtEtaPhiJet%s",suffix[itype].Data()));
     if(! h3PtEtaPhiJet[itype]) continue;

     //jet pt spectra
     Int_t binMin = 1;
     if(ptmin>0.) binMin = h3PtEtaPhiJet[itype]->GetXaxis()->FindBin(ptmin+0.00001);
     h3PtEtaPhiJet[itype]->GetXaxis()->SetRange(binMin, h3PtEtaPhiJet[itype]->GetNbinsX());

     hPtJet[itype] = (TH1F*) h3PtEtaPhiJet[itype]->Project3D("x");
     if(nEvents>0)
       hPtJet[itype]->Scale(1./nEvents,"width");
     SetHist((TH1F*) hPtJet[itype],"p_{T,corr}^{jet} (GeV)","1/N_{evt} dN/dp_{T,corr}^{jet} (GeV^{-1})");
     hPtJet[itype]->SetName(Form("hPtJet%s",suffix[itype].Data())); 

     //eta versus phi
     hEtaPhiJet[itype] = (TH2F*) h3PtEtaPhiJet[itype]->Project3D("yz");
     SetHist((TH1F*) hEtaPhiJet[itype],"#varphi^{jet} (rad)","#eta^{jet}");
     hEtaPhiJet[itype]->SetName(Form("hEtaPhiJet%s",suffix[itype].Data()));
  }

  //rho versus centrality
  hRhoCent = (TH2F*) histListRho->FindObject("fHistRhovsCent");
  if(hRhoCent) {
     SetHist((TH1F*) hRhoCent, "Centrality (%)", "#rho (GeV/c*rad^{-1})");
     hRhoCent->SetName("hRhoCent");
  }

  //______________
  //Draw histograms

  TCanvas *c[100];
  TH1F *frame[100];
  Int_t nCan = 0;
  TLegend *leg; 
 
  for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets

     //draw pt spectrum
     if(!hPtJet[itype]) continue; 
     c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: Pt %s jets",nCan,suffix[itype].Data()),600,450);
     SetCanvas((TCanvas*) c[nCan]);
     c[nCan]->SetLogy();

     frame[nCan] = gPad->DrawFrame(hPtJet[itype]->GetBinLowEdge(1),  
                                   1e-7,
                                   hPtJet[itype]->GetBinLowEdge(hPtJet[itype]->GetNbinsX()+1), 
                                   hPtJet[itype]->GetBinContent(hPtJet[itype]->GetMaximumBin())*2.);

     SetHist((TH1F*) frame[nCan],hPtJet[itype]->GetXaxis()->GetTitle(),hPtJet[itype]->GetYaxis()->GetTitle());
 

     hPtJet[itype]->DrawCopy("same");
 
     leg = new TLegend(0.35,0.5,0.88,0.88);
     SetLeg(leg);
     
     TString txt =  Form("%s jets AKT R=%.1f",suffix[itype].Data(),jetR);
     if(run>0) txt += Form(" run:%d",run); 
                              ;
     leg->AddEntry((TObject*) 0, txt.Data(),"");
     leg->AddEntry((TObject*) 0, Form("p_{T,trk}> %d MeV",TMath::Nint(minTrkPT*1000)),"");
     if(itype==1) leg->AddEntry((TObject*) 0, Form("E_{T}>%d MeV",TMath::Nint(minClusterET*1000)),""); 
     leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
     leg->Draw();

     c[nCan]->SaveAs(Form("%s_Pt_AKT%02d_pT%04d_ET%04d_Run%d_%s.%s",prefix.Data(),TMath::Nint(jetR*10),
			  TMath::Nint(minTrkPT*1000),TMath::Nint(minClusterET*1000), run, suffix[itype].Data(), suftype.Data()));

     nCan++;

     //_________________
     //draw eta versus phi

     if(!hEtaPhiJet[itype]) continue; 
     c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: eta-phi %s jets",nCan,suffix[itype].Data()),600,450);
     SetCanvas((TCanvas*) c[nCan]);
     c[nCan]->SetRightMargin(0.15);

     frame[nCan] = gPad->DrawFrame(hEtaPhiJet[itype]->GetXaxis()->GetBinLowEdge(1),  
                                   hEtaPhiJet[itype]->GetYaxis()->GetBinLowEdge(1),
                                   hEtaPhiJet[itype]->GetXaxis()->GetBinLowEdge(hEtaPhiJet[itype]->GetNbinsX()),
                                   hEtaPhiJet[itype]->GetYaxis()->GetBinLowEdge(hEtaPhiJet[itype]->GetNbinsY()));

     SetHist((TH1F*) frame[nCan],hEtaPhiJet[itype]->GetXaxis()->GetTitle(),hEtaPhiJet[itype]->GetYaxis()->GetTitle());
 

     hEtaPhiJet[itype]->DrawCopy("same,colz");
 
     leg = new TLegend(0.35,0.5,0.88,0.88);
     SetLeg(leg);

     leg->AddEntry((TObject*) 0, txt.Data(),"");
     leg->AddEntry((TObject*) 0, Form("p_{T,trk}> %d MeV",TMath::Nint(minTrkPT*1000)),"");
     if(itype==1) leg->AddEntry((TObject*) 0, Form("E_{T}>%d MeV",TMath::Nint(minClusterET*1000)),"");
     leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
     leg->AddEntry((TObject*) 0, Form("p_{T,corr}^{jet} > %.1f GeV", ptmin),"");
     leg->Draw();

     c[nCan]->SaveAs(Form("%s_EtaPhi_AKT%02d_pT%04d_ET%04d_Run%d_%s.%s",prefix.Data(),TMath::Nint(jetR*10),
			  TMath::Nint(minTrkPT*1000),TMath::Nint(minClusterET*1000), run, suffix[itype].Data(), suftype.Data()));

     nCan++;

  }//end of the loop over charged and full jets

  //draw rho versus centrality
  if(hRhoCent) {
     c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: Rho-Cent",nCan),600,450);
     SetCanvas((TCanvas*) c[nCan]);
     c[nCan]->SetRightMargin(0.15);
     c[nCan]->SetLogz();

     frame[nCan] = gPad->DrawFrame(hRhoCent->GetXaxis()->GetBinLowEdge(1),
                                   hRhoCent->GetYaxis()->GetBinLowEdge(1),
                                   hRhoCent->GetXaxis()->GetBinLowEdge(hRhoCent->GetNbinsX()),
                                   hRhoCent->GetYaxis()->GetBinLowEdge(hRhoCent->GetNbinsY()));

     SetHist((TH1F*) frame[nCan],hRhoCent->GetXaxis()->GetTitle(),hRhoCent->GetYaxis()->GetTitle());

     hRhoCent->DrawCopy("colz");

     leg = new TLegend(0.35,0.5,0.88,0.88);
     SetLeg(leg);

     TString txt2 = Form("Jets KT R=%.1f", jetR);
     if(run>0) txt2 += Form(" run:%d",run);

     leg->AddEntry((TObject*) 0, txt2.Data(),"");
     leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
     leg->Draw();

     c[nCan]->SaveAs(Form("%s_RhoCent_AKT%02d_pT%04d_ET%04d_Run%d.%s",prefix.Data(),TMath::Nint(jetR*10),
			  TMath::Nint(minTrkPT*1000),TMath::Nint(minClusterET*1000), run, suftype.Data()));

     nCan++;
  }

  //---------------------------------------------------------------------------------------------------
  //                       WRITE OUTPUT TO ROOT FILE
  //---------------------------------------------------------------------------------------------------


  /* Standalone output
    TFile *histOut = new TFile(Form("%s_AKT%02d_pT%04d_ET%04d_jetPtMin%.1f_Run%d.root",prefix.Data(),
    TMath::Nint(jetR*10), TMath::Nint(minTrkPT*1000), TMath::Nint(minClusterET*1000),ptmin,run),"RECREATE");
    
    for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets
     
    if(hPtJet[itype])     hPtJet[itype]->Write();
    if(hEtaPhiJet[itype]) hEtaPhiJet[itype]->Write();
    if(hRhoCent)          hRhoCent->Write();
    }
    
    histOut->Close();
  
  */

  // Common output - 
  // Added by sjena

  TFile *fout = TFile::Open(outfile,"UPDATE");
  fout->ls();
  
  TDirectoryFile *cdd = NULL;
  cdd = (TDirectoryFile*)fout->Get("JE");
  if(!cdd) {
    Printf("Warning: JE <dir> doesn't exist, creating a new one");
    cdd = (TDirectoryFile*)fout->mkdir("JE");
  }
  cdd->cd();
  cdd->ls();
  
  for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets
    
    if(hPtJet[itype])     hPtJet[itype]->Write(Form("%s%d_%s",prefix.Data(), itype, hPtJet[itype]->GetName()));
    if(hEtaPhiJet[itype]) hEtaPhiJet[itype]->Write(Form("%s%d_%s",prefix.Data(), itype, hEtaPhiJet[itype]->GetName()));
    if(hRhoCent)          hRhoCent->Write(Form("%s%d_%s",prefix.Data(), itype, hRhoCent->GetName()));
  }
  
  fout->Close();


}
//__________________________________________________________

void SetHist(TH1* h,TString titx, TString tity){

   h->GetXaxis()->SetTitle(titx.Data());
   h->GetYaxis()->SetTitle(tity.Data());
   h->GetXaxis()->SetTitleSize(0.06);
   h->GetYaxis()->SetTitleSize(0.06);
   h->GetYaxis()->SetTitleOffset(1.);
   h->GetXaxis()->SetTitleOffset(1.);
   h->SetLineWidth(3);

}
//_____________________________________________________________________

void SetCanvas(TCanvas* c){
   c->SetLeftMargin(0.15);
   c->SetBottomMargin(0.15);
   c->SetRightMargin(0.05);
   c->SetTopMargin(0.05);
   c->SetTickx();
   c->SetTicky();
}
//_____________________________________________________________________

void SetLeg(TLegend* le){
   le->SetFillColor(10);
   le->SetBorderSize(0);
   le->SetFillStyle(0);
   le->SetTextSize(0.05);
}
Commit	Line	Data
948b38b2	1	/***************************************************
	2	processJETriggerQA:
	3	To procees JE Triggere QA wagon's output
	4
	5	****************************************************/
	6
	7
	8
	9	const Float_t ptmin = 0. ; //lower cutoff of jet pt spectrum
	10
	11	void processJETriggerQA(TString strFileIn = "AnalysisResults.root",
	12	TString suftype="eps",
	13	Float_t jetR = 0.2,
	14	Float_t minTrkPT = 0.15,
	15	Float_t minClusterET = 0.3,
	16	Int_t run = 0,
	17	const char* outfile = "JETriggerQA_outfile.root"
	18	){
	19
	20	gStyle->SetOptStat(0);
	21	gStyle->SetOptTitle(0);
	22
	23	TString prefix = "fig_je_TriggerQA_";
	24
	25	TFile * f1 = TFile::Open(strFileIn.Data());
	26
	27	//Load histogram list
	28	TList *histList = 0x0;
	29	histList = (TList*)f1->Get(Form("TriggerQA_Jet_AKTFullR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_Jet_AKTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TC/TriggerQA_Jet_AKTFullR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_Jet_AKTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TC",
	30	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000),
	31	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000),
	32	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000),
	33	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000)));
	34
	35	//Load a second histogram list in order to incorporate results from the Rho Task as well
	36	TList *histListRho = 0x0;
	37	histListRho = (TList*)f1->Get(Form("Rho_Jet_KTChargedR%02d0_PicoTracks_pT%04d_CaloClustersCorr_ET%04d_pt_scheme_TPC_histos",
	38	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000)));
	39
	40	//---------------------------------------------------------------------------------------------------
	41	// jet histograms
	42	//---------------------------------------------------------------------------------------------------
	43
	44	const Int_t kJetType = 2;
	45	TString suffix [kJetType] = {"Charged","Full"};
	46	TH3F *h3PtEtaPhiJet[kJetType];
	47	TH1F *hPtJet[kJetType];
	48	TH2F *hEtaPhiJet[kJetType];
	49	TH2F *hRhoCent;
	50
	51	TH1F *hNEventSel = histList->FindObject("fhNEvents");
	52	Float_t nEvents = hNEventSel->GetBinContent(2);
	53
	54	for(Int_t itype = 0; itype < kJetType; itype++){
	55
	56	h3PtEtaPhiJet[itype] = (TH3F*) histList->FindObject(Form("fh3PtEtaPhiJet%s",suffix[itype].Data()));
	57	if(! h3PtEtaPhiJet[itype]) continue;
	58
	59	//jet pt spectra
	60	Int_t binMin = 1;
	61	if(ptmin>0.) binMin = h3PtEtaPhiJet[itype]->GetXaxis()->FindBin(ptmin+0.00001);
	62	h3PtEtaPhiJet[itype]->GetXaxis()->SetRange(binMin, h3PtEtaPhiJet[itype]->GetNbinsX());
	63
	64	hPtJet[itype] = (TH1F*) h3PtEtaPhiJet[itype]->Project3D("x");
805ebfa6	65	if(nEvents>0)
805ebfa6	66	hPtJet[itype]->Scale(1./nEvents,"width");
948b38b2	67	SetHist((TH1F*) hPtJet[itype],"p_{T,corr}^{jet} (GeV)","1/N_{evt} dN/dp_{T,corr}^{jet} (GeV^{-1})");
	68	hPtJet[itype]->SetName(Form("hPtJet%s",suffix[itype].Data()));
	69
	70	//eta versus phi
	71	hEtaPhiJet[itype] = (TH2F*) h3PtEtaPhiJet[itype]->Project3D("yz");
	72	SetHist((TH1F*) hEtaPhiJet[itype],"#varphi^{jet} (rad)","#eta^{jet}");
	73	hEtaPhiJet[itype]->SetName(Form("hEtaPhiJet%s",suffix[itype].Data()));
	74	}
	75
	76	//rho versus centrality
	77	hRhoCent = (TH2F*) histListRho->FindObject("fHistRhovsCent");
	78	if(hRhoCent) {
	79	SetHist((TH1F) hRhoCent, "Centrality (%)", "#rho (GeV/crad^{-1})");
	80	hRhoCent->SetName("hRhoCent");
	81	}
	82
	83	//______________
	84	//Draw histograms
	85
	86	TCanvas *c[100];
	87	TH1F *frame[100];
	88	Int_t nCan = 0;
	89	TLegend *leg;
	90
	91	for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets
	92
	93	//draw pt spectrum
	94	if(!hPtJet[itype]) continue;
	95	c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: Pt %s jets",nCan,suffix[itype].Data()),600,450);
	96	SetCanvas((TCanvas*) c[nCan]);
	97	c[nCan]->SetLogy();
	98
	99	frame[nCan] = gPad->DrawFrame(hPtJet[itype]->GetBinLowEdge(1),
	100	1e-7,
	101	hPtJet[itype]->GetBinLowEdge(hPtJet[itype]->GetNbinsX()+1),
	102	hPtJet[itype]->GetBinContent(hPtJet[itype]->GetMaximumBin())*2.);
	103
	104	SetHist((TH1F*) frame[nCan],hPtJet[itype]->GetXaxis()->GetTitle(),hPtJet[itype]->GetYaxis()->GetTitle());
	105
	106
	107	hPtJet[itype]->DrawCopy("same");
	108
	109	leg = new TLegend(0.35,0.5,0.88,0.88);
	110	SetLeg(leg);
	111
	112	TString txt = Form("%s jets AKT R=%.1f",suffix[itype].Data(),jetR);
	113	if(run>0) txt += Form(" run:%d",run);
	114	;
	115	leg->AddEntry((TObject*) 0, txt.Data(),"");
	116	leg->AddEntry((TObject) 0, Form("p_{T,trk}> %d MeV",TMath::Nint(minTrkPT1000)),"");
	117	if(itype==1) leg->AddEntry((TObject) 0, Form("E_{T}>%d MeV",TMath::Nint(minClusterET1000)),"");
	118	leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
	119	leg->Draw();
	120
	121	c[nCan]->SaveAs(Form("%s_Pt_AKT%02d_pT%04d_ET%04d_Run%d_%s.%s",prefix.Data(),TMath::Nint(jetR*10),
	122	TMath::Nint(minTrkPT1000),TMath::Nint(minClusterET1000), run, suffix[itype].Data(), suftype.Data()));
	123
	124	nCan++;
	125
	126	//_________________
	127	//draw eta versus phi
	128
	129	if(!hEtaPhiJet[itype]) continue;
	130	c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: eta-phi %s jets",nCan,suffix[itype].Data()),600,450);
131	SetCanvas((TCanvas*) c[nCan]);
132	c[nCan]->SetRightMargin(0.15);
133
134	frame[nCan] = gPad->DrawFrame(hEtaPhiJet[itype]->GetXaxis()->GetBinLowEdge(1),
135	hEtaPhiJet[itype]->GetYaxis()->GetBinLowEdge(1),
136	hEtaPhiJet[itype]->GetXaxis()->GetBinLowEdge(hEtaPhiJet[itype]->GetNbinsX()),
137	hEtaPhiJet[itype]->GetYaxis()->GetBinLowEdge(hEtaPhiJet[itype]->GetNbinsY()));
138
139	SetHist((TH1F*) frame[nCan],hEtaPhiJet[itype]->GetXaxis()->GetTitle(),hEtaPhiJet[itype]->GetYaxis()->GetTitle());
140
141
142	hEtaPhiJet[itype]->DrawCopy("same,colz");
143
144	leg = new TLegend(0.35,0.5,0.88,0.88);
145	SetLeg(leg);
146
147	leg->AddEntry((TObject*) 0, txt.Data(),"");
148	leg->AddEntry((TObject) 0, Form("p_{T,trk}> %d MeV",TMath::Nint(minTrkPT1000)),"");
149	if(itype==1) leg->AddEntry((TObject) 0, Form("E_{T}>%d MeV",TMath::Nint(minClusterET1000)),"");
150	leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
151	leg->AddEntry((TObject*) 0, Form("p_{T,corr}^{jet} > %.1f GeV", ptmin),"");
152	leg->Draw();
153
154	c[nCan]->SaveAs(Form("%s_EtaPhi_AKT%02d_pT%04d_ET%04d_Run%d_%s.%s",prefix.Data(),TMath::Nint(jetR*10),
155	TMath::Nint(minTrkPT1000),TMath::Nint(minClusterET1000), run, suffix[itype].Data(), suftype.Data()));
156
157	nCan++;
158
159	}//end of the loop over charged and full jets
160
161	//draw rho versus centrality
162	if(hRhoCent) {
163	c[nCan] = new TCanvas(Form("c%d",nCan),Form("c%d: Rho-Cent",nCan),600,450);
164	SetCanvas((TCanvas*) c[nCan]);
165	c[nCan]->SetRightMargin(0.15);
166	c[nCan]->SetLogz();
167
168	frame[nCan] = gPad->DrawFrame(hRhoCent->GetXaxis()->GetBinLowEdge(1),
169	hRhoCent->GetYaxis()->GetBinLowEdge(1),
170	hRhoCent->GetXaxis()->GetBinLowEdge(hRhoCent->GetNbinsX()),
171	hRhoCent->GetYaxis()->GetBinLowEdge(hRhoCent->GetNbinsY()));
172
173	SetHist((TH1F*) frame[nCan],hRhoCent->GetXaxis()->GetTitle(),hRhoCent->GetYaxis()->GetTitle());
174
175	hRhoCent->DrawCopy("colz");
176
177	leg = new TLegend(0.35,0.5,0.88,0.88);
178	SetLeg(leg);
179
180	TString txt2 = Form("Jets KT R=%.1f", jetR);
181	if(run>0) txt2 += Form(" run:%d",run);
182
183	leg->AddEntry((TObject*) 0, txt2.Data(),"");
184	leg->AddEntry((TObject*) 0, Form("#it{N}_{events} = %.0f",nEvents),"");
185	leg->Draw();
186
187	c[nCan]->SaveAs(Form("%s_RhoCent_AKT%02d_pT%04d_ET%04d_Run%d.%s",prefix.Data(),TMath::Nint(jetR*10),
188	TMath::Nint(minTrkPT1000),TMath::Nint(minClusterET1000), run, suftype.Data()));
189
190	nCan++;
191	}
192
193	//---------------------------------------------------------------------------------------------------
194	// WRITE OUTPUT TO ROOT FILE
195	//---------------------------------------------------------------------------------------------------
196
197
198	/* Standalone output
199	TFile *histOut = new TFile(Form("%s_AKT%02d_pT%04d_ET%04d_jetPtMin%.1f_Run%d.root",prefix.Data(),
200	TMath::Nint(jetR10), TMath::Nint(minTrkPT1000), TMath::Nint(minClusterET*1000),ptmin,run),"RECREATE");
201
202	for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets
203
204	if(hPtJet[itype]) hPtJet[itype]->Write();
205	if(hEtaPhiJet[itype]) hEtaPhiJet[itype]->Write();
206	if(hRhoCent) hRhoCent->Write();
207	}
208
209	histOut->Close();
210
211	*/
212
213	// Common output -
214	// Added by sjena
215
216	TFile *fout = TFile::Open(outfile,"UPDATE");
217	fout->ls();
218
219	TDirectoryFile *cdd = NULL;
220	cdd = (TDirectoryFile*)fout->Get("JE");
221	if(!cdd) {
222	Printf("Warning: JE <dir> doesn't exist, creating a new one");
223	cdd = (TDirectoryFile*)fout->mkdir("JE");
224	}
225	cdd->cd();
226	cdd->ls();
227
228	for(Int_t itype = 0; itype < kJetType; itype++){ //loop over charged and full jets
229
230	if(hPtJet[itype]) hPtJet[itype]->Write(Form("%s%d_%s",prefix.Data(), itype, hPtJet[itype]->GetName()));
231	if(hEtaPhiJet[itype]) hEtaPhiJet[itype]->Write(Form("%s%d_%s",prefix.Data(), itype, hEtaPhiJet[itype]->GetName()));
232	if(hRhoCent) hRhoCent->Write(Form("%s%d_%s",prefix.Data(), itype, hRhoCent->GetName()));
233	}
234
235	fout->Close();
236
237
238
239	}
240	//__________________________________________________________
241
242	void SetHist(TH1* h,TString titx, TString tity){
243
244	h->GetXaxis()->SetTitle(titx.Data());
245	h->GetYaxis()->SetTitle(tity.Data());
246	h->GetXaxis()->SetTitleSize(0.06);
247	h->GetYaxis()->SetTitleSize(0.06);
248	h->GetYaxis()->SetTitleOffset(1.);
249	h->GetXaxis()->SetTitleOffset(1.);
250	h->SetLineWidth(3);
251
252	}
253	//_____________________________________________________________________
254
255	void SetCanvas(TCanvas* c){
256	c->SetLeftMargin(0.15);
257	c->SetBottomMargin(0.15);
258	c->SetRightMargin(0.05);
259	c->SetTopMargin(0.05);
260	c->SetTickx();
261	c->SetTicky();
262	}
263	//_____________________________________________________________________
264
265	void SetLeg(TLegend* le){
266	le->SetFillColor(10);
267	le->SetBorderSize(0);
268	le->SetFillStyle(0);
269	le->SetTextSize(0.05);
270	}