C++ TChain::Draw方法代码示例

double correlation_coef( string source1, string source2 ){

	TChain * c = new TChain( "rcpPicoDst" );

	c->Add( "/Users/danielbrandenburg/bnl/local/data/RcpAnalysis/data/RcpPicoDst/tuple_93477DF676758D2254762F71F52E3D44_14*.root" );

	TCanvas * c1 = new TCanvas( "c1", source1.c_str() );
	c->Draw( (source1 + " >> hSource1 ").c_str() );

	TCanvas * c2 = new TCanvas( "c2", source2.c_str() );
	c->Draw( (source2 + " >> hSource2 ").c_str() );

	TH1 * hSource1 = (TH1*)gDirectory->Get( "hSource1" );
	TH1 * hSource2 = (TH1*)gDirectory->Get( "hSource2" );

	double m1 = hSource1->GetMean();
	double m2 = hSource2->GetMean();

	double s1 = hSource1->GetRMS();
	double s2 = hSource2->GetRMS();

	TCanvas * c3 = new TCanvas( "c3", "correlation" );
	string expr = "( " + source1 + " - " + dts( m1 ) + " ) * ( " + source2 + " - " + dts( m2 ) + " ) / ( " + dts( s1 * s2 ) + " ) >> hCorr";
	c->Draw( expr.c_str() );

	TH1 * hCorr = (TH1*)gDirectory->Get( "hCorr" );


	return m1;

}

void paint(TString dir, TString a, TString b)
{
	TChain *t = new TChain(a);
	t->Add(dir + a + "_" + b +".root");

	std::cout << "painting the histograms for: " << dir + a + "_" + b +".root" << std::endl;

	TDirectory *direc = new TDirectory("dir", "dir");
	direc->cd();

	TH2D * histo = new TH2D("jets", "jets", 50, -2.5, 2.5, 40, 4.17438727, 6.95654544315);//pt starting from 15 and until 1000
	histo->SetDirectory(direc);

	t->Draw("log(jetPt+50):jetEta >> +jets", "", "Lego goff");

  TH2D * histo_lin = new TH2D("jets_lin", "jets_lin", 50, -2.5, 2.5, 40, 15, 1000);//pt starting from 15 and until 1000
	t->Draw("jetPt:jetEta >> +jets_lin", "", "Lego goff");

	std::cout << "saving the histograms: " << a + "_" + b +"_histo.root" << std::endl;
	TFile g(a + "_" + b +"_histo.root", "RECREATE");
	histo->SetDirectory(&g);
	histo_lin->SetDirectory(&g);
	delete direc;

	g.cd();
	histo->Write();
	histo_lin->Write();
	g.Close();

}

void resol() {
  TChain* trajTree = new TChain("trajTree");
  trajTree->Add("/data/jkarancs/CMSSW/TimingStudy/CMSSW_7_1_0_pre9/src/DPGAnalysis/PixelTimingStudy/Ntuple.root/trajTree");
  
  TH1D* dx_l1 = new TH1D("dx_l1", "PXB2 residuals #Deltax, p_{t} > 12, lever 1;PXB2 #Deltax [#mum];hits", 100, -150, 150); // h420_1
  TH1D* dx_l2 = new TH1D("dx_l2", "PXB1 residuals #Deltax, p_{t} > 12, lever 3;PXB1 #Deltax [#mum];hits", 100, -150, 150); // h520_3
  TH1D* dx_l3 = new TH1D("dx_l3", "PXB3 residuals #Deltax, p_{t} > 12, lever 3;PXB3 #Deltax [#mum];hits", 100, -150, 150); // i520_3
  
  TCanvas c("c");
  trajTree->Draw("res_dx>>dx_l1","layer==1&&pt>12&&lev>1.43&&lev<1.57");
  trajTree->Draw("res_dx>>dx_l2","layer==2&&pt>12&&lev>0&&lev<0.015","SAME");
  trajTree->Draw("res_dx>>dx_l3","layer==3&&pt>12&&lev>1.43&&lev<1.57","SAME");
  c.Close();
  
  double sigma_dx_l1 = fittp0(dx_l1);
  double sigma_dx_l2 = fittp0(dx_l2);
  double sigma_dx_l3 = fittp0(dx_l3);
  
  double res_dx_l1=sqrt((-sigma_dx_l1*sigma_dx_l1-10*sigma_dx_l2*sigma_dx_l2+5*sigma_dx_l3*sigma_dx_l3))/3;
  double res_dx_l2=sqrt((-sigma_dx_l1*sigma_dx_l1+14*sigma_dx_l2*sigma_dx_l2-sigma_dx_l3*sigma_dx_l3))/3;
  double res_dx_l3=sqrt((5*sigma_dx_l1*sigma_dx_l1-10*sigma_dx_l2*sigma_dx_l2-sigma_dx_l3*sigma_dx_l3))/3;
  
  std::cout<<"Resolution dx: Lay1 = "<<res_dx_l1<<" Lay2 = "<<res_dx_l2<<" Lay3 = "<<res_dx_l3<<std::endl;
  
  // Comparison with PixelTriplets:
  // PixelTimingStudy Resolution dx: Lay1 = 10.4437 Lay2 = 10.7449 Lay3 = 10.2808
  // PixelTriplets    Resolution dx: Lay1 = 10.5472 Lay2 = 10.7889 Lay3 = 10.4725
}

void drawThesis(){

  // TFile* f = TFile::Open("testing/comboStudies/generateAngularResolutionTreePlots_352_2016-03-02_11-57-40.root");
  // TFile* f = TFile::Open("testing/comboStudies/generateAngularResolutionTreePlots_352_2016-03-02_13-05-06.root");
  // TTree* t = (TTree*) f->Get("angResTree");

  TChain* t = new TChain("angResTree");
  t->Add("testing/comboStudies/generateAngularResolutionTreePlots_*");

  

  const int numBins = 64; //128;
  const double maxDeg = 2;
  TH1D* hTheta = new TH1D("hTheta", "#delta#theta Distribution of WAIS pulses; #delta#theta (Degrees); Events per bin", numBins, -maxDeg, maxDeg);

  TH1D* hPhi = new TH1D("hPhi", "#delta#phi Distribution of WAIS pulses; #delta#phi (Degrees); Events per bin", numBins, -maxDeg, maxDeg); 

  // t->Draw("deltaPhiDeg>>h2", "TMath::Abs(deltaPhiDeg) < 5", "colz");
  t->Draw("deltaThetaDeg3>>hTheta", "TMath::Abs(deltaPhiDeg3) < 5", "goff");
  t->Draw("deltaPhiDeg3>>hPhi", "TMath::Abs(deltaPhiDeg3) < 5", "goff");

  


  TH1D* hs[2] = {hTheta, hPhi};
  hPhi->SetLineColor(kRed);
  hTheta->SetLineColor(kBlue);  
  auto c1 = RootTools::drawHistsWithStatsBoxes(2, hs, "", "mre");

  TLegend* l1 = new TLegend(0.8, 0.8, 1, 1);
  l1->AddEntry(hPhi, "#delta#phi (Degrees)", "l");
  l1->AddEntry(hTheta, "#delta#phi (Degrees)", "l");
  hs[0]->SetTitle("Angular resolution of WAIS pulses");
}

void pixelHistos(const char *filenames = "/export/data1/dambach/hardware/ntpls/gradeC/*.root", const char *cuts =  "defects == 0")
{

  Init();      // do not forget to load the Utilities.C

  TH1D *gainH = new TH1D("gainH", "Gain", 1000, 0., 10.);
  TH1D *pedH = new TH1D("pedH", "Pedestal", 1000, -1000., 1000.);
  TH1D *noiseH = new TH1D("noiseH", "Noise", 1100, -100., 1000.);
  TH1D *tthrH = new TH1D("tthrH", "Trimmed Threshold", 100, 0., 100.);
  
  TChain *c = new TChain("mod");
  c->Add(filenames);
  
  canvas->Clear();
  canvas->Divide(2,2);

  canvas->cd(1);
  InitPad("log");
  gainH->GetXaxis()->SetTitle("ADC units/DAC units");
  gainH->GetYaxis()->SetTitle("# pixels");
  gainH->GetXaxis()->SetTitleSize(0.055);
  gainH->GetYaxis()->SetTitleSize(0.055);
  c->Draw("gain>>gainH", cuts);

  TBox *box = new TBox();

  canvas->cd(2);
  InitPad("log");
  pedH->GetXaxis()->SetTitle("DAC units");
  pedH->GetYaxis()->SetTitle("# pixels");
  pedH->GetXaxis()->SetTitleSize(0.055);
  pedH->GetYaxis()->SetTitleSize(0.055);
  c->Draw("ped>>pedH", cuts);

  canvas->cd(3);
  InitPad("log");
  noiseH->GetXaxis()->SetTitle("e^{-}");
  noiseH->GetYaxis()->SetTitle("# pixels");
  noiseH->GetYaxis()->SetRangeUser(0.1,400000);
  noiseH->GetXaxis()->SetTitleSize(0.055);
  noiseH->GetYaxis()->SetTitleSize(0.055);
  c->Draw("noise>>noiseH", cuts);
  noiseH->Fit("gaus");
//   box->SetFillStyle(3003);
//   box->SetFillColor(kBlue);
//   box->DrawBox(0.1,0.,500.,400000.);
//   box->SetFillColor(kRed);
//   box->DrawBox(500.,0.,1000.,400000.);

  canvas->cd(4);
  tthrH->GetXaxis()->SetTitle("DAC units");
  tthrH->GetYaxis()->SetTitle("# pixels");
  tthrH->GetXaxis()->SetTitleSize(0.055);
  tthrH->GetYaxis()->SetTitleSize(0.055);
  InitPad("log");
  c->Draw("tthr>>tthrH", cuts);

}

plot_BH2Raw()
{
  gROOT->Reset();
  TCanvas *c1 = new TCanvas("c1","c1",1600,900);
  TCanvas *c2 = new TCanvas("c2","c2",1600,900);

  int n = 7;
  TCanvas* c11[24];
  TCanvas* c22[24];

  c1->Divide(7,2);
  c2->Divide(7,2);

  char title1[255];
  char title2[255];

  TChain *chain = new TChain("tree");
  //chain->Add("test0.root");
  chain->Add("./root/lithium/old/run04520_Hodoscope.root");
  //chain->Add("./root/tdcshift/run04554_Hodoscope.root");

  //ADC
  for( int i=0; i<n; ++i){
    c11[i] = (TCanvas*) c1->GetPad(i+1);
    c1->cd(i+1);
    c11[i]->SetLogy();
    sprintf(title1,"bh2ua[%d]>>h10%d(500,0,2000)",i,i);
    //   sprintf(title2,"bh2ut[%d]>0",i);
    //     chain->Draw( title1, title2 );
    chain->Draw( title1 );
  }
  for( int i=0; i<n; ++i){
    c11[i+n] = (TCanvas*) c1->GetPad(i+1+n);
    c1->cd(i+1+n);
    c11[i+n]->SetLogy();
    sprintf(title1,"bh2da[%d]>>h20%d(500,0,1500)",i,i);
    //     sprintf(title2,"bh2dt[%d]>0",i);
    //     chain->Draw( title1, title2 );
    chain->Draw( title1 );
  }

  //TDC
  for( int i=0; i<n; ++i){
    c22[i] = (TCanvas*) c2->GetPad(i+1);
    c2->cd(i+1);
    c22[i]->SetLogy();
    sprintf(title1,"bh2ut[%d]>>h30%d(1000,500,1500)",i,i);
    chain->Draw( title1 );
  }
  for( int i=0; i<n; ++i){
    c22[i+n] = (TCanvas*) c2->GetPad(i+1+n);
    c2->cd(i+1+n);
    c22[i+n]->SetLogy();
    sprintf(title1,"bh2dt[%d]>>h40%d(1000,500,1500)",i,i);
    chain->Draw( title1 );
  }
}

void MakePlot(TString yvar = "mht30", int nbinsy = 100, double ylow=0., double yhigh=1000., TString xvar = "true_mu_ptW", int nbinsx = 100, double xlow=0., double xhigh=1000., TString axistitles="p_{T}^{W} [GeV];H_{T}^{miss} [GeV]", const TCut cuts="true_mu_pt>0&&num_true_mus==1&&num_true_els+num_true_had_taus==0", const TString treestring = "reduced_trees/13TeV/TTJets_MSDecaysCKM_central_Tune4C_13TeV-madgraph-tauola_Phys14DR-PU20bx25_PHYS14_25_V1-v1_MINIAODSIM_UCSB2406_v78/*.root", const TString plotNote="ttbar_1mu", bool logz = true, TString xvar_name="", TString yvar_name="") 
{
  set_plot_style();
  TH2D* h2 = new TH2D("h2",";"+axistitles,nbinsx,xlow,xhigh,nbinsy,ylow,yhigh);
  h2->SetStats(0);
  h2->GetZaxis()->SetLabelSize(0.02);
  h2->GetYaxis()->SetLabelSize(0.03);
  h2->GetYaxis()->SetTitleOffset(1.3);
  h2->GetXaxis()->SetLabelSize(0.03);
  TCut weighted_selection((cuts)*"(weightppb*10000.)");

  TChain * pypu = new TChain("reduced_tree");
  pypu->Add(treestring);// treestring is passed as an argument
  TCanvas * thecanvas= new TCanvas("thecanvas","the canvas",850,800);
  pypu->Draw(yvar+":"+xvar+">>h2",weighted_selection,"colz");

  h2->Scale(1/h2->Integral(1,nbinsx+1,1,nbinsy+1));

  thecanvas->SetLogz(logz);
  if (xvar_name.Length()>0) xvar=xvar_name;
  if (yvar_name.Length()>0) yvar=yvar_name;
  thecanvas->Print("macros/junews/lostleptons/conceptual_plots/"+yvar+"_vs_"+xvar+"_"+plotNote+"_scatter.pdf");
  
    
  delete h2;
  delete thecanvas;

}

void plotTowers_dj(
    TString inFile0Name="~/scratch01/ana/Hydj_BSC_HF_L1Emul/oh0928_v2/openhlt_hiReco_RAW2DIGI_RECO_*_djana.root",
    TString outdir = "out/Hydj/BSC_HF_L1Emul/devHLT381")
{
  TChain * djTree = new TChain("dijetAna_data_calojet_tower/djTree");
  djTree->Add(inFile0Name);
  aliases_dj(djTree);
  cout << " # entries: " << djTree->GetEntries() << endl;

  gSystem->mkdir(outdir.Data(),kTRUE);
  Int_t NTowerMax=70;


  pair<TString,TString> t = make_pair("name","cut");
  vector<pair<TString,TCut> > trigs;
  trigs.push_back(make_pair("All","1==1"));
  trigs.push_back(make_pair("HfAnyHit","hlt[2]"));
  trigs.push_back(make_pair("HfCoic1","hlt[3]"));
  trigs.push_back(make_pair("HfCoic2","hlt[4]"));

  for (unsigned int i=0; i<trigs.size(); ++i) {
    TCanvas * cHfNvsP = new TCanvas("HfNvsP_"+trigs[i].first,"HfNvsP_"+trigs[i].first,500,500);
    //cHfNvsP->SetLogx();
    //cHfNvsP->SetLogy();
    TH2D * hHfNvsP = new TH2D("hHfNvsP"+trigs[i].first,"hHfNvsP_",NTowerMax,0,NTowerMax,NTowerMax,0,NTowerMax);
    hHfNvsP->SetAxisRange(0.1,400,"X");
    hHfNvsP->SetAxisRange(0.1,400,"Y");
    djTree->Draw("NHfN:NHfP>>hHfNvsP"+trigs[i].first,trigs[i].second,"colz");
    hHfNvsP->SetTitle("HF Towers;# HF+ Towers (E>3GeV);# HF- Towers (E>3GeV)");
    cHfNvsP->Print(Form("%s/trHfNvsP_%s.gif",outdir.Data(),trigs[i].first.Data()));
  }
}

void paint(TString a, TString b)
{
	TChain *t = new TChain(a);
	t->Add(a + "_" + b +".root");

	std::cout << "painting the histograms for: " << a + "_" + b +".root" << std::endl;

	TDirectory *dir = new TDirectory("dir", "dir");
	dir->cd();

//	histo = new TH2D("jets", "jets", 50, -2.5, 2.5, 60, 10, 610);
//	histo = new TH2D("jets", "jets", 50, -2.5, 2.5, 40, 4.7004803657924166, 6.404803657924166);
//	TH2D * histo = new TH2D("jets", "jets", 25, -2.5, 2.5, 20, 4.0943445622221004, 6.1943445622221004); 
	TH2D * histo = new TH2D("jets", "jets", 50, -2.5, 2.5, 40, 4.0943445622221004, 6.1943445622221004);//original
//	TH2D * histo = new TH2D("jets", "jets", 50, -2.5, 2.5, 40, 4.0943445622221004, 7.8);
	histo->SetDirectory(dir);

	//the varexp part of the draw synatx means: draw log(jetPt+50) versus jetEta and append the existing ("+" -> avoid recreation) histogram called "jets"
	//selection is an expression with a combination of the Tree variables -> no selection applied in this case ""
	//option is the drawing option -> if option contains the string "goff", no graphics is generated.
	//fourth and fifth arguments are: Int_t nevents, Int_t firstevent
	t->Draw("log(jetPt+50):jetEta >> +jets", "", "Lego goff");
	//std::cout <<"jetPt " << log(jetPt+50) << " and jetEta " << jetEta << std::endl;


	std::cout << "saving the histograms: " << a + "_" + b +"_histo.root" << std::endl;
	TFile g(a + "_" + b +"_histo.root", "RECREATE");
	histo->SetDirectory(&g);
	delete dir;

	g.cd();
	histo->Write();
	g.Close();
}

void make_psi2s_mupt(float ptmin=0.0, float ptmax=30.0, float ymin=0.0, float ymax=2.4, bool absRapidity=true, bool saveFile=false)
{
  TCanvas *c1 = new TCanvas("c1","c1");
  c1->SetLogy();
  TCut defaultCut = "Reco_QQ_sign==0&&Reco_QQ_4mom.M()>3.35&&Reco_QQ_4mom.M()<4.0&&Reco_QQ_ctauTrue>-10";
  TCut ptCut = Form("Reco_QQ_4mom.Pt()>%4.1f&&Reco_QQ_4mom.Pt()<%4.1f",ptmin,ptmax);
  TCut rapCut;
  if (absRapidity)
    rapCut = Form("abs(Reco_QQ_4mom.Rapidity())>%3.1f&&abs(Reco_QQ_4mom.Rapidity())<%3.1f",ymin,ymax);
  else
    rapCut = Form("Reco_QQ_4mom.Rapidity()>%3.1f&&Reco_QQ_4mom.Rapidity()<%3.1f",ymin,ymax);

  unsigned int trigBit=2; // DoubleMu0_HighQ
  TCut trigCut = Form("(HLTriggers&%u)==%u&&(Reco_QQ_trig&%u)==%u",trigBit,trigBit,trigBit,trigBit);

  TString fname;
  fname = Form("20140324/MC_psi2s_pp_mupt_Rap_%3.1f-%3.1f_Pt_%3.1f-%3.1f.root",ymin,ymax,ptmin,ptmax);

  std::cout << fname << std::endl;


  std::cout << "default: " << defaultCut.GetTitle() << std::endl;
  std::cout << "pt cut: " << ptCut.GetTitle() << std::endl;
  std::cout << "rapidity cut: " << rapCut.GetTitle() << std::endl;
  std::cout << "trigger bit: " << trigCut.GetTitle() << std::endl;

  TH1F *hMuPlPtRec = new TH1F("hMuPlPtRec","hMuPlPtRec;p_{T} (#mu^{+}) (GeV/c);Events",200,0,20);
  TH1F *hMuMiPtRec = new TH1F("hMuMiPtRec","hMuMiPtRec;p_{T} (#mu^{-}) (GeV/c);Events",200,0,20);
  TH1F *hMuPtRec = new TH1F("hMuPtRec","hMuPtRec;p_{T} (#mu^{#pm}) (GeV/c);Events",200,0,20);
  hMuPlPtRec->Sumw2();
  hMuMiPtRec->Sumw2();
  hMuPtRec->Sumw2();

  hMuPlPtRec->SetMarkerColor(kRed);
  hMuMiPtRec->SetMarkerColor(kBlue);

  TChain *myTree = new TChain("myTree");
  myTree->Add("../root_files/PRpsi2SMC_Histos_2013pp_GlbGlb_STARTHI53_V28-v1_muLessPV.root");

  myTree->Draw("Reco_QQ_mupl_4mom.Pt()>>hMuPlPtRec",defaultCut&&ptCut&&rapCut&&trigCut,"e");
  myTree->Draw("Reco_QQ_mumi_4mom.Pt()>>hMuMiPtRec",defaultCut&&ptCut&&rapCut&&trigCut,"e");

  hMuPtRec->Add(hMuPlPtRec,hMuMiPtRec);

  hMuPtRec->Draw();
  hMuPlPtRec->Draw("same");
  hMuMiPtRec->Draw("same");

  TFile *outf = NULL;
  if (saveFile) {
    outf = new TFile(fname,"RECREATE");
    hMuPtRec->Write();
    outf->Close();
  }

  return;
}

void anaFragShape(TString infile="dj_HCPR-GoodTrkAndPixel_CleanEvt1130.root")
{
  TChain * djcalo = new TChain("djcalo/djTree");
  djcalo->Add(infile);
  aliases_dijet(djcalo);
  TChain * djcalopix = new TChain("djcalo_pxl/djTree");
  djcalopix->Add(infile);
  aliases_dijet(djcalopix);
  TString evtSel("(cent<10 && nljet>120 && abs(nljeta)<2 && aljet>50 && abs(aljeta)<2 && jdphi>2.5 && Aj>0.24)");
  TH1::SetDefaultSumw2();

  TCanvas * c0 = new TCanvas("c0","c0",500,500);
  djcalo->Draw("Aj>>hAj(20,0,1)",evtSel);
  djcalopix->Draw("Aj",evtSel,"Esame");

  TCanvas * c2 = new TCanvas("c2","c2",500,500);
  TH1D * hJDPhi = new TH1D("hJDPhi","hJDPhi",50,0,3.1416);
  TH1D * hJDPhi2 = new TH1D("hJDPhi2","hJDPhi",50,0,3.1416);
  Float_t numDJ = djcalo->Draw("jdphi>>hJDPhi",evtSel);
  Float_t numDJ2 = djcalopix->Draw("jdphi>>hJDPhi2",evtSel);
  cout << "num djs: " << numDJ << "  djs(pix)" << numDJ2 << endl;

  TH1D * hPNDRTrk = new TH1D("hPNDRTrk",";R(trk,jet);1/(N_{DJ} 2#piR) dp_{T}^{Trk}/dR",10,0,TMath::PiOver2());
  TH1D * hPNDRPix = new TH1D("hPNDRPix",";R(trk,jet);1/(N_{DJ} 2#piR) dp_{T}^{Trk}/dR",10,0,TMath::PiOver2());
  TH1D * hPADRTrk = new TH1D("hPADRTrk",";R(trk,jet);1/(N_{DJ} 2#piR) dp_{T}^{Trk}/dR",10,0,TMath::PiOver2());
  TH1D * hPADRPix = new TH1D("hPADRPix",";R(trk,jet);1/(N_{DJ} 2#piR) dp_{T}^{Trk}/dR",10,0,TMath::PiOver2());
  djcalo->Draw("pndr>>hPNDRTrk",Form("(%s&&ppt>=1.2)*(ppt/(TMath::TwoPi()*pndr))",evtSel.Data()));
  djcalopix->Draw("pndr>>hPNDRPix",Form("(%s&&ppt>=0.3&&ppt<1.2)*(ppt/(TMath::TwoPi()*pndr))",evtSel.Data()));
  djcalo->Draw("padr>>hPADRTrk",Form("(%s&&ppt>=1.2)*(ppt/(TMath::TwoPi()*padr))",evtSel.Data()));
  djcalopix->Draw("padr>>hPADRPix",Form("(%s&&ppt>=0.3&&ppt<1.2)*(ppt/(TMath::TwoPi()*padr))",evtSel.Data()));
  hPNDRTrk->Scale(1./(numDJ*hPNDRTrk->GetBinWidth(1)));
  hPNDRPix->Scale(1./(numDJ*hPNDRPix->GetBinWidth(1)));
  hPADRTrk->Scale(1./(numDJ*hPADRTrk->GetBinWidth(1)));
  hPADRPix->Scale(1./(numDJ*hPADRPix->GetBinWidth(1)));
  hPNDRPix->SetMinimum(50);
  hPNDRPix->SetMaximum(100);
  hPNDRTrk->SetMinimum(10);
  hPNDRTrk->SetMaximum(2000);
  hPNDRPix->SetMarkerStyle(kFullCircle);
  hPNDRTrk->SetMarkerStyle(kFullSquare);
  hPNDRPix->SetMarkerColor(kRed);
  hPNDRTrk->SetMarkerColor(kRed);
  hPADRPix->SetMarkerStyle(kOpenCircle);
  hPADRTrk->SetMarkerStyle(kOpenSquare);
  hPADRPix->SetMarkerColor(kBlue);
  hPADRTrk->SetMarkerColor(kBlue);
  TCanvas * c3 = new TCanvas("c3","c3",500,500);
  hPNDRPix->Draw("E");
  hPADRPix->Draw("Esame");
  TCanvas * c3_2 = new TCanvas("c3_2","c3_2",500,500);
  c3_2->SetLogy();
  hPNDRTrk->Draw("E");
  hPADRTrk->Draw("Esame");
}

void quickMakeHisto() {

  
  new TCanvas();
  TChain *nocontainerTree = new TChain("pcaTree");
  nocontainerTree->Add("/home/rjn/creamtea/data/nocontainer/pca_nocontainer_0_1000files.root");
  nocontainerTree->Add("/home/rjn/creamtea/data/nocontainer/pca_nocontainer_1000_2000files.root");
  nocontainerTree->Add("/home/rjn/creamtea/data/nocontainer/pca_nocontainer_3000_4000files.root");

  TFile *fpNocontainer = new TFile("nocontainer.root","RECREATE");
  TH3F *histNocontainer = new TH3F("histNoContainer","histNoContainer",100,-6000,+6000,100,-6000,+6000,100,-6000,+6000);
  
  nocontainerTree->Draw("xPosTrue:yPosTrue:zPosTrue>>histNoContainer","abs(zPosTrue)<3000 && abs(xPosTrue)<6000 && abs(yPosTrue)<6000 && thetaTrue>0.00","");
  histNocontainer->Write();


  new TCanvas();
  TChain *containerTree = new TChain("pcaTree");
  containerTree->Add("/home/rjn/creamtea/data/container13m/pca_container_1_1000files.root");
  containerTree->Add("/home/rjn/creamtea/data/container13m/pca_container_1001_2000files.root");
  containerTree->Add("/home/rjn/creamtea/data/container13m/pca_container_2001_3000files.root");

  TFile *fpcontainer = new TFile("container.root","RECREATE");
  TH3F *histcontainer = new TH3F("histContainer","histContainer",100,-6000,+6000,100,-6000,+6000,100,-6000,+6000);
  
  containerTree->Draw("xPosTrue:yPosTrue:zPosTrue>>histContainer","abs(zPosTrue)<3000 && abs(xPosTrue)<6000 && abs(yPosTrue)<6000 && thetaTrue>0.00","");
  histcontainer->Write();


  new TCanvas();
  TChain *targetTree = new TChain("pcaTree");
  targetTree->Add("/home/rjn/creamtea/data/target/pca_container_1_1000files.root");
  targetTree->Add("/home/rjn/creamtea/data/target/pca_container_1001_2000files.root");
  targetTree->Add("/home/rjn/creamtea/data/target/pca_container_2001_3000files.root");

  TFile *fptarget = new TFile("target.root","RECREATE");
  TH3F *histtarget = new TH3F("histTarget","histTarget",100,-6000,+6000,100,-6000,+6000,100,-6000,+6000);
  
  targetTree->Draw("xPosTrue:yPosTrue:zPosTrue>>histTarget","abs(zPosTrue)<3000 && abs(xPosTrue)<6000 && abs(yPosTrue)<6000 && thetaTrue>0.00","");
  histtarget->Write();



}

void Sums(const char* savename,unsigned batchmode=0){
  
  if (batchmode ==1){
    gROOT->SetBatch(kTRUE); //set batch mode to suppress canvas output
     }
  //TChain* V;
  // gROOT->ProcessLine(".x dual-seg-chainload.C");
  TChain *P = new TChain("Tab");
 
   for (int i=1;i<11;i++){
       stringstream ss;
       if (i <10){
  	ss<<"094rb_14-seg-0"<<i<<".root";
       }
       else if (i >= 10){
  	ss<<"094rb_14-seg-"<<i<<".root";
       }
       P->Add(ss.str().c_str());
   }
  

   cout<<"chain loaded"<<endl;
      
  TFile *f = new TFile(savename,"RECREATE");

  
  TH2D* NvG = new TH2D("NvG"," NaI (x) vs Ge (y) Paddbacks",8000.,0.,8000.,8000.,0.,8000.);
  TH2D* LvG = new TH2D("LvG"," LaBr3 (x) vs Ge (y) Paddbacks",8000.,0.,8000.,8000.,0.,8000.);
  TH2D* NvL = new TH2D("NvL"," NaI (x) vs LaBr3 (y) Paddbacks",8000.,0.,8000.,8000.,0.,8000.);

  cout<<"Hist + files made"<<endl;
    
  NvG->SetMinimum(1);
  LvG->SetMinimum(1);
  NvL->SetMinimum(1);
  cout<<"minimum set"<<endl;
  P->Draw("pab.GabE:pab.NabE>>NvG");
  cout<<"NvG done"<<endl; 
  P->Draw("pab.GabE:pab.LabE>>LvG");
  cout<<"LvG done"<<endl; 
  P->Draw("pab.LabE:pab.NabE>>NvL");
  
  f->Write();
}

void DrawdEdxResolExample(){
  /// Example analysis to make an space point resolution study

  TChain * chain  = AliXRDPROOFtoolkit::MakeChain("trackerSimul.list", "simulTrack",0,100); 
  chain->SetCacheSize(10000000000);

  //
  // 1.) Qmax/Qtot as function of the input ionization density
  //
  chain->Draw("tr.CookdEdxDmax(0,0.6,1,0,1,0)/tr.CookdEdxDtot(0,0.6,1,0,1,0):tr.fMNprim>>hisQtotMax(10,10,50)","","prof",10000);
  //
  // 2.) Non linearity due to the truncation Qtot_{60%}/Qtot 100% 
  //
  chain->Draw("tr.CookdEdxDtot(0,0.6,1,0,1,0)/tr.CookdEdxDtot(0,0.99,1,0,1,0):tr.fMNprim>>hisQtot60100(10,10,50)","","prof",10000);
  //
  // 3.) 
  //
  chain->Draw("tr.CookdEdxDtot(0,0.6,1,0,1,0)/tr.fMNprim:tr.fMNprim>>profQtot60(10,10,50)","","prof",10000);
  chain->Draw("tr.CookdEdxDtot(0,0.99,1,0,1,0)/tr.fMNprim:tr.fMNprim>>profQtot100(10,10,50)","","profsame",10000);

}

TChain * compJOO(
    TString infile="hltana_Exp_151058.root",
    TString header="HLT_HIJet90U"
    )
{
  TChain * dj = new TChain("hltanalysis/HltTree");

  dj->Add(infile);
  dj->AddFriend("onl = hltanalysisOnl/HltTree",infile);
  cout << "Trig Total: " << dj->GetEntries() << endl;

  TCut evtSel("HLT_HIJet90U && hiBin*2.5<100");

  TH1D * hJEtIc5pu = plotJEt(dj,evtSel,"recoJetCalPt[0]","hJEtIc5pu",0);
  TH1D * hJEtIc5puOnl = plotJEt(dj,evtSel,"onl.recoJetCalPt[0]","hJEtIc5puOnl",0);

  TCut matEvtSel(evtSel);
  bool doMat = true;
  if (doMat) matEvtSel = evtSel && "abs(onl.recoJetCalEta[0]-recoJetCalEta[0])<0.2";
  TH2D * hJEtOnlIc5pu = plotJEtCorr(dj,matEvtSel,"onl.recoJetCalPt[0]/recoJetCalPt[0]:recoJetCalPt[0]","hJEtOnlIc5pu",";leading E_{T}^{Ic5pu} [GeV];E_{T}^{Ic5pu Online}/E_{T}^{Ic5pu}");

  TCanvas * cJEt = new TCanvas("cJEt","cJEt",500,500);
  CPlot cpJEt("JEt","JEt","leading E_{T}^{Jet}","#");
  cpJEt.SetLogy();
  cpJEt.AddHist1D(hJEtIc5puOnl,"Online Reco","hist",kGreen+2,kOpenSquare);
  cpJEt.AddHist1D(hJEtIc5pu,"Offline Reco","E",kBlack,kFullCircle);
  cpJEt.SetLegend(0.62,0.78,0.94,0.9);
  cpJEt.SetLegendHeader(header);
  cpJEt.Draw(cJEt,false);

  TCanvas * cDEta = new TCanvas("cDEta","cDEta",500,500);
  cDEta->SetLogy();
  dj->Draw("onl.recoJetCalEta[0]-recoJetCalEta[0]>>hDEta",evtSel,"E");
  TCanvas * cDPhi = new TCanvas("cDPhi","cDPhi",500,500);
  cDPhi->SetLogy();
  dj->Draw("onl.recoJetCalPhi[0]-recoJetCalPhi[0]>>hDPhi",evtSel,"E");

  return dj;
}