AliPhysics  9b6b435 (9b6b435)
AliTrackletdNdetaPost.C
Go to the documentation of this file.
1 #ifndef ALITRACKLETDNDETAPOST_H
2 #define ALITRACKLETDNDETAPOST_H
4 #ifndef __CINT__
5 #include <THStack.h>
6 #include <TFile.h>
7 #include <TLatex.h>
8 #include <TCanvas.h>
9 #include <TLegend.h>
10 #include <TError.h>
11 #include <TParameter.h>
12 #include <TMath.h>
13 #include <TLine.h>
14 #include <TLegendEntry.h>
15 #include <TStyle.h>
16 #include <TBrowser.h>
17 #include <TGraphAsymmErrors.h>
18 #include <TProfile.h>
19 #include <TFitResult.h>
20 #include <TF1.h>
21 #else
22 class TPad;
23 class TLatex;
24 class TObject;
25 class TSeqCollection;
26 class TH1;
27 class TH2;
28 class TF1;
29 class TFitResultPtr;
30 class THStack;
31 class TCanvas;
32 class TVirtualPad;
33 class TFile;
34 class TAxis;
35 class TLegend;
36 class TDirectory;
37 #endif
38 
39 //====================================================================
41 {
46  enum {
48  kRealComb = 0x00001,
50  kSimComb = 0x00002,
52  kdNdetas = 0x00004,
54  kGeneral = 0x00008,
56  kParameters = 0x00010,
58  kSimInfo = 0x00020,
60  kSpecies = 0x00040,
62  kDeltas = 0x00080,
64  kBackgrounds = 0x00100,
66  kAlphas = 0x00200,
68  kClosure = 0x00400,
70  kPDF = 0x01000,
72  kPause = 0x02000,
74  kLandscape = 0x04000,
76  kFixedScale = 0x08000,
78  kDynCentScale = 0x10000,
80  kAltMarker = 0x20000,
82  kDefault = 0x133ff
83  };
88  };
89  enum {
90  kTopBackground = kAzure-8
91  };
92  //==================================================================
94  TCanvas* fCanvas;
96  TPad* fTop;
98  TPad* fBody;
110  TLatex* fHeader;
137  void ClearCanvas();
143  void CreateCanvas(const TString& outputName);
148  void CloseCanvas();
155  void PrintCanvas(const TString& title, Float_t size=.7);
156  //____________________________________________________________________
167  TLegend* DrawInPad(TVirtualPad* c, Int_t pad, TObject* o, Option_t* opt);
168  /* @} */
169  //____________________________________________________________________
181  TFile* OpenFile(const char* filename);
195  TH1* SetAttr(TH1* h,
196  Color_t color,
197  Style_t marker=20,
198  Double_t size=1.,
199  Style_t fill=0,
200  Style_t line=1,
201  Width_t width=1);
202  /* @} */
203  //____________________________________________________________________
215  void DrawParam(const char* name, Double_t& y, const char* val);
223  void DrawParam(const char* name, Double_t& y, Double_t val);
231  void DrawParam(const char* name, Double_t& y, Int_t val);
239  void DrawParam(const char* name, Double_t& y, Bool_t val);
247  void DrawParams(Container* pars, const char* title, Bool_t comb);
256  void DrawParams(Container* realSums, Bool_t realComb,
257  Container* simSums, Bool_t simComb);
258  /* @} */
259  //____________________________________________________________________
276  THStack* Make2Stack(const char* name,
277  const char* title,
278  Container* realList,
279  Container* simList,
280  Option_t* dataOpt="",
281  Option_t* simOpt="");
290  void DrawClusters(Int_t pad, const char* name,
291  Container* realList, Container* simList);
298  void DrawGeneral(Container* realList, Container* simList);
299  /* @} */
300  //____________________________________________________________________
313  TH1* FindDelta(Container* ress, const char* sub, Bool_t scaled=false);
323  TH1* FindSub(Container* ress, const char* sub, const char* name);
329  void DrawSpecies(Container* mcList);
336  void DrawDeltas(Container* realList, Container* simList);
348  Double_t DrawBackground(TVirtualPad* c,
349  Int_t pad,
350  Container* ress,
351  const char* name,
352  const char* pretty);
359  void DrawBackground(Container* realList, Container* simList);
367  TH2* CutAlpha(TH2* alpha);
377  void DrawAlpha(TVirtualPad* c,
378  Int_t pad,
379  Container* ress,
380  const char* name,
381  const char* pretty);
388  void DrawAlpha(Container* realList, Container* simList);
394  void DrawSim(Container* simList);
397  //__________________________________________________________________
407  const char* ObsTitle() const { return "d#it{N}_{ch}/d#eta"; }
426  Style_t MS(Int_t what, Bool_t sim, Bool_t alt) const;
434  TH1* GetDeltaScaleDist(Container* c) const;
442  TH1* GetDeltaIntDist(Container* c) const;
520  void DrawdNdeta(Container* realList,
521  Bool_t realComb,
522  Container* simList,
523  Bool_t simComb,
524  Color_t color=kBlack,
525  THStack* stack=0,
526  TDirectory* out=0,
527  Bool_t alt=false);
528  //====================================================================
541  void ProcessBin(UInt_t what, Int_t bin,
542  Double_t c1, Double_t c2,
543  Container* realList, Container* simList,
544  THStack* stack=0, TDirectory* out=0);
545  //====================================================================
554  void Run(UInt_t what = kDefault,
555  UShort_t maxBins = 9,
556  const char* dataName = "data.root",
557  const char* simName = "sim.root");
558 };
559 
560 
561 //____________________________________________________________________
564  fCanvas(0),
565  fTop(0),
566  fBody(0),
567  fLandscape(false),
568  fPDF(false),
569  fPause(false),
570  fLastTitle(""),
571  fLastBin(""),
572  fHeader(0),
573  fDeltaCut(0),
574  fTailDelta(0),
575  fMaxDelta(0),
576  fCombinatoricsScale(1.3),
578  fAlphaMin(0),
579  fAlphaMax(2.5)
580 {}
581 //====================================================================
582 struct SuppressGuard
583 {
584  Int_t save = 0;
586  {
587  save = gErrorIgnoreLevel;
588  gErrorIgnoreLevel = lvl;
589  }
591  {
592  gErrorIgnoreLevel = save;
593  }
594 };
595 
596 //____________________________________________________________________
598 {
599  fTop->Clear();
600  fTop->SetNumber(1);
601  fTop->SetFillColor(kTopBackground);
602  fTop->SetFillStyle(1001);
603  fTop->SetBorderSize(0);
604  fTop->SetBorderMode(0);
605 
606  fBody->Clear();
607  fBody->SetNumber(2);
608  fBody->SetFillColor(0);
609  fBody->SetFillStyle(0);
610  fBody->SetBorderSize(0);
611  fBody->SetBorderMode(0);
612  fBody->SetTopMargin(0.01);
613  fBody->SetRightMargin(0.01);
614  fBody->SetBottomMargin(0.10);
615  fBody->SetLeftMargin(0.10);
616  fBody->SetTicks();
617 
618  fCanvas->cd();
619 }
620 //____________________________________________________________________
622 {
623  gStyle->SetOptStat(0);
624  gStyle->SetOptTitle(0);
625 
626  Int_t h = 1000;
627  Int_t w = h / TMath::Sqrt(2);
628  if (fLandscape) {
629  Int_t t = h;
630  h = w;
631  w = t;
632  }
633 
634  fCanvas = new TCanvas("c",outputName,w,h);
635  fCanvas->SetFillColor(0);
636  fCanvas->SetBorderSize(0);
637  fCanvas->SetBorderMode(0);
638 
639  if (fPDF) {
640  SuppressGuard g;
641  fCanvas->Print(Form("%s[", outputName.Data()),
642  Form("pdf %s", fLandscape ? "Landscape" : ""));
643  }
644  fCanvas->SetLeftMargin (0.10);
645  fCanvas->SetRightMargin (0.05);
646  fCanvas->SetTopMargin (0.05);
647  fCanvas->SetBottomMargin(0.10);
648 
649  Float_t dy = 0.05;
650  fTop = new TPad("top","Top",0,1-dy,1,1,0,0);
651  // fTop->SetNumber(1);
652  // fTop->SetFillColor(kTopBackground);
653  // fTop->SetBorderSize(0);
654  // fTop->SetBorderMode(0);
655  fCanvas->cd();
656  fTop->Draw();
657 
658  fBody = new TPad("body","Body",0,0,1,1-dy,0,0);
659  fBody->SetNumber(2);
660  // fBody->SetFillColor(0);
661  // fBody->SetFillStyle(0);
662  // fBody->SetBorderSize(0);
663  // fBody->SetBorderMode(0);
664  fCanvas->cd();
665  fBody->Draw();
666 
667  ClearCanvas();
668 
669  fHeader = new TLatex(.5, .5, "Title");
670  fHeader->SetNDC();
671  fHeader->SetTextAlign(22);
672  // fHeader->SetTextColor(kWhite);
673  fHeader->SetTextFont(62);
674  fHeader->SetTextSize(0.7);
675 
676  fCanvas->cd();
677 }
678 //____________________________________________________________________
680 {
681  if (fPDF && fCanvas) {
682  SuppressGuard g;
683  fCanvas->Print(Form("%s]", fCanvas->GetTitle()),
684  Form("pdf %s Title:%s",
685  fLandscape ? "Landscape" : "",
686  fLastTitle.Data()));
687  Printf("PDF %s written", fCanvas->GetTitle());
688  }
689  if (fCanvas) fCanvas->Close();
690  fCanvas = 0;
691 }
692 
693 //____________________________________________________________________
695 {
696  if (fTop) {
697  fTop->cd();
698  fHeader->SetTextSize(size);
699  fHeader->DrawLatex(.5,.5,title);
700  }
701  fCanvas->Modified();
702  fCanvas->Update();
703  fCanvas->cd();
704 
705  if (fPDF) {
706  TString tit;
707  tit.Form("pdf %s Title:%s",
708  fLandscape ? "Landscape" : "",
709  title.Data());
710  // Suppress prints
711  SuppressGuard g;
712  fCanvas->Print(fCanvas->GetTitle(), tit);
713  }
714  fLastTitle = title;
715  if (fPause) fCanvas->WaitPrimitive();
716 }
717 //====================================================================
719 {
720  TFile* file = TFile::Open(filename, "READ");
721  if (!file) {
722  Warning("OpenFile", "Failed to open \"%s\"", filename);
723  return 0;
724  }
725  return file;
726 }
727 //____________________________________________________________________
729  Color_t color,
730  Style_t marker,
731  Double_t size,
732  Style_t fill,
733  Style_t line,
734  Width_t width)
735 {
736  h->SetMarkerColor(color);
737  h->SetMarkerStyle(marker);
738  h->SetMarkerSize(size);
739  h->SetFillColor(color);
740  h->SetFillStyle(fill);
741  h->SetLineColor(color);
742  h->SetLineStyle(line);
743  h->SetLineWidth(width);
744  h->GetXaxis()->SetNdivisions(210);
745  h->GetYaxis()->SetNdivisions(210);
746 }
747 
748 //____________________________________________________________________
749 TLegend* AliTrackletdNdetaPost::DrawInPad(TVirtualPad* c,
750  Int_t pad,
751  TObject* o,
752  Option_t* opt)
753 {
754  if (!o) {
755  Warning("", "Nothing to draw in pad %d", pad);
756  return 0;
757  }
758  TLegend* l = 0;
759  TVirtualPad* p = c->cd(pad);
760  TString option(opt);
761  option.ToLower();
762  if (option.Contains("logx")) { p->SetLogx(); option.ReplaceAll("logx",""); }
763  if (option.Contains("logy")) { p->SetLogy(); option.ReplaceAll("logy",""); }
764  if (option.Contains("logz")) { p->SetLogz(); option.ReplaceAll("logz",""); }
765  if (option.Contains("grid")) { p->SetGridx(); p->SetGridy();
766  option.ReplaceAll("grid",""); }
767  Int_t leg = 0;
768  if (option.Contains("leg2")) { leg = 2; option.ReplaceAll("leg2",""); }
769  if (option.Contains("leg")) { leg = 1; option.ReplaceAll("leg",""); }
770  o->Draw(opt);
771  if (leg) {
772  l = p->BuildLegend(0.5, 0.73, .98, .98);
773  l->SetNColumns(leg);
774  TObject* frame = 0;
775  TIter next(l->GetListOfPrimitives());
776  TLegendEntry* ent = 0;
777  while ((ent = static_cast<TLegendEntry*>(next()))) {
778  if (TString(ent->GetLabel()).EqualTo("frame")) frame = ent;
779  }
780  if (frame) l->GetListOfPrimitives()->Remove(frame);
781  // l->GetListOfPrimitives()->Print();
782  }
783  p->Modified();
784  // p->Update();
785  // p->cd();
786  return l;
787 }
788 //====================================================================
789 void AliTrackletdNdetaPost::DrawParam(const char* name,
790  Double_t& y,
791  const char* val)
792 {
793  TLatex* ln = new TLatex(.49, y, name);
794  ln->SetTextAlign(31);
795  ln->SetTextSize(0.02/gPad->GetHNDC());
796  ln->SetNDC();
797  ln->Draw();
798  TLatex* lv = new TLatex(.51, y, val);
799  lv->SetTextAlign(11);
800  lv->SetTextSize(0.02/gPad->GetHNDC());
801  lv->SetNDC();
802  lv->Draw();
803  y -= 0.025/gPad->GetHNDC();
804 }
805 //____________________________________________________________________
806 void AliTrackletdNdetaPost::DrawParam(const char* name,
807  Double_t& y,
808  Double_t val)
809 {
810  DrawParam(name, y, Form("%f", val));
811 }
812 //____________________________________________________________________
813 void AliTrackletdNdetaPost::DrawParam(const char* name,
814  Double_t& y,
815  Int_t val)
816 {
817  DrawParam(name, y, Form("%d", val));
818 }
819 //____________________________________________________________________
820 void AliTrackletdNdetaPost::DrawParam(const char* name,
821  Double_t& y,
822  Bool_t val)
823 {
824  DrawParam(name, y, val ? "yes" : "no");
825 }
826 //____________________________________________________________________
828  const char* title,
829  Bool_t comb)
830 {
831  // c->Clear();
832  Double_t y = .9;
833  TLatex* latex = new TLatex(.5, y, title);
834  latex->SetTextAlign(21);
835  latex->SetTextSize(0.023/gPad->GetHNDC());
836  latex->SetNDC();
837  latex->Draw();
838  y -= 0.028/gPad->GetHNDC();
839  if (!pars) return;
840  Int_t mode = GetI(pars, "RecMode");
841  DrawParam("Reconstruction mode",y,
842  Form("%s%s%s",
843  (mode & 0x2 ? "recon" : ""),
844  (mode & 0x4 ? " inj" : ""),
845  (mode & 0x8 ? " rot" : "")));
846  DrawParam("Scale by sin^{2}(#it{#theta})",y,
847  bool(GetB(pars, "ScaleDTheta")));
848  DrawParam("#delta#phi shift", y, GetD(pars, "DPhiShift"));
849  DrawParam("Shifted #delta#phi cut",
850  y, GetD(pars, "ShiftedDPhiCut"));
851  Double_t tmp = GetD(pars, "ScaledDThetaCut");
852  DrawParam("Scaled #delta#theta cut",
853  y, (tmp < 0 ? "n/a" : Form("%f", tmp)));
854  DrawParam("#Delta cut", y, GetD(pars, "DeltaCut"));
855  DrawParam("max#Delta", y, GetD(pars, "MaxDelta"));
856  DrawParam("min#Delta_{tail}", y,
857  GetD(pars,"TailDelta"));
858  DrawParam("d#theta window", y, GetD(pars, "DThetaWindow"));
859  DrawParam("d#phi window", y, GetD(pars, "DPhiWindow"));
860  DrawParam("#phi overlap cut", y, GetD(pars, "PhiOverlapCut"));
861  DrawParam("#it{z}-#eta overlap cut",y, GetD(pars, "ZEtaOverlapCut"));
862  DrawParam("#phi rotation", y, TMath::RadToDeg()*
863  GetD(pars, "PhiRotation"));
864  y -= 0.045;
865  DrawParam("Background method", y, comb ? "Combinatorics" : "Injection");
866 }
867 //____________________________________________________________________
869  Container* simSums, Bool_t simComb)
870 {
871  ClearCanvas();
872  fBody->Divide(1,3,0,0);
873 
874  Double_t yr = .1;
875  TVirtualPad* p1 = fBody->GetPad(1); p1->SetPad(0,1-yr,1,1);
876  TVirtualPad* p2 = fBody->GetPad(2); p2->SetPad(0,(1-yr)/2,1,1-yr);
877  TVirtualPad* p3 = fBody->GetPad(3); p3->SetPad(0,0,1,(1-yr)/2);
878 
879  // Post-processing stuff
880  fBody->cd(1);
881  Double_t y = .80;
882  TLatex* latex = new TLatex(.5, y, "Post-processing");
883  latex->SetTextAlign(21);
884  latex->SetTextSize(0.023/p1->GetHNDC());
885  latex->SetNDC();
886  latex->Draw();
887  y -= 0.028/p1->GetHNDC();
888  DrawParam("Scaling of comb. bg.", y,
890  Form("%f", fCombinatoricsScale) :
892  "per centrality" :
893  "per eta,centrality"));
894  DrawParam("min#alpha", y, fAlphaMin);
895  DrawParam("max#alpha", y, fAlphaMax);
896 
897  // From tasks
898  fBody->cd(2);
899  DrawParams(GetC(realSums, "parameters"), "Real data", realComb);
900  fBody->cd(3);
901  DrawParams(GetC(simSums, "parameters"), "Simulated data", simComb);
902  PrintCanvas("Parameters");
903 }
904 //====================================================================
905 THStack* AliTrackletdNdetaPost::Make2Stack(const char* name,
906  const char* title,
907  Container* realList,
908  Container* simList,
909  Option_t* dataOpt,
910  Option_t* simOpt)
911 {
912  TString nme(name);
913  THStack* stack = new THStack(name, title);
914  TH1* data = CopyH1(realList, name, Form("real%s",name));
915  TH1* sim = CopyH1(simList, name, Form("sim%s",name));
916  TString dtit(data->GetTitle());
917  if (dtit.Contains("\\")) dtit.Form("%s\\hbox{ - real}", data->GetTitle());
918  else dtit.Form("%s - real", data->GetTitle());
919  data->SetTitle(dtit);
920  TString stit(sim->GetTitle());
921  if (stit.Contains("\\")) stit.Form("%s\\hbox{ - sim.}", sim->GetTitle());
922  else stit.Form("%s - sim.", sim->GetTitle());
923  sim->SetTitle(stit);
924  stack->Add(data, dataOpt);
925  stack->Add(sim, simOpt);
926  return stack;
927 }
928 //____________________________________________________________________
930  const char* name,
931  Container* realList,
932  Container* simList)
933 {
934  fBody->cd(pad);
935  TH2* real = GetH2(realList, name);
936  TH2* sim = GetH2(simList, name);
937  real->Draw("colz");
938  sim->Draw("box same");
939  fBody->GetPad(pad)->Modified();
940 }
941 
942 //____________________________________________________________________
943 void ModLegend(TVirtualPad* p, TLegend* l,
944  Double_t x1, Double_t y1,
945  Double_t x2, Double_t y2)
946 {
947  Double_t px1 = p->GetX1();
948  Double_t px2 = p->GetX2();
949  Double_t py1 = p->GetY1();
950  Double_t py2 = p->GetY2();
951  l->SetX1(px1+(px2-px1)*x1);
952  l->SetX2(px1+(px2-px1)*x2);
953  l->SetY1(py1+(py2-py1)*y1);
954  l->SetY2(py1+(py2-py1)*y2);
955  p->Modified();
956 }
957 //____________________________________________________________________
959  Container* simList)
960 {
961  THStack* ipz = Make2Stack("ipz", "IP_{#it{z}}", realList,simList);
962  THStack* cent = Make2Stack("cent", "Centrality [%]", realList,simList);
963  THStack* status = Make2Stack("status","Task status", realList,simList,
964  "B text90", "B text90");
965  TH1* mcStatus = GetH1(simList, "statusMC");
966  mcStatus->SetMaximum(1.2*mcStatus->GetMaximum());
967  ClearCanvas();
968  fBody->Divide(2,4);
969  for (Int_t i = 1; i <= 8; i++) {
970  fBody->GetPad(i)->SetRightMargin(0.01);
971  fBody->GetPad(i)->SetTopMargin(0.01);
972  }
973  TLegend* l = 0;
974  l = DrawInPad(fBody,1,ipz, "nostack leg");
975  ModLegend(fBody->GetPad(1),l,.4,.1,.75,.4);
976  l = DrawInPad(fBody,2,cent, "nostack leg");
977  ModLegend(fBody->GetPad(2),l,.6,.1,.99,.4);
978  l = DrawInPad(fBody,3,status, "nostack hist text90 leg");
979  ModLegend(fBody->GetPad(3),l,.5,.7,.99,.99);
980  DrawInPad(fBody,4,mcStatus,"hist text90");
981  DrawClusters(5, "allClusters0", realList, simList);
982  DrawClusters(6, "allClusters1", realList, simList);
983  DrawClusters(7, "usedClusters0", realList, simList);
984  DrawClusters(8, "usedClusters1", realList, simList);
985  PrintCanvas("General information");
986 }
987 //====================================================================
989  const char* sub,
990  const char* name)
991 {
992  if (!ress) return 0;
993  TObjArray* tokens = TString(sub).Tokenize("/");
994  Container* subCont = ress;
995  TIter next(tokens);
996  TObjString* dir = 0;
997  while ((dir = static_cast<TObjString*>(next()))) {
998  subCont = GetC(ress, dir->GetName());
999  if (!subCont) break;
1000  }
1001  TH1* ret = 0;
1002  if (subCont) ret = GetH1(subCont, name);
1003  tokens->Delete();
1004  return ret;
1005 }
1006 //____________________________________________________________________
1008  const char* sub,
1009  Bool_t scaled)
1010 {
1011  if (!ress) return 0;
1012  Container* subCont = GetC(ress, sub);
1013  if (!subCont) return 0;
1014  TH1* h = GetH1(subCont, Form("delta%s", scaled ? "Scaled" : ""));
1015  if (!h) return 0;
1016  if (!scaled) return h;
1017  TString tit(h->GetTitle());
1018  if (tit.Contains("\\")) tit.Append(Form("\\hbox{ - %s}",sub));
1019  else tit.Append(Form(" - %s", sub));
1020  return h;
1021 }
1022 //____________________________________________________________________
1024 {
1025  THStack* particlePdgs = new THStack("particlePdg","Particle species");
1026  particlePdgs->Add(GetH1(mcList, "allPrimaryPdg"));
1027  particlePdgs->Add(GetH1(mcList, "cutPrimaryPdg"));
1028  particlePdgs->Add(GetH1(mcList, "allSecondaryPdg"));
1029  particlePdgs->Add(GetH1(mcList, "cutSecondaryPdg"));
1030  particlePdgs->SetMaximum(particlePdgs->GetMaximum("nostack")*50);
1031  particlePdgs->SetMinimum(8*1e-7);
1032  THStack* parentPdgs = new THStack("parentPdg","Parent species");
1033  parentPdgs->Add(GetH1(mcList, "allPrimaryParentPdg"));
1034  parentPdgs->Add(GetH1(mcList, "cutPrimaryParentPdg"));
1035  parentPdgs->Add(GetH1(mcList, "allSecondaryParentPdg"));
1036  parentPdgs->Add(GetH1(mcList, "cutSecondaryParentPdg"));
1037  parentPdgs->SetMaximum(parentPdgs->GetMaximum("nostack")*50);
1038  parentPdgs->SetMinimum(8*1e-5);
1039 
1040  ClearCanvas();
1041  fBody->Divide(1,2,0,0);
1042  fBody->GetPad(1)->SetRightMargin(0.01);
1043  fBody->GetPad(2)->SetRightMargin(0.01);
1044 
1045  TLegend* p1Leg = DrawInPad(fBody, 1, particlePdgs, "nostack logy leg2 grid");
1046  TLegend* p2Leg = DrawInPad(fBody, 2, parentPdgs, "nostack logy leg2 grid");
1047  p1Leg->SetHeader(particlePdgs->GetTitle());
1048  p2Leg->SetHeader(parentPdgs ->GetTitle());
1049  ModLegend(fBody->GetPad(1),p1Leg,.2,.8,.99,1);
1050  ModLegend(fBody->GetPad(2),p2Leg,.2,.8,.99,1);
1051 
1052  particlePdgs->GetHistogram()->GetXaxis()->LabelsOption("v");
1053  parentPdgs ->GetHistogram()->GetXaxis()->LabelsOption("v");
1054  // particlePdgs->SetMaximum(1.4*particlePdgs->GetMaximum("nostack"));
1055  // parentPdgs ->SetMaximum(1.4*parentPdgs ->GetMaximum("nostack"));
1056 
1057  fBody->GetPad(1)->Modified();
1058  fBody->GetPad(2)->Modified();
1059 
1060  PrintCanvas(Form("%s - species", fLastBin.Data()));
1061 }
1062 //____________________________________________________________________
1064 {
1065 
1066  ClearCanvas();
1067  fBody->Divide(1,2,0,0);
1068  fBody->GetPad(1)->SetRightMargin(0.01);
1069  fBody->GetPad(2)->SetRightMargin(0.01);
1070 
1071 
1072  THStack* orig = new THStack("orig", "#Delta");
1073  orig->Add(FindDelta(realList, "data", false));
1074  orig->Add(FindDelta(realList, "injection", false));
1075  orig->Add(FindDelta(simList, "data", false));
1076  orig->Add(FindDelta(simList, "injection", false));
1077  orig->Add(FindDelta(simList, "primaries", false));
1078  orig->Add(FindDelta(simList, "combinatorics", false));
1079  orig->Add(FindDelta(simList, "secondaries", false));
1080  orig->Add(FindDelta(simList, "uncorrelatedCombinatorics",false));
1081 
1082  THStack* scaled = new THStack("scaled", "#Delta (scaled)");
1083  scaled->Add(FindDelta(realList, "data", false));
1084  scaled->Add(FindDelta(realList, "injection", true));
1085  scaled->Add(FindDelta(simList, "data", false));
1086  scaled->Add(FindDelta(simList, "injection", true));
1087  scaled->Add(FindDelta(simList, "primaries", false));
1088  scaled->Add(FindDelta(simList, "combinatorics", true));
1089  scaled->Add(FindDelta(simList, "secondaries", false));
1090  scaled->Add(FindDelta(simList, "uncorrelatedCombinatorics",true));
1091 
1092  Double_t max = /*5*/1.2*orig->GetMaximum("nostack");
1093  TH1* templ = FindDelta(realList, "data", false);
1094  TH2* frame = new TH2D("frame","",
1095  templ->GetNbinsX(),
1096  templ->GetXaxis()->GetXmin(),
1097  templ->GetXaxis()->GetXmax(),
1098  500, 0.001, max);
1099  frame->SetXTitle(templ->GetXaxis()->GetTitle());
1100  frame->SetYTitle(templ->GetYaxis()->GetTitle());
1101  frame->SetStats(0);
1102  frame->SetDirectory(0);
1103  TLine* cuts = new TLine(fDeltaCut, 0.001, fDeltaCut, .99*max);
1104  cuts->SetLineColor(kGreen+1);
1105  frame->GetListOfFunctions()->Add(cuts);
1106  TLine* cutt = new TLine(fTailDelta, 0.001, fTailDelta, .99*max);
1107  cutt->SetLineColor(kRed+1);
1108  frame->GetListOfFunctions()->Add(cutt);
1109 
1110  TLegend* l = 0;
1111  TH1* f1 = static_cast<TH1*>(frame->Clone());
1112  f1->SetDirectory(0);
1113  DrawInPad(fBody, 1, f1, "axis logx logy grid");
1114  l = DrawInPad(fBody, 1, orig, "nostack same leg2");
1115  // ModLegend(fBody->GetPad(1), l, .2, .8, .99, .99);
1116  ModLegend(fBody->GetPad(1), l,
1117  fBody->GetPad(1)->GetLeftMargin(),
1118  fBody->GetPad(1)->GetBottomMargin(),
1119  .6, .35);
1120 
1121 
1122  TH1* f2 = static_cast<TH1*>(frame->Clone());
1123  f2->SetDirectory(0);
1124  DrawInPad(fBody, 2, f2, "axis logx logy grid");
1125  l = DrawInPad(fBody, 2, scaled, "nostack same leg2");
1126  // ModLegend(fBody->GetPad(1), l, .2, .8, .99, 1);
1127  ModLegend(fBody->GetPad(2), l,
1128  fBody->GetPad(2)->GetLeftMargin(),
1129  fBody->GetPad(2)->GetBottomMargin(),
1130  .6, .45);
1131 
1132  PrintCanvas(Form("%s - #Delta", fLastBin.Data()));
1133 }
1134 //____________________________________________________________________
1136  Int_t pad,
1137  Container* ress,
1138  const char* name,
1139  const char* pretty)
1140 {
1141  SuppressGuard g(kError);
1142  Container* sc = GetC(ress, name);
1143  TH2* mask = GetH2(sc, "fiducial");
1144  TH2* data = CopyH2(sc, "etaVsIPz", "raw");
1145  TH2* bgEst = CopyH2(sc, "backgroundEst", "rawBg");
1146  TH2* mBeta = CopyH2(sc, "beta"); // oneMinusBeta
1147  TH2* sigEst= CopyH2(sc, "signalEst", "raw");
1148  if (data && mask) data ->Multiply(mask);
1149  if (bgEst && mask) bgEst ->Multiply(mask);
1150  if (mBeta && mask) mBeta ->Multiply(mask);
1151  if (sigEst && mask) sigEst->Multiply(mask);
1152  mBeta->SetMinimum(0);
1153  mBeta->SetMaximum(.7);
1154  Double_t max = TMath::Max(data->GetMaximum(), sigEst->GetMaximum());
1155  data ->SetMinimum(0);
1156  sigEst->SetMinimum(0);
1157  bgEst ->SetMinimum(0);
1158  data ->SetMaximum(max);
1159  sigEst->SetMaximum(max);
1160 
1161  TVirtualPad* q = c->cd(pad);
1162  q->SetRightMargin(0.10);
1163  q->Divide(4,1,0,0);
1164  DrawInPad(q,1,data, "colz");
1165  DrawInPad(q,2,bgEst, "colz");
1166  DrawInPad(q,3,mBeta, "colz");
1167  DrawInPad(q,4,sigEst, "colz");
1168  TVirtualPad* r = q->GetPad(4);
1169  r->SetRightMargin(0.12);
1170  r->SetPad(r->GetXlowNDC(), r->GetYlowNDC(),
1171  .99, r->GetYlowNDC()+r->GetHNDC());
1172 
1173  TLatex* ltx = new TLatex(.01,.5, pretty);
1174  ltx->SetNDC();
1175  ltx->SetTextFont(62);
1176  ltx->SetTextAngle(90);
1177  ltx->SetTextAlign(23);
1178  ltx->SetTextSize(0.07);
1179  DrawInPad(q,1,ltx,"");
1180 
1181  return max;
1182 }
1183 //____________________________________________________________________
1185  Container* simList)
1186 {
1187  ClearCanvas();
1188 
1189  fBody->SetLeftMargin(0.2);
1190  fBody->SetTopMargin(0.10);
1191  fBody->Divide(1,4,0,0);
1192  Double_t max = 0;
1193  max = TMath::Max(max, DrawBackground(fBody, 1, realList,
1194  "injection", "Injection - Real"));
1195  max = TMath::Max(max, DrawBackground(fBody, 2, simList,
1196  "injection", "Injection - Sim."));
1197  max = TMath::Max(max, DrawBackground(fBody, 3, simList,
1198  "combinatorics", "MC Labels"));
1199  max = TMath::Max(max, DrawBackground(fBody, 4, simList,
1200  "uncorrelatedCombinatorics",
1201  "MC Labels - Uncorr."));
1202  Int_t i = 0;
1203  for (i = 1; i <= 4; i++) {
1204  TVirtualPad* p = fBody->GetPad(i);
1205  if (!p) continue;
1206  for (Int_t j = 1; j <= 4; j++) {
1207  TVirtualPad* q = p->GetPad(j);
1208  TH1* h = static_cast<TH1*>(q->FindObject("raw"));
1209  if (h) {
1210  h->SetMaximum(max);
1211  continue;
1212  }
1213  h = static_cast<TH1*>(q->FindObject("rawBg"));
1214  if (h) h->SetMaximum(max/5);
1215  }
1216  }
1217  const char* headings[] = { "Measured",
1218  "Background est.",
1219  "#beta",
1220  "Signal", 0 };
1221  const char** ptr = headings;
1222  i = 0;
1223  while (*ptr) {
1224  TLatex* head = new TLatex(fBody->GetRightMargin()+
1225  (1-fBody->GetRightMargin())/4*(i+.5),
1226  .99, *ptr);
1227  head->SetNDC();
1228  head->SetTextAlign(23);
1229  head->SetTextFont(62);
1230  head->SetTextSize(0.025);
1231  DrawInPad(fBody, 0, head, "");
1232  ptr++;
1233  i++;
1234  }
1235  PrintCanvas(Form("%s - backgrounds", fLastBin.Data()));
1236 }
1237 //____________________________________________________________________
1238 void PrintH(TH2* h, Int_t prec=2)
1239 {
1240  Printf("Content of %s - %s", h->GetName(), h->GetTitle());
1241  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
1242  printf("%3d: ", i);
1243  for (Int_t j = 1; j <= h->GetNbinsY(); j++) {
1244  Double_t c = h->GetBinContent(i,j);
1245  Double_t e = h->GetBinError (i,j);
1246  if (TMath::IsNaN(c) || TMath::IsNaN(e))
1247  printf("*** NAN ***");
1248  else
1249  printf("%.*f+/-%.*f ", prec, c, prec, e);
1250  }
1251  printf("\n");
1252  }
1253 }
1254 //____________________________________________________________________
1256 {
1257  TH2* r = static_cast<TH2*>(h->Clone("fiducial"));
1258  r->SetDirectory(0);
1259  r->Reset();
1260  r->SetTitle("Fiducial cut");
1261  for (Int_t i = 1; i <= r->GetNbinsX(); i++) {
1262  for (Int_t j = 1; j <= r->GetNbinsY(); j++) {
1263  Double_t c = h->GetBinContent(i,j);
1264  if (c <= fAlphaMin || c > fAlphaMax) continue;
1265  r->SetBinContent(i,j,1);
1266  r->SetBinError (i,j,0);
1267  }
1268  }
1269  // PrintH(h, 1);
1270  // PrintH(r, 0);
1271  return r;
1272 }
1273 
1274 //____________________________________________________________________
1276  Int_t pad,
1277  Container* ress,
1278  const char* name,
1279  const char* pretty)
1280 {
1281  Container* sc = GetC(ress, name);
1282  TH2* mask = GetH2(sc, "fiducial");
1283  TH2* al = CopyH2(sc, "alpha");
1284  TH2* als = CopyH2(sc, "alphaSel");
1285  if (al && mask) al ->Multiply(mask);
1286  if (als && mask) als->Multiply(mask);
1287  al ->SetMinimum(fAlphaMin); al ->SetMaximum(fAlphaMax);
1288  als->SetMinimum(fAlphaMin); als->SetMaximum(fAlphaMax);
1289 
1290  TVirtualPad* q = c->cd(pad);
1291  q->SetRightMargin(0.10);
1292  q->Divide(2,1,0,0);
1293 
1294  DrawInPad(q,1,al, "colz");
1295  DrawInPad(q,2,als, "colz");
1296  TVirtualPad* r = q->GetPad(2);
1297  r->SetRightMargin(0.12);
1298  r->SetPad(r->GetXlowNDC(), r->GetYlowNDC(),
1299  .99, r->GetYlowNDC()+r->GetHNDC());
1300 
1301  TLatex* ltx = new TLatex(.01,.5, pretty);
1302  ltx->SetNDC();
1303  ltx->SetTextFont(62);
1304  ltx->SetTextAngle(90);
1305  ltx->SetTextAlign(23);
1306  ltx->SetTextSize(0.07);
1307  DrawInPad(q,1,ltx,"");
1308 }
1309 //____________________________________________________________________
1311 {
1312  ClearCanvas();
1313 
1314  fBody->SetLeftMargin(0.2);
1315  fBody->SetTopMargin(0.10);
1316  fBody->Divide(1,3,0,0);
1317  DrawAlpha(fBody, 1, simList, "injection", "Injection - Sim.");
1318  DrawAlpha(fBody, 2, simList, "combinatorics", "MC Labels");
1319  DrawAlpha(fBody, 3, simList, "uncorrelatedCombinatorics",
1320  "MC Labels - Uncorr.");
1321 
1322  const char* headings[] = { "#alpha all MC",
1323  "#alpha selected MC",
1324  0 };
1325  const char** ptr = headings;
1326  Int_t i = 0;
1327  while (*ptr) {
1328  TLatex* head = new TLatex(fBody->GetRightMargin()+
1329  (1-fBody->GetRightMargin())/2*(i+.5),
1330  .99, *ptr);
1331  head->SetNDC();
1332  head->SetTextAlign(23);
1333  head->SetTextFont(62);
1334  head->SetTextSize(0.025);
1335  DrawInPad(fBody, 0, head, "");
1336  ptr++;
1337  i++;
1338  }
1339  PrintCanvas(Form("%s - #alpha", fLastBin.Data()));
1340 }
1341 //____________________________________________________________________
1343 {
1344  ClearCanvas();
1345 
1346  fBody->Divide(1,2);
1347  THStack* dNdeta = new THStack("dndeta",ObsTitle());
1348  TH1* all = GetH1(simList, "truedNdeta");
1349  TH1* sel = GetH1(simList, "truedNdetaSel");
1350  dNdeta->Add(all);
1351  dNdeta->Add(sel);
1352  fBody->GetPad(1)->SetTopMargin(0.01);
1353  fBody->GetPad(1)->SetRightMargin(0.01);
1354  TLegend* l = DrawInPad(fBody, 1, dNdeta, "nostack grid leg");
1355  ModLegend(fBody->GetPad(1), l, .4, .1, .99, .4);
1356  // dNdeta->SetMaximum(1.5*dNdeta->GetMaximum("nostack"));
1357 
1358  TVirtualPad* p = fBody->cd(2);
1359  p->SetRightMargin(0.0);
1360  p->SetBottomMargin(0.15);
1361  p->Divide(1,3,0,0);
1362  p->GetPad(1)->SetRightMargin(0.01);
1363  p->GetPad(2)->SetRightMargin(0.01);
1364  p->GetPad(3)->SetRightMargin(0.01);
1365  THStack* ipz = new THStack("ipz", "IP_{#it{z}}");
1366  ipz->Add(GetH1(simList, "ipz"));
1367  ipz->Add(GetH1(simList, "ipzGen"));
1368  ipz->Add(GetH1(simList, "ipzSel"));
1369  l = DrawInPad(p, 1, ipz, "nostack grid leg");
1370  ipz->GetHistogram()->SetYTitle("Fraction");
1371  ipz->GetHistogram()->GetYaxis()->SetTitleSize(0.12);
1372  ipz->GetHistogram()->GetYaxis()->SetTitleOffset(0.4);
1373  ipz->GetHistogram()->GetYaxis()->SetLabelSize(0.09);
1374  ModLegend(p->GetPad(1), l, .4, 0, .75, .4);
1375  TGraphAsymmErrors* ipzEff =
1376  static_cast<TGraphAsymmErrors*>(GetO(simList, "ipzEff",
1377  TGraphAsymmErrors::Class()));
1378  DrawInPad(p, 2, ipzEff, "ap5 grid");
1379  ipzEff->GetHistogram()->SetYTitle("#epsilon_{IP_{#it{z}}}");
1380  ipzEff->GetHistogram()->GetYaxis()->SetTitleSize(0.12);
1381  ipzEff->GetHistogram()->GetYaxis()->SetLabelSize(0.09);
1382  ipzEff->GetHistogram()->GetYaxis()->SetTitleOffset(0.4);
1383  ipzEff->GetHistogram()->GetXaxis()->SetRangeUser(ipz->GetXaxis()->GetXmin(),
1384  ipz->GetXaxis()->GetXmax());
1385  p->GetPad(2)->Modified();
1386  TH1* ipzRes = 0;
1387  DrawInPad(p, 3, ipzRes = GetH1(simList, "ipzVsGenIPz"), "grid");
1388  ipzRes->GetYaxis()->SetTitleSize(0.12);
1389  ipzRes->GetYaxis()->SetTitleOffset(0.4);
1390  ipzRes->GetYaxis()->SetLabelSize(0.09);
1391  ipzRes->GetXaxis()->SetTitleSize(0.09);
1392  ipzRes->GetXaxis()->SetLabelSize(0.09);
1393  ipzRes->GetXaxis()->SetTitleOffset(0.4);
1394 
1395  PrintCanvas(Form("%s - MC information",fLastBin.Data()));
1396 }
1397 //____________________________________________________________________
1399 {
1400  Container* data = GetC(c, "data");
1401  Container* inj = GetC(c, "injection");
1402  // Get integral and error from data
1403  Double_t dataIntg = GetD(data, "deltaTailIntg");
1404  Double_t dataInte = GetD(data, "deltaTailIntE");
1405  // Get integral and error from injection
1406  Double_t injIntg = GetD(inj, "deltaTailIntg");
1407  Double_t injInte = GetD(inj, "deltaTailIntE");
1408  // Return the ratio of data to injection
1409  Double_t r = RatioE(dataIntg, dataInte, injIntg, injInte, e);
1410  // Printf("%20s: data=%9g+/-%9g inj=%9g+/-%9g -> %9g+/-%9g",
1411  // c->GetName(), dataIntg, dataInte, injIntg, injInte, r, e);
1412  return r;
1413 }
1414 
1415 //____________________________________________________________________
1417 {
1418  TH1* deltaIntg = 0;
1419  {
1420  // SuppressGuard g(kError);
1421  deltaIntg = GetH1(c, "deltaTailInt");
1422  }
1423  if (deltaIntg) return deltaIntg;
1424 
1425  TH2* etaVsDelta = GetH2(c, "etaVsDelta");
1426  Double_t maxDelta = etaVsDelta->GetYaxis()->GetXmax();
1427  deltaIntg = etaVsDelta->ProjectionX("deltaTailInt");
1428  deltaIntg->SetDirectory(0);
1429  deltaIntg->Reset();
1430  deltaIntg->SetYTitle(Form("#int_{%3.1f}^{%4.1f}d#Delta",
1431  fTailDelta, maxDelta));
1432  for (Int_t i = 1; i <= deltaIntg->GetNbinsX(); i++) {
1433  TH1* tmp = etaVsDelta->ProjectionY("tmp",i,i);
1434  Double_t eintg;
1435  Double_t intg = Integrate(tmp, fTailDelta, maxDelta, eintg);
1436  deltaIntg->SetBinContent(i, intg);
1437  deltaIntg->SetBinError (i, eintg);
1438  delete tmp;
1439  }
1440  return deltaIntg;
1441 }
1442 //____________________________________________________________________
1444 {
1445  Container* data = GetC(c, "data");
1446  Container* inj = GetC(c, "injection");
1447  TH1* dataInt = GetDeltaIntDist(data);
1448  TH1* injInt = GetDeltaIntDist(inj);
1449  TH1* ret = static_cast<TH1*>(dataInt->Clone("deltaTailRatio"));
1450  ret->SetDirectory(0);
1451  ret->SetYTitle("Ratio of tails (data/injection)");
1452  ret->Divide(injInt);
1453  return ret;
1454 }
1455 
1456 //____________________________________________________________________
1457 Style_t AliTrackletdNdetaPost::MS(Int_t what, Bool_t sim, Bool_t alt) const
1458 {
1459  if (!alt) {
1460  switch (what) {
1461  case 0: /* dndeta */ return (sim ? 24 : 20); // cirlce
1462  case 1: /* M */ return (sim ? 24 : 20); // circle
1463  case 2: /* B */ return (sim ? 25 : 21); // square
1464  case 3: /* S */ return (sim ? 26 : 22); // up-triangle
1465  default: return (sim ? 5 : 2); //
1466  }
1467  }
1468  switch (what) {
1469  case 0: /* dndeta */ return (sim ? 25 : 21); // square
1470  case 1: /* M */ return (sim ? 28 : 34); // cross
1471  case 2: /* B */ return (sim ? 27 : 33); // diamond
1472  case 3: /* S */ return (sim ? 30 : 29); // star
1473  }
1474  return (sim ? 5 : 2); //
1475 }
1476 //____________________________________________________________________
1478  Bool_t realComb,
1479  Container* simList,
1480  Bool_t simComb,
1481  Color_t color,
1482  THStack* stack,
1483  TDirectory* out,
1484  Bool_t alt)
1485 {
1486  Double_t k = (realList == simList ? 1 : fCombinatoricsScale);
1487  Container* realData = GetC(realList, "data");
1488  Container* simData = GetC(simList, "data");
1489  Container* simCombC = GetC(simList, "combinatorics");
1490  Container* simBgC = simComb ? simCombC : GetC(simList, "injection");
1491  Container* realBgC = realComb ? simCombC : GetC(realList, "injection");
1492 
1493 
1494  TH2* realMeas = CopyH2(realData, "etaVsIPz", "realMeas");
1495  TH1* realIPz = CopyH1(realList, "ipz", "realIPz");
1496  TH2* realBg = CopyH2(realBgC, "backgroundEst", "realBg");
1497  TH2* realSig = 0;
1498  Double_t scale = 0;
1499  Double_t scaleE = 0;
1500  TH1* scaleH = 0;
1501  if (realComb) {
1502  // Scale combinatorial background
1504  // fixed scalar
1505  Scale(realBg, fCombinatoricsScale, 0);
1506  scale = k;
1507  }
1508  else if (fCombinatoricsScaleMode == kPerCent) {
1509  // Printf("");
1510  // Get scalar and error from data
1511  Double_t eReal = 0;
1512  Double_t rReal = GetDeltaScale(realList, eReal);
1513  // Get integral and error from simulations
1514  Double_t eSim = 0;
1515  Double_t rSim = GetDeltaScale(simList, eSim);
1516  // Calculate the scalar of the combinatorial background
1517  scale = RatioE(rReal, eReal, rSim, eSim, scaleE);
1518  Scale(realBg, scale, scaleE);
1519  }
1520  else {
1521  TH1* realTail = GetDeltaScaleDist(realList);
1522  TH1* simTail = GetDeltaScaleDist(simList);
1523  scaleH = static_cast<TH1*>(realTail->Clone("scaleBg"));
1524  scaleH->SetDirectory(0);
1525  scaleH->SetYTitle("k_{#eta}");
1526  scaleH->SetTitle("Scalar");
1527  SetAttr(scaleH, scaleH->GetMarkerColor(), 34);
1528  scaleH->Divide(simTail);
1529  TF1* scaleF = new TF1("ff","pol0",
1530  scaleH->GetXaxis()->GetXmin(),
1531  scaleH->GetXaxis()->GetXmax());
1532  TFitResultPtr scaleR = scaleH->Fit(scaleF, "QN0SR");
1533  scale = scaleR->Parameter(0);
1534  scaleE = scaleR->ParError(0);
1535  Scale(realBg, scaleH);
1536  }
1537  // Get the background beta, and scale it
1538  TH2* beta = CopyH2(realBgC, "beta");
1539  Scale(beta, scale, scaleE);
1540  // Create unit histogram
1541  TH2* one = CopyH2(realBgC, "beta", "one");
1542  one->Reset();
1543  for (Int_t i = 1; i <= one->GetNbinsX(); i++)
1544  for (Int_t j = 1; j <= one->GetNbinsY(); j++)
1545  one->SetBinContent(i,j,1);
1546  // Subtract k times beta
1547  one->Add(beta,-1);
1548  // Multiply on to the measured distribution
1549  realSig = CopyH2(realData, "etaVsIPz", "realSig");
1550  realSig->SetMarkerStyle(beta->GetMarkerStyle());
1551  realSig->Multiply(one);
1552  // Clean-up temporary histograms
1553  delete beta;
1554  delete one;
1555  }
1556  else {
1557  // We can get the signal estimate directory from the background
1558  // container.
1559  realSig = CopyH2(realBgC, "signalEst");
1560  }
1561  if (!realMeas || !realSig || !realBg || !realIPz) {
1562  Warning("DrawdNdeta", "One or more real data histograms missing: "
1563  "meas=%p sig=%p bg=%p ipz=%p", realMeas, realSig, realBg, realIPz);
1564  return;
1565  }
1566  if (realBg ->GetMarkerStyle() == 30) realBg ->SetMarkerStyle(29);
1567  if (realBg ->GetMarkerStyle() == 27) realBg ->SetMarkerStyle(33);
1568  if (realSig->GetMarkerStyle() == 30) realSig->SetMarkerStyle(29);
1569  if (realSig->GetMarkerStyle() == 27) realSig->SetMarkerStyle(33);
1570 
1571  TH2* simMeas = CopyH2(simData, "etaVsIPz", "simMeas");
1572  TH1* simIPz = CopyH1(simList, "ipz", "simIPz");
1573  TH2* simBg = CopyH2(simBgC, "backgroundEst", "simBg");
1574  TH2* simSig = CopyH2(simBgC, "signalEst", "simSig");
1575  TH2* alpha = CopyH2(simBgC, "alpha");
1576  TH2* fiducial = CutAlpha(alpha);
1577  if (!simMeas || !simSig || !simBg || !simIPz || !alpha) {
1578  Warning("DrawdNdeta", "One or more simuluated data histograms missing: "
1579  "meas=%p sig=%p bg=%p ipz=%p alpha=%p",
1580  simMeas, simSig, simBg, simIPz, alpha);
1581  return;
1582  }
1583  alpha->Multiply(fiducial);
1584  TH2* trueGen = CopyH2(simList, "etaVsIPzMC", "trueGen");
1585  TH1* trueIPz = CopyH1(simList, "ipzGen", "trueIPz");
1586  if (!trueGen || !trueIPz) {
1587  Warning("DrawdNdeta", "One or more generator data histograms missing: "
1588  "gen=%p ipz=%p", trueGen, trueIPz);
1589  return;
1590  }
1591 
1592  realMeas->SetTitle("Measured (real)");
1593  realBg ->SetTitle("Background (real)");
1594  realSig ->SetTitle("Signal (real)");
1595  realIPz ->SetTitle("IP_{z} (real)");
1596  simMeas ->SetTitle("Measured (simulated)");
1597  simBg ->SetTitle("Background (simulated)");
1598  simSig ->SetTitle("Signal (simulated)");
1599  simIPz ->SetTitle("IP_{z} (simulated)");
1600  trueGen ->SetTitle("Generated");
1601  trueIPz ->SetTitle("IP_{z} (generated)");
1602 
1603  // Create result as alpha times real signal
1604  TH2* result = static_cast<TH2*>(realSig->Clone("result"));
1605  result->SetTitle("Result");
1606  result->SetDirectory(0);
1607  result->Multiply(alpha);
1608 
1609  UShort_t mode = 1;
1610  // Calculate the primary dN/deta
1611  TH1* truth = AverageOverIPz(trueGen, "truth", mode, trueIPz, 0);
1612  truth->SetYTitle(ObsTitle());
1613  SetAttr(truth, color, MS(0,true,alt), 1.5, 0, 7, 1);
1614 
1615  // Calculate the real dN/deta
1616  TH1* dndeta = AverageOverIPz(result, "dndeta", mode, realIPz, 0);
1617  dndeta->SetYTitle(ObsTitle());
1618  SetAttr(dndeta, color, MS(0,false,alt), 1.2, 0, 1, 1);
1619 
1620  // Multiply histograms by fiducial cut
1621  realMeas->Multiply(fiducial);
1622  simMeas ->Multiply(fiducial);
1623  realSig ->Multiply(fiducial);
1624  simSig ->Multiply(fiducial);
1625  realBg ->Multiply(fiducial);
1626  simBg ->Multiply(fiducial);
1627  realSig ->Multiply(fiducial);
1628  simSig ->Multiply(fiducial);
1629  // Get some step histograms. Note, we use the result histogram to
1630  // select which bins we project so that we get the used sub-sample.
1631  TH1* realAvgMeas = AverageOverIPz(realMeas,"realAvgMeas",mode,realIPz,result);
1632  TH1* realAvgSig = AverageOverIPz(realSig, "realAvgSig", mode,realIPz,result);
1633  TH1* realAvgBg = AverageOverIPz(realBg, "realAvgBg", mode,realIPz,result);
1634  TH1* simAvgMeas = AverageOverIPz(simMeas, "simAvgMeas", mode,simIPz, result);
1635  TH1* simAvgSig = AverageOverIPz(simSig, "simAvgSig", mode,simIPz, result);
1636  TH1* simAvgBg = AverageOverIPz(simBg, "simAvgBg", mode,simIPz, result);
1637  SetAttr(realAvgMeas,realAvgMeas->GetMarkerColor(), MS(1,false,alt));
1638  SetAttr(realAvgSig, realAvgSig ->GetMarkerColor(), MS(2,false,alt));
1639  SetAttr(realAvgBg, realAvgSig ->GetMarkerColor(), MS(3,false,alt));
1640  SetAttr(simAvgMeas, simAvgSig ->GetMarkerColor(), MS(1,true, alt));
1641  SetAttr(simAvgSig, simAvgSig ->GetMarkerColor(), MS(2,true, alt));
1642  SetAttr(simAvgBg, simAvgBg ->GetMarkerColor(), MS(3,true, alt));
1643  realAvgMeas->SetYTitle("#LT M#GT");
1644  realAvgSig ->SetYTitle("#LT S#GT");
1645  realAvgBg ->SetYTitle("#LT B#GT");
1646  simAvgMeas ->SetYTitle("#LT M'#GT");
1647  simAvgSig ->SetYTitle("#LT S'#GT");
1648  simAvgBg ->SetYTitle("#LT B'#GT");
1649 
1650  // Put everything together
1651  THStack* summary = new THStack("summary", dndeta->GetYaxis()->GetTitle());
1652  summary->Add(dndeta, "e2");
1653  summary->Add(truth, "e");
1654  summary->Add(realAvgMeas);
1655  summary->Add(simAvgMeas);
1656  summary->Add(realAvgSig);
1657  summary->Add(simAvgSig);
1658  summary->Add(realAvgBg);
1659  summary->Add(simAvgBg);
1660  if (scaleH) {
1661  TH1* tmp = static_cast<TH1*>(scaleH->Clone("scaleBg"));
1662  Double_t sf = TMath::Power(10,Int_t(TMath::Log10(realBg->GetMaximum())+.3));
1663  tmp->SetTitle(Form("%s#times%3.0f",scaleH->GetTitle(),sf));
1664  tmp->Scale(sf);
1665  tmp->SetDirectory(0);
1666  summary->Add(tmp);
1667  }
1668 
1669  ClearCanvas();
1670  fBody->SetLeftMargin(0.15);
1671  TLegend* l = DrawInPad(fBody, 0, summary, "nostack leg2 grid");
1672  summary->GetHistogram()->SetXTitle("#eta");
1673  summary->GetHistogram()->SetYTitle(ObsTitle());
1674  summary->GetHistogram()->GetYaxis()->SetTitleOffset(1.7);
1675  fBody->Modified();
1676  fBody->Update();
1677  fBody->cd();
1678  l->SetName("legend");
1679  l->SetFillStyle(0);
1680  l->SetBorderSize(0);
1681  ModLegend(fBody, l, fBody->GetLeftMargin(), .7, .99, .99);
1682  Double_t lim = .5;
1683  TF1* f = new TF1("ff", "pol0", -lim, lim);
1684  TFitResultPtr r = dndeta->Fit(f, "N0QSR", "", -lim, lim);
1685  // r->Print();
1686  TLatex* ltx = new TLatex(fBody->GetLeftMargin()+
1687  (1-fBody->GetLeftMargin()-fBody->GetRightMargin())/2,
1688  .3,
1689  Form("%s|_{|#eta|<%3.1f}"
1690  "=%6.1f#pm%6.1f (#chi^{2}/#nu=%5.2f)",
1691  ObsTitle(), lim, r->Parameter(0),
1692  r->ParError(0), r->Chi2()/r->Ndf()));
1693  ltx->SetTextAlign(22);
1694  ltx->SetNDC();
1695  ltx->SetTextFont(42);
1696  ltx->SetTextSize(0.035);
1697  ltx->Draw();
1698  printf("%6.1f +/- %6.1f (%5.2f - %5.3f+/-%5.3f)",
1699  r->Parameter(0), r->ParError(0),
1700  r->Chi2()/r->Ndf(), scale, scaleE);
1701  if (scale > 0)
1702  ltx->DrawLatex(fBody->GetLeftMargin()+
1703  (1-fBody->GetLeftMargin()-fBody->GetRightMargin())/2,
1704  .25,
1705  Form("%sk%s=%4.2f #pm %4.2f",
1706  scaleH ? "#LT" : "",
1707  scaleH ? "#GT" : "",
1708  scale, scaleE));
1709  summary->GetHistogram()->GetListOfFunctions()->Add(l);
1710  summary->GetHistogram()->GetListOfFunctions()->Add(ltx);
1711 
1712  fBody->Modified();
1713  if (summary->GetHistogram()->GetMinimum() < 1)
1714  summary->SetMinimum(1);
1715  summary->SetMaximum(1.5*summary->GetMaximum("nostack"));
1716 
1717  if (stack) {
1718  dndeta->SetTitle(fLastBin);
1719  stack->Add(dndeta, "e2");
1720  stack->Add(truth, "e");
1721  }
1722 
1723  if (out) {
1724  out->cd();
1725  result ->Write();
1726  dndeta ->Write();
1727  summary ->Write();
1728  if (scaleH) scaleH->Write();
1729  (new TParameter<double>("scale", scale)) ->Write();
1730  (new TParameter<double>("scaleE", scaleE))->Write();
1731 
1732  TDirectory* details = out->mkdir("details");
1733  details->cd();
1734  realMeas->Write();
1735  realBg ->Write();
1736  realSig ->Write();
1737  simMeas ->Write();
1738  simBg ->Write();
1739  simSig ->Write();
1740  trueGen ->Write();
1741  alpha ->Write();
1742  fiducial->Write();
1743  GetH1(realData, "delta") ->Write("realDataDelta");
1744  GetH1(GetC(realList, "injection"), "delta")->Write("realInjDelta");
1745  GetH1(simData, "delta") ->Write("simDataDelta");
1746  GetH1(GetC(simList, "injection"), "delta") ->Write("simInjDelta");
1747 
1748  TDirectory* averages = out->mkdir("averages");
1749  averages->cd();
1750  realAvgMeas ->Write();
1751  realAvgSig ->Write();
1752  realAvgBg ->Write();
1753  simAvgMeas ->Write();
1754  simAvgSig ->Write();
1755  simAvgBg ->Write();
1756  truth ->Write();
1757  realIPz ->Write();
1758  simIPz ->Write();
1759  trueIPz ->Write();
1760  }
1761  PrintCanvas(Form("%s - %s", fLastBin.Data(),ObsTitle()));
1762 }
1763 //____________________________________________________________________
1765  Double_t c1, Double_t c2,
1766  Container* realList, Container* simList,
1767  THStack* stack, TDirectory* out)
1768 {
1769  fLastBin.Form("%5.1f%% - %5.1f%%", c1, c2);
1770  printf("Centrality bin %5.1f%% - %5.1f%%: ", c1, c2);
1771  TString name;
1772  name.Form("cent%03dd%02d_%03dd%02d",
1773  Int_t(c1), Int_t(c1*100)%100,
1774  Int_t(c2), Int_t(c2*100)%100);
1775  Container* realListBin = GetC(realList, name);
1776  Container* simListBin = GetC(simList, name);
1777  if (!realListBin || !simListBin) {
1778  Warning("PostBin", "Missing bin for %5.1f%% - %5.1f%% (%s): %p %p",
1779  c1, c2, name.Data(), realListBin, simListBin);
1780  return;
1781  }
1782  const Color_t cc[] = { kMagenta+2, // 0
1783  kBlue+2, // 1
1784  kAzure-1, // 2 // 10,
1785  kCyan+2, // 3
1786  kGreen+1, // 4
1787  kSpring+5, // 5 //+10,
1788  kYellow+1, // 6
1789  kOrange+5, // 7 //+10,
1790  kRed+1, // 8
1791  kPink+5, // 9 //+10,
1792  kBlack }; // 10
1793  Color_t color = (bin <= 0 ? cc[10] : cc[(bin-1)%11]);
1794  TDirectory* binOut = 0;
1795  Bool_t alt = (what & kAltMarker);
1796  if (out) binOut = out->mkdir(name);
1797 
1798  if (what & kSimInfo) DrawSim (simListBin);
1799  if (what & kSpecies) DrawSpecies (simListBin);
1800  if (what & kDeltas) DrawDeltas (realListBin, simListBin);
1801  if (what & kBackgrounds) DrawBackground(realListBin, simListBin);
1802  if (what & kAlphas) DrawAlpha (realListBin, simListBin);
1803  if (what & kdNdetas) DrawdNdeta (realListBin, what & kRealComb,
1804  simListBin, what & kSimComb,
1805  color, stack, binOut, alt);
1806  printf(" done \n");
1807 }
1808 //====================================================================
1810  UShort_t maxBins,
1811  const char* dataName,
1812  const char* simName)
1813 {
1814  fLandscape = what & kLandscape;
1815  fPause = what & kPause;
1816  fPDF = what & kPDF;
1818  else if (what & kDynCentScale) fCombinatoricsScaleMode = kPerCent;
1820 
1821  TFile* dataFile = 0;
1822  TFile* simFile = 0;
1823  if (!(dataFile = OpenFile(dataName))) return;
1824  if (!(simFile = OpenFile(simName))) return;
1825 
1826  Container* realSums = GetC(dataFile, "MiddNdetaSums");
1827  Container* realRess = GetC(dataFile, "MiddNdetaResults");
1828  Container* simSums = GetC(simFile, "MidMCdNdetaSums");
1829  Container* simRess = GetC(simFile, "MidMCdNdetaResults");
1830  if (!realSums || !realRess || !simSums || !simRess) return;
1831 
1832  Container* params = GetC(realSums, "parameters");
1833  fDeltaCut = GetD(params, "DeltaCut");
1834  fTailDelta = GetD(params, "TailDelta");
1835  fMaxDelta = GetD(params, "MaxDelta");
1836 
1837  CreateCanvas(Form("MiddNdeta_0x%04x.pdf", what));
1838  if (what & kGeneral) DrawGeneral(realRess, simRess);
1839  if (what & kParameters) DrawParams(realSums, what & kRealComb,
1840  simSums, what & kSimComb);
1841 
1842  TH1* realCent = CopyH1(realSums, "cent", "realCent");
1843  TH1* simCent = CopyH1(simSums, "cent", "simCent");
1844  TH1* realIPz = GetH1(realSums, "ipz");
1845  TH1* simIPz = GetH1(simSums, "ipz");
1846 
1847  if (!CheckConsistency(realCent, simCent)) {
1848  Warning("Post", "Centrality bins are incompatible, giving up");
1849  return;
1850  }
1851  if (!CheckConsistency(realIPz, simIPz)) {
1852  Warning("Post", "Centrality bins are incompatible, giving up");
1853  return;
1854  }
1855  if (what & kClosure) realRess = simRess;
1856 
1857  THStack* stack = new THStack("all", ObsTitle());
1858  // Draw "min-bias" bin
1859  ProcessBin(what, 0, 0, 100, realRess, simRess, stack);
1860 
1861  TFile* out = TFile::Open(Form("MiddNdeta_0x%1x.root", 0x3 & what),
1862  "RECREATE");
1863  realCent->Write();
1864  simCent ->Write();
1865  for (Int_t i = 1; i <= realCent->GetNbinsX() && i <= maxBins ; i++) {
1866  Double_t c1 = realCent->GetXaxis()->GetBinLowEdge(i);
1867  Double_t c2 = realCent->GetXaxis()->GetBinUpEdge (i);
1868 
1869  ProcessBin(what, i, c1, c2, realRess, simRess, stack, out);
1870  }
1871 
1872  if (stack->GetHists() && stack->GetHists()->GetEntries() > 0) {
1873  ClearCanvas();
1874  DrawInPad(fBody, 0, stack, "nostack grid e2");
1875  if (stack->GetHistogram()->GetMinimum() < 1)
1876  stack->SetMinimum(1);
1877  out->cd();
1878  stack->Write();
1879  TFile* res = 0;
1880  {
1881  SuppressGuard g(kFatal);
1882  res = TFile::Open("result.root","READ");
1883  }
1884  if (res) {
1885  THStack* other = static_cast<THStack*>(res->Get("result"));
1886  if (other)
1887  DrawInPad(fBody, 0, other, "nostack same");
1888  out->cd();
1889  other->Write();
1890  }
1891  PrintCanvas(ObsTitle());
1892  }
1893  Printf("Results stored in %s", out->GetName());
1894  out->Write();
1895  CloseCanvas();
1896 }
1897 
1898 #endif
1899 // Local Variables:
1900 // mode: C++
1901 // End:
1902 
Int_t color[]
print message on plot with ok/not ok
const char * filename
Definition: TestFCM.C:1
static TH2 * CopyH2(Container *parent, const char *name, const char *newName=0)
const Color_t cc[]
Definition: DrawKs.C:1
THStack * Make2Stack(const char *name, const char *title, Container *realList, Container *simList, Option_t *dataOpt="", Option_t *simOpt="")
TH1 * FindDelta(Container *ress, const char *sub, Bool_t scaled=false)
void Run(UInt_t what=kDefault, UShort_t maxBins=9, const char *dataName="data.root", const char *simName="sim.root")
const char * ObsTitle() const
TLegend * DrawInPad(TVirtualPad *c, Int_t pad, TObject *o, Option_t *opt)
double Double_t
Definition: External.C:58
TCanvas * fCanvas
TPad * fTop
Double_t fAlphaMax
void DrawSim(Container *simList)
const char * title
Definition: MakeQAPdf.C:27
TPad * fBody
static TH1 * AverageOverIPz(TH2 *h, const char *name, UShort_t mode, TH1 *ipz, TH2 *mask=0)
TH1 * SetAttr(TH1 *h, Color_t color, Style_t marker=20, Double_t size=1., Style_t fill=0, Style_t line=1, Width_t width=1)
ECombinatoricsScaleMode
AliTrackletdNdetaPost()
Double_t GetDeltaScale(Container *c, Double_t &e) const
static Int_t GetI(Container *parent, const char *name, Int_t def=-1)
void PrintH(TH2 *h, Int_t prec=2)
void DrawClusters(Int_t pad, const char *name, Container *realList, Container *simList)
static Double_t RatioE(Double_t n, Double_t en, Double_t d, Double_t ed, Double_t &er)
TH1 * GetDeltaIntDist(Container *c) const
TCanvas * c
Definition: TestFitELoss.C:172
AliStack * stack
void ClearCanvas()
static Double_t Integrate(TH1 *h, Double_t min, Double_t max, Double_t &err)
Bool_t scaled
Double_t fDeltaCut
void DrawAlpha(TVirtualPad *c, Int_t pad, Container *ress, const char *name, const char *pretty)
void DrawDeltas(Container *realList, Container *simList)
static Container * GetC(Container *parent, const char *name)
static TH1 * CopyH1(Container *parent, const char *name, const char *newName=0)
TLatex * fHeader
static Int_t GetB(Container *parent, const char *name, Bool_t def=false)
int Int_t
Definition: External.C:63
Double_t DrawBackground(TVirtualPad *c, Int_t pad, Container *ress, const char *name, const char *pretty)
void PrintCanvas(const TString &title, Float_t size=.7)
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
Bool_t fPause
Style_t MS(Int_t what, Bool_t sim, Bool_t alt) const
Definition: External.C:228
Double_t fAlphaMin
void ModLegend(TVirtualPad *p, TLegend *l, Double_t x1, Double_t y1, Double_t x2, Double_t y2)
void CloseCanvas()
void CreateCanvas(const TString &outputName)
void DrawParam(const char *name, Double_t &y, const char *val)
static TH1 * Scale(TH1 *h, Double_t x, Double_t xe)
Int_t mode
Definition: anaM.C:41
static Double_t GetD(Container *parent, const char *name, Double_t def=-1)
Bool_t fLandscape
void ProcessBin(UInt_t what, Int_t bin, Double_t c1, Double_t c2, Container *realList, Container *simList, THStack *stack=0, TDirectory *out=0)
Double_t fMaxDelta
SuppressGuard(Int_t lvl=2000)
Definition: External.C:220
static TH1 * GetH1(Container *parent, const char *name)
TFile * file
TList with histograms for a given trigger.
TH2 * CutAlpha(TH2 *alpha)
unsigned short UShort_t
Definition: External.C:28
void DrawSpecies(Container *mcList)
Double_t fTailDelta
void DrawParams(Container *pars, const char *title, Bool_t comb)
const char Option_t
Definition: External.C:48
void DrawdNdeta(Container *realList, Bool_t realComb, Container *simList, Bool_t simComb, Color_t color=kBlack, THStack *stack=0, TDirectory *out=0, Bool_t alt=false)
bool Bool_t
Definition: External.C:53
ECombinatoricsScaleMode fCombinatoricsScaleMode
static Bool_t CheckConsistency(const TH1 *h1, const TH1 *h2)
TString fLastBin
TFile * OpenFile(const char *filename)
Bool_t fPDF
Double_t fCombinatoricsScale
static TObject * GetO(Container *parent, const char *name, TClass *cls)
TH1 * FindSub(Container *ress, const char *sub, const char *name)
TString fLastTitle
Definition: External.C:196
void DrawGeneral(Container *realList, Container *simList)
static TH2 * GetH2(Container *parent, const char *name)
TDirectoryFile * dir
TH1 * GetDeltaScaleDist(Container *c) const