AliPhysics  35e5fca (35e5fca)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTrackletAODdNdeta.C
Go to the documentation of this file.
1 
12 #include <AliAnalysisTaskSE.h>
13 #include <AliTrackletAODUtils.C>
14 #include <AliAODTracklet.C>
15 #ifndef __CINT__
16 #include <cctype>
17 #include "AliAODSimpleHeader.C"
18 #include <AliVVertex.h>
19 #include <AliVertex.h>
20 #include <AliAnalysisManager.h>
21 #include <AliVEventHandler.h>
22 #include <AliInputEventHandler.h>
23 #include <AliMultSelection.h>
24 #include <AliCentrality.h>
25 #include <AliLog.h>
26 #include <TClonesArray.h>
27 #include "AliTrackletWeights.C"
28 #include <TUrl.h>
29 #include <TFile.h>
30 #include <TGraphErrors.h>
31 #else
32 // class AliAODTracklet;
34 class AliVEvent;
35 class AliMultSelection; // Auto-load
36 class TClonesArray;
37 #endif
38 
39 //====================================================================
46  public AliTrackletAODUtils
47 {
48 public:
52  enum {
62  kPrimaryMask = 0x0,
71  };
72  // -----------------------------------------------------------------
74  enum {
75  kAll=1, // Count all events
76  kEvent, // Have event
77  kTracklets, // Have tracklets
78  kTrigger, // Have trigger
79  kIP, // Have IP
80  kCentrality, // Have centrality
81  kCompleted // Have completed
82  };
83  enum EStatus {
84  kOKEvent = ((1<<(kAll -1))|
85  (1<<(kEvent -1))|
86  (1<<(kTracklets-1))),
87  kOKTrigger = (kOKEvent|(1<<(kTrigger -1))),
88  kOKIPz = (kOKTrigger|(1<<kIP-1)),
90  };
92  typedef TList Container;
100  AliTrackletAODdNdeta(const char* name);
128  virtual void Print(Option_t* option="") const;
135  Bool_t Connect(const char* sumFile=0,const char* resFile=0);
146  static AliTrackletAODdNdeta* Create(Bool_t mc=false,
147  const char* weights=0,
148  const char* sumFile=0,
149  const char* resFile=0);
150  /* @} */
151 
152  // -----------------------------------------------------------------
165  void UserExec(Option_t*);
170  void FinishTaskOutput() { /*WorkerFinalize();*/ }
175  void Terminate(Option_t*);
176  /* @} */
177 
178  // -----------------------------------------------------------------
188  void SetCentralityMethod(const TString& name) { fCentMethod = name; }
196  {
197  SetAxis(fCentAxis,n,bins);
198  }
204  void SetCentralityAxis(const TString& spec)
205  {
206  SetAxis(fCentAxis, spec, "-:,");
207  }
215  void SetIPzAxis(Int_t n, Double_t min, Double_t max)
216  {
217  SetAxis(fIPzAxis, n, min, max);
218  }
225  void SetIPzAxis(Int_t n, Double_t max)
226  {
227  SetAxis(fIPzAxis, n, max);
228  }
234  void SetIPzAxis(const TString& spec)
235  {
236  SetAxis(fIPzAxis, spec);
237  }
245  void SetPhiAxis(Int_t n, Double_t min, Double_t max)
246  {
247  SetAxis(fPhiAxis, n, min, max);
248  }
255  void SetPhiAxis(Int_t n, Double_t max)
256  {
257  SetAxis(fPhiAxis, n, max);
258  }
266  void SetEtaAxis(Int_t n, Double_t min, Double_t max)
267  {
268  SetAxis(fEtaAxis, n, min, max);
269  }
276  void SetEtaAxis(Int_t n, Double_t max)
277  {
278  SetAxis(fEtaAxis, n, max);
279  }
285  void SetEtaAxis(const TString& spec)
286  {
287  SetAxis(fEtaAxis, spec);
288  }
294  void SetDPhiShift(Double_t x=0.0045) { fDPhiShift = x; }
306  void SetDeltaCut(Double_t x=1.5) { fDeltaCut = x; }
312  void SetMaxDelta(Double_t x=25) { fMaxDelta = x; }
318  void SetTailDelta(Double_t x=5) { fTailDelta = x; }
324  void SetTailMaximum(Double_t x=-1) { fTailMax = x; }
338  void SetAbsMinCent(Double_t x=-1) { fAbsMinCent = x; }
344  virtual void SetWeights(AliTrackletBaseWeights* w) {};
350  virtual void SetWeightCalc(UChar_t mode=0) {}
356  virtual void SetWeightMask(UChar_t mask=0xFF) {}
364  virtual void SetWeightInverse(Bool_t inv) {}
370  virtual void SetWeightVeto(UChar_t veto=0xFF) {}
371  /* @} */
372  //__________________________________________________________________
376  struct Sub : public TObject
377  {
382  Sub(const char* name="") : fName(name), fContainer(0), fDebug(0) {}
388  Sub(const Sub& o) : TObject(o), fName(o.fName), fContainer(0), fDebug(0) {}
392  virtual ~Sub() {}
398  Sub& operator=(const Sub&) { return *this; }
402  const char* GetName() const { return fName.Data(); }
413  virtual Bool_t WorkerInit(Container* parent,
414  const TAxis& etaAxis,
415  const TAxis& ipzAxis,
416  const TAxis& deltaAxis)
417  {
418  fContainer = new Container;
419  fContainer->SetName(fName);
420  fContainer->SetOwner();
421  if (parent) parent->Add(fContainer);
422  return true;
423  }
434  virtual Bool_t ProcessTracklet(AliAODTracklet* tracklet,
435  Double_t ipz,
436  UShort_t signal,
437  Double_t weight) = 0;
447  virtual Bool_t FinalizeInit(Container* parent) = 0;
458  virtual Bool_t MasterFinalize(Container* parent,
459  TH1* ipz,
460  Double_t tailCut,
461  Double_t tailMax) = 0;
467  virtual void SetDebug(UShort_t lvl) { fDebug = lvl; }
468  protected:
472  Container* fContainer;
475  ClassDef(Sub,1);
476  };
477 
478  //__________________________________________________________________
488  struct Histos : public Sub
489  {
494  Histos(const char* name="", Color_t color=kBlack, Style_t style=1,
495  UChar_t mask=0, UChar_t veto=0)
496  : Sub(name),
497  fMask(mask),
498  fVeto(veto),
499  fEtaIPz(0),
500  fEtaDeltaIPz(0),
501  fEtaDeltaPdg(0),
502  fEtaPdgIPz(0),
503  fEtaPdg(0),
504  fEtaPt(0),
505  fColor(color),
506  fStyle(style)
507  {}
513  Histos(const Histos& o)
514  : Sub(o),
515  fMask(o.fMask),
516  fVeto(o.fVeto),
517  fEtaIPz(0),
518  fEtaDeltaIPz(0),
519  fEtaDeltaPdg(0),
520  fEtaPdgIPz(0),
521  fEtaPdg(0),
522  fEtaPt(0),
523  fColor(o.fColor),
524  fStyle(o.fStyle)
525  {}
529  virtual ~Histos() {}
535  Histos& operator=(const Histos&) { return *this; }
540  {
541  return fMask == kMeasuredMask && fVeto == kMeasuredVeto;
542  }
547  {
548  return fMask == kInjectedMask && fVeto == kInjectedVeto;
549  }
554  {
556  }
562  {
563  return fMask == kDistinctMask && fVeto == kDistinctVeto;
564  }
569  {
570  return fMask == kPrimaryMask && fVeto == kPrimaryVeto;
571  }
576  {
577  return fMask == kSecondaryMask && fVeto == kSecondaryVeto;
578  }
583  {
584  return fMask == kGeneratedMask && fVeto == kGeneratedVeto;
585  }
596  Bool_t WorkerInit(Container* parent,
597  const TAxis& etaAxis,
598  const TAxis& ipzAxis,
599  const TAxis& deltaAxis);
606  void SetAttr(Color_t c, Style_t m);
618  Double_t ipz,
619  UShort_t signal,
620  Double_t weight);
630  Bool_t FinalizeInit(Container* parent);
641  Bool_t MasterFinalize(Container* parent,
642  TH1* ipz,
643  Double_t tailCut,
644  Double_t tailMax);
645 
646  UChar_t GetMask() const { return fMask; }
647  UChar_t GetVeto() const { return fVeto; }
653  void Print(Option_t* option="") const;
654  protected:
655  Bool_t ProjectEtaPdg(Container* result, Int_t nEvents);
656  Bool_t ProjectEtaDeltaPdgPart(Container* result,
657  Int_t nEvents,
658  Double_t tailDelta,
659  Double_t tailMax,
660  const TString& pre,
661  const TString& tit);
662  Bool_t ProjectEtaDeltaPdg(Container* result,
663  Int_t nEvents,
664  Double_t tailCut,
665  Double_t tailMax);
666  Bool_t ProjectEtaPdgIPz(Container* result,
667  TH1* ipz,
668  const TString& shn);
669  const UChar_t fMask;
670  const UChar_t fVeto;
677  public:
678  Color_t fColor;
679  Style_t fStyle;
680 
681  ClassDef(Histos,1);
682  };
683 
684  //__________________________________________________________________
700  struct CentBin : public Sub
701  {
706  : Sub(""),
707  fSubs(0),
708  fLow(0),
709  fHigh(0),
710  fStatus(0),
711  fIPz(0),
712  fCent(0),
713  fCentIPz(0),
714  fMeasured(0),
715  fInjection(0)
716  {
717  }
724  CentBin(Double_t c1, Double_t c2);
730  CentBin(const CentBin& o)
731  : Sub(o),
732  fSubs(0),
733  fLow(o.fLow),
734  fStatus(0),
735  fIPz(0),
736  fCent(0),
737  fCentIPz(0),
738  fHigh(o.fHigh),
739  fMeasured(0),
740  fInjection(0)
741  {}
745  virtual ~CentBin() {}
751  CentBin& operator=(const CentBin&) { return *this; }
763  Histos* MakeHistos(const char* name,
764  Color_t color,
765  Style_t style,
766  UShort_t mask,
767  UShort_t veto);
778  Bool_t WorkerInit(Container* parent,
779  const TAxis& etaAxis,
780  const TAxis& ipzAxis,
781  const TAxis& deltaAxis);
787  Bool_t IsAllBin() const;
797  Bool_t Accept(UInt_t status, Double_t cent, Double_t ipz);
809  Double_t ipz,
810  UShort_t signal,
811  Double_t weight);
816  void Completed();
826  Bool_t FinalizeInit(Container* parent);
837  Bool_t MasterFinalize(Container* parent,
838  TH1* ipz,
839  Double_t tailCut,
840  Double_t tailMax);
853  Bool_t EstimateBackground(Container* result,
854  Container* measCont,
855  Container* genCont,
856  Histos* h,
857  Double_t tailCut,
858  Double_t tailMax);
864  void Print(Option_t* option="") const;
870  virtual void SetDebug(UShort_t lvl);
871  protected:
872  Container* fSubs;
878  TProfile* fCentIPz;
881 
882  ClassDef(CentBin,2);
883  };
884 
885 protected:
886  // -----------------------------------------------------------------
896  virtual Bool_t WorkerInit();
903  Bool_t InitCentBins(Container* existing);
913  {
914  return new CentBin(c1, c2);
915  }
927  virtual CentBin* InitCentBin(Int_t bin,
928  Float_t c1,
929  Float_t c2,
930  Container* existing,
931  const TAxis& deltaAxis);
932  /* @} */
933 
934  // -----------------------------------------------------------------
939  virtual const char* GetBranchName() const { return "AliAODTracklets"; }
947  static UInt_t Bin2Flag(UInt_t bin);
957  UInt_t CheckEvent(Double_t& cent,
958  const AliVVertex*& ip,
959  TClonesArray*& tracklets);
967  TClonesArray* FindTracklets(AliVEvent* event);
977  const AliVVertex* FindSimpleIP(AliVEvent* event,
978  Double_t maxDispersion=0.04,
979  Double_t maxZError=0.25);
989  const AliVVertex* FindRealIP(AliVEvent* event,
990  Double_t maxDispersion=0.04,
991  Double_t maxZError=0.25);
1001  const AliVVertex* CheckIP(const AliVVertex* ip,
1002  Double_t maxDispersion=0.04,
1003  Double_t maxZError=0.25);
1004 
1014  const AliVVertex* FindIP(AliVEvent* event,
1015  Double_t maxDispersion=0.04,
1016  Double_t maxZError=0.25);
1027  Double_t FindMultCentrality(AliVEvent* event, Int_t& nTracklets);
1037  Double_t FindCompatCentrality(AliVEvent* event);
1049  Double_t FindCentrality(AliVEvent* event, Int_t& nTracklets);
1055  Bool_t FindTrigger();
1056  /* @} */
1057 
1058  // -----------------------------------------------------------------
1072  virtual Double_t LookupWeight(AliAODTracklet* tracklet,
1073  Double_t cent,
1074  Double_t ipz);
1084  virtual UShort_t CheckTracklet(AliAODTracklet* tracklet) const;
1085  /* @} */
1086 
1087  // -----------------------------------------------------------------
1100  void ProcessEvent(UInt_t status,
1101  Double_t cent,
1102  const AliVVertex* ip,
1103  TClonesArray* tracklets);
1104  /* @} */
1105 
1106  // -----------------------------------------------------------------
1118  virtual Bool_t MasterFinalize(Container* results);
1119  /* @} */
1120 
1128  static TH1* MakeStatus(Container* container);
1129 
1130  // -----------------------------------------------------------------
1132  Container* fContainer;
1134  Container* fCentBins;
1137 
1139 
1141 
1143 
1144  TProfile* fNBareVsGood;
1146  TProfile* fNBareVsFake;
1148  TProfile* fNTrackletVsGood;
1150  TProfile* fNTrackletVsFake;
1154  TProfile* fNGeneratedVsFake;
1156  TProfile* fCentTracklets;
1157 
1158  TProfile* fCentEst;
1159 
1187 
1189 };
1190 //====================================================================
1197 {
1198 public:
1200  typedef TList Container;
1206  {}
1210  AliTrackletAODMCdNdeta(const char* name)
1211  : AliTrackletAODdNdeta(name)
1212  {}
1220  {}
1225  //__________________________________________________________________
1242  {
1248  fCombinatorics(0),
1249  fDistinct(0),
1250  fPrimaries(0),
1251  fSecondaries(0),
1252  fGenerated(0)
1253  {
1254  }
1261  CentBin(Double_t c1, Double_t c2);
1267  CentBin(const CentBin& o)
1269  fCombinatorics(0),
1270  fPrimaries(0),
1271  fSecondaries(0),
1272  fGenerated(0)
1273  {}
1277  virtual ~CentBin() {}
1283  CentBin& operator=(const CentBin&) { return *this; }
1284  protected:
1290 
1291  ClassDef(CentBin,1);
1292  };
1293 
1294 protected:
1295  // -----------------------------------------------------------------
1309  {
1310  return new CentBin(c1, c2);
1311  }
1312  /* @} */
1313 
1314  // -----------------------------------------------------------------
1319  virtual const char* GetBranchName() const { return "AliAODMCTracklets"; }
1320  /* @} */
1321 
1323 };
1324 //====================================================================
1331 {
1332 public:
1334  typedef TList Container;
1340  fWeights(0),
1341  fEtaWeight(0),
1342  fWeightCorr(0)
1343  {}
1348  : AliTrackletAODMCdNdeta(name),
1349  fWeights(0),
1350  fEtaWeight(0),
1351  fWeightCorr(0)
1352  {}
1360  fWeights(0),
1361  fEtaWeight(0),
1362  fWeightCorr(0)
1363  {}
1382  void Print(Option_t* option="") const;
1383 
1384  // -----------------------------------------------------------------
1394  void SetWeights(AliTrackletBaseWeights* w) { fWeights = w; }
1400  void SetWeightCalc(UChar_t mode=0) {
1401  Info("SetWeightCalc", "Use=%d", mode);
1402  if (fWeights) fWeights->SetCalc(mode);}
1408  void SetWeightMask(UChar_t mask=0xFF) {
1409  Info("SetWeightMask", "mask=0x%x", mask);
1410  if (fWeights) fWeights->SetMask(mask); }
1416  void SetWeightVeto(UChar_t veto=0x0) {
1417  Info("SetWeightVeto", "veto=0x%x", veto);
1418  if (fWeights) fWeights->SetVeto(veto); }
1426  void SetWeightInverse(Bool_t inv) { if (fWeights) fWeights->SetInverse(inv); }
1427  /* @} */
1428 protected:
1429  // -----------------------------------------------------------------
1439  Bool_t WorkerInit();
1449  virtual Double_t LookupWeight(AliAODTracklet* tracklet,
1450  Double_t cent,
1451  Double_t ipz);
1452  /* @} */
1453  // -----------------------------------------------------------------
1465  Bool_t MasterFinalize(Container* results);
1466  /* @} */
1467  // -----------------------------------------------------------------
1469  TProfile2D* fEtaWeight;
1472 };
1473 
1474 //____________________________________________________________________
1476  : AliAnalysisTaskSE(),
1477  fContainer(0),
1478  fCentBins(0),
1479  fIPz(0),
1480  fCent(0),
1481  fStatus(0),
1482  fEtaPhi(0),
1483  fNBareVsGood(0),
1484  fNBareVsFake(0),
1485  fNTrackletVsGood(0),
1486  fNTrackletVsFake(0),
1487  fNGeneratedVsGood(0),
1488  fNGeneratedVsFake(0),
1489  fCentTracklets(0),
1490  fCentEst(0),
1491  fCentMethod(""),
1492  fCentIdx(-1),
1493  fTrkIdx(-1),
1494  fCentAxis(1,0,0),
1495  fIPzAxis(1,0,0),
1496  fEtaAxis(1,0,0),
1497  fPhiAxis(1,0,0),
1498  fMaxDelta(0),
1499  fTailDelta(0),
1500  fTailMax(-1),
1501  fDPhiShift(0),
1502  fShiftedDPhiCut(0),
1503  fDeltaCut(0),
1504  fAbsMinCent(-1)
1505 {}
1506 //____________________________________________________________________
1508  : AliAnalysisTaskSE(name),
1509  fContainer(0),
1510  fCentBins(0),
1511  fIPz(0),
1512  fCent(0),
1513  fStatus(0),
1514  fEtaPhi(0),
1515  fNBareVsGood(0),
1516  fNBareVsFake(0),
1517  fNTrackletVsGood(0),
1518  fNTrackletVsFake(0),
1519  fNGeneratedVsGood(0),
1520  fNGeneratedVsFake(0),
1521  fCentTracklets(0),
1522  fCentEst(0),
1523  fCentMethod("V0M"),
1524  fCentIdx(-1),
1525  fTrkIdx(-1),
1526  fCentAxis(10,0,100),
1527  fIPzAxis(30,-15,+15),
1528  fEtaAxis(16,-2,+2),
1529  fPhiAxis(100,0,TMath::TwoPi()),
1530  fMaxDelta(25),
1531  fTailDelta(5),
1532  fTailMax(-1),
1533  fDPhiShift(0.0045),
1534  fShiftedDPhiCut(-1),
1535  fDeltaCut(1.5),
1536  fAbsMinCent(-1)
1537 {
1538  FixAxis(fCentAxis, "Centrality [%]");
1539  FixAxis(fIPzAxis, "IP_{#it{z}} [cm]");
1540  FixAxis(fEtaAxis, "#eta");
1541  FixAxis(fPhiAxis, "#phi");
1542 
1543  DefineOutput(1, Container::Class());
1544  DefineOutput(2, Container::Class());
1545 }
1546 //____________________________________________________________________
1548  : AliAnalysisTaskSE(o),
1549  fContainer(0),
1550  fCentBins(0),
1551  fIPz(0),
1552  fCent(0),
1553  fStatus(0),
1554  fEtaPhi(0),
1555  fNBareVsGood(0),
1556  fNBareVsFake(0),
1557  fNTrackletVsGood(0),
1558  fNTrackletVsFake(0),
1559  fNGeneratedVsGood(0),
1560  fNGeneratedVsFake(0),
1561  fCentTracklets(0),
1562  fCentEst(0),
1563  fCentMethod(o.fCentMethod),
1564  fCentIdx(o.fCentIdx),
1565  fTrkIdx(o.fTrkIdx),
1566  fCentAxis(o.fCentAxis),
1567  fIPzAxis(o.fIPzAxis),
1568  fEtaAxis(o.fEtaAxis),
1569  fPhiAxis(o.fPhiAxis),
1570  fMaxDelta(o.fMaxDelta),
1571  fTailDelta(o.fTailDelta),
1572  fTailMax(o.fTailMax),
1573  fDPhiShift(o.fDPhiShift),
1574  fShiftedDPhiCut(o.fShiftedDPhiCut),
1575  fDeltaCut(o.fDeltaCut),
1576  fAbsMinCent(o.fAbsMinCent)
1577 {}
1578 //____________________________________________________________________
1581 {
1582  if (&o == this) return *this;
1583  if (fContainer) {
1584  delete fContainer;
1585  fContainer = 0;
1586  }
1587  if (fCentBins) {
1588  delete fCentBins;
1589  fCentBins = 0;
1590  }
1591  fIPz = 0;
1592  fCent = 0;
1594  fCentIdx = o.fCentIdx;
1595  fTrkIdx = o.fTrkIdx;
1596  fCentAxis = o.fCentAxis;
1597  fIPzAxis = o.fIPzAxis;
1598  fEtaAxis = o.fEtaAxis;
1599  fPhiAxis = o.fPhiAxis;
1600  fMaxDelta = o.fMaxDelta;
1601  fTailDelta = o.fTailDelta;
1602  fTailMax = o.fTailMax;
1603  fDPhiShift = o.fDPhiShift;
1605  fDeltaCut = o.fDeltaCut;
1607  return *this;
1608 }
1609 //____________________________________________________________________
1612  o)
1613 {
1614  if (&o == this) return *this;
1616  return *this;
1617 }
1618 //____________________________________________________________________
1620  const char* resFile)
1621 {
1622  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1623  if (!mgr) {
1624  AliError("No analysis manager to connect to.");
1625  return false;
1626  }
1627 
1628  // Add to the manager
1629  mgr->AddTask(this);
1630 
1631  // Create and connect output containers
1632  TString sumOut;
1633  TString resOut;
1634  if (sumFile && sumFile[0] != '\0') sumOut = sumFile;
1635  if (resFile && resFile[0] != '\0') resOut = resFile;
1636  else if (sumFile && sumFile[0] != '\0') resOut = sumFile;
1637 
1638  // If the string is null or 'default' connect to standard output file
1639  if (sumOut.IsNull() || sumOut.EqualTo("default", TString::kIgnoreCase))
1640  sumOut = AliAnalysisManager::GetCommonFileName();
1641  // If the string is null or 'default' connect to standard output file
1642  if (resOut.IsNull() || resOut.EqualTo("default", TString::kIgnoreCase))
1643  resOut = AliAnalysisManager::GetCommonFileName();
1644 
1645  // Always connect input
1646  mgr->ConnectInput(this, 0, mgr->GetCommonInputContainer());
1647 
1648  // Connect sum list unless the output 'none' is specified
1649  if (!sumOut.EqualTo("none", TString::kIgnoreCase)) {
1650  TString sumName(Form("%sSums", GetName()));
1651  AliAnalysisDataContainer* sumCon =
1652  mgr->CreateContainer(sumName, TList::Class(),
1653  AliAnalysisManager::kOutputContainer, sumOut);
1654  mgr->ConnectOutput(this, 1, sumCon);
1655  }
1656 
1657  // Connect the result list unless the output 'none' is specified
1658  if (!resOut.EqualTo("none", TString::kIgnoreCase)) {
1659  TString resName(Form("%sResults", GetName()));
1660  AliAnalysisDataContainer* resCon =
1661  mgr->CreateContainer(resName, TList::Class(),
1662  AliAnalysisManager::kParamContainer, resOut);
1663  mgr->ConnectOutput(this, 2, resCon);
1664  }
1665  return true;
1666 }
1667 //____________________________________________________________________
1669 {
1670  Double_t shiftedDPhiCut = fShiftedDPhiCut;
1671  if (shiftedDPhiCut < 0) shiftedDPhiCut = TMath::Sqrt(fDeltaCut)*0.06;
1672  Double_t tailMax = fTailMax;
1673  if (tailMax < 0) tailMax = fMaxDelta;
1674 
1675  Printf("%s: %s", ClassName(), GetName());
1676  Printf(" %22s: 0x%08x", "Off-line trigger mask", fOfflineTriggerMask);
1677  Printf(" %22s: %f", "Delta phi shift", fDPhiShift);
1678  Printf(" %22s: %f", "Shifted Delta phi cut", shiftedDPhiCut);
1679  Printf(" %22s: %f", "Delta cut", fDeltaCut);
1680  Printf(" %22s: %f", "max Delta", fMaxDelta);
1681  Printf(" %22s: %f", "tail Delta", fTailDelta);
1682  Printf(" %22s: %f", "tail maximum", tailMax);
1683  Printf(" %22s: %f%%", "Absolute least c", fAbsMinCent);
1684  Printf(" %22s: %s (%d)", "Centrality estimator", fCentMethod.Data(),fCentIdx);
1687  PrintAxis(fIPzAxis,1,"IPz");
1688  PrintAxis(fCentAxis,0);
1689 
1690  if (!fCentBins) return;
1691 
1692  Printf("--- Centrality bins");
1693  TIter next(fCentBins);
1694  CentBin* bin = 0;
1695  while ((bin = static_cast<CentBin*>(next()))) {
1696  bin->Print(option);
1697  }
1698 }
1699 //____________________________________________________________________
1701 {
1703  if (!fWeights) return;
1704  fWeights->Print(option);
1705 }
1706 
1707 //____________________________________________________________________
1709 {
1710  Printf(" Centrality bin: %s", fName.Data());
1711  Printf(" Low cut: %5.1f", fLow);
1712  Printf(" High cut: %5.1f", fHigh);
1713 
1714  if (!fSubs) return;
1715  Printf(" --- Histogram sets");
1716  TIter next(fSubs);
1717  Histos* h = 0;
1718  while ((h = static_cast<Histos*>(next()))) {
1719  h->Print(option);
1720  }
1721 }
1722 
1723 //____________________________________________________________________
1725 {
1726  Sub::SetDebug(lvl);
1727  TIter next(fSubs);
1728  Histos* h = 0;
1729  while ((h = static_cast<Histos*>(next()))) {
1730  h->SetDebug(lvl);
1731  }
1732 }
1733 namespace {
1734  void Bits2String(UChar_t m, char out[7])
1735  {
1736  if (m & AliAODTracklet::kInjection) out[0] = 'I'; else out[0] = '-';
1737  if (m & AliAODTracklet::kCombinatorics) out[1] = 'C'; else out[1] = '-';
1738  if (m & AliAODTracklet::kSecondary) out[2] = 'S'; else out[2] = '-';
1739  if (m & AliAODTracklet::kDistinct) out[3] = 'D'; else out[3] = '-';
1740  if (m & AliAODTracklet::kSimulated) out[4] = 'X'; else out[4] = '-';
1741  if (m & AliAODTracklet::kGenerated) out[5] = 'G'; else out[5] = '-';
1742  out[6] = '\0';
1743  }
1744 }
1745 
1746 //____________________________________________________________________
1748 {
1749  char cMask[7]; Bits2String(fMask, cMask);
1750  char cVeto[7]; Bits2String(fVeto, cVeto);
1751  Printf(" Histograms: %s", fName.Data());
1752  Printf(" Mask: 0x%02x (%s)", fMask, cMask);
1753  Printf(" Veto: 0x%02x (%s)", fVeto, cVeto);
1754  Printf(" Delta: %s", fEtaDeltaIPz ? "yes" : "no");
1755  Printf(" Delta per PDG: %s", fEtaDeltaPdg ? "yes" : "no");
1756 }
1757 
1758 
1759 //____________________________________________________________________
1760 void
1762 {
1763  if (!WorkerInit()) {
1764  AliWarning("Failed to initialize on worker");
1765  return;
1766  }
1767  PostData(1,fContainer);
1768 }
1769 
1770 //____________________________________________________________________
1772 {
1773  if (DebugLevel() > 1) Printf("Initialising on worker");
1774  if (fShiftedDPhiCut < 0) fShiftedDPhiCut = TMath::Sqrt(fDeltaCut)*0.06;
1775  if (fTailMax < 0) fTailMax = fMaxDelta;
1776 
1777  fContainer = new Container;
1778  fContainer->SetName(Form("%sSums", GetName()));
1779  fContainer->SetOwner();
1780 
1781  Bool_t save = TH1::GetDefaultSumw2();
1782  TH1::SetDefaultSumw2();
1783  fIPz = Make1D(fContainer, "ipz", "", kMagenta+2, 20, fIPzAxis);
1784  fCent = Make1D(fContainer, "cent", "", kMagenta+2, 20, fCentAxis);
1785  fEtaPhi = Make2D(fContainer, "etaPhi","",kMagenta+2, 20, fEtaAxis,fPhiAxis);
1786 
1787  TAxis trackletAxis(1000, 0, 10000);
1788  FixAxis(trackletAxis, "#it{N}_{tracklets}");
1789  fNBareVsGood = Make1P(fContainer, "nBareVsGood", "Good",
1790  kGreen+2, 20, trackletAxis);
1791  fNBareVsFake = Make1P(fContainer, "nBareVsFake", "Fake",
1792  kRed+2, 20, trackletAxis);
1793  fNBareVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1794  fNBareVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1795  fNBareVsGood->SetStats(0);
1796  fNBareVsFake->SetStats(0);
1797 
1798  fNTrackletVsGood = Make1P(fContainer, "nTrackletVsGood", "Good",
1799  kGreen+2, 20, trackletAxis);
1800  fNTrackletVsFake = Make1P(fContainer, "nTrackletVsFake", "Fake",
1801  kRed+2, 20, trackletAxis);
1802  fNTrackletVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1803  fNTrackletVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1804  fNTrackletVsGood->SetStats(0);
1805  fNTrackletVsFake->SetStats(0);
1806 
1807  TAxis generatedAxis(1000, 0, 15000);
1808  FixAxis(generatedAxis, "#it{N}_{generated,|#eta|<2}");
1809  fNGeneratedVsGood = Make1P(fContainer, "nGeneratedVsGood", "Good",
1810  kGreen+2, 24, generatedAxis);
1811  fNGeneratedVsFake = Make1P(fContainer, "nGeneratedVsFake", "Fake",
1812  kRed+2, 24, generatedAxis);
1813  fNGeneratedVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1814  fNGeneratedVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1815  fNGeneratedVsGood->SetStats(0);
1816  fNGeneratedVsFake->SetStats(0);
1817 
1818  fCentTracklets = new TProfile("centTracklets",
1819  "Mean number of tracklets per centrality",
1820  1050, 0, 105);
1821  fCentTracklets->SetXTitle("Centrality [%]");
1822  fCentTracklets->SetYTitle("#LTtracklets#GT");
1823  fCentTracklets->SetDirectory(0);
1824  fCentTracklets->SetMarkerStyle(20);
1825  fCentTracklets->SetMarkerColor(kMagenta+2);
1826  fCentTracklets->SetLineColor(kMagenta+2);
1827  fContainer->Add(fCentTracklets);
1828  fCentEst = Make1P(fContainer,"centEstimator","",kMagenta+2, 20, fCentAxis);
1829  fCentEst->SetYTitle(Form("#LT%s%GT",fCentMethod.Data()));
1830  TH1::SetDefaultSumw2(save);
1831 
1833 
1834  typedef TParameter<double> DP;
1835  typedef TParameter<bool> BP;
1836  typedef TParameter<int> IP;
1837  Container* params = new Container;
1838  params->SetName("parameters");
1839  params->SetOwner();
1840  fContainer->Add(params);
1841  params->Add(new DP("DPhiShift", fDPhiShift, 'f'));
1842  params->Add(new DP("ShiftedDPhiCut", fShiftedDPhiCut, 'f'));
1843  params->Add(new DP("DeltaCut", fDeltaCut, 'f'));
1844  params->Add(new DP("MaxDelta", fMaxDelta, 'f'));
1845  params->Add(new DP("TailDelta", fTailDelta, 'f'));
1846  params->Add(new DP("TailMax", fTailMax, 'f'));
1847  params->Add(new DP("AbsMinCent", fAbsMinCent, 'f'));
1848  // Create our centrality bins
1849  if (!InitCentBins(0)) {
1850  AliWarning("Failed to initialize centrality bins");
1851  return false;
1852  }
1853 
1854  // Print information to log
1855  Print();
1856  return true;
1857 }
1858 //____________________________________________________________________
1860 {
1861  if (!AliTrackletAODMCdNdeta::WorkerInit()) return false;
1862  fEtaWeight = 0;
1863  if (!fWeights) {
1864  AliFatal("No weights set!");
1865  return false;
1866  }
1867 
1868  TAxis wAxis(100,0,10);
1869  FixAxis(wAxis,"#it{w}_{i}");
1870  fEtaWeight = Make2P(fContainer, "etaWeight", "#LTw#GT", kYellow+2, 24,
1871  fEtaAxis, fCentAxis);
1872  fWeightCorr = Make2D(fContainer, "weightCorr", "w_{1} vs w_{2}",
1873  kCyan+2, 24, wAxis, wAxis);
1874  fWeights->Store(fContainer);
1875  return true;
1876 }
1877 //____________________________________________________________________
1880  Float_t c1,
1881  Float_t c2,
1882  Container* existing,
1883  const TAxis& deltaAxis)
1884 {
1885  CentBin* ret = MakeCentBin(c1, c2);
1886  Bool_t ok = true;
1887  if (!existing)
1888  ok = ret->WorkerInit(fContainer,fEtaAxis,fIPzAxis,deltaAxis);
1889  else
1890  ok = ret->FinalizeInit(existing);
1891  if (!ok) {
1892  AliWarningF("Failed to initialize bin %s", ret->GetName());
1893  delete ret;
1894  return 0;
1895  }
1896  ret->SetDebug(DebugLevel());
1897  fCentBins->AddAt(ret, bin);
1898  return ret;
1899 
1900 }
1901 //____________________________________________________________________
1903 {
1904  if (DebugLevel() > 1)
1905  Printf("Initialising on centrality bins on %s",
1906  existing ? "master" : "worker");
1907  if (fCentBins) return true;
1908  fCentBins = new Container;
1909  fCentBins->SetName("centralityBins");
1910  fCentBins->SetOwner();
1911 
1912  TAxis deltaAxis (Int_t(5*fMaxDelta),0,fMaxDelta);
1913  FixAxis(deltaAxis,
1914  "#Delta=[(#Delta#phi-#delta#phi)/#sigma_{#phi}]^{2}+"
1915  "[#Delta#thetasin^{-2}(#theta)/#sigma_{#theta}]^{2}");
1916 
1917  // Add min-bias bin
1918  if (!InitCentBin(0,0,0,existing,deltaAxis)) return false;
1919 
1920  // Add other bins
1921  Int_t nCentBins = fCentAxis.GetNbins();
1922  for (Int_t i = 1; i <= nCentBins; i++) {
1923  Float_t c1 = fCentAxis.GetBinLowEdge(i);
1924  Float_t c2 = fCentAxis.GetBinUpEdge(i);
1925  if (!InitCentBin(i, c1, c2, existing, deltaAxis)) return false;
1926  }
1927  if (!InitCentBin(nCentBins+1, 0, 100, existing, deltaAxis)) return false;
1928  return true;
1929 }
1930 
1931 //____________________________________________________________________
1933  : Sub(""),
1934  fSubs(0),
1935  fLow(c1),
1936  fHigh(c2),
1937  fStatus(0),
1938  fIPz(0),
1939  fCent(0),
1940  fCentIPz(0),
1941  fMeasured(0),
1942  fInjection(0)
1943 {
1944  if (c1 >= c2)
1945  fName = "all";
1946  else
1948  fMeasured = MakeHistos("measured", kRed+2, 20,
1949  kMeasuredMask, // No requirements, just veto
1950  kMeasuredVeto);
1951  fInjection = MakeHistos("injected", kOrange+2, 21,
1952  kInjectedMask,
1953  kInjectedVeto);
1954  fSubs = new Container;
1955  fSubs->SetOwner(true);
1956  fSubs->Add(fMeasured);
1957  fSubs->Add(fInjection);
1958 }
1959 
1960 //____________________________________________________________________
1962  : AliTrackletAODdNdeta::CentBin(c1, c2),
1963  fCombinatorics(0),
1964  fDistinct(0),
1965  fPrimaries(0),
1966  fSecondaries(0),
1967  fGenerated(0)
1968 {
1969  // Combinatorics is everything that has the combinatorics bit
1970  // Primaries all with simulated bit, but not secondary or
1971  // combinatorics bit.
1972  // Secondaries are all those with the secondary bit set
1973  fCombinatorics = MakeHistos("combinatorics", kMagenta+2, 30,
1976  fDistinct = MakeHistos("distinct", kCyan+2, 27,
1977  kDistinctMask,
1978  kDistinctVeto);
1979  fPrimaries = MakeHistos("primaries", kGreen+2, 26,
1980  kPrimaryMask,
1981  kPrimaryVeto);
1982  fSecondaries = MakeHistos("secondaries", kBlue+2, 32,
1984  kSecondaryVeto);
1985  fGenerated = MakeHistos("generated", kGray+1, 28,
1987  kGeneratedVeto);
1988  fMeasured->fStyle = 24;
1989  fInjection->fStyle = 25;
1990  fSubs->Add(fCombinatorics);
1991  fSubs->Add(fDistinct);
1992  fSubs->Add(fPrimaries);
1993  fSubs->Add(fSecondaries);
1994  fSubs->AddAfter(fMeasured, fGenerated);
1995 }
1996 
1997 //____________________________________________________________________
2000  Color_t color,
2001  Style_t style,
2002  UShort_t mask,
2003  UShort_t veto)
2004 {
2005  return new Histos(name, color, style, mask, veto);
2006 }
2007 
2008 //____________________________________________________________________
2010  const TAxis& etaAxis,
2011  const TAxis& ipzAxis,
2012  const TAxis& deltaAxis)
2013 {
2014  if (fDebug > 0)
2015  Printf("Initializing centrality bin %s", fName.Data());
2016  if (!Sub::WorkerInit(parent, etaAxis, ipzAxis, deltaAxis)) return false;
2017 
2018  TAxis centAxis(20, fLow, fHigh);
2019  FixAxis(centAxis, "Centrality [%]");
2021  fCent = Make1D(fContainer,"cent","Centrality [%]",
2022  kMagenta+2,20,centAxis);
2023  fIPz = Make1D(fContainer,"ipz","IP_{#it{z}} [cm]",kRed+2,20,ipzAxis);
2024  fCentIPz = Make1P(fContainer,"centIpz","#LTc#GT vs IP_{#it{z}}",kPink+2,
2025  20,ipzAxis);
2026  fCentIPz->SetYTitle("#LTc#GT");
2027  TIter next(fSubs);
2028  Histos* h = 0;
2029  while ((h = static_cast<Histos*>(next()))) {
2030  if (!h ->WorkerInit(fContainer, etaAxis,ipzAxis,deltaAxis))
2031  return false;
2032  }
2033 
2034  return true;
2035 }
2036 //____________________________________________________________________
2038  const TAxis& etaAxis,
2039  const TAxis& ipzAxis,
2040  const TAxis& deltaAxis)
2041 {
2042  if (!Sub::WorkerInit(parent, etaAxis, ipzAxis, deltaAxis)) return false;
2043  TString shrt(Form("%c", GetName()[0]));
2044  shrt.ToUpper();
2045  if (GetC(parent, "generated", false) != 0) shrt.Append("'");
2046 
2047  // Do not make eta vs IPz for secondaries and primaries
2048  if (!IsPrimary() && !IsSecondary())
2049  fEtaIPz = Make2D(fContainer, "etaIPz", shrt.Data(),
2050  kRed+2, 20, etaAxis, ipzAxis);
2051  // Always make eta vs Delta distribution, except for MC truth
2052  if (!IsGenerated())
2053  fEtaDeltaIPz = Make3D(fContainer, "etaDeltaIPz",
2054  Form("#Delta_{%s}",shrt.Data()),
2055  kBlue+2, 21, etaAxis, deltaAxis, ipzAxis);
2056  if (!IsGenerated() &&
2057  // !IsPrimary() &&
2058  // !IsSecondary() &&
2059  !IsInjected() &&
2060  fStyle != 20) // Last condition to not make this for real data
2061  fEtaPdgIPz = Make3D(fContainer, "etaPdgIPz",
2062  "Parent particle type",
2063  kGreen+2, 22, etaAxis, PdgAxis(), ipzAxis);
2064 
2065  if (IsGenerated()) {
2066  fEtaPdg = Make2D(fContainer, "etaPdg",
2067  "Primary particle type",
2068  kYellow+2, 30, etaAxis, PdgAxis());
2069  TAxis ptAxis(100, 0, 5);
2070  ptAxis.SetTitle("#it{p}_{T}");
2071  FixAxis(ptAxis);
2072  fEtaPt = Make2D(fContainer, "etaPt",
2073  "Primary transverse momentum",
2074  kCyan+2, 30, etaAxis, ptAxis);
2075  }
2076  if (IsPrimary() || IsSecondary() || IsCombinatoric()) {
2077  fEtaDeltaPdg = Make3D(fContainer, "etaDeltaPdg",
2078  "#Delta by primary particle type",
2079  kMagenta, 22, etaAxis, deltaAxis, PdgAxis());
2080  }
2081 
2082  SetAttr(fColor, fStyle);
2083  return true;
2084 }
2085 //____________________________________________________________________
2087 {
2088  if (fEtaIPz) {
2089  fEtaIPz->SetMarkerStyle(m);
2090  fEtaIPz->SetMarkerColor(c);
2091  fEtaIPz->SetLineColor(c);
2092  fEtaIPz->SetFillColor(c);
2093  }
2094  if (fEtaDeltaIPz) {
2095  fEtaDeltaIPz->SetMarkerStyle(m);
2096  fEtaDeltaIPz->SetMarkerColor(c);
2097  fEtaDeltaIPz->SetLineColor(c);
2098  fEtaDeltaIPz->SetFillColor(c);
2099  }
2100  if (fEtaDeltaPdg) {
2101  fEtaDeltaPdg->SetMarkerStyle(m);
2102  fEtaDeltaPdg->SetMarkerColor(c);
2103  fEtaDeltaPdg->SetLineColor(c);
2104  fEtaDeltaPdg->SetFillColor(c);
2105  }
2106  if (fEtaPdgIPz) {
2107  fEtaPdgIPz->SetMarkerStyle(m);
2108  fEtaPdgIPz->SetMarkerColor(c);
2109  fEtaPdgIPz->SetLineColor(c);
2110  fEtaPdgIPz->SetFillColor(c);
2111  }
2112  if (fEtaPdg) {
2113  fEtaPdg->SetMarkerStyle(m);
2114  fEtaPdg->SetMarkerColor(c);
2115  fEtaPdg->SetLineColor(c);
2116  fEtaPdg->SetFillColor(c);
2117  }
2118  if (fEtaPt) {
2119  fEtaPt->SetMarkerStyle(m);
2120  fEtaPt->SetMarkerColor(c);
2121  fEtaPt->SetLineColor(c);
2122  fEtaPt->SetFillColor(c);
2123  }
2124 }
2125 
2126 //____________________________________________________________________
2127 void
2129 {
2130  if (DebugLevel() > 0) Printf("In user exec");
2131  Double_t cent = -1;
2132  const AliVVertex* ip = 0;
2133  TClonesArray* tracklets = 0;
2134  UInt_t status = CheckEvent(cent, ip, tracklets);
2135  if ((status & kOKEvent) != kOKEvent) {
2136  AliWarningF("Event didn't pass %f, %p, %p", cent, ip, tracklets);
2137  return;
2138  }
2139  if (DebugLevel() > 0) Printf("Got centrality=%f ipZ=%f %d tracklets",
2140  cent, ip->GetZ(), tracklets->GetEntriesFast());
2141  ProcessEvent(status, cent, ip, tracklets);
2142 
2143  PostData(1,fContainer);
2144  fStatus->Fill(kCompleted);
2145 }
2146 
2147 //____________________________________________________________________
2149 {
2150  return (1 << (bin-1));
2151 }
2152 
2153 //____________________________________________________________________
2155 {
2156  TH1* status = new TH1F("status", "Status of task",
2157  kCompleted, .5, kCompleted+.5);
2158  status->SetMarkerSize(2);
2159  status->SetMarkerColor(kMagenta+2);
2160  status->SetLineColor(kMagenta+2);
2161  status->SetFillColor(kMagenta+2);
2162  status->SetFillStyle(1001);
2163  status->SetBarOffset(0.1);
2164  status->SetBarWidth(0.4);
2165  status->SetDirectory(0);
2166  status->SetStats(0);
2167  status->SetXTitle("Event have");
2168  status->SetYTitle("# Events");
2169  status->GetXaxis()->SetBinLabel(kAll, "Been seen");
2170  status->GetXaxis()->SetBinLabel(kEvent, "Event data");
2171  status->GetXaxis()->SetBinLabel(kTracklets, "Tracklets");
2172  status->GetXaxis()->SetBinLabel(kTrigger, "Trigger");
2173  status->GetXaxis()->SetBinLabel(kIP, "IP");
2174  status->GetXaxis()->SetBinLabel(kCentrality, "Centrality");
2175  status->GetXaxis()->SetBinLabel(kCompleted, "Completed");
2176  container->Add(status);
2177  return status;
2178 }
2179 
2180 //____________________________________________________________________
2182  const AliVVertex*& ip,
2183  TClonesArray*& tracklets)
2184 {
2185  UInt_t ret = Bin2Flag(kAll);
2186  // Count all events
2187  fStatus->Fill(kAll);
2188 
2189  // Check for event
2190  AliVEvent* event = InputEvent();
2191  if (!event) {
2192  AliWarning("No event");
2193  return ret;
2194  }
2195  ret |= Bin2Flag(kEvent);
2196  fStatus->Fill(kEvent);
2197 
2198  // Check if we have the tracklets
2199  tracklets = FindTracklets(event);
2200  if (!tracklets) return ret;
2201  ret |= Bin2Flag(kTracklets);
2202  fStatus->Fill(kTracklets);
2203 
2204  // Check if event was triggered
2205  Bool_t trg = FindTrigger();
2206  if (!trg) return ret;
2207  ret |= Bin2Flag(kTrigger);
2208  fStatus->Fill(kTrigger);
2209 
2210  // Check the interaction point
2211  ip = FindIP(event);
2212  if (!ip) return ret;
2213  ret |= Bin2Flag(kIP);
2214  fStatus->Fill(kIP);
2215 
2216  // Check the centrality
2217  Int_t nTracklets = 0;
2218  cent = FindCentrality(event, nTracklets);
2219  // Do not fail on missing centrality
2220  if (cent >= 0) {
2221  ret |= Bin2Flag(kCentrality);
2222  fStatus->Fill(kCentrality);
2223  }
2224 
2225  fIPz->Fill(ip->GetZ());
2226  fCent->Fill(cent);
2227  fCentTracklets->Fill(cent, nTracklets);
2228 
2229  return ret;
2230 }
2231 
2232 //____________________________________________________________________
2233 TClonesArray* AliTrackletAODdNdeta::FindTracklets(AliVEvent* event)
2234 {
2235  // Check the multiplicity
2236  TObject* obj = event->FindListObject(GetBranchName());
2237  if (!obj) {
2238  AliWarningF("Couldn't get object %s", GetBranchName());
2239  // event->GetList()->Print();
2240  return 0;
2241  }
2242  if (!obj->IsA()->InheritsFrom(TClonesArray::Class())) {
2243  AliWarningF("Object %s is not a TClonesArray but a %s",
2244  obj->GetName(), obj->ClassName());
2245  return 0;
2246  }
2247  return static_cast<TClonesArray*>(obj);
2248 }
2249 
2250 //____________________________________________________________________
2251 const AliVVertex* AliTrackletAODdNdeta::CheckIP(const AliVVertex* ip,
2252  Double_t maxDispersion,
2253  Double_t maxZError)
2254 {
2255  if (!ip) return 0;
2256  if (ip->GetNContributors() <= 0) {
2257  AliWarning("Not enough contributors for IP");
2258  return 0;
2259  }
2260  // If this is from the Z vertexer, do some checks
2261  if (ip->IsFromVertexerZ()) {
2262  // Get covariance matrix
2263  Double_t covar[6];
2264  ip->GetCovarianceMatrix(covar);
2265  Double_t sigmaZ = TMath::Sqrt(covar[5]);
2266  if (sigmaZ >= maxZError) {
2267  AliWarningF("IPz resolution = %f >= %f", sigmaZ, maxZError);
2268  return 0;
2269  }
2270 
2271  // If this IP doesn not derive from AliVertex, don't check dispersion.
2272  if (ip->IsA()->InheritsFrom(AliVertex::Class())) {
2273  const AliVertex* ipv = static_cast<const AliVertex*>(ip);
2274  // Dispersion is the parameter used by the vertexer for finding the IP.
2275  if (ipv->GetDispersion() >= maxDispersion) {
2276  AliWarningF("IP dispersion = %f >= %f",
2277  ipv->GetDispersion(), maxDispersion);
2278  return 0;
2279  }
2280  }
2281  }
2282 
2283  // If we get here, we either have a full 3D vertex or track
2284  // vertex, and we should check if it is in range
2285  if (ip->GetZ() < fIPzAxis.GetXmin() || ip->GetZ() > fIPzAxis.GetXmax()) {
2286  AliWarningF("IPz = %fcm out of range [%f,%f]cm",
2287  ip->GetZ(), fIPzAxis.GetXmin(), fIPzAxis.GetXmax());
2288  return 0;
2289  }
2290  return ip;
2291 }
2292 
2293 //____________________________________________________________________
2294 const AliVVertex* AliTrackletAODdNdeta::FindRealIP(AliVEvent* event,
2295  Double_t maxDispersion,
2296  Double_t maxZError)
2297 {
2298  const AliVVertex* ip = event->GetPrimaryVertex();
2299  if (!ip) {
2300  if (DebugLevel() > 1) AliWarning("No real IP for this event found!");
2301  return 0;
2302  }
2303  return CheckIP(ip, maxDispersion, maxZError);
2304 }
2305 
2306 //____________________________________________________________________
2307 const AliVVertex* AliTrackletAODdNdeta::FindSimpleIP(AliVEvent* event,
2308  Double_t maxDispersion,
2309  Double_t maxZError)
2310 {
2311  static AliVertex* ip;
2312  AliAODSimpleHeader* head =
2313  static_cast<AliAODSimpleHeader*>(event
2314  ->FindListObject("AliAODSimpleHeader"));
2315  if (!head) {
2316  if (DebugLevel() > 1) AliWarning("No simple header");
2317  return 0;
2318  }
2319  if (!ip)
2320  ip = new AliVertex;
2321  ip->SetXv(head->fRecIP.X());
2322  ip->SetYv(head->fRecIP.Y());
2323  ip->SetZv(head->fRecIP.Z());
2324  ip->SetNContributors(10);
2325  ip->SetTitle("");
2326  return CheckIP(ip, maxDispersion, maxZError);
2327 }
2328 
2329 //____________________________________________________________________
2330 const AliVVertex* AliTrackletAODdNdeta::FindIP(AliVEvent* event,
2331  Double_t maxDispersion,
2332  Double_t maxZError)
2333 {
2334  const AliVVertex* ip = FindRealIP(event,maxDispersion,maxZError);
2335  if (!ip) {
2336  ip = FindSimpleIP(event,maxDispersion,maxZError);
2337  if (!ip) {
2338  if (DebugLevel() > 1) AliWarning("No IP for this event found!");
2339  return 0;
2340  }
2341  }
2342  // Good vertex, return it
2343  return ip;
2344 }
2345 //____________________________________________________________________
2347 {
2348  UInt_t evBits = fInputHandler->IsEventSelected();
2349  Bool_t trgOK = (evBits & fOfflineTriggerMask || fOfflineTriggerMask == 0);
2350  if (DebugLevel())
2351  Printf("Trigger bits=0x%08x mask=0x%08x masked=0x%08x -> %s",
2352  evBits, fOfflineTriggerMask, evBits & fOfflineTriggerMask,
2353  trgOK ? "selected" : "rejected");
2354  return trgOK;
2355 }
2356 //____________________________________________________________________
2358 {
2359  static TString centMeth;
2360  if (centMeth.IsNull()) {
2361  centMeth = fCentMethod(3,fCentMethod.Length()-3);
2362  }
2363  AliCentrality* cent = event->GetCentrality();
2364  if (!cent) return -1;
2365 
2366  Double_t centPer = cent->GetCentralityPercentileUnchecked(centMeth);
2367  return centPer;
2368 }
2369 //____________________________________________________________________
2371  Int_t& nTracklets)
2372 {
2373  AliMultSelection* cent =
2374  static_cast<AliMultSelection*>(event->FindListObject("MultSelection"));
2375  if (!cent) {
2376  AliWarning("No centrality in event");
2377  event->GetList()->Print();
2378  return -1;
2379  }
2380  if (fCentIdx < 0) {
2381  TString trkName("SPDTracklets");
2382  for (Int_t i = 0; i < cent->GetNEstimators(); i++) {
2383  AliMultEstimator* e = cent->GetEstimator(i);
2384  if (!e) continue;
2385  if (fCentMethod.EqualTo(e->GetName())) fCentIdx = i;
2386  if (trkName.EqualTo(e->GetName())) fTrkIdx = i;
2387  }
2388  }
2389  if (fCentIdx < 0) {
2390  AliWarningF("Centrality estimator %s not found", fCentMethod.Data());
2391  return -1;
2392  }
2393 
2394  // Use cached index to look up the estimator
2395  AliMultEstimator* estCent = cent->GetEstimator(fCentIdx);
2396  if (!estCent) {
2397  AliWarningF("Centrality estimator %s not available for this event",
2398  fCentMethod.Data());
2399  return -1;
2400  }
2401  fCentEst->Fill(estCent->GetPercentile(), estCent->GetValue());
2402 
2403  // Now look up tracklet value using cached index
2404  if (fTrkIdx < 0) return estCent->GetPercentile();
2405  AliMultEstimator* estTracklets = cent->GetEstimator(fTrkIdx);
2406  if (estTracklets) nTracklets = estTracklets->GetValue();
2407 
2408  return estCent->GetPercentile();
2409 }
2410 
2411 //____________________________________________________________________
2413  Int_t& nTracklets)
2414 {
2415  if (fCentMethod.EqualTo("MB", TString::kIgnoreCase)) {
2416  Printf("MB centrality - not checked");
2417  return 0;
2418  }
2419  Double_t centPer = -1;
2420  if (fCentMethod.BeginsWith("OLD"))
2421  centPer = FindCompatCentrality(event);
2422  else
2423  centPer = FindMultCentrality(event, nTracklets);
2424  const Double_t safety = 1e-3;
2425  if (DebugLevel() > 1) Printf("Read centrality: %f%%", centPer);
2426  if (centPer < -safety) return -2;
2427  if (centPer < +safety) centPer = safety;
2428  else if (centPer > 100-safety) centPer = 100-safety;
2429 
2430  if (fAbsMinCent >= 0 && centPer < fAbsMinCent) {
2431  AliWarningF("Centrality = %f lower than absolute minimum [%f]",
2432  centPer, fAbsMinCent);
2433  return -3;
2434  }
2435  if (centPer < fCentAxis.GetXmin() || centPer > fCentAxis.GetXmax()) {
2436  AliWarningF("Centrality = %f out of range [%f,%f]",
2437  centPer, fCentAxis.GetXmin(), fCentAxis.GetXmax());
2438  return -3;
2439  }
2440  return centPer;
2441 }
2442 
2443 //____________________________________________________________________
2445  Double_t cent,
2446  Double_t ipz)
2447 {
2448  return 1;
2449 }
2450 //____________________________________________________________________
2452  Double_t cent,
2453  Double_t ipz)
2454 {
2455  // We don't check for weights, as we must have them to come this far
2456  // if (!fWeights) {
2457  // AliWarning("No weights defined");
2458  // return 1;
2459  // }
2460  Double_t w = fWeights->LookupWeight(tracklet, cent, ipz, fWeightCorr);
2461  if (tracklet->IsMeasured())
2462  fEtaWeight->Fill(tracklet->GetEta(), cent, w);
2463 
2464  // printf("Looking up weight of tracklet -> %f ", w);
2465  // tracklet->Print();
2466  return w;
2467 }
2468 //____________________________________________________________________
2470 {
2471  Double_t dPhiS = (tracklet->GetDPhi() -
2472  TMath::Sign(fDPhiShift,Double_t(tracklet->GetDPhi())));
2473  UShort_t ret = 0;
2474  ret |= ((tracklet->GetDelta() <= fDeltaCut) ? 0x1 : 0);
2475  ret |= ((dPhiS < fShiftedDPhiCut) ? 0x2 : 0);
2476  // Printf("dPhiS=%f (%f) Delta=%f (%f) -> 0x%x",
2477  // dPhiS, fShiftedDPhiCut,
2478  // tracklet->GetDelta(), fDeltaCut,
2479  // ret);
2480  return ret;
2481 }
2482 
2483 
2484 //____________________________________________________________________
2486  Double_t cent,
2487  const AliVVertex* ip,
2488  TClonesArray* tracklets)
2489 {
2490  // Figure out which centrality bins to fill
2491  Double_t ipz = (ip ? ip->GetZ() : -1000);
2492  Int_t nAcc = 0;
2493  TIter nextAcc(fCentBins);
2494  CentBin* bin = 0;
2495  TList toRun;
2496  while ((bin = static_cast<CentBin*>(nextAcc()))) {
2497  // Not in range for this bin
2498  if (!bin->Accept(status, cent, ipz)) continue;
2499  toRun.Add(bin);
2500  nAcc++;
2501  }
2502  // If we have no centrality bins to fill, we return immediately
2503  if (nAcc <= 0) return;
2504 
2505  Double_t nBare = 0;
2506  Double_t nMeasured = 0;
2507  Double_t nGenerated = 0;
2508  Double_t nGood = 0;
2509  Double_t nFake = 0;
2510  AliAODTracklet* tracklet = 0;
2511  TIter nextTracklet(tracklets);
2512  while ((tracklet = static_cast<AliAODTracklet*>(nextTracklet()))) {
2513  Double_t weight = LookupWeight(tracklet, cent, ipz);
2514  UShort_t signal = CheckTracklet(tracklet);
2515  if (signal) fEtaPhi->Fill(tracklet->GetEta(), tracklet->GetPhi());
2516  if (tracklet->IsMeasured() && TMath::Abs(tracklet->GetEta()) < 0.7) {
2517  nMeasured += weight;
2518  nBare ++;
2519  if (tracklet->IsCombinatorics()) nFake += weight;
2520  else nGood += weight;
2521  }
2522  else if (tracklet->IsGenerated() &&
2523  !tracklet->IsNeutral() &&
2524  TMath::Abs(tracklet->GetEta()) <= 1) nGenerated ++; // += weight;
2525  TIter nextBin(&toRun);
2526  while ((bin = static_cast<CentBin*>(nextBin()))) {
2527  bin->ProcessTracklet(tracklet, ip->GetZ(), signal, weight);
2528  }
2529  }
2530  TIter nextBin(&toRun);
2531  while ((bin = static_cast<CentBin*>(nextBin()))) {
2532  bin->Completed();
2533  }
2534  fNBareVsFake ->Fill(nBare, nFake);
2535  fNBareVsGood ->Fill(nBare, nGood);
2536  fNTrackletVsFake ->Fill(nMeasured, nFake);
2537  fNTrackletVsGood ->Fill(nMeasured, nGood);
2538  fNGeneratedVsFake->Fill(nGenerated, nFake);
2539  fNGeneratedVsGood->Fill(nGenerated, nGood);
2540 }
2541 
2542 //____________________________________________________________________
2544 {
2545  return fLow >= fHigh;
2546 }
2547 //____________________________________________________________________
2549  Double_t cent,
2550  Double_t ipz)
2551 {
2552  Bool_t centOK = (IsAllBin() ||
2553  (status & Bin2Flag(kCentrality) &&
2554  cent >= fLow && cent < fHigh));
2555  // We get out here already as we're not in this centrality bin
2556  if (!centOK) return false;
2557  fStatus->Fill(kAll);
2558  if (status & Bin2Flag(kEvent)) fStatus->Fill(kEvent);
2559  if (status & Bin2Flag(kTracklets)) fStatus->Fill(kTracklets);
2560  if (status & Bin2Flag(kTrigger))
2561  fStatus->Fill(kTrigger);
2562  else
2563  return false; // In case we have no trigger
2564  if (status & Bin2Flag(kCentrality)) {
2565  fStatus->Fill(kCentrality);
2566  fCent->Fill(cent);
2567  }
2568  // We explicity check the IP status here, since we need that to
2569  // calculate the per-bin vertex efficiency
2570  if (status & Bin2Flag(kIP)) {
2571  fStatus ->Fill(kIP);
2572  fIPz ->Fill(ipz);
2573  fCentIPz->Fill(ipz, cent);
2574  return true;
2575  }
2576  return false;
2577 }
2578 
2579 //____________________________________________________________________
2581  Double_t ipZ,
2582  UShort_t signal,
2583  Double_t weight)
2584 {
2585  TIter next(fSubs);
2586  Histos* h = 0;
2587  if (fDebug > 3) tracklet->Print();
2588  while ((h = static_cast<Histos*>(next())))
2589  h->ProcessTracklet(tracklet, ipZ, signal, weight);
2590 
2591  return true;
2592 }
2593 //____________________________________________________________________
2595  Double_t ipZ,
2596  UShort_t signal,
2597  Double_t weight)
2598 {
2599  if (!fEtaIPz && !fEtaDeltaIPz) return true;
2600  // Get tracklet info
2601  Double_t eta = tracklet->GetEta();
2602  Double_t delta = tracklet->GetDelta();
2603  UChar_t flags = tracklet->GetFlags();
2604  Int_t pdgBin = -1;
2605  Int_t pdgBin2 = -1;
2606 
2607  // For debugging
2608  char m[7];
2609  if (fDebug > 3) Bits2String(flags, m);
2610 
2611  // Check the filter mask
2612  if (fMask != 0 && (flags & fMask) == 0) {
2613  if (fDebug > 3) Printf("%14s (0x%02x,----) %6s %7s (0x%02x) ",
2614  GetName(),fMask,"reject",m,flags);
2615  return false;
2616  }
2617 
2618  // If we need the PDG code, get that here
2619  if (fEtaPdgIPz || fEtaPdg || fEtaDeltaPdg) {
2620  pdgBin = PdgBin(tracklet->GetParentPdg());
2621  if (tracklet->IsCombinatorics())
2622  pdgBin2 = PdgBin(tracklet->GetParentPdg(true));
2623  }
2624 
2625  // If we have the eta-vs-pdg histogram, we should fill before the
2626  // veto, which filters out the neutral particles.
2627  if (fEtaPdg) fEtaPdg->Fill(eta, pdgBin, weight);
2628 
2629  // Check the veto mask
2630  if (fVeto != 0 && (flags & fVeto) != 0) {
2631  if (fDebug > 3) Printf("%14s (----,0x%02x) %6s %7s (0x%02x) ",
2632  GetName(), fVeto, "veto", m, flags);
2633  return false;
2634  }
2635  // Debug message that we accepted tracklet
2636  if (fDebug > 3) Printf("%14s (0x%02x,0x%02x) %6s %7s (0x%02x) ",
2637  GetName(), fMask, fVeto, "accept", m, flags);
2638 
2639  // Fill PDG,eta dependent Delta distributions
2640  if (fEtaDeltaPdg) {
2641  fEtaDeltaPdg ->Fill(eta, delta, pdgBin, weight);
2642  // Fill both particles
2643  if (pdgBin2 >= 0) fEtaPdgIPz->Fill(eta, delta, pdgBin2, weight);
2644  }
2645 
2646  // Both reguirements (Delta < cut, dPhi < cut)
2647  if (fEtaIPz && (signal == 0x3)) fEtaIPz->Fill(eta, ipZ, weight);
2648  // Just check dPhi
2649  if (fEtaDeltaIPz && (signal & 0x2)) fEtaDeltaIPz->Fill(eta,delta,ipZ,weight);
2650 
2651  // If we do not need the eta-PDG-IPz or eta-Pt, return now
2652  if (!fEtaPdgIPz && !fEtaPt) return true;
2653 
2654  if (fEtaPdgIPz) {
2655  fEtaPdgIPz ->Fill(eta, pdgBin, ipZ, weight);
2656  // Fill both particles
2657  if (pdgBin2 >= 0) fEtaPdgIPz->Fill(eta, pdgBin2, ipZ, weight);
2658  }
2659 
2660  if (fEtaPt) fEtaPt ->Fill(eta, tracklet->GetParentPt(), weight);
2661 
2662  return true;
2663 }
2664 
2665 //____________________________________________________________________
2667 {
2668  fStatus->Fill(kCompleted);
2669 }
2670 
2671 //____________________________________________________________________
2672 void
2674 {
2675  TString resName; resName.Form("%sResults",GetName());
2676  if (GetOutputData(2)) {
2677  Warning("Terminate", "Already have a result container, making a new one");
2678  resName.Append("New");
2679  }
2680  Container* results = new Container;
2681  results->SetName(resName);
2682  results->SetOwner();
2683 
2684  Print("");
2685  fContainer = static_cast<Container*>(GetOutputData(1));
2686  if (!fContainer) {
2687  AliWarning("No sum container found!");
2688  return;
2689  }
2690 
2691  if (!InitCentBins(fContainer)) {
2692  AliWarningF("Failed to initialize centrality bins from %s",
2693  fContainer->GetName());
2694  return;
2695  }
2696 
2697  if (!MasterFinalize(results)) {
2698  AliWarning("Failed to finalize results");
2699  return;
2700  }
2701 
2702  PostData(2, results);
2703 }
2704 
2705 
2706 //____________________________________________________________________
2708 {
2709  fContainer = GetC(parent, fName);
2710  fStatus = GetH1(fContainer, "status");
2711  fCent = GetH1(fContainer, "cent");
2712  fIPz = GetH1(fContainer, "ipz");
2713  fCentIPz = GetP1(fContainer, "centIpz");
2714  if (!fContainer || !fStatus || !fCent || !fIPz) return false;
2715  TIter next(fSubs);
2716  Histos* h = 0;
2717  while ((h = static_cast<Histos*>(next())))
2718  if (!h->FinalizeInit(fContainer)) return false;
2719  return true;
2720 }
2721 //____________________________________________________________________
2723 {
2724  fContainer = GetC(parent, fName);
2725  fEtaIPz = GetH2(fContainer, "etaIPz", false); // No complaints
2726  fEtaDeltaIPz = GetH3(fContainer, "etaDeltaIPz", false);
2727  fEtaDeltaPdg = GetH3(fContainer, "etaDeltaPdg", false);
2728  fEtaPdgIPz = GetH3(fContainer, "etaPdgIPz", false);
2729  fEtaPdg = GetH2(fContainer, "etaPdg", false);
2730  fEtaPt = GetH2(fContainer, "etaPt", false);
2731  if (GetName()[0] == 'm' && GetC(parent,"generated")) {
2732  // Fix up titles
2733  if (fEtaIPz)
2734  fEtaIPz->SetTitle(Form("%s'", fEtaIPz->GetTitle()));
2735  if (fEtaDeltaIPz)
2736  fEtaDeltaIPz->SetTitle(Form("#Delta_{%s}", fEtaIPz->GetTitle()));
2737  }
2738  return (fContainer != 0); // && fEtaIPz != 0);
2739 }
2740 
2741 //____________________________________________________________________
2743 {
2744  // Make copies of histograms and store
2745  fIPz = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"ipz")));
2746  fCent = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"cent")));
2747  fStatus = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"status")));
2748  fEtaPhi = static_cast<TH2*>(CloneAndAdd(results,GetH2(fContainer,"etaPhi")));
2749  fCentTracklets =
2750  static_cast<TProfile*>(CloneAndAdd(results,GetP(fContainer,
2751  "centTracklets")));
2752  fCentEst =
2753  static_cast<TProfile*>(CloneAndAdd(results,GetP(fContainer,
2754  "centEstimator")));
2755  Double_t nEvents = fIPz->GetEntries();
2756  Printf("Event summary:");
2757  for (Int_t i = 1; i <= fStatus->GetNbinsX(); i++)
2758  Printf(" %10d %s",
2759  Int_t(fStatus->GetBinContent(i)),
2760  fStatus->GetXaxis()->GetBinLabel(i));
2761  for (Int_t i = 1; i <= fCent->GetNbinsX(); i++)
2762  Printf(" %6.2f-%6.2f%%: %d",
2763  fCent->GetXaxis()->GetBinLowEdge(i),
2764  fCent->GetXaxis()->GetBinUpEdge(i),
2765  Int_t(fCent->GetBinContent(i)));
2766 
2767 
2768  fIPz ->Scale(1./nEvents);
2769  fCent ->Scale(1./fCent->GetEntries());
2770  fStatus->Scale(1./fStatus->GetBinContent(1));
2771  fEtaPhi->Scale(1./nEvents);
2772 
2773  if (fTailMax < 0) fTailMax = fMaxDelta;
2774  TIter next(fCentBins);
2775  CentBin* bin = 0;
2776  while ((bin = static_cast<CentBin*>(next()))) {
2777  if (!bin->MasterFinalize(results, 0, fTailDelta, fTailMax)) {
2778  AliWarningF("Failed to finalize %s", bin->GetName());
2779  return false;
2780  }
2781  }
2782  return true;
2783 }
2784 //____________________________________________________________________
2786 {
2787  if (!AliTrackletAODMCdNdeta::MasterFinalize(results)) return false;
2788 
2789  TObject* o = fContainer->FindObject("etaWeight");
2790  if (o && o->IsA()->InheritsFrom(TH1::Class())) {
2791  TH1* h = static_cast<TH1*>(o->Clone());
2792  h->SetDirectory(0);
2793  results->Add(h);
2794  }
2795  else {
2796  AliWarningF("Object %p (etaWeight) is not a TH1 or not found",o);
2797  }
2798  o = fContainer->FindObject("weightCorr");
2799  if (o && o->IsA()->InheritsFrom(TH1::Class())) {
2800  TH1* h = static_cast<TH1*>(o->Clone());
2801  h->SetDirectory(0);
2802  results->Add(h);
2803  }
2804  else {
2805  AliWarningF("Object %p (weightCorr) is not a TH1 or not found",o);
2806  }
2807  if (!fWeights) return true;
2808 
2809  if (!fWeights->Retrieve(fContainer)) return false;
2810  fWeights->Store(results);
2811 
2812  return true;
2813 }
2814 
2815 
2816 //____________________________________________________________________
2818  TH1* ,
2819  Double_t tailDelta,
2820  Double_t tailMax)
2821 {
2822  Container* result = new Container;
2823  result->SetName(fName);
2824  result->SetOwner(true);
2825  parent->Add(result);
2826 
2827  TH1* status = static_cast<TH1*>(CloneAndAdd(result, fStatus));
2828  status->Scale(1./status->GetBinContent(kAll));
2829 
2830  Double_t trigEff = status->GetBinContent(kTrigger);
2831  Double_t vtxEff = status->GetBinContent(kIP) / trigEff;
2832  typedef TParameter<double> DP;
2833 
2834  result->Add(new DP("triggerEfficiency", trigEff));
2835  result->Add(new DP("ipEfficiency", vtxEff));
2836 
2837  Double_t nEvents = fIPz->GetEntries();
2838  // Copy ipZ histogram and scale by number of events
2839  TH1* ipZ = static_cast<TH1*>(CloneAndAdd(result, fIPz));
2840  ipZ->Scale(1./nEvents);
2841 
2842  TH1* cent = static_cast<TH1*>(CloneAndAdd(result, fCent));
2843  cent->Scale(1./nEvents);
2844 
2845 
2846  CloneAndAdd(result, fCentIPz);
2847 
2848  Container* measCont = 0;
2849  Container* genCont = 0;
2850  TIter next(fSubs);
2851  Histos* h = 0;
2852  while ((h = static_cast<Histos*>(next()))) {
2853  if (h->GetMask() != AliAODTracklet::kInjection &&
2855  // For the anything but injection or MC-labels, we just finalize
2856  if (!h->MasterFinalize(result, fIPz, tailDelta, tailMax)) {
2857  AliWarningF("Failed to finalize %s/%s", GetName(), h->GetName());
2858  return false;
2859  }
2860  if (h->GetMask() == AliAODTracklet::kGenerated)
2861  // Set MC truth container
2862  genCont = GetC(result, h->GetName());
2863  if (h == fMeasured)
2864  measCont = GetC(result, h->GetName());
2865  continue;
2866  }
2867  if (!EstimateBackground(result, measCont, genCont, h, tailDelta, tailMax)) {
2868  AliWarningF("Failed to estimate Bg in %s/%s", GetName(), h->GetName());
2869  return false;
2870  }
2871  }
2872  return true;
2873 }
2874 
2875 //____________________________________________________________________
2876 Bool_t
2878  Int_t nEvents)
2879 {
2880  // Scale distribution of PDGs to number of events and bin width
2881  if (!fEtaPdg) return true;
2882 
2883 
2884  TH2* etaPdg = static_cast<TH2*>(fEtaPdg->Clone());
2885  etaPdg->SetDirectory(0);
2886  etaPdg->Scale(1. / nEvents, "width");
2887  result->Add(etaPdg);
2888 
2889  // Loop over PDG types and create 2D and 1D distributions
2890  TAxis* yaxis = etaPdg->GetYaxis();
2891  Int_t first = yaxis->GetFirst();
2892  Int_t last = yaxis->GetLast();
2893  THStack* pdgs = new THStack("all","");
2894  THStack* toPion = new THStack("toPion", "");
2895  THStack* toAll = new THStack("toAll", "");
2896  TH1* pion = 0;
2897  Container* pdgOut = new Container();
2898  pdgOut->SetName("mix");
2899  result->Add(pdgOut);
2900  pdgOut->Add(pdgs);
2901  pdgOut->Add(toPion);
2902  pdgOut->Add(toAll);
2903 
2904  TH1* all = static_cast<TH1*>(etaPdg->ProjectionX("total",
2905  1,etaPdg->GetNbinsY()));
2906  all->SetDirectory(0);
2907  all->SetName("total");
2908  all->SetTitle("All");
2909  all->SetYTitle(Form("d#it{N}_{%s}/d#eta", "all"));
2910  all->SetFillColor(kBlack);
2911  all->SetMarkerColor(kBlack);
2912  all->SetLineColor(kBlack);
2913  all->SetMarkerStyle(20);
2914  all->SetFillStyle(0);
2915  all->SetFillColor(0);
2916  all->Reset();
2917  pdgs->Add(all);
2918  for (Int_t i = 1; i <= etaPdg->GetNbinsY(); i++) {
2919  Int_t pdg = TString(yaxis->GetBinLabel(i)).Atoi();
2920  TString nme;
2921  Style_t sty;
2922  Color_t col;
2923  PdgAttr(pdg, nme, col, sty);
2924  if (pdg < 0) pdg = 0;
2925  if (pdg == 22) continue; // Ignore photons
2926 
2927  TH1* h1 = static_cast<TH1*>(etaPdg->ProjectionX(Form("h%d", pdg),i,i));
2928  if (h1->GetEntries() <= 0) continue; // Do not store if empty
2929  h1->SetDirectory(0);
2930  h1->SetName(Form("eta_%d", pdg));
2931  h1->SetTitle(nme);
2932  h1->SetYTitle(Form("d#it{N}_{%s}/d#eta", nme.Data()));
2933  h1->SetFillColor(col);
2934  h1->SetMarkerColor(col);
2935  h1->SetLineColor(col);
2936  h1->SetMarkerStyle(sty);
2937  h1->SetFillStyle(0);
2938  h1->SetFillColor(0);
2939  h1->SetBinContent(0,0);
2940  all->Add(h1);
2941  pdgs->Add(h1);
2942  switch (pdg) {
2943  case 321: h1->SetBinContent(0, 0.15); break; // PRC88,044910
2944  case 2212: h1->SetBinContent(0, 0.05); break; // PRC88,044910
2945  case 310: h1->SetBinContent(0, 0.075); break; // PRL111,222301
2946  case 3122: h1->SetBinContent(0, 0.018); break; // PRL111,222301
2947  case 3212: h1->SetBinContent(0, 0.0055); break; // NPA904,539
2948  case 3322: h1->SetBinContent(0, 0.005); break; // PLB734,409
2949  case 211: h1->SetBinContent(0, 1); break; // it self
2950  default: h1->SetBinContent(0, -1); break; // Unknown
2951  }
2952  pdgOut->Add(h1);
2953 
2954  if (pdg == 211) pion = h1;
2955  }
2956  if (!pdgs->GetHists()) return true;
2957 
2958  TIter next(pdgs->GetHists());
2959  TH1* tmp = 0;
2960  Double_t rmin = +1e9;
2961  Double_t rmax = -1e9;
2962  while ((tmp = static_cast<TH1*>(next()))) {
2963  if (tmp == all) continue;
2964  // Calculate ratio to all
2965  TH1* rat = static_cast<TH1*>(tmp->Clone());
2966  rat->Divide(all);
2967  rat->SetDirectory(0);
2968  rat->SetTitle(Form("%s / all", tmp->GetTitle()));
2969  toAll->Add(rat);
2970 
2971  if (tmp == pion) continue;
2972  Double_t r276 = tmp->GetBinContent(0);
2973  if (r276 < 0 || r276 >= 1) continue;
2974 
2975  // Calulate ratio to pions
2976  rat = static_cast<TH1*>(tmp->Clone());
2977  rat->Divide(pion);
2978  rat->SetTitle(Form("%s / %s", tmp->GetTitle(), pion->GetTitle()));
2979  rat->SetDirectory(0);
2980 
2981  TGraphErrors* g = new TGraphErrors(1);
2982  g->SetName(Form("%s_2760", rat->GetName()));
2983  g->SetTitle(Form("%s in #sqrt{s_{NN}}=2.76TeV", rat->GetTitle()));
2984  g->SetPoint(0,0,r276);
2985  g->SetPointError(0,.5,0);
2986  g->SetLineColor(rat->GetLineColor());
2987  g->SetLineStyle(rat->GetLineStyle());
2988  g->SetMarkerColor(rat->GetMarkerColor());
2989  g->SetMarkerStyle(rat->GetMarkerStyle());
2990  g->SetMarkerSize(1.5*rat->GetMarkerSize());
2991  rat->GetListOfFunctions()->Add(g,"p");
2992  rat->SetMaximum(TMath::Max(rat->GetMaximum(),r276));
2993  rat->SetMinimum(TMath::Min(rat->GetMinimum(),r276));
2994  rmin = TMath::Min(rat->GetMinimum(),rmin);
2995  rmax = TMath::Max(rat->GetMaximum(),rmax);
2996 
2997  toPion->Add(rat);
2998  }
2999  // toPion->SetMinimum(rmin);
3000  toPion->SetMaximum(1.1*rmax);
3001 
3002  return true;
3003 }
3004 
3005 //____________________________________________________________________
3006 Bool_t
3008  Int_t nEvents,
3009  Double_t tailDelta,
3010  Double_t tailMax,
3011  const TString& pre,
3012  const TString& tit)
3013 {
3014  Container* pdgOut = new Container();
3015  pdgOut->SetName(pre);
3016  result->Add(pdgOut);
3017  TAxis* zaxis = fEtaDeltaPdg->GetZaxis();
3018  Int_t first = zaxis->GetFirst();
3019  Int_t last = zaxis->GetLast();
3020  Bool_t isMid = TString(pre).EqualTo("mid");
3021  Int_t nX = fEtaDeltaPdg->GetNbinsX();
3022  Int_t nY = fEtaDeltaPdg->GetNbinsY();
3023  Int_t nZ = fEtaDeltaPdg->GetNbinsZ();
3024  Int_t l1 = fEtaDeltaPdg->GetXaxis()->FindBin(-1);
3025  Int_t r1 = fEtaDeltaPdg->GetXaxis()->FindBin(+1);
3026  Int_t b1 = (isMid ? l1 : 1);
3027  Int_t b2 = (isMid ? r1 : l1-1);
3028  Int_t b3 = (isMid ? -1 : r1+1);
3029  Int_t b4 = (isMid ? -1 : nX);
3030  THStack* stack = new THStack("all", tit);
3031  pdgOut->Add(stack);
3032 
3033  TH1* total = fEtaDeltaPdg->ProjectionY("totalMid",1, 1,1, nZ);
3034  total->SetDirectory(0);
3035  total->SetName("total");
3036  total->SetTitle("Total");
3037  total->SetYTitle(Form("d#it{N}_{%s}/d#Delta", "total"));
3038  total->SetFillColor(kBlack);
3039  total->SetMarkerColor(kBlack);
3040  total->SetLineColor(kBlack);
3041  total->SetMarkerStyle(24);
3042  total->SetFillStyle(0);
3043  total->SetFillColor(0);
3044  total->Reset();
3045  stack->Add(total);
3046 
3047  THStack* ratios = new THStack("ratios","");
3048  TH1* ratioSig = Make1D(0, "ratioSig", "Ratio to all", kGreen+1, 24, *zaxis);
3049  TH1* ratioBg = Make1D(0, "ratioBg", "Ratio to all", kRed+1, 25, *zaxis);
3050  ratioSig->GetXaxis()->LabelsOption("v");
3051  ratioBg ->GetXaxis()->LabelsOption("v");
3052  ratioSig->SetFillColor(kGreen+1);
3053  ratioBg ->SetFillColor(kRed +1);
3054  ratioSig->SetFillStyle(1001);
3055  ratioBg ->SetFillStyle(1001);
3056  ratioSig->SetBarWidth(0.4);
3057  ratioBg ->SetBarWidth(0.4);
3058  ratioSig->SetBarOffset(0.05);
3059  ratioBg ->SetBarOffset(0.55);
3060  ratios->Add(ratioSig, "bar0 e text30");
3061  ratios->Add(ratioBg, "bar0 e text30");
3062  pdgOut->Add(ratios);
3063 
3064  THStack* rows = new THStack("rows","");
3065  TH1* rowSig = new TH1D("rowSig","row",6,0,6);
3066  rowSig->SetDirectory(0);
3067  rowSig->SetYTitle("Fraction of signal");
3068  rowSig->GetXaxis()->SetBinLabel(1, "K^{0}_{S}");
3069  rowSig->GetXaxis()->SetBinLabel(2, "K^{#pm}");
3070  rowSig->GetXaxis()->SetBinLabel(3, "#Lambda");
3071  rowSig->GetXaxis()->SetBinLabel(4, "#Xi");
3072  rowSig->GetXaxis()->SetBinLabel(5, "#Sigma");
3073  rowSig->GetXaxis()->SetBinLabel(6, "Other");
3074  rowSig->GetXaxis()->LabelsOption("v");
3075  rowSig->SetMaximum(100);
3076  rowSig->SetMinimum(0);
3077  rowSig->SetLineColor(kGreen+1);
3078  rowSig->SetMarkerColor(kGreen+1);
3079  rowSig->SetMarkerStyle(24);
3080  TH1* rowBg = static_cast<TH1*>(rowSig->Clone("rowBg"));
3081  rowBg->SetDirectory(0);
3082  rowBg->SetYTitle("Fraction of background");
3083  rowBg->SetLineColor(kRed+1);
3084  rowBg->SetMarkerColor(kRed+1);
3085  rowBg->SetMarkerStyle(25);
3086  rowSig->SetFillStyle(1001);
3087  rowBg ->SetFillColor(1001);
3088  rowSig->SetFillColor(kGreen+1);
3089  rowBg ->SetFillColor(kRed +1);
3090  rowSig->SetBarWidth(0.4);
3091  rowBg ->SetBarWidth(0.4);
3092  rowSig->SetBarOffset(0.05);
3093  rowBg ->SetBarOffset(0.55);
3094  rows->Add(rowSig, "bar0 e text30");
3095  rows->Add(rowBg, "bar0 e text30");
3096  pdgOut->Add(rows);
3097 
3098  for (Int_t i = 1; i <= zaxis->GetNbins(); i++) {
3099  Int_t pdg = TString(zaxis->GetBinLabel(i)).Atoi();
3100  TString nme;
3101  Style_t sty;
3102  Color_t col;
3103  Int_t ipdg = pdg;
3104  PdgAttr(pdg, nme, col, sty);
3105  if (pdg < 0) pdg = 0;
3106  ratioSig->GetXaxis()->SetBinLabel(i, nme);
3107  ratioBg ->GetXaxis()->SetBinLabel(i, nme);
3108  if (pdg == 22) continue; // Ignore photons
3109 
3110 
3111  TH1* h1 = fEtaDeltaPdg->ProjectionY(Form("%d", pdg), b1, b2, i,i);
3112  h1->SetDirectory(0);
3113  if (b3 < b4) {
3114  TH1* h2 = fEtaDeltaPdg->ProjectionY(Form("t%d", pdg), b3, b4, i,i);
3115  h2->SetDirectory(0);
3116  h1->Add(h2);
3117  delete h2;
3118  }
3119  h1->SetUniqueID(i); // Store bin number
3120  if (h1->GetEntries() <= 0) continue; // Do not store if empty
3121  h1->SetName(Form("%d", pdg));
3122  h1->SetTitle(nme);
3123  h1->SetYTitle(Form("d#it{N}_{%s}/d#Delta", nme.Data()));
3124  h1->SetFillColor(col);
3125  h1->SetMarkerColor(col);
3126  h1->SetLineColor(col);
3127  h1->SetMarkerStyle(sty);
3128  h1->SetFillStyle(0);
3129  h1->SetFillColor(0);
3130  h1->SetBinContent(0,ipdg);
3131  total->Add(h1);
3132  stack->Add(h1);
3133  }
3134  total->SetBinContent(0,0);
3135  total->Scale(1./nEvents);
3136  Double_t deltaCut = 1.5;
3137  Double_t eTot, iTot = Integrate(total, 0, tailMax, eTot);
3138  Double_t eSig, iSig = Integrate(total, 0, deltaCut, eSig);
3139  Double_t eBg, iBg = Integrate(total, tailDelta, tailMax, eBg);
3140 
3141  TH1* totInt = new TH1D("totalIntegrals", "Integrals of total:", 3, 0, 3);
3142  totInt->SetDirectory(0);
3143  totInt->GetXaxis()->SetBinLabel(1,Form("0#minus%4.1f", tailMax));
3144  totInt->GetXaxis()->SetBinLabel(2,Form("0#minus%3.1f",deltaCut));
3145  totInt->GetXaxis()->SetBinLabel(3,Form("%3.1f#minus%4.1f",tailDelta,tailMax));
3146  totInt->SetBinContent(1,iTot); totInt->SetBinError(1, eTot);
3147  totInt->SetBinContent(2,iSig); totInt->SetBinError(2, eSig);
3148  totInt->SetBinContent(3,iBg); totInt->SetBinError(3, eBg);
3149  pdgOut->Add(totInt);
3150 
3151  if (!stack->GetHists()) return true;
3152 
3153  Double_t sigOth = 0, eSigOth = 0, bgOth = 0, eBgOth = 0;
3154  Double_t sigK0s = 0, eSigK0s = 0, bgK0s = 0, eBgK0s = 0;
3155  Double_t sigKpm = 0, eSigKpm = 0, bgKpm = 0, eBgKpm = 0;
3156  Double_t sigLam = 0, eSigLam = 0, bgLam = 0, eBgLam = 0;
3157  Double_t sigXi0 = 0, eSigXi0 = 0, bgXi0 = 0, eBgXi0 = 0;
3158  Double_t sigSig = 0, eSigSig = 0, bgSig = 0, eBgSig = 0;
3159  TIter next(stack->GetHists());
3160  TH1* tmp = 0;
3161  while ((tmp = static_cast<TH1*>(next()))) {
3162  // Calculate ratio to all
3163  // Scale(tmp, iTot, eTot);
3164  if (tmp == total) continue;
3165  Int_t pdg = tmp->GetBinContent(0); tmp->SetBinContent(0,0);
3166  tmp->Scale(1./ nEvents);
3167  Double_t eiSig, iiSig = Integrate(tmp, 0, deltaCut, eiSig);
3168  Double_t eiBg, iiBg = Integrate(tmp, tailDelta, tailMax, eiBg);
3169  Double_t erS, rS = RatioE(iiSig, eiSig, iSig, eSig, erS);
3170  Double_t erB, rB = RatioE(iiBg, eiBg, iBg, eBg, erB);
3171  Int_t bin = tmp->GetUniqueID();
3172 #if 0
3173  Printf(" %10s(%4d,%3d) "
3174  "signal=%6.1f/%6.1f=%5.3f+/-%5.3f bg=%6.1f/%6.1f=%5.3f+/-%5.3f",
3175  tmp->GetTitle(), pdg, bin, iiSig, iSig, rS, erS, iiBg, iBg, rB, erB);
3176 #endif
3177  ratioSig->SetBinContent(bin, rS);
3178  ratioSig->SetBinError (bin, erS);
3179  ratioBg ->SetBinContent(bin, rB);
3180  ratioBg ->SetBinError (bin, erB);
3181 
3182  switch (pdg) {
3183  case 321: // K^{+}
3184  sigKpm += rS; eSigKpm += erS*erS; bgKpm += rB; eBgKpm += erB*erB; break;
3185  case 310: // K^{0}_{S}
3186  sigK0s += rS; eSigK0s += erS*erS; bgK0s += rB; eBgK0s += erB*erB; break;
3187  case 3112: // #Sigma^{-}
3188  case 3222: // #Sigma^{+}
3189  // case 3114: // #Sigma^{*-}
3190  // case 3224: // #Sigma^{*+}
3191  // case 4214: // #Sigma^{*+}_{c}
3192  // case 4224: // #Sigma^{*++}_{c}
3193  // case 3212: // #Sigma^{0}
3194  // case 4114: // #Sigma^{*0}_{c}
3195  // case 3214: // #Sigma^{*0}
3196  sigSig += rS; eSigSig += erS*erS; bgSig += rB; eBgSig += erB*erB; break;
3197  case 3312: // #Xi^{-}
3198  // case 3314: // #Xi^{*-}
3199  // case 3324: // #Xi^{*0}
3200  // case 4132: // #Xi^{0}_{c}
3201  // case 4314: // #Xi^{*0}_{c}
3202  sigXi0 += rS; eSigXi0 += erS*erS; bgXi0 += rB; eBgXi0 += erB*erB; break;
3203  // case 4122: // #Lambda^{+}_{c}
3204  case 3122: // #Lambda
3205  sigLam += rS; eSigLam += erS*erS; bgLam += rB; eBgLam += erB*erB; break;
3206  default:
3207  sigOth += rS; eSigOth += erS*erS; bgOth += rB; eBgOth += erB*erB; break;
3208  }
3209  } // while(tmp)
3210  rowSig->SetBinContent(1,100*sigK0s);
3211  rowSig->SetBinError (1,100*TMath::Sqrt(eSigK0s));
3212  rowSig->SetBinContent(2,100*sigKpm);
3213  rowSig->SetBinError (2,100*TMath::Sqrt(eSigKpm));
3214  rowSig->SetBinContent(3,100*sigLam);
3215  rowSig->SetBinError (3,100*TMath::Sqrt(eSigLam));
3216  rowSig->SetBinContent(4,100*sigXi0);
3217  rowSig->SetBinError (4,100*TMath::Sqrt(eSigXi0));
3218  rowSig->SetBinContent(5,100*sigSig);
3219  rowSig->SetBinError (5,100*TMath::Sqrt(eSigSig));
3220  rowSig->SetBinContent(6,100*sigOth);
3221  rowSig->SetBinError (6,100*TMath::Sqrt(eSigOth));
3222  rowBg ->SetBinContent(1,100* bgK0s);
3223  rowBg ->SetBinError (1,100*TMath::Sqrt( eBgK0s));
3224  rowBg ->SetBinContent(2,100* bgKpm);
3225  rowBg ->SetBinError (2,100*TMath::Sqrt( eBgKpm));
3226  rowBg ->SetBinContent(3,100* bgLam);
3227  rowBg ->SetBinError (3,100*TMath::Sqrt( eBgLam));
3228  rowBg ->SetBinContent(4,100* bgXi0);
3229  rowBg ->SetBinError (4,100*TMath::Sqrt( eBgXi0));
3230  rowBg ->SetBinContent(5,100* bgSig);
3231  rowBg ->SetBinError (5,100*TMath::Sqrt( eBgSig));
3232  rowBg ->SetBinContent(6,100* bgOth);
3233  rowBg ->SetBinError (6,100*TMath::Sqrt( eBgOth));
3234 
3235  return true;
3236 }
3237 
3238 //____________________________________________________________________
3239 Bool_t
3241  Int_t nEvents,
3242  Double_t tailDelta,
3243  Double_t tailMax)
3244 {
3245  // Scale distribution of PDGs to number of events and bin width
3246  if (!fEtaDeltaPdg) return true;
3247 
3248  Container* pdgOut = new Container();
3249  pdgOut->SetName("specie");
3250  result->Add(pdgOut);
3251 
3252  ProjectEtaDeltaPdgPart(pdgOut,
3253  nEvents,
3254  tailDelta,
3255  tailMax,
3256  "mid",
3257  "|#eta|<1");
3258  ProjectEtaDeltaPdgPart(pdgOut,
3259  nEvents,
3260  tailDelta,
3261  tailMax,
3262  "fwd",
3263  "|#eta|>1");
3264 }
3265 
3266 //____________________________________________________________________
3267 Bool_t
3269  TH1* ipz,
3270  const TString& shn)
3271 {
3272  // If we have the PDG distributions, we project on eta,IPz, and then
3273  // average over IPz to get dNpdg/deta.
3274  if (!fEtaPdgIPz) return true;
3275 
3276  // Scale each vertex range by number of events in that range
3277  TH3* etaPdgIPz = ScaleToIPz(fEtaPdgIPz, ipz, false);
3278  result->Add(etaPdgIPz);
3279 
3280  // Loop over PDG types and create 2D and 1D distributions
3281  TAxis* yaxis = etaPdgIPz->GetYaxis();
3282  Int_t first = yaxis->GetFirst();
3283  Int_t last = yaxis->GetLast();
3284  THStack* pdgs = new THStack("all","");
3285  THStack* ratios = new THStack("toPion", "");
3286  TH1* pion = 0;
3287  TH2* dtfs = 0;
3288  Container* pdgOut = new Container();
3289  pdgOut->SetName("types");
3290  result->Add(pdgOut);
3291  pdgOut->Add(pdgs);
3292  pdgOut->Add(ratios);
3293  for (Int_t i = 1; i <= etaPdgIPz->GetNbinsY(); i++) {
3294  yaxis->SetRange(i,i);
3295 
3296  Int_t pdg = TString(yaxis->GetBinLabel(i)).Atoi();
3297  TString nme;
3298  Style_t sty;
3299  Color_t col;
3300  PdgAttr(pdg, nme, col, sty);
3301  if (pdg < 0) pdg = 0;
3302 
3303  TH2* h2 = static_cast<TH2*>(etaPdgIPz->Project3D("zx e"));
3304  if (h2->GetEntries() <= 0) continue; // Do not store if empty
3305  h2->SetDirectory(0);
3306  h2->SetName(Form("etaIPz_%d", pdg));
3307  h2->SetTitle(Form("%s#rightarrowX_{%s}", nme.Data(), shn.Data()));
3308  h2->SetYTitle(Form("d#it{N}^{2}_{%s#rightarrow%s}/"
3309  "(d#etadIP_{#it{z}})",
3310  nme.Data(), shn.Data()));
3311  h2->SetFillColor(col);
3312  h2->SetMarkerColor(col);
3313  h2->SetLineColor(col);
3314  pdgOut->Add(h2);
3315 
3316  TH1* h1 = AverageOverIPz(h2, Form("eta_%d",pdg), 1, ipz, 0, false);
3317  h1->SetDirectory(0);
3318  h1->SetYTitle(Form("d#it{N}_{%s#rightarrow%s}/d#eta",
3319  nme.Data(), shn.Data()));
3320  h1->SetMarkerStyle(sty);
3321  pdgs->Add(h1);
3322 
3323  if (pdg == 211) pion = h1;
3324  TH2* tmp = 0;
3325  switch (pdg) {
3326  case 321: // Strange meson K^{+} break;
3327  case 323: // Strange meson K^{*+} break;
3328  case 310: // Strange meson K^{0}_{S} break;
3329  case 130: // Strange meson K^{0}_{L} break;
3330  case 311: // Strange meson K^{0} break;
3331  case 313: // Strange meson K^{*} break;
3332  case 221: // Strange meson #eta break;
3333  case 333: // Strange meson #varphi break;
3334  case 331: // Strange meson #eta' break;
3335  case 3112: // Strange baryon #Sigma^{-} break;
3336  case 3222: // Strange baryon #Sigma^{+} break;
3337  case 3114: // Strange baryon #Sigma^{*-} break;
3338  case 3224: // Strange baryon #Sigma^{*+} break;
3339  case 3312: // Strange baryon #Xi^{-} break;
3340  case 3314: // Strange baryon #Xi^{*-} break;
3341  case 3122: // Strange baryon #Lambda break;
3342  case 3212: // Strange baryon #Sigma^{0} break;
3343  case 3214: // Strange baryon #Sigma^{*0} break;
3344  case 3322: // Strange baryon #Xi^{0} break;
3345  case 3324: // Strange baryon #Xi^{*0} break;
3346  tmp = h2;
3347  break;
3348  default: break;
3349  }
3350  if (tmp) {
3351  if (!dtfs) {
3352  dtfs = static_cast<TH2*>(tmp->Clone("dtfs"));
3353  dtfs->SetDirectory(0);
3354  dtfs->SetTitle("X_{s}#rightarrowX");
3355  dtfs->SetYTitle("IP_{#it{z}} [cm]");
3356  dtfs->Reset();
3357  result->Add(dtfs);
3358  }
3359  // Add up contributions from strange particles
3360  dtfs->Add(tmp);
3361  }
3362  }
3363  if (!pdgs->GetHists()) {
3364  yaxis->SetRange(first, last);
3365  return true;
3366  }
3367 
3368  TIter next(pdgs->GetHists());
3369  TH1* tmp = 0;
3370  while ((tmp = static_cast<TH1*>(next()))) {
3371  if (tmp == pion) continue;
3372  TH1* rat = static_cast<TH1*>(tmp->Clone());
3373  rat->Divide(pion);
3374  rat->SetDirectory(0);
3375  ratios->Add(rat);
3376  }
3377 
3378  yaxis->SetRange(first, last);
3379 }
3380 
3381 //____________________________________________________________________
3382 Bool_t
3384  TH1* ipz,
3385  Double_t tailDelta,
3386  Double_t tailMax)
3387 {
3388  Container* result = new Container;
3389  result->SetName(fName);
3390  result->SetOwner(true);
3391  parent->Add(result);
3392 
3393  // Get the number of events
3394  Double_t nEvents = ipz->GetEntries();
3395 
3396  // Scale each vertex range by number of events in that range
3397  TH2* etaIPz = 0;
3398  if (fEtaIPz) {
3399  etaIPz = ScaleToIPz(fEtaIPz, ipz);
3400  result->Add(etaIPz);
3401  }
3402  // Scale distribution of Pt to number of events and bin width
3403  if (fEtaPt) {
3404  TH2* etaPt = static_cast<TH2*>(fEtaPt->Clone());
3405  etaPt->SetDirectory(0);
3406  etaPt->Scale(1. / nEvents, "width");
3407  result->Add(etaPt);
3408  }
3409  // Short-hand-name
3410  TString shn(etaIPz ? etaIPz->GetTitle() : "X");
3411 
3412  ProjectEtaPdg (result, nEvents);
3413  ProjectEtaDeltaPdg(result, nEvents, tailDelta, tailMax);
3414  ProjectEtaPdgIPz (result, ipz, shn);
3415 
3416  // If we do not have eta vs Delta, just return
3417  if (!fEtaDeltaIPz) return true;
3418 
3419  // Normalize delta distribution to integral number of events
3420  // static_cast<TH3*>(CloneAndAdd(result, fEtaDeltaIPz));
3421  TH3* etaDeltaIPz = ScaleToIPz(fEtaDeltaIPz, ipz, false);
3422  // ipz->GetEntries()>1000);
3423  result->Add(etaDeltaIPz);
3424 
3425  // Make 2D projection to eta,Delta
3426  TH2* etaDelta = ProjectEtaDelta(etaDeltaIPz);
3427  result->Add(etaDelta);
3428 
3429  // Make projection of delta
3430  TH1* delta = ProjectDelta(etaDelta);
3431  result->Add(delta);
3432 
3433  // PArameters of integrals
3434  Double_t maxDelta = etaDeltaIPz->GetYaxis()->GetXmax();
3435  Int_t lowBin = etaDeltaIPz->GetYaxis()->FindBin(tailDelta);
3436  Int_t sigBin = etaDeltaIPz->GetYaxis()->FindBin(1.5);
3437  Int_t highBin = TMath::Min(etaDeltaIPz->GetYaxis()->FindBin(tailMax),
3438  etaDeltaIPz->GetYaxis()->GetNbins());
3439 
3440 
3441  TH1* etaDeltaTail = etaDelta->ProjectionX("etaDeltaTail");
3442  etaDeltaTail ->SetDirectory(0);
3443  etaDeltaTail ->Reset();
3444  etaDeltaTail ->SetTitle(Form("#scale[.7]{#int}_{%3.1f}^{%4.1f}"
3445  "d%s d#it{N}/d%s",
3446  tailDelta,maxDelta,
3447  etaDeltaIPz->GetTitle(),
3448  etaDeltaIPz->GetTitle()));
3449  etaDeltaTail ->SetYTitle("#scale[.7]{#int}_{tail}d#Delta d#it{N}/d#Delta");
3450 
3451  TH2* etaIPzDeltaTail = static_cast<TH2*>(etaDeltaIPz->Project3D("zx e"));
3452  etaIPzDeltaTail->SetName("etaIPzDeltaTail");
3453  etaIPzDeltaTail->SetDirectory(0);
3454  etaIPzDeltaTail->Reset();
3455  etaIPzDeltaTail->SetTitle(etaDeltaTail->GetTitle());
3456  etaIPzDeltaTail->SetZTitle(etaDelta->GetYaxis()->GetTitle());
3457  TH2* etaIPzDeltaMain = static_cast<TH2*>(etaDeltaIPz->Project3D("zx e"));
3458  etaIPzDeltaMain->SetName("etaIPzDeltaMain");
3459  etaIPzDeltaMain->SetDirectory(0);
3460  etaIPzDeltaMain->Reset();
3461  etaIPzDeltaMain->SetTitle(etaDeltaTail->GetTitle());
3462  etaIPzDeltaMain->SetZTitle(etaDelta->GetYaxis()->GetTitle());
3463  // Loop over eta
3464  Double_t intg = 0, eintg = 0;
3465  for (Int_t i = 1; i <= etaDeltaTail->GetNbinsX(); i++) {
3466  // Integrate over Delta
3467  intg = etaDelta->IntegralAndError(i, i, lowBin, highBin, eintg);
3468  etaDeltaTail->SetBinContent(i, intg);
3469  etaDeltaTail->SetBinError (i, eintg);
3470  // Loop over IPz
3471  for (Int_t j = 1; j <= etaIPzDeltaTail->GetNbinsY(); j++) {
3472  // Integrate over Delta
3473  intg = etaDeltaIPz->IntegralAndError(i,i,lowBin,highBin,j,j,eintg);
3474  etaIPzDeltaTail->SetBinContent(i, j, intg);
3475  etaIPzDeltaTail->SetBinError (i, j, eintg);
3476 
3477  intg = etaDeltaIPz->IntegralAndError(i,i,1,sigBin,j,j,eintg);
3478  etaIPzDeltaMain->SetBinContent(i, j, intg);
3479  etaIPzDeltaMain->SetBinError (i, j, eintg);
3480  }
3481  }
3482  result->Add(etaIPzDeltaTail);
3483  result->Add(etaIPzDeltaMain);
3484  result->Add(etaDeltaTail);
3485 
3486  // Integrate full tail
3487  intg = etaDeltaIPz->IntegralAndError(1,etaDeltaIPz->GetNbinsX(),
3488  lowBin, highBin,
3489  1,etaDeltaIPz->GetNbinsZ(),
3490  eintg);
3491  result->Add(new TParameter<double>("deltaTail", intg));
3492  result->Add(new TParameter<double>("deltaTailError", eintg));
3493 
3494  // Some consistency checks:
3495  if (fDebug > 1) {
3496  Printf("%10s: Integral over eta,IPz: %9.4f +/- %9.4f",
3497  GetName(), intg, eintg);
3498  intg = etaDelta->IntegralAndError(1,etaDeltaIPz->GetNbinsX(),
3499  lowBin, highBin,
3500  eintg);
3501  Printf("%10s: Integral over eta: %9.4f +/- %9.4f",
3502  GetName(), intg, eintg);
3503  intg = delta->IntegralAndError(lowBin, highBin, eintg);
3504  Printf("%10s: Integral: %9.4f +/- %9.4f",
3505  GetName(), intg, eintg);
3506  }
3507 
3508  return true;
3509 }
3510 
3511 //____________________________________________________________________
3513  Container* measCont,
3514  Container* genCont,
3515  Histos* h,
3516  Double_t tailCut,
3517  Double_t tailMax)
3518 {
3519  if (!h || !measCont) {
3520  AliWarningF("No sub-histos or measured container in %s", GetName());
3521  return false;
3522  }
3523 
3524  if (!h->MasterFinalize(result, fIPz, tailCut,tailMax)) {
3525  AliWarningF("Failed to finalize %s/%s", GetName(), h->GetName());
3526  return false;
3527  }
3528 
3529  Container* bgCont = GetC(result, h->GetName());
3530  if (!bgCont) {
3531  AliWarningF("%s/%s didn't put a container on output",
3532  GetName(), h->GetName());
3533  return false;
3534  }
3535  const char* sub = h->GetName();
3536  TString shrt(Form("%c", sub[0]));
3537  shrt.ToUpper();
3538  if (genCont) shrt.Append("'");
3539 
3540  TH2* background = 0;
3541  TH2* etaIPzScale = 0;
3542  if (h->GetMask() == AliAODTracklet::kInjection) {
3543  // Get the tail ratio per eta, ipZ
3544  // We get the integral of the measured eta,IPz,Delta dist
3545  // We get the integral of the injected eta,IPz,Delta dist
3546  // Then we take the ratio
3547  etaIPzScale = CopyH2(measCont, "etaIPzDeltaTail",
3548  "etaIPzScale");
3549  etaIPzScale->Divide(GetH2(bgCont, "etaIPzDeltaTail"));
3550  etaIPzScale->SetZTitle("k_{#eta,IP_{#it{z}}}");
3551  etaIPzScale->SetTitle(Form("k_{%s,#eta,IP_{#it{z}}}",shrt.Data()));
3552  bgCont->Add(etaIPzScale);
3553 
3554  // Get the tail ratio per eta
3555  // We get the integral of the measured eta,Delta dist
3556  // We get the integral of the injected eta,Detla dist
3557  // Then we take the ratio
3558  TH1* etaScale = CopyH1(measCont, "etaDeltaTail",
3559  "etaScale");
3560  etaScale->Divide(GetH1(bgCont, "etaDeltaTail"));
3561  etaScale->SetYTitle("k_{#eta}");
3562  etaScale->SetTitle(Form("k_{%s,#eta}", shrt.Data()));
3563  bgCont->Add(etaScale);
3564 
3565  // Get the integrated tail ratio.
3566  // We get the integrated tail of the measured delta dist
3567  // We get the integrated tail of the ianjection delta dist
3568  // We then get the ratio.
3569  Double_t measIntg = GetD(measCont, "deltaTail", -1);
3570  Double_t measIntE = GetD(measCont, "deltaTailError", -1);
3571  Double_t bgIntg = GetD(bgCont, "deltaTail", -1);
3572  Double_t bgIntE = GetD(bgCont, "deltaTailError", -1);
3573  Double_t scaleE = 0;
3574  Double_t scale = RatioE(measIntg, measIntE, bgIntg, bgIntE, scaleE);
3575  bgCont->Add(new TParameter<double>("scale", scale));
3576  bgCont->Add(new TParameter<double>("scaleError", scaleE));
3577 
3578  // Get the fully differential Delta distribution and scale by the
3579  // eta,IPz scalar.
3580  TH3* scaledEtaDeltaIPz = ScaleDelta(CopyH3(bgCont, "etaDeltaIPz",
3581  "scaleEtaDeltaIPz"),
3582  etaIPzScale);
3583  scaledEtaDeltaIPz->SetTitle(Form("%5.3f#times%s",
3584  scale, scaledEtaDeltaIPz->GetTitle()));
3585  scaledEtaDeltaIPz->SetDirectory(0);
3586  scaledEtaDeltaIPz->SetYTitle("k#timesd^{3}#it{N}/"
3587  "(d#Deltad#etadIP_{#it{z}})");
3588  bgCont->Add(scaledEtaDeltaIPz);
3589 #if 0
3590  // scale by derived scalars, rather than by taking the scaled full
3591  // distribution.
3592  TH2* scaledEtaDelta = CopyH2(bgCont, "etaDelta", "scaledEtaDelta");
3593  scaledEtaDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3594  scaledEtaDelta->SetZTitle("k#timesd^{2}#it{N}/(d#Deltad#eta)");
3595  Scale(scaledEtaDelta, etaScale);
3596  bgCont->Add(scaledEtaDelta);
3597 
3598  TH1* scaledDelta = CopyH1(bgCont, "delta", "scaledDelta");
3599  scaledDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3600  scaledDelta->SetYTitle("k#timesd#it{N}/d#Delta");
3601  Scale(scaledDelta,scale,scaleE);
3602  bgCont->Add(scaledDelta);
3603 #else
3604  // Make 2D projection to eta,Delta
3605  TH2* scaledEtaDelta = ProjectEtaDelta(scaledEtaDeltaIPz);
3606  scaledEtaDelta->SetName("scaledEtaDelta");
3607  scaledEtaDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3608  scaledEtaDelta->SetYTitle("k#timesd^{2}#it{N}/(d#Deltad#eta)");
3609  bgCont->Add(scaledEtaDelta);
3610 
3611  // Make projection of delta
3612  TH1* scaledDelta = ProjectDelta(scaledEtaDelta);
3613  scaledDelta->SetName("scaledDelta");
3614  scaledDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3615  scaledDelta->SetYTitle("k#timesd#it{N}/d#Delta");
3616  bgCont->Add(scaledDelta);
3617 #endif
3618  // Make background scaled by full tail
3619  background = CopyH2(bgCont, "etaIPz", "background");
3620  if (!background) AliWarningF("Didn't get background in %s", sub);
3621  else background->Multiply(etaIPzScale);
3622  // else Scale(background, scale, scaleE);
3623  }
3624  else {
3625  // For non-injection sets (i.e., combinatorics) we cannot form the
3626  // scaled background until we have the real data, so we calculate
3627  // beta instead.
3628  background = CopyH2(bgCont, "etaIPz", "background");
3629  TH2* beta = CopyH2(bgCont, "etaIPz", "beta");
3630  if (!background || !beta)
3631  AliWarningF("Didn't get background or beta in %s", sub);
3632  else {
3633  beta->Divide(GetH1(measCont, "etaIPz"));
3634  beta->SetTitle(Form("#beta_{%s}", shrt.Data()));
3635  bgCont->Add(beta);
3636  }
3637  }
3638  if (!background) {
3639  AliWarningF("Didn't get background in %s", sub);
3640  return false;
3641  }
3642  background->SetTitle(Form("#it{B}_{%s}", shrt.Data()));
3643  bgCont->Add(background);
3644 
3645  TH2* signal = CopyH2(measCont, "etaIPz", "signal");
3646  if (!signal) {
3647  AliWarningF("Didn't get signal in %s", sub);
3648  return false;
3649  }
3650  else {
3651  signal->SetTitle(Form("#it{S}_{%s}", shrt.Data()));
3652  signal->Add(background,-1);
3653  // Zero small bins
3654  for (Int_t i = 1; i <= signal->GetNbinsX(); i++) {
3655  for (Int_t j = 1; j <= signal->GetNbinsX(); j++) {
3656  if (signal->GetBinContent(i,j)<1e-6) {
3657  signal->SetBinContent(i,j,0);
3658  signal->SetBinError (i,j,0);
3659  }
3660  }
3661  }
3662  CopyAttr(background, signal);
3663  bgCont->Add(signal);
3664  }
3665 
3666  TH1* alpha = 0;
3667  if (genCont) {
3668  alpha = CopyH2(genCont, "etaIPz", "alpha");
3669  if (alpha && signal) {
3670  alpha->Divide(signal);
3671  alpha->SetTitle(Form("#alpha_{%s}", shrt.Data()));
3672  CopyAttr(signal, alpha);
3673  bgCont->Add(alpha);
3674  }
3675  }
3676 
3677  return true;
3678 }
3679 
3680 //====================================================================
3683  const char* weights,
3684  const char* sumFile,
3685  const char* resFile)
3686 {
3687  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
3688  if (!mgr) {
3689  ::Error("Create","No analysis manager to connect to.");
3690  return 0;
3691  }
3692  AliTrackletAODdNdeta* ret = 0;
3693  if (mc) {
3694  if (weights && weights[0] != '\0') {
3696  new AliTrackletAODWeightedMCdNdeta("MidRapidityMC");
3697  TUrl wurl(weights);
3698  TFile* wfile = TFile::Open(wurl.GetFile());
3699  if (!wfile) {
3700  ::Warning("Create", "Failed to open weights file: %s",
3701  wurl.GetUrl());
3702  return 0;
3703  }
3704  TString wnam(wurl.GetAnchor());
3705  if (wnam.IsNull()) wnam = "weights";
3706  TObject* wobj = wfile->Get(wnam);
3707  if (!wobj) {
3708  ::Warning("Create", "Failed to get weights %s from file %s",
3709  wnam.Data(), wfile->GetName());
3710  return 0;
3711  }
3712  if (!wobj->IsA()->InheritsFrom(AliTrackletBaseWeights::Class())) {
3713  ::Warning("Create", "Object %s from file %s not an "
3714  "AliTrackletBaseWeights but a %s",
3715  wnam.Data(), wfile->GetName(), wobj->ClassName());
3716  return 0;
3717  }
3718  wret->SetWeights(static_cast<AliTrackletBaseWeights*>(wobj));
3719  ret = wret;
3720  }
3721  else
3722  ret = new AliTrackletAODMCdNdeta("MidRapidityMC");
3723  }
3724  else ret = new AliTrackletAODdNdeta("MidRapidity");
3725  if (ret) ret->Connect();
3726 
3727  return ret;
3728 
3729 }
3730 
3731 
3732 
3733 
3734 //____________________________________________________________________
virtual Bool_t ProcessTracklet(AliAODTracklet *tracklet, Double_t ipz, UShort_t signal, Double_t weight)=0
Int_t color[]
static TProfile * GetP1(Container *parent, const char *name, Bool_t verb=true)
virtual void SetWeights(AliTrackletBaseWeights *w)
Int_t pdg
Histos & operator=(const Histos &)
static const TAxis & PdgAxis()
static void PdgAttr(Int_t pdg, TString &nme, Color_t &c, Style_t &s)
static TH2 * Make2D(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis, const TAxis &yAxis)
virtual void SetWeightCalc(UChar_t mode=0)
void Print(Option_t *option="") const
Real_t GetEta() const
double Double_t
Definition: External.C:58
Bool_t Connect(const char *sumFile=0, const char *resFile=0)
const AliVVertex * CheckIP(const AliVVertex *ip, Double_t maxDispersion=0.04, Double_t maxZError=0.25)
virtual Bool_t FinalizeInit(Container *parent)=0
virtual void SetWeightMask(UChar_t mask=0xFF)
AliTrackletAODWeightedMCdNdeta & operator=(const AliTrackletAODWeightedMCdNdeta &o)
Bool_t MasterFinalize(Container *parent, TH1 *ipz, Double_t tailCut, Double_t tailMax)
Bool_t FinalizeInit(Container *parent)
Definition: External.C:244
static void PrintAxis(const TAxis &axis, Int_t nSig=2, const char *alt=0)
Double_t FindMultCentrality(AliVEvent *event, Int_t &nTracklets)
static void SetAxis(TAxis &axis, Int_t n, Double_t *borders)
virtual Double_t LookupWeight(AliAODTracklet *tracklet, Double_t cent, Double_t ipz)
Bool_t ProjectEtaDeltaPdg(Container *result, Int_t nEvents, Double_t tailCut, Double_t tailMax)
Double_t FindCompatCentrality(AliVEvent *event)
virtual void SetWeightInverse(Bool_t inv)
static TProfile2D * Make2P(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis, const TAxis &yAxis)
void ProcessEvent(UInt_t status, Double_t cent, const AliVVertex *ip, TClonesArray *tracklets)
static TH1 * GetH1(Container *parent, const char *name, Bool_t verb=true)
static const char * CentName(Double_t c1, Double_t c2)
void SetEtaAxis(Int_t n, Double_t max)
Bool_t ProjectEtaDeltaPdgPart(Container *result, Int_t nEvents, Double_t tailDelta, Double_t tailMax, const TString &pre, const TString &tit)
static TH3 * ScaleDelta(TH3 *h, TH2 *scale)
Bool_t WorkerInit(Container *parent, const TAxis &etaAxis, const TAxis &ipzAxis, const TAxis &deltaAxis)
TCanvas * c
Definition: TestFitELoss.C:172
virtual void SetWeightVeto(UChar_t veto=0xFF)
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)
void SetMaxDelta(Double_t x=25)
ClassDef(AliTrackletAODdNdeta, 1)
static TH2 * ScaleToIPz(TH2 *h, TH1 *ipZ, Bool_t full=false)
Int_t nCentBins
static TH2 * ProjectEtaDelta(TH3 *h)
void SetCentralityAxis(const TString &spec)
static TH1 * AverageOverIPz(TH2 *h, const char *name, UShort_t mode, TH1 *ipz, TH2 *mask=0, Bool_t verb=true)
static TH1 * ProjectDelta(TH2 *h)
void SetWeightVeto(UChar_t veto=0x0)
static Double_t GetD(Container *parent, const char *name, Double_t def=-1, Bool_t verb=true)
virtual Real_t GetParentPt(Bool_t second=false) const
void SetTailDelta(Double_t x=5)
static TH3 * Make3D(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis, const TAxis &yAxis, const TAxis &zAxis)
static TH1 * Scale(TH1 *h, Double_t x, Double_t xe)
Definition: External.C:92
Bool_t IsCombinatorics() const
CentBin & operator=(const CentBin &)
virtual Bool_t WorkerInit(Container *parent, const TAxis &etaAxis, const TAxis &ipzAxis, const TAxis &deltaAxis)
void SetIPzAxis(Int_t n, Double_t max)
Utilities for midrapidity analysis.
void SetAbsMinCent(Double_t x=-1)
Tracklet AOD object.
static TH1 * MakeStatus(Container *container)
CentBin & operator=(const CentBin &)
void Print(Option_t *option="") const
TClonesArray * FindTracklets(AliVEvent *event)
Bool_t MasterFinalize(Container *parent, TH1 *ipz, Double_t tailCut, Double_t tailMax)
void SetTailMaximum(Double_t x=-1)
void SetEtaAxis(Int_t n, Double_t min, Double_t max)
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
Bool_t InitCentBins(Container *existing)
float Float_t
Definition: External.C:68
virtual void SetDebug(UShort_t lvl)
Real_t GetDPhi() const
static Double_t Integrate(TH1 *h, Double_t min, Double_t max, Double_t &err)
virtual Bool_t MasterFinalize(Container *parent, TH1 *ipz, Double_t tailCut, Double_t tailMax)=0
void SetCentralityMethod(const TString &name)
Bool_t FinalizeInit(Container *parent)
Real_t GetPhi() const
Encode simulation weights for 2nd pass.
const AliVVertex * FindIP(AliVEvent *event, Double_t maxDispersion=0.04, Double_t maxZError=0.25)
Definition: External.C:212
const AliVVertex * FindSimpleIP(AliVEvent *event, Double_t maxDispersion=0.04, Double_t maxZError=0.25)
void SetDPhiShift(Double_t x=0.0045)
static UInt_t Bin2Flag(UInt_t bin)
AliTrackletAODMCdNdeta(const AliTrackletAODdNdeta &o)
Double_t FindCentrality(AliVEvent *event, Int_t &nTracklets)
virtual AliTrackletAODdNdeta::CentBin * MakeCentBin(Float_t c1, Float_t c2)
void SetShiftedDPhiCut(Double_t x=-1)
void Print(Option_t *option="") const
void SetIPzAxis(Int_t n, Double_t min, Double_t max)
const AliVVertex * FindRealIP(AliVEvent *event, Double_t maxDispersion=0.04, Double_t maxZError=0.25)
Bool_t IsMeasured() const
static TH2 * CopyH2(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
virtual Short_t GetParentPdg(Bool_t second=false) const
AliTrackletAODWeightedMCdNdeta(const char *name)
virtual void Print(Option_t *option="") const
Bool_t EstimateBackground(Container *result, Container *measCont, Container *genCont, Histos *h, Double_t tailCut, Double_t tailMax)
void SetDeltaCut(Double_t x=1.5)
Int_t mode
Definition: anaM.C:40
static TH1 * CopyH1(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
virtual const char * GetBranchName() const
static TH3 * GetH3(Container *parent, const char *name, Bool_t verb=true)
const char * GetName() const
static Double_t RatioE(Double_t n, Double_t en, Double_t d, Double_t ed, Double_t &er)
static TH2 * GetH2(Container *parent, const char *name, Bool_t verb=true)
static AliTrackletAODdNdeta * Create(Bool_t mc=false, const char *weights=0, const char *sumFile=0, const char *resFile=0)
UChar_t GetFlags() const
Histos * MakeHistos(const char *name, Color_t color, Style_t style, UShort_t mask, UShort_t veto)
Bool_t IsNeutral() const
AliTrackletAODMCdNdeta(const char *name)
Bool_t MasterFinalize(Container *results)
void SetPhiAxis(Int_t n, Double_t min, Double_t max)
Bool_t IsGenerated() const
Bool_t ProjectEtaPdgIPz(Container *result, TH1 *ipz, const TString &shn)
Float_t nEvents[nProd]
void SetAttr(Color_t c, Style_t m)
Definition: External.C:220
void SetEtaAxis(const TString &spec)
virtual Double_t LookupWeight(AliAODTracklet *tracklet, Double_t cent, Double_t ipz)
virtual void Print(Option_t *option="") const
static Container * GetC(Container *parent, const char *name, Bool_t verb=true)
static Int_t PdgBin(Int_t pdg)
AliTrackletAODWeightedMCdNdeta(const AliTrackletAODWeightedMCdNdeta &o)
Bool_t ProjectEtaPdg(Container *result, Int_t nEvents)
Bool_t ProcessTracklet(AliAODTracklet *tracklet, Double_t ipz, UShort_t signal, Double_t weight)
unsigned short UShort_t
Definition: External.C:28
void Print(Option_t *option="") const
virtual UShort_t CheckTracklet(AliAODTracklet *tracklet) const
void SetIPzAxis(const TString &spec)
const char Option_t
Definition: External.C:48
virtual Bool_t MasterFinalize(Container *results)
Real_t GetDelta() const
void SetWeightMask(UChar_t mask=0xFF)
AliTrackletBaseWeights * fWeights
UInt_t CheckEvent(Double_t &cent, const AliVVertex *&ip, TClonesArray *&tracklets)
static TH3 * CopyH3(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
Histos(const char *name="", Color_t color=kBlack, Style_t style=1, UChar_t mask=0, UChar_t veto=0)
bool Bool_t
Definition: External.C:53
static TObject * CloneAndAdd(Container *c, TObject *o)
virtual void SetDebug(UShort_t lvl)
virtual const char * GetBranchName() const
Bool_t WorkerInit(Container *parent, const TAxis &etaAxis, const TAxis &ipzAxis, const TAxis &deltaAxis)
static TH1 * Make1D(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis)
Bool_t Accept(UInt_t status, Double_t cent, Double_t ipz)
virtual CentBin * MakeCentBin(Float_t c1, Float_t c2)
static TProfile * Make1P(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis)
void SetPhiAxis(Int_t n, Double_t max)
static void FixAxis(TAxis &axis, const char *title=0)
A simplified AOD header.
virtual CentBin * InitCentBin(Int_t bin, Float_t c1, Float_t c2, Container *existing, const TAxis &deltaAxis)
Definition: External.C:196
void SetCentralityAxis(Int_t n, Double_t *bins)
void SetWeights(AliTrackletBaseWeights *w)
AliTrackletAODdNdeta & operator=(const AliTrackletAODdNdeta &o)
Bool_t ProcessTracklet(AliAODTracklet *tracklet, Double_t ipz, UShort_t signal, Double_t weight)
static TProfile * GetP(Container *parent, const char *name, Bool_t verb=true)
static void CopyAttr(const TH1 *src, TH1 *dest)