[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / scripts / SimpledNdeta.C

#ifndef __CINT__
# include <TStyle.h>
# include <TFile.h>
# include <TList.h>
# include <TH1.h>
# include <THStack.h>
# include <TLegend.h>
# include <TLegendEntry.h>
# include <TLatex.h>
# include <TLine.h>
# include <TString.h>
# include <TCanvas.h>
# include <TError.h>
# include <TColor.h>
#else
class THStack;
class TAxis;
#endif

/**
 * A simple script to draw results from MakedNdeta.C (or similar)
 * 
 * @ingroup pwglf_forward_scripts
 */
/** 
 * Get a stack from the passed list 
 * 
 * @param list   List to get the stack from 
 * @param name   Name of stack 
 * @param rebin  Optional rebinning - must exists in list 
 * 
 * @return Stack or null
 * 
 * @ingroup pwglf_forward_scripts
 */
THStack*
GetStack(const TList* list, const char* name, Int_t rebin)
{
  if (!list) { 
    Warning("GetStack", "List is null");
    return 0;
  }

  TString n(name);
  if (rebin > 1) n.Append(Form("_rebin%02d", rebin));
  TObject* o = list->FindObject(n);
  if (!o) { 
    Warning("GetStack", "No %s object found in %s", n.Data(), list->GetName());
    return 0;
  }
  return static_cast<THStack*>(o);
}

/** 
 * Get a histogram from a list 
 * 
 * @param list   List 
 * @param name   Name of histogram
 * @param rebin  Rebinning factor
 * 
 * @return Histogram or null
 * 
 * @ingroup pwglf_forward_scripts
 */
TH1*
GetHist(const TList* list, const char* name, Int_t rebin)
{
  if (!list) { 
    Warning("GetStack", "List is null");
    return 0;
  }
  TList* all = static_cast<TList*>(list->FindObject("all"));
  if (!all) { 
    Warning("GetHist", "List all not found in %s", list->GetName());
    return 0;
  }

  TString n(name);
  if (rebin > 1) n.Append(Form("_rebin%02d", rebin));
  TObject* o = all->FindObject(n);
  if (!o) { 
    Warning("GetStack", "No %s object found in %s", n.Data(), all->GetName());
    return 0;
  }
  return static_cast<TH1*>(o);
}

/** 
 * Add histograms from one stack to another 
 * 
 * @param p      Parent stack to add to 
 * @param list   List to look for the stack in 
 * @param name   Name of stack to add 
 * @param rebin  Optional rebinning - must exists in list 
 * 
 * @return Added stack or null
 * 
 * @ingroup pwglf_forward_scripts
 */
THStack*
AddStack(THStack* p, const TList* list, const char* name, Int_t rebin)
{
  THStack* s = GetStack(list, name, rebin);
  if (!s) return 0;
  
  TIter next(s->GetHists());
  TH1*  hist = 0;
  while ((hist = static_cast<TH1*>(next()))) 
    p->Add(hist);
  return s;
}

/** 
 * Build up a centrality legend 
 * 
 * @param c Centrality axis 
 * 
 * @ingroup pwglf_forward_scripts
 */
void
BuildCentLegend(const TAxis* c)
{
  if (!c) return;
  TLegend* l = new TLegend(1.2*gPad->GetLeftMargin(), 
			   .55, .35, 1-gPad->GetTopMargin());
  l->SetFillColor(0);
  l->SetFillStyle(0);
  l->SetBorderSize(0);
  l->SetTextFont(22);
  l->SetHeader("Centralities");
  l->SetTextFont(132);
  
  Int_t    nCol     = gStyle->GetNumberOfColors();
  for (Int_t i = 0; i < c->GetNbins(); i++) {
    UShort_t low    = c->GetBinLowEdge(i+1);
    UShort_t high   = c->GetBinUpEdge(i+1);
    Float_t  fc     = (low+double(high-low)/2) / 100;
    Int_t    icol   = TMath::Min(nCol-1,int(fc * nCol + .5));
    Int_t    col    = gStyle->GetColorPalette(icol);
    TLegendEntry* e = l->AddEntry(Form("dummy%02d", i), 
				  Form("%3d%% - %3d%%", low, high), "p");
    e->SetMarkerColor(col);
				  
  }
  l->Draw();
}

/** 
 * Build a legend.  Histograms a filtered for the same title 
 * 
 * @param stack Stack of histograms 
 * @param c     Centrality axis.  If present, markers are black 
 * 
 * @ingroup pwglf_forward_scripts
 */
void
BuildLegend(const THStack* stack, const TAxis* c)
{
  Double_t x1 = .75, x2 = 1-gPad->GetRightMargin();
  Double_t y1 = .8,  y2 = 1-gPad->GetTopMargin();
  if (!c) { 
    // PP 
    x1 = .4; 
    y1 = .4;
    x2 = .8;
    y2 = .6;
  }
  TLegend* l = new TLegend(x1, y1, x2, y2, "");
  l->SetFillColor(0);
  l->SetFillStyle(0);
  l->SetBorderSize(0);
  l->SetTextFont(132);
  
  TObjArray seen;
  seen.SetOwner();
  TIter next(stack->GetHists());
  TH1* hist = 0;
  while ((hist = static_cast<TH1*>(next()))) { 
    TString n(hist->GetTitle());
    if (n.Contains("mirrored")) continue;
    if (seen.FindObject(n.Data())) continue;
    TObjString* ns = new TObjString(n);
    ns->SetUniqueID(((hist->GetMarkerStyle() & 0xFFFF) << 16) |
		    ((hist->GetMarkerColor() & 0xFFFF) <<  0));
    seen.Add(ns);
  }
  seen.Sort();
  TIter nextu(&seen);
  TObject* s = 0;
  Int_t i = 0;
  while ((s = nextu())) { 
    TLegendEntry* dd = l->AddEntry(Form("data%2d", i++), 
				   s->GetName(), "p");
    Int_t style = (s->GetUniqueID() >> 16) & 0xFFFF;
    Int_t color = (s->GetUniqueID() >>  0) & 0xFFFF;
    dd->SetLineColor(kBlack);
    if (c) dd->SetMarkerColor(kBlack);
    else   dd->SetMarkerColor(color);
    dd->SetMarkerStyle(style);
    Int_t st = dd->GetMarkerStyle();
    if (st == 27 || st == 28 || st == 29 || st == 30 || st == 33 || st == 34) 
      dd->SetMarkerSize(1.4*dd->GetMarkerSize());
  }
  // TLegendEntry* sep = l->AddEntry("s", "_", "h");
  // sep->SetTextSize(0.01);
  TLine* sep = new TLine(0,0,1,1);
  sep->SetLineWidth(1);
  sep->DrawLineNDC(x1+.02, y1-.005, x2-.03, y1-.005);
  // sep->Draw();
  
  TLegend* l2 = new TLegend(x1, y1-.005, x2, y1-.15, "");
  l2->SetFillColor(0);
  l2->SetFillStyle(0);
  l2->SetBorderSize(0);
  l2->SetTextFont(132);
  
  TLegendEntry* d1 = l2->AddEntry("d1", "Data", "p");
  d1->SetLineColor(kBlack);
  d1->SetMarkerColor(kBlack);
  d1->SetMarkerStyle(20);
  TLegendEntry* d2 = l2->AddEntry("d2", "Mirrored data", "p");
  d2->SetLineColor(kBlack);
  d2->SetMarkerColor(kBlack);
  d2->SetMarkerStyle(24);
  
  l->Draw();
  l2->Draw();
}

/** 
 * Add additional information
 *  
 * @param forward  List of info
 * @param prelim   Preliminary mark 
 */
void
AddInformation(TList* forward, bool prelim=true)
{
  Double_t x  = .12;
  Double_t y  = .95;
  Double_t ts = 0.05;
  if (prelim) {
    TLatex* wip = new TLatex(x, y, "Work in progress");
    wip->SetNDC();
    wip->SetTextColor(TColor::GetColor(234,26,46));
    wip->SetTextAlign(13);
    wip->SetTextFont(132);
    wip->SetTextSize(ts);
    wip->Draw();
    y -= ts;
  }

  TObject* sNN = forward->FindObject("sNN");
  TObject* sys = forward->FindObject("sys");
  TObject* trg = forward->FindObject("trigger");
  TObject* vtx = forward->FindObject("vtxAxis");
  TObject* sch = forward->FindObject("scheme");

  TString t(sys->GetTitle());
  Bool_t isPP = t == "pp";

  
  TString s = Form("%s @ #sqrt{s%s}=", 
		   sys->GetTitle(),
		   (isPP ? "" : "_{NN}"));
  Int_t cms = sNN->GetUniqueID();
  if (cms > 1000) s.Append(Form("%5.2fTeV", float(cms)/1000));
  else            s.Append(Form("%3dGeV", cms));
  s.Append(Form(", %s", trg->GetTitle()));

  // if (isPP) { x = .3; y = .3; }
  if (isPP) { 
    x  = 1-gPad->GetRightMargin()-.01;
    y  = 1-gPad->GetTopMargin()-.01;
    ts = .04;
  }
  TLatex* ltx = new TLatex(x, y, s.Data());
  ltx->SetNDC();
  ltx->SetTextColor(TColor::GetColor(41,73,156));
  ltx->SetTextAlign((isPP ? 33 : 13));
  ltx->SetTextFont(132);
  ltx->SetTextSize(ts);
  ltx->Draw();
  y -= ts;
  ltx->DrawLatex(x, y, vtx->GetTitle());
  y -= ts;
  ltx->DrawLatex(x, y, sch->GetTitle());
}  

/** 
 * A function (double Gaussian)
 * 
 * @param xp Independent variables
 * @param pp Parameters 
 * 
 * @return Value of function
 * 
 * @ingroup pwglf_forward_scripts
 */
Double_t myFunc(Double_t* xp, Double_t* pp)
{
  Double_t x  = xp[0];
  Double_t a1 = pp[0];
  Double_t a2 = pp[1];
  Double_t s1 = pp[2];
  Double_t s2 = pp[3];
  return a1*(TMath::Gaus(x, 0, s1) - a2 * TMath::Gaus(x, 0, s2));
}

/** 
 * Make systematic error band 
 * 
 * @param cen     Central result
 * @param fwd     Forward result
 * @param sysErr  Systematic error (fractional)
 * 
 * @return 
 * 
 * @ingroup pwglf_forward_scripts
 */
TH1* 
MakeSysError(const TH1* cen, const TH1* fwd, Double_t sysErr=0.7)
{
  TH1* tmp = static_cast<TH1*>(fwd->Clone("dndetaFitted"));
  tmp->SetMarkerStyle(0);
  tmp->SetFillColor(kGray);
  tmp->SetFillStyle(3001);
  tmp->SetDirectory(0);
  Double_t xlow  = 100;
  Double_t xhigh = 0;
  for (Int_t i = 1; i <= tmp->GetNbinsX(); i++) {
    Double_t cc = cen->GetBinContent(i);
    Double_t cf = fwd->GetBinContent(i);
    Double_t ec = fwd->GetBinError(i);
    Double_t ef = fwd->GetBinError(i);
    Double_t nc = cf;
    Double_t ne = ef;
    if (cc < 0.001 && cf < 0.01) continue;
    xlow  = TMath::Min(tmp->GetXaxis()->GetBinLowEdge(i),xlow);
    xhigh = TMath::Max(tmp->GetXaxis()->GetBinUpEdge(i),xhigh);
    if (cc > 0.001) {
      nc = cc;
      ne = ec;
      if (cf > 0.001) {
	nc  = (cf + cc) / 2;
	ne  = TMath::Sqrt(ec*ec + ef*ef);
      }
    }
    tmp->SetBinContent(i, nc);
    tmp->SetBinError(i, ne);
  }
  TF1* tmpf  = new TF1("tmpf", "gaus", xlow, xhigh);
  tmp->Fit(tmpf, "NQ", "");
  
  TF1* fit = new TF1("f", myFunc, xlow, xhigh, 4);
  fit->SetParNames("a_{1}", "a_{2}", "#sigma_{1}", "#sigma_{2}");
  fit->SetParameters(tmpf->GetParameter(0), 
		     .2, 
		     tmpf->GetParameter(2), 
		     tmpf->GetParameter(2)/4);
  tmp->Fit(fit,"0WQ","");
  for (Int_t i = 1; i <= tmp->GetNbinsX(); i++) {
    Double_t tc = tmp->GetBinContent(i);
    if (tc < 0.01) continue;
    Double_t x = tmp->GetXaxis()->GetBinCenter(i);
    Double_t y = fit->Eval(x);
    tmp->SetBinContent(i, y);
    tmp->SetBinError(i,sysErr*y);
  }
  return tmp;
}

/** 
 * Function to draw the results from forward_dndeta.root file 
 * 
 * @param what     What to draw 
 * @param rebin    Rebinnig.  Note, the data must be present in the file
 * @param filename File to open and draw stuff from >
 * 
 * @ingroup pwglf_forward_scripts
 */
void
SimpledNdeta(Int_t what=0x5, 
	     Int_t rebin=5, const char* filename="forward_dndeta.root")
{
  TFile* file = TFile::Open(filename, "READ");
  if (!file) { 
    Error("SimpledNdeta", "File %s not found", filename);
    return;
  }

  TList* forward = static_cast<TList*>(file->Get("ForwardResults"));
  TList* central = static_cast<TList*>(file->Get("CentralResults"));
  TList* mctruth = static_cast<TList*>(file->Get("MCTruthResults"));

  THStack* all = new THStack("all", "All");
  if (what & 0x1) AddStack(all, forward, "dndeta",      rebin);
  if (what & 0x2) AddStack(all, forward, "dndetaMC",    rebin);
  if (what & 0x1) AddStack(all, central, "dndeta",      rebin);
  if (what & 0x2) AddStack(all, central, "dndetaMC",    rebin);
  if (what & 0x4) AddStack(all, mctruth, "dndeta",      rebin);
  
  TH1* tmp = 0;
  if (what & 0x1) { 
    Double_t sysErr = 0.07;
    TH1* fwd = GetHist(forward, "dndetaForward", rebin);
    TH1* cen = GetHist(central, "dndetaCentral", rebin);
    tmp = MakeSysError(cen, fwd, sysErr);
    // fit = tmpf;
  }
  
  TAxis* centAxis = static_cast<TAxis*>(forward->FindObject("centAxis"));
  Bool_t isPP = centAxis == 0;

  gStyle->SetPalette(1);
  gStyle->SetOptTitle(0);
  TCanvas* c = new TCanvas("dndeta", "dN/deta results", 900, 600);
  c->SetFillColor(0);
  c->SetFillStyle(0);
  c->SetBorderSize(0);
  c->SetBorderMode(0);
  c->SetTopMargin(0.03);
  c->SetRightMargin(0.03);
  
  all->Draw("nostack");
  TAxis* xa = all->GetHistogram()->GetXaxis();
  xa->SetTitleFont(132);
  xa->SetLabelFont(132);
  xa->SetTitle("#eta");
  TAxis* ya = all->GetHistogram()->GetYaxis();
  ya->SetTitleFont(132);
  ya->SetLabelFont(132);
  ya->SetTitle("dN_{ch}/d#eta");

  if (tmp) { 
    tmp->Draw("e5 same");
    all->Draw("same nostack");
  }
  

  // if (fit) fit->Draw("same");
  // if (tmp) tmp->Draw("same");

  BuildCentLegend(centAxis);
  BuildLegend(all, centAxis);

  AddInformation(forward);

  c->SaveAs("dndeta_simple.C");
  c->SaveAs("dndeta_simple.png");
  c->SaveAs("dndeta_simple.root");
}
//
// EOF
//
Commit	Line	Data
0919aa56	1	#ifndef __CINT__
	2	# include <TStyle.h>
	3	# include <TFile.h>
	4	# include <TList.h>
	5	# include <TH1.h>
	6	# include <THStack.h>
	7	# include <TLegend.h>
	8	# include <TLegendEntry.h>
797161e8	9	# include <TLatex.h>
0919aa56	10	# include <TLine.h>
	11	# include <TString.h>
	12	# include <TCanvas.h>
	13	# include <TError.h>
797161e8	14	# include <TColor.h>
0919aa56	15	#else
	16	class THStack;
	17	class TAxis;
	18	#endif
	19
	20	/**
	21	* A simple script to draw results from MakedNdeta.C (or similar)
	22	*
bd6f5206	23	* @ingroup pwglf_forward_scripts
0919aa56	24	*/
	25	/**
	26	* Get a stack from the passed list
	27	*
	28	* @param list List to get the stack from
	29	* @param name Name of stack
	30	* @param rebin Optional rebinning - must exists in list
	31	*
	32	* @return Stack or null
56199f2b	33	*
bd6f5206	34	* @ingroup pwglf_forward_scripts
0919aa56	35	*/
	36	THStack*
	37	GetStack(const TList* list, const char* name, Int_t rebin)
	38	{
	39	if (!list) {
	40	Warning("GetStack", "List is null");
	41	return 0;
	42	}
	43
	44	TString n(name);
	45	if (rebin > 1) n.Append(Form("_rebin%02d", rebin));
	46	TObject* o = list->FindObject(n);
	47	if (!o) {
	48	Warning("GetStack", "No %s object found in %s", n.Data(), list->GetName());
	49	return 0;
	50	}
	51	return static_cast<THStack*>(o);
	52	}
	53
56199f2b	54	/**
	55	* Get a histogram from a list
	56	*
	57	* @param list List
	58	* @param name Name of histogram
	59	* @param rebin Rebinning factor
	60	*
	61	* @return Histogram or null
	62	*
bd6f5206	63	* @ingroup pwglf_forward_scripts
56199f2b	64	*/
797161e8	65	TH1*
	66	GetHist(const TList* list, const char* name, Int_t rebin)
	67	{
	68	if (!list) {
	69	Warning("GetStack", "List is null");
	70	return 0;
	71	}
	72	TList* all = static_cast<TList*>(list->FindObject("all"));
	73	if (!all) {
	74	Warning("GetHist", "List all not found in %s", list->GetName());
	75	return 0;
	76	}
	77
	78	TString n(name);
	79	if (rebin > 1) n.Append(Form("_rebin%02d", rebin));
	80	TObject* o = all->FindObject(n);
	81	if (!o) {
	82	Warning("GetStack", "No %s object found in %s", n.Data(), all->GetName());
	83	return 0;
	84	}
	85	return static_cast<TH1*>(o);
	86	}
	87
0919aa56	88	/**
	89	* Add histograms from one stack to another
	90	*
	91	* @param p Parent stack to add to
	92	* @param list List to look for the stack in
	93	* @param name Name of stack to add
	94	* @param rebin Optional rebinning - must exists in list
	95	*
	96	* @return Added stack or null
56199f2b	97	*
bd6f5206	98	* @ingroup pwglf_forward_scripts
0919aa56	99	*/
	100	THStack*
	101	AddStack(THStack* p, const TList* list, const char* name, Int_t rebin)
	102	{
	103	THStack* s = GetStack(list, name, rebin);
	104	if (!s) return 0;
	105
	106	TIter next(s->GetHists());
	107	TH1* hist = 0;
	108	while ((hist = static_cast<TH1*>(next())))
	109	p->Add(hist);
	110	return s;
	111	}
	112
	113	/**
	114	* Build up a centrality legend
	115	*
	116	* @param c Centrality axis
56199f2b	117	*
bd6f5206	118	* @ingroup pwglf_forward_scripts
0919aa56	119	*/
	120	void
	121	BuildCentLegend(const TAxis* c)
	122	{
	123	if (!c) return;
	124	TLegend* l = new TLegend(1.2*gPad->GetLeftMargin(),
	125	.55, .35, 1-gPad->GetTopMargin());
	126	l->SetFillColor(0);
	127	l->SetFillStyle(0);
	128	l->SetBorderSize(0);
	129	l->SetTextFont(22);
	130	l->SetHeader("Centralities");
	131	l->SetTextFont(132);
	132
	133	Int_t nCol = gStyle->GetNumberOfColors();
	134	for (Int_t i = 0; i < c->GetNbins(); i++) {
	135	UShort_t low = c->GetBinLowEdge(i+1);
	136	UShort_t high = c->GetBinUpEdge(i+1);
	137	Float_t fc = (low+double(high-low)/2) / 100;
	138	Int_t icol = TMath::Min(nCol-1,int(fc * nCol + .5));
	139	Int_t col = gStyle->GetColorPalette(icol);
	140	TLegendEntry* e = l->AddEntry(Form("dummy%02d", i),
	141	Form("%3d%% - %3d%%", low, high), "p");
	142	e->SetMarkerColor(col);
	143
	144	}
	145	l->Draw();
	146	}
	147
	148	/**
	149	* Build a legend. Histograms a filtered for the same title
	150	*
	151	* @param stack Stack of histograms
	152	* @param c Centrality axis. If present, markers are black
56199f2b	153	*
bd6f5206	154	* @ingroup pwglf_forward_scripts
0919aa56	155	*/
	156	void
	157	BuildLegend(const THStack* stack, const TAxis* c)
	158	{
797161e8	159	Double_t x1 = .75, x2 = 1-gPad->GetRightMargin();
	160	Double_t y1 = .8, y2 = 1-gPad->GetTopMargin();
	161	if (!c) {
	162	// PP
	163	x1 = .4;
	164	y1 = .4;
	165	x2 = .8;
	166	y2 = .6;
	167	}
	168	TLegend* l = new TLegend(x1, y1, x2, y2, "");
0919aa56	169	l->SetFillColor(0);
	170	l->SetFillStyle(0);
	171	l->SetBorderSize(0);
	172	l->SetTextFont(132);
	173
	174	TObjArray seen;
	175	seen.SetOwner();
	176	TIter next(stack->GetHists());
	177	TH1* hist = 0;
	178	while ((hist = static_cast<TH1*>(next()))) {
	179	TString n(hist->GetTitle());
	180	if (n.Contains("mirrored")) continue;
	181	if (seen.FindObject(n.Data())) continue;
	182	TObjString* ns = new TObjString(n);
	183	ns->SetUniqueID(((hist->GetMarkerStyle() & 0xFFFF) << 16) \|
	184	((hist->GetMarkerColor() & 0xFFFF) << 0));
	185	seen.Add(ns);
	186	}
	187	seen.Sort();
	188	TIter nextu(&seen);
	189	TObject* s = 0;
	190	Int_t i = 0;
	191	while ((s = nextu())) {
	192	TLegendEntry* dd = l->AddEntry(Form("data%2d", i++),
	193	s->GetName(), "p");
	194	Int_t style = (s->GetUniqueID() >> 16) & 0xFFFF;
	195	Int_t color = (s->GetUniqueID() >> 0) & 0xFFFF;
	196	dd->SetLineColor(kBlack);
	197	if (c) dd->SetMarkerColor(kBlack);
	198	else dd->SetMarkerColor(color);
	199	dd->SetMarkerStyle(style);
797161e8	200	Int_t st = dd->GetMarkerStyle();
	201	if (st == 27 \|\| st == 28 \|\| st == 29 \|\| st == 30 \|\| st == 33 \|\| st == 34)
	202	dd->SetMarkerSize(1.4*dd->GetMarkerSize());
0919aa56	203	}
	204	// TLegendEntry* sep = l->AddEntry("s", "_", "h");
	205	// sep->SetTextSize(0.01);
	206	TLine* sep = new TLine(0,0,1,1);
	207	sep->SetLineWidth(1);
797161e8	208	sep->DrawLineNDC(x1+.02, y1-.005, x2-.03, y1-.005);
0919aa56	209	// sep->Draw();
0919aa56	210
797161e8	211	TLegend* l2 = new TLegend(x1, y1-.005, x2, y1-.15, "");
0919aa56	212	l2->SetFillColor(0);
	213	l2->SetFillStyle(0);
	214	l2->SetBorderSize(0);
	215	l2->SetTextFont(132);
	216
	217	TLegendEntry* d1 = l2->AddEntry("d1", "Data", "p");
	218	d1->SetLineColor(kBlack);
	219	d1->SetMarkerColor(kBlack);
	220	d1->SetMarkerStyle(20);
	221	TLegendEntry* d2 = l2->AddEntry("d2", "Mirrored data", "p");
	222	d2->SetLineColor(kBlack);
	223	d2->SetMarkerColor(kBlack);
	224	d2->SetMarkerStyle(24);
797161e8	225
0919aa56	226	l->Draw();
	227	l2->Draw();
	228	}
797161e8	229
56199f2b	230	/**
	231	* Add additional information
	232	*
	233	* @param forward List of info
	234	* @param prelim Preliminary mark
	235	*/
797161e8	236	void
	237	AddInformation(TList* forward, bool prelim=true)
	238	{
	239	Double_t x = .12;
	240	Double_t y = .95;
	241	Double_t ts = 0.05;
	242	if (prelim) {
	243	TLatex* wip = new TLatex(x, y, "Work in progress");
	244	wip->SetNDC();
	245	wip->SetTextColor(TColor::GetColor(234,26,46));
	246	wip->SetTextAlign(13);
	247	wip->SetTextFont(132);
	248	wip->SetTextSize(ts);
	249	wip->Draw();
	250	y -= ts;
	251	}
	252
	253	TObject* sNN = forward->FindObject("sNN");
	254	TObject* sys = forward->FindObject("sys");
	255	TObject* trg = forward->FindObject("trigger");
	256	TObject* vtx = forward->FindObject("vtxAxis");
	257	TObject* sch = forward->FindObject("scheme");
	258
	259	TString t(sys->GetTitle());
	260	Bool_t isPP = t == "pp";
	261
0919aa56	262
797161e8	263	TString s = Form("%s @ #sqrt{s%s}=",
	264	sys->GetTitle(),
	265	(isPP ? "" : "_{NN}"));
	266	Int_t cms = sNN->GetUniqueID();
	267	if (cms > 1000) s.Append(Form("%5.2fTeV", float(cms)/1000));
	268	else s.Append(Form("%3dGeV", cms));
	269	s.Append(Form(", %s", trg->GetTitle()));
	270
	271	// if (isPP) { x = .3; y = .3; }
	272	if (isPP) {
	273	x = 1-gPad->GetRightMargin()-.01;
	274	y = 1-gPad->GetTopMargin()-.01;
	275	ts = .04;
	276	}
	277	TLatex* ltx = new TLatex(x, y, s.Data());
	278	ltx->SetNDC();
	279	ltx->SetTextColor(TColor::GetColor(41,73,156));
	280	ltx->SetTextAlign((isPP ? 33 : 13));
	281	ltx->SetTextFont(132);
	282	ltx->SetTextSize(ts);
	283	ltx->Draw();
	284	y -= ts;
	285	ltx->DrawLatex(x, y, vtx->GetTitle());
	286	y -= ts;
	287	ltx->DrawLatex(x, y, sch->GetTitle());
	288	}
	289
56199f2b	290	/**
	291	* A function (double Gaussian)
	292	*
	293	* @param xp Independent variables
	294	* @param pp Parameters
	295	*
	296	* @return Value of function
	297	*
bd6f5206	298	* @ingroup pwglf_forward_scripts
56199f2b	299	*/
797161e8	300	Double_t myFunc(Double_t* xp, Double_t* pp)
	301	{
	302	Double_t x = xp[0];
	303	Double_t a1 = pp[0];
	304	Double_t a2 = pp[1];
	305	Double_t s1 = pp[2];
	306	Double_t s2 = pp[3];
	307	return a1(TMath::Gaus(x, 0, s1) - a2 TMath::Gaus(x, 0, s2));
	308	}
	309
56199f2b	310	/**
	311	* Make systematic error band
	312	*
	313	* @param cen Central result
	314	* @param fwd Forward result
	315	* @param sysErr Systematic error (fractional)
	316	*
	317	* @return
	318	*
bd6f5206	319	* @ingroup pwglf_forward_scripts
56199f2b	320	*/
797161e8	321	TH1*
	322	MakeSysError(const TH1* cen, const TH1* fwd, Double_t sysErr=0.7)
	323	{
	324	TH1* tmp = static_cast<TH1*>(fwd->Clone("dndetaFitted"));
	325	tmp->SetMarkerStyle(0);
	326	tmp->SetFillColor(kGray);
	327	tmp->SetFillStyle(3001);
	328	tmp->SetDirectory(0);
	329	Double_t xlow = 100;
	330	Double_t xhigh = 0;
	331	for (Int_t i = 1; i <= tmp->GetNbinsX(); i++) {
	332	Double_t cc = cen->GetBinContent(i);
	333	Double_t cf = fwd->GetBinContent(i);
	334	Double_t ec = fwd->GetBinError(i);
	335	Double_t ef = fwd->GetBinError(i);
	336	Double_t nc = cf;
	337	Double_t ne = ef;
	338	if (cc < 0.001 && cf < 0.01) continue;
	339	xlow = TMath::Min(tmp->GetXaxis()->GetBinLowEdge(i),xlow);
	340	xhigh = TMath::Max(tmp->GetXaxis()->GetBinUpEdge(i),xhigh);
	341	if (cc > 0.001) {
	342	nc = cc;
	343	ne = ec;
	344	if (cf > 0.001) {
	345	nc = (cf + cc) / 2;
	346	ne = TMath::Sqrt(ecec + efef);
	347	}
	348	}
	349	tmp->SetBinContent(i, nc);
	350	tmp->SetBinError(i, ne);
	351	}
	352	TF1* tmpf = new TF1("tmpf", "gaus", xlow, xhigh);
	353	tmp->Fit(tmpf, "NQ", "");
	354
	355	TF1* fit = new TF1("f", myFunc, xlow, xhigh, 4);
	356	fit->SetParNames("a_{1}", "a_{2}", "#sigma_{1}", "#sigma_{2}");
	357	fit->SetParameters(tmpf->GetParameter(0),
	358	.2,
	359	tmpf->GetParameter(2),
	360	tmpf->GetParameter(2)/4);
	361	tmp->Fit(fit,"0WQ","");
	362	for (Int_t i = 1; i <= tmp->GetNbinsX(); i++) {
	363	Double_t tc = tmp->GetBinContent(i);
	364	if (tc < 0.01) continue;
	365	Double_t x = tmp->GetXaxis()->GetBinCenter(i);
	366	Double_t y = fit->Eval(x);
	367	tmp->SetBinContent(i, y);
	368	tmp->SetBinError(i,sysErr*y);
	369	}
	370	return tmp;
	371	}
0919aa56	372
	373	/**
	374	* Function to draw the results from forward_dndeta.root file
	375	*
970b1a8a	376	* @param what What to draw
0919aa56	377	* @param rebin Rebinnig. Note, the data must be present in the file
0919aa56	378	* @param filename File to open and draw stuff from >
56199f2b	379	*
bd6f5206	380	* @ingroup pwglf_forward_scripts
0919aa56	381	*/
0919aa56	382	void
797161e8	383	SimpledNdeta(Int_t what=0x5,
797161e8	384	Int_t rebin=5, const char* filename="forward_dndeta.root")
0919aa56	385	{
	386	TFile* file = TFile::Open(filename, "READ");
	387	if (!file) {
	388	Error("SimpledNdeta", "File %s not found", filename);
	389	return;
	390	}
	391
	392	TList* forward = static_cast<TList*>(file->Get("ForwardResults"));
	393	TList* central = static_cast<TList*>(file->Get("CentralResults"));
797161e8	394	TList* mctruth = static_cast<TList*>(file->Get("MCTruthResults"));
0919aa56	395
0919aa56	396	THStack* all = new THStack("all", "All");
797161e8	397	if (what & 0x1) AddStack(all, forward, "dndeta", rebin);
	398	if (what & 0x2) AddStack(all, forward, "dndetaMC", rebin);
	399	if (what & 0x1) AddStack(all, central, "dndeta", rebin);
	400	if (what & 0x2) AddStack(all, central, "dndetaMC", rebin);
	401	if (what & 0x4) AddStack(all, mctruth, "dndeta", rebin);
	402
	403	TH1* tmp = 0;
	404	if (what & 0x1) {
	405	Double_t sysErr = 0.07;
	406	TH1* fwd = GetHist(forward, "dndetaForward", rebin);
	407	TH1* cen = GetHist(central, "dndetaCentral", rebin);
	408	tmp = MakeSysError(cen, fwd, sysErr);
	409	// fit = tmpf;
	410	}
0919aa56	411
0919aa56	412	TAxis* centAxis = static_cast<TAxis*>(forward->FindObject("centAxis"));
797161e8	413	Bool_t isPP = centAxis == 0;
0919aa56	414
	415	gStyle->SetPalette(1);
	416	gStyle->SetOptTitle(0);
	417	TCanvas* c = new TCanvas("dndeta", "dN/deta results", 900, 600);
	418	c->SetFillColor(0);
	419	c->SetFillStyle(0);
	420	c->SetBorderSize(0);
	421	c->SetBorderMode(0);
	422	c->SetTopMargin(0.03);
	423	c->SetRightMargin(0.03);
	424
	425	all->Draw("nostack");
	426	TAxis* xa = all->GetHistogram()->GetXaxis();
	427	xa->SetTitleFont(132);
	428	xa->SetLabelFont(132);
	429	xa->SetTitle("#eta");
	430	TAxis* ya = all->GetHistogram()->GetYaxis();
	431	ya->SetTitleFont(132);
	432	ya->SetLabelFont(132);
	433	ya->SetTitle("dN_{ch}/d#eta");
	434
797161e8	435	if (tmp) {
	436	tmp->Draw("e5 same");
	437	all->Draw("same nostack");
	438	}
	439
	440
	441	// if (fit) fit->Draw("same");
	442	// if (tmp) tmp->Draw("same");
	443
0919aa56	444	BuildCentLegend(centAxis);
	445	BuildLegend(all, centAxis);
	446
797161e8	447	AddInformation(forward);
	448
	449	c->SaveAs("dndeta_simple.C");
	450	c->SaveAs("dndeta_simple.png");
	451	c->SaveAs("dndeta_simple.root");
0919aa56	452	}
	453	//
	454	// EOF
	455	//