C++ TH1D::SetTitle方法代码示例

	void plotUstatistic( IPDF * pdf, IDataSet * data, PhaseSpaceBoundary * phase, string plot )
	{
		unsigned int nD = (unsigned) data->GetDataNumber();
		int bins = 30;
		double level = (double)nD/(double)bins;
		TH1D * distances = new TH1D( "distances", "U", bins, 0., 1.);
		distances->Sumw2();
		calculateUstatisticNum(pdf, data, phase, distances);

		char buff[10];
		sprintf( buff, "%f", level );
		TF1 * line = new TF1("line", buff, 0, 1);
		line->SetLineColor(kBlue);

		TCanvas* ca = new TCanvas("ca","canvas", 3000, 3000);
		gStyle->SetOptStat(0);
		distances->SetTitle("");
		distances->GetXaxis()->SetNdivisions(508);
		distances->GetXaxis()->SetTitle("U");
		distances->Draw();
		line->Draw("same");
		ca->Update();
		ca->SaveAs(plot.c_str());
		//delete distances;
		//delete ca;
	}

void PrintHist(const char *filename = "test.root", const char *histname = "histname") {

  Init();

  TFile *f = new TFile(filename);
  TH1D *histo;

  histo = (TH1D*)(f->Get(histname));
  histo->GetXaxis()->SetTitleSize(0.055);
  histo->GetYaxis()->SetTitleSize(0.055);
  histo->GetXaxis()->SetLabelSize(0.04);
  histo->GetYaxis()->SetLabelSize(0.05);
  histo->GetXaxis()->SetTitleOffset(1.15);
  histo->GetYaxis()->SetTitleOffset(1.1);
  histo->SetTitle("");
  histo->SetLineStyle(1);
  histo->SetLineWidth(2);
//   histo->GetXaxis()->SetTitle("Vcal [low range DAC units]");
//   histo->GetYaxis()->SetTitle("PH [ADC units]");
     histo->GetXaxis()->SetTitle("par1");
     histo->GetYaxis()->SetTitle("# pixels");
//    histo->GetXaxis()->SetTitle("position resolution [#mum]");
//    histo->GetYaxis()->SetTitle("# pixels");
//     histo->GetXaxis()->SetTitle("Vana [DAC units]");
//     histo->GetYaxis()->SetTitle("par1");
//      histo->GetXaxis()->SetTitle("nChip");
//      histo->GetYaxis()->SetTitle("par1");


     histo->Draw();
     //histo->Draw("colz");

}

void mc2ibd_NE(TChain *tMC, TFile *fBgnd, TCanvas *cv)
{
	char str[1024];
	TLatex *txt = new TLatex();

	TH1D *hExp = (TH1D *) fBgnd->Get("hNEA-diff");
	if (!hExp) {
		printf("Histogram hNEA-diff not found in %s\n", fBgnd->GetName());
		return;
	}
	hExp->SetTitle("Delayed event energy;MeV;Events/200 keV");
	hExp->SetLineColor(kBlack);
	hExp->SetLineWidth(3);

	gROOT->cd();
	TH1D *hMC = new TH1D("hNEMC", "Delayed event energy (MC);MeV;Events/200 keV", 45, 3, 12);
	hMC->SetLineColor(kBlue);
	tMC->Project(hMC->GetName(), "NeutronEnergy", cX && cY && cZ && cR && c20 && cGamma && cGammaMax && cPe);
	hMC->Sumw2();
	hMC->Scale(hExp->Integral(15, 45) / hMC->Integral(15,45));
	
	cv->Clear();
	hMC->Draw("hist");
	hExp->DrawCopy("same");
	TLegend *lg = new TLegend(0.65, 0.8, 0.89, 0.89);
	lg->AddEntry(hExp, "IBD", "LE");
	lg->AddEntry(hMC, "MC", "L");
	lg->Draw();
}

void Plot1DimFoams(TList& foam_list, TMVA::ECellValue cell_value,
                   const TString& cell_value_description,
                   TMVA::PDEFoamKernelBase* kernel)
{
   // visualize a 1 dimensional PDEFoam via a histogram
   TCanvas* canvas = NULL;
   TH1D* projection = NULL;

   // loop over all foams and draw the histogram
   TListIter it(&foam_list);
   TPair* fm_pair = NULL;    // the (foam, caption) pair
   while (fm_pair = (TPair*) it()) {
      TMVA::PDEFoam* foam = (TMVA::PDEFoam*) fm_pair->Key();
      if (!foam) continue;
      TString foam_caption(((TObjString*) fm_pair->Value())->String());
      TString variable_name(foam->GetVariableName(0)->String());

      canvas = new TCanvas(Form("canvas_%u",foam),
                           "1-dimensional PDEFoam", 400, 400);

      projection = foam->Draw1Dim(cell_value, 100, kernel);
      projection->SetTitle(cell_value_description + " of " + foam_caption
                           + ";" + variable_name);
      projection->Draw();
      projection->SetDirectory(0);

      canvas->Update();
   }
}

void Drawtrackeff(){
    gStyle->SetOptStat(kFALSE);
    TString type="eta";
    TCanvas *c1 = new TCanvas("c1","c1",600,600);
    TH1D* heff = Draw(type,1,20);
    TH1D* hFrame = new TH1D("","",2030,-3,200);
    hFrame->SetTitle("");
    if(type=="pt"){
        c1->SetLogx();
   //   c1->SetLogy();
        hFrame->GetXaxis()->SetTitle("p_{T}");
        hFrame->GetXaxis()->SetRangeUser(0,200);
    }
    else if(type=="eta"){
        hFrame->GetXaxis()->SetTitle("#eta");
        hFrame->GetXaxis()->SetRangeUser(-3,3);
    }
    hFrame->GetYaxis()->SetRangeUser(0,1);
    hFrame->GetYaxis()->SetTitle("tracking efficiency");
    hFrame->GetYaxis()->SetTitleOffset(1.1);
    TLegend *leg = new TLegend(0.3,0.75,0.85,0.88);
    leg->SetBorderSize(0);
    leg->SetFillColor(0);
    leg->SetTextSize(0.025);
    TLine *l = new TLine(hFrame->GetXaxis()->GetXmin(),1,hFrame->GetXaxis()->GetXmax(),1);
    l->SetLineStyle(2);
    l->SetLineWidth(1.2);
    hFrame->Draw();
    heff->Draw("Psame");
    l->Draw("same");
    //leg->Draw("same");
    
    c1->Print(Form("trackeff_%s.png",type.Data()));
}

void plotter::draw_purity(TH1D* numerator_, TH1D* denominator_, TString file_name){
  TH1D* numerator = (TH1D*) numerator_->Clone("numerator");
  TH1D* denominator = (TH1D*) denominator_->Clone("denominator");

  TH1D* purity = numerator; // just to set correct binning
  purity->Divide(numerator, denominator, 1., 1., "B");

  purity->SetTitle(" ");
  purity->GetXaxis()->SetTitle("m_{gen}");
  purity->GetYaxis()->SetTitle("purity");
  purity->GetYaxis()->SetRangeUser(0,1);

  purity->GetXaxis()->SetTitleSize(0.05);
  purity->GetYaxis()->SetTitleSize(0.05);
  purity->GetXaxis()->SetTitleOffset(0.9);
  purity->GetYaxis()->SetTitleOffset(0.8);
  purity->GetXaxis()->SetNdivisions(505);
  purity->GetYaxis()->SetNdivisions(505);

  purity->SetMarkerStyle(20);
  purity->SetMarkerSize(0.8);
  purity->SetLineColor(1);

  TCanvas *c= new TCanvas("Purity","",600,600);
  gPad->SetLeftMargin(0.15);
  TGaxis::SetMaxDigits(3);
  purity->Draw("E1");
  gPad->RedrawAxis();
  c->SaveAs(directory + file_name + ".pdf");
  delete c;
}

void mc2ibd_R2(TChain *tMC, TFile *fBgnd, TCanvas *cv)
{
	char str[1024];
	TLatex *txt = new TLatex();

	TH1D *hExp = (TH1D *) fBgnd->Get("hR2A-diff");
	if (!hExp) {
		printf("Histogram hR2A-diff not found in %s\n", fBgnd->GetName());
		return;
	}
	hExp->SetTitle("Distance between positron and neutron, 3D case;cm;Events/4cm");
	hExp->SetLineColor(kBlack);
	hExp->SetLineWidth(3);

	gROOT->cd();
	TH1D *hMC = new TH1D("hR2MC", "Distance between positron and neutron, 3D case (MC);cm;Events/4cm", 40, 0, 160);
	hMC->SetLineColor(kBlue);
	tMC->Project(hMC->GetName(), "Distance", cX && cY && cZ && cRXY && c20 && cGamma && cGammaMax && cPe && cN);
	hMC->Sumw2();
	hMC->Scale(hExp->Integral() / hMC->Integral());
	
	cv->Clear();
	hExp->DrawCopy();
	hMC->Draw("hist,same");
	TLegend *lg = new TLegend(0.65, 0.8, 0.89, 0.89);
	lg->AddEntry(hExp, "IBD", "LE");
	lg->AddEntry(hMC, "MC", "L");
	lg->Draw();
}

TH1D* Analysis(TFile *f)
{
    
    const int firstRoc = 1;
    const int lastRoc = 2;
    const int firstCol = 0;
    const int lastCol = 52;
    const int firstRow = 0;
    const int lastRow = 80;
        
    const int numRoc = lastRoc - firstRoc + 1;
    const int numCol = lastCol - firstCol + 1;
    const int numRow = lastRow - firstRow + 1;
        
    f->cd();


    TH1D *hist[numRoc][numCol][numRow];
    TH1D *qualityHist = new TH1D("quality","quality",1000,-0.5,0.5);
    
    qualityHist->SetLineColor(color);
    qualityHist->GetYaxis()->SetTitle("# pixels");


    qualityHist->GetXaxis()->SetTitle("a/b");
    qualityHist->SetTitle("a over b");
    qualityHist->GetXaxis()->SetRangeUser(-0.05,0.05);
    qualityHist->GetYaxis()->SetRangeUser(0,2500);
    
    int numRows = 0; int numCols = 0; int numRocs = 0;
    double quality;
    for (int roc = firstRoc; roc < lastRoc; roc++)
    {
	for (int col = firstCol; col < lastCol; col++)
	{ 
	    for (int row = firstRow; row < lastRow; row++)
	    {
		if (debug) cout << "roc = " << roc << "   col = " << col << "   row = " << row << endl;
		if (debug) cout << "numRocs = " << numRocs << "   numCols = " << numCols << "   numRows = " << numRows << endl;

		hist[numRocs][numCols][numRows] = (TH1D*)(f->Get(Form("PhVcal_c%ir%i_C%i;1",col,row,roc)))->Clone(Form("PhVcal_c%ir%i_C%i;1",col,row,roc)); 
 		quality = QualityLowRange(hist[numRocs][numCols][numRows]);

		if (debug) cout << "quality = " << quality << endl;

		if (-99 != quality) qualityHist->Fill(quality);		
		numRows++;
	    }
	    numCols++;
	    numRows = 0;
	}
	numRocs++;
	numCols = 0;
    }
    numRocs = 0;

    return qualityHist;    

}

void  makeGaussianSignals(SigData_t& m_sigdata)
{
  //std::vector<TH1D *> vgsh(NUMCHAN);
  std::vector<TH1D *> vcdfh(NUMCHAN);

  if( m_sigdata.find("gs") == m_sigdata.end() ) {
    cerr << "Gaussian signal data not found, not making CDF signal!" << endl;
    return;
  }

  for (int ichan=0; ichan<NUMCHAN; ichan++) {
    TH1D *cdfh;

    TString channame(channames[ichan]);
    TString name;

    TH1D * gsh = m_sigdata["gs"].at(ichan);

    assert(gsh) ;

#if 0
    name = "Signalgs_"+channame;
    gsh = (TH1D *)tch->Clone(name.Data());

    assert(gsh);

    gsh->SetTitle("Gaussian signal");
    
    gsh->Reset();

    TF1 *g = (TF1 *)gROOT->GetFunction("gaus");
    g->SetParameters(1,gaussian_mean_mass_gev,gaussian_mass_sigma_gev);
    gsh->FillRandom("gaus",100000);

    // norm to 1pb signal with 1/fb integrated luminosity
    double norm = 1000 * gseffxacc[ichan]/gsh->Integral(0,gsh->GetNbinsX()+1,"width");

    //gsh->Scale(norm/eff_fudgefactor); // kludge: pre-undo the fudge in the next module
    gsh->Scale(norm);

    vgsh[ichan] = gsh;
#endif

    // New CDF bump, same as Gauss but set to CDF (obs/theor)*(LHC theor) = 3.43pb
    cdfh = (TH1D *)gsh->Clone("CDFbump");

    cdfh->Scale(3.43);

    vcdfh[ichan] = cdfh;

    cdfh->Draw();

    gsh->Draw("same");

  } // channel loop

  //m_sigdata["gs"]  = vgsh;
  m_sigdata["cdf"] = vcdfh;
}                                                           // makeGaussianSignals

void peakAnalysis( string filename2 , string filenamelist ) {
   
    TFile * file0 = new TFile(filename2.c_str());
    TFile * file1 = new TFile(filenamelist.c_str());

    TTree * tree = (TTree*)file0->Get("tree");
    TEventList * listsel = (TEventList*)file1->Get("listofselected");

    TCanvas * can = new TCanvas( "canPeaks" , "Peak Analysis" , 7500 , 5500 , 900 , 600 );
    //can->SetLogy();
    can->SetGrid();
    can->cd();

    tree->SetEventList(listsel);
    tree->Draw("GEMDEnergyGauss_1.energy[0]>>htemp(10000,10000,10000)");
    TH1D * spectrum = (TH1D*)gDirectory->Get("htemp");

    spectrum->SetTitle("Energy Spectrum");
    spectrum->SetXTitle("Energy [keV]");
    spectrum->SetYTitle("Counts");
    //spectrum->SetStats(kFALSE);
    
    TSpectrum analyzer( 8 , 2 );
        // arg1: max number of peaks
        // arg2: resolution between peaks
    analyzer.Search( spectrum , 2 , "" , 0.0025);
        // arg2: sigma of the searched peaks
        // arg3: options
        // arg4: peaks with amplitude less than threshold*highest_peak are discarded

    vector<double> peaks;
    for ( int i = 0 ; i < analyzer.GetNPeaks() ; i++ ) {
        peaks.push_back(analyzer.GetPositionX()[i]);
    }
    sort(peaks.begin(),peaks.end());

    cout << endl << "########## PEAKS ##########" << endl;
    for ( int i = 0 ; i < peaks.size() ; i++ ) cout << i+1 << "\t" << peaks.at(i) << endl;

    cout << "###########################" << endl << endl;
   
	if ( peaks.size() < 6 ) { 
		cout << "ERROR: Not enough peaks found, try to modify TSpectrum parameters.\n\n";
		return;
	}

    ofstream file("calib.txt");

    file << "// calibration"                                     << endl
         << "// ch\tkeV"                                         << endl
         << peaks.at(0) << "\t1460.882 // 40K"                   << endl
         << peaks.at(1) << "\t1512.700 // 212Bi"                 << endl
         << peaks.at(2) << "\t1592.515 // 208Tl (double escape)" << endl
         << peaks.at(3) << "\t1620.738 // 212Bi"                 << endl
         << peaks.at(5) << "\t2103.513 // 208Tl (single escape)" << endl
         << peaks.at(6) << "\t2614.511 // 208Tl";

    return;
}

///////////////////////////////////////////////////////////////////
//////// Go4 GUI example script fft.C
//          J.Adamczewski, gsi, 30 May 2012
// NOTE: to be run in Go4 GUI local command line only!
//       NEVER call this script in remote analysis process!!!
/////// Functionality:
// perfroms fft on histogram of name1 using the option as explained in root TVirtualFFT:FFT
/////// Usage:
// The draw flag switches if the results are displayed each time this macro is called
// if display is switched off, the result histogram is just updated in browser and existing displays
///////
Bool_t fft(const char* name1, Option_t*  opt = "R2C M", Bool_t draw=kTRUE)
{
   if(TGo4AbstractInterface::Instance()==0 || go4!=TGo4AbstractInterface::Instance()) {
      std::cout <<"FATAL: Go4 gui macro executed outside Go4 GUI!! returning." << std::endl;
      return kFALSE;
   }
   TString newname;
   TString fullname1 = go4->FindItem(name1);
   TObject* ob1 = go4->GetObject(fullname1,1000); // 1000=timeout to get object from analysis in ms

   if ((ob1==0) || !ob1->InheritsFrom("TH1")) {
     std::cout <<"fft could not get histogram "<<fullname1 << std::endl;
     return kFALSE;
   }

   if(ob1->InheritsFrom("TH2") || ob1->InheritsFrom("TH3")){  // 2d
      std::cout <<"fft does not support 2d/3d histogram "<<fullname1 << std::endl;
      return kFALSE;
   }

   TH1* his1=(TH1*)ob1;
   TString n1=his1->GetName();
   TString t1=his1->GetTitle();
   newname.Form("_fft_%s",opt);
   TString finalname = n1+newname;
   TString finaltitle = t1+newname;


   // do fft here:
   Int_t N = his1->GetNbinsX();
   TH1D* result = new TH1D(finalname, finaltitle,N,0,N);
   result->SetName(finalname);
   result->SetTitle(finaltitle);
   result->Reset("");
   Double_t *in = new Double_t[N];
   // since we do not know type of input histo, we copy contents to Double array:
   for(Int_t ix=0; ix<N;++ix)
   {
      in[ix]=his1->GetBinContent(ix+1);
   }
   TVirtualFFT *thefft = TVirtualFFT::FFT(1, &N, opt);
   thefft->SetPoints(in);
   thefft->Transform();
   Double_t re, im;
   for (Int_t i=0; i<N; i++) {
      thefft->GetPointComplex(i, re, im);
      result->SetBinContent(i+1,TMath::Sqrt(re*re + im*im));
   }
   result->SetDirectory(0);

   TString rname = go4->SaveToMemory("FFT", result, kTRUE);
   std::cout<< "Saved result histogram to " << rname.Data() <<std::endl;
   if(draw){
      ViewPanelHandle vpanel = go4->StartViewPanel();
      go4->DrawItem(rname, vpanel);
   }
   return kTRUE;
}

void ProcYields::Proc_hPhi(hel_t hel/*=UNPOL*/){
	Info(TString::Format("Proc_hPhi(%s)",helTitle[hel].Data()), "");

	TDirectory* dir_phi=NULL;
	if (hel==UNPOL)	{
		dir_phi = _dirQ2W->mkdir("hPhi");
	}else if(hel==POS||hel==NEG){
		dir_phi = _dirQ2W->mkdir(TString::Format("hPhi_%s",helTitle[hel].Data()));
	}

	TDirectory* dir_varset=NULL;
	for (Int_t iVarset=0;iVarset<nVARSET;iVarset++){
		dir_varset = dir_phi->mkdir(TString::Format("Varset%d", iVarset+1));
		dir_varset->cd();
		//!Get relevent h5D
		THnSparse* hY5D = NULL;
		if (hel==UNPOL)	{
			//! To make hphi from ACC_CORR evts
			hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D/Varset%d/hY5D_ACC_CORR",iVarset+1));
			//! To make hphi from FULL evts
			//hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D/Varset%d/hY5D_FULL",iVarset+1));
			//! To make hphi from TH evts
			//hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D/Varset%d/hY5D_TH",iVarset+1));
		}else if (hel==POS||hel==NEG){
			//! To make hphi from ACC_CORR evts
			hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D_%s/Varset%d/hY5D_ACC_CORR",helTitle[hel].Data(),iVarset+1));
			//! To make hphi from FULL evts
			//hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D_%s/Varset%d/hY5D_FULL",helTitle[hel].Data(),iVarset+1));
			//! To make hphi from TH evts
			//hY5D = (THnSparseF*)_dirQ2W->Get(TString::Format("hY5D_%s/Varset%d/hY5D_TH",helTitle[hel].Data(),iVarset+1));
		}

		//! PROJECT phi bins for each VAR
		TDirectory* dir_var=NULL;
		for(Int_t iVar=0; iVar<nVAR; iVar++){
			if (iVar==ALPHA || iVar==M1 || iVar==M2 || iVar==PHI) continue; 
			dir_var = dir_varset->mkdir(varName[iVar].Data());
			dir_var->cd();
			int nvarbins = hY5D->GetAxis(iVar)->GetNbins();
			//! Loop over number of bins in Var and make phi projections for each
			for (int ivarbin=0; ivarbin<nvarbins; ivarbin++){
				Float_t varbin_lowedge = hY5D->GetAxis(iVar)->GetBinLowEdge(ivarbin+1);
				Float_t varbin_highedge = varbin_lowedge + hY5D->GetAxis(iVar)->GetBinWidth(ivarbin+1);
				//! Create Phi projection histogram
				hY5D->GetAxis(iVar)->SetRange(ivarbin+1, ivarbin+1);
				TH1D* hphi = (TH1D*)hY5D->Projection(PHI,"E");
				hphi->Sumw2();
				TString name = TString::Format("hphi_proj_%02d",ivarbin+1);
				TString title = TString::Format("#phi projection for %s = [%.2f,%.2f) | h=%s | top=%s | q2w = %s", 
					                            varTitle[iVarset][iVar].Data(), varbin_lowedge, varbin_highedge,
					                            helTitle[hel].Data(), _topList.Data(), _Q2Wdirname);
				hphi->SetName(name);
				hphi->SetTitle(title);
			}//end nvarbins loop
		}//end nVAR loop
	}//end nVARSET loop
	Info(TString::Format("Proc_hPhi(%s)",helTitle[hel].Data()), "done\n");
}

//--------------------------------------------------------------------------------------------------
void plotFrame(Double_t xMin, Double_t xMax, Double_t max)
{
  // Make a good frame
  TH1D * h = new TH1D("tmp","Time Series of Failures",1,xMin,xMax);
  MitStyle::InitHistWide(h,"","",kBlack);
  h->SetTitle("; Epoch Time [sec]; SmartCache Failures / 90 sec");
  h->SetMaximum(max*1.2);
  h->Draw("hist");
}

//--------------------------------------------------------------------------------------------------
void plotFrame(Double_t xMin, Double_t xMax, Double_t maxRate)
{
  // Make a good frame
  TH1D * h = new TH1D("tmp","Time Series of Rates",1,xMin,xMax);
  MitRootStyle::InitHistWide(h,"","",kBlack);  
  h->SetTitle("; Epoch Time [sec]; SmartCache transfer rate [MB/sec]");
  h->SetMaximum(maxRate*1.2);
  h->Draw("hist");
}

void plotter::draw_output_pseudo(TH1* output_, TH1D* pseudotruth_, TH1D* mctruth_, bool norm, TString file_name){

  TH1* output = (TH1*) output_->Clone("output");
  TH1D* pseudotruth = (TH1D*) pseudotruth_->Clone("pseudotruth");
  TH1D* mctruth = (TH1D*) mctruth_->Clone("mctruth");

  double ymax_temp = 0;
  if(pseudotruth->GetMaximum() > ymax_temp) ymax_temp = pseudotruth->GetMaximum();
  if(mctruth->GetMaximum() > ymax_temp) ymax_temp = mctruth->GetMaximum();
  if(output->GetMaximum() > ymax_temp) ymax_temp = output->GetMaximum();
  double ymax = 1.5 * ymax_temp;


  pseudotruth->SetTitle(" ");
  pseudotruth->GetYaxis()->SetRangeUser(0., ymax);
  pseudotruth->GetXaxis()->SetTitle("Leading-jet mass [GeV]");
  if(norm) pseudotruth->GetYaxis()->SetTitle("#frac{1}{#sigma} #frac{d#sigma}{dm_{jet}} [#frac{1}{GeV}]");
  else     pseudotruth->GetYaxis()->SetTitle("events");
  pseudotruth->GetYaxis()->SetTitleOffset(1.1);
  pseudotruth->GetXaxis()->SetTitleOffset(0.9);
  pseudotruth->GetYaxis()->SetTitleSize(0.05);
  pseudotruth->GetXaxis()->SetTitleSize(0.05);
  pseudotruth->GetYaxis()->SetNdivisions(505);

  pseudotruth->SetLineWidth(4);
  pseudotruth->SetLineColor(kRed);
  mctruth->SetLineWidth(3);
  mctruth->SetLineStyle(2);
  mctruth->SetLineColor(kBlue);

  output->SetLineColor(kBlack);
  output->SetMarkerColor(kBlack);
  output->SetMarkerStyle(8);
  output->SetMarkerSize(1);

  TCanvas *c= new TCanvas("Particle Level","",600,600);
  gPad->SetLeftMargin(0.15);
  TGaxis::SetMaxDigits(3);
  pseudotruth->Draw("HIST SAME");
  mctruth->Draw("HIST SAME");
  output->Draw("E1 SAME");
  TLegend *l;
  if(mctruth->GetSize() > 20) l=new TLegend(0.2,0.6,0.4,0.85);
  else                       l=new TLegend(0.55,0.6,0.85,0.85);
  l->SetBorderSize(0);
  l->SetFillStyle(0);
  l->AddEntry(output,"pseudo data","pl");
  l->AddEntry(pseudotruth,"pseudo data truth","pl");
  l->AddEntry(mctruth,"MC truth","pl");
  l->SetTextSize(0.04);

  l->Draw();
  gPad->RedrawAxis();
  c->SaveAs(directory + file_name + ".pdf");
  delete c;
}