AliPhysics  a5cd6b6 (a5cd6b6)
 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;
33 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  // -----------------------------------------------------------------
192  void SetCentralityMethod(const TString& name) { fCentMethod = name; }
200  {
201  SetAxis(fCentAxis,n,bins);
202  }
208  void SetCentralityAxis(const TString& spec)
209  {
210  SetAxis(fCentAxis, spec, "-:,");
211  }
225  void SetAbsMinCent(Double_t x=-1) { fAbsMinCent = x; }
233  void SetMaxNTracklet(Double_t maxN) { fMaxNTracklet = maxN; }
234  /* @} */
235  //------------------------------------------------------------------
247  void SetIPzAxis(Int_t n, Double_t min, Double_t max)
248  {
249  SetAxis(fIPzAxis, n, min, max);
250  }
257  void SetIPzAxis(Int_t n, Double_t max)
258  {
259  SetAxis(fIPzAxis, n, max);
260  }
266  void SetIPzAxis(const TString& spec)
267  {
268  SetAxis(fIPzAxis, spec);
269  }
270  /* @} */
271  //------------------------------------------------------------------
283  void SetPhiAxis(Int_t n, Double_t min, Double_t max)
284  {
285  SetAxis(fPhiAxis, n, min, max);
286  }
293  void SetPhiAxis(Int_t n, Double_t max)
294  {
295  SetAxis(fPhiAxis, n, max);
296  }
297  /* @} */
298  //------------------------------------------------------------------
310  void SetEtaAxis(Int_t n, Double_t min, Double_t max)
311  {
312  SetAxis(fEtaAxis, n, min, max);
313  }
320  void SetEtaAxis(Int_t n, Double_t max)
321  {
322  SetAxis(fEtaAxis, n, max);
323  }
329  void SetEtaAxis(const TString& spec)
330  {
331  SetAxis(fEtaAxis, spec);
332  }
333  /* @} */
334  //------------------------------------------------------------------
344  void SetDPhiShift(Double_t x=0.0045) { fDPhiShift = x; }
356  void SetDeltaCut(Double_t x=1.5) { fDeltaCut = x; }
362  void SetMaxDelta(Double_t x=25) { fMaxDelta = x; }
368  void SetTailDelta(Double_t x=5) { fTailDelta = x; }
374  void SetTailMaximum(Double_t x=-1) { fTailMax = x; }
375  /* @} */
376  //------------------------------------------------------------------
386  virtual void SetWeights(AliTrackletBaseWeights* w) {};
392  virtual void SetWeightCalc(UChar_t mode=0) {}
398  virtual void SetWeightMask(UChar_t mask=0xFF) {}
406  virtual void SetWeightInverse(Bool_t inv) {}
412  virtual void SetWeightVeto(UChar_t veto=0xFF) {}
413  /* @} */
414  virtual void SetTriggerEfficiency(Double_t eff) { fTriggerEff = eff; }
415  virtual void SetMinEta1(Int_t n) { fMinEta1 = n; }
416  /* @} */
417  //__________________________________________________________________
421  struct Sub : public TObject
422  {
427  Sub(const char* name="") : fName(name), fContainer(0), fDebug(0) {}
433  Sub(const Sub& o) : TObject(o), fName(o.fName), fContainer(0), fDebug(0) {}
437  virtual ~Sub() {}
443  Sub& operator=(const Sub&) { return *this; }
447  const char* GetName() const { return fName.Data(); }
458  virtual Bool_t WorkerInit(Container* parent,
459  const TAxis& etaAxis,
460  const TAxis& ipzAxis,
461  const TAxis& deltaAxis)
462  {
463  fContainer = new Container;
464  fContainer->SetName(fName);
465  fContainer->SetOwner();
466  if (parent) parent->Add(fContainer);
467  return true;
468  }
479  virtual Bool_t ProcessTracklet(AliAODTracklet* tracklet,
480  Double_t ipz,
481  UShort_t signal,
482  Double_t weight) = 0;
492  virtual Bool_t FinalizeInit(Container* parent) = 0;
503  virtual Bool_t MasterFinalize(Container* parent,
504  TH1* ipz,
505  Double_t tailCut,
506  Double_t tailMax) = 0;
512  virtual void SetDebug(UShort_t lvl) { fDebug = lvl; }
513  protected:
517  Container* fContainer;
520  ClassDef(Sub,1);
521  };
522 
523  //__________________________________________________________________
533  struct Histos : public Sub
534  {
539  Histos(const char* name="", Color_t color=kBlack, Style_t style=1,
540  UChar_t mask=0, UChar_t veto=0)
541  : Sub(name),
542  fMask(mask),
543  fVeto(veto),
544  fEtaIPz(0),
545  fEtaDeltaIPz(0),
546  fEtaDeltaPdg(0),
547  fEtaPdgIPz(0),
548  fEtaPdg(0),
549  fEtaPt(0),
550  fColor(color),
551  fStyle(style)
552  {}
558  Histos(const Histos& o)
559  : Sub(o),
560  fMask(o.fMask),
561  fVeto(o.fVeto),
562  fEtaIPz(0),
563  fEtaDeltaIPz(0),
564  fEtaDeltaPdg(0),
565  fEtaPdgIPz(0),
566  fEtaPdg(0),
567  fEtaPt(0),
568  fColor(o.fColor),
569  fStyle(o.fStyle)
570  {}
574  virtual ~Histos() {}
580  Histos& operator=(const Histos&) { return *this; }
585  {
586  return fMask == kMeasuredMask && fVeto == kMeasuredVeto;
587  }
592  {
593  return fMask == kInjectedMask && fVeto == kInjectedVeto;
594  }
599  {
601  }
607  {
608  return fMask == kDistinctMask && fVeto == kDistinctVeto;
609  }
614  {
615  return fMask == kPrimaryMask && fVeto == kPrimaryVeto;
616  }
621  {
622  return fMask == kSecondaryMask && fVeto == kSecondaryVeto;
623  }
628  {
629  return fMask == kGeneratedMask && fVeto == kGeneratedVeto;
630  }
641  Bool_t WorkerInit(Container* parent,
642  const TAxis& etaAxis,
643  const TAxis& ipzAxis,
644  const TAxis& deltaAxis);
651  void SetAttr(Color_t c, Style_t m);
663  Double_t ipz,
664  UShort_t signal,
665  Double_t weight);
675  Bool_t FinalizeInit(Container* parent);
686  Bool_t MasterFinalize(Container* parent,
687  TH1* ipz,
688  Double_t tailCut,
689  Double_t tailMax);
690 
691  UChar_t GetMask() const { return fMask; }
692  UChar_t GetVeto() const { return fVeto; }
698  void Print(Option_t* option="") const;
699  protected:
700  Bool_t ProjectEtaPdg(Container* result, Int_t nEvents);
701  Bool_t ProjectEtaDeltaPdgPart(Container* result,
702  Int_t nEvents,
703  Double_t tailDelta,
704  Double_t tailMax,
705  const TString& pre,
706  const TString& tit);
707  Bool_t ProjectEtaDeltaPdg(Container* result,
708  Int_t nEvents,
709  Double_t tailCut,
710  Double_t tailMax);
711  Bool_t ProjectEtaPdgIPz(Container* result,
712  TH1* ipz,
713  const TString& shn);
714  const UChar_t fMask;
715  const UChar_t fVeto;
722  public:
723  Color_t fColor;
724  Style_t fStyle;
725 
726  ClassDef(Histos,1);
727  };
728 
729  //__________________________________________________________________
745  struct CentBin : public Sub
746  {
751  : Sub(""),
752  fSubs(0),
753  fLow(0),
754  fHigh(0),
755  fStatus(0),
756  fIPz(0),
757  fCent(0),
758  fCentIPz(0),
759  fMeasured(0),
760  fInjection(0)
761  {
762  }
769  CentBin(Double_t c1, Double_t c2);
775  CentBin(const CentBin& o)
776  : Sub(o),
777  fSubs(0),
778  fLow(o.fLow),
779  fStatus(0),
780  fIPz(0),
781  fCent(0),
782  fCentIPz(0),
783  fHigh(o.fHigh),
784  fMeasured(0),
785  fInjection(0)
786  {}
790  virtual ~CentBin() {}
796  CentBin& operator=(const CentBin&) { return *this; }
808  Histos* MakeHistos(const char* name,
809  Color_t color,
810  Style_t style,
811  UShort_t mask,
812  UShort_t veto);
823  Bool_t WorkerInit(Container* parent,
824  const TAxis& etaAxis,
825  const TAxis& ipzAxis,
826  const TAxis& deltaAxis);
832  Bool_t IsAllBin() const;
842  Bool_t Accept(UInt_t status, Double_t cent, Double_t ipz);
854  Double_t ipz,
855  UShort_t signal,
856  Double_t weight);
861  void Completed();
871  Bool_t FinalizeInit(Container* parent);
882  Bool_t MasterFinalize(Container* parent,
883  TH1* ipz,
884  Double_t tailCut,
885  Double_t tailMax);
898  Bool_t EstimateBackground(Container* result,
899  Container* measCont,
900  Container* genCont,
901  Histos* h,
902  Double_t tailCut,
903  Double_t tailMax);
909  void Print(Option_t* option="") const;
915  virtual void SetDebug(UShort_t lvl);
916  protected:
917  Container* fSubs;
923  TProfile* fCentIPz;
926 
927  ClassDef(CentBin,2);
928  };
929 
930 protected:
931  // -----------------------------------------------------------------
941  virtual Bool_t WorkerInit();
948  Bool_t InitCentBins(Container* existing);
958  {
959  return new CentBin(c1, c2);
960  }
972  virtual CentBin* InitCentBin(Int_t bin,
973  Float_t c1,
974  Float_t c2,
975  Container* existing,
976  const TAxis& deltaAxis);
977  /* @} */
978 
979  // -----------------------------------------------------------------
984  virtual const char* GetBranchName() const { return "AliAODTracklets"; }
992  static UInt_t Bin2Flag(UInt_t bin);
1002  UInt_t CheckEvent(Double_t& cent,
1003  const AliVVertex*& ip,
1004  TClonesArray*& tracklets);
1012  TClonesArray* FindTracklets(AliVEvent* event);
1022  const AliVVertex* FindSimpleIP(AliVEvent* event,
1023  Double_t maxDispersion=0.04,
1024  Double_t maxZError=0.25);
1034  const AliVVertex* FindRealIP(AliVEvent* event,
1035  Double_t maxDispersion=0.04,
1036  Double_t maxZError=0.25);
1046  const AliVVertex* CheckIP(const AliVVertex* ip,
1047  Double_t maxDispersion=0.04,
1048  Double_t maxZError=0.25);
1049 
1059  const AliVVertex* FindIP(AliVEvent* event,
1060  Double_t maxDispersion=0.04,
1061  Double_t maxZError=0.25);
1075  Double_t FindMultCentrality(AliVEvent* event,
1076  Double_t& value,
1077  Int_t& nTracklets,
1078  Int_t& nCl0,
1079  Int_t& nCl1);
1080  Double_t FindMultCentrality(AliVEvent* event,
1081  Int_t& nTracklets);
1092  Double_t FindFakeCentrality(AliVEvent* event, Int_t& nTracklets);
1102  Double_t FindCompatCentrality(AliVEvent* event);
1114  Double_t FindCentrality(AliVEvent* event, Int_t& nTracklets);
1120  Bool_t FindTrigger();
1121  /* @} */
1122 
1123  // -----------------------------------------------------------------
1137  virtual Double_t LookupWeight(AliAODTracklet* tracklet,
1138  Double_t cent,
1139  Double_t ipz);
1149  virtual UShort_t CheckTracklet(AliAODTracklet* tracklet) const;
1150  /* @} */
1151 
1152  // -----------------------------------------------------------------
1165  void ProcessEvent(UInt_t status,
1166  Double_t cent,
1167  const AliVVertex* ip,
1168  TClonesArray* tracklets);
1169  /* @} */
1170 
1171  // -----------------------------------------------------------------
1183  virtual Bool_t MasterFinalize(Container* results);
1184  /* @} */
1185 
1193  static TH1* MakeStatus(Container* container);
1194 
1195  // -----------------------------------------------------------------
1197  Container* fContainer;
1199  Container* fCentBins;
1202 
1204 
1206 
1208 
1209  TProfile* fNBareVsGood;
1211  TProfile* fNBareVsFake;
1213  TProfile* fNTrackletVsGood;
1215  TProfile* fNTrackletVsFake;
1219  TProfile* fNGeneratedVsFake;
1221  TProfile* fCentTracklets;
1222 
1223  TProfile* fCentEst;
1224 
1225  TProfile2D* fCl0Cl1Tracklets;
1264 
1266 };
1267 //====================================================================
1274 {
1275 public:
1277  typedef TList Container;
1283  {}
1287  AliTrackletAODMCdNdeta(const char* name)
1288  : AliTrackletAODdNdeta(name)
1289  {}
1297  {}
1302  //__________________________________________________________________
1319  {
1325  fCombinatorics(0),
1326  fDistinct(0),
1327  fPrimaries(0),
1328  fSecondaries(0),
1329  fGenerated(0)
1330  {
1331  }
1338  CentBin(Double_t c1, Double_t c2);
1344  CentBin(const CentBin& o)
1346  fCombinatorics(0),
1347  fPrimaries(0),
1348  fSecondaries(0),
1349  fGenerated(0)
1350  {}
1354  virtual ~CentBin() {}
1360  CentBin& operator=(const CentBin&) { return *this; }
1361  protected:
1367 
1368  ClassDef(CentBin,1);
1369  };
1370 
1371 protected:
1372  // -----------------------------------------------------------------
1386  {
1387  return new CentBin(c1, c2);
1388  }
1389  /* @} */
1390 
1391  // -----------------------------------------------------------------
1396  virtual const char* GetBranchName() const { return "AliAODMCTracklets"; }
1397  /* @} */
1398 
1400 };
1401 //====================================================================
1408 {
1409 public:
1411  typedef TList Container;
1417  fWeights(0),
1418  fEtaWeight(0),
1419  fWeightCorr(0)
1420  {}
1425  : AliTrackletAODMCdNdeta(name),
1426  fWeights(0),
1427  fEtaWeight(0),
1428  fWeightCorr(0)
1429  {}
1437  fWeights(0),
1438  fEtaWeight(0),
1439  fWeightCorr(0)
1440  {}
1459  void Print(Option_t* option="") const;
1460 
1461  // -----------------------------------------------------------------
1471  void SetWeights(AliTrackletBaseWeights* w) { fWeights = w; }
1477  void SetWeightCalc(UChar_t mode=0) {
1478  Info("SetWeightCalc", "Use=%d", mode);
1479  if (fWeights) fWeights->SetCalc(mode);}
1485  void SetWeightMask(UChar_t mask=0xFF) {
1486  Info("SetWeightMask", "mask=0x%x", mask);
1487  if (fWeights) fWeights->SetMask(mask); }
1493  void SetWeightVeto(UChar_t veto=0x0) {
1494  Info("SetWeightVeto", "veto=0x%x", veto);
1495  if (fWeights) fWeights->SetVeto(veto); }
1503  void SetWeightInverse(Bool_t inv) { if (fWeights) fWeights->SetInverse(inv); }
1504  /* @} */
1505 protected:
1506  // -----------------------------------------------------------------
1516  Bool_t WorkerInit();
1526  virtual Double_t LookupWeight(AliAODTracklet* tracklet,
1527  Double_t cent,
1528  Double_t ipz);
1529  /* @} */
1530  // -----------------------------------------------------------------
1542  Bool_t MasterFinalize(Container* results);
1543  /* @} */
1544  // -----------------------------------------------------------------
1546  TProfile2D* fEtaWeight;
1549 };
1550 
1551 //____________________________________________________________________
1553  : AliAnalysisTaskSE(),
1554  fContainer(0),
1555  fCentBins(0),
1556  fIPz(0),
1557  fCent(0),
1558  fStatus(0),
1559  fEtaPhi(0),
1560  fNBareVsGood(0),
1561  fNBareVsFake(0),
1562  fNTrackletVsGood(0),
1563  fNTrackletVsFake(0),
1564  fNGeneratedVsGood(0),
1565  fNGeneratedVsFake(0),
1566  fCentTracklets(0),
1567  fCentEst(0),
1568  fCl0Cl1Tracklets(0),
1569  fCentMethod(""),
1570  fCentIdx(-1),
1571  fTrkIdx(-1),
1572  fCl0Idx(-1),
1573  fCl1Idx(-1),
1574  fCentAxis(1,0,0),
1575  fIPzAxis(1,0,0),
1576  fEtaAxis(1,0,0),
1577  fPhiAxis(1,0,0),
1578  fMaxDelta(0),
1579  fTailDelta(0),
1580  fTailMax(-1),
1581  fDPhiShift(0),
1582  fShiftedDPhiCut(0),
1583  fDeltaCut(0),
1584  fAbsMinCent(-1),
1585  fMaxNTracklet(6000),
1586  fTriggerEff(0),
1587  fMinEta1(0)
1588 {}
1589 //____________________________________________________________________
1591  : AliAnalysisTaskSE(name),
1592  fContainer(0),
1593  fCentBins(0),
1594  fIPz(0),
1595  fCent(0),
1596  fStatus(0),
1597  fEtaPhi(0),
1598  fNBareVsGood(0),
1599  fNBareVsFake(0),
1600  fNTrackletVsGood(0),
1601  fNTrackletVsFake(0),
1602  fNGeneratedVsGood(0),
1603  fNGeneratedVsFake(0),
1604  fCentTracklets(0),
1605  fCentEst(0),
1606  fCl0Cl1Tracklets(0),
1607  fCentMethod("V0M"),
1608  fCentIdx(-1),
1609  fTrkIdx(-1),
1610  fCl0Idx(-1),
1611  fCl1Idx(-1),
1612  fCentAxis(10,0,100),
1613  fIPzAxis(30,-15,+15),
1614  fEtaAxis(16,-2,+2),
1615  fPhiAxis(100,0,TMath::TwoPi()),
1616  fMaxDelta(25),
1617  fTailDelta(5),
1618  fTailMax(-1),
1619  fDPhiShift(0.0045),
1620  fShiftedDPhiCut(-1),
1621  fDeltaCut(1.5),
1622  fAbsMinCent(-1),
1623  fMaxNTracklet(6000),
1624  fTriggerEff(0),
1625  fMinEta1(0)
1626 {
1627  FixAxis(fCentAxis, "Centrality [%]");
1628  FixAxis(fIPzAxis, "IP_{#it{z}} [cm]");
1629  FixAxis(fEtaAxis, "#eta");
1630  FixAxis(fPhiAxis, "#phi");
1631 
1632  DefineOutput(1, Container::Class());
1633  DefineOutput(2, Container::Class());
1634 }
1635 //____________________________________________________________________
1637  : AliAnalysisTaskSE(o),
1638  fContainer(0),
1639  fCentBins(0),
1640  fIPz(0),
1641  fCent(0),
1642  fStatus(0),
1643  fEtaPhi(0),
1644  fNBareVsGood(0),
1645  fNBareVsFake(0),
1646  fNTrackletVsGood(0),
1647  fNTrackletVsFake(0),
1648  fNGeneratedVsGood(0),
1649  fNGeneratedVsFake(0),
1650  fCentTracklets(0),
1651  fCentEst(0),
1652  fCl0Cl1Tracklets(0),
1653  fCentMethod(o.fCentMethod),
1654  fCentIdx(o.fCentIdx),
1655  fTrkIdx(o.fTrkIdx),
1656  fCl0Idx(o.fCl0Idx),
1657  fCl1Idx(o.fCl1Idx),
1658  fCentAxis(o.fCentAxis),
1659  fIPzAxis(o.fIPzAxis),
1660  fEtaAxis(o.fEtaAxis),
1661  fPhiAxis(o.fPhiAxis),
1662  fMaxDelta(o.fMaxDelta),
1663  fTailDelta(o.fTailDelta),
1664  fTailMax(o.fTailMax),
1665  fDPhiShift(o.fDPhiShift),
1666  fShiftedDPhiCut(o.fShiftedDPhiCut),
1667  fDeltaCut(o.fDeltaCut),
1668  fAbsMinCent(o.fAbsMinCent),
1669  fMaxNTracklet(o.fMaxNTracklet),
1670  fTriggerEff(o.fTriggerEff),
1671  fMinEta1(o.fMinEta1)
1672 {}
1673 //____________________________________________________________________
1676 {
1677  if (&o == this) return *this;
1678  if (fContainer) {
1679  delete fContainer;
1680  fContainer = 0;
1681  }
1682  if (fCentBins) {
1683  delete fCentBins;
1684  fCentBins = 0;
1685  }
1686  fIPz = 0;
1687  fCent = 0;
1689  fCentIdx = o.fCentIdx;
1690  fTrkIdx = o.fTrkIdx;
1691  fCl0Idx = o.fCl0Idx;
1692  fCl1Idx = o.fCl1Idx;
1693  fCentAxis = o.fCentAxis;
1694  fIPzAxis = o.fIPzAxis;
1695  fEtaAxis = o.fEtaAxis;
1696  fPhiAxis = o.fPhiAxis;
1697  fMaxDelta = o.fMaxDelta;
1698  fTailDelta = o.fTailDelta;
1699  fTailMax = o.fTailMax;
1700  fDPhiShift = o.fDPhiShift;
1702  fDeltaCut = o.fDeltaCut;
1706  fMinEta1 = o.fMinEta1;
1707  return *this;
1708 }
1709 //____________________________________________________________________
1712  o)
1713 {
1714  if (&o == this) return *this;
1716  return *this;
1717 }
1718 //____________________________________________________________________
1720  const char* resFile)
1721 {
1722  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1723  if (!mgr) {
1724  AliError("No analysis manager to connect to.");
1725  return false;
1726  }
1727 
1728  // Add to the manager
1729  mgr->AddTask(this);
1730 
1731  // Create and connect output containers
1732  TString sumOut;
1733  TString resOut;
1734  if (sumFile && sumFile[0] != '\0') sumOut = sumFile;
1735  if (resFile && resFile[0] != '\0') resOut = resFile;
1736  else if (sumFile && sumFile[0] != '\0') resOut = sumFile;
1737 
1738  // If the string is null or 'default' connect to standard output file
1739  if (sumOut.IsNull() || sumOut.EqualTo("default", TString::kIgnoreCase))
1740  sumOut = AliAnalysisManager::GetCommonFileName();
1741  // If the string is null or 'default' connect to standard output file
1742  if (resOut.IsNull() || resOut.EqualTo("default", TString::kIgnoreCase))
1743  resOut = AliAnalysisManager::GetCommonFileName();
1744 
1745  // Always connect input
1746  mgr->ConnectInput(this, 0, mgr->GetCommonInputContainer());
1747 
1748  // Connect sum list unless the output 'none' is specified
1749  if (!sumOut.EqualTo("none", TString::kIgnoreCase)) {
1750  TString sumName(Form("%sSums", GetName()));
1751  AliAnalysisDataContainer* sumCon =
1752  mgr->CreateContainer(sumName, TList::Class(),
1753  AliAnalysisManager::kOutputContainer, sumOut);
1754  mgr->ConnectOutput(this, 1, sumCon);
1755  }
1756 
1757  // Connect the result list unless the output 'none' is specified
1758  if (!resOut.EqualTo("none", TString::kIgnoreCase)) {
1759  TString resName(Form("%sResults", GetName()));
1760  AliAnalysisDataContainer* resCon =
1761  mgr->CreateContainer(resName, TList::Class(),
1762  AliAnalysisManager::kParamContainer, resOut);
1763  mgr->ConnectOutput(this, 2, resCon);
1764  }
1765  return true;
1766 }
1767 //____________________________________________________________________
1769 {
1770  Double_t shiftedDPhiCut = fShiftedDPhiCut;
1771  if (shiftedDPhiCut < 0) shiftedDPhiCut = TMath::Sqrt(fDeltaCut)*0.06;
1772  Double_t tailMax = fTailMax;
1773  if (tailMax < 0) tailMax = fMaxDelta;
1774 
1775  Printf("%s: %s", ClassName(), GetName());
1776  Printf(" %22s: 0x%08x", "Off-line trigger mask", fOfflineTriggerMask);
1777  Printf(" %22s: %f", "Delta phi shift", fDPhiShift);
1778  Printf(" %22s: %f", "Shifted Delta phi cut", shiftedDPhiCut);
1779  Printf(" %22s: %f", "Delta cut", fDeltaCut);
1780  Printf(" %22s: %f", "max Delta", fMaxDelta);
1781  Printf(" %22s: %f", "tail Delta", fTailDelta);
1782  Printf(" %22s: %f", "tail maximum", tailMax);
1783  Printf(" %22s: %f%%", "Absolute least c", fAbsMinCent);
1784  Printf(" %22s: %f", "Max(Ntracklet)", fMaxNTracklet);
1785  Printf(" %22s: %f", "Trigger eff)", fTriggerEff);
1786  Printf(" %22s: %d", "Min(Ntracklet)", fMinEta1);
1787  Printf(" %22s: %s (%d)", "Centrality estimator", fCentMethod.Data(),fCentIdx);
1790  PrintAxis(fIPzAxis,1,"IPz");
1791  PrintAxis(fCentAxis,0);
1792 
1793  if (!fCentBins) return;
1794 
1795  Printf("--- Centrality bins");
1796  TIter next(fCentBins);
1797  CentBin* bin = 0;
1798  while ((bin = static_cast<CentBin*>(next()))) {
1799  bin->Print(option);
1800  }
1801 }
1802 //____________________________________________________________________
1804 {
1806  if (!fWeights) return;
1807  fWeights->Print(option);
1808 }
1809 
1810 //____________________________________________________________________
1812 {
1813  Printf(" Centrality bin: %s", fName.Data());
1814  Printf(" Low cut: %5.1f", fLow);
1815  Printf(" High cut: %5.1f", fHigh);
1816 
1817  if (!fSubs) return;
1818  Printf(" --- Histogram sets");
1819  TIter next(fSubs);
1820  Histos* h = 0;
1821  while ((h = static_cast<Histos*>(next()))) {
1822  h->Print(option);
1823  }
1824 }
1825 
1826 //____________________________________________________________________
1828 {
1829  Sub::SetDebug(lvl);
1830  TIter next(fSubs);
1831  Histos* h = 0;
1832  while ((h = static_cast<Histos*>(next()))) {
1833  h->SetDebug(lvl);
1834  }
1835 }
1836 namespace {
1837  void Bits2String(UChar_t m, char out[7])
1838  {
1839  if (m & AliAODTracklet::kInjection) out[0] = 'I'; else out[0] = '-';
1840  if (m & AliAODTracklet::kCombinatorics) out[1] = 'C'; else out[1] = '-';
1841  if (m & AliAODTracklet::kSecondary) out[2] = 'S'; else out[2] = '-';
1842  if (m & AliAODTracklet::kDistinct) out[3] = 'D'; else out[3] = '-';
1843  if (m & AliAODTracklet::kSimulated) out[4] = 'X'; else out[4] = '-';
1844  if (m & AliAODTracklet::kGenerated) out[5] = 'G'; else out[5] = '-';
1845  out[6] = '\0';
1846  }
1847 }
1848 
1849 //____________________________________________________________________
1851 {
1852  char cMask[7]; Bits2String(fMask, cMask);
1853  char cVeto[7]; Bits2String(fVeto, cVeto);
1854  Printf(" Histograms: %s", fName.Data());
1855  Printf(" Mask: 0x%02x (%s)", fMask, cMask);
1856  Printf(" Veto: 0x%02x (%s)", fVeto, cVeto);
1857  Printf(" Delta: %s", fEtaDeltaIPz ? "yes" : "no");
1858  Printf(" Delta per PDG: %s", fEtaDeltaPdg ? "yes" : "no");
1859 }
1860 
1861 
1862 //____________________________________________________________________
1863 void
1865 {
1866  if (!WorkerInit()) {
1867  AliWarning("Failed to initialize on worker");
1868  return;
1869  }
1870  PostData(1,fContainer);
1871 }
1872 
1873 //____________________________________________________________________
1875 {
1876  if (DebugLevel() > 1) Printf("Initialising on worker");
1877  if (fShiftedDPhiCut < 0) fShiftedDPhiCut = TMath::Sqrt(fDeltaCut)*0.06;
1878  if (fTailMax < 0) fTailMax = fMaxDelta;
1879 
1880  fContainer = new Container;
1881  fContainer->SetName(Form("%sSums", GetName()));
1882  fContainer->SetOwner();
1883 
1884  Bool_t save = TH1::GetDefaultSumw2();
1885  TH1::SetDefaultSumw2();
1886  fIPz = Make1D(fContainer, "ipz", "", kMagenta+2, 20, fIPzAxis);
1887  fCent = Make1D(fContainer, "cent", "", kMagenta+2, 20, fCentAxis);
1888  fEtaPhi = Make2D(fContainer, "etaPhi","",kMagenta+2, 20, fEtaAxis,fPhiAxis);
1889 
1890  TAxis trackletAxis(1000, 0, 10000);
1891  FixAxis(trackletAxis, "#it{N}_{tracklets}");
1892  fNBareVsGood = Make1P(fContainer, "nBareVsGood", "Good",
1893  kGreen+2, 20, trackletAxis);
1894  fNBareVsFake = Make1P(fContainer, "nBareVsFake", "Fake",
1895  kRed+2, 20, trackletAxis);
1896  fNBareVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1897  fNBareVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1898  fNBareVsGood->SetStats(0);
1899  fNBareVsFake->SetStats(0);
1900 
1901  fNTrackletVsGood = Make1P(fContainer, "nTrackletVsGood", "Good",
1902  kGreen+2, 20, trackletAxis);
1903  fNTrackletVsFake = Make1P(fContainer, "nTrackletVsFake", "Fake",
1904  kRed+2, 20, trackletAxis);
1905  fNTrackletVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1906  fNTrackletVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1907  fNTrackletVsGood->SetStats(0);
1908  fNTrackletVsFake->SetStats(0);
1909 
1910  TAxis generatedAxis(1000, 0, 15000);
1911  FixAxis(generatedAxis, "#it{N}_{generated,|#eta|<2}");
1912  fNGeneratedVsGood = Make1P(fContainer, "nGeneratedVsGood", "Good",
1913  kGreen+2, 24, generatedAxis);
1914  fNGeneratedVsFake = Make1P(fContainer, "nGeneratedVsFake", "Fake",
1915  kRed+2, 24, generatedAxis);
1916  fNGeneratedVsGood->SetYTitle("#LT#it{N}_{tracklet}#GT");
1917  fNGeneratedVsFake->SetYTitle("#LT#it{N}_{tracklet}#GT");
1918  fNGeneratedVsGood->SetStats(0);
1919  fNGeneratedVsFake->SetStats(0);
1920 
1921  fCentTracklets = new TProfile("centTracklets",
1922  "Mean number of tracklets per centrality",
1923  1050, 0, 105);
1924  fCentTracklets->SetXTitle("Centrality [%]");
1925  fCentTracklets->SetYTitle("#LTtracklets#GT");
1926  fCentTracklets->SetDirectory(0);
1927  fCentTracklets->SetMarkerStyle(20);
1928  fCentTracklets->SetMarkerColor(kMagenta+2);
1929  fCentTracklets->SetLineColor(kMagenta+2);
1930  fContainer->Add(fCentTracklets);
1931  fCentEst = Make1P(fContainer,"centEstimator","",kMagenta+2, 20, fCentAxis);
1932  fCentEst->SetYTitle(Form("#LT%s#GT",fCentMethod.Data()));
1933 
1934  TAxis clAxis(150, 0, 15000);
1935  FixAxis(clAxis, "#it{N}_{cluster}");
1936  fCl0Cl1Tracklets = Make2P(fContainer, "cl0cl1Tracklets","",kMagenta+2, 20,
1937  clAxis, clAxis);
1938  fCl0Cl1Tracklets->SetXTitle(Form("Layer 0 %s", clAxis.GetTitle()));
1939  fCl0Cl1Tracklets->SetYTitle(Form("Layer 1 %s", clAxis.GetTitle()));
1940  fCl0Cl1Tracklets->SetZTitle("#LT#it{N}_{tracklet}#GT");
1941 
1942  TH1::SetDefaultSumw2(save);
1943 
1945 
1946  typedef TParameter<double> DP;
1947  typedef TParameter<bool> BP;
1948  typedef TParameter<int> IP;
1949  Container* params = new Container;
1950  params->SetName("parameters");
1951  params->SetOwner();
1952  fContainer->Add(params);
1953  params->Add(new DP("DPhiShift", fDPhiShift, 'f'));
1954  params->Add(new DP("ShiftedDPhiCut", fShiftedDPhiCut, 'f'));
1955  params->Add(new DP("DeltaCut", fDeltaCut, 'f'));
1956  params->Add(new DP("MaxDelta", fMaxDelta, 'f'));
1957  params->Add(new DP("TailDelta", fTailDelta, 'f'));
1958  params->Add(new DP("TailMax", fTailMax, 'f'));
1959  params->Add(new DP("AbsMinCent", fAbsMinCent, 'f'));
1960  // Create our centrality bins
1961  if (!InitCentBins(0)) {
1962  AliWarning("Failed to initialize centrality bins");
1963  return false;
1964  }
1965 
1966  // Print information to log
1967  Print();
1968  return true;
1969 }
1970 //____________________________________________________________________
1972 {
1973  if (!AliTrackletAODMCdNdeta::WorkerInit()) return false;
1974  fEtaWeight = 0;
1975  if (!fWeights) {
1976  AliFatal("No weights set!");
1977  return false;
1978  }
1979  fWeights->SetDebug(fDebug);
1980 
1981  TAxis wAxis(100,0,10);
1982  FixAxis(wAxis,"#it{w}_{i}");
1983  fEtaWeight = Make2P(fContainer, "etaWeight", "#LTw#GT", kYellow+2, 24,
1984  fEtaAxis, fCentAxis);
1985  fWeightCorr = Make2D(fContainer, "weightCorr", "w_{1} vs w_{2}",
1986  kCyan+2, 24, wAxis, wAxis);
1987  fWeights->Store(fContainer);
1988  return true;
1989 }
1990 //____________________________________________________________________
1993  Float_t c1,
1994  Float_t c2,
1995  Container* existing,
1996  const TAxis& deltaAxis)
1997 {
1998  CentBin* ret = MakeCentBin(c1, c2);
1999  Bool_t ok = true;
2000  if (!existing)
2001  ok = ret->WorkerInit(fContainer,fEtaAxis,fIPzAxis,deltaAxis);
2002  else
2003  ok = ret->FinalizeInit(existing);
2004  if (!ok) {
2005  AliWarningF("Failed to initialize bin %s", ret->GetName());
2006  delete ret;
2007  return 0;
2008  }
2009  ret->SetDebug(DebugLevel());
2010  fCentBins->AddAt(ret, bin);
2011  return ret;
2012 
2013 }
2014 //____________________________________________________________________
2016 {
2017  if (DebugLevel() > 1)
2018  Printf("Initialising on centrality bins on %s",
2019  existing ? "master" : "worker");
2020  if (fCentBins) return true;
2021  fCentBins = new Container;
2022  fCentBins->SetName("centralityBins");
2023  fCentBins->SetOwner();
2024 
2025  TAxis deltaAxis (Int_t(5*fMaxDelta),0,fMaxDelta);
2026  FixAxis(deltaAxis,
2027  "#Delta=[(#Delta#phi-#delta#phi)/#sigma_{#phi}]^{2}+"
2028  "[#Delta#thetasin^{-2}(#theta)/#sigma_{#theta}]^{2}");
2029 
2030  // Add min-bias bin
2031  if (!InitCentBin(0,0,0,existing,deltaAxis)) return false;
2032 
2033  // Add other bins
2034  Int_t nCentBins = fCentAxis.GetNbins();
2035  for (Int_t i = 1; i <= nCentBins; i++) {
2036  Float_t c1 = fCentAxis.GetBinLowEdge(i);
2037  Float_t c2 = fCentAxis.GetBinUpEdge(i);
2038  if (!InitCentBin(i, c1, c2, existing, deltaAxis)) return false;
2039  }
2040  if (!InitCentBin(nCentBins+1, 0, 100, existing, deltaAxis)) return false;
2041  return true;
2042 }
2043 
2044 //____________________________________________________________________
2046  : Sub(""),
2047  fSubs(0),
2048  fLow(c1),
2049  fHigh(c2),
2050  fStatus(0),
2051  fIPz(0),
2052  fCent(0),
2053  fCentIPz(0),
2054  fMeasured(0),
2055  fInjection(0)
2056 {
2057  if (c1 >= c2)
2058  fName = "all";
2059  else
2061  fMeasured = MakeHistos("measured", kRed+2, 20,
2062  kMeasuredMask, // No requirements, just veto
2063  kMeasuredVeto);
2064  fInjection = MakeHistos("injected", kOrange+2, 21,
2065  kInjectedMask,
2066  kInjectedVeto);
2067  fSubs = new Container;
2068  fSubs->SetOwner(true);
2069  fSubs->Add(fMeasured);
2070  fSubs->Add(fInjection);
2071 }
2072 
2073 //____________________________________________________________________
2075  : AliTrackletAODdNdeta::CentBin(c1, c2),
2076  fCombinatorics(0),
2077  fDistinct(0),
2078  fPrimaries(0),
2079  fSecondaries(0),
2080  fGenerated(0)
2081 {
2082  // Combinatorics is everything that has the combinatorics bit
2083  // Primaries all with simulated bit, but not secondary or
2084  // combinatorics bit.
2085  // Secondaries are all those with the secondary bit set
2086  fCombinatorics = MakeHistos("combinatorics", kMagenta+2, 30,
2089  fDistinct = MakeHistos("distinct", kCyan+2, 27,
2090  kDistinctMask,
2091  kDistinctVeto);
2092  fPrimaries = MakeHistos("primaries", kGreen+2, 26,
2093  kPrimaryMask,
2094  kPrimaryVeto);
2095  fSecondaries = MakeHistos("secondaries", kBlue+2, 32,
2097  kSecondaryVeto);
2098  fGenerated = MakeHistos("generated", kGray+1, 28,
2100  kGeneratedVeto);
2101  fMeasured->fStyle = 24;
2102  fInjection->fStyle = 25;
2103  fSubs->Add(fCombinatorics);
2104  fSubs->Add(fDistinct);
2105  fSubs->Add(fPrimaries);
2106  fSubs->Add(fSecondaries);
2107  fSubs->AddAfter(fMeasured, fGenerated);
2108 }
2109 
2110 //____________________________________________________________________
2113  Color_t color,
2114  Style_t style,
2115  UShort_t mask,
2116  UShort_t veto)
2117 {
2118  return new Histos(name, color, style, mask, veto);
2119 }
2120 
2121 //____________________________________________________________________
2123  const TAxis& etaAxis,
2124  const TAxis& ipzAxis,
2125  const TAxis& deltaAxis)
2126 {
2127  if (fDebug > 0)
2128  Printf("Initializing centrality bin %s", fName.Data());
2129  if (!Sub::WorkerInit(parent, etaAxis, ipzAxis, deltaAxis)) return false;
2130 
2131  TAxis centAxis(20, fLow, fHigh);
2132  FixAxis(centAxis, "Centrality [%]");
2134  fCent = Make1D(fContainer,"cent","Centrality [%]",
2135  kMagenta+2,20,centAxis);
2136  fIPz = Make1D(fContainer,"ipz","IP_{#it{z}} [cm]",kRed+2,20,ipzAxis);
2137  fCentIPz = Make1P(fContainer,"centIpz","#LTc#GT vs IP_{#it{z}}",kPink+2,
2138  20,ipzAxis);
2139  fCentIPz->SetYTitle("#LTc#GT");
2140  TIter next(fSubs);
2141  Histos* h = 0;
2142  while ((h = static_cast<Histos*>(next()))) {
2143  if (!h ->WorkerInit(fContainer, etaAxis,ipzAxis,deltaAxis))
2144  return false;
2145  }
2146 
2147  return true;
2148 }
2149 //____________________________________________________________________
2151  const TAxis& etaAxis,
2152  const TAxis& ipzAxis,
2153  const TAxis& deltaAxis)
2154 {
2155  if (!Sub::WorkerInit(parent, etaAxis, ipzAxis, deltaAxis)) return false;
2156  TString shrt(Form("%c", GetName()[0]));
2157  shrt.ToUpper();
2158  if (GetC(parent, "generated", false) != 0) shrt.Append("'");
2159 
2160  // Do not make eta vs IPz for secondaries and primaries
2161  if (!IsPrimary() && !IsSecondary())
2162  fEtaIPz = Make2D(fContainer, "etaIPz", shrt.Data(),
2163  kRed+2, 20, etaAxis, ipzAxis);
2164  // Always make eta vs Delta distribution, except for MC truth
2165  if (!IsGenerated())
2166  fEtaDeltaIPz = Make3D(fContainer, "etaDeltaIPz",
2167  Form("#Delta_{%s}",shrt.Data()),
2168  kBlue+2, 21, etaAxis, deltaAxis, ipzAxis);
2169  if (!IsGenerated() &&
2170  // !IsPrimary() &&
2171  // !IsSecondary() &&
2172  !IsInjected() &&
2173  fStyle != 20) // Last condition to not make this for real data
2174  fEtaPdgIPz = Make3D(fContainer, "etaPdgIPz",
2175  "Parent particle type",
2176  kGreen+2, 22, etaAxis, PdgAxis(), ipzAxis);
2177 
2178  if (IsGenerated()) {
2179  fEtaPdg = Make2D(fContainer, "etaPdg",
2180  "Primary particle type",
2181  kYellow+2, 30, etaAxis, PdgAxis());
2182  TAxis ptAxis(100, 0, 5);
2183  ptAxis.SetTitle("#it{p}_{T}");
2184  FixAxis(ptAxis);
2185  fEtaPt = Make2D(fContainer, "etaPt",
2186  "Primary transverse momentum",
2187  kCyan+2, 30, etaAxis, ptAxis);
2188  }
2189  if (IsPrimary() || IsSecondary() || IsCombinatoric()) {
2190  fEtaDeltaPdg = Make3D(fContainer, "etaDeltaPdg",
2191  "#Delta by primary particle type",
2192  kMagenta, 22, etaAxis, deltaAxis, PdgAxis());
2193  }
2194 
2195  SetAttr(fColor, fStyle);
2196  return true;
2197 }
2198 //____________________________________________________________________
2200 {
2201  if (fEtaIPz) {
2202  fEtaIPz->SetMarkerStyle(m);
2203  fEtaIPz->SetMarkerColor(c);
2204  fEtaIPz->SetLineColor(c);
2205  fEtaIPz->SetFillColor(c);
2206  }
2207  if (fEtaDeltaIPz) {
2208  fEtaDeltaIPz->SetMarkerStyle(m);
2209  fEtaDeltaIPz->SetMarkerColor(c);
2210  fEtaDeltaIPz->SetLineColor(c);
2211  fEtaDeltaIPz->SetFillColor(c);
2212  }
2213  if (fEtaDeltaPdg) {
2214  fEtaDeltaPdg->SetMarkerStyle(m);
2215  fEtaDeltaPdg->SetMarkerColor(c);
2216  fEtaDeltaPdg->SetLineColor(c);
2217  fEtaDeltaPdg->SetFillColor(c);
2218  }
2219  if (fEtaPdgIPz) {
2220  fEtaPdgIPz->SetMarkerStyle(m);
2221  fEtaPdgIPz->SetMarkerColor(c);
2222  fEtaPdgIPz->SetLineColor(c);
2223  fEtaPdgIPz->SetFillColor(c);
2224  }
2225  if (fEtaPdg) {
2226  fEtaPdg->SetMarkerStyle(m);
2227  fEtaPdg->SetMarkerColor(c);
2228  fEtaPdg->SetLineColor(c);
2229  fEtaPdg->SetFillColor(c);
2230  }
2231  if (fEtaPt) {
2232  fEtaPt->SetMarkerStyle(m);
2233  fEtaPt->SetMarkerColor(c);
2234  fEtaPt->SetLineColor(c);
2235  fEtaPt->SetFillColor(c);
2236  }
2237 }
2238 
2239 //____________________________________________________________________
2240 void
2242 {
2243  if (DebugLevel() > 0) Printf("In user exec");
2244  Double_t cent = -1;
2245  const AliVVertex* ip = 0;
2246  TClonesArray* tracklets = 0;
2247  UInt_t status = CheckEvent(cent, ip, tracklets);
2248  if ((status & kOKEvent) != kOKEvent) {
2249  AliWarningF("Event didn't pass %f, %p, %p", cent, ip, tracklets);
2250  return;
2251  }
2252  if (DebugLevel() > 0) Printf("Got centrality=%f ipZ=%f %d tracklets",
2253  cent,
2254  (ip ? ip->GetZ() : -1000),
2255  tracklets->GetEntriesFast());
2256  ProcessEvent(status, cent, ip, tracklets);
2257 
2258  PostData(1,fContainer);
2259  fStatus->Fill(kCompleted);
2260 }
2261 
2262 //____________________________________________________________________
2264 {
2265  return (1 << (bin-1));
2266 }
2267 
2268 //____________________________________________________________________
2270 {
2271  TH1* status = new TH1F("status", "Status of task",
2272  kCompleted, .5, kCompleted+.5);
2273  status->SetMarkerSize(2);
2274  status->SetMarkerColor(kMagenta+2);
2275  status->SetLineColor(kMagenta+2);
2276  status->SetFillColor(kMagenta+2);
2277  status->SetFillStyle(1001);
2278  status->SetBarOffset(0.1);
2279  status->SetBarWidth(0.4);
2280  status->SetDirectory(0);
2281  status->SetStats(0);
2282  status->SetXTitle("Event have");
2283  status->SetYTitle("# Events");
2284  status->GetXaxis()->SetBinLabel(kAll, "Been seen");
2285  status->GetXaxis()->SetBinLabel(kEvent, "Event data");
2286  status->GetXaxis()->SetBinLabel(kTracklets, "Tracklets");
2287  status->GetXaxis()->SetBinLabel(kTrigger, "Trigger");
2288  status->GetXaxis()->SetBinLabel(kIP, "IP");
2289  status->GetXaxis()->SetBinLabel(kCentrality, "Centrality");
2290  status->GetXaxis()->SetBinLabel(kCompleted, "Completed");
2291  container->Add(status);
2292  return status;
2293 }
2294 
2295 //____________________________________________________________________
2297  const AliVVertex*& ip,
2298  TClonesArray*& tracklets)
2299 {
2300  UInt_t ret = Bin2Flag(kAll);
2301  // Count all events
2302  fStatus->Fill(kAll);
2303 
2304  // Check for event
2305  AliVEvent* event = InputEvent();
2306  if (!event) {
2307  AliWarning("No event");
2308  return ret;
2309  }
2310  ret |= Bin2Flag(kEvent);
2311  fStatus->Fill(kEvent);
2312 
2313  // Check if we have the tracklets
2314  tracklets = FindTracklets(event);
2315  if (!tracklets) return ret;
2316  ret |= Bin2Flag(kTracklets);
2317  fStatus->Fill(kTracklets);
2318 
2319  // Check for least number of tracklets in |Eta|<1
2320  if (fMinEta1 > 0) {
2321  Int_t nOK = 0;
2322  AliAODTracklet* tracklet = 0;
2323  TIter nextTracklet(tracklets);
2324  while ((tracklet = static_cast<AliAODTracklet*>(nextTracklet()))) {
2325  UShort_t signal = CheckTracklet(tracklet);
2326  if (signal && tracklet->IsMeasured() && TMath::Abs(tracklet->GetEta()) < 1) nOK++;
2327  }
2328  if (nOK < fMinEta1) return ret;
2329  }
2330  // Check if event was triggered
2331  Bool_t trg = FindTrigger();
2332  if (!trg) return ret;
2333  ret |= Bin2Flag(kTrigger);
2334  fStatus->Fill(kTrigger);
2335 
2336  // Check the interaction point
2337  ip = FindIP(event);
2338  if (!ip) return ret;
2339  ret |= Bin2Flag(kIP);
2340  fStatus->Fill(kIP);
2341 
2342  // Check the centrality
2343  Int_t nTracklets = 0;
2344  cent = FindCentrality(event, nTracklets);
2345  // Do not fail on missing centrality
2346  if (cent >= 0) {
2347  ret |= Bin2Flag(kCentrality);
2348  fStatus->Fill(kCentrality);
2349  }
2350 
2351  fIPz->Fill(ip->GetZ());
2352  fCent->Fill(cent);
2353  fCentTracklets->Fill(cent, nTracklets);
2354 
2355  return ret;
2356 }
2357 
2358 //____________________________________________________________________
2359 TClonesArray* AliTrackletAODdNdeta::FindTracklets(AliVEvent* event)
2360 {
2361  // Check the multiplicity
2362  TObject* obj = event->FindListObject(GetBranchName());
2363  if (!obj) {
2364  AliWarningF("Couldn't get object %s", GetBranchName());
2365  // event->GetList()->Print();
2366  return 0;
2367  }
2368  if (!obj->IsA()->InheritsFrom(TClonesArray::Class())) {
2369  AliWarningF("Object %s is not a TClonesArray but a %s",
2370  obj->GetName(), obj->ClassName());
2371  return 0;
2372  }
2373  return static_cast<TClonesArray*>(obj);
2374 }
2375 
2376 //____________________________________________________________________
2377 const AliVVertex* AliTrackletAODdNdeta::CheckIP(const AliVVertex* ip,
2378  Double_t maxDispersion,
2379  Double_t maxZError)
2380 {
2381  if (!ip) return 0;
2382  if (ip->GetNContributors() <= 0) {
2383  AliWarning("Not enough contributors for IP");
2384  return 0;
2385  }
2386  // If this is from the Z vertexer, do some checks
2387  if (ip->IsFromVertexerZ()) {
2388  // Get covariance matrix
2389  Double_t covar[6];
2390  ip->GetCovarianceMatrix(covar);
2391  Double_t sigmaZ = TMath::Sqrt(covar[5]);
2392  if (sigmaZ >= maxZError) {
2393  AliWarningF("IPz resolution = %f >= %f", sigmaZ, maxZError);
2394  return 0;
2395  }
2396 
2397  // If this IP doesn not derive from AliVertex, don't check dispersion.
2398  if (ip->IsA()->InheritsFrom(AliVertex::Class())) {
2399  const AliVertex* ipv = static_cast<const AliVertex*>(ip);
2400  // Dispersion is the parameter used by the vertexer for finding the IP.
2401  if (ipv->GetDispersion() >= maxDispersion) {
2402  AliWarningF("IP dispersion = %f >= %f",
2403  ipv->GetDispersion(), maxDispersion);
2404  return 0;
2405  }
2406  }
2407  }
2408 
2409  // If we get here, we either have a full 3D vertex or track
2410  // vertex, and we should check if it is in range
2411  if (ip->GetZ() < fIPzAxis.GetXmin() || ip->GetZ() > fIPzAxis.GetXmax()) {
2412  AliWarningF("IPz = %fcm out of range [%f,%f]cm",
2413  ip->GetZ(), fIPzAxis.GetXmin(), fIPzAxis.GetXmax());
2414  return 0;
2415  }
2416  return ip;
2417 }
2418 
2419 //____________________________________________________________________
2420 const AliVVertex* AliTrackletAODdNdeta::FindRealIP(AliVEvent* event,
2421  Double_t maxDispersion,
2422  Double_t maxZError)
2423 {
2424  const AliVVertex* ip = event->GetPrimaryVertex();
2425  if (!ip) {
2426  if (DebugLevel() > 1) AliWarning("No real IP for this event found!");
2427  return 0;
2428  }
2429  return CheckIP(ip, maxDispersion, maxZError);
2430 }
2431 
2432 //____________________________________________________________________
2433 const AliVVertex* AliTrackletAODdNdeta::FindSimpleIP(AliVEvent* event,
2434  Double_t maxDispersion,
2435  Double_t maxZError)
2436 {
2437  static AliVertex* ip;
2438  AliAODSimpleHeader* head =
2439  static_cast<AliAODSimpleHeader*>(event
2440  ->FindListObject("AliAODSimpleHeader"));
2441  if (!head) {
2442  if (DebugLevel() > 1) AliWarning("No simple header");
2443  return 0;
2444  }
2445  if (!ip)
2446  ip = new AliVertex;
2447  ip->SetXv(head->fRecIP.X());
2448  ip->SetYv(head->fRecIP.Y());
2449  ip->SetZv(head->fRecIP.Z());
2450  ip->SetNContributors(10);
2451  ip->SetTitle("");
2452  return CheckIP(ip, maxDispersion, maxZError);
2453 }
2454 
2455 //____________________________________________________________________
2456 const AliVVertex* AliTrackletAODdNdeta::FindIP(AliVEvent* event,
2457  Double_t maxDispersion,
2458  Double_t maxZError)
2459 {
2460  const AliVVertex* ip = FindRealIP(event,maxDispersion,maxZError);
2461  if (!ip) {
2462  ip = FindSimpleIP(event,maxDispersion,maxZError);
2463  if (!ip) {
2464  if (DebugLevel() > 1) AliWarning("No IP for this event found!");
2465  return 0;
2466  }
2467  }
2468  // Good vertex, return it
2469  return ip;
2470 }
2471 //____________________________________________________________________
2473 {
2474  UInt_t evBits = fInputHandler->IsEventSelected();
2475  Bool_t trgOK = (evBits & fOfflineTriggerMask || fOfflineTriggerMask == 0);
2476  if (DebugLevel())
2477  Printf("Trigger bits=0x%08x mask=0x%08x masked=0x%08x -> %s",
2478  evBits, fOfflineTriggerMask, evBits & fOfflineTriggerMask,
2479  trgOK ? "selected" : "rejected");
2480  return trgOK;
2481 }
2482 //____________________________________________________________________
2484 {
2485  static TString centMeth;
2486  if (centMeth.IsNull()) {
2487  centMeth = fCentMethod(3,fCentMethod.Length()-3);
2488  }
2489  AliCentrality* cent = event->GetCentrality();
2490  if (!cent) return -1;
2491 
2492  Double_t centPer = cent->GetCentralityPercentileUnchecked(centMeth);
2493  return centPer;
2494 }
2495 //____________________________________________________________________
2497  Double_t& value,
2498  Int_t& nTracklets,
2499  Int_t& nCl0,
2500  Int_t& nCl1)
2501 {
2502  AliMultSelection* cent =
2503  static_cast<AliMultSelection*>(event->FindListObject("MultSelection"));
2504  if (!cent) {
2505  AliWarning("No centrality in event");
2506  event->GetList()->Print();
2507  return -1;
2508  }
2509  if (fCentIdx < 0) {
2510  TString estName(fCentMethod);
2511  TString trkName("SPDTracklets");
2512  TString cl0Name("CL0");
2513  TString cl1Name("CL1");
2514  if (estName.EqualTo("fake", TString::kIgnoreCase)) estName = "V0M";
2515  for (Int_t i = 0; i < cent->GetNEstimators(); i++) {
2516  AliMultEstimator* e = cent->GetEstimator(i);
2517  if (!e) continue;
2518  if (estName.EqualTo(e->GetName())) fCentIdx = i;
2519  if (trkName.EqualTo(e->GetName())) fTrkIdx = i;
2520  else if (cl0Name.EqualTo(e->GetName())) fCl0Idx = i;
2521  else if (cl1Name.EqualTo(e->GetName())) fCl1Idx = i;
2522  }
2523  if (fTrkIdx < 0)
2524  AliWarning("Index for tracklet estimator not found!");
2525  if (fCl0Idx < 0)
2526  AliWarning("Index for layer 0 cluster estimator not found!");
2527  if (fCl1Idx < 0)
2528  AliWarning("Index for layer 1 cluster estimator not found!");
2529  }
2530  if (fCentIdx < 0) {
2531  AliWarningF("Centrality estimator %s not found", fCentMethod.Data());
2532  return -1;
2533  }
2534 
2535  // Use cached index to look up the estimator
2536  AliMultEstimator* estCent = cent->GetEstimator(fCentIdx);
2537  if (!estCent) {
2538  AliWarningF("Centrality estimator %s not available for this event",
2539  fCentMethod.Data());
2540  return -1;
2541  }
2542  value = estCent->GetValue();
2543 
2544  // Now look up tracklet value using cached index
2545  if (fTrkIdx >= 0) {
2546  AliMultEstimator* est = cent->GetEstimator(fTrkIdx);
2547  if (est) nTracklets = est->GetValue();
2548  }
2549  if (fCl0Idx >= 0) {
2550  AliMultEstimator* est = cent->GetEstimator(fCl0Idx);
2551  if (est) nCl0 = est->GetValue();
2552  }
2553  if (fCl1Idx >= 0) {
2554  AliMultEstimator* est = cent->GetEstimator(fCl1Idx);
2555  if (est) nCl1 = est->GetValue();
2556  }
2557 
2558  return estCent->GetPercentile();
2559 }
2560 //____________________________________________________________________
2562  Int_t& nTracklets)
2563 {
2564  Int_t nCl0 = 0, nCl1 = 0;
2565  Double_t value = 0;
2566  Double_t cent = FindMultCentrality(event, value, nTracklets, nCl0, nCl1);
2567  fCentEst ->Fill(cent, value);
2568  fCl0Cl1Tracklets->Fill(nCl0, nCl1, nTracklets);
2569  return cent;
2570 }
2571 //____________________________________________________________________
2573  Int_t& nTracklets)
2574 
2575 {
2576  Int_t nCl0 = 0, nCl1 = 0;
2577  Double_t value = 0;
2578  FindMultCentrality(event, value, nTracklets, nCl0, nCl1);
2579 
2580  AliAODTracklet* tracklet = 0;
2581  TClonesArray* tracklets = FindTracklets(event);
2582  TIter nextTracklet(tracklets);
2583  Int_t nTracklet = 0;
2584  while ((tracklet = static_cast<AliAODTracklet*>(nextTracklet()))) {
2585  if (tracklet->IsMeasured())
2586  nTracklet += LookupWeight(tracklet,0,0);
2587  }
2588  nTracklets = nTracklet;
2589  fCl0Cl1Tracklets->Fill(nCl0, nCl1, nTracklets);
2590 
2591  if (nTracklet > fMaxNTracklet) return -1;
2592 
2593  Double_t c = 100*(1-Float_t(nTracklet)/fMaxNTracklet);
2594 
2595  fCentEst->Fill(c, nTracklet);
2596  return c;
2597 }
2598 
2599 //____________________________________________________________________
2601  Int_t& nTracklets)
2602 {
2603  if (fCentMethod.EqualTo("MB", TString::kIgnoreCase)) {
2604  Printf("MB centrality - not checked");
2605  return 0;
2606  }
2607  Double_t centPer = -1;
2608 
2609  if (fCentMethod.EqualTo("fake",TString::kIgnoreCase))
2610  centPer = FindFakeCentrality(event, nTracklets);
2611  else if (fCentMethod.BeginsWith("OLD"))
2612  centPer = FindCompatCentrality(event);
2613  else
2614  centPer = FindMultCentrality(event, nTracklets);
2615  const Double_t safety = 1e-3;
2616  if (DebugLevel() > 1) Printf("Read centrality: %f%%", centPer);
2617  if (centPer < -safety) return -2;
2618  if (centPer < +safety) centPer = safety;
2619  else if (centPer > 100-safety) centPer = 100-safety;
2620 
2621  if (fAbsMinCent >= 0 && centPer < fAbsMinCent) {
2622  AliWarningF("Centrality = %f lower than absolute minimum [%f]",
2623  centPer, fAbsMinCent);
2624  return -3;
2625  }
2626  if (centPer < fCentAxis.GetXmin() || centPer > fCentAxis.GetXmax()) {
2627  AliWarningF("Centrality = %f out of range [%f,%f]",
2628  centPer, fCentAxis.GetXmin(), fCentAxis.GetXmax());
2629  return -3;
2630  }
2631  return centPer;
2632 }
2633 
2634 //____________________________________________________________________
2636  Double_t cent,
2637  Double_t ipz)
2638 {
2639  return 1;
2640 }
2641 //____________________________________________________________________
2643  Double_t cent,
2644  Double_t ipz)
2645 {
2646  // We don't check for weights, as we must have them to come this far
2647  // if (!fWeights) {
2648  // AliWarning("No weights defined");
2649  // return 1;
2650  // }
2651  Double_t w = fWeights->LookupWeight(tracklet, cent, ipz, fWeightCorr);
2652  if (tracklet->IsMeasured())
2653  fEtaWeight->Fill(tracklet->GetEta(), cent, w);
2654 
2655  if (fDebug > 4) {
2656  Bool_t ok=true;
2657  switch (TMath::Abs(tracklet->GetParentPdg())) {
2658  case 310: printf("K0s "); break;
2659  case 321: printf("K+- "); break;
2660  case 3122: printf("Lambda"); break;
2661  case 3112: printf("Sigma-"); break;
2662  // case 3212: printf("Sigma0"); break;// Old, wrong code
2663  case 3222: printf("Sigma+"); break;
2664  case 3312: printf("Xi- "); break;
2665  // case 3322: printf("Xi0 "); break;// Old, wrong code
2666  default: ok = false; break;
2667  }
2668  if (ok) {
2669  switch (TMath::Abs(tracklet->GetParentPdg(true))) {
2670  case 0: printf(" "); break;
2671  case 310: printf("-K0s "); break;
2672  case 321: printf("-K+- "); break;
2673  case 3122: printf("-Lambda"); break;
2674  case 3112: printf("-Sigma-"); break;
2675  // case 3212: printf("-Sigma0"); break;// Old, wrong code
2676  case 3222: printf("-Sigma+"); break;
2677  case 3312: printf("-Xi- "); break;
2678  // case 3322: printf("-Xi0 "); break;// Old, wrong code
2679  default: printf("-Other ");
2680  }
2681  printf(": %f", w); tracklet->Print();
2682  }
2683  }
2684 
2685  if (fDebug > 3) {
2686  printf("Looking up weight of tracklet -> %f ", w);
2687  tracklet->Print();
2688  }
2689  return w;
2690 }
2691 //____________________________________________________________________
2693 {
2694  Double_t dPhiS = (tracklet->GetDPhi() -
2695  TMath::Sign(fDPhiShift,Double_t(tracklet->GetDPhi())));
2696  UShort_t ret = 0;
2697  ret |= ((tracklet->GetDelta() <= fDeltaCut) ? 0x1 : 0);
2698  ret |= ((dPhiS < fShiftedDPhiCut) ? 0x2 : 0);
2699  // Printf("dPhiS=%f (%f) Delta=%f (%f) -> 0x%x",
2700  // dPhiS, fShiftedDPhiCut,
2701  // tracklet->GetDelta(), fDeltaCut,
2702  // ret);
2703  return ret;
2704 }
2705 
2706 
2707 //____________________________________________________________________
2709  Double_t cent,
2710  const AliVVertex* ip,
2711  TClonesArray* tracklets)
2712 {
2713  // Figure out which centrality bins to fill
2714  Double_t ipz = (ip ? ip->GetZ() : -1000);
2715  Int_t nAcc = 0;
2716  TIter nextAcc(fCentBins);
2717  CentBin* bin = 0;
2718  TList toRun;
2719  while ((bin = static_cast<CentBin*>(nextAcc()))) {
2720  // Not in range for this bin
2721  if (!bin->Accept(status, cent, ipz)) continue;
2722  toRun.Add(bin);
2723  nAcc++;
2724  }
2725  // If we have no centrality bins to fill, we return immediately
2726  if (nAcc <= 0) return;
2727 
2728  Double_t nBare = 0;
2729  Double_t nMeasured = 0;
2730  Double_t nGenerated = 0;
2731  Double_t nGood = 0;
2732  Double_t nFake = 0;
2733  AliAODTracklet* tracklet = 0;
2734  TIter nextTracklet(tracklets);
2735  while ((tracklet = static_cast<AliAODTracklet*>(nextTracklet()))) {
2736  Double_t weight = LookupWeight(tracklet, cent, ipz);
2737  UShort_t signal = CheckTracklet(tracklet);
2738  if (signal) fEtaPhi->Fill(tracklet->GetEta(), tracklet->GetPhi());
2739  if (tracklet->IsMeasured() && TMath::Abs(tracklet->GetEta()) < 0.7) {
2740  nMeasured += weight;
2741  nBare ++;
2742  if (tracklet->IsCombinatorics()) nFake += weight;
2743  else nGood += weight;
2744  }
2745  else if (tracklet->IsGenerated() &&
2746  !tracklet->IsNeutral() &&
2747  TMath::Abs(tracklet->GetEta()) <= 1) nGenerated ++; // += weight;
2748  TIter nextBin(&toRun);
2749  while ((bin = static_cast<CentBin*>(nextBin()))) {
2750  bin->ProcessTracklet(tracklet, ip->GetZ(), signal, weight);
2751  }
2752  }
2753  TIter nextBin(&toRun);
2754  while ((bin = static_cast<CentBin*>(nextBin()))) {
2755  bin->Completed();
2756  }
2757  fNBareVsFake ->Fill(nBare, nFake);
2758  fNBareVsGood ->Fill(nBare, nGood);
2759  fNTrackletVsFake ->Fill(nMeasured, nFake);
2760  fNTrackletVsGood ->Fill(nMeasured, nGood);
2761  fNGeneratedVsFake->Fill(nGenerated, nFake);
2762  fNGeneratedVsGood->Fill(nGenerated, nGood);
2763 }
2764 
2765 //____________________________________________________________________
2767 {
2768  return fLow >= fHigh;
2769 }
2770 //____________________________________________________________________
2772  Double_t cent,
2773  Double_t ipz)
2774 {
2775  Bool_t centOK = (IsAllBin() ||
2776  (status & Bin2Flag(kCentrality) &&
2777  cent >= fLow && cent < fHigh));
2778  // We get out here already as we're not in this centrality bin
2779  if (!centOK) return false;
2780  fStatus->Fill(kAll);
2781  if (status & Bin2Flag(kEvent)) fStatus->Fill(kEvent);
2782  if (status & Bin2Flag(kTracklets)) fStatus->Fill(kTracklets);
2783  if (status & Bin2Flag(kTrigger))
2784  fStatus->Fill(kTrigger);
2785  else
2786  return false; // In case we have no trigger
2787  if (status & Bin2Flag(kCentrality)) {
2788  fStatus->Fill(kCentrality);
2789  fCent->Fill(cent);
2790  }
2791  // We explicity check the IP status here, since we need that to
2792  // calculate the per-bin vertex efficiency
2793  if (status & Bin2Flag(kIP)) {
2794  fStatus ->Fill(kIP);
2795  fIPz ->Fill(ipz);
2796  fCentIPz->Fill(ipz, cent);
2797  return true;
2798  }
2799  return false;
2800 }
2801 
2802 //____________________________________________________________________
2804  Double_t ipZ,
2805  UShort_t signal,
2806  Double_t weight)
2807 {
2808  TIter next(fSubs);
2809  Histos* h = 0;
2810  if (fDebug > 3) tracklet->Print();
2811  while ((h = static_cast<Histos*>(next())))
2812  h->ProcessTracklet(tracklet, ipZ, signal, weight);
2813 
2814  return true;
2815 }
2816 //____________________________________________________________________
2818  Double_t ipZ,
2819  UShort_t signal,
2820  Double_t weight)
2821 {
2822  if (!fEtaIPz && !fEtaDeltaIPz) return true;
2823  // Get tracklet info
2824  Double_t eta = tracklet->GetEta();
2825  Double_t delta = tracklet->GetDelta();
2826  UChar_t flags = tracklet->GetFlags();
2827  Int_t pdgBin = -1;
2828  Int_t pdgBin2 = -1;
2829 
2830  // For debugging
2831  char m[7];
2832  if (fDebug > 3) Bits2String(flags, m);
2833 
2834  // Check the filter mask
2835  if (fMask != 0 && (flags & fMask) == 0) {
2836  if (fDebug > 3) Printf("%14s (0x%02x,----) %6s %7s (0x%02x) ",
2837  GetName(),fMask,"reject",m,flags);
2838  return false;
2839  }
2840 
2841  // If we need the PDG code, get that here
2842  if (fEtaPdgIPz || fEtaPdg || fEtaDeltaPdg) {
2843  pdgBin = PdgBin(tracklet->GetParentPdg());
2844  if (tracklet->IsCombinatorics())
2845  pdgBin2 = PdgBin(tracklet->GetParentPdg(true));
2846  }
2847 
2848  // If we have the eta-vs-pdg histogram, we should fill before the
2849  // veto, which filters out the neutral particles. We do however,
2850  // check if the particle was suppressed
2851  if (fEtaPdg && !(tracklet->IsSuppressed()))
2852  fEtaPdg->Fill(eta, pdgBin, weight);
2853 
2854  // Check the veto mask
2855  if (fVeto != 0 && (flags & fVeto) != 0) {
2856  if (fDebug > 3) Printf("%14s (----,0x%02x) %6s %7s (0x%02x) ",
2857  GetName(), fVeto, "veto", m, flags);
2858  return false;
2859  }
2860  // Debug message that we accepted tracklet
2861  if (fDebug > 3) Printf("%14s (0x%02x,0x%02x) %6s %7s (0x%02x) ",
2862  GetName(), fMask, fVeto, "accept", m, flags);
2863 
2864  // Fill PDG,eta dependent Delta distributions
2865  if (fEtaDeltaPdg) {
2866  fEtaDeltaPdg ->Fill(eta, delta, pdgBin, weight);
2867  // Fill both particles
2868  if (pdgBin2 >= 0) fEtaPdgIPz->Fill(eta, delta, pdgBin2, weight);
2869  }
2870 
2871  // Both reguirements (Delta < cut, dPhi < cut)
2872  if (fEtaIPz && (signal == 0x3)) fEtaIPz->Fill(eta, ipZ, weight);
2873  // Just check dPhi
2874  if (fEtaDeltaIPz && (signal & 0x2)) fEtaDeltaIPz->Fill(eta,delta,ipZ,weight);
2875 
2876  // If we do not need the eta-PDG-IPz or eta-Pt, return now
2877  if (!fEtaPdgIPz && !fEtaPt) return true;
2878 
2879  if (fEtaPdgIPz) {
2880  fEtaPdgIPz ->Fill(eta, pdgBin, ipZ, weight);
2881  // Fill both particles
2882  if (pdgBin2 >= 0) fEtaPdgIPz->Fill(eta, pdgBin2, ipZ, weight);
2883  }
2884 
2885  if (fEtaPt) fEtaPt ->Fill(eta, tracklet->GetParentPt(), weight);
2886 
2887  return true;
2888 }
2889 
2890 //____________________________________________________________________
2892 {
2893  fStatus->Fill(kCompleted);
2894 }
2895 
2896 //____________________________________________________________________
2897 void
2899 {
2900  Bool_t save = TH1::GetDefaultSumw2();
2901  TH1::SetDefaultSumw2();
2902 
2903  TString resName; resName.Form("%sResults",GetName());
2904  if (GetOutputData(2)) {
2905  Warning("Terminate", "Already have a result container, making a new one");
2906  resName.Append("New");
2907  }
2908  Container* results = new Container;
2909  results->SetName(resName);
2910  results->SetOwner();
2911 
2912  Print("");
2913  fContainer = static_cast<Container*>(GetOutputData(1));
2914  if (!fContainer) {
2915  AliWarning("No sum container found!");
2916  return;
2917  }
2918 
2919  if (!InitCentBins(fContainer)) {
2920  AliWarningF("Failed to initialize centrality bins from %s",
2921  fContainer->GetName());
2922  return;
2923  }
2924 
2925  if (!MasterFinalize(results)) {
2926  AliWarning("Failed to finalize results");
2927  return;
2928  }
2929 
2930  PostData(2, results);
2931  TH1::SetDefaultSumw2(save);
2932 }
2933 
2934 
2935 //____________________________________________________________________
2937 {
2938  fContainer = GetC(parent, fName);
2939  fStatus = GetH1(fContainer, "status");
2940  fCent = GetH1(fContainer, "cent");
2941  fIPz = GetH1(fContainer, "ipz");
2942  fCentIPz = GetP1(fContainer, "centIpz");
2943  if (!fContainer || !fStatus || !fCent || !fIPz) return false;
2944  TIter next(fSubs);
2945  Histos* h = 0;
2946  while ((h = static_cast<Histos*>(next())))
2947  if (!h->FinalizeInit(fContainer)) return false;
2948  return true;
2949 }
2950 //____________________________________________________________________
2952 {
2953  fContainer = GetC(parent, fName);
2954  fEtaIPz = GetH2(fContainer, "etaIPz", false); // No complaints
2955  fEtaDeltaIPz = GetH3(fContainer, "etaDeltaIPz", false);
2956  fEtaDeltaPdg = GetH3(fContainer, "etaDeltaPdg", false);
2957  fEtaPdgIPz = GetH3(fContainer, "etaPdgIPz", false);
2958  fEtaPdg = GetH2(fContainer, "etaPdg", false);
2959  fEtaPt = GetH2(fContainer, "etaPt", false);
2960  if (GetName()[0] == 'm' && GetC(parent,"generated")) {
2961  // Fix up titles
2962  if (fEtaIPz)
2963  fEtaIPz->SetTitle(Form("%s'", fEtaIPz->GetTitle()));
2964  if (fEtaDeltaIPz)
2965  fEtaDeltaIPz->SetTitle(Form("#Delta_{%s}", fEtaIPz->GetTitle()));
2966  }
2967  return (fContainer != 0); // && fEtaIPz != 0);
2968 }
2969 
2970 //____________________________________________________________________
2972 {
2973  // Make copies of histograms and store
2974  fIPz = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"ipz")));
2975  fCent = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"cent")));
2976  fStatus = static_cast<TH1*>(CloneAndAdd(results,GetH1(fContainer,"status")));
2977  fEtaPhi = static_cast<TH2*>(CloneAndAdd(results,GetH2(fContainer,"etaPhi")));
2978  fCentTracklets =
2979  static_cast<TProfile*>(CloneAndAdd(results,GetP(fContainer,
2980  "centTracklets")));
2981  fCentEst =
2982  static_cast<TProfile*>(CloneAndAdd(results,GetP(fContainer,
2983  "centEstimator")));
2984  typedef TParameter<double> DP;
2985  results->Add(new DP("triggerEfficiency", fTriggerEff));
2986 
2987  TString tstr("UNKNOWN");
2988  if (fOfflineTriggerMask & AliVEvent::kMB) {
2989  tstr = "MBOR ";
2990  if (fTriggerEff > 0 && fTriggerEff < 1) {
2991  tstr = "INEL";
2992  if (fMinEta1 > 0) tstr.Append(Form(">%d",fMinEta1-1));
2993  }
2994  }
2995  else if (fOfflineTriggerMask & AliVEvent::kINT7) {
2996  tstr = "V0AND ";
2997  if (fTriggerEff > 0 && fTriggerEff < 1) tstr = "NSD";
2998  }
2999  else if (fOfflineTriggerMask & AliVEvent::kINT5)
3000  tstr = "MBAND";
3001  else if (fOfflineTriggerMask == AliVEvent::kAny)
3002  tstr = "OFFLINE";
3003  results->Add(new TNamed("trigger", tstr.Data()));
3004 
3005 
3006  Double_t nEvents = fIPz->GetEntries();
3007  Printf("Event summary:");
3008  for (Int_t i = 1; i <= fStatus->GetNbinsX(); i++)
3009  Printf(" %10d %s",
3010  Int_t(fStatus->GetBinContent(i)),
3011  fStatus->GetXaxis()->GetBinLabel(i));
3012  for (Int_t i = 1; i <= fCent->GetNbinsX(); i++)
3013  Printf(" %6.2f-%6.2f%%: %d",
3014  fCent->GetXaxis()->GetBinLowEdge(i),
3015  fCent->GetXaxis()->GetBinUpEdge(i),
3016  Int_t(fCent->GetBinContent(i)));
3017 
3018 
3019  fIPz ->Scale(1./nEvents);
3020  fCent ->Scale(1./fCent->GetEntries());
3021  fStatus->Scale(1./fStatus->GetBinContent(1));
3022  fEtaPhi->Scale(1./nEvents);
3023 
3024  if (fTailMax < 0) fTailMax = fMaxDelta;
3025  TIter next(fCentBins);
3026  CentBin* bin = 0;
3027  while ((bin = static_cast<CentBin*>(next()))) {
3028  if (!bin->MasterFinalize(results, 0, fTailDelta, fTailMax)) {
3029  AliWarningF("Failed to finalize %s", bin->GetName());
3030  return false;
3031  }
3032  }
3033  return true;
3034 }
3035 //____________________________________________________________________
3037 {
3038  if (!AliTrackletAODMCdNdeta::MasterFinalize(results)) return false;
3039 
3040  TObject* o = fContainer->FindObject("etaWeight");
3041  if (o && o->IsA()->InheritsFrom(TH1::Class())) {
3042  TH1* h = static_cast<TH1*>(o->Clone());
3043  h->SetDirectory(0);
3044  results->Add(h);
3045  }
3046  else {
3047  AliWarningF("Object %p (etaWeight) is not a TH1 or not found",o);
3048  }
3049  o = fContainer->FindObject("weightCorr");
3050  if (o && o->IsA()->InheritsFrom(TH1::Class())) {
3051  TH1* h = static_cast<TH1*>(o->Clone());
3052  h->SetDirectory(0);
3053  results->Add(h);
3054  }
3055  else {
3056  AliWarningF("Object %p (weightCorr) is not a TH1 or not found",o);
3057  }
3058  if (!fWeights) return true;
3059 
3060  if (!fWeights->Retrieve(fContainer)) return false;
3061  fWeights->Store(results);
3062 
3063  return true;
3064 }
3065 
3066 
3067 //____________________________________________________________________
3069  TH1* ,
3070  Double_t tailDelta,
3071  Double_t tailMax)
3072 {
3073  Container* result = new Container;
3074  result->SetName(fName);
3075  result->SetOwner(true);
3076  parent->Add(result);
3077 
3078  TH1* status = static_cast<TH1*>(CloneAndAdd(result, fStatus));
3079  status->Scale(1./status->GetBinContent(kAll));
3080 
3081  Double_t trigEff = status->GetBinContent(kTrigger);
3082  Double_t vtxEff = status->GetBinContent(kIP) / trigEff;
3083  typedef TParameter<double> DP;
3084 
3085  result->Add(new DP("triggerEfficiency", trigEff));
3086  result->Add(new DP("ipEfficiency", vtxEff));
3087 
3088  Double_t nEvents = fIPz->GetEntries();
3089  // Copy ipZ histogram and scale by number of events
3090  TH1* ipZ = static_cast<TH1*>(CloneAndAdd(result, fIPz));
3091  ipZ->Scale(1./nEvents);
3092 
3093  TH1* cent = static_cast<TH1*>(CloneAndAdd(result, fCent));
3094  cent->Scale(1./nEvents);
3095 
3096 
3097  CloneAndAdd(result, fCentIPz);
3098 
3099  Container* measCont = 0;
3100  Container* genCont = 0;
3101  TIter next(fSubs);
3102  Histos* h = 0;
3103  while ((h = static_cast<Histos*>(next()))) {
3104  if (h->GetMask() != AliAODTracklet::kInjection &&
3106  // For the anything but injection or MC-labels, we just finalize
3107  if (!h->MasterFinalize(result, fIPz, tailDelta, tailMax)) {
3108  AliWarningF("Failed to finalize %s/%s", GetName(), h->GetName());
3109  return false;
3110  }
3111  if (h->GetMask() == AliAODTracklet::kGenerated)
3112  // Set MC truth container
3113  genCont = GetC(result, h->GetName());
3114  if (h == fMeasured)
3115  measCont = GetC(result, h->GetName());
3116  continue;
3117  }
3118  if (!EstimateBackground(result, measCont, genCont, h, tailDelta, tailMax)) {
3119  AliWarningF("Failed to estimate Bg in %s/%s", GetName(), h->GetName());
3120  return false;
3121  }
3122  }
3123  return true;
3124 }
3125 
3126 //____________________________________________________________________
3127 Bool_t
3129  Int_t nEvents)
3130 {
3131  // Scale distribution of PDGs to number of events and bin width
3132  if (!fEtaPdg) return true;
3133 
3134 
3135  TH2* etaPdg = static_cast<TH2*>(fEtaPdg->Clone());
3136  etaPdg->SetDirectory(0);
3137  etaPdg->Scale(1. / nEvents, "width");
3138  result->Add(etaPdg);
3139 
3140  // Loop over PDG types and create 2D and 1D distributions
3141  TAxis* yaxis = etaPdg->GetYaxis();
3142  Int_t first = yaxis->GetFirst();
3143  Int_t last = yaxis->GetLast();
3144  THStack* pdgs = new THStack("all","");
3145  THStack* toPion = new THStack("toPion", "");
3146  THStack* toAll = new THStack("toAll", "");
3147  TH1* pion = 0;
3148  Container* pdgOut = new Container();
3149  pdgOut->SetName("mix");
3150  result->Add(pdgOut);
3151  pdgOut->Add(pdgs);
3152  pdgOut->Add(toPion);
3153  pdgOut->Add(toAll);
3154 
3155  TH1* all = static_cast<TH1*>(etaPdg->ProjectionX("total",
3156  1,etaPdg->GetNbinsY()));
3157  all->SetDirectory(0);
3158  all->SetName("total");
3159  all->SetTitle("All");
3160  all->SetYTitle(Form("d#it{N}_{%s}/d#eta", "all"));
3161  all->SetFillColor(kBlack);
3162  all->SetMarkerColor(kBlack);
3163  all->SetLineColor(kBlack);
3164  all->SetMarkerStyle(20);
3165  all->SetFillStyle(0);
3166  all->SetFillColor(0);
3167  all->Reset();
3168  pdgs->Add(all);
3169  for (Int_t i = 1; i <= etaPdg->GetNbinsY(); i++) {
3170  Int_t pdg = TString(yaxis->GetBinLabel(i)).Atoi();
3171  TString nme;
3172  Style_t sty;
3173  Color_t col;
3174  PdgAttr(pdg, nme, col, sty);
3175  if (pdg < 0) pdg = 0;
3176  if (pdg == 22) continue; // Ignore photons
3177 
3178  TH1* h1 = static_cast<TH1*>(etaPdg->ProjectionX(Form("h%d", pdg),i,i));
3179  if (h1->GetEntries() <= 0) continue; // Do not store if empty
3180  h1->SetDirectory(0);
3181  h1->SetName(Form("eta_%d", pdg));
3182  h1->SetTitle(nme);
3183  h1->SetYTitle(Form("d#it{N}_{%s}/d#eta", nme.Data()));
3184  h1->SetFillColor(col);
3185  h1->SetMarkerColor(col);
3186  h1->SetLineColor(col);
3187  h1->SetMarkerStyle(sty);
3188  h1->SetFillStyle(0);
3189  h1->SetFillColor(0);
3190  h1->SetBinContent(0,0);
3191  all->Add(h1);
3192  pdgs->Add(h1);
3193  switch (pdg) {
3194  case 321: h1->SetBinContent(0, 0.15); break; // PRC88,044910
3195  case 2212: h1->SetBinContent(0, 0.05); break; // PRC88,044910
3196  case 310: h1->SetBinContent(0, 0.075); break; // PRL111,222301
3197  case 3122: h1->SetBinContent(0, 0.018); break; // PRL111,222301
3198  case 3212: h1->SetBinContent(0, 0.0055); break; // NPA904,539
3199  case 3322: h1->SetBinContent(0, 0.005); break; // PLB734,409
3200  case 211: h1->SetBinContent(0, 1); break; // it self
3201  default: h1->SetBinContent(0, -1); break; // Unknown
3202  }
3203  pdgOut->Add(h1);
3204 
3205  if (pdg == 211) pion = h1;
3206  }
3207  if (!pdgs->GetHists()) return true;
3208 
3209  TIter next(pdgs->GetHists());
3210  TH1* tmp = 0;
3211  Double_t rmin = +1e9;
3212  Double_t rmax = -1e9;
3213  while ((tmp = static_cast<TH1*>(next()))) {
3214  if (tmp == all) continue;
3215  // Calculate ratio to all
3216  TH1* rat = static_cast<TH1*>(tmp->Clone());
3217  rat->Divide(all);
3218  rat->SetDirectory(0);
3219  rat->SetTitle(Form("%s / all", tmp->GetTitle()));
3220  toAll->Add(rat);
3221 
3222  if (tmp == pion) continue;
3223  Double_t r276 = tmp->GetBinContent(0);
3224  if (r276 < 0 || r276 >= 1) continue;
3225 
3226  // Calulate ratio to pions
3227  rat = static_cast<TH1*>(tmp->Clone());
3228  rat->Divide(pion);
3229  rat->SetTitle(Form("%s / %s", tmp->GetTitle(), pion->GetTitle()));
3230  rat->SetDirectory(0);
3231 
3232  TGraphErrors* g = new TGraphErrors(1);
3233  g->SetName(Form("%s_2760", rat->GetName()));
3234  g->SetTitle(Form("%s in #sqrt{s_{NN}}=2.76TeV", rat->GetTitle()));
3235  g->SetPoint(0,0,r276);
3236  g->SetPointError(0,.5,0);
3237  g->SetLineColor(rat->GetLineColor());
3238  g->SetLineStyle(rat->GetLineStyle());
3239  g->SetMarkerColor(rat->GetMarkerColor());
3240  g->SetMarkerStyle(rat->GetMarkerStyle());
3241  g->SetMarkerSize(1.5*rat->GetMarkerSize());
3242  rat->GetListOfFunctions()->Add(g,"p");
3243  rat->SetMaximum(TMath::Max(rat->GetMaximum(),r276));
3244  rat->SetMinimum(TMath::Min(rat->GetMinimum(),r276));
3245  rmin = TMath::Min(rat->GetMinimum(),rmin);
3246  rmax = TMath::Max(rat->GetMaximum(),rmax);
3247 
3248  toPion->Add(rat);
3249  }
3250  // toPion->SetMinimum(rmin);
3251  toPion->SetMaximum(1.1*rmax);
3252 
3253  return true;
3254 }
3255 
3256 //____________________________________________________________________
3257 Bool_t
3259  Int_t nEvents,
3260  Double_t tailDelta,
3261  Double_t tailMax,
3262  const TString& pre,
3263  const TString& tit)
3264 {
3265  Container* pdgOut = new Container();
3266  pdgOut->SetName(pre);
3267  result->Add(pdgOut);
3268  TAxis* zaxis = fEtaDeltaPdg->GetZaxis();
3269  Int_t first = zaxis->GetFirst();
3270  Int_t last = zaxis->GetLast();
3271  Bool_t isMid = TString(pre).EqualTo("mid");
3272  Int_t nX = fEtaDeltaPdg->GetNbinsX();
3273  Int_t nY = fEtaDeltaPdg->GetNbinsY();
3274  Int_t nZ = fEtaDeltaPdg->GetNbinsZ();
3275  Int_t l1 = fEtaDeltaPdg->GetXaxis()->FindBin(-1);
3276  Int_t r1 = fEtaDeltaPdg->GetXaxis()->FindBin(+1);
3277  Int_t b1 = (isMid ? l1 : 1);
3278  Int_t b2 = (isMid ? r1 : l1-1);
3279  Int_t b3 = (isMid ? -1 : r1+1);
3280  Int_t b4 = (isMid ? -1 : nX);
3281  THStack* stack = new THStack("all", tit);
3282  pdgOut->Add(stack);
3283 
3284  TH1* total = fEtaDeltaPdg->ProjectionY("totalMid",1, 1,1, nZ);
3285  total->SetDirectory(0);
3286  total->SetName("total");
3287  total->SetTitle("Total");
3288  total->SetYTitle(Form("d#it{N}_{%s}/d#Delta", "total"));
3289  total->SetFillColor(kBlack);
3290  total->SetMarkerColor(kBlack);
3291  total->SetLineColor(kBlack);
3292  total->SetMarkerStyle(24);
3293  total->SetFillStyle(0);
3294  total->SetFillColor(0);
3295  total->Reset();
3296  stack->Add(total);
3297 
3298  THStack* ratios = new THStack("ratios","");
3299  TH1* ratioSig = Make1D(0, "ratioSig", "Ratio to all", kGreen+1, 24, *zaxis);
3300  TH1* ratioBg = Make1D(0, "ratioBg", "Ratio to all", kRed+1, 25, *zaxis);
3301  ratioSig->GetXaxis()->LabelsOption("v");
3302  ratioBg ->GetXaxis()->LabelsOption("v");
3303  ratioSig->SetFillColor(kGreen+1);
3304  ratioBg ->SetFillColor(kRed +1);
3305  ratioSig->SetFillStyle(1001);
3306  ratioBg ->SetFillStyle(1001);
3307  ratioSig->SetBarWidth(0.4);
3308  ratioBg ->SetBarWidth(0.4);
3309  ratioSig->SetBarOffset(0.05);
3310  ratioBg ->SetBarOffset(0.55);
3311  ratios->Add(ratioSig, "bar0 e text30");
3312  ratios->Add(ratioBg, "bar0 e text30");
3313  pdgOut->Add(ratios);
3314 
3315  THStack* rows = new THStack("rows","");
3316  TH1* rowSig = new TH1D("rowSig","row",6,0,6);
3317  rowSig->SetDirectory(0);
3318  rowSig->SetYTitle("Fraction of signal");
3319  rowSig->GetXaxis()->SetBinLabel(1, "K^{0}_{S}");
3320  rowSig->GetXaxis()->SetBinLabel(2, "K^{#pm}");
3321  rowSig->GetXaxis()->SetBinLabel(3, "#Lambda");
3322  rowSig->GetXaxis()->SetBinLabel(4, "#Xi");
3323  rowSig->GetXaxis()->SetBinLabel(5, "#Sigma");
3324  rowSig->GetXaxis()->SetBinLabel(6, "Other");
3325  rowSig->GetXaxis()->LabelsOption("v");
3326  rowSig->SetMaximum(100);
3327  rowSig->SetMinimum(0);
3328  rowSig->SetLineColor(kGreen+1);
3329  rowSig->SetMarkerColor(kGreen+1);
3330  rowSig->SetMarkerStyle(24);
3331  TH1* rowBg = static_cast<TH1*>(rowSig->Clone("rowBg"));
3332  rowBg->SetDirectory(0);
3333  rowBg->SetYTitle("Fraction of background");
3334  rowBg->SetLineColor(kRed+1);
3335  rowBg->SetMarkerColor(kRed+1);
3336  rowBg->SetMarkerStyle(25);
3337  rowSig->SetFillStyle(1001);
3338  rowBg ->SetFillColor(1001);
3339  rowSig->SetFillColor(kGreen+1);
3340  rowBg ->SetFillColor(kRed +1);
3341  rowSig->SetBarWidth(0.4);
3342  rowBg ->SetBarWidth(0.4);
3343  rowSig->SetBarOffset(0.05);
3344  rowBg ->SetBarOffset(0.55);
3345  rows->Add(rowSig, "bar0 e text30");
3346  rows->Add(rowBg, "bar0 e text30");
3347  pdgOut->Add(rows);
3348 
3349  TH1* hK0S = 0;
3350  TH1* hKpm = 0;
3351  TH1* hLambda = 0;
3352  TH1* hXi = 0;
3353  TH1* hSigma = 0;
3354  TH1* hOther = 0;
3355  for (Int_t i = 1; i <= zaxis->GetNbins(); i++) {
3356  Int_t pdg = TString(zaxis->GetBinLabel(i)).Atoi();
3357  TString nme;
3358  Style_t sty;
3359  Color_t col;
3360  Int_t ipdg = pdg;
3361  PdgAttr(pdg, nme, col, sty);
3362  if (pdg < 0) pdg = 0;
3363  ratioSig->GetXaxis()->SetBinLabel(i, nme);
3364  ratioBg ->GetXaxis()->SetBinLabel(i, nme);
3365  if (pdg == 22) continue; // Ignore photons
3366 
3367 
3368  TH1* h1 = fEtaDeltaPdg->ProjectionY(Form("%d", pdg), b1, b2, i,i);
3369  h1->SetDirectory(0);
3370  if (b3 < b4) {
3371  TH1* h2 = fEtaDeltaPdg->ProjectionY(Form("t%d", pdg), b3, b4, i,i);
3372  h2->SetDirectory(0);
3373  h1->Add(h2);
3374  delete h2;
3375  }
3376  h1->SetUniqueID(i); // Store bin number
3377  if (h1->GetEntries() <= 0) continue; // Do not store if empty
3378  h1->SetName(Form("%d", pdg));
3379  h1->SetTitle(nme);
3380  h1->SetYTitle(Form("d#it{N}_{%s}/d#Delta", nme.Data()));
3381  h1->SetFillColor(col);
3382  h1->SetMarkerColor(col);
3383  h1->SetLineColor(col);
3384  h1->SetMarkerStyle(sty);
3385  h1->SetFillStyle(0);
3386  h1->SetFillColor(0);
3387  h1->SetBinContent(0,ipdg);
3388  total->Add(h1);
3389  stack->Add(h1);
3390 
3391  TH1** ptr = 0;
3392  switch (pdg) {
3393  case 321: ptr = &hKpm; break; // K^{+}
3394  case 310: ptr = &hK0S; break; // K^{0}_{S}
3395  // case 3212: // #Sigma^{0} // Old, wrong code
3396  case 3112: // #Sigma^{-}
3397  case 3222: ptr = &hSigma; break; // #Sigma^{+} // New, right code
3398 
3399  // case 3322: // #Xi^{0} // Old, wrong code
3400  case 3312: ptr = &hXi; break; // #Xi^{-} // New, right code
3401  case 3122: ptr = &hLambda;break; // #Lambda
3402  default: ptr = &hOther; break;
3403  }
3404  if (!*ptr) {
3405  *ptr = static_cast<TH1*>(h1->Clone());
3406  (*ptr)->Reset();
3407  (*ptr)->SetDirectory(0);
3408  if (ptr == &hOther) {
3409  (*ptr)->SetTitle("Other");
3410  (*ptr)->SetName("0");
3411  }
3412  }
3413  // Printf("Adding %s to %s", h1->GetName(), (*ptr)->GetName());
3414  (*ptr)->Add(h1);
3415  }
3416  total->SetBinContent(0,0);
3417  total->Scale(1./nEvents);
3418  Double_t deltaCut = 1.5;
3419  Double_t eTot, iTot = Integrate(total, 0, tailMax, eTot);
3420  Double_t eSig, iSig = Integrate(total, 0, deltaCut, eSig);
3421  Double_t eBg, iBg = Integrate(total, tailDelta, tailMax, eBg);
3422 
3423  TH1* totInt = new TH1D("totalIntegrals", "Integrals of total:", 3, 0, 3);
3424  totInt->SetDirectory(0);
3425  totInt->GetXaxis()->SetBinLabel(1,Form("0#minus%4.1f", tailMax));
3426  totInt->GetXaxis()->SetBinLabel(2,Form("0#minus%3.1f",deltaCut));
3427  totInt->GetXaxis()->SetBinLabel(3,Form("%3.1f#minus%4.1f",tailDelta,tailMax));
3428  totInt->SetBinContent(1,iTot); totInt->SetBinError(1, eTot);
3429  totInt->SetBinContent(2,iSig); totInt->SetBinError(2, eSig);
3430  totInt->SetBinContent(3,iBg); totInt->SetBinError(3, eBg);
3431  pdgOut->Add(totInt);
3432 
3433  if (!stack->GetHists()) {
3434  AliWarningF("No histograms in the stack %s", stack->GetName());
3435  return true;
3436  }
3437  TH1* spec[] = { hK0S, hKpm, hLambda, hXi, hSigma, hOther, 0 };
3438  TH1** ptr = spec;
3439  Int_t pbin = 1;
3440  for (Int_t pbin=1; pbin<=6; pbin++) {
3441  TH1* h = spec[pbin-1];
3442  if (!h) continue;
3443  h->Scale(1./ nEvents);
3444  Double_t eiSig, iiSig = Integrate(h, 0, deltaCut, eiSig);
3445  Double_t eiBg, iiBg = Integrate(h, tailDelta, tailMax, eiBg);
3446  Double_t erS, rS = RatioE(iiSig, eiSig, iSig, eSig, erS);
3447  Double_t erB, rB = RatioE(iiBg, eiBg, iBg, eBg, erB);
3448 #if 0
3449  Printf("%10s (%d) signal=%6.4f+/-%6.4f background=%6.4f+/-%6.4f "
3450  "ratios %6.4f+/-%6.4f %6.4f+/-%6.4f",
3451  h->GetTitle(), pbin, iiSig, eiSig, iiBg, eiBg, rS, erS, rB, erB);
3452 #endif
3453  rowSig->SetBinContent(pbin,100*rS);
3454  rowSig->SetBinError (pbin,100*TMath::Sqrt(erS));
3455  rowBg ->SetBinContent(pbin,100*rB);
3456  rowBg ->SetBinError (pbin,100*TMath::Sqrt(erB));
3457  pdgOut->Add(h);
3458  }
3459  TIter next(stack->GetHists());
3460  TH1* tmp = 0;
3461  while ((tmp = static_cast<TH1*>(next()))) {
3462  // Calculate ratio to all
3463  // Scale(tmp, iTot, eTot);
3464  if (tmp == total) continue;
3465  Int_t pdg = tmp->GetBinContent(0); tmp->SetBinContent(0,0);
3466  tmp->Scale(1./ nEvents);
3467  Double_t eiSig, iiSig = Integrate(tmp, 0, deltaCut, eiSig);
3468  Double_t eiBg, iiBg = Integrate(tmp, tailDelta, tailMax, eiBg);
3469  Double_t erS, rS = RatioE(iiSig, eiSig, iSig, eSig, erS);
3470  Double_t erB, rB = RatioE(iiBg, eiBg, iBg, eBg, erB);
3471  Int_t bin = tmp->GetUniqueID();
3472 #if 0
3473  Printf(" %10s(%4d,%3d) "
3474  "signal=%6.1f/%6.1f=%5.3f+/-%5.3f bg=%6.1f/%6.1f=%5.3f+/-%5.3f",
3475  tmp->GetTitle(), pdg, bin, iiSig, iSig, rS, erS, iiBg, iBg, rB, erB);
3476 #endif
3477  ratioSig->SetBinContent(bin, rS);
3478  ratioSig->SetBinError (bin, erS);
3479  ratioBg ->SetBinContent(bin, rB);
3480  ratioBg ->SetBinError (bin, erB);
3481  } // while(tmp)
3482 
3483  return true;
3484 }
3485 
3486 //____________________________________________________________________
3487 Bool_t
3489  Int_t nEvents,
3490  Double_t tailDelta,
3491  Double_t tailMax)
3492 {
3493  // Scale distribution of PDGs to number of events and bin width
3494  if (!fEtaDeltaPdg) return true;
3495 
3496  Container* pdgOut = new Container();
3497  pdgOut->SetName("specie");
3498  result->Add(pdgOut);
3499 
3500  ProjectEtaDeltaPdgPart(pdgOut,
3501  nEvents,
3502  tailDelta,
3503  tailMax,
3504  "mid",
3505  "|#eta|<1");
3506  ProjectEtaDeltaPdgPart(pdgOut,
3507  nEvents,
3508  tailDelta,
3509  tailMax,
3510  "fwd",
3511  "|#eta|>1");
3512 }
3513 
3514 //____________________________________________________________________
3515 Bool_t
3517  TH1* ipz,
3518  const TString& shn)
3519 {
3520  // If we have the PDG distributions, we project on eta,IPz, and then
3521  // average over IPz to get dNpdg/deta.
3522  if (!fEtaPdgIPz) return true;
3523 
3524  // Scale each vertex range by number of events in that range
3525  TH3* etaPdgIPz = ScaleToIPz(fEtaPdgIPz, ipz, false);
3526  result->Add(etaPdgIPz);
3527 
3528  // Loop over PDG types and create 2D and 1D distributions
3529  TAxis* yaxis = etaPdgIPz->GetYaxis();
3530  Int_t first = yaxis->GetFirst();
3531  Int_t last = yaxis->GetLast();
3532  THStack* pdgs = new THStack("all","");
3533  THStack* ratios = new THStack("toPion", "");
3534  TH1* pion = 0;
3535  TH2* dtfs = 0;
3536  Container* pdgOut = new Container();
3537  pdgOut->SetName("types");
3538  result->Add(pdgOut);
3539  pdgOut->Add(pdgs);
3540  pdgOut->Add(ratios);
3541  for (Int_t i = 1; i <= etaPdgIPz->GetNbinsY(); i++) {
3542  yaxis->SetRange(i,i);
3543 
3544  Int_t pdg = TString(yaxis->GetBinLabel(i)).Atoi();
3545  TString nme;
3546  Style_t sty;
3547  Color_t col;
3548  PdgAttr(pdg, nme, col, sty);
3549  if (pdg < 0) pdg = 0;
3550 
3551  TH2* h2 = static_cast<TH2*>(etaPdgIPz->Project3D("zx e"));
3552  if (h2->GetEntries() <= 0) continue; // Do not store if empty
3553  h2->SetDirectory(0);
3554  h2->SetName(Form("etaIPz_%d", pdg));
3555  h2->SetTitle(Form("%s#rightarrowX_{%s}", nme.Data(), shn.Data()));
3556  h2->SetYTitle(Form("d#it{N}^{2}_{%s#rightarrow%s}/"
3557  "(d#etadIP_{#it{z}})",
3558  nme.Data(), shn.Data()));
3559  h2->SetFillColor(col);
3560  h2->SetMarkerColor(col);
3561  h2->SetLineColor(col);
3562  pdgOut->Add(h2);
3563 
3564  TH1* h1 = AverageOverIPz(h2, Form("eta_%d",pdg), 1, ipz, 0, 0, false);
3565  h1->SetDirectory(0);
3566  h1->SetYTitle(Form("d#it{N}_{%s#rightarrow%s}/d#eta",
3567  nme.Data(), shn.Data()));
3568  h1->SetMarkerStyle(sty);
3569  pdgs->Add(h1);
3570 
3571  if (pdg == 211) pion = h1;
3572  TH2* tmp = 0;
3573  switch (pdg) {
3574  case 321: // Strange meson K^{+} break;
3575  case 323: // Strange meson K^{*+} break;
3576  case 310: // Strange meson K^{0}_{S} break;
3577  case 130: // Strange meson K^{0}_{L} break;
3578  case 311: // Strange meson K^{0} break;
3579  case 313: // Strange meson K^{*} break;
3580  case 221: // Strange meson #eta break;
3581  case 333: // Strange meson #varphi break;
3582  case 331: // Strange meson #eta' break;
3583  case 3112: // Strange baryon #Sigma^{-} break;
3584  case 3222: // Strange baryon #Sigma^{+} break;
3585  case 3114: // Strange baryon #Sigma^{*-} break;
3586  case 3224: // Strange baryon #Sigma^{*+} break;
3587  case 3312: // Strange baryon #Xi^{-} break;
3588  case 3314: // Strange baryon #Xi^{*-} break;
3589  case 3122: // Strange baryon #Lambda break;
3590  case 3212: // Strange baryon #Sigma^{0} break;
3591  case 3214: // Strange baryon #Sigma^{*0} break;
3592  case 3322: // Strange baryon #Xi^{0} break;
3593  case 3324: // Strange baryon #Xi^{*0} break;
3594  tmp = h2;
3595  break;
3596  default: break;
3597  }
3598  if (tmp) {
3599  if (!dtfs) {
3600  dtfs = static_cast<TH2*>(tmp->Clone("dtfs"));
3601  dtfs->SetDirectory(0);
3602  dtfs->SetTitle("X_{s}#rightarrowX");
3603  dtfs->SetYTitle("IP_{#it{z}} [cm]");
3604  dtfs->Reset();
3605  result->Add(dtfs);
3606  }
3607  // Add up contributions from strange particles
3608  dtfs->Add(tmp);
3609  }
3610  }
3611  if (!pdgs->GetHists()) {
3612  yaxis->SetRange(first, last);
3613  return true;
3614  }
3615 
3616  TIter next(pdgs->GetHists());
3617  TH1* tmp = 0;
3618  while ((tmp = static_cast<TH1*>(next()))) {
3619  if (tmp == pion) continue;
3620  TH1* rat = static_cast<TH1*>(tmp->Clone());
3621  rat->Divide(pion);
3622  rat->SetDirectory(0);
3623  ratios->Add(rat);
3624  }
3625 
3626  yaxis->SetRange(first, last);
3627 }
3628 
3629 //____________________________________________________________________
3630 Bool_t
3632  TH1* ipz,
3633  Double_t tailDelta,
3634  Double_t tailMax)
3635 {
3636  Container* result = new Container;
3637  result->SetName(fName);
3638  result->SetOwner(true);
3639  parent->Add(result);
3640 
3641  // Get the number of events
3642  Double_t nEvents = ipz->GetEntries();
3643 
3644  // Scale each vertex range by number of events in that range
3645  TH2* etaIPz = 0;
3646  if (fEtaIPz) {
3647  etaIPz = ScaleToIPz(fEtaIPz, ipz);
3648  result->Add(etaIPz);
3649  }
3650  // Scale distribution of Pt to number of events and bin width
3651  if (fEtaPt) {
3652  TH2* etaPt = static_cast<TH2*>(fEtaPt->Clone());
3653  etaPt->SetDirectory(0);
3654  etaPt->Scale(1. / nEvents, "width");
3655  result->Add(etaPt);
3656  }
3657  // Short-hand-name
3658  TString shn(etaIPz ? etaIPz->GetTitle() : "X");
3659 
3660  ProjectEtaPdg (result, nEvents);
3661  ProjectEtaDeltaPdg(result, nEvents, tailDelta, tailMax);
3662  ProjectEtaPdgIPz (result, ipz, shn);
3663 
3664  // If we do not have eta vs Delta, just return
3665  if (!fEtaDeltaIPz) return true;
3666 
3667  // Normalize delta distribution to integral number of events
3668  // static_cast<TH3*>(CloneAndAdd(result, fEtaDeltaIPz));
3669  TH3* etaDeltaIPz = ScaleToIPz(fEtaDeltaIPz, ipz, false);
3670  // ipz->GetEntries()>1000);
3671  result->Add(etaDeltaIPz);
3672 
3673  // Make 2D projection to eta,Delta
3674  TH2* etaDelta = ProjectEtaDelta(etaDeltaIPz);
3675  result->Add(etaDelta);
3676 
3677  // Make projection of delta
3678  TH1* delta = ProjectDelta(etaDelta);
3679  result->Add(delta);
3680 
3681  // PArameters of integrals
3682  Double_t maxDelta = etaDeltaIPz->GetYaxis()->GetXmax();
3683  Int_t lowBin = etaDeltaIPz->GetYaxis()->FindBin(tailDelta);
3684  Int_t sigBin = etaDeltaIPz->GetYaxis()->FindBin(1.5);
3685  Int_t highBin = TMath::Min(etaDeltaIPz->GetYaxis()->FindBin(tailMax),
3686  etaDeltaIPz->GetYaxis()->GetNbins());
3687 
3688 
3689  TH1* etaDeltaTail = etaDelta->ProjectionX("etaDeltaTail");
3690  etaDeltaTail ->SetDirectory(0);
3691  etaDeltaTail ->Reset();
3692  etaDeltaTail ->SetTitle(Form("#scale[.7]{#int}_{%3.1f}^{%4.1f}"
3693  "d%s d#it{N}/d%s",
3694  tailDelta,maxDelta,
3695  etaDeltaIPz->GetTitle(),
3696  etaDeltaIPz->GetTitle()));
3697  etaDeltaTail ->SetYTitle("#scale[.7]{#int}_{tail}d#Delta d#it{N}/d#Delta");
3698 
3699  TH2* etaIPzDeltaTail = static_cast<TH2*>(etaDeltaIPz->Project3D("zx e"));
3700  etaIPzDeltaTail->SetName("etaIPzDeltaTail");
3701  etaIPzDeltaTail->SetDirectory(0);
3702  etaIPzDeltaTail->Reset();
3703  etaIPzDeltaTail->SetTitle(etaDeltaTail->GetTitle());
3704  etaIPzDeltaTail->SetZTitle(etaDelta->GetYaxis()->GetTitle());
3705  TH2* etaIPzDeltaMain = static_cast<TH2*>(etaDeltaIPz->Project3D("zx e"));
3706  etaIPzDeltaMain->SetName("etaIPzDeltaMain");
3707  etaIPzDeltaMain->SetDirectory(0);
3708  etaIPzDeltaMain->Reset();
3709  etaIPzDeltaMain->SetTitle(etaDeltaTail->GetTitle());
3710  etaIPzDeltaMain->SetZTitle(etaDelta->GetYaxis()->GetTitle());
3711  // Loop over eta
3712  Double_t intg = 0, eintg = 0;
3713  for (Int_t i = 1; i <= etaDeltaTail->GetNbinsX(); i++) {
3714  // Integrate over Delta
3715  intg = etaDelta->IntegralAndError(i, i, lowBin, highBin, eintg);
3716  etaDeltaTail->SetBinContent(i, intg);
3717  etaDeltaTail->SetBinError (i, eintg);
3718  // Loop over IPz
3719  for (Int_t j = 1; j <= etaIPzDeltaTail->GetNbinsY(); j++) {
3720  // Integrate over Delta
3721  intg = etaDeltaIPz->IntegralAndError(i,i,lowBin,highBin,j,j,eintg);
3722  etaIPzDeltaTail->SetBinContent(i, j, intg);
3723  etaIPzDeltaTail->SetBinError (i, j, eintg);
3724 
3725  intg = etaDeltaIPz->IntegralAndError(i,i,1,sigBin,j,j,eintg);
3726  etaIPzDeltaMain->SetBinContent(i, j, intg);
3727  etaIPzDeltaMain->SetBinError (i, j, eintg);
3728  }
3729  }
3730  result->Add(etaIPzDeltaTail);
3731  result->Add(etaIPzDeltaMain);
3732  result->Add(etaDeltaTail);
3733 
3734  // Integrate full tail
3735  intg = etaDeltaIPz->IntegralAndError(1,etaDeltaIPz->GetNbinsX(),
3736  lowBin, highBin,
3737  1,etaDeltaIPz->GetNbinsZ(),
3738  eintg);
3739  result->Add(new TParameter<double>("deltaTail", intg));
3740  result->Add(new TParameter<double>("deltaTailError", eintg));
3741 
3742  // Some consistency checks:
3743  if (fDebug > 1) {
3744  Printf("%10s: Integral over eta,IPz: %9.4f +/- %9.4f",
3745  GetName(), intg, eintg);
3746  intg = etaDelta->IntegralAndError(1,etaDeltaIPz->GetNbinsX(),
3747  lowBin, highBin,
3748  eintg);
3749  Printf("%10s: Integral over eta: %9.4f +/- %9.4f",
3750  GetName(), intg, eintg);
3751  intg = delta->IntegralAndError(lowBin, highBin, eintg);
3752  Printf("%10s: Integral: %9.4f +/- %9.4f",
3753  GetName(), intg, eintg);
3754  }
3755 
3756  return true;
3757 }
3758 
3759 //____________________________________________________________________
3761  Container* measCont,
3762  Container* genCont,
3763  Histos* h,
3764  Double_t tailCut,
3765  Double_t tailMax)
3766 {
3767  if (!h || !measCont) {
3768  AliWarningF("No sub-histos or measured container in %s", GetName());
3769  return false;
3770  }
3771 
3772  if (!h->MasterFinalize(result, fIPz, tailCut,tailMax)) {
3773  AliWarningF("Failed to finalize %s/%s", GetName(), h->GetName());
3774  return false;
3775  }
3776 
3777  Container* bgCont = GetC(result, h->GetName());
3778  if (!bgCont) {
3779  AliWarningF("%s/%s didn't put a container on output",
3780  GetName(), h->GetName());
3781  return false;
3782  }
3783  const char* sub = h->GetName();
3784  TString shrt(Form("%c", sub[0]));
3785  shrt.ToUpper();
3786  if (genCont) shrt.Append("'");
3787 
3788  TH2* background = 0;
3789  TH2* etaIPzScale = 0;
3790  if (h->GetMask() == AliAODTracklet::kInjection) {
3791  // Get the tail ratio per eta, ipZ
3792  // We get the integral of the measured eta,IPz,Delta dist
3793  // We get the integral of the injected eta,IPz,Delta dist
3794  // Then we take the ratio
3795  etaIPzScale = CopyH2(measCont, "etaIPzDeltaTail",
3796  "etaIPzScale");
3797  etaIPzScale->Divide(GetH2(bgCont, "etaIPzDeltaTail"));
3798  etaIPzScale->SetZTitle("k_{#eta,IP_{#it{z}}}");
3799  etaIPzScale->SetTitle(Form("k_{%s,#eta,IP_{#it{z}}}",shrt.Data()));
3800  bgCont->Add(etaIPzScale);
3801 
3802  // Get the tail ratio per eta
3803  // We get the integral of the measured eta,Delta dist
3804  // We get the integral of the injected eta,Detla dist
3805  // Then we take the ratio
3806  TH1* etaScale = CopyH1(measCont, "etaDeltaTail",
3807  "etaScale");
3808  etaScale->Divide(GetH1(bgCont, "etaDeltaTail"));
3809  etaScale->SetYTitle("k_{#eta}");
3810  etaScale->SetTitle(Form("k_{%s,#eta}", shrt.Data()));
3811  bgCont->Add(etaScale);
3812 
3813  // Get the integrated tail ratio.
3814  // We get the integrated tail of the measured delta dist
3815  // We get the integrated tail of the ianjection delta dist
3816  // We then get the ratio.
3817  Double_t measIntg = GetD(measCont, "deltaTail", -1);
3818  Double_t measIntE = GetD(measCont, "deltaTailError", -1);
3819  Double_t bgIntg = GetD(bgCont, "deltaTail", -1);
3820  Double_t bgIntE = GetD(bgCont, "deltaTailError", -1);
3821  Double_t scaleE = 0;
3822  Double_t scale = RatioE(measIntg, measIntE, bgIntg, bgIntE, scaleE);
3823  bgCont->Add(new TParameter<double>("scale", scale));
3824  bgCont->Add(new TParameter<double>("scaleError", scaleE));
3825 
3826  // Get the fully differential Delta distribution and scale by the
3827  // eta,IPz scalar.
3828  TH3* scaledEtaDeltaIPz = ScaleDelta(CopyH3(bgCont, "etaDeltaIPz",
3829  "scaleEtaDeltaIPz"),
3830  etaIPzScale);
3831  scaledEtaDeltaIPz->SetTitle(Form("%5.3f#times%s",
3832  scale, scaledEtaDeltaIPz->GetTitle()));
3833  scaledEtaDeltaIPz->SetDirectory(0);
3834  scaledEtaDeltaIPz->SetYTitle("k#timesd^{3}#it{N}/"
3835  "(d#Deltad#etadIP_{#it{z}})");
3836  bgCont->Add(scaledEtaDeltaIPz);
3837 #if 0
3838  // scale by derived scalars, rather than by taking the scaled full
3839  // distribution.
3840  TH2* scaledEtaDelta = CopyH2(bgCont, "etaDelta", "scaledEtaDelta");
3841  scaledEtaDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3842  scaledEtaDelta->SetZTitle("k#timesd^{2}#it{N}/(d#Deltad#eta)");
3843  Scale(scaledEtaDelta, etaScale);
3844  bgCont->Add(scaledEtaDelta);
3845 
3846  TH1* scaledDelta = CopyH1(bgCont, "delta", "scaledDelta");
3847  scaledDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3848  scaledDelta->SetYTitle("k#timesd#it{N}/d#Delta");
3849  Scale(scaledDelta,scale,scaleE);
3850  bgCont->Add(scaledDelta);
3851 #else
3852  // Make 2D projection to eta,Delta
3853  TH2* scaledEtaDelta = ProjectEtaDelta(scaledEtaDeltaIPz);
3854  scaledEtaDelta->SetName("scaledEtaDelta");
3855  scaledEtaDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3856  scaledEtaDelta->SetYTitle("k#timesd^{2}#it{N}/(d#Deltad#eta)");
3857  bgCont->Add(scaledEtaDelta);
3858 
3859  // Make projection of delta
3860  TH1* scaledDelta = ProjectDelta(scaledEtaDelta);
3861  scaledDelta->SetName("scaledDelta");
3862  scaledDelta->SetTitle(scaledEtaDeltaIPz->GetTitle());
3863  scaledDelta->SetYTitle("k#timesd#it{N}/d#Delta");
3864  bgCont->Add(scaledDelta);
3865 #endif
3866  // Make background scaled by full tail
3867  background = CopyH2(bgCont, "etaIPz", "background");
3868  if (!background) AliWarningF("Didn't get background in %s", sub);
3869  else background->Multiply(etaIPzScale);
3870  // else Scale(background, scale, scaleE);
3871  }
3872  else {
3873  // For non-injection sets (i.e., combinatorics) we cannot form the
3874  // scaled background until we have the real data, so we calculate
3875  // beta instead.
3876  background = CopyH2(bgCont, "etaIPz", "background");
3877  TH2* beta = CopyH2(bgCont, "etaIPz", "beta");
3878  if (!background || !beta)
3879  AliWarningF("Didn't get background or beta in %s", sub);
3880  else {
3881  beta->Divide(GetH1(measCont, "etaIPz"));
3882  beta->SetTitle(Form("#beta_{%s}", shrt.Data()));
3883  bgCont->Add(beta);
3884  }
3885  }
3886  if (!background) {
3887  AliWarningF("Didn't get background in %s", sub);
3888  return false;
3889  }
3890  background->SetTitle(Form("#it{B}_{%s}", shrt.Data()));
3891  bgCont->Add(background);
3892 
3893  TH2* signal = CopyH2(measCont, "etaIPz", "signal");
3894  if (!signal) {
3895  AliWarningF("Didn't get signal in %s", sub);
3896  return false;
3897  }
3898  else {
3899  signal->SetTitle(Form("#it{S}_{%s}", shrt.Data()));
3900  signal->Add(background,-1);
3901  // Zero small bins
3902  for (Int_t i = 1; i <= signal->GetNbinsX(); i++) {
3903  for (Int_t j = 1; j <= signal->GetNbinsX(); j++) {
3904  if (signal->GetBinContent(i,j)<1e-6) {
3905  signal->SetBinContent(i,j,0);
3906  signal->SetBinError (i,j,0);
3907  }
3908  }
3909  }
3910  CopyAttr(background, signal);
3911  bgCont->Add(signal);
3912  }
3913 
3914  TH1* alpha = 0;
3915  if (genCont) {
3916  alpha = CopyH2(genCont, "etaIPz", "alpha");
3917  if (alpha && signal) {
3918  alpha->Divide(signal);
3919  alpha->SetTitle(Form("#alpha_{%s}", shrt.Data()));
3920  CopyAttr(signal, alpha);
3921  bgCont->Add(alpha);
3922  }
3923  }
3924 
3925  return true;
3926 }
3927 
3928 //====================================================================
3931  const char* weights,
3932  const char* sumFile,
3933  const char* resFile)
3934 {
3935  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
3936  if (!mgr) {
3937  ::Error("Create","No analysis manager to connect to.");
3938  return 0;
3939  }
3940  AliTrackletAODdNdeta* ret = 0;
3941  if (mc) {
3942  if (weights && weights[0] != '\0') {
3944  new AliTrackletAODWeightedMCdNdeta("MidRapidityMC");
3945  TUrl wurl(weights);
3946  TFile* wfile = TFile::Open(wurl.GetFile());
3947  if (!wfile) {
3948  ::Warning("Create", "Failed to open weights file: %s",
3949  wurl.GetUrl());
3950  return 0;
3951  }
3952  TString wnam(wurl.GetAnchor());
3953  if (wnam.IsNull()) wnam = "weights";
3954  TObject* wobj = wfile->Get(wnam);
3955  if (!wobj) {
3956  ::Warning("Create", "Failed to get weights %s from file %s",
3957  wnam.Data(), wfile->GetName());
3958  return 0;
3959  }
3960  if (!wobj->IsA()->InheritsFrom(AliTrackletBaseWeights::Class())) {
3961  ::Warning("Create", "Object %s from file %s not an "
3962  "AliTrackletBaseWeights but a %s",
3963  wnam.Data(), wfile->GetName(), wobj->ClassName());
3964  return 0;
3965  }
3966  wret->SetWeights(static_cast<AliTrackletBaseWeights*>(wobj));
3967  ret = wret;
3968  }
3969  else
3970  ret = new AliTrackletAODMCdNdeta("MidRapidityMC");
3971  }
3972  else ret = new AliTrackletAODdNdeta("MidRapidity");
3973  if (ret) ret->Connect();
3974 
3975  return ret;
3976 
3977 }
3978 
3979 
3980 
3981 
3982 //____________________________________________________________________
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)
Double_t FindMultCentrality(AliVEvent *event, Double_t &value, Int_t &nTracklets, Int_t &nCl0, Int_t &nCl1)
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
TH1D * hSigma
static void PrintAxis(const TAxis &axis, Int_t nSig=2, const char *alt=0)
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)
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 * 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)
Bool_t IsSuppressed() const
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)
ClassDef(AliTrackletAODdNdeta, 2)
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 void SetTriggerEfficiency(Double_t eff)
static TH1 * AverageOverIPz(TH2 *h, const char *name, UShort_t mode, TH1 *ipz, Double_t ipEff=0, TH2 *mask=0, Bool_t verb=false)
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)
void SetMaxNTracklet(Double_t maxN)
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
Double_t FindFakeCentrality(AliVEvent *event, Int_t &nTracklets)
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)
virtual void SetMinEta1(Int_t n)
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)