THnSparse* target = GetGrid(i)->GetGrid();
Int_t* binIdx = new Int_t[fNVars];
+ Int_t* nBins = new Int_t[fNVars];
for (Int_t j=0; j<fNVars; j++)
+ {
binIdx[j] = 1;
+ nBins[j] = target->GetAxis(j)->GetNbins();
+ }
Long64_t count = 0;
for (Int_t j=fNVars-1; j>0; j--)
{
- if (binIdx[j] > target->GetAxis(j)->GetNbins())
+ if (binIdx[j] > nBins[j])
{
binIdx[j] = 1;
binIdx[j-1]++;
}
}
- if (binIdx[0] > target->GetAxis(0)->GetNbins())
+ if (binIdx[0] > nBins[0])
break;
}
AliInfo(Form("Step %d: copied %lld entries out of %lld bins", i, count, GetGlobalBinIndex(binIdx)));
+
+ delete[] binIdx;
+ delete[] nBins;
+ }
+}
+
+void AliTHn::ReduceAxis()
+{
+ // "removes" one axis by summing over the axis and putting the entry to bin 1
+ // TODO presently only implemented for the last axis
+
+ Int_t axis = fNVars-1;
+
+ for (Int_t i=0; i<fNSteps; i++)
+ {
+ if (!fValues[i])
+ continue;
+
+ Float_t* source = fValues[i]->GetArray();
+ Float_t* sourceSumw2 = 0;
+ if (fSumw2[i])
+ sourceSumw2 = fSumw2[i]->GetArray();
+
+ THnSparse* target = GetGrid(i)->GetGrid();
+
+ Int_t* binIdx = new Int_t[fNVars];
+ Int_t* nBins = new Int_t[fNVars];
+ for (Int_t j=0; j<fNVars; j++)
+ {
+ binIdx[j] = 1;
+ nBins[j] = target->GetAxis(j)->GetNbins();
+ }
+
+ Long64_t count = 0;
+
+ while (1)
+ {
+ // sum over axis <axis>
+ Float_t sumValues = 0;
+ Float_t sumSumw2 = 0;
+ for (Int_t j=1; j<=nBins[axis]; j++)
+ {
+ binIdx[axis] = j;
+ Long64_t globalBin = GetGlobalBinIndex(binIdx);
+ sumValues += source[globalBin];
+ source[globalBin] = 0;
+
+ if (sourceSumw2)
+ {
+ sumSumw2 += sourceSumw2[globalBin];
+ sourceSumw2[globalBin] = 0;
+ }
+ }
+ binIdx[axis] = 1;
+
+ Long64_t globalBin = GetGlobalBinIndex(binIdx);
+ source[globalBin] = sumValues;
+ if (sourceSumw2)
+ sourceSumw2[globalBin] = sumSumw2;
+
+ count++;
+
+ // next bin
+ binIdx[fNVars-2]++;
+
+ for (Int_t j=fNVars-2; j>0; j--)
+ {
+ if (binIdx[j] > nBins[j])
+ {
+ binIdx[j] = 1;
+ binIdx[j-1]++;
+ }
+ }
+
+ if (binIdx[0] > nBins[0])
+ break;
+ }
+
+ AliInfo(Form("Step %d: reduced %lld bins to %lld entries", i, GetGlobalBinIndex(binIdx), count));
+
+ delete[] binIdx;
+ delete[] nBins;
}
}
return tracks;
}
-//____________________________________________________________________
-void AliUEHist::CorrectTracks(CFStep step1, CFStep step2, TH1* trackCorrection, Int_t var1, Int_t var2)
+void AliUEHist::MultiplyHistograms(THnSparse* grid, THnSparse* target, TH1* histogram, Int_t var1, Int_t var2)
{
- // corrects from step1 to step2 by multiplying the tracks with trackCorrection
- // trackCorrection can be a function of eta (var1 == 0), pT (var1 == 1), leading pT (var1 == 2), multiplicity (var1 == 3), delta phi (var1 == 4)
- // if var2 >= 0 a two dimension correction is assumed in trackCorrection
+ // multiplies <grid> with <histogram> and put the result in <target>
+ // <grid> has usually more dimensions than histogram. The axis which are used to choose the value
+ // from <histogram> are given in <var1> and <var2>
//
- // if trackCorrection is 0, just copies content from step1 to step2
-
- for (UInt_t region=0; region<fkRegions; region++)
- CorrectTracks(step1, step2, region, trackCorrection, var1, var2);
-}
-
-//____________________________________________________________________
-void AliUEHist::CorrectTracks(CFStep step1, CFStep step2, Int_t region, TH1* trackCorrection, Int_t var1, Int_t var2)
-{
- //
- // see documentation of CorrectTracks above
- //
-
- if (!fTrackHist[region])
- return;
-
- THnSparse* grid = fTrackHist[region]->GetGrid(step1)->GetGrid();
- THnSparse* target = fTrackHist[region]->GetGrid(step2)->GetGrid();
+ // if <histogram> is 0, just copies content from step1 to step2
// clear target histogram
target->Reset();
- if (trackCorrection != 0)
+ if (histogram != 0)
{
- if (grid->GetAxis(var1)->GetNbins() != trackCorrection->GetNbinsX())
- AliFatal(Form("Invalid binning (var1): %d %d", grid->GetAxis(var1)->GetNbins(), trackCorrection->GetNbinsX()));
+ if (grid->GetAxis(var1)->GetNbins() != histogram->GetNbinsX())
+ AliFatal(Form("Invalid binning (var1): %d %d", grid->GetAxis(var1)->GetNbins(), histogram->GetNbinsX()));
- if (var2 >= 0 && grid->GetAxis(var2)->GetNbins() != trackCorrection->GetNbinsY())
- AliFatal(Form("Invalid binning (var2): %d %d", grid->GetAxis(var2)->GetNbins(), trackCorrection->GetNbinsY()));
+ if (var2 >= 0 && grid->GetAxis(var2)->GetNbins() != histogram->GetNbinsY())
+ AliFatal(Form("Invalid binning (var2): %d %d", grid->GetAxis(var2)->GetNbins(), histogram->GetNbinsY()));
}
+
+ if (grid->GetNdimensions() > 6)
+ AliFatal("Too many dimensions in THnSparse");
+
+ Int_t bins[6];
// optimized implementation
for (Int_t binIdx = 0; binIdx < grid->GetNbins(); binIdx++)
{
- Int_t bins[5];
Double_t value = grid->GetBinContent(binIdx, bins);
Double_t error = grid->GetBinError(binIdx);
- if (trackCorrection != 0)
+ if (histogram != 0)
{
if (var2 < 0)
{
- value *= trackCorrection->GetBinContent(bins[var1]);
- error *= trackCorrection->GetBinContent(bins[var1]);
+ value *= histogram->GetBinContent(bins[var1]);
+ error *= histogram->GetBinContent(bins[var1]);
}
else
{
- value *= trackCorrection->GetBinContent(bins[var1], bins[var2]);
- error *= trackCorrection->GetBinContent(bins[var1], bins[var2]);
+ value *= histogram->GetBinContent(bins[var1], bins[var2]);
+ error *= histogram->GetBinContent(bins[var1], bins[var2]);
}
}
target->SetBinContent(bins, value);
target->SetBinError(bins, error);
}
-
+}
+
+//____________________________________________________________________
+void AliUEHist::CorrectTracks(CFStep step1, CFStep step2, TH1* trackCorrection, Int_t var1, Int_t var2)
+{
+ // corrects from step1 to step2 by multiplying the tracks with trackCorrection
+ // trackCorrection can be a function of eta (var1 == 0), pT (var1 == 1), leading pT (var1 == 2), multiplicity (var1 == 3), delta phi (var1 == 4)
+ // if var2 >= 0 a two dimension correction is assumed in trackCorrection
+ //
+ // if trackCorrection is 0, just copies content from step1 to step2
+
+ for (UInt_t region=0; region<fkRegions; region++)
+ CorrectTracks(step1, step2, region, trackCorrection, var1, var2);
+}
+
+//____________________________________________________________________
+void AliUEHist::CorrectTracks(CFStep step1, CFStep step2, Int_t region, TH1* trackCorrection, Int_t var1, Int_t var2)
+{
+ //
+ // see documentation of CorrectTracks above
+ //
+
+ if (!fTrackHist[region])
+ return;
+
+ THnSparse* grid = fTrackHist[region]->GetGrid(step1)->GetGrid();
+ THnSparse* target = fTrackHist[region]->GetGrid(step2)->GetGrid();
+
+ 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);
}
AliCFGridSparse* grid = fEventHist->GetGrid(step1);
AliCFGridSparse* target = fEventHist->GetGrid(step2);
- // clear target histogram
- target->GetGrid()->Reset();
-
- if (eventCorrection != 0 && grid->GetNBins(var1) != eventCorrection->GetNbinsX())
- AliFatal(Form("Invalid binning: %d %d", grid->GetNBins(var1), eventCorrection->GetNbinsX()));
-
- if (eventCorrection != 0 && var2 != -1 && grid->GetNBins(var2) != eventCorrection->GetNbinsY())
- AliFatal(Form("Invalid binning: %d %d", grid->GetNBins(var2), eventCorrection->GetNbinsY()));
-
- Int_t bins[2];
- for (Int_t x = 1; x <= grid->GetNBins(0); x++)
- {
- bins[0] = x;
- for (Int_t y = 1; y <= grid->GetNBins(1); y++)
- {
- bins[1] = y;
-
- Double_t value = grid->GetElement(bins);
- if (value != 0)
- {
- Double_t error = grid->GetElementError(bins);
-
- if (eventCorrection != 0)
- {
- if (var2 == -1)
- {
- value *= eventCorrection->GetBinContent(bins[var1]);
- error *= eventCorrection->GetBinContent(bins[var1]);
- }
- else
- {
- value *= eventCorrection->GetBinContent(bins[var1], bins[var2]);
- error *= eventCorrection->GetBinContent(bins[var1], bins[var2]);
- }
- }
-
- target->SetElement(bins, value);
- target->SetElementError(bins, error);
- }
- }
- }
-
+ 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);
}
}
}
- //new TCanvas; correlatedContamination->DrawCopy("COLZ");
+// new TCanvas; correlatedContamination->DrawCopy("COLZ");
CorrectCorrelatedContamination(kCFStepTrackedOnlyPrim, 0, correlatedContamination);
+// Printf("\n\n\nWARNING ---> SKIPPING CorrectCorrelatedContamination\n\n\n");
delete correlatedContamination;
// optimized implementation
for (Int_t binIdx = 0; binIdx < grid->GetNbins(); binIdx++)
{
- Int_t bins[5];
+ Int_t bins[6];
Double_t value = grid->GetBinContent(binIdx, bins);
Double_t error = grid->GetBinError(binIdx);
// optimized implementation
for (Int_t binIdx = 0; binIdx < grid->GetNbins(); binIdx++)
{
- Int_t bins[5];
+ Int_t bins[6];
Double_t value = grid->GetBinContent(binIdx, bins);
Double_t error = grid->GetBinError(binIdx);