#include "TCanvas.h"
#include "TF1.h"
#include "AliTHn.h"
+#include "THn.h"
ClassImp(AliUEHist)
-const Int_t AliUEHist::fgkCFSteps = 10;
+const Int_t AliUEHist::fgkCFSteps = 11;
-AliUEHist::AliUEHist(const char* reqHist) :
+AliUEHist::AliUEHist(const char* reqHist, const char* binning) :
TObject(),
fkRegions(4),
fEventHist(0),
fTrackHistEfficiency(0),
+ fFakePt(0),
fEtaMin(0),
fEtaMax(0),
fPtMin(0),
fPtMax(0),
+ fPartSpecies(-1),
fCentralityMin(0),
fCentralityMax(0),
fZVtxMin(0),
fZVtxMax(0),
+ fPt2Min(0),
fContaminationEnhancement(0),
fCombineMinMax(0),
+ fTrackEtaCut(0),
+ fWeightPerEvent(0),
+ fSkipScaleMixedEvent(kFALSE),
fCache(0),
+ fGetMultCacheOn(kFALSE),
+ fGetMultCache(0),
fHistogramType(reqHist)
{
// Constructor
if (strlen(reqHist) == 0)
return;
+
+ Printf("Creating AliUEHist with %s (binning: %s)", reqHist, binning);
AliLog::SetClassDebugLevel("AliCFContainer", -1);
AliLog::SetClassDebugLevel("AliCFGridSparse", -3);
// track level
Int_t nTrackVars = 4; // eta vs pT vs pT,lead (vs delta phi) vs multiplicity
- Int_t iTrackBin[6];
- Double_t* trackBins[6];
- const char* trackAxisTitle[6];
+ Int_t iTrackBin[7];
+ Double_t* trackBins[7];
+ const char* trackAxisTitle[7];
// eta
- iTrackBin[0] = 20;
- Double_t etaBins[20+1];
- for (Int_t i=0; i<=iTrackBin[0]; i++)
- etaBins[i] = -1.0 + 0.1 * i;
+ Int_t nEtaBins = -1;
+ Double_t* etaBins = GetBinning(binning, "eta", nEtaBins);
+ const char* etaTitle = "#eta";
+
+ iTrackBin[0] = nEtaBins;
trackBins[0] = etaBins;
- trackAxisTitle[0] = "#eta";
+ trackAxisTitle[0] = etaTitle;
// delta eta
- const Int_t kNDeltaEtaBins = 40+4;
- Double_t deltaEtaBins[kNDeltaEtaBins+1] = { -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1,
- -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1,
- -0.05, -0.025, 0, 0.025, 0.05,
- 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0,
- 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0 };
-
+ Int_t nDeltaEtaBins = -1;
+ Double_t* deltaEtaBins = GetBinning(binning, "delta_eta", nDeltaEtaBins);
+
// pT
- iTrackBin[1] = 22;
- Double_t pTBins[] = {0.15, 0.2, 0.3, 0.4, 0.5, 0.75, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0, 12.0, 15.0, 20.0};
- trackBins[1] = pTBins;
+ trackBins[1] = GetBinning(binning, "p_t_assoc", iTrackBin[1]);
trackAxisTitle[1] = "p_{T} (GeV/c)";
+ // pT, fine
+ Int_t npTBinsFine = -1;
+ Double_t* pTBinsFine = GetBinning(binning, "p_t_eff", npTBinsFine);
+
// pT,lead binning 1
- const Int_t kNLeadingpTBins = 100;
- Double_t leadingpTBins[kNLeadingpTBins+1];
- for (Int_t i=0; i<=kNLeadingpTBins; i++)
- leadingpTBins[i] = 0.5 * i;
-
+ Int_t nLeadingpTBins = -1;
+ Double_t* leadingpTBins = GetBinning(binning, "p_t_leading", nLeadingpTBins);
+
// pT,lead binning 2
- const Int_t kNLeadingpTBins2 = 9;
- Double_t leadingpTBins2[] = { 0.5, 1.0, 2.0, 3.0, 4.0, 6.0, 8.0, 10.0, 15.0, 20.0 };
-
- // phi,lead; this binning starts at -pi/2 and is modulo 3
- const Int_t kNLeadingPhiSpacing = 72;
- const Int_t kNLeadingPhiBins = kNLeadingPhiSpacing+4;
- Double_t leadingPhiBins[kNLeadingPhiBins+4+1];
- for (Int_t i=0; i<=17; i++)
- leadingPhiBins[i] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * i;
- leadingPhiBins[18] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * 17 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() / 2;
- leadingPhiBins[19] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * 17 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() / 2 * 1.5;
- leadingPhiBins[20] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * 17 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() / 2 * 2; // = 0
- leadingPhiBins[21] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * 17 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() / 2 * 2.5;
- leadingPhiBins[22] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * 17 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() / 2 * 3;
- for (Int_t i=19; i<=kNLeadingPhiSpacing; i++)
- leadingPhiBins[i+4] = -TMath::Pi() / 2 + 1.0 / kNLeadingPhiSpacing * TMath::TwoPi() * i;
-
- // multiplicity
- const Int_t kNMultiplicityBins = 15;
- Double_t multiplicityBins[kNMultiplicityBins+1];
- for (Int_t i=0; i<=kNMultiplicityBins; i++)
- multiplicityBins[i] = -0.5 + i;
- multiplicityBins[kNMultiplicityBins] = 200;
-
- trackBins[3] = multiplicityBins;
- iTrackBin[3] = kNMultiplicityBins;
- trackAxisTitle[3] = "multiplicity";
+ Int_t nLeadingpTBins2 = -1;
+ Double_t* leadingpTBins2 = GetBinning(binning, "p_t_leading_course", nLeadingpTBins2);
+
+ // phi,lead
+ Int_t nLeadingPhiBins = -1;
+ Double_t* leadingPhiBins = GetBinning(binning, "delta_phi", nLeadingPhiBins);
- // centrality (in %)
- const Int_t kNCentralityBins = 5 + 1 + 9;
- Double_t centralityBins[] = { 0, 1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100.1 };
+ trackBins[3] = GetBinning(binning, "multiplicity", iTrackBin[3]);
+ trackAxisTitle[3] = "multiplicity";
// particle species
const Int_t kNSpeciesBins = 4; // pi, K, p, rest
Double_t speciesBins[] = { -0.5, 0.5, 1.5, 2.5, 3.5 };
// vtx-z axis
- const Int_t kNVertexBins = 7;
- Double_t vertexBins[] = { -7, -5, -3, -1, 1, 3, 5, 7 };
+ const char* vertexTitle = "z-vtx (cm)";
+ Int_t nVertexBins = -1;
+ Double_t* vertexBins = GetBinning(binning, "vertex", nVertexBins);
+ Int_t nVertexBinsEff = -1;
+ Double_t* vertexBinsEff = GetBinning(binning, "vertex_eff", nVertexBinsEff);
- Bool_t useVtxAxis = kFALSE;
+ Int_t useVtxAxis = 0;
// selection depending on requested histogram
Int_t axis = -1; // 0 = pT,lead, 1 = phi,lead
axis = 1;
title = "d^{2}N_{ch}/d#varphid#eta";
}
- else if (strcmp(reqHist, "NumberDensityPhiCentrality") == 0 || strcmp(reqHist, "NumberDensityPhiCentralityVtx") == 0)
+ else if (TString(reqHist).BeginsWith("NumberDensityPhiCentrality"))
{
- if (strcmp(reqHist, "NumberDensityPhiCentralityVtx") == 0)
- {
- reqHist = "NumberDensityPhiCentrality";
- fHistogramType = reqHist;
- useVtxAxis = kTRUE;
- }
+ if (TString(reqHist).Contains("Vtx"))
+ useVtxAxis = 1;
+
+ reqHist = "NumberDensityPhiCentrality";
+ fHistogramType = reqHist;
axis = 2;
title = "d^{2}N_{ch}/d#varphid#eta";
}
axis = 0;
title = "d^{2}#Sigma p_{T}/d#varphid#eta";
}
+ else if (TString(reqHist).BeginsWith("TwoPlusOne"))
+ {
+ useVtxAxis = 1;
+
+ reqHist = "TwoPlusOne";
+ fHistogramType = reqHist;
+ axis = 3;
+ title = "d^{2}N_{ch}/d#varphid#eta";
+ }
else
AliFatal(Form("Invalid histogram requested: %s", reqHist));
if (axis == 0)
{
trackBins[2] = leadingpTBins;
- iTrackBin[2] = kNLeadingpTBins;
+ iTrackBin[2] = nLeadingpTBins;
trackAxisTitle[2] = "leading p_{T} (GeV/c)";
}
nTrackVars = 5;
initRegions = 1;
- iTrackBin[2] = kNLeadingpTBins2;
+ iTrackBin[2] = nLeadingpTBins2;
trackBins[2] = leadingpTBins2;
trackAxisTitle[2] = "leading p_{T} (GeV/c)";
- iTrackBin[4] = kNLeadingPhiBins;
+ iTrackBin[4] = nLeadingPhiBins;
trackBins[4] = leadingPhiBins;
trackAxisTitle[4] = "#Delta #varphi w.r.t. leading track";
}
nTrackVars = 5;
initRegions = 1;
- iTrackBin[0] = kNDeltaEtaBins;
+ iTrackBin[0] = nDeltaEtaBins;
trackBins[0] = deltaEtaBins;
trackAxisTitle[0] = "#Delta#eta";
- iTrackBin[2] = kNLeadingpTBins2;
+ iTrackBin[2] = nLeadingpTBins2;
trackBins[2] = leadingpTBins2;
trackAxisTitle[2] = "leading p_{T} (GeV/c)";
- trackBins[3] = centralityBins;
- iTrackBin[3] = kNCentralityBins;
trackAxisTitle[3] = "centrality";
- iTrackBin[4] = kNLeadingPhiBins;
+ iTrackBin[4] = nLeadingPhiBins;
trackBins[4] = leadingPhiBins;
- trackAxisTitle[4] = "#Delta#varphi (rad.)";
+ trackAxisTitle[4] = "#Delta#varphi (rad)";
- if (useVtxAxis)
+ if (useVtxAxis > 0)
{
nTrackVars = 6;
- iTrackBin[5] = kNVertexBins;
+ iTrackBin[5] = nVertexBins;
trackBins[5] = vertexBins;
- trackAxisTitle[5] = "z-vtx (cm)";
+ trackAxisTitle[5] = vertexTitle;
}
}
+ else if (axis == 3)
+ {
+ nTrackVars = 7;
+ initRegions = 1;
+
+ iTrackBin[0] = nDeltaEtaBins;
+ trackBins[0] = deltaEtaBins;
+ trackAxisTitle[0] = "#Delta#eta";
+
+ iTrackBin[2] = nLeadingpTBins;
+ trackBins[2] = leadingpTBins;
+ trackAxisTitle[2] = "Trigger 1 p_{T} (GeV/c)";
+
+ trackAxisTitle[3] = "centrality";
+
+ iTrackBin[4] = nLeadingPhiBins;
+ trackBins[4] = leadingPhiBins;
+ trackAxisTitle[4] = "#Delta#varphi (rad)";
+
+ iTrackBin[5] = nVertexBins;
+ trackBins[5] = vertexBins;
+ trackAxisTitle[5] = vertexTitle;
+
+ iTrackBin[6] = nLeadingpTBins2;
+ trackBins[6] = leadingpTBins2;
+ trackAxisTitle[6] = "Trigger 2 p_{T} (GeV/c)";
+ }
for (UInt_t i=0; i<initRegions; i++)
{
- if (strcmp(reqHist, "NumberDensityPhiCentrality") == 0)
+ if (axis >= 2)
fTrackHist[i] = new AliTHn(Form("fTrackHist_%d", i), title, fgkCFSteps, nTrackVars, iTrackBin);
else
fTrackHist[i] = new AliCFContainer(Form("fTrackHist_%d", i), title, fgkCFSteps, nTrackVars, iTrackBin);
// event level
Int_t nEventVars = 2;
- Int_t iEventBin[3];
+ Int_t iEventBin[4];
// track 3rd and 4th axis --> event 1st and 2nd axis
iEventBin[0] = iTrackBin[2];
iEventBin[1] = iTrackBin[3];
// plus track 5th axis (in certain cases)
- if (axis == 2 && useVtxAxis)
+ if (axis >= 2 && useVtxAxis)
{
nEventVars = 3;
iEventBin[2] = iTrackBin[5];
}
-
+ if (axis >= 3)
+ {
+ nEventVars = 4;
+ iEventBin[3] = iTrackBin[6];
+ }
fEventHist = new AliCFContainer("fEventHist", title, fgkCFSteps, nEventVars, iEventBin);
fEventHist->SetBinLimits(0, trackBins[2]);
fEventHist->SetBinLimits(1, trackBins[3]);
fEventHist->SetVarTitle(1, trackAxisTitle[3]);
- if (axis == 2 && useVtxAxis)
+ if (axis >= 2 && useVtxAxis)
{
fEventHist->SetBinLimits(2, trackBins[5]);
fEventHist->SetVarTitle(2, trackAxisTitle[5]);
}
+ if (axis >= 3)
+ {
+ fEventHist->SetBinLimits(3, trackBins[6]);
+ fEventHist->SetVarTitle(3, trackAxisTitle[6]);
+ }
SetStepNames(fEventHist);
+ iTrackBin[0] = nEtaBins;
+ iTrackBin[1] = npTBinsFine;
iTrackBin[2] = kNSpeciesBins;
+ iTrackBin[4] = nVertexBinsEff;
- fTrackHistEfficiency = new AliCFContainer("fTrackHistEfficiency", "Tracking efficiency", 3, 4, iTrackBin);
- fTrackHistEfficiency->SetBinLimits(0, trackBins[0]);
- fTrackHistEfficiency->SetVarTitle(0, trackAxisTitle[0]);
- fTrackHistEfficiency->SetBinLimits(1, trackBins[1]);
+ fTrackHistEfficiency = new AliCFContainer("fTrackHistEfficiency", "Tracking efficiency", 6, 5, iTrackBin);
+ fTrackHistEfficiency->SetBinLimits(0, etaBins);
+ fTrackHistEfficiency->SetVarTitle(0, etaTitle);
+ fTrackHistEfficiency->SetBinLimits(1, pTBinsFine);
fTrackHistEfficiency->SetVarTitle(1, trackAxisTitle[1]);
fTrackHistEfficiency->SetBinLimits(2, speciesBins);
fTrackHistEfficiency->SetVarTitle(2, "particle species");
fTrackHistEfficiency->SetBinLimits(3, trackBins[3]);
fTrackHistEfficiency->SetVarTitle(3, trackAxisTitle[3]);
+ fTrackHistEfficiency->SetBinLimits(4, vertexBinsEff);
+ fTrackHistEfficiency->SetVarTitle(4, vertexTitle);
+
+ fFakePt = new TH3F("fFakePt","fFakePt;p_{T,rec};p_{T};centrality", 200, 0, 20, 200, 0, 20, 20, 0, 100);
+
+ delete[] deltaEtaBins;
+ delete[] pTBinsFine;
+ delete[] leadingpTBins;
+ delete[] leadingpTBins2;
+ delete[] leadingPhiBins;
+ delete[] vertexBins;
+ delete[] vertexBinsEff;
+}
+
+TString AliUEHist::CombineBinning(TString defaultBinning, TString customBinning)
+{
+ // combine default binning with custom binning
+ // replaces binnings in default binning if it is defined in custom binning
+
+ TString binningStr;
+
+ TObjArray* lines = defaultBinning.Tokenize("\n");
+ for (Int_t i=0; i<lines->GetEntriesFast(); i++)
+ {
+ TString line(lines->At(i)->GetName());
+ TString tag = line(0, line.Index(":")+1);
+ if (!customBinning.BeginsWith(tag) && !customBinning.Contains(TString("\n") + tag))
+ binningStr += line + "\n";
+ else
+ Printf("Using custom binning for %s", tag.Data());
+ }
+ delete lines;
+ binningStr += customBinning;
+
+ return binningStr;
+}
+
+Double_t* AliUEHist::GetBinning(const char* configuration, const char* tag, Int_t& nBins)
+{
+ // takes the binning from <configuration> identified by <tag>
+ // configuration syntax example:
+ // eta: 2.4, -2.3, -2.2, -2.1, -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.1, 2.2, 2.3, 2.4
+ // phi: .....
+ //
+ // returns bin edges which have to be deleted by the caller
+
+ TString config(configuration);
+ TObjArray* lines = config.Tokenize("\n");
+ for (Int_t i=0; i<lines->GetEntriesFast(); i++)
+ {
+ TString line(lines->At(i)->GetName());
+ if (line.BeginsWith(TString(tag) + ":"))
+ {
+ line.Remove(0, strlen(tag) + 1);
+ line.ReplaceAll(" ", "");
+ TObjArray* binning = line.Tokenize(",");
+ Double_t* bins = new Double_t[binning->GetEntriesFast()];
+ for (Int_t j=0; j<binning->GetEntriesFast(); j++)
+ bins[j] = TString(binning->At(j)->GetName()).Atof();
+
+ nBins = binning->GetEntriesFast() - 1;
+
+ delete binning;
+ delete lines;
+ return bins;
+ }
+ }
+
+ delete lines;
+ AliFatal(Form("Tag %s not found in %s", tag, configuration));
+ return 0;
}
//_____________________________________________________________________________
fkRegions(4),
fEventHist(0),
fTrackHistEfficiency(0),
+ fFakePt(0),
fEtaMin(0),
fEtaMax(0),
fPtMin(0),
fPtMax(0),
+ fPartSpecies(-1),
fCentralityMin(0),
fCentralityMax(0),
fZVtxMin(0),
fZVtxMax(0),
+ fPt2Min(0),
fContaminationEnhancement(0),
fCombineMinMax(0),
+ fTrackEtaCut(0),
+ fWeightPerEvent(0),
+ fSkipScaleMixedEvent(kFALSE),
fCache(0),
+ fGetMultCacheOn(kFALSE),
+ fGetMultCache(0),
fHistogramType()
{
//
fTrackHistEfficiency = 0;
}
+ if (fFakePt)
+ {
+ delete fFakePt;
+ fFakePt = 0;
+ }
+
if (fCache)
{
delete fCache;
if (fTrackHistEfficiency)
target.fTrackHistEfficiency = dynamic_cast<AliCFContainer*> (fTrackHistEfficiency->Clone());
+ if (fFakePt)
+ target.fFakePt = dynamic_cast<TH3F*> (fFakePt->Clone());
+
target.fEtaMin = fEtaMin;
target.fEtaMax = fEtaMax;
target.fPtMin = fPtMin;
target.fPtMax = fPtMax;
+ target.fPartSpecies = fPartSpecies;
target.fCentralityMin = fCentralityMin;
target.fCentralityMax = fCentralityMax;
target.fZVtxMin = fZVtxMin;
target.fZVtxMax = fZVtxMax;
+ target.fPt2Min = fPt2Min;
if (fContaminationEnhancement)
target.fContaminationEnhancement = dynamic_cast<TH1F*> (fContaminationEnhancement->Clone());
target.fCombineMinMax = fCombineMinMax;
+ target.fTrackEtaCut = fTrackEtaCut;
+ target.fWeightPerEvent = fWeightPerEvent;
+ target.fSkipScaleMixedEvent = fSkipScaleMixedEvent;
target.fHistogramType = fHistogramType;
}
TObject* obj;
// collections of objects
- const UInt_t kMaxLists = fkRegions+2;
+ const UInt_t kMaxLists = fkRegions+3;
TList** lists = new TList*[kMaxLists];
for (UInt_t i=0; i<kMaxLists; i++)
lists[fkRegions]->Add(entry->fEventHist);
lists[fkRegions+1]->Add(entry->fTrackHistEfficiency);
+ if (entry->fFakePt)
+ lists[fkRegions+2]->Add(entry->fFakePt);
count++;
}
fEventHist->Merge(lists[fkRegions]);
fTrackHistEfficiency->Merge(lists[fkRegions+1]);
+ if (fFakePt)
+ fFakePt->Merge(lists[fkRegions+2]);
for (UInt_t i=0; i<kMaxLists; i++)
delete lists[i];
}
//____________________________________________________________________
-void AliUEHist::SetBinLimits(AliCFGridSparse* grid)
+void AliUEHist::SetBinLimits(THnBase* grid)
{
// sets the bin limits in eta and pT defined by fEtaMin/Max, fPtMin/Max
if (fEtaMax > fEtaMin)
- grid->SetRangeUser(0, fEtaMin, fEtaMax);
+ grid->GetAxis(0)->SetRangeUser(fEtaMin, fEtaMax);
if (fPtMax > fPtMin)
- grid->SetRangeUser(1, fPtMin, fPtMax);
+ grid->GetAxis(1)->SetRangeUser(fPtMin, fPtMax);
+ if (fPt2Min > 0)
+ grid->GetAxis(6)->SetRangeUser(fPt2Min, grid->GetAxis(6)->GetXmax() - 0.01);
}
+//____________________________________________________________________
+void AliUEHist::SetBinLimits(AliCFGridSparse* grid)
+{
+ // sets the bin limits in eta and pT defined by fEtaMin/Max, fPtMin/Max
+
+ SetBinLimits(grid->GetGrid());
+}
+
+//____________________________________________________________________
+void AliUEHist::ResetBinLimits(THnBase* grid)
+{
+ // resets all bin limits
+
+ for (Int_t i=0; i<grid->GetNdimensions(); i++)
+ if (grid->GetAxis(i)->TestBit(TAxis::kAxisRange))
+ grid->GetAxis(i)->SetRangeUser(0, -1);
+}
+
//____________________________________________________________________
void AliUEHist::ResetBinLimits(AliCFGridSparse* grid)
{
// resets all bin limits
- for (Int_t i=0; i<grid->GetNVar(); i++)
- if (grid->GetGrid()->GetAxis(i)->TestBit(TAxis::kAxisRange))
- grid->SetRangeUser(i, 0, -1);
+ ResetBinLimits(grid->GetGrid());
}
-
+
//____________________________________________________________________
void AliUEHist::CountEmptyBins(AliUEHist::CFStep step, Float_t ptLeadMin, Float_t ptLeadMax)
{
}
//____________________________________________________________________
-TH1* AliUEHist::GetUEHist(AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd, Int_t twoD, Bool_t etaNorm, Int_t* normEvents)
+TH1* AliUEHist::GetUEHist(AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd, Int_t twoD, Bool_t etaNorm, Long64_t* normEvents)
{
// Extracts the UE histogram at the given step and in the given region by projection and dividing tracks by events
//
// NOTE fEventHist contains the number of events for the underlying event analysis and the number of trigger particles for the azimuthal correlation analysis. In the latter case the naming is therefore somewhat misleading!
TH1D* events = fEventHist->ShowProjection(0, step);
- Int_t nEvents = (Int_t) events->Integral(firstBin, lastBin);
- Printf("Calculated histogram --> %d events", nEvents);
+ Long64_t nEvents = (Long64_t) events->Integral(firstBin, lastBin);
+ Printf("Calculated histogram --> %lld events", nEvents);
if (normEvents)
*normEvents = nEvents;
}
//____________________________________________________________________
-TH2* AliUEHist::GetSumOfRatios2(AliUEHist* mixed, AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd)
+void AliUEHist::GetHistsZVtxMult(AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, THnBase** trackHist, TH2** eventHist)
+{
+ // Calculates a 4d histogram with deltaphi, deltaeta, zvtx, multiplicity on track level and
+ // a 2d histogram on event level (as fct of zvtx, multiplicity)
+ // Histograms has to be deleted by the caller of the function
+
+ if (!fTrackHist[kToward]->GetGrid(0)->GetGrid()->GetAxis(5))
+ AliFatal("Histogram without vertex axis provided");
+
+ THnBase* sparse = fTrackHist[region]->GetGrid(step)->GetGrid();
+ if (fGetMultCacheOn)
+ {
+ if (!fGetMultCache)
+ {
+ fGetMultCache = ChangeToThn(sparse);
+ // should work but causes SEGV in ProjectionND below
+ }
+ sparse = fGetMultCache;
+ }
+
+ // unzoom all axes
+ ResetBinLimits(sparse);
+ ResetBinLimits(fEventHist->GetGrid(step));
+
+ SetBinLimits(sparse);
+
+ Int_t firstBin = sparse->GetAxis(2)->FindBin(ptLeadMin);
+ Int_t lastBin = sparse->GetAxis(2)->FindBin(ptLeadMax);
+ Printf("Using leading pT range %d --> %d", firstBin, lastBin);
+ sparse->GetAxis(2)->SetRange(firstBin, lastBin);
+ fEventHist->GetGrid(step)->GetGrid()->GetAxis(0)->SetRange(firstBin, lastBin);
+
+ // cut on the second trigger particle if there is a minimum set
+ if (fPt2Min > 0)
+ {
+ Int_t firstBinPt2 = sparse->GetAxis(6)->FindBin(fPt2Min);
+ Int_t lastBinPt2 = sparse->GetAxis(6)->GetNbins();
+
+ fEventHist->GetGrid(step)->GetGrid()->GetAxis(3)->SetRange(firstBinPt2, lastBinPt2);
+ }
+
+ Int_t dimensions[] = { 4, 0, 5, 3 };
+ THnBase* tmpTrackHist = sparse->ProjectionND(4, dimensions, "E");
+ *eventHist = (TH2*) fEventHist->GetGrid(step)->Project(2, 1);
+
+ ResetBinLimits(sparse);
+ ResetBinLimits(fEventHist->GetGrid(step));
+
+ // convert to THn
+ *trackHist = ChangeToThn(tmpTrackHist);
+ delete tmpTrackHist;
+}
+
+//____________________________________________________________________
+TH2* AliUEHist::GetSumOfRatios2(AliUEHist* mixed, AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd, Bool_t normalizePerTrigger, Int_t stepForMixed)
{
- // Calls GetUEHist(...) for *each* vertex bin and performs a sum of ratios:
+ // Calls GetUEHist(...) for *each* vertex bin and multiplicity bin and performs a sum of ratios:
// 1_N [ (same/mixed)_1 + (same/mixed)_2 + (same/mixed)_3 + ... ]
- // where N is the total number of events/trigger particles and the subscript is the vertex bin
+ // where N is the total number of events/trigger particles and the subscript is the vertex/multiplicity bin
// where mixed is normalized such that the information about the number of pairs in same is kept
//
// returns a 2D histogram: deltaphi, deltaeta
//
// Parameters:
- // mixed: AliUEHist containing mixed event corresponding to this object
+ // mixed: AliUEHist containing mixed event corresponding to this object (the histograms are taken from step <stepForMixed> if defined otherwise from step <step>)
// <other parameters> : check documentation of AliUEHist::GetUEHist
+ // normalizePerTrigger: divide through number of triggers
+
+ // do not add this hists to the directory
+ Bool_t oldStatus = TH1::AddDirectoryStatus();
+ TH1::AddDirectory(kFALSE);
TH2* totalTracks = 0;
- TH3* trackSameAll = 0;
- TH3* trackMixedAll = 0;
- TH1* eventSameAll = 0;
- TH1* eventMixedAll = 0;
+ THnBase* trackSameAll = 0;
+ THnBase* trackMixedAll = 0;
+ THnBase* trackMixedAllStep6 = 0;
+ TH2* eventSameAll = 0;
+ TH2* eventMixedAll = 0;
+ TH2* eventMixedAllStep6 = 0;
- Int_t totalEvents = 0;
+ Long64_t totalEvents = 0;
+ Int_t nCorrelationFunctions = 0;
- GetHistsZVtx(step, region, ptLeadMin, ptLeadMax, multBinBegin, multBinEnd, &trackSameAll, &eventSameAll);
- mixed->GetHistsZVtx(step, region, ptLeadMin, ptLeadMax, multBinBegin, multBinEnd, &trackMixedAll, &eventMixedAll);
+ GetHistsZVtxMult(step, region, ptLeadMin, ptLeadMax, &trackSameAll, &eventSameAll);
+ mixed->GetHistsZVtxMult((stepForMixed == -1) ? step : (CFStep) stepForMixed, region, ptLeadMin, ptLeadMax, &trackMixedAll, &eventMixedAll);
- TAxis* vertexAxis = trackSameAll->GetZaxis();
- for (Int_t vertexBin = 1; vertexBin <= vertexAxis->GetNbins(); vertexBin++)
+ // If we ask for histograms from step8 (TTR cut applied) there is a hole at 0,0; so this cannot be used for the
+ // mixed-event normalization. If step6 is available, the normalization factor is read out from that one.
+ // If step6 is not available we fallback to taking the normalization along all delta phi (WARNING requires a
+ // flat delta phi distribution)
+ if (stepForMixed == -1 && step == kCFStepBiasStudy && mixed->fEventHist->GetGrid(kCFStepReconstructed)->GetEntries() > 0 && !fSkipScaleMixedEvent)
{
- trackSameAll->GetZaxis()->SetRange(vertexBin, vertexBin);
- trackMixedAll->GetZaxis()->SetRange(vertexBin, vertexBin);
+ Printf("Using mixed-event normalization factors from step %d", kCFStepReconstructed);
+ mixed->GetHistsZVtxMult(kCFStepReconstructed, region, ptLeadMin, ptLeadMax, &trackMixedAllStep6, &eventMixedAllStep6);
+ }
+
+// Printf("%f %f %f %f", trackSameAll->GetEntries(), eventSameAll->GetEntries(), trackMixedAll->GetEntries(), eventMixedAll->GetEntries());
+
+// TH1* normParameters = new TH1F("normParameters", "", 100, 0, 2);
+
+// trackSameAll->Dump();
- TH2* tracksSame = (TH2*) trackSameAll->Project3D("yx1");
- TH2* tracksMixed = (TH2*) trackMixedAll->Project3D("yx2");
-
- // asssume flat in dphi, gain in statistics
-// TH1* histMixedproj = mixedTwoD->ProjectionY();
-// histMixedproj->Scale(1.0 / mixedTwoD->GetNbinsX());
-//
-// for (Int_t x=1; x<=mixedTwoD->GetNbinsX(); x++)
-// for (Int_t y=1; y<=mixedTwoD->GetNbinsY(); y++)
-// mixedTwoD->SetBinContent(x, y, histMixedproj->GetBinContent(y));
+ TAxis* multAxis = trackSameAll->GetAxis(3);
+
+ if (multBinEnd < multBinBegin)
+ {
+ multBinBegin = 1;
+ multBinEnd = multAxis->GetNbins();
+ }
+
+ for (Int_t multBin = TMath::Max(1, multBinBegin); multBin <= TMath::Min(multAxis->GetNbins(), multBinEnd); multBin++)
+ {
+ trackSameAll->GetAxis(3)->SetRange(multBin, multBin);
+ trackMixedAll->GetAxis(3)->SetRange(multBin, multBin);
+ if (trackMixedAllStep6)
+ trackMixedAllStep6->GetAxis(3)->SetRange(multBin, multBin);
-// delete histMixedproj;
+ Double_t mixedNorm = 1;
+ Double_t mixedNormError = 0;
- // get mixed event normalization by assuming full acceptance at deta of 0
- Double_t mixedNorm = tracksMixed->Integral(tracksMixed->GetXaxis()->FindBin(-0.01), tracksMixed->GetXaxis()->FindBin(0.01), tracksMixed->GetYaxis()->FindBin(-0.01), tracksMixed->GetYaxis()->FindBin(0.01));
- mixedNorm /= (tracksMixed->GetXaxis()->FindBin(0.01) - tracksMixed->GetXaxis()->FindBin(-0.01) + 1) * (tracksMixed->GetYaxis()->FindBin(0.01) - tracksMixed->GetYaxis()->FindBin(-0.01) + 1);
- if (mixedNorm <= 0)
+ if (!fSkipScaleMixedEvent)
{
- Printf("ERROR: Skipping vertex bin %d because mixed event is empty at (0,0)", vertexBin);
+ // get mixed normalization correction factor: is independent of vertex bin if scaled with number of triggers
+ TH2* tracksMixed = 0;
+ if (trackMixedAllStep6)
+ {
+ trackMixedAllStep6->GetAxis(2)->SetRange(0, -1);
+ tracksMixed = trackMixedAllStep6->Projection(1, 0, "E");
+ }
+ else
+ {
+ trackMixedAll->GetAxis(2)->SetRange(0, -1);
+ tracksMixed = trackMixedAll->Projection(1, 0, "E");
+ }
+ // Printf("%f", tracksMixed->Integral());
+ Float_t binWidthEta = tracksMixed->GetYaxis()->GetBinWidth(1);
+
+ if (step == kCFStepBiasStudy && !trackMixedAllStep6)
+ {
+ // get mixed event normalization by assuming full acceptance at deta at 0 (integrate over dphi), excluding (0, 0)
+ Float_t phiExclude = 0.41;
+ mixedNorm = tracksMixed->IntegralAndError(1, tracksMixed->GetXaxis()->FindBin(-phiExclude)-1, tracksMixed->GetYaxis()->FindBin(-0.01), tracksMixed->GetYaxis()->FindBin(0.01), mixedNormError);
+ Double_t mixedNormError2 = 0;
+ Double_t mixedNorm2 = tracksMixed->IntegralAndError(tracksMixed->GetXaxis()->FindBin(phiExclude)+1, tracksMixed->GetNbinsX(), tracksMixed->GetYaxis()->FindBin(-0.01), tracksMixed->GetYaxis()->FindBin(0.01), mixedNormError2);
+
+ if (mixedNormError == 0 || mixedNormError2 == 0)
+ {
+ Printf("ERROR: Skipping multiplicity %d because mixed event is empty %f %f %f %f", multBin, mixedNorm, mixedNormError, mixedNorm2, mixedNormError2);
+ continue;
+ }
+
+ Int_t nBinsMixedNorm = (tracksMixed->GetXaxis()->FindBin(-phiExclude) - 1 - 1 + 1) * (tracksMixed->GetYaxis()->FindBin(0.01) - tracksMixed->GetYaxis()->FindBin(-0.01) + 1);
+ mixedNorm /= nBinsMixedNorm;
+ mixedNormError /= nBinsMixedNorm;
+
+ Int_t nBinsMixedNorm2 = (tracksMixed->GetNbinsX() - tracksMixed->GetXaxis()->FindBin(phiExclude) - 1 + 1) * (tracksMixed->GetYaxis()->FindBin(0.01) - tracksMixed->GetYaxis()->FindBin(-0.01) + 1);
+ mixedNorm2 /= nBinsMixedNorm2;
+ mixedNormError2 /= nBinsMixedNorm2;
+
+ mixedNorm = mixedNorm / mixedNormError / mixedNormError + mixedNorm2 / mixedNormError2 / mixedNormError2;
+ mixedNormError = TMath::Sqrt(1.0 / (1.0 / mixedNormError / mixedNormError + 1.0 / mixedNormError2 / mixedNormError2));
+ mixedNorm *= mixedNormError * mixedNormError;
+ }
+ else
+ {
+ // get mixed event normalization at (0,0)
+
+ mixedNorm = tracksMixed->IntegralAndError(tracksMixed->GetXaxis()->FindBin(-0.01), tracksMixed->GetXaxis()->FindBin(0.01), tracksMixed->GetYaxis()->FindBin(-0.01), tracksMixed->GetYaxis()->FindBin(0.01), mixedNormError);
+ Int_t nBinsMixedNorm = (tracksMixed->GetXaxis()->FindBin(0.01) - tracksMixed->GetXaxis()->FindBin(-0.01) + 1) * (tracksMixed->GetYaxis()->FindBin(0.01) - tracksMixed->GetYaxis()->FindBin(-0.01) + 1);
+ mixedNorm /= nBinsMixedNorm;
+ mixedNormError /= nBinsMixedNorm;
+
+ if (mixedNormError == 0)
+ {
+ Printf("ERROR: Skipping multiplicity %d because mixed event is empty %f %f", multBin, mixedNorm, mixedNormError);
+ continue;
+ }
+ }
+
+ // finite bin correction
+ if (fTrackEtaCut > 0)
+ {
+ Double_t finiteBinCorrection = -1.0 / (2*fTrackEtaCut) * binWidthEta / 2 + 1;
+ Printf("Finite bin correction: %f", finiteBinCorrection);
+ mixedNorm /= finiteBinCorrection;
+ mixedNormError /= finiteBinCorrection;
+ }
+ else
+ {
+ Printf("ERROR: fTrackEtaCut not set. Finite bin correction cannot be applied. Continuing anyway...");
+ }
+
+ Float_t triggers = eventMixedAll->Integral(1, eventMixedAll->GetNbinsX(), multBin, multBin);
+// Printf("%f +- %f | %f | %f", mixedNorm, mixedNormError, triggers, mixedNorm / triggers);
+ if (triggers <= 0)
+ {
+ Printf("ERROR: Skipping multiplicity %d because mixed event is empty", multBin);
+ continue;
+ }
+
+ mixedNorm /= triggers;
+ mixedNormError /= triggers;
+
+ delete tracksMixed;
}
else
+ Printf("WARNING: Skipping mixed-event scaling! fSkipScaleMixedEvent IS set!");
+
+ if (mixedNorm <= 0)
{
- tracksMixed->Scale(1.0 / mixedNorm);
+ Printf("ERROR: Skipping multiplicity %d because mixed event is empty at (0,0)", multBin);
+ continue;
+ }
+
+// Printf("Norm: %f +- %f", mixedNorm, mixedNormError);
- // tracksSame->Scale(tracksMixed->Integral() / tracksSame->Integral());
+// normParameters->Fill(mixedNorm);
- tracksSame->Divide(tracksMixed);
+ TAxis* vertexAxis = trackSameAll->GetAxis(2);
+ Int_t vertexBinBegin = 1;
+ Int_t vertexBinEnd = vertexAxis->GetNbins();
+
+ if (fZVtxMax > fZVtxMin)
+ {
+ vertexBinBegin = vertexAxis->FindBin(fZVtxMin);
+ vertexBinEnd = vertexAxis->FindBin(fZVtxMax);
+ }
+
+ for (Int_t vertexBin = vertexBinBegin; vertexBin <= vertexBinEnd; vertexBin++)
+ {
+ trackSameAll->GetAxis(2)->SetRange(vertexBin, vertexBin);
+ trackMixedAll->GetAxis(2)->SetRange(vertexBin, vertexBin);
+
+ TH2* tracksSame = trackSameAll->Projection(1, 0, "E");
+ TH2* tracksMixed = trackMixedAll->Projection(1, 0, "E");
- // code to draw contributions
- /*
- TH1* proj = tracksSame->ProjectionX("projx", tracksSame->GetYaxis()->FindBin(-1.59), tracksSame->GetYaxis()->FindBin(1.59));
- proj->SetTitle(Form("Bin %d", vertexBin));
- proj->SetLineColor(vertexBin);
- proj->DrawCopy((vertexBin > 1) ? "SAME" : "");
- */
+ // asssume flat in dphi, gain in statistics
+ // TH1* histMixedproj = mixedTwoD->ProjectionY();
+ // histMixedproj->Scale(1.0 / mixedTwoD->GetNbinsX());
+ //
+ // for (Int_t x=1; x<=mixedTwoD->GetNbinsX(); x++)
+ // for (Int_t y=1; y<=mixedTwoD->GetNbinsY(); y++)
+ // mixedTwoD->SetBinContent(x, y, histMixedproj->GetBinContent(y));
+
+ // delete histMixedproj;
- if (!totalTracks)
- totalTracks = (TH2*) tracksSame->Clone("totalTracks");
+ Float_t triggers2 = eventMixedAll->Integral(vertexBin, vertexBin, multBin, multBin);
+ if (triggers2 <= 0)
+ {
+ Printf("ERROR: Skipping multiplicity %d vertex bin %d because mixed event is empty", multBin, vertexBin);
+ }
else
- totalTracks->Add(tracksSame);
+ {
+ if (!fSkipScaleMixedEvent)
+ tracksMixed->Scale(1.0 / triggers2 / mixedNorm);
+ else if (tracksMixed->Integral() > 0)
+ tracksMixed->Scale(1.0 / tracksMixed->Integral());
+ // tracksSame->Scale(tracksMixed->Integral() / tracksSame->Integral());
+
+// new TCanvas; tracksSame->DrawClone("SURF1");
+// new TCanvas; tracksMixed->DrawClone("SURF1");
+
+ // some code to judge the relative contribution of the different correlation functions to the overall uncertainty
+ Double_t sums[] = { 0, 0, 0 };
+ Double_t errors[] = { 0, 0, 0 };
+
+ for (Int_t x=1; x<=tracksSame->GetNbinsX(); x++)
+ for (Int_t y=1; y<=tracksSame->GetNbinsY(); y++)
+ {
+ sums[0] += tracksSame->GetBinContent(x, y);
+ errors[0] += tracksSame->GetBinError(x, y);
+ sums[1] += tracksMixed->GetBinContent(x, y);
+ errors[1] += tracksMixed->GetBinError(x, y);
+ }
+
+ tracksSame->Divide(tracksMixed);
+
+ for (Int_t x=1; x<=tracksSame->GetNbinsX(); x++)
+ for (Int_t y=1; y<=tracksSame->GetNbinsY(); y++)
+ {
+ sums[2] += tracksSame->GetBinContent(x, y);
+ errors[2] += tracksSame->GetBinError(x, y);
+ }
+
+ for (Int_t x=0; x<3; x++)
+ if (sums[x] > 0)
+ errors[x] /= sums[x];
+
+ Printf("The correlation function %d %d has uncertainties %f %f %f (Ratio S/M %f)", multBin, vertexBin, errors[0], errors[1], errors[2], (errors[1] > 0) ? errors[0] / errors[1] : -1);
+ // code to draw contributions
+ /*
+ TH1* proj = tracksSame->ProjectionX("projx", tracksSame->GetYaxis()->FindBin(-1.59), tracksSame->GetYaxis()->FindBin(1.59));
+ proj->SetTitle(Form("Bin %d", vertexBin));
+ proj->SetLineColor(vertexBin);
+ proj->DrawCopy((vertexBin > 1) ? "SAME" : "");
+ */
+
+ if (!totalTracks)
+ totalTracks = (TH2*) tracksSame->Clone("totalTracks");
+ else
+ totalTracks->Add(tracksSame);
+
+ totalEvents += eventSameAll->GetBinContent(vertexBin, multBin);
+
+// new TCanvas; tracksMixed->DrawCopy("SURF1");
+ }
delete tracksSame;
delete tracksMixed;
- totalEvents += eventSameAll->GetBinContent(vertexBin);
+ nCorrelationFunctions++;
}
}
if (totalTracks) {
- Printf("Dividing %f tracks by %d events", totalTracks->Integral(), totalEvents);
- if (totalEvents > 0)
- totalTracks->Scale(1.0 / totalEvents);
+ Double_t sums[] = { 0, 0, 0 };
+ Double_t errors[] = { 0, 0, 0 };
+
+ for (Int_t x=1; x<=totalTracks->GetNbinsX(); x++)
+ for (Int_t y=1; y<=totalTracks->GetNbinsY(); y++)
+ {
+ sums[0] += totalTracks->GetBinContent(x, y);
+ errors[0] += totalTracks->GetBinError(x, y);
+ }
+ if (sums[0] > 0)
+ errors[0] /= sums[0];
+
+ if (normalizePerTrigger)
+ {
+ Printf("Dividing %f tracks by %lld events (%d correlation function(s)) (error %f)", totalTracks->Integral(), totalEvents, nCorrelationFunctions, errors[0]);
+ if (totalEvents > 0)
+ totalTracks->Scale(1.0 / totalEvents);
+ }
// normalizate to dphi width
Float_t normalization = totalTracks->GetXaxis()->GetBinWidth(1);
delete trackSameAll;
delete trackMixedAll;
+ delete trackMixedAllStep6;
delete eventSameAll;
delete eventMixedAll;
+ delete eventMixedAllStep6;
+
+// new TCanvas; normParameters->Draw();
+
+ TH1::AddDirectory(oldStatus);
return totalTracks;
}
+TH1* AliUEHist::GetTriggersAsFunctionOfMultiplicity(AliUEHist::CFStep step, Float_t ptLeadMin, Float_t ptLeadMax)
+{
+ // returns the distribution of triggers as function of centrality/multiplicity
+
+ ResetBinLimits(fEventHist->GetGrid(step));
+
+ Int_t firstBin = fEventHist->GetGrid(step)->GetGrid()->GetAxis(0)->FindBin(ptLeadMin);
+ Int_t lastBin = fEventHist->GetGrid(step)->GetGrid()->GetAxis(0)->FindBin(ptLeadMax);
+ Printf("Using pT range %d --> %d", firstBin, lastBin);
+ fEventHist->GetGrid(step)->GetGrid()->GetAxis(0)->SetRange(firstBin, lastBin);
+
+ if (fZVtxMax > fZVtxMin)
+ {
+ fEventHist->GetGrid(step)->GetGrid()->GetAxis(2)->SetRangeUser(fZVtxMin, fZVtxMax);
+ Printf("Restricting z-vtx: %f-%f", fZVtxMin, fZVtxMax);
+ }
+
+ TH1* eventHist = fEventHist->GetGrid(step)->Project(1);
+
+ ResetBinLimits(fEventHist->GetGrid(step));
+
+ return eventHist;
+}
+
+/*
+TH2* AliUEHist::GetPtDistInPhiRegion(AliUEHist* mixed, AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd, Float_t phiBegin, Float_t phiEnd)
+{
+ // Returns pT,assoc distribution in the given pt,trig, multiplicity, phi region
+ // Does not use sum of ratios for mixed event correction (TODO to be improved)
+ // returns a 2D histogram: deltaphi, deltaeta
+ //
+ // Parameters:
+ // mixed: AliUEHist containing mixed event corresponding to this object
+*/
+
//____________________________________________________________________
TH2* AliUEHist::GetSumOfRatios(AliUEHist* mixed, AliUEHist::CFStep step, AliUEHist::Region region, Float_t ptLeadMin, Float_t ptLeadMax, Int_t multBinBegin, Int_t multBinEnd, Bool_t etaNorm, Bool_t useVertexBins)
{
multIter++;
}
- Int_t nEvents = 0;
+ Long64_t nEvents = 0;
TH2* tracks = (TH2*) GetUEHist(step, region, ptLeadMin, ptLeadMax, multBinBeginLocal, multBinEndLocal, 1, etaNorm, &nEvents);
// undo normalization
tracks->Scale(nEvents);
// if <histogram> is 0, just copies content from step1 to step2
// clear target histogram
- target->Reset();
+ if (grid != target)
+ target->Reset();
if (histogram != 0)
{
THnSparse* grid = fTrackHist[region]->GetGrid(step1)->GetGrid();
THnSparse* target = fTrackHist[region]->GetGrid(step2)->GetGrid();
+ Float_t entriesBefore = grid->GetEntries();
+
MultiplyHistograms(grid, target, trackCorrection, var1, var2);
- Printf("AliUEHist::CorrectTracks: Corrected from %f to %f entries. Correction histogram: %f entries (integral: %f)", grid->GetEntries(), target->GetEntries(), (trackCorrection) ? trackCorrection->GetEntries() : -1.0, (trackCorrection) ? trackCorrection->Integral() : -1.0);
+ Printf("AliUEHist::CorrectTracks: Corrected from %f (step %d) to %f (step %d) entries. Correction histogram: %f entries (integral: %f)", entriesBefore, step1, target->GetEntries(), step2, (trackCorrection) ? trackCorrection->GetEntries() : -1.0, (trackCorrection) ? trackCorrection->Integral() : -1.0);
}
//____________________________________________________________________
AliCFGridSparse* grid = fEventHist->GetGrid(step1);
AliCFGridSparse* target = fEventHist->GetGrid(step2);
+ Float_t entriesBefore = grid->GetEntries();
+
MultiplyHistograms(grid->GetGrid(), target->GetGrid(), eventCorrection, var1, var2);
- Printf("AliUEHist::CorrectEvents: Corrected from %f to %f entries. Correction histogram: %f entries (integral: %f)", grid->GetEntries(), target->GetEntries(), (eventCorrection) ? eventCorrection->GetEntries() : -1.0, (eventCorrection) ? eventCorrection->Integral() : -1.0);
+ Printf("AliUEHist::CorrectEvents: Corrected from %f (step %d) to %f (step %d) entries. Correction histogram: %f entries (integral: %f)", entriesBefore, step1, target->GetEntries(), step2, (eventCorrection) ? eventCorrection->GetEntries() : -1.0, (eventCorrection) ? eventCorrection->Integral() : -1.0);
}
//____________________________________________________________________
}
else if (strcmp(fHistogramType, "NumberDensityPhiCentrality") == 0)
{
- // copy
- CorrectTracks(kCFStepReconstructed, kCFStepTracked, 0, -1);
- CorrectEvents(kCFStepReconstructed, kCFStepTracked, 0, -1);
-
- // Dont use eta in the following, because it is a Delta-eta axis
-
- // contamination correction
- // correct single-particle contamination for associated particles
-
- TH1* contamination = corrections->GetTrackingContamination(1);
-
- if (0)
+ if (fTrackHist[0]->GetNVar() <= 5)
{
- Printf("Applying contamination enhancement");
+ // do corrections copying between steps
+// CFStep step = kCFStepReconstructed;
+ CFStep step = kCFStepBiasStudy;
+
+ // copy
+ CorrectTracks(step, kCFStepTracked, 0, -1);
+ CorrectEvents(step, kCFStepTracked, 0, -1);
+
+ // Dont use eta in the following, because it is a Delta-eta axis
- for (Int_t bin = 1; bin <= contamination->GetNbinsX(); bin++)
+ // contamination correction
+ // correct single-particle contamination for associated particles
+
+ TH1* contamination = corrections->GetTrackingContamination(1);
+
+ if (0)
{
- printf("%f", contamination->GetBinContent(bin));
- if (contamination->GetBinContent(bin) > 0)
- contamination->SetBinContent(bin, 1.0 + 1.1 * (contamination->GetBinContent(bin) - 1.0));
- printf(" --> %f\n", contamination->GetBinContent(bin));
+ Printf("Applying contamination enhancement");
+
+ for (Int_t bin = 1; bin <= contamination->GetNbinsX(); bin++)
+ {
+ printf("%f", contamination->GetBinContent(bin));
+ if (contamination->GetBinContent(bin) > 0)
+ contamination->SetBinContent(bin, 1.0 + 1.1 * (contamination->GetBinContent(bin) - 1.0));
+ printf(" --> %f\n", contamination->GetBinContent(bin));
+ }
}
- }
+
+ CorrectTracks(kCFStepTracked, kCFStepTrackedOnlyPrim, contamination, 1);
+ delete contamination;
- CorrectTracks(kCFStepTracked, kCFStepTrackedOnlyPrim, contamination, 1);
- delete contamination;
-
- // correct for additional contamination due to trigger particle around phi ~ 0
- TH2* correlatedContamination = corrections->GetCorrelatedContamination();
- if (0)
+ // correct for additional contamination due to trigger particle around phi ~ 0
+ TH2* correlatedContamination = corrections->GetCorrelatedContamination();
+ if (0)
+ {
+ Printf("Applying contamination enhancement");
+
+ for (Int_t bin = 1; bin <= correlatedContamination->GetNbinsX(); bin++)
+ for (Int_t bin2 = 1; bin2 <= correlatedContamination->GetNbinsY(); bin2++)
+ {
+ printf("%f", correlatedContamination->GetBinContent(bin, bin2));
+ if (correlatedContamination->GetBinContent(bin, bin2) > 0)
+ correlatedContamination->SetBinContent(bin, bin2, 1.0 + 1.1 * (correlatedContamination->GetBinContent(bin, bin2) - 1.0));
+ printf(" --> %f\n", correlatedContamination->GetBinContent(bin, bin2));
+ }
+ }
+
+// new TCanvas; correlatedContamination->DrawCopy("COLZ");
+ CorrectCorrelatedContamination(kCFStepTrackedOnlyPrim, 0, correlatedContamination);
+// Printf("\n\n\nWARNING ---> SKIPPING CorrectCorrelatedContamination\n\n\n");
+
+ delete correlatedContamination;
+
+ // TODO correct for contamination of trigger particles (for tracks AND events)
+ CorrectEvents(kCFStepTracked, kCFStepTrackedOnlyPrim, 0, 0);
+
+ // --- efficiency correction ---
+ // correct single-particle efficiency for associated particles
+ // in addition correct for efficiency on trigger particles (tracks AND events)
+
+ // in bins of pT and centrality
+ TH1* efficiencyCorrection = corrections->GetTrackingEfficiencyCorrectionCentrality();
+// new TCanvas; efficiencyCorrection->DrawCopy("COLZ");
+ // use kCFStepAnaTopology as a temporary step
+ CorrectTracks(kCFStepTrackedOnlyPrim, kCFStepAnaTopology, efficiencyCorrection, 1, 3);
+ delete efficiencyCorrection;
+
+ // correct pT,T in bins of pT and centrality
+ efficiencyCorrection = corrections->GetTrackEfficiency(kCFStepTrackedOnlyPrim, kCFStepAnaTopology, 1, 3, 2);
+ CorrectEvents(kCFStepTrackedOnlyPrim, kCFStepVertex, efficiencyCorrection, 0, 1);
+ CorrectTracks(kCFStepAnaTopology, kCFStepVertex, efficiencyCorrection, 2, 3);
+ delete efficiencyCorrection;
+
+ // no correction for vertex finding efficiency and trigger efficiency needed in PbPb
+ // copy
+ CorrectTracks(kCFStepVertex, kCFStepAll, 0, -1);
+ CorrectEvents(kCFStepVertex, kCFStepAll, 0, -1);
+ }
+ else
{
- Printf("Applying contamination enhancement");
+ // with 6 axes there is not enough memory, do the corrections in-place
+ Printf(">>>>>>>> Applying corrections in place to reduce memory consumption");
+ CFStep step = kCFStepBiasStudy;
+// CFStep step = kCFStepReconstructed;
+
+ // Dont use eta in the following, because it is a Delta-eta axis
- for (Int_t bin = 1; bin <= correlatedContamination->GetNbinsX(); bin++)
- for (Int_t bin2 = 1; bin2 <= correlatedContamination->GetNbinsY(); bin2++)
- {
- printf("%f", correlatedContamination->GetBinContent(bin, bin2));
- if (correlatedContamination->GetBinContent(bin, bin2) > 0)
- correlatedContamination->SetBinContent(bin, bin2, 1.0 + 1.1 * (correlatedContamination->GetBinContent(bin, bin2) - 1.0));
- printf(" --> %f\n", correlatedContamination->GetBinContent(bin, bin2));
- }
+ // --- contamination correction ---
+ // correct single-particle contamination for associated particles
+ // correct contamination for trigger particles (tracks AND events)
+
+ // in bins of p,T
+ TH1* contamination = corrections->GetTrackingContamination(1);
+
+ if (0)
+ {
+ Printf("Applying contamination enhancement");
+
+ for (Int_t bin = 1; bin <= contamination->GetNbinsX(); bin++)
+ {
+ printf("%f", contamination->GetBinContent(bin));
+ if (contamination->GetBinContent(bin) > 0)
+ contamination->SetBinContent(bin, 1.0 + 1.1 * (contamination->GetBinContent(bin) - 1.0));
+ printf(" --> %f\n", contamination->GetBinContent(bin));
+ }
+ }
+
+ // correct pT,A in bins of pT
+ CorrectTracks(step, step, contamination, 1);
+ delete contamination;
+
+ // correct pT,T in bins of pT (for tracks AND events)
+ contamination = corrections->GetTrackEfficiency(kCFStepTracked, kCFStepTrackedOnlyPrim, 1, -1, 2);
+ new TCanvas; contamination->DrawCopy();
+ CorrectEvents(step, step, contamination, 0);
+ CorrectTracks(step, step, contamination, 2);
+ delete contamination;
+
+ // correct for additional contamination due to trigger particle around phi ~ 0
+ if (0)
+ {
+ TH2* correlatedContamination = corrections->GetCorrelatedContamination();
+ if (1)
+ {
+ Printf("Applying contamination enhancement");
+
+ for (Int_t bin = 1; bin <= correlatedContamination->GetNbinsX(); bin++)
+ for (Int_t bin2 = 1; bin2 <= correlatedContamination->GetNbinsY(); bin2++)
+ {
+ printf("%f", correlatedContamination->GetBinContent(bin, bin2));
+ if (correlatedContamination->GetBinContent(bin, bin2) > 0)
+ correlatedContamination->SetBinContent(bin, bin2, 1.0 + 1.1 * (correlatedContamination->GetBinContent(bin, bin2) - 1.0));
+ printf(" --> %f\n", correlatedContamination->GetBinContent(bin, bin2));
+ }
+ }
+
+ // new TCanvas; correlatedContamination->DrawCopy("COLZ");
+ CorrectCorrelatedContamination(step, 0, correlatedContamination);
+
+ delete correlatedContamination;
+ }
+ else
+ Printf("\n\n\nWARNING ---> SKIPPING CorrectCorrelatedContamination\n\n\n");
+
+ // --- tracking efficiency correction ---
+ // correct single-particle efficiency for associated particles
+ // correct for efficiency on trigger particles (tracks AND events)
+
+ // in bins of pT and centrality
+ TH1* efficiencyCorrection = corrections->GetTrackingEfficiencyCorrectionCentrality();
+ // new TCanvas; efficiencyCorrection->DrawCopy("COLZ");
+
+ // correct pT,A in bins of pT and centrality
+ CorrectTracks(step, step, efficiencyCorrection, 1, 3);
+ delete efficiencyCorrection;
+
+ // correct pT,T in bins of pT and centrality
+ efficiencyCorrection = corrections->GetTrackEfficiency(kCFStepTrackedOnlyPrim, kCFStepAnaTopology, 1, 3, 2);
+ CorrectEvents(step, step, efficiencyCorrection, 0, 1);
+ CorrectTracks(step, step, efficiencyCorrection, 2, 3);
+
+ delete efficiencyCorrection;
}
-
- new TCanvas; correlatedContamination->DrawCopy("COLZ");
-// CorrectCorrelatedContamination(kCFStepTrackedOnlyPrim, 0, correlatedContamination);
- Printf("\n\n\nWARNING ---> SKIPPING CorrectCorrelatedContamination\n\n\n");
-
- delete correlatedContamination;
-
- // TODO correct for contamination of trigger particles (for tracks AND events)
- CorrectEvents(kCFStepTracked, kCFStepTrackedOnlyPrim, 0, 0);
-
- // --- efficiency correction ---
- // correct single-particle efficiency for associated particles
- // in addition correct for efficiency on trigger particles (tracks AND events)
-
- // in bins of pT and centrality
- TH1* efficiencyCorrection = corrections->GetTrackingEfficiencyCorrectionCentrality();
- new TCanvas; efficiencyCorrection->DrawCopy("COLZ");
- // use kCFStepAnaTopology as a temporary step
- CorrectTracks(kCFStepTrackedOnlyPrim, kCFStepAnaTopology, efficiencyCorrection, 1, 3);
- delete efficiencyCorrection;
-
- // correct pT,T in bins of pT and centrality
- efficiencyCorrection = corrections->GetTrackEfficiency(kCFStepTrackedOnlyPrim, kCFStepAnaTopology, 1, 3, 2);
- CorrectEvents(kCFStepTrackedOnlyPrim, kCFStepVertex, efficiencyCorrection, 0, 1);
- CorrectTracks(kCFStepAnaTopology, kCFStepVertex, efficiencyCorrection, 2, 3);
- delete efficiencyCorrection;
-
- // no correction for vertex finding efficiency and trigger efficiency needed in PbPb
- // copy
- CorrectTracks(kCFStepVertex, kCFStepAll, 0, -1);
- CorrectEvents(kCFStepVertex, kCFStepAll, 0, -1);
}
else
AliFatal(Form("Unknown histogram for correction: %s", GetTitle()));
}
+//____________________________________________________________________
+THnBase* AliUEHist::GetTrackEfficiencyND(CFStep step1, CFStep step2)
+{
+ // creates a track-level efficiency by dividing step2 by step1
+ // in all dimensions but the particle species one
+
+ AliCFContainer* sourceContainer = fTrackHistEfficiency;
+ // step offset because we start with kCFStepAnaTopology
+ step1 = (CFStep) ((Int_t) step1 - (Int_t) kCFStepAnaTopology);
+ step2 = (CFStep) ((Int_t) step2 - (Int_t) kCFStepAnaTopology);
+
+ ResetBinLimits(sourceContainer->GetGrid(step1));
+ ResetBinLimits(sourceContainer->GetGrid(step2));
+
+ if (fEtaMax > fEtaMin)
+ {
+ Printf("Restricted eta-range to %f %f", fEtaMin, fEtaMax);
+ sourceContainer->GetGrid(step1)->SetRangeUser(0, fEtaMin, fEtaMax);
+ sourceContainer->GetGrid(step2)->SetRangeUser(0, fEtaMin, fEtaMax);
+ }
+
+ Int_t dimensions[] = { 0, 1, 3, 4 };
+ THnBase* generated = sourceContainer->GetGrid(step1)->GetGrid()->ProjectionND(4, dimensions);
+ THnBase* measured = sourceContainer->GetGrid(step2)->GetGrid()->ProjectionND(4, dimensions);
+
+// Printf("%d %d %f %f", step1, step2, generated->GetEntries(), measured->GetEntries());
+
+ ResetBinLimits(sourceContainer->GetGrid(step1));
+ ResetBinLimits(sourceContainer->GetGrid(step2));
+
+ measured->Divide(measured, generated, 1, 1, "B");
+
+ delete generated;
+
+ return measured;
+}
+
//____________________________________________________________________
TH1* AliUEHist::GetTrackEfficiency(CFStep step1, CFStep step2, Int_t axis1, Int_t axis2, Int_t source, Int_t axis3)
{
}
if (fPtMax > fPtMin && axis1 != 1 && axis2 != 1 && axis3 != 1)
{
+ Printf("Restricted pt-range to %f %f", fPtMin, fPtMax);
sourceContainer->GetGrid(step1)->SetRangeUser(1, fPtMin, fPtMax);
sourceContainer->GetGrid(step2)->SetRangeUser(1, fPtMin, fPtMax);
}
+ if (fPartSpecies != -1 && axis1 != 2 && axis2 != 2 && axis3 != 2)
+ {
+ Printf("Restricted to particle species %d", fPartSpecies);
+ sourceContainer->GetGrid(step1)->SetRangeUser(2, fPartSpecies, fPartSpecies);
+ sourceContainer->GetGrid(step2)->SetRangeUser(2, fPartSpecies, fPartSpecies);
+ }
if (fCentralityMax > fCentralityMin && axis1 != 3 && axis2 != 3 && axis3 != 3)
{
+ Printf("Restricted centrality range to %f %f", fCentralityMin, fCentralityMax);
sourceContainer->GetGrid(step1)->SetRangeUser(3, fCentralityMin, fCentralityMax);
sourceContainer->GetGrid(step2)->SetRangeUser(3, fCentralityMin, fCentralityMax);
}
+ if (fZVtxMax > fZVtxMin && axis1 != 4 && axis2 != 4 && axis3 != 4)
+ {
+ Printf("Restricted z-vtx range to %f %f", fZVtxMin, fZVtxMax);
+ sourceContainer->GetGrid(step1)->SetRangeUser(4, fZVtxMin, fZVtxMax);
+ sourceContainer->GetGrid(step2)->SetRangeUser(4, fZVtxMin, fZVtxMax);
+ }
TH1* measured = 0;
TH1* generated = 0;
return dynamic_cast<TH2D*> (GetTrackEfficiency(kCFStepAnaTopology, kCFStepTrackedOnlyPrim, 0, 1));
}
+//____________________________________________________________________
+TH2D* AliUEHist::GetFakeRate()
+{
+ return dynamic_cast<TH2D*> (GetTrackEfficiency(kCFStepTracked, (CFStep) (kCFStepTracked+3), 0, 1));
+}
+
//____________________________________________________________________
TH2D* AliUEHist::GetTrackingEfficiencyCentrality()
{
return dynamic_cast<TH1D*> (GetTrackEfficiency(kCFStepAnaTopology, kCFStepTrackedOnlyPrim, axis));
}
+//____________________________________________________________________
+TH1D* AliUEHist::GetFakeRate(Int_t axis)
+{
+ return dynamic_cast<TH1D*> (GetTrackEfficiency(kCFStepTracked, (CFStep) (kCFStepTracked+3), axis));
+}
//____________________________________________________________________
TH2D* AliUEHist::GetTrackingCorrection()
{
return "Bias study applying tracking efficiency";
case kCFStepBiasStudy2:
return "Bias study applying tracking efficiency in two steps";
+ case kCFStepCorrected:
+ return "Corrected for efficiency on-the-fly";
}
return 0;
fEventHist->SetGrid(AliUEHist::kCFStepBiasStudy, from->fEventHist->GetGrid(AliUEHist::kCFStepBiasStudy));
}
+void AliUEHist::DeepCopy(AliUEHist* from)
+{
+ // copies the entries of this object's members from the object <from> to this object
+ // fills using the fill function and thus allows that the objects have different binning
+
+ for (Int_t region=0; region<4; region++)
+ {
+ if (!fTrackHist[region] || !from->fTrackHist[region])
+ continue;
+
+ for (Int_t step=0; step<fTrackHist[region]->GetNStep(); step++)
+ {
+ Printf("Copying region %d step %d", region, step);
+ THnSparse* target = fTrackHist[region]->GetGrid(step)->GetGrid();
+ THnSparse* source = from->fTrackHist[region]->GetGrid(step)->GetGrid();
+
+ target->Reset();
+ target->RebinnedAdd(source);
+ }
+ }
+
+ for (Int_t step=0; step<fEventHist->GetNStep(); step++)
+ {
+ Printf("Ev: Copying step %d", step);
+ THnSparse* target = fEventHist->GetGrid(step)->GetGrid();
+ THnSparse* source = from->fEventHist->GetGrid(step)->GetGrid();
+
+ target->Reset();
+ target->RebinnedAdd(source);
+ }
+
+ for (Int_t step=0; step<TMath::Min(fTrackHistEfficiency->GetNStep(), from->fTrackHistEfficiency->GetNStep()); step++)
+ {
+ if (!from->fTrackHistEfficiency->GetGrid(step))
+ continue;
+
+ Printf("Eff: Copying step %d", step);
+ THnSparse* target = fTrackHistEfficiency->GetGrid(step)->GetGrid();
+ THnSparse* source = from->fTrackHistEfficiency->GetGrid(step)->GetGrid();
+
+ target->Reset();
+ target->RebinnedAdd(source);
+ }
+}
+
+void AliUEHist::SymmetrizepTBins()
+{
+ // copy pt,a < pt,t bins to pt,a > pt,t (inverting deltaphi and delta eta as it should be) including symmetric bins
+
+ for (Int_t region=0; region<4; region++)
+ {
+ if (!fTrackHist[region])
+ continue;
+
+ for (Int_t step=0; step<fTrackHist[region]->GetNStep(); step++)
+ {
+ Printf("Copying region %d step %d", region, step);
+ THnSparse* target = fTrackHist[region]->GetGrid(step)->GetGrid();
+ if (target->GetEntries() == 0)
+ continue;
+
+ // for symmetric bins
+ THnSparse* source = (THnSparse*) target->Clone();
+
+ Int_t zVtxBins = 1;
+ if (target->GetNdimensions() > 5)
+ zVtxBins = target->GetAxis(5)->GetNbins();
+
+ // axes: 0 delta eta; 1 pT,a; 2 pT,t; 3 centrality; 4 delta phi; 5 vtx-z
+ for (Int_t i3 = 1; i3 <= target->GetAxis(3)->GetNbins(); i3++)
+ for (Int_t i5 = 1; i5 <= zVtxBins; i5++)
+ {
+ for (Int_t i1 = 1; i1 <= target->GetAxis(1)->GetNbins(); i1++)
+ for (Int_t i2 = 1; i2 <= target->GetAxis(2)->GetNbins(); i2++)
+ {
+ // find source bin
+ Int_t binA = target->GetAxis(1)->FindBin(target->GetAxis(2)->GetBinCenter(i2));
+ Int_t binT = target->GetAxis(2)->FindBin(target->GetAxis(1)->GetBinCenter(i1));
+
+ Printf("(%d %d) Copying from %d %d to %d %d", i3, i5, binA, binT, i1, i2);
+
+ for (Int_t i0 = 1; i0 <= target->GetAxis(0)->GetNbins(); i0++)
+ for (Int_t i4 = 1; i4 <= target->GetAxis(4)->GetNbins(); i4++)
+ {
+ Int_t binEta = target->GetAxis(0)->FindBin(-target->GetAxis(0)->GetBinCenter(i0));
+ Double_t phi = -target->GetAxis(4)->GetBinCenter(i4);
+ if (phi < -TMath::Pi()/2)
+ phi += TMath::TwoPi();
+ Int_t binPhi = target->GetAxis(4)->FindBin(phi);
+
+ Int_t binSource[] = { binEta, binA, binT, i3, binPhi, i5 };
+ Int_t binTarget[] = { i0, i1, i2, i3, i4, i5 };
+
+ Double_t value = source->GetBinContent(binSource);
+ Double_t error = source->GetBinError(binSource);
+
+ if (error == 0)
+ continue;
+
+ Double_t value2 = target->GetBinContent(binTarget);
+ Double_t error2 = target->GetBinError(binTarget);
+
+ Double_t sum = value;
+ Double_t err = error;
+
+ if (error2 > 0)
+ {
+ sum = value + value2;
+ err = TMath::Sqrt(error * error + error2 * error2);
+ }
+
+ // Printf(" Values: %f +- %f; %f +- %f --> %f +- %f", value, error, value2, error2, sum, err);
+
+ target->SetBinContent(binTarget, sum);
+ target->SetBinError(binTarget, err);
+ }
+ }
+ }
+
+ delete source;
+ }
+ }
+}
+
//____________________________________________________________________
void AliUEHist::ExtendTrackingEfficiency(Bool_t verbose)
{
for (Int_t step=0; step<fTrackHistEfficiency->GetNStep(); step++)
fTrackHistEfficiency->GetGrid(step)->GetGrid()->Reset();
}
+
+THnBase* AliUEHist::ChangeToThn(THnBase* sparse)
+{
+ // change the object to THn for faster processing
+
+ // convert to THn (SEGV's for some strange reason...)
+ // x = THn::CreateHn("a", "a", sparse);
+
+ // own implementation
+ Int_t nBins[10];
+ for (Int_t i=0; i<sparse->GetNdimensions(); i++)
+ nBins[i] = sparse->GetAxis(i)->GetNbins();
+ THn* tmpTHn = new THnF(Form("%s_thn", sparse->GetName()), sparse->GetTitle(), sparse->GetNdimensions(), nBins, 0, 0);
+ for (Int_t i=0; i<sparse->GetNdimensions(); i++)
+ {
+ tmpTHn->SetBinEdges(i, sparse->GetAxis(i)->GetXbins()->GetArray());
+ tmpTHn->GetAxis(i)->SetTitle(sparse->GetAxis(i)->GetTitle());
+ }
+ tmpTHn->RebinnedAdd(sparse);
+
+ return tmpTHn;
+}
+
+void AliUEHist::CondenseBin(THnSparse* grid, THnSparse* target, Int_t axis, Float_t targetValue, Float_t from, Float_t to)
+{
+ //
+ // loops through the histogram and moves all entries to a single point <targetValue> on the axis <axis>
+ // if <from> and <to> are set, then moving only occurs if value on <axis> is betweem <from> and <to>
+ //
+
+ if (grid->GetNdimensions() > 6)
+ AliFatal("Too many dimensions in THnSparse");
+
+ Int_t targetBin = grid->GetAxis(axis)->FindBin(targetValue);
+ AliInfo(Form("Target bin on axis %d with value %f is %d", axis, targetValue, targetBin));
+
+ Int_t fromBin = 1;
+ Int_t toBin = grid->GetAxis(axis)->GetNbins();
+ if (to > from)
+ {
+ fromBin = grid->GetAxis(axis)->FindBin(from);
+ toBin = grid->GetAxis(axis)->FindBin(to);
+ AliInfo(Form("Only condensing between bin %d and %d", fromBin, toBin));
+ }
+
+ Int_t bins[6];
+ for (Int_t binIdx = 0; binIdx < grid->GetNbins(); binIdx++)
+ {
+ Double_t value = grid->GetBinContent(binIdx, bins);
+ Double_t error = grid->GetBinError(binIdx);
+
+ if (bins[axis] >= fromBin && bins[axis] <= toBin)
+ bins[axis] = targetBin;
+
+ value += target->GetBinContent(bins);
+ error = TMath::Sqrt(error * error + target->GetBinError(bins) * target->GetBinError(bins));
+
+ target->SetBinContent(bins, value);
+ target->SetBinError(bins, error);
+ }
+}
+
+void AliUEHist::CondenseBin(CFStep step, Int_t trackAxis, Int_t eventAxis, Float_t targetValue, Float_t from, Float_t to, CFStep tmpStep)
+{
+ // loops through the histogram at <step> and moves all entries to a single point <targetValue> on the axes
+ // <trackAxis> and <eventAxis>. <tmpStep> is used to temporary store the data
+ // This is useful e.g. to move bin content around for MC productions where the centrality selection did
+ // not yield the desired result
+
+ // reset tmpStep
+ fEventHist->GetGrid(tmpStep)->GetGrid()->Reset();
+ for (UInt_t i=0; i<fkRegions; i++)
+ if (fTrackHist[i])
+ fTrackHist[i]->GetGrid(tmpStep)->GetGrid()->Reset();
+
+ // copy to tmpStep
+ CorrectTracks(step, tmpStep, 0, -1);
+ CorrectEvents(step, tmpStep, 0, -1);
+
+ // reset step
+ fEventHist->GetGrid(step)->GetGrid()->Reset();
+ for (UInt_t i=0; i<fkRegions; i++)
+ if (fTrackHist[i])
+ fTrackHist[i]->GetGrid(step)->GetGrid()->Reset();
+
+ // rewriting
+ for (UInt_t i=0; i<fkRegions; i++)
+ {
+ if (!fTrackHist[i])
+ continue;
+
+ THnSparse* grid = fTrackHist[i]->GetGrid(tmpStep)->GetGrid();
+ THnSparse* target = fTrackHist[i]->GetGrid(step)->GetGrid();
+
+ CondenseBin(grid, target, trackAxis, targetValue, from, to);
+ }
+ CondenseBin(fEventHist->GetGrid(tmpStep)->GetGrid(), fEventHist->GetGrid(step)->GetGrid(), eventAxis, targetValue, from, to);
+
+ // reset tmpStep
+ fEventHist->GetGrid(tmpStep)->GetGrid()->Reset();
+ for (UInt_t i=0; i<fkRegions; i++)
+ if (fTrackHist[i])
+ fTrackHist[i]->GetGrid(tmpStep)->GetGrid()->Reset();
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff