AliPhysics  ff1d528 (ff1d528)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTrackletAODUtils.C
Go to the documentation of this file.
1 
12 #ifndef ALITRACKLETAODUTILS_H
13 #define ALITRACKLETAODUTILS_H
14 #ifndef __CINT__
15 # include <TH1.h>
16 # include <TH2.h>
17 # include <TH3.h>
18 # include <TList.h>
19 # include <TParameter.h>
20 # include <TError.h>
21 # include <TMath.h>
22 # include <TFile.h>
23 # include <TDirectory.h>
24 # include <THashList.h>
25 # include <TProfile.h>
26 # include <TProfile2D.h>
27 # include <TDatabasePDG.h>
28 # include <THStack.h>
29 # include <cstdarg>
30 #else
31 class TList;
32 class TH1;
33 class TH2;
34 class TH3;
35 class TProfile;
36 class TProfile2D;
37 class TAxis;
38 class TDirectory;
39 class TFile;
40 class TDatabasePDG; // Auto-load
41 class THStack;
42 #endif
43 
57 {
58 public:
59  typedef TList Container;
67  virtual ~AliTrackletAODUtils() {}
68  //_________________________________________________________________________
82  static Bool_t CheckAxisNBins(const char* which,
83  const TAxis* a1,
84  const TAxis* a2);
94  static Bool_t CheckAxisLimits(const char* which,
95  const TAxis* a1,
96  const TAxis* a2);
107  static Bool_t CheckAxisBins(const char* which,
108  const TAxis* a1,
109  const TAxis* a2);
119  static Bool_t CheckAxisLabels(const char* which,
120  const TAxis* a1,
121  const TAxis* a2);
132  static Bool_t CheckAxis(const char* which,
133  const TAxis* a1,
134  const TAxis* a2,
135  Bool_t alsoLbls);
145  static Bool_t CheckConsistency(const TH1* h1, const TH1* h2);
146  /* @} */
147  //__________________________________________________________________
162  static TObject* GetO(Container* parent,
163  const char* name,
164  TClass* cls=0,
165  Bool_t verb=true);
176  static TObject* GetO(TDirectory* parent,
177  const char* name,
178  TClass* cls=0,
179  Bool_t verb=true);
189  static TH1* GetH1(Container* parent, const char* name, Bool_t verb=true);
199  static TH1* GetH1(TDirectory* parent, const char* name, Bool_t verb=true);
209  static TProfile* GetP1(Container* parent, const char* name, Bool_t verb=true);
219  static TH2* GetH2(Container* parent, const char* name, Bool_t verb=true);
229  static TH2* GetH2(TDirectory* parent, const char* name, Bool_t verb=true);
239  static TH3* GetH3(Container* parent, const char* name, Bool_t verb=true);
249  static TH3* GetH3(TDirectory* parent, const char* name, Bool_t verb=true);
259  static TProfile2D* GetP2(Container* parent, const char* name,
260  Bool_t verb=true);
270  static TProfile* GetP(Container* parent, const char* name, Bool_t verb=true)
271  {
272  return GetP1(parent, name, verb);
273  }
283  static THStack* GetHS(Container* parent, const char* name,Bool_t verb=true);
293  static THStack* GetHS(TDirectory* parent, const char* name,Bool_t verb=true);
303  static Container* GetC(Container* parent, const char* name, Bool_t verb=true);
313  static Container* GetC(TDirectory* parent, const char* name,Bool_t verb=true);
323  static TDirectory* GetT(TDirectory* parent,const char* name,Bool_t verb=true);
334  static Double_t GetD(Container* parent, const char* name,
335  Double_t def=-1, Bool_t verb=true);
346  static Int_t GetI(Container* parent, const char* name,
347  Int_t def=-1, Bool_t verb=true);
358  static Int_t GetB(Container* parent, const char* name,
359  Bool_t def=false, Bool_t verb=true);
370  static TH1* CopyH1(Container* parent,
371  const char* name,
372  const char* newName=0,
373  Bool_t verb=true);
384  static TH2* CopyH2(Container* parent,
385  const char* name,
386  const char* newName=0,
387  Bool_t verb=true);
398  static TH3* CopyH3(Container* parent,
399  const char* name,
400  const char* newName=0,
401  Bool_t verb=true);
408  static void CopyAttr(const TH1* src, TH1* dest);
421  static TH1* Make1D(Container* c,
422  const TString& name,
423  const TString& title,
424  Color_t color,
425  Style_t style,
426  const TAxis& xAxis);
439  static TProfile* Make1P(Container* c,
440  const TString& name,
441  const TString& title,
442  Color_t color,
443  Style_t style,
444  const TAxis& xAxis);
458  static TH2* Make2D(Container* c,
459  const TString& name,
460  const TString& title,
461  Color_t color,
462  Style_t style,
463  const TAxis& xAxis,
464  const TAxis& yAxis);
479  static TH3* Make3D(Container* c,
480  const TString& name,
481  const TString& title,
482  Color_t color,
483  Style_t style,
484  const TAxis& xAxis,
485  const TAxis& yAxis,
486  const TAxis& zAxis);
500  static TProfile2D* Make2P(Container* c,
501  const TString& name,
502  const TString& title,
503  Color_t color,
504  Style_t style,
505  const TAxis& xAxis,
506  const TAxis& yAxis);
507 
517  static TH1* Scale(TH1* h, Double_t x, Double_t xe);
527  static TH2* Scale(TH2* h, Double_t x, Double_t xe);
537  static TH2* Scale(TH2* h, TH1* s);
544  static void FixAxis(TAxis& axis, const char* title=0);
552  static void ScaleAxis(TAxis& ret, Double_t fact=1);
560  static void SetAxis(TAxis& axis, Int_t n, Double_t* borders);
569  static void SetAxis(TAxis& axis, const TString& spec, const char* sep=":,");
578  static void SetAxis(TAxis& axis, Int_t n, Double_t l, Double_t h);
586  static void SetAxis(TAxis& axis, Int_t n, Double_t m);
594  static void PrintAxis(const TAxis& axis, Int_t nSig=2, const char* alt=0);
605  static TH2* ScaleToIPz(TH2* h, TH1* ipZ, Bool_t full=false);
606  static TH3* ScaleToIPz(TH3* h, TH1* ipZ, Bool_t full=false);
607  static TH3* ScaleDelta(TH3* h, TH2* scale);
632  static TH2* ProjectEtaDelta(TH3* h);
648  static TH1* ProjectDelta(TH2* h);
664  static TH1* ProjectDeltaFull(TH3* h);
690  static TH1* AverageOverIPz(TH2* h,
691  const char* name,
692  UShort_t mode,
693  TH1* ipz,
694  Double_t ipEff=0,
695  TH2* mask=0,
696  Bool_t verb=false);
705  static TObject* CloneAndAdd(Container* c, TObject* o);
718  static Double_t Integrate(TH1* h, Double_t min, Double_t max, Double_t& err);
730  static Double_t RatioE(Double_t n, Double_t en,
731  Double_t d, Double_t ed,
732  Double_t& er);
744  static Double_t RatioE2(Double_t n, Double_t e2n,
745  Double_t d, Double_t e2d,
746  Double_t& e2r);
756  static TH1* RatioH(const TH1* num, const TH1* denom, const char* name=0);
763  static void FixMinMax(TH1* h, Bool_t ignoreZero=true);
770  static void FixMinMax(THStack* h, Bool_t ignoreZero=true);
771 
772  /* @} */
781  static void PdgAttr(Int_t pdg, TString& nme, Color_t& c, Style_t& s);
786  static Int_t* PdgArray(Int_t& size);
792  static const TAxis& PdgAxis();
800  static Int_t PdgBin(Int_t pdg);
801  /* @} */
810  static const char* CentName(Double_t c1, Double_t c2);
820  static Color_t CentColor(const TAxis& axis,
821  Double_t c1,
822  Double_t c2);
830  static Color_t CentColor(Int_t bin);
838  static TFile* OpenFile(const char* filename);
839 
840  //__________________________________________________________________
845  struct DebugGuard
846  {
856  DebugGuard(Int_t lvl, Int_t msgLvl, const char* format, ...);
860  ~DebugGuard();
868  static void Message(Int_t lvl, Int_t msgLvl, const char* format, ...);
869  private:
876  static void Output(int in, TString& msg);
884  static void Format(TString& out, const char* format, va_list ap);
886  };
887 protected:
888  ClassDef(AliTrackletAODUtils,1); // Utilities
889 };
890 
891 //====================================================================
892 // Utilities
893 //____________________________________________________________________
895  const TAxis* a1,
896  const TAxis* a2)
897 {
898  if (a1->GetNbins() != a2->GetNbins()) {
899  ::Warning("CheckAxisNBins", "Incompatible number %s bins: %d vs %d",
900  which, a1->GetNbins(), a2->GetNbins());
901  return false;
902  }
903  return true;
904 }
905 //____________________________________________________________________
907  const TAxis* a1,
908  const TAxis* a2)
909 {
910  if (TMath::AreEqualAbs(a1->GetXmin(), a2->GetXmin(), 1.e-3) &&
911  TMath::AreEqualAbs(a1->GetXmax(), a2->GetXmax(), 1.e-3)) return true;
912  if (!TMath::AreEqualRel(a1->GetXmin(), a2->GetXmin(),1.E-12) ||
913  !TMath::AreEqualRel(a1->GetXmax(), a2->GetXmax(),1.E-12)) {
914  Warning("CheckAxisLimits",
915  "Limits of %s axis incompatible [%f,%f] vs [%f,%f]", which,
916  a1->GetXmin(), a1->GetXmax(), a2->GetXmin(), a2->GetXmax());
917  return false;
918  }
919  return true;
920 }
921 //____________________________________________________________________
923  const TAxis* a1,
924  const TAxis* a2)
925 {
926  const TArrayD * h1Array = a1->GetXbins();
927  const TArrayD * h2Array = a2->GetXbins();
928  Int_t fN = h1Array->fN;
929  if ( fN == 0 ) return true;
930  if (h2Array->fN != fN) {
931  // Redundant?
932  Warning("CheckAxisBins", "Not equal number of %s bin limits: %d vs %d",
933  which, fN, h2Array->fN);
934  return false;
935  }
936  else {
937  for (int i = 0; i < fN; ++i) {
938  if (TMath::AreEqualAbs(h1Array->GetAt(i),
939  h2Array->GetAt(i), 1.e-3)) continue;
940  if (!TMath::AreEqualRel(h1Array->GetAt(i),h2Array->GetAt(i),1E-10)) {
941  Warning("CheckAxisBins",
942  "%s limit # %3d incompatible: %f vs %f",
943  which, i, h1Array->GetAt(i),h2Array->GetAt(i));
944  return false;
945  }
946  }
947  }
948  return true;
949 }
950 //____________________________________________________________________
952  const TAxis* a1,
953  const TAxis* a2)
954 {
955  // check that axis have same labels
956  THashList *l1 = (const_cast<TAxis*>(a1))->GetLabels();
957  THashList *l2 = (const_cast<TAxis*>(a2))->GetLabels();
958 
959  if (!l1 && !l2) return true;
960  if (!l1 || !l2) {
961  Warning("CheckAxisLabels", "Difference in %s labels: %p vs %p",
962  which, l1, l2);
963  return false;
964  }
965  // check now labels sizes are the same
966  if (l1->GetSize() != l2->GetSize()) {
967  Warning("CheckAxisLabels", "Different number of %s labels: %d vs %d",
968  which, l1->GetSize(), l2->GetSize());
969  return false;
970  }
971  for (int i = 1; i <= a1->GetNbins(); ++i) {
972  TString label1 = a1->GetBinLabel(i);
973  TString label2 = a2->GetBinLabel(i);
974  if (label1 != label2) {
975  Warning("CheckAxisLabels", "%s label # %d not the same: '%s' vs '%s'",
976  which, i, label1.Data(), label2.Data());
977  return false;
978  }
979  }
980 
981  return true;
982 }
983 //____________________________________________________________________
985  const TAxis* a1,
986  const TAxis* a2,
987  Bool_t alsoLbls)
988 {
989  if (a1->GetNbins() == 0 && a2->GetNbins() == 0) return true;
990  if (!CheckAxisNBins (which, a1, a2)) return false;
991  if (!CheckAxisLimits(which, a1, a2)) return false;
992  if (!CheckAxisBins (which, a1, a2)) return false;
993  if (alsoLbls && !CheckAxisLabels(which, a1, a2)) return false;
994  return true;
995 }
996 //____________________________________________________________________
998 {
999  // Check histogram compatibility
1000  if (h1 == h2) return true;
1001 
1002  if (h1->GetDimension() != h2->GetDimension() ) {
1003  Warning("CheckConsistency",
1004  "%s and %s have different dimensions %d vs %d",
1005  h1->GetName(), h2->GetName(),
1006  h1->GetDimension(), h2->GetDimension());
1007  return false;
1008  }
1009  Int_t dim = h1->GetDimension();
1010 
1011  Bool_t alsoLbls = (h1->GetEntries() != 0 && h2->GetEntries() != 0);
1012  if (!CheckAxis("X", h1->GetXaxis(), h2->GetXaxis(), alsoLbls)) return false;
1013  if (dim > 1 &&
1014  !CheckAxis("Y", h1->GetYaxis(), h2->GetYaxis(), alsoLbls)) return false;
1015  if (dim > 2 &&
1016  !CheckAxis("Z", h1->GetZaxis(), h2->GetZaxis(), alsoLbls)) return false;
1017 
1018  return true;
1019 }
1020 
1021 //____________________________________________________________________
1023 {
1024  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
1025  Double_t dc = s->GetBinContent(i);
1026  Double_t de = s->GetBinError (i);
1027  Double_t dr = (dc > 1e-6 ? de/dc : 0);
1028  for (Int_t j = 1; j <= h->GetNbinsY(); j++) {
1029  Double_t nc = h->GetBinContent(i,j);
1030  Double_t ne = h->GetBinError (i,j);
1031  Double_t ns = (nc > 0 ? ne/nc : 0);
1032  Double_t sc = dc * nc;
1033  Double_t se = sc*TMath::Sqrt(ns*ns+dr*dr);
1034  // Printf("Setting bin %2d,%2d to %f +/- %f", sc, se);
1035  h->SetBinContent(i,j,sc);
1036  h->SetBinError (i,j,se);
1037  }
1038  }
1039  return h;
1040 }
1041 
1042 //____________________________________________________________________
1044 {
1045  Double_t rr = xe/x;
1046  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
1047  for (Int_t j = 1; j <= h->GetNbinsY(); j++) {
1048  Double_t c = h->GetBinContent(i,j);
1049  Double_t e = h->GetBinError (i,j);
1050  Double_t s = (c > 0 ? e/c : 0);
1051  Double_t sc = x * c;
1052  Double_t se = sc*TMath::Sqrt(s*s+rr*rr);
1053  // Printf("Setting bin %2d,%2d to %f +/- %f", sc, se);
1054  h->SetBinContent(i,j,sc);
1055  h->SetBinError (i,j,se);
1056  }
1057  }
1058  return h;
1059 }
1060 
1061 //____________________________________________________________________
1063 {
1064  Double_t rr = xe/x;
1065  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
1066  Double_t c = h->GetBinContent(i);
1067  Double_t e = h->GetBinError (i);
1068  Double_t s = (c > 0 ? e/c : 0);
1069  Double_t sc = x * c;
1070  Double_t se = sc*TMath::Sqrt(s*s+rr*rr);
1071  // Printf("Setting bin %2d,%2d to %f +/- %f", sc, se);
1072  h->SetBinContent(i,sc);
1073  h->SetBinError (i,se);
1074  }
1075  return h;
1076 }
1077 
1078 //____________________________________________________________________
1080  const char* name,
1081  TClass* cls,
1082  Bool_t verb)
1083 {
1084  if (!parent) {
1085  if (verb) ::Warning("GetO", "No parent container passed");
1086  return 0;
1087  }
1088  TObject* o = parent->FindObject(name);
1089  if (!o) {
1090  if (verb) ::Warning("GetO", "Object \"%s\" not found in \"%s\"",
1091  name, parent->GetName());
1092  // parent->ls();
1093  return 0;
1094  }
1095  if (!cls) return o;
1096  if (!o->IsA()->InheritsFrom(cls)) {
1097  if (verb) ::Warning("GetO", "\"%s\" is an object of class %s, not %s",
1098  name, o->ClassName(), cls->GetName());
1099  return 0;
1100  }
1101  return o;
1102 }
1103 //____________________________________________________________________
1105  const char* name,
1106  TClass* cls,
1107  Bool_t verb)
1108 {
1109  if (!parent) {
1110  if (!verb) ::Warning("GetO", "No parent directory passed");
1111  return 0;
1112  }
1113  TObject* o = parent->Get(name);
1114  if (!o) {
1115  if (verb) ::Warning("GetO", "Object \"%s\" not found in \"%s\"",
1116  name, parent->GetName());
1117  // parent->ls();
1118  return 0;
1119  }
1120  if (!cls) return o;
1121  if (!o->IsA()->InheritsFrom(cls)) {
1122  if (verb) ::Warning("GetO", "\"%s\" is an object of class %s, not %s",
1123  name, o->ClassName(), cls->GetName());
1124  return 0;
1125  }
1126  return o;
1127 }
1128 //____________________________________________________________________
1129 TH1* AliTrackletAODUtils::GetH1(Container* parent, const char* name, Bool_t v)
1130 {
1131  return static_cast<TH1*>(GetO(parent, name, TH1::Class(), v));
1132 }
1133 //____________________________________________________________________
1134 TH1* AliTrackletAODUtils::GetH1(TDirectory* parent, const char* name, Bool_t v)
1135 {
1136  return static_cast<TH1*>(GetO(parent, name, TH1::Class(), v));
1137 }
1138 //____________________________________________________________________
1139 TH2* AliTrackletAODUtils::GetH2(Container* parent, const char* name, Bool_t v)
1140 {
1141  return static_cast<TH2*>(GetO(parent, name, TH2::Class(), v));
1142 }
1143 //____________________________________________________________________
1144 TH2* AliTrackletAODUtils::GetH2(TDirectory* parent, const char* name, Bool_t v)
1145 {
1146  return static_cast<TH2*>(GetO(parent, name, TH2::Class(), v));
1147 }
1148 //____________________________________________________________________
1149 TH3* AliTrackletAODUtils::GetH3(Container* parent, const char* name, Bool_t v)
1150 {
1151  return static_cast<TH3*>(GetO(parent, name, TH3::Class(), v));
1152 }
1153 //____________________________________________________________________
1154 TH3* AliTrackletAODUtils::GetH3(TDirectory* parent, const char* name, Bool_t v)
1155 {
1156  return static_cast<TH3*>(GetO(parent, name, TH3::Class(), v));
1157 }
1158 //____________________________________________________________________
1159 TProfile* AliTrackletAODUtils::GetP1(Container* parent,const char* name,
1160  Bool_t v)
1161 {
1162  return static_cast<TProfile*>(GetO(parent, name, TProfile::Class(), v));
1163 }
1164 //____________________________________________________________________
1165 TProfile2D* AliTrackletAODUtils::GetP2(Container* parent, const char* name,
1166  Bool_t v)
1167 {
1168  return static_cast<TProfile2D*>(GetO(parent, name, TProfile2D::Class(), v));
1169 }
1170 //____________________________________________________________________
1171 THStack*
1172 AliTrackletAODUtils::GetHS(TDirectory* parent, const char* name, Bool_t v)
1173 {
1174  return static_cast<THStack*>(GetO(parent, name, THStack::Class(), v));
1175 }
1176 //____________________________________________________________________
1177 THStack*
1178 AliTrackletAODUtils::GetHS(Container* parent, const char* name, Bool_t v)
1179 {
1180  return static_cast<THStack*>(GetO(parent, name, THStack::Class(), v));
1181 }
1182 //____________________________________________________________________
1184 AliTrackletAODUtils::GetC(Container* parent, const char* name, Bool_t v)
1185 {
1186  return static_cast<Container*>(GetO(parent, name, Container::Class(), v));
1187 }
1188 //____________________________________________________________________
1190 AliTrackletAODUtils::GetC(TDirectory* parent, const char* name, Bool_t v)
1191 {
1192  return static_cast<Container*>(GetO(parent, name, Container::Class(), v));
1193 }
1194 //____________________________________________________________________
1195 TDirectory*
1196 AliTrackletAODUtils::GetT(TDirectory* parent, const char* name, Bool_t v)
1197 {
1198  if (!parent) {
1199  if (v) ::Warning("GetT", "No parent directory passed");
1200  return 0;
1201  }
1202  TDirectory* d = parent->GetDirectory(name);
1203  if (!d) {
1204  if (v) ::Warning("GetT", "Directory \"%s\" not found in \"%s\"",
1205  name, parent->GetName());
1206  return 0;
1207  }
1208  return d;
1209 }
1210 //____________________________________________________________________
1212  const char* name,
1213  Double_t def,
1214  Bool_t v)
1215 {
1216  TParameter<double>* p =
1217  static_cast<TParameter<double>*>(GetO(parent, name,
1219  if (!p) return def;
1220  return p->GetVal();
1221 }
1222 //____________________________________________________________________
1224  const char* name,
1225  Int_t def,
1226  Bool_t v)
1227 {
1228  TParameter<int>* p =
1229  static_cast<TParameter<int>*>(GetO(parent, name,
1230  TParameter<int>::Class(),v));
1231  if (!p) return def;
1232  return p->GetVal();
1233 }
1234 //____________________________________________________________________
1236  const char* name,
1237  Bool_t def,
1238  Bool_t v)
1239 {
1240  TParameter<bool>* p =
1241  static_cast<TParameter<bool>*>(GetO(parent, name,
1242  TParameter<bool>::Class(), v));
1243  if (!p) return def;
1244  return p->GetVal();
1245 }
1246 //____________________________________________________________________
1248  const char* name,
1249  const char* newName,
1250  Bool_t v)
1251 {
1252  TH1* orig = GetH1(parent, name, v);
1253  if (!orig) return 0;
1254  TH1* ret = static_cast<TH1*>(orig->Clone(newName ? newName : name));
1255  ret->SetDirectory(0); // Release from file container
1256  return ret;
1257 }
1258 //____________________________________________________________________
1260  const char* name,
1261  const char* newName,
1262  Bool_t v)
1263 {
1264  TH2* orig = GetH2(parent, name, v);
1265  if (!orig) return 0;
1266  TH2* ret = static_cast<TH2*>(orig->Clone(newName ? newName : name));
1267  ret->SetDirectory(0); // Release from file container
1268  return ret;
1269 }
1270 //____________________________________________________________________
1272  const char* name,
1273  const char* newName,
1274  Bool_t v)
1275 {
1276  TH3* orig = GetH3(parent, name, v);
1277  if (!orig) return 0;
1278  TH3* ret = static_cast<TH3*>(orig->Clone(newName ? newName : name));
1279  ret->SetDirectory(0); // Release from file container
1280  return ret;
1281 }
1282 
1283 //____________________________________________________________________
1284 void AliTrackletAODUtils::CopyAttr(const TH1* src, TH1* dest)
1285 {
1286  if (!src || !dest) return;
1287  dest->SetMarkerStyle(src->GetMarkerStyle());
1288  dest->SetMarkerColor(src->GetMarkerColor());
1289  dest->SetMarkerSize (src->GetMarkerSize());
1290  dest->SetLineStyle (src->GetLineStyle());
1291  dest->SetLineColor (src->GetLineColor());
1292  dest->SetLineWidth (src->GetLineWidth());
1293  dest->SetFillStyle (src->GetFillStyle());
1294  dest->SetFillColor (src->GetFillColor());
1295 }
1296 
1297 
1298 //____________________________________________________________________
1300  const TString& name,
1301  const TString& title,
1302  Color_t color,
1303  Style_t style,
1304  const TAxis& xAxis)
1305 {
1306  TString n = name;
1307  TString t = title;
1308  TH1* ret = 0;
1309  Int_t nx = xAxis.GetNbins();
1310  if (t.IsNull()) t = xAxis.GetTitle();
1311  if (xAxis.GetXbins() && xAxis.GetXbins()->GetArray())
1312  ret = new TH1D(n,t,nx,xAxis.GetXbins()->GetArray());
1313  else
1314  ret = new TH1D(n,t,nx,xAxis.GetXmin(),xAxis.GetXmax());
1315  ret->Sumw2();
1316  ret->SetXTitle(xAxis.GetTitle());
1317  static_cast<const TAttAxis&>(xAxis).Copy(*(ret->GetXaxis()));
1318  ret->SetDirectory(0);
1319  ret->SetLineColor(color);
1320  ret->SetMarkerColor(color);
1321  ret->SetFillColor(kWhite);// color);
1322  ret->SetFillStyle(0);
1323  ret->SetMarkerStyle(style);
1324  if (const_cast<TAxis&>(xAxis).GetLabels()) {
1325  for (Int_t i = 1; i <= xAxis.GetNbins(); i++)
1326  ret->GetXaxis()->SetBinLabel(i, xAxis.GetBinLabel(i));
1327  }
1328  switch (style) {
1329  case 27:
1330  case 33:
1331  ret->SetMarkerSize(1.4);
1332  break;
1333  case 29:
1334  case 30:
1335  ret->SetMarkerSize(1.2);
1336  break;
1337  }
1338  if (c) c->Add(ret);
1339  return ret;
1340 }
1341 //____________________________________________________________________
1343  const TString& name,
1344  const TString& title,
1345  Color_t color,
1346  Style_t style,
1347  const TAxis& xAxis)
1348 {
1349  TString n = name;
1350  TString t = title;
1351  TProfile* ret = 0;
1352  Int_t nx = xAxis.GetNbins();
1353  if (t.IsNull()) t = xAxis.GetTitle();
1354  if (xAxis.GetXbins() && xAxis.GetXbins()->GetArray())
1355  ret = new TProfile(n,t,nx,xAxis.GetXbins()->GetArray());
1356  else
1357  ret = new TProfile(n,t,nx,xAxis.GetXmin(),xAxis.GetXmax());
1358  ret->Sumw2();
1359  ret->SetXTitle(xAxis.GetTitle());
1360  static_cast<const TAttAxis&>(xAxis).Copy(*(ret->GetXaxis()));
1361  ret->SetDirectory(0);
1362  ret->SetLineColor(color);
1363  ret->SetMarkerColor(color);
1364  ret->SetFillColor(kWhite);// color);
1365  ret->SetFillStyle(0);
1366  ret->SetMarkerStyle(style);
1367  if (const_cast<TAxis&>(xAxis).GetLabels()) {
1368  for (Int_t i = 1; i <= xAxis.GetNbins(); i++)
1369  ret->GetXaxis()->SetBinLabel(i, xAxis.GetBinLabel(i));
1370  }
1371  switch (style) {
1372  case 27:
1373  case 33:
1374  ret->SetMarkerSize(1.4);
1375  break;
1376  case 29:
1377  case 30:
1378  ret->SetMarkerSize(1.2);
1379  break;
1380  }
1381  if (c) c->Add(ret);
1382  return ret;
1383 }
1384 //____________________________________________________________________
1386  const TString& name,
1387  const TString& title,
1388  Color_t color,
1389  Style_t style,
1390  const TAxis& xAxis,
1391  const TAxis& yAxis)
1392 {
1393  TString n = name;
1394  TString t = title;
1395  TH2* ret = 0;
1396  Int_t nx = xAxis.GetNbins();
1397  Int_t ny = yAxis.GetNbins();
1398  const Double_t* xb = (xAxis.GetXbins() && xAxis.GetXbins()->GetArray() ?
1399  xAxis.GetXbins()->GetArray() : 0);
1400  const Double_t* yb = (yAxis.GetXbins() && yAxis.GetXbins()->GetArray() ?
1401  yAxis.GetXbins()->GetArray() : 0);
1402  if (t.IsNull())
1403  t.Form("%s vs %s", yAxis.GetTitle(), xAxis.GetTitle());
1404  if (xb) {
1405  if (yb) ret = new TH2D(n,t,nx,xb,ny,yb);
1406  else ret = new TH2D(n,t,
1407  nx,xb,
1408  ny,yAxis.GetXmin(),yAxis.GetXmax());
1409  }
1410  else {
1411  if (yb) ret = new TH2D(n,t,
1412  nx,xAxis.GetXmin(),xAxis.GetXmax(),
1413  ny,yb);
1414  else ret = new TH2D(n,t,
1415  nx,xAxis.GetXmin(),xAxis.GetXmax(),
1416  ny,yAxis.GetXmin(),yAxis.GetXmax());
1417  }
1418  ret->Sumw2();
1419  ret->SetXTitle(xAxis.GetTitle());
1420  ret->SetYTitle(yAxis.GetTitle());
1421  ret->SetLineColor(color);
1422  ret->SetMarkerColor(color);
1423  ret->SetFillColor(color);
1424  ret->SetMarkerStyle(style);
1425  static_cast<const TAttAxis&>(xAxis).Copy(*(ret->GetXaxis()));
1426  static_cast<const TAttAxis&>(yAxis).Copy(*(ret->GetYaxis()));
1427  ret->SetDirectory(0);
1428  if (const_cast<TAxis&>(xAxis).GetLabels()) {
1429  for (Int_t i = 1; i <= xAxis.GetNbins(); i++)
1430  ret->GetXaxis()->SetBinLabel(i, xAxis.GetBinLabel(i));
1431  }
1432  if (const_cast<TAxis&>(yAxis).GetLabels()) {
1433  for (Int_t i = 1; i <= yAxis.GetNbins(); i++)
1434  ret->GetYaxis()->SetBinLabel(i, yAxis.GetBinLabel(i));
1435  }
1436  if (c) c->Add(ret);
1437  return ret;
1438 }
1439 //____________________________________________________________________
1441  const TString& name,
1442  const TString& title,
1443  Color_t color,
1444  Style_t style,
1445  const TAxis& xAxis,
1446  const TAxis& yAxis,
1447  const TAxis& zAxis)
1448 {
1449  TString n = name;
1450  TString t = title;
1451  TH3* ret = 0;
1452  Int_t nx = xAxis.GetNbins();
1453  Int_t ny = yAxis.GetNbins();
1454  Int_t nz = zAxis.GetNbins();
1455  Double_t* xb = (xAxis.GetXbins() && xAxis.GetXbins()->GetArray() ?
1456  const_cast<Double_t*>(xAxis.GetXbins()->GetArray()) : 0);
1457  Double_t* yb = (yAxis.GetXbins() && yAxis.GetXbins()->GetArray() ?
1458  const_cast<Double_t*>(yAxis.GetXbins()->GetArray()) : 0);
1459  Double_t* zb = (zAxis.GetXbins() && zAxis.GetXbins()->GetArray() ?
1460  const_cast<Double_t*>(zAxis.GetXbins()->GetArray()) : 0);
1461  if (t.IsNull())
1462  t.Form("%s vs %s vs %s",
1463  zAxis.GetTitle(), yAxis.GetTitle(), xAxis.GetTitle());
1464  if (xb || yb || zb) {
1465  // One or more axis are defined as arrays. Make sure the rest are
1466  // also arrays.
1467  if (xb) {
1468  xb = new Double_t[nx+1];
1469  Double_t dx = (xAxis.GetXmax()-xAxis.GetXmin())/nx;
1470  xb[0] = xAxis.GetXmin();
1471  for (Int_t i = 1; i <= nx; i++) xb[i] = xb[i-1]+dx;
1472  }
1473  if (yb) {
1474  yb = new Double_t[ny+1];
1475  Double_t dy = (yAxis.GetXmax()-yAxis.GetXmin())/ny;
1476  yb[0] = yAxis.GetXmin();
1477  for (Int_t i = 1; i <= ny; i++) yb[i] = yb[i-1]+dy;
1478  }
1479  if (zb) {
1480  zb = new Double_t[nz+1];
1481  Double_t dz = (zAxis.GetXmax()-zAxis.GetXmin())/nz;
1482  zb[0] = zAxis.GetXmin();
1483  for (Int_t i = 1; i <= nz; i++) zb[i] = zb[i-1]+dz;
1484  }
1485  ret = new TH3D(n,t,nx,xb,ny,yb,nz,zb);
1486  }
1487  else {
1488  ret = new TH3D(n,t,
1489  nx, xAxis.GetXmin(), xAxis.GetXmax(),
1490  ny, yAxis.GetXmin(), yAxis.GetXmax(),
1491  nz, zAxis.GetXmin(), zAxis.GetXmax());
1492  }
1493  ret->Sumw2();
1494  ret->SetXTitle(xAxis.GetTitle());
1495  ret->SetYTitle(yAxis.GetTitle());
1496  ret->SetZTitle(zAxis.GetTitle());
1497  ret->SetLineColor(color);
1498  ret->SetMarkerColor(color);
1499  ret->SetFillColor(color);
1500  ret->SetMarkerStyle(style);
1501  static_cast<const TAttAxis&>(xAxis).Copy(*(ret->GetXaxis()));
1502  static_cast<const TAttAxis&>(yAxis).Copy(*(ret->GetYaxis()));
1503  static_cast<const TAttAxis&>(zAxis).Copy(*(ret->GetZaxis()));
1504  ret->SetDirectory(0);
1505  if (const_cast<TAxis&>(xAxis).GetLabels()) {
1506  for (Int_t i = 1; i <= xAxis.GetNbins(); i++)
1507  ret->GetXaxis()->SetBinLabel(i, xAxis.GetBinLabel(i));
1508  }
1509  if (const_cast<TAxis&>(yAxis).GetLabels()) {
1510  for (Int_t i = 1; i <= yAxis.GetNbins(); i++)
1511  ret->GetYaxis()->SetBinLabel(i, yAxis.GetBinLabel(i));
1512  }
1513  if (const_cast<TAxis&>(zAxis).GetLabels()) {
1514  for (Int_t i = 1; i <= zAxis.GetNbins(); i++)
1515  ret->GetZaxis()->SetBinLabel(i, zAxis.GetBinLabel(i));
1516  }
1517  if (c) c->Add(ret);
1518  return ret;
1519 }
1520 //____________________________________________________________________
1522  const TString& name,
1523  const TString& title,
1524  Color_t color,
1525  Style_t style,
1526  const TAxis& xAxis,
1527  const TAxis& yAxis)
1528 {
1529  TString n = name;
1530  TString t = title;
1531  TProfile2D* ret = 0;
1532  Int_t nx = xAxis.GetNbins();
1533  Int_t ny = yAxis.GetNbins();
1534  const Double_t* xb = (xAxis.GetXbins() && xAxis.GetXbins()->GetArray() ?
1535  xAxis.GetXbins()->GetArray() : 0);
1536  const Double_t* yb = (yAxis.GetXbins() && yAxis.GetXbins()->GetArray() ?
1537  yAxis.GetXbins()->GetArray() : 0);
1538  if (t.IsNull())
1539  t.Form("%s vs %s", yAxis.GetTitle(), xAxis.GetTitle());
1540  if (xb) {
1541  if (yb) ret = new TProfile2D(n,t,nx,xb,ny,yb);
1542  else ret = new TProfile2D(n,t,
1543  nx,xb,
1544  ny,yAxis.GetXmin(),yAxis.GetXmax());
1545  }
1546  else {
1547  if (yb) ret = new TProfile2D(n,t,
1548  nx,xAxis.GetXmin(),xAxis.GetXmax(),
1549  ny,yb);
1550  else ret = new TProfile2D(n,t,
1551  nx,xAxis.GetXmin(),xAxis.GetXmax(),
1552  ny,yAxis.GetXmin(),yAxis.GetXmax());
1553  }
1554  ret->Sumw2();
1555  ret->SetXTitle(xAxis.GetTitle());
1556  ret->SetYTitle(yAxis.GetTitle());
1557  ret->SetLineColor(color);
1558  ret->SetMarkerColor(color);
1559  ret->SetFillColor(color);
1560  ret->SetMarkerStyle(style);
1561  static_cast<const TAttAxis&>(xAxis).Copy(*(ret->GetXaxis()));
1562  static_cast<const TAttAxis&>(yAxis).Copy(*(ret->GetYaxis()));
1563  ret->SetDirectory(0);
1564  if (const_cast<TAxis&>(xAxis).GetLabels()) {
1565  for (Int_t i = 1; i <= xAxis.GetNbins(); i++)
1566  ret->GetXaxis()->SetBinLabel(i, xAxis.GetBinLabel(i));
1567  }
1568  if (const_cast<TAxis&>(yAxis).GetLabels()) {
1569  for (Int_t i = 1; i <= yAxis.GetNbins(); i++)
1570  ret->GetYaxis()->SetBinLabel(i, yAxis.GetBinLabel(i));
1571  }
1572  if (c) c->Add(ret);
1573  return ret;
1574 }
1575 //____________________________________________________________________
1576 void AliTrackletAODUtils::FixAxis(TAxis& axis, const char* title)
1577 {
1578  if (title && title[0] != '\0') axis.SetTitle(title);
1579  axis. SetNdivisions(210);
1580  axis. SetLabelFont(42);
1581  axis. SetLabelSize(0.03);
1582  axis. SetLabelOffset(0.005);
1583  axis. SetLabelColor(kBlack);
1584  axis. SetTitleOffset(1);
1585  axis. SetTitleFont(42);
1586  axis. SetTitleSize(0.04);
1587  axis. SetTitleColor(kBlack);
1588  axis. SetTickLength(0.03);
1589  axis. SetAxisColor(kBlack);
1590 }
1591 //____________________________________________________________________
1593  Double_t fact)
1594 {
1595  if (ret.GetXbins()->GetArray()) {
1596  TArrayD bins(*ret.GetXbins());
1597  for (Int_t i = 0; i < bins.GetSize(); i++) bins[i] *= fact;
1598  ret.Set(ret.GetNbins(), bins.GetArray());
1599  }
1600  else {
1601  ret.Set(ret.GetNbins(), fact*ret.GetXmin(), fact*ret.GetXmax());
1602  }
1603  FixAxis(ret);
1604 }
1605 //____________________________________________________________________
1607 {
1608  axis.Set(n, borders);
1609  FixAxis(axis);
1610 }
1611 //____________________________________________________________________
1613  const TString& spec,
1614  const char* sep)
1615 {
1616  TString s(spec);
1617  Bool_t isRange = false, isUnit = false;
1618  if (s[0] == 'r' || s[0] == 'R') {
1619  isRange = true;
1620  s.Remove(0,1);
1621  }
1622  if (s[0] == 'u' || s[0] == 'U') {
1623  isUnit = true;
1624  s.Remove(0,1);
1625  }
1626  TObjArray* tokens = s.Tokenize(sep);
1627  TArrayD bins(tokens->GetEntries());
1628  TObjString* token = 0;
1629  TIter next(tokens);
1630  Int_t i = 0;
1631  while ((token = static_cast<TObjString*>(next()))) {
1632  Double_t v = token->String().Atof();
1633  bins[i] = v;
1634  i++;
1635  }
1636  delete tokens;
1637  if (isUnit) {
1638  if (bins.GetSize() > 1)
1639  SetAxis(axis, Int_t(bins[1]-bins[0]), bins[0], bins[1]);
1640  else
1641  SetAxis(axis, 2*Int_t(bins[0]), bins[0]);
1642  }
1643  else if (isRange) {
1644  Int_t nBins = Int_t(bins[0]);
1645  if (bins.GetSize() > 2)
1646  SetAxis(axis, nBins, bins[1], bins[2]);
1647  else
1648  SetAxis(axis, nBins, bins[1]);
1649  }
1650  else
1651  SetAxis(axis, bins.GetSize()-1,bins.GetArray());
1652 }
1653 //____________________________________________________________________
1655  Int_t n,
1656  Double_t l,
1657  Double_t h)
1658 {
1659  axis.Set(n, l, h);
1660  FixAxis(axis);
1661 }
1662 //____________________________________________________________________
1664 {
1665  SetAxis(axis, n, -TMath::Abs(m), +TMath::Abs(m));
1666 }
1667 //____________________________________________________________________
1669  Int_t nSig,
1670  const char* alt)
1671 {
1672  printf(" %17s axis: ", alt ? alt : axis.GetTitle());
1673  if (axis.GetXbins() && axis.GetXbins()->GetArray()) {
1674  printf("%.*f", nSig, axis.GetBinLowEdge(1));
1675  for (Int_t i = 1; i <= axis.GetNbins(); i++)
1676  printf(":%.*f", nSig, axis.GetBinUpEdge(i));
1677  }
1678  else
1679  printf("%d bins between %.*f and %.*f",
1680  axis.GetNbins(), nSig, axis.GetXmin(),nSig,axis.GetXmax());
1681  printf("\n");
1682 }
1683 //____________________________________________________________________
1685 {
1686  if (!h) {
1687  ::Warning("ScaleToIPz","Nothing to scale");
1688  return 0;
1689  }
1690  if (!ipZ) {
1691  ::Warning("ScaleToIPz","Nothing to scale by");
1692  return 0;
1693  }
1694  TH2* ret = static_cast<TH2*>(h->Clone());
1695  ret->SetDirectory(0);
1696  if (!ipZ) return ret;
1697  for (Int_t iy = 1; iy <= ret->GetNbinsY(); iy++) {
1698  Double_t z = ret->GetYaxis()->GetBinCenter(iy);
1699  Int_t bin = ipZ->GetXaxis()->FindBin(z);
1700  Double_t nEv = ipZ->GetBinContent(bin);
1701  Double_t eEv = ipZ->GetBinError (bin);
1702  Double_t esc = (nEv > 0 ? 1./nEv : 0);
1703  Double_t rE2 = esc*esc*eEv*eEv;
1704  for (Int_t ix = 1; ix <= ret->GetNbinsX(); ix++) {
1705  Double_t c = ret->GetBinContent(ix,iy);
1706  Double_t e = ret->GetBinError (ix,iy);
1707  Double_t r = (c > 0 ? e/c : 0);
1708  // Scale by number of events, and error propagate
1709  Double_t sc = c * esc;
1710  Double_t se = 0;
1711  if (full) se = sc * TMath::Sqrt(r*r+rE2);
1712  else se = e * esc;
1713  Double_t scl = 1 / ret->GetXaxis()->GetBinWidth(ix);
1714  ret->SetBinContent(ix, iy, scl*sc);
1715  ret->SetBinError (ix, iy, scl*se);
1716  }
1717  }
1718  return ret;
1719 }
1720 //____________________________________________________________________
1722 {
1723  if (!h) {
1724  ::Warning("ScaleToIPz","Nothing to scale");
1725  return 0;
1726  }
1727  if (!ipZ) {
1728  ::Warning("ScaleToIPz","Nothing to scale by");
1729  return 0;
1730  }
1731  TH3* ret = static_cast<TH3*>(h->Clone());
1732  ret->SetDirectory(0);
1733  if (!ipZ) return ret;
1734  for (Int_t iz = 1; iz <= ret->GetNbinsZ(); iz++) {
1735  Double_t z = ret->GetZaxis()->GetBinCenter(iz);
1736  Int_t bin = ipZ->GetXaxis()->FindBin(z);
1737  Double_t nEv = ipZ->GetBinContent(bin);
1738  Double_t eEv = ipZ->GetBinError (bin);
1739  Double_t esc = (nEv > 0 ? 1./nEv : 0);
1740  Double_t rE2 = esc*esc*eEv*eEv;
1741  for (Int_t iy = 1; iy <= ret->GetNbinsY(); iy++) {
1742  for (Int_t ix = 1; ix <= ret->GetNbinsX(); ix++) {
1743  Double_t c = ret->GetBinContent(ix,iy,iz);
1744  Double_t e = ret->GetBinError (ix,iy,iz);
1745  Double_t r = (c > 0 ? e/c : 0);
1746  // Scale by number of events, and error propagate
1747  Double_t sc = c * esc;
1748  Double_t se = 0;
1749  if (full) se = sc * TMath::Sqrt(r*r+rE2);
1750  else se = e * esc;
1751  Double_t scl = 1 / ret->GetXaxis()->GetBinWidth(ix);
1752  ret->SetBinContent(ix, iy, iz, scl*sc);
1753  ret->SetBinError (ix, iy, iz, scl*se);
1754  }
1755  }
1756  }
1757  return ret;
1758 }
1759 //____________________________________________________________________
1761 {
1762  for (Int_t i = 1; i <= h->GetNbinsX(); i++) { // eta
1763  for (Int_t j = 1; j <= h->GetNbinsZ(); j++) { // IPz
1764  Double_t scale = etaIPzScale->GetBinContent(i,j);
1765  Double_t scaleE = etaIPzScale->GetBinError (i,j);
1766  Double_t q = (scale > 0 ? scaleE / scale : 0);
1767  for (Int_t k = 1; k <= h->GetNbinsY(); k++) { // delta
1768  Double_t c = h->GetBinContent(i,k,j);
1769  Double_t e = h->GetBinError (i,k,j);
1770  Double_t r = (c > 0 ? e / c : 0);
1771  Double_t w = c * scale;
1772  Double_t v = w * TMath::Sqrt(r*r + q*q);
1773  h->SetBinContent(i,k,j,w);
1774  h->SetBinError (i,k,j,v);
1775 #if 0
1776  Printf("%2d,%3d,%2d=(%9g+/-%9g)*(%9g+/-%9g)=(%9g+/-%9g)",
1777  i,k,j,c,e,scale,scaleE,w,v);
1778 #endif
1779  }
1780  }
1781  }
1782  return h;
1783 }
1784 //____________________________________________________________________
1785 // @cond
1787 {
1788  TH2* etaDelta = static_cast<TH2*>(h->Project3D("yx e"));
1789  etaDelta->SetName("etaDelta");
1790  etaDelta->SetTitle(h->GetTitle());
1791  etaDelta->SetDirectory(0);
1792  etaDelta->SetZTitle("d^{2}#it{N}/(d#etad#Delta)");
1793 #if 1
1794  // Reset content of projected histogram and calculate averages
1795  etaDelta->Reset();
1796  for (Int_t i = 1; i <= h->GetNbinsX(); i++) { // Loop over eta
1797  for (Int_t j = 1; j <= h->GetNbinsY(); j++) { // Loop over Delta
1798  Double_t sum = 0;
1799  Double_t sumw = 0;
1800  Int_t cnt = 0;
1801  for (Int_t k = 1; k <= h->GetNbinsZ(); k++) { // Loop over IPz
1802  Double_t c = h->GetBinContent(i,j,k);
1803  Double_t e = h->GetBinError (i,j,k);
1804  if (c < 1e-6 || e/c > 1) continue;
1805 #if 0
1806  Double_t w = 1/e/e;
1807  sum += w * c;
1808  sumw += w;
1809 #else
1810  sum += c;
1811  sumw += e*e;
1812 #endif
1813  cnt += 1;
1814  }
1815  if (sumw < 1e-6) continue;
1816 #if 0
1817  etaDelta->SetBinContent(i,j,sum/sumw);
1818  etaDelta->SetBinError (i,j,TMath::Sqrt(1/sumw));
1819 #else
1820  etaDelta->SetBinContent(i,j,sum/cnt);
1821  etaDelta->SetBinError (i,j,TMath::Sqrt(sumw)/cnt);
1822 #endif
1823  }
1824  }
1825 #endif
1826  return etaDelta;
1827 }
1828 // @endcond
1829 //____________________________________________________________________
1831 {
1832  if (!h) {
1833  Warning("ProjectDelta", "No histogram passed");
1834  return 0;
1835  }
1836  TH1* delta = h->ProjectionY("delta");
1837  delta->SetDirectory(0);
1838  delta->SetTitle(h->GetTitle());
1839  delta->SetYTitle("d#it{N}/d#Delta");
1840  delta->Scale(1./h->GetNbinsX());
1841  return delta;
1842 }
1843 //____________________________________________________________________
1845 {
1846  if (!h) {
1847  Warning("ProjectDelta", "No histogram passed");
1848  return 0;
1849  }
1850  TH2* tmp = ProjectEtaDelta(h);
1851  TH1* delta = ProjectDelta(tmp);
1852  delta->SetDirectory(0);
1853  tmp->SetDirectory(0);
1854  delete tmp;
1855  return delta;
1856 }
1857 
1858 //____________________________________________________________________
1860  const char* name,
1861  UShort_t mode,
1862  TH1* ipz,
1863  Double_t ipEff,
1864  TH2* other,
1865  Bool_t verb)
1866 {
1867  if (!h) return 0;
1868 
1869  Int_t nIPz = h->GetNbinsY();
1870  Int_t nEta = h->GetNbinsX();
1871  TH1* p = h->ProjectionX(name,1,nIPz,"e");
1872  TH2* mask = (other ? other : h);
1873  p->SetDirectory(0);
1874  p->SetFillColor(0);
1875  p->SetFillStyle(0);
1876  p->SetYTitle(Form("#LT%s#GT", h->GetYaxis()->GetTitle()));
1877  p->Reset();
1878  if (ipz && mode < 2 && ipEff > 1e-6) ipz->Scale(1./ipEff);
1879 
1880  for (Int_t etaBin = 1; etaBin <= nEta; etaBin++) {
1881  TArrayD hv(nIPz);
1882  TArrayD he(nIPz);
1883  TArrayD hr(nIPz);
1884  TArrayI hb(nIPz);
1885  Int_t j = 0;
1886  for (Int_t ipzBin = 1; ipzBin <= nIPz; ipzBin++) {
1887  Double_t bc = mask->GetBinContent(etaBin, ipzBin);
1888  if (bc < 1e-9) continue; // Low value
1889  Double_t be = mask->GetBinError (etaBin, ipzBin);
1890  if (TMath::IsNaN(be)) continue; // Bad error value
1891  hv[j] = bc;
1892  he[j] = be;
1893  hr[j] = be/bc;
1894  hb[j] = ipzBin;
1895  j++;
1896  }
1897  // Now we have all interesting values. Sort the relative error
1898  // array to get the most significant first. Note, we sort on the
1899  // mask errors which may not be the same as the histogram errors.
1900  TArrayI idx(nIPz);
1901  TMath::Sort(j, hr.fArray, idx.fArray, false);
1902  Double_t nsum = 0; // Running weighted sum
1903  Double_t nsume = 0; // Running weighted sum error
1904  Double_t dsum = 0;
1905  Double_t dsume = 0;
1906  Int_t n = 0;
1907  Double_t rat = 1e99;
1908 
1909  Int_t k = 0;
1910  for (k = 0; k < j; k++) {
1911  Int_t l = idx[k]; // Sorted index - ipzBin
1912  Int_t ipzBin = hb[l];
1913  Double_t hvv = hv[l];
1914  Double_t hee = he[l];
1915  Double_t x = TMath::Sqrt(nsume+hee*hee)/(nsum+hvv);
1916  if (x > rat) {
1917  continue; // Ignore - does not help
1918  }
1919 
1920  Double_t by = mask->GetYaxis()->GetBinCenter(ipzBin);
1921  Int_t ib = ipz ? ipz->FindBin(by) : 0;
1922  rat = x;
1923  nsum += h->GetBinContent(etaBin, ipzBin);
1924  nsume += TMath::Power(h->GetBinError(etaBin, ipzBin), 2);
1925  // If we do not have the vertex distribution, then just count
1926  // number of observations.
1927  dsum += !ipz ? 1 : ipz->GetBinContent(ib);
1928  dsume += !ipz ? 0 : TMath::Power(ipz->GetBinError(ib),2);
1929  n++;
1930  }
1931  if (k == 0 || n == 0) {
1932  if (verb)
1933  ::Warning("Average", "Eta bin # %3d of %10s has no data",
1934  etaBin, h->GetName());
1935  continue; // This eta bin empty!
1936  }
1937  Double_t norm = (mode > 0 ? n : dsum);
1938  Double_t rne = nsume/nsum/nsum;
1939  Double_t rde = dsume/dsum/dsum;
1940  Double_t av = nsum/norm;
1941  Double_t ave = 0;
1942  if (mode==2) ave = TMath::Sqrt(nsume)/n;
1943  else ave = av*TMath::Sqrt(rne+rde);
1944  p->SetBinContent(etaBin, av);
1945  p->SetBinError (etaBin, ave);
1946  }
1947  if (mode == 0) p->Scale(1, "width");
1948  return p;
1949 }
1950 
1951 //____________________________________________________________________
1953 {
1954  if (!o) {
1955  ::Warning("CloneAndAdd", "Nothing to clone");
1956  return 0;
1957  }
1958  TObject* copy = o->Clone();
1959  if (copy->IsA()->InheritsFrom(TH1::Class()))
1960  // Release from underlying directory
1961  static_cast<TH1*>(copy)->SetDirectory(0);
1962  if (c)
1963  // Add to output container
1964  c->Add(copy);
1965  return copy;
1966 }
1967 //____________________________________________________________________
1969  Double_t min,
1970  Double_t max,
1971  Double_t& err)
1972 {
1973  const Double_t epsilon = 1e-6;
1974  Int_t bMin = h->FindBin(min+epsilon);
1975  Int_t bMax = h->FindBin(max-epsilon);
1976  if (bMin < 1) bMin = 0; // Include underflow bin
1977  Double_t val = h->IntegralAndError(bMin, bMax, err);
1978  // For a-symmetric histograms, return
1979  if (TMath::Abs(h->GetXaxis()->GetXmin()+h->GetXaxis()->GetXmax())>=epsilon)
1980  return val;
1981 
1982  // Otherwise, also integrate negative side
1983  Double_t err2;
1984  bMin = h->FindBin(-min+epsilon);
1985  bMax = h->FindBin(-max-epsilon);
1986  val += h->IntegralAndError(bMin, bMax, err2);
1987  err = TMath::Sqrt(err*err+err2*err2);
1988  return val;
1989 }
1990 //____________________________________________________________________
1992  Double_t d, Double_t ed,
1993  Double_t& er)
1994 {
1995  Double_t r = 0;
1996  er = 0;
1997  if (TMath::Abs(n) < 1e-16 || TMath::Abs(d) < 1e-16) return 0;
1998  r = n/d;
1999  er = TMath::Sqrt(en*en/n/n + ed*ed/d/d);
2000  return r;
2001 }
2002 //____________________________________________________________________
2004  Double_t d, Double_t e2d,
2005  Double_t& e2r)
2006 {
2007  Double_t r = 0;
2008  e2r = 0;
2009  if (TMath::Abs(n) < 1e-16 || TMath::Abs(d) < 1e-16) return 0;
2010  r = n/d;
2011  e2r = (e2n/n/n + e2d/d/d);
2012  return r;
2013 }
2014 //____________________________________________________________________
2016  const TH1* denom,
2017  const char* name)
2018 {
2019  if (!num || !denom) {
2020  ::Warning("RatioH", "Numerator (%p) or denominator (%p) missing",
2021  num, denom);
2022  return 0;
2023  }
2024  if (!CheckConsistency(num, denom)) return 0;
2025 
2026  TH1* ratio = static_cast<TH1*>(num->Clone(name ? name : num->GetName()));
2027  ratio->SetDirectory(0);
2028  ratio->SetMinimum(-1111);
2029  ratio->SetMaximum(-1111);
2030  ratio->GetListOfFunctions()->Clear();
2031  ratio->Divide(denom);
2032  ratio->SetYTitle(Form("%s / %s",
2033  num ->GetYaxis()->GetTitle(),
2034  denom->GetYaxis()->GetTitle()));
2035  ratio->SetTitle(Form("%s / %s", num ->GetTitle(), denom->GetTitle()));
2036  FixMinMax(ratio, true);
2037 
2038  return ratio;
2039 }
2040 //____________________________________________________________________
2042 {
2043  Double_t min = +1e99;
2044  Double_t max = -1e99;
2045  for (Int_t ix = 1; ix <= h->GetNbinsX(); ix++) {
2046  for (Int_t iy = 1; iy <= h->GetNbinsY(); iy++) {
2047  for (Int_t iz = 1; iz <= h->GetNbinsZ(); iz++) {
2048  Int_t bin = h->GetBin(ix,iy,iz);
2049  Double_t c = h->GetBinContent(bin);
2050  Double_t e = h->GetBinError (bin);
2051  if (c == 0 && ignoreZero) continue;
2052  min = TMath::Min(min, c-e);
2053  max = TMath::Max(max, c+e);
2054  }
2055  }
2056  }
2057  h->SetMinimum(min);
2058  h->SetMaximum(max);
2059 }
2060 //____________________________________________________________________
2061 void AliTrackletAODUtils::FixMinMax(THStack* stack, Bool_t ignoreZero)
2062 {
2063  Double_t min = +1e99;
2064  Double_t max = -1e99;
2065  TIter next(stack->GetHists());
2066  TH1* h = 0;
2067  while ((h = static_cast<TH1*>(next()))) {
2068  FixMinMax(h, ignoreZero);
2069  min = TMath::Min(min, h->GetMinimum());
2070  max = TMath::Max(max, h->GetMaximum());
2071  }
2072  stack->SetMinimum(min);
2073  stack->SetMaximum(max);
2074 }
2075 
2076 //____________________________________________________________________
2078 {
2079  static Int_t codes[] = {
2080  211, // #pi^{+}
2081  2212, // p
2082  321, // K^{+}
2083  323, // K^{*+}
2084  11, // e^{-}
2085  13, // #mu^{-}
2086  213, // #rho^{+}
2087  411, // D^{+}
2088  413, // D^{*+}
2089  431, // D_{s}^{+}
2090  433, // D_{s}^{*+}
2091  1114, // #Delta^{-}
2092  2214, // #Delta^{+}
2093  2224, // #Delta^{++}
2094  3112, // #Sigma^{-}
2095  3222, // #Sigma^{+}
2096  3114, // #Sigma^{*-}
2097  3224, // #Sigma^{*+}
2098  4214, // #Sigma^{*+}_{c}
2099  4224, // #Sigma^{*++}_{c}
2100  3312, // #Xi^{-}
2101  3314, // #Xi^{*-}
2102  4122, // #Lambda^{+}_{c}
2103  2112, // n
2104  2114, // #Delta^{0}
2105  22, // #gamma
2106  310, // K^{0}_{S}
2107  130, // K^{0}_{L}
2108  311, // K^{0}
2109  313, // K^{*}
2110  221, // #eta
2111  111, // #pi^{0}
2112  113, // #rho^{0}
2113  333, // #varphi
2114  331, // #eta'
2115  223, // #omega
2116  3122, // #Lambda
2117  3212, // #Sigma^{0}
2118  4114, // #Sigma^{*0}_{c}
2119  3214, // #Sigma^{*0}
2120  421, // D^{0}
2121  423, // D^{*0}
2122  3322, // #Xi_{0}
2123  3324, // #Xi^{*0}
2124  4132, // #Xi^{0}_{c}
2125  4314, // #Xi^{*0}_{c}
2126  1000000000 // Nuclei
2127  };
2128  static Int_t nCodes = sizeof(codes) / sizeof(Int_t);
2129  static Int_t* sorted = 0;
2130  size = nCodes;
2131  if (sorted) return sorted;
2132 
2133  sorted = new Int_t[nCodes];
2134  Int_t* idx = new Int_t[nCodes];
2135  TMath::Sort(nCodes, codes, idx, false);
2136  for (Int_t i = 0; i < nCodes; i++) sorted[i] = codes[idx[i]];
2137  delete [] idx;
2138  return sorted;
2139 }
2140 
2141 //____________________________________________________________________
2143  Color_t& c, Style_t& s)
2144 {
2145  // Leptons are black
2146  // Non-strange, non-charm meson are red
2147  // Non-strange, non-charm baryons are magenta
2148  // Strange-mesons are blue
2149  // Strange-baryons are cyan
2150  // Charmed mesons are green
2151  // Charmed baryons are yellow
2152  switch (pdg) {
2153  case -1: c=kPink +7; s=20; nme="?"; break;
2154  case 11: c=kBlack +0; s=20; nme="e^{-}"; break;
2155  case 13: c=kBlack +0; s=21; nme="#mu^{-}"; break;
2156  case 22: c=kBlack +0; s=22; nme="#gamma"; break;
2157  case 111: c=kRed +2; s=20; nme="#pi^{0}"; break;
2158  case 211: c=kRed +2; s=21; nme="#pi^{+}"; break;
2159  case 213: c=kRed +2; s=22; nme="#rho^{+}"; break;
2160  case 113: c=kRed +2; s=23; nme="#rho^{0}"; break;
2161  case 223: c=kRed +2; s=24; nme="#omega"; break;
2162  case 321: c=kBlue +2; s=20; nme="K^{+}"; break;
2163  case 323: c=kBlue +2; s=21; nme="K^{*+}"; break;
2164  case 310: c=kBlue +2; s=22; nme="K^{0}_{S}"; break;
2165  case 130: c=kBlue +2; s=23; nme="K^{0}_{L}"; break;
2166  case 311: c=kBlue +2; s=24; nme="K^{0}"; break;
2167  case 313: c=kBlue +2; s=25; nme="K^{*}"; break;
2168  case 221: c=kBlue +2; s=26; nme="#eta"; break;
2169  case 333: c=kBlue +2; s=27; nme="#varphi"; break;
2170  case 331: c=kBlue +2; s=28; nme="#eta'"; break;
2171  case 411: c=kGreen +2; s=20; nme="D^{+}"; break;
2172  case 413: c=kGreen +2; s=21; nme="D^{*+}"; break;
2173  case 421: c=kGreen +2; s=22; nme="D^{0}"; break;
2174  case 423: c=kGreen +2; s=23; nme="D^{*0}"; break;
2175  case 431: c=kGreen +2; s=24; nme="D_{s}^{+}"; break;
2176  case 433: c=kGreen +2; s=25; nme="D_{s}^{*+}"; break;
2177  case 2212: c=kMagenta+2; s=20; nme="p"; break;
2178  case 2112: c=kMagenta+2; s=21; nme="n"; break;
2179  case 2114: c=kMagenta+2; s=22; nme="#Delta^{0}"; break;
2180  case 1114: c=kMagenta+2; s=23; nme="#Delta^{-}"; break;
2181  case 2214: c=kMagenta+2; s=24; nme="#Delta^{+}"; break;
2182  case 2224: c=kMagenta+2; s=25; nme="#Delta^{++}"; break;
2183  case 3112: c=kCyan +2; s=20; nme="#Sigma^{-}"; break;
2184  case 3222: c=kCyan +2; s=21; nme="#Sigma^{+}"; break;
2185  case 3114: c=kCyan +2; s=22; nme="#Sigma^{*-}"; break;
2186  case 3224: c=kCyan +2; s=23; nme="#Sigma^{*+}"; break;
2187  case 3312: c=kCyan +2; s=24; nme="#Xi^{-}"; break;
2188  case 3314: c=kCyan +2; s=25; nme="#Xi^{*-}"; break;
2189  case 3122: c=kCyan +2; s=26; nme="#Lambda"; break;
2190  case 3212: c=kCyan +2; s=27; nme="#Sigma^{0}"; break;
2191  case 3214: c=kCyan +2; s=28; nme="#Sigma^{*0}"; break;
2192  case 3322: c=kCyan +2; s=29; nme="#Xi^{0}"; break;
2193  case 3324: c=kCyan +2; s=30; nme="#Xi^{*0}"; break;
2194  case 4214: c=kYellow +2; s=20; nme="#Sigma^{*+}_{c}"; break;
2195  case 4224: c=kYellow +2; s=21; nme="#Sigma^{*++}_{c}";break;
2196  case 4122: c=kYellow +2; s=22; nme="#Lambda^{+}_{c}"; break;
2197  case 4114: c=kYellow +2; s=23; nme="#Sigma^{*0}_{c}"; break;
2198  case 4132: c=kYellow +2; s=24; nme="#Xi^{0}_{c}"; break;
2199  case 4314: c=kYellow +2; s=25; nme="#Xi^{*0}_{c}"; break;
2200  case 1000000000:c=kPink +2; s=20; nme="Nuclei"; break;
2201  default: c=kGray; s=1; nme.Form("%d",pdg); break;
2202  };
2203 }
2204 
2205 //____________________________________________________________________
2207 {
2208  Int_t size;
2209  Int_t* array = PdgArray(size);
2210  Int_t apdg = TMath::Abs(pdg);
2211  Int_t idx = TMath::BinarySearch(size, array, apdg);
2212  if (idx == size - 1) return idx+1;
2213  if (array[idx] != apdg) return size+1;
2214  return idx+1;
2215 }
2216 
2217 //____________________________________________________________________
2219 {
2220  static TAxis* axis = 0;
2221  if (axis) return *axis;
2222 
2223  Int_t size;
2224  Int_t* array = PdgArray(size);
2225  axis = new TAxis(size+1, +.5, size+1.5);
2226  // TDatabasePDG* pdgDb = TDatabasePDG::Instance();
2227  for (Int_t i = 1; i <= size; i++) {
2228  Int_t pdg = array[i-1];
2229  axis->SetBinLabel(i, Form("%d", pdg));
2230  // TString name = "?";
2231  // TParticlePDG* type = pdgDb->GetParticle(pdg);
2232  // if (type) name = type->GetName();
2233  // axis->SetBinLabel(i, name.Data());
2234  }
2235  axis->SetBinLabel(size+1, "-1");
2236  FixAxis(*axis);
2237  return *axis;
2238 }
2239 //____________________________________________________________________
2241 {
2242  TFile* file = TFile::Open(filename, "READ");
2243  if (!file) {
2244  ::Warning("OpenFile", "Failed to open \"%s\"", filename);
2245  return 0;
2246  }
2247  return file;
2248 }
2249 //____________________________________________________________________
2251 {
2252  static TString tmp;
2253  tmp = "all";
2254  if (c1 < c2)
2255  tmp.Form("cent%03dd%02d_%03dd%02d",
2256  Int_t(c1), Int_t(c1*100)%100,
2257  Int_t(c2), Int_t(c2*100)%100);
2258  return tmp.Data();
2259 }
2260 //____________________________________________________________________
2262 {
2263  const Color_t cc[] = { kMagenta+2, // 0
2264  kBlue+2, // 1
2265  kAzure-1, // 2 // 10,
2266  kCyan+2, // 3
2267  kGreen+1, // 4
2268  kSpring+5, // 5 //+10,
2269  kYellow+1, // 6
2270  kOrange+5, // 7 //+10,
2271  kRed+1, // 8
2272  kPink+5, // 9 //+10,
2273  kBlack }; // 10
2274  Color_t tc = cc[bin % 10];
2275  return tc;
2276 }
2277 //____________________________________________________________________
2279  Double_t c1,
2280  Double_t c2)
2281 {
2282  Double_t avgC = (c1+c2)/2;
2283  Int_t binC = const_cast<TAxis&>(axis).FindBin(avgC);
2284  return CentColor(binC);
2285 }
2286 
2287 //====================================================================
2289  const char* format, ...)
2290  : fMsg("")
2291 {
2292  if (lvl < msgLvl) return;
2293  va_list ap;
2294  va_start(ap, format);
2295  Format(fMsg, format, ap);
2296  va_end(ap);
2297  Output(+1, fMsg);
2298 }
2299 //____________________________________________________________________
2301 {
2302  if (fMsg.IsNull()) return;
2303  Output(-1, fMsg);
2304 }
2305 //____________________________________________________________________
2306 void
2308  const char* format, ...)
2309 {
2310  if (lvl < msgLvl) return;
2311  TString msg;
2312  va_list ap;
2313  va_start(ap, format);
2314  Format(msg, format, ap);
2315  va_end(ap);
2316  Output(0, msg);
2317 }
2318 
2319 //____________________________________________________________________
2320 void
2322  va_list ap)
2323 {
2324  static char buf[512];
2325  Int_t n = gROOT->GetDirLevel() + 2;
2326  for (Int_t i = 0; i < n; i++) buf[i] = ' ';
2327  vsnprintf(&(buf[n]), 511-n, format, ap);
2328  buf[511] = '\0';
2329  out = buf;
2330 }
2331 //____________________________________________________________________
2332 void
2334 {
2335  msg[0] = (in > 0 ? '>' : in < 0 ? '<' : '=');
2336  if (in < 0) gROOT->DecreaseDirLevel();
2337  // gROOT->IndentLevel();
2338  Printf(msg);
2339  if (in > 0) gROOT->IncreaseDirLevel();
2340 }
2341 
2342 
2343 #endif
2344 // Local Variables:
2345 // mode: C++
2346 // End:
2347 
Int_t color[]
static void Message(Int_t lvl, Int_t msgLvl, const char *format,...)
static TProfile * GetP1(Container *parent, const char *name, Bool_t verb=true)
static Double_t RatioE2(Double_t n, Double_t e2n, Double_t d, Double_t e2d, Double_t &e2r)
Int_t pdg
const char * filename
Definition: TestFCM.C:1
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)
static TDirectory * GetT(TDirectory *parent, const char *name, Bool_t verb=true)
const Color_t cc[]
Definition: DrawKs.C:1
double Double_t
Definition: External.C:58
static void ScaleAxis(TAxis &ret, Double_t fact=1)
const char * title
Definition: MakeQAPdf.C:26
Definition: External.C:244
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)
static TFile * OpenFile(const char *filename)
static TProfile2D * Make2P(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis, const TAxis &yAxis)
static Bool_t CheckAxisBins(const char *which, const TAxis *a1, const TAxis *a2)
static TH1 * GetH1(Container *parent, const char *name, Bool_t verb=true)
static const char * CentName(Double_t c1, Double_t c2)
static Bool_t CheckAxisLabels(const char *which, const TAxis *a1, const TAxis *a2)
static void Format(TString &out, const char *format, va_list ap)
static TH3 * ScaleDelta(TH3 *h, TH2 *scale)
TCanvas * c
Definition: TestFitELoss.C:172
static Bool_t CheckAxis(const char *which, const TAxis *a1, const TAxis *a2, Bool_t alsoLbls)
static Int_t GetI(Container *parent, const char *name, Int_t def=-1, Bool_t verb=true)
static TH2 * ScaleToIPz(TH2 *h, TH1 *ipZ, Bool_t full=false)
static TH2 * ProjectEtaDelta(TH3 *h)
static TH1 * ProjectDelta(TH2 *h)
static Double_t GetD(Container *parent, const char *name, Double_t def=-1, Bool_t verb=true)
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 Bool_t CheckConsistency(const TH1 *h1, const TH1 *h2)
static TH1 * Scale(TH1 *h, Double_t x, Double_t xe)
static Int_t GetB(Container *parent, const char *name, Bool_t def=false, Bool_t verb=true)
static void FixMinMax(TH1 *h, Bool_t ignoreZero=true)
int Int_t
Definition: External.C:63
static TProfile2D * GetP2(Container *parent, const char *name, Bool_t verb=true)
static TH1 * ProjectDeltaFull(TH3 *h)
static Int_t * PdgArray(Int_t &size)
static Double_t Integrate(TH1 *h, Double_t min, Double_t max, Double_t &err)
Definition: External.C:252
static Bool_t CheckAxisNBins(const char *which, const TAxis *a1, const TAxis *a2)
DebugGuard(Int_t lvl, Int_t msgLvl, const char *format,...)
Definition: External.C:228
Definition: External.C:212
ClassDef(AliTrackletAODUtils, 1)
static Bool_t CheckAxisLimits(const char *which, const TAxis *a1, const TAxis *a2)
static Color_t CentColor(const TAxis &axis, Double_t c1, Double_t c2)
static TH2 * CopyH2(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
Int_t mode
Definition: anaM.C:41
static TH1 * CopyH1(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
static TH1 * AverageOverIPz(TH2 *h, const char *name, UShort_t mode, TH1 *ipz, Double_t ipEff=0, TH2 *mask=0, Bool_t verb=false)
static TH3 * GetH3(Container *parent, const char *name, Bool_t verb=true)
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 TObject * GetO(Container *parent, const char *name, TClass *cls=0, Bool_t verb=true)
static THStack * GetHS(Container *parent, const char *name, Bool_t verb=true)
Definition: External.C:220
TFile * file
static Container * GetC(Container *parent, const char *name, Bool_t verb=true)
static void Output(int in, TString &msg)
static Int_t PdgBin(Int_t pdg)
unsigned short UShort_t
Definition: External.C:28
static TH3 * CopyH3(Container *parent, const char *name, const char *newName=0, Bool_t verb=true)
bool Bool_t
Definition: External.C:53
static TObject * CloneAndAdd(Container *c, TObject *o)
static TH1 * RatioH(const TH1 *num, const TH1 *denom, const char *name=0)
static TH1 * Make1D(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis)
static TProfile * Make1P(Container *c, const TString &name, const TString &title, Color_t color, Style_t style, const TAxis &xAxis)
static void FixAxis(TAxis &axis, const char *title=0)
Definition: External.C:196
static TProfile * GetP(Container *parent, const char *name, Bool_t verb=true)
static void CopyAttr(const TH1 *src, TH1 *dest)