#include <TH1D.h>
#include <THStack.h>
#include <TList.h>
-#include <TParameter.h>
#include <AliAnalysisManager.h>
#include <AliAODEvent.h>
#include <AliAODHandler.h>
#include "AliAODForwardMult.h"
#include <TFile.h>
#include <TStyle.h>
+#include <TROOT.h>
//____________________________________________________________________
AliBasedNdetaTask::AliBasedNdetaTask()
- : AliAnalysisTaskSE(),
- fSums(0), // Container of sums
- fOutput(0), // Container of output
- fVtxMin(0), // Minimum v_z
- fVtxMax(0), // Maximum v_z
- fTriggerMask(0), // Trigger mask
+ : AliBaseAODTask(),
fRebin(0), // Rebinning factor
fCutEdges(false),
fSymmetrice(true),
fCorrEmpty(true),
fUseROOTProj(false),
fTriggerEff(1),
+ fTriggerEff0(1),
fShapeCorr(0),
fListOfCentralities(0),
- fSNNString(0),
- fSysString(0),
- fCent(0),
- fCentAxis(0),
fNormalizationScheme(kFull),
- fSchemeString(0),
- fTriggerString(0),
- fFinalMCCorrFile(""),
- fIsESD(false)
+ fFinalMCCorrFile(""),
+ fSatelliteVertices(0),
+ fglobalempiricalcorrection(0),
+ fmeabsignalvscentr(0)
{
//
// Constructor
//
- DGUARD(0,0,"Default construction of AliBasedNdetaTask");
+ DGUARD(fDebug,3,"Default CTOR of AliBasedNdetaTask");
}
//____________________________________________________________________
AliBasedNdetaTask::AliBasedNdetaTask(const char* name)
- : AliAnalysisTaskSE(name),
- fSums(0), // Container of sums
- fOutput(0), // Container of output
- fVtxMin(-10), // Minimum v_z
- fVtxMax(10), // Maximum v_z
- fTriggerMask(AliAODForwardMult::kInel),
+ : AliBaseAODTask(Form("%sdNdeta", name)),
fRebin(5), // Rebinning factor
fCutEdges(false),
fSymmetrice(true),
fCorrEmpty(true),
fUseROOTProj(false),
fTriggerEff(1),
+ fTriggerEff0(1),
fShapeCorr(0),
fListOfCentralities(0),
- fSNNString(0),
- fSysString(0),
- fCent(0),
- fCentAxis(0),
fNormalizationScheme(kFull),
- fSchemeString(0),
- fTriggerString(0),
- fFinalMCCorrFile(""),
- fIsESD(false)
+ fFinalMCCorrFile(""),
+ fSatelliteVertices(0),
+ fglobalempiricalcorrection(0),
+ fmeabsignalvscentr(0)
{
//
// Constructor
//
- DGUARD(0,0,"Named construction of AliBasedNdetaTask: %s", name);
+ DGUARD(fDebug, 3,"Named CTOR of AliBasedNdetaTask: %s", name);
+
+ fTriggerMask = AliAODForwardMult::kInel;
+ fMinIpZ = -10;
+ fMaxIpZ = +10;
fListOfCentralities = new TObjArray(1);
// Set the normalisation scheme
SetNormalizationScheme(kFull);
- // Set the trigger mask
- SetTriggerMask(AliAODForwardMult::kInel);
-
- // Output slot #1 writes into a TH1 container
- DefineOutput(1, TList::Class());
- DefineOutput(2, TList::Class());
}
-//____________________________________________________________________
-AliBasedNdetaTask::AliBasedNdetaTask(const AliBasedNdetaTask& o)
- : AliAnalysisTaskSE(o),
- fSums(o.fSums), // TList* - Container of sums
- fOutput(o.fOutput), // Container of output
- fVtxMin(o.fVtxMin), // Double_t - Minimum v_z
- fVtxMax(o.fVtxMax), // Double_t - Maximum v_z
- fTriggerMask(o.fTriggerMask),// Int_t - Trigger mask
- fRebin(o.fRebin), // Int_t - Rebinning factor
- fCutEdges(o.fCutEdges), // Bool_t - Whether to cut edges when rebinning
- fSymmetrice(o.fSymmetrice),
- fCorrEmpty(o.fCorrEmpty),
- fUseROOTProj(o.fUseROOTProj),
- fTriggerEff(o.fTriggerEff),
- fShapeCorr(o.fShapeCorr),
- fListOfCentralities(o.fListOfCentralities),
- fSNNString(o.fSNNString),
- fSysString(o.fSysString),
- fCent(o.fCent),
- fCentAxis(o.fCentAxis),
- fNormalizationScheme(o.fNormalizationScheme),
- fSchemeString(o.fSchemeString),
- fTriggerString(o.fTriggerString),
- fFinalMCCorrFile(o.fFinalMCCorrFile),
- fIsESD(o.fIsESD)
-{
- DGUARD(0,0,"Copy construction of AliBasedNdetaTask");
-}
//____________________________________________________________________
AliBasedNdetaTask::~AliBasedNdetaTask()
// Destructor
//
DGUARD(fDebug,3,"Destruction of AliBasedNdetaTask");
- if (fSums) {
- fSums->Delete();
- delete fSums;
- fSums = 0;
- }
- if (fOutput) {
- fOutput->Delete();
- delete fOutput;
- fOutput = 0;
- }
}
//________________________________________________________________________
}
}
-//________________________________________________________________________
-void
-AliBasedNdetaTask::SetCentralityAxis(UShort_t n, Short_t* bins)
-{
- DGUARD(fDebug,3,"Set centrality axis, %d bins", n);
- if (!fCentAxis) {
- fCentAxis = new TAxis();
- fCentAxis->SetName("centAxis");
- fCentAxis->SetTitle("Centrality [%]");
- }
- TArrayD dbins(n+1);
- for (UShort_t i = 0; i <= n; i++)
- dbins[i] = (bins[i] == 100 ? 100.1 : bins[i]);
- fCentAxis->Set(n, dbins.GetArray());
-}
-
//________________________________________________________________________
void
AliBasedNdetaTask::AddCentralityBin(UShort_t at, Short_t low, Short_t high)
// low Low cut
// high High cut
//
- DGUARD(fDebug,3,"Add a centrality bin [%f,%f] @ %d", low, high, at);
+ DGUARD(fDebug,3,"Add a centrality bin [%d,%d] @ %d", low, high, at);
CentralityBin* bin = MakeCentralityBin(GetName(), low, high);
+ if (!bin) {
+ Error("AddCentralityBin",
+ "Failed to create centrality bin for %s [%d,%d] @ %d",
+ GetName(), low, high, at);
+ return;
+ }
+ bin->SetSatelliteVertices(fSatelliteVertices);
bin->SetDebugLevel(fDebug);
fListOfCentralities->AddAtAndExpand(bin, at);
}
// Return:
// A newly created centrality bin
//
- DGUARD(fDebug,3,"Make a centrality bin [%f,%f]: %s", low, high, name);
+ DGUARD(fDebug,3,"Make a centrality bin %s [%d,%d]", name, low, high);
return new CentralityBin(name, low, high);
}
+
+#define TESTAPPEND(SCHEME,BIT,STRING) \
+ do { if (!(SCHEME & BIT)) break; \
+ if (!s.IsNull()) s.Append(","); s.Append(STRING); } while(false)
+
//________________________________________________________________________
-void
-AliBasedNdetaTask::SetNormalizationScheme(const char* what)
+const Char_t*
+AliBasedNdetaTask::NormalizationSchemeString(UShort_t scheme)
+{
+ // Create a string from normalization scheme bits
+ static TString s;
+ s = "";
+
+ if (scheme == kNone)
+ return s.Data();
+ if (scheme == kFull) {
+ s = "FULL";
+ return s.Data();
+ }
+ TESTAPPEND(scheme, kEventLevel, "EVENT");
+ TESTAPPEND(scheme, kShape, "SHAPE");
+ TESTAPPEND(scheme, kBackground, "BACKGROUND");
+ TESTAPPEND(scheme, kTriggerEfficiency, "TRIGGER");
+ TESTAPPEND(scheme, kZeroBin, "ZEROBIN");
+
+ return s.Data();
+}
+//________________________________________________________________________
+UShort_t
+AliBasedNdetaTask::ParseNormalizationScheme(const char* what)
{
- //
- // Set normalisation scheme
- //
- DGUARD(fDebug,3,"Set the normalization scheme: %s", what);
UShort_t scheme = 0;
TString twhat(what);
twhat.ToUpper();
TObjString* opt;
- TIter next(twhat.Tokenize(" ,|"));
+ TObjArray* token = twhat.Tokenize(" ,|");
+ TIter next(token);
while ((opt = static_cast<TObjString*>(next()))) {
TString s(opt->GetString());
if (s.IsNull()) continue;
else if (s.CompareTo("NONE") == 0) bit = kNone;
else if (s.CompareTo("ZEROBIN") == 0) bit = kZeroBin;
else
- Warning("SetNormalizationScheme", "Unknown option %s", s.Data());
+ ::Warning("SetNormalizationScheme", "Unknown option %s", s.Data());
if (add) scheme |= bit;
else scheme ^= bit;
}
- SetNormalizationScheme(scheme);
-}
-//________________________________________________________________________
-void
-AliBasedNdetaTask::SetNormalizationScheme(UShort_t scheme)
-{
- DGUARD(fDebug,3,"Set the normalization scheme: 0x%x", scheme);
- fNormalizationScheme = scheme;
- TString tit = "";
- if (scheme == kFull) tit = "FULL";
- else {
- if (scheme & kEventLevel) tit.Append("EVENT ");
- if (scheme & kShape) tit.Append("SHAPE ");
- if (scheme & kBackground) tit.Append("BACKGROUND ");
- if (scheme & kTriggerEfficiency) tit.Append("TRIGGER ");
- if (scheme & kZeroBin) tit.Append("ZEROBIN ");
- }
- tit = tit.Strip(TString::kBoth);
- if (!fSchemeString) fSchemeString = new TNamed("scheme", "");
- fSchemeString->SetTitle(tit);
- fSchemeString->SetUniqueID(fNormalizationScheme);
-}
+ delete token;
+ return scheme;
+}
//________________________________________________________________________
void
-AliBasedNdetaTask::SetTriggerMask(const char* mask)
+AliBasedNdetaTask::SetNormalizationScheme(const char* what)
{
//
- // Set the trigger maskl
+ // Set normalisation scheme
//
- // Parameters:
- // mask Trigger mask
- //
- DGUARD(fDebug,3,"Set the trigger mask: %s", mask);
- SetTriggerMask(AliAODForwardMult::MakeTriggerMask(mask));
+ DGUARD(fDebug,3,"Set the normalization scheme: %s", what);
+ SetNormalizationScheme(ParseNormalizationScheme(what));
}
//________________________________________________________________________
void
-AliBasedNdetaTask::SetTriggerMask(UShort_t mask)
-{
- DGUARD(fDebug,3,"Set the trigger mask: 0x%0x", mask);
- fTriggerMask = mask;
- TString tit(AliAODForwardMult::GetTriggerString(mask));
- tit = tit.Strip(TString::kBoth);
- if (!fTriggerString) fTriggerString = new TNamed("trigger", "");
- fTriggerString->SetTitle(tit);
- fTriggerString->SetUniqueID(fTriggerMask);
+AliBasedNdetaTask::SetNormalizationScheme(UShort_t scheme)
+{
+ DGUARD(fDebug,3,"Set the normalization scheme: 0x%x", scheme);
+ fNormalizationScheme = scheme;
}
-
//________________________________________________________________________
void
-AliBasedNdetaTask::SetShapeCorrection(const TH1* c)
+AliBasedNdetaTask::SetShapeCorrection(const TH2F* c)
{
//
// Set the shape correction (a.k.a., track correction) for selected
//
DGUARD(fDebug,3,"Set the shape correction: %p", c);
if (!c) return;
- fShapeCorr = static_cast<TH1*>(c->Clone());
+ fShapeCorr = static_cast<TH2F*>(c->Clone());
fShapeCorr->SetDirectory(0);
}
-
//________________________________________________________________________
void
-AliBasedNdetaTask::UserCreateOutputObjects()
+AliBasedNdetaTask::InitializeCentBins()
{
- //
- // Create output objects.
- //
- // This is called once per slave process
- //
- DGUARD(fDebug,1,"Create user ouput object");
- fSums = new TList;
- fSums->SetName(Form("%s_sums", GetName()));
- fSums->SetOwner();
+ if (fListOfCentralities->GetEntries() > 0) return;
// Automatically add 'all' centrality bin if nothing has been defined.
AddCentralityBin(0, 0, 0);
- if (fCentAxis && fCentAxis->GetNbins() > 0 && fCentAxis->GetXbins()) {
- const TArrayD* bins = fCentAxis->GetXbins();
- Int_t nbin = fCentAxis->GetNbins();
+ if (HasCentrality() && fCentAxis.GetXbins()) {
+ const TArrayD* bins = fCentAxis.GetXbins();
+ Int_t nbin = fCentAxis.GetNbins();
for (Int_t i = 0; i < nbin; i++)
AddCentralityBin(i+1, Short_t((*bins)[i]), Short_t((*bins)[i+1]));
}
- if (fCentAxis) fSums->Add(fCentAxis);
+}
+
+//________________________________________________________________________
+Bool_t
+AliBasedNdetaTask::Book()
+{
+ //
+ // Create output objects.
+ //
+ // This is called once per slave process
+ //
+ DGUARD(fDebug,1,"Create user ouput object");
+ fSums->Add(AliForwardUtil::MakeParameter("empirical",
+ fglobalempiricalcorrection != 0));
+ fSums->Add(AliForwardUtil::MakeParameter("scheme", fNormalizationScheme));
- // Centrality histogram
- fCent = new TH1D("cent", "Centrality", 100, 0, 100);
- fCent->SetDirectory(0);
- fCent->SetXTitle(0);
- fSums->Add(fCent);
+ // Make our centrality bins
+ InitializeCentBins();
// Loop over centrality bins
TIter next(fListOfCentralities);
CentralityBin* bin = 0;
while ((bin = static_cast<CentralityBin*>(next())))
- bin->CreateOutputObjects(fSums);
-
- // Check that we have an AOD input handler
- AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
- AliAODInputHandler* ah =
- dynamic_cast<AliAODInputHandler*>(am->GetInputEventHandler());
- if (!ah) {
- fIsESD = true;
- const char* depName = GetTitle();
- AliAnalysisTask* dep = am->GetTask(depName); // "Central"
- if (!dep) {
- AliErrorF("The multiplicity task %s wasn't added to the train", depName);
- fIsESD = false;
- }
- }
- if (!fIsESD && !ah) AliFatal("No AOD input handler set in analysis manager");
+ bin->CreateOutputObjects(fSums, fTriggerMask);
+
+ fmeabsignalvscentr=new TH2D("meanAbsSignalVsCentr",
+ "Mean absolute signal versus centrality",
+ 400, 0, 20, 100, 0, 100);
+ fSums->Add(fmeabsignalvscentr);
- // Post data for ALL output slots >0 here, to get at least an empty histogram
- PostData(1, fSums);
+ return true;
}
+
//____________________________________________________________________
-AliAODEvent*
-AliBasedNdetaTask::GetAODEvent(bool isESD)
+Bool_t
+AliBasedNdetaTask::CheckEvent(const AliAODForwardMult& fwd)
{
- // Get the AOD event
- //
- // Parameters:
- // isESD Whether we have ESD input or not
- //
- // Return:
- // Pointer to AOD event or null
- DGUARD(fDebug,1,"Get the AOD event (%s analysis)", isESD ? "ESD" : "AOD");
- AliAODEvent* aod = 0;
-
- if (!isESD) aod = dynamic_cast<AliAODEvent*>(InputEvent());
- else aod = AODEvent();
- if (!aod) AliError("Cannot get the AOD event");
-
- return aod;
+ AliBaseAODTask::CheckEvent(fwd);
+ // Here, we always return true, as the centrality bins will do
+ // their own checks on the events - this is needed for event
+ // normalization etc.
+ return true;
}
-
+
//____________________________________________________________________
-void
-AliBasedNdetaTask::UserExec(Option_t *)
+Bool_t
+AliBasedNdetaTask::Event(AliAODEvent& aod)
{
//
// Process a single event
//
// Main loop
DGUARD(fDebug,1,"Analyse the AOD event");
- AliAODEvent* aod = GetAODEvent(fIsESD);
- if (!aod) return;
-
- TObject* obj = aod->FindListObject("Forward");
- if (!obj) {
- AliWarning("No forward object found");
- return;
- }
- AliAODForwardMult* forward = static_cast<AliAODForwardMult*>(obj);
+ AliAODForwardMult* forward = GetForward(aod);
+ if (!forward) return false;
// Fill centrality histogram
- Float_t cent = forward->GetCentrality();
- fCent->Fill(cent);
+
+ Double_t vtx = forward->GetIpZ();
+ TH2D* data = GetHistogram(aod, false);
+ TH2D* dataMC = GetHistogram(aod, true);
+ if (!data) return false;
+
+ CheckEventData(vtx, data, dataMC);
+
+ if (!ApplyEmpiricalCorrection(forward,data))
+ return false;
- // Get the histogram(s)
- TH2D* data = GetHistogram(aod, false);
- TH2D* dataMC = GetHistogram(aod, true);
Bool_t isZero = ((fNormalizationScheme & kZeroBin) &&
!forward->IsTriggerBits(AliAODForwardMult::kNClusterGt0));
-
-
+ Bool_t taken = false;
+
// Loop over centrality bins
CentralityBin* allBin =
static_cast<CentralityBin*>(fListOfCentralities->At(0));
- allBin->ProcessEvent(forward, fTriggerMask, isZero,
- fVtxMin, fVtxMax, data, dataMC);
+ if (allBin->ProcessEvent(forward, fTriggerMask, isZero,
+ fMinIpZ, fMaxIpZ, data, dataMC)) taken = true;
// Find this centrality bin
- if (fCentAxis && fCentAxis->GetNbins() > 0) {
- Int_t icent = fCentAxis->FindBin(cent);
+ if (HasCentrality()) {
+ Double_t cent = forward->GetCentrality();
+ Int_t icent = fCentAxis.FindBin(cent);
CentralityBin* thisBin = 0;
- if (icent >= 1 && icent <= fCentAxis->GetNbins())
+ if (icent >= 1 && icent <= fCentAxis.GetNbins())
thisBin = static_cast<CentralityBin*>(fListOfCentralities->At(icent));
if (thisBin)
- thisBin->ProcessEvent(forward, fTriggerMask, isZero, fVtxMin, fVtxMax,
- data, dataMC);
+ if (thisBin->ProcessEvent(forward, fTriggerMask, isZero, fMinIpZ,
+ fMaxIpZ, data, dataMC)) taken = true;
}
-
- // Here, we get the update
- if (!fSNNString) {
- UShort_t sNN = forward->GetSNN();
- fSNNString = new TNamed("sNN", "");
- fSNNString->SetTitle(AliForwardUtil::CenterOfMassEnergyString(sNN));
- fSNNString->SetUniqueID(sNN);
- fSums->Add(fSNNString);
- UShort_t sys = forward->GetSystem();
- fSysString = new TNamed("sys", "");
- fSysString->SetTitle(AliForwardUtil::CollisionSystemString(sys));
- fSysString->SetUniqueID(sys);
- fSums->Add(fSysString);
-
- fSums->Add(fSchemeString);
- fSums->Add(fTriggerString);
+ return taken;
+}
- // Print();
- }
-
- PostData(1, fSums);
+//________________________________________________________________________
+void AliBasedNdetaTask::CheckEventData(Double_t,
+ TH2*,
+ TH2*)
+{
}
//________________________________________________________________________
h->SetMarkerStyle(marker);
h->SetMarkerSize(marker == 29 || marker == 30 ? 1.2 : 1);
h->SetFillStyle(0);
- h->SetYTitle(ytitle);
+ TString ytit;
+ if (ytitle && ytitle[0] != '\0') ytit = ytitle;
+ ytit = "#frac{1}{#it{N}}#frac{d#it{N}_{ch}}{d#it{#eta}}";
+ h->SetYTitle(ytit);
h->GetXaxis()->SetTitleFont(132);
h->GetXaxis()->SetLabelFont(132);
h->GetXaxis()->SetNdivisions(10);
}
}
}
+//________________________________________________________________________
+void
+AliBasedNdetaTask::ScaleToCoverage(TH1D* copy, const TH1D* norm)
+{
+ // Normalize to the acceptance -
+ // dndeta->Divide(accNorm);
+ for (Int_t i = 1; i <= copy->GetNbinsX(); i++) {
+ Double_t a = norm->GetBinContent(i);
+ if (a <= 0) {
+ copy->SetBinContent(i,0);
+ copy->SetBinError(i,0);
+ continue;
+ }
+ Double_t c = copy->GetBinContent(i);
+ Double_t e = copy->GetBinError(i);
+ copy->SetBinContent(i, c / a);
+ copy->SetBinError(i, e / a);
+ }
+}
//________________________________________________________________________
TH1D*
Int_t first = firstbin;
Int_t last = lastbin;
- if (first < 0) first = 0;
- else if (first >= yaxis->GetNbins()+1) first = yaxis->GetNbins();
+ if (first < 0) first = 1;
+ else if (first >= yaxis->GetNbins()+2) first = yaxis->GetNbins()+1;
if (last < 0) last = yaxis->GetNbins();
- else if (last > yaxis->GetNbins()+1) last = yaxis->GetNbins();
+ else if (last >= yaxis->GetNbins()+2) last = yaxis->GetNbins()+1;
if (last-first < 0) {
AliWarningGeneral("AliBasedNdetaTask",
Form("Nothing to project [%d,%d]", first, last));
}
// Loop over X bins
- // AliInfo(Form("Projecting bins [%d,%d] of %s", first, last, h->GetName()));
+ //DMSG(fDebug,3,"Projecting bins [%d,%d] of %s", first, last, h->GetName()));
Int_t ybins = (last-first+1);
for (Int_t xbin = 0; xbin <= xaxis->GetNbins()+1; xbin++) {
Double_t content = 0;
for (Int_t ybin = first; ybin <= last; ybin++) {
- Double_t c1 = h->GetCellContent(xbin, ybin);
- Double_t e1 = h->GetCellError(xbin, ybin);
+ Double_t c1 = h->GetBinContent(h->GetBin(xbin, ybin));
+ Double_t e1 = h->GetBinError(h->GetBin(xbin, ybin));
// Ignore empty bins
if (c1 < 1e-12) continue;
return ret;
}
-
+
//________________________________________________________________________
-void
-AliBasedNdetaTask::Terminate(Option_t *)
+Bool_t
+AliBasedNdetaTask::Finalize()
{
//
// Called at end of event processing..
// once at the end of the query
DGUARD(fDebug,1,"Process final merged results");
- fSums = dynamic_cast<TList*> (GetOutputData(1));
- if(!fSums) {
- AliError("Could not retrieve TList fSums");
- return;
- }
-
- fOutput = new TList;
- fOutput->SetName(Form("%s_result", GetName()));
- fOutput->SetOwner();
-
- fSNNString = static_cast<TNamed*>(fSums->FindObject("sNN"));
- fSysString = static_cast<TNamed*>(fSums->FindObject("sys"));
- fCentAxis = static_cast<TAxis*>(fSums->FindObject("centAxis"));
- fSchemeString = static_cast<TNamed*>(fSums->FindObject("scheme"));
- fTriggerString = static_cast<TNamed*>(fSums->FindObject("trigger"));
-
- if(fSysString && fSNNString &&
- fSysString->GetUniqueID() == AliForwardUtil::kPP)
- LoadNormalizationData(fSysString->GetUniqueID(),
- fSNNString->GetUniqueID());
+ UShort_t sNN;
+ UShort_t sys;
+ ULong_t trig;
+ UShort_t scheme;
+ AliForwardUtil::GetParameter(fSums->FindObject("sNN"), sNN);
+ AliForwardUtil::GetParameter(fSums->FindObject("sys"), sys);
+ AliForwardUtil::GetParameter(fSums->FindObject("trigger"), trig);
+ AliForwardUtil::GetParameter(fSums->FindObject("scheme"), scheme);
+
+ TAxis* cA = static_cast<TAxis*>(fSums->FindObject("centAxis"));
+ if (cA) cA->Copy(fCentAxis);
+
+ if(sNN > 0 && sys == AliForwardUtil::kPP)
+ LoadNormalizationData(sys, sNN);
+
+ InitializeCentBins();
+
// Print before we loop
Print();
Int_t style = GetMarker();
Int_t color = GetColor();
- AliInfo(Form("Marker style=%d, color=%d", style, color));
+ DMSG(fDebug,3,"Marker style=%d, color=%d", style, color);
while ((bin = static_cast<CentralityBin*>(next()))) {
-
- bin->End(fSums, fOutput, fNormalizationScheme, fShapeCorr, fTriggerEff,
+ bin->End(fSums, fResults, fNormalizationScheme, fShapeCorr,
+ fTriggerEff, fTriggerEff0,
fSymmetrice, fRebin, fUseROOTProj, fCorrEmpty, fCutEdges,
fTriggerMask, style, color, mclist, truthlist);
- if (fCentAxis && bin->IsAllBin()) continue;
- TH1* dndeta = bin->GetResult(0, false, "");
- TH1* dndetaSym = bin->GetResult(0, true, "");
- TH1* dndetaMC = bin->GetResult(0, false, "MC");
- TH1* dndetaMCSym = bin->GetResult(0, true, "MC");
+ if (HasCentrality() && bin->IsAllBin())
+ // If we have centrality bins, do not add the min-bias
+ // distribution to the output stack.
+ continue;
+ TH1* dndeta = bin->GetResult(0, false, "");
+ TH1* dndetaSym = fSymmetrice ? bin->GetResult(0, true, "") : 0;
+ TH1* dndetaMC = bin->GetResult(0, false, "MC", false);
+ TH1* dndetaMCSym = fSymmetrice ? bin->GetResult(0, true, "MC", false) : 0;
+ DMSG(fDebug,2,"Results: bare=%p sym=%p mcbare=%p mcsym=%p",
+ dndeta, dndetaSym, dndetaMC, dndetaMCSym);
if (dndeta) dndetaStack->Add(dndeta);
if (dndetaSym) dndetaStack->Add(dndetaSym);
if (dndetaMC) dndetaMCStack->Add(dndetaMC);
if (dndetaMCSym) dndetaMCStack->Add(dndetaMCSym);
if (fRebin > 1) {
- dndeta = bin->GetResult(fRebin, false, "");
- dndetaSym = bin->GetResult(fRebin, true, "");
- dndetaMC = bin->GetResult(fRebin, false, "MC");
- dndetaMCSym = bin->GetResult(fRebin, true, "MC");
+ dndeta = bin->GetResult(fRebin, false, "");
+ dndetaSym = fSymmetrice ? bin->GetResult(fRebin, true, "") : 0;
+ dndetaMC = bin->GetResult(fRebin, false, "MC", false);
+ dndetaMCSym = fSymmetrice ? bin->GetResult(fRebin, true, "MC", false): 0;
if (dndeta) dndetaStackRebin->Add(dndeta);
if (dndetaSym) dndetaStackRebin->Add(dndetaSym);
if (dndetaMC) dndetaMCStackRebin->Add(dndetaMC);
}
}
// Output the stack
- fOutput->Add(dndetaStack);
+ fResults->Add(dndetaStack);
// If available output rebinned stack
- if (!dndetaStackRebin->GetHists() ||
- dndetaStackRebin->GetHists()->GetEntries() <= 0) {
+ if (fRebin > 0 && (!dndetaStackRebin->GetHists() ||
+ dndetaStackRebin->GetHists()->GetEntries() <= 0)) {
AliWarning("No rebinned histograms found");
delete dndetaStackRebin;
dndetaStackRebin = 0;
}
- if (dndetaStackRebin) fOutput->Add(dndetaStackRebin);
+ if (dndetaStackRebin) fResults->Add(dndetaStackRebin);
// If available, output track-ref stack
if (!dndetaMCStack->GetHists() ||
dndetaMCStack->GetHists()->GetEntries() <= 0) {
- AliWarning("No MC histograms found");
+ // AliWarning("No MC histograms found");
delete dndetaMCStack;
dndetaMCStack = 0;
}
- if (dndetaMCStack) fOutput->Add(dndetaMCStack);
+ if (dndetaMCStack) fResults->Add(dndetaMCStack);
// If available, output rebinned track-ref stack
- if (!dndetaMCStackRebin->GetHists() ||
- dndetaMCStackRebin->GetHists()->GetEntries() <= 0) {
+ if ((fRebin > 0 && dndetaMCStack)
+ && (!dndetaMCStackRebin->GetHists() ||
+ dndetaMCStackRebin->GetHists()->GetEntries() <= 0)) {
AliWarning("No rebinned MC histograms found");
delete dndetaMCStackRebin;
dndetaMCStackRebin = 0;
}
- if (dndetaMCStackRebin) fOutput->Add(dndetaMCStackRebin);
+ if (dndetaMCStackRebin) fResults->Add(dndetaMCStackRebin);
// Output collision energy string
- if (fSNNString) fOutput->Add(fSNNString->Clone());
+ if (sNN > 0) {
+ TNamed* sNNObj = new TNamed("sNN",
+ AliForwardUtil::CenterOfMassEnergyString(sNN));
+ sNNObj->SetUniqueID(sNN);
+ fResults->Add(sNNObj); // fSNNString->Clone());
+ }
// Output collision system string
- if (fSysString) fOutput->Add(fSysString->Clone());
+ if (sys > 0) {
+ TNamed* sysObj = new TNamed("sys",
+ AliForwardUtil::CollisionSystemString(sys));
+ sysObj->SetUniqueID(sys);
+ fResults->Add(sysObj); // fSysString->Clone());
+ }
// Output centrality axis
- if (fCentAxis) fOutput->Add(fCentAxis);
+ fResults->Add(&fCentAxis);
// Output trigger string
- if (fTriggerString) fOutput->Add(fTriggerString->Clone());
+ if (trig) {
+ TNamed* maskObj = new TNamed("trigger",
+ AliAODForwardMult::GetTriggerString(trig));
+ maskObj->SetUniqueID(trig);
+ fResults->Add(maskObj); // fTriggerString->Clone());
+ }
// Normalization string
- if (fSchemeString) fOutput->Add(fSchemeString->Clone());
+ if (scheme > 0) {
+ TNamed* schemeObj = new TNamed("scheme",
+ NormalizationSchemeString(scheme));
+ schemeObj->SetUniqueID(scheme);
+ fResults->Add(schemeObj); // fSchemeString->Clone());
+ }
// Output vertex axis
- TAxis* vtxAxis = new TAxis(1,fVtxMin,fVtxMax);
+ TAxis* vtxAxis = new TAxis(1,fMinIpZ,fMaxIpZ);
vtxAxis->SetName("vtxAxis");
- vtxAxis->SetTitle(Form("v_{z}#in[%+5.1f,%+5.1f]cm", fVtxMin,fVtxMax));
- fOutput->Add(vtxAxis);
+ vtxAxis->SetTitle(Form("v_{z}#in[%+5.1f,%+5.1f]cm", fMinIpZ,fMaxIpZ));
+ fResults->Add(vtxAxis);
- // Output trigger efficiency and shape correction
- fOutput->Add(new TParameter<Double_t>("triggerEff", fTriggerEff));
- if (fShapeCorr) fOutput->Add(fShapeCorr);
+ // Output trigger efficiency and shape correction
+ fResults->Add(AliForwardUtil::MakeParameter("triggerEff", fTriggerEff));
+ fResults->Add(AliForwardUtil::MakeParameter("triggerEff0", fTriggerEff0));
+ if (fShapeCorr) fResults->Add(fShapeCorr);
TNamed* options = new TNamed("options","");
TString str;
str.Append(Form("Empty bins %scorrected for, ", fCorrEmpty ? "" : "not "));
str.Append(Form("TH2::ProjectionX %sused", fUseROOTProj ? "" : "not "));
options->SetTitle(str);
- fOutput->Add(options);
+ fResults->Add(options);
- PostData(2, fOutput);
+ return true;
}
//________________________________________________________________________
void
if(energy == 7000) snn.Form("7000");
if(energy == 2750) snn.Form("2750");
- if(fShapeCorr && (fTriggerEff != 1)) {
- AliInfo("Objects already set for normalization - no action taken");
+ // Check if shape correction/trigger efficiency was requsted and not
+ // already set
+ Bool_t needShape = ((fNormalizationScheme & kShape) && !fShapeCorr);
+ Bool_t needEff = ((fNormalizationScheme & kTriggerEfficiency) &&
+ ((1 - fTriggerEff) < 1e-6) && fTriggerEff > 0);
+ if (needShape) DMSG(fDebug, 0, "Will load shape correction");
+ if (needEff) DMSG(fDebug, 0, "Will load trigger efficiency, was=%f, %f",
+ fTriggerEff, fTriggerEff0);
+ if(!needShape) { // && !needEff) {
+ DMSG(fDebug,1,"Objects already set for normalization - no action taken");
return;
}
-
- TFile* fin = TFile::Open(Form("$ALICE_ROOT/PWGLF/FORWARD/corrections/"
- "Normalization/normalizationHists_%s_%s.root",
- type.Data(),snn.Data()));
+
+ TString fname(Form("$ALICE_ROOT/PWGLF/FORWARD/corrections/"
+ "Normalization/normalizationHists_%s_%s.root",
+ type.Data(),snn.Data()));
+ AliWarningF("Using old-style corrections from %s", fname.Data());
+ TFile* fin = TFile::Open(fname, "READ");
if(!fin) {
- AliWarning(Form("no file for normalization of %d/%d", sys, energy));
+ AliWarningF("no file for normalization of %d/%d (%s)",
+ sys, energy, fname.Data());
return;
}
// Shape correction
- if ((fNormalizationScheme & kShape) && !fShapeCorr) {
+ if (needShape) {
TString trigName("All");
if (fTriggerMask == AliAODForwardMult::kInel ||
fTriggerMask == AliAODForwardMult::kNClusterGt0)
}
// Trigger efficiency
- TString effName(Form("%sTriggerEff",
- fTriggerMask == AliAODForwardMult::kInel ? "inel" :
- fTriggerMask == AliAODForwardMult::kNSD ? "nsd" :
- fTriggerMask == AliAODForwardMult::kInelGt0 ?
- "inelgt0" : "all"));
- TParameter<float>* eff = 0;
- if (fNormalizationScheme & kTriggerEfficiency)
- eff = static_cast<TParameter<float>*>(fin->Get(effName));
- Double_t trigEff = eff ? eff->GetVal() : 1;
- if (fTriggerEff != 1) SetTriggerEff(trigEff);
- if (fTriggerEff < 0) fTriggerEff = 1;
-
+ if (needEff) {
+ TString effName(Form("%sTriggerEff",
+ fTriggerMask == AliAODForwardMult::kInel ? "inel" :
+ fTriggerMask == AliAODForwardMult::kNSD ? "nsd" :
+ fTriggerMask == AliAODForwardMult::kInelGt0 ?
+ "inelgt0" : "all"));
+ Double_t trigEff = 1;
+ TObject* eff = fin->Get(effName);
+ if (eff) AliForwardUtil::GetParameter(eff, trigEff);
+
+ if (trigEff <= 0)
+ AliWarningF("Retrieved trigger efficiency %s is %f<=0, ignoring",
+ effName.Data(), trigEff);
+ else
+ SetTriggerEff(trigEff);
+
+ // Trigger efficiency
+ TString eff0Name(effName);
+ eff0Name.Append("0");
+
+ Double_t trigEff0 = 1;
+ TObject* eff0 = fin->Get(eff0Name);
+ if (eff0) AliForwardUtil::GetParameter(eff, trigEff0);
+ if (trigEff0 < 0)
+ AliWarningF("Retrieved trigger efficiency %s is %f<0, ignoring",
+ eff0Name.Data(), trigEff0);
+ else
+ SetTriggerEff0(trigEff0);
+ }
+
// TEMPORARY FIX
// Rescale the shape correction by the trigger efficiency
if (fShapeCorr) {
AliWarning(Form("Rescaling shape correction by trigger efficiency: "
- "1/E_X=1/%f", trigEff));
- fShapeCorr->Scale(1. / trigEff);
+ "1/E_X=1/%f", fTriggerEff));
+ fShapeCorr->Scale(1. / fTriggerEff);
}
+ if (fin) fin->Close();
// Print - out
- if (fShapeCorr && fTriggerEff) AliInfo("Loaded objects for normalization.");
+ if (fDebug > 1 && fShapeCorr && fTriggerEff)
+ DMSG(fDebug, 1, "Loaded objects for normalization.");
}
+#define PF(N,V,...) \
+ AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
+#define PFB(N,FLAG) \
+ do { \
+ AliForwardUtil::PrintName(N); \
+ std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
+ } while(false)
+#define PFV(N,VALUE) \
+ do { \
+ AliForwardUtil::PrintName(N); \
+ std::cout << (VALUE) << std::endl; } while(false)
+
//________________________________________________________________________
void
-AliBasedNdetaTask::Print(Option_t*) const
+AliBasedNdetaTask::Print(Option_t* option) const
{
//
// Print information
//
- std::cout << this->ClassName() << ": " << this->GetName() << "\n"
- << std::boolalpha
- << " Trigger: " << (fTriggerString ?
- fTriggerString->GetTitle() :
- "none") << "\n"
- << " Vertex range: [" << fVtxMin << ":" << fVtxMax << "]\n"
- << " Rebin factor: " << fRebin << "\n"
- << " Cut edges: " << fCutEdges << "\n"
- << " Symmertrice: " << fSymmetrice << "\n"
- << " Use TH2::ProjectionX: " << fUseROOTProj << "\n"
- << " Correct for empty: " << fCorrEmpty << "\n"
- << " Normalization scheme: " << (fSchemeString ?
- fSchemeString->GetTitle() :
- "none") <<"\n"
- << " Trigger efficiency: " << fTriggerEff << "\n"
- << " Shape correction: " << (fShapeCorr ?
- fShapeCorr->GetName() :
- "none") << "\n"
- << " sqrt(s_NN): " << (fSNNString ?
- fSNNString->GetTitle() :
- "unknown") << "\n"
- << " Collision system: " << (fSysString ?
- fSysString->GetTitle() :
- "unknown") << "\n"
- << " Centrality bins: " << (fCentAxis ? "" : "none");
- if (fCentAxis) {
- Int_t nBins = fCentAxis->GetNbins();
- const Double_t* bins = fCentAxis->GetXbins()->GetArray();
- for (Int_t i = 0; i <= nBins; i++)
- std::cout << (i==0 ? "" : "-") << bins[i];
- }
- std::cout << std::noboolalpha << std::endl;
-
+ AliBaseAODTask::Print(option);
+ TString schemeString(NormalizationSchemeString(fNormalizationScheme));
+
+ gROOT->IncreaseDirLevel();
+ PFV("Rebin factor", fRebin);
+ PFB("Cut edges", fCutEdges);
+ PFB("Symmertrice", fSymmetrice);
+ PFB("Use TH2::ProjectionX", fUseROOTProj);
+ PFB("Correct for empty", fCorrEmpty);
+ PFV("Normalization scheme", schemeString );
+ PFV("Trigger efficiency", fTriggerEff);
+ PFV("Bin-0 Trigger efficiency", fTriggerEff0);
+ PFV("Shape correction", (fShapeCorr?fShapeCorr->GetName():"none"));;
+ gROOT->DecreaseDirLevel();
}
//________________________________________________________________________
TH1D* tmp = static_cast<TH1D*>(h->Clone(Form("%s_rebin%02d",
h->GetName(), rebin)));
tmp->Rebin(rebin);
+ // Hist should be reset, as it otherwise messes with the cutEdges option
+ tmp->Reset();
tmp->SetDirectory(0);
// The new number of bins
for(Int_t j = 1; j<=rebin;j++) {
Int_t bin = (i-1)*rebin + j;
Double_t c = h->GetBinContent(bin);
- if (c <= 0) continue;
+ if (c <= 0) {
+ continue; // old TODO: check
+ //content = -1; // new
+ //break; // also new
+ }
if (cutEdges) {
if (h->GetBinContent(bin+1)<=0 ||
//____________________________________________________________________
TString
-AliBasedNdetaTask::Sum::GetHistName(Int_t what) const
+AliBasedNdetaTask::Sum::GetHistName(const char* /*name*/,
+ Int_t what, const char* post)
{
- TString n(GetName());
- if (what == 1) n.Append("0");
- else if (what == 2) n.Append("Events");
- const char* postfix = GetTitle();
- if (postfix && postfix[0] != '\0') n.Append(postfix);
+ TString n;
+ switch (what) {
+ case 0: n = "sum"; break;
+ case 1: n = "sum0"; break;
+ case 2: n = "events"; break;
+ }
+ if (post && post[0] != '\0') n.Append(post);
return n;
}
+//____________________________________________________________________
+TString
+AliBasedNdetaTask::Sum::GetHistName(Int_t what) const
+{
+ return GetHistName(GetName(), what, GetTitle());
+}
+
//____________________________________________________________________
void
AliBasedNdetaTask::Sum::Add(const TH2D* data, Bool_t isZero)
//____________________________________________________________________
TH2D*
-AliBasedNdetaTask::Sum::GetSum(const TList* input,
- TList* output,
- Double_t& ntotal,
- Double_t epsilon0,
- Double_t epsilon,
- Int_t marker,
- Bool_t rootProj,
- Bool_t corrEmpty) const
+AliBasedNdetaTask::Sum::CalcSum(TList* output,
+ Double_t& ntotal,
+ Double_t epsilon0,
+ Double_t epsilon,
+ Int_t marker,
+ Bool_t rootProj,
+ Bool_t corrEmpty) const
{
- DGUARD(fDebug,2,"Get sums from input list %s", input->GetName());
- TH2D* sum = static_cast<TH2D*>(input->FindObject(GetHistName(0)));
- TH2D* sum0 = static_cast<TH2D*>(input->FindObject(GetHistName(1)));
- TH1I* events = static_cast<TH1I*>(input->FindObject(GetHistName(2)));
- if (!sum || !sum0 || !events) {
- AliWarning(Form("Failed to find one or more histograms: "
- "%s (%p) %s (%p) %s (%p)",
- GetHistName(0).Data(), sum,
- GetHistName(1).Data(), sum0,
- GetHistName(2).Data(), events));
- return 0;
- }
+ DGUARD(fDebug,2,"Calculating final summed histogram %s", fSum->GetName());
- TH2D* ret = static_cast<TH2D*>(sum->Clone(sum->GetName()));
+ // The return value `ret' is not scaled in anyway
+ TH2D* ret = static_cast<TH2D*>(fSum->Clone(fSum->GetName()));
ret->SetDirectory(0);
- ret->Reset();
- Int_t n = Int_t(events->GetBinContent(1));
- Int_t n0 = Int_t(events->GetBinContent(2));
-
- AliInfo(Form("Adding histograms %s and %s with weights %f and %f resp.",
- sum0->GetName(), sum->GetName(), 1./epsilon, 1./epsilon0));
- // Generate merged histogram
- ret->Add(sum0, sum, 1. / epsilon0, 1. / epsilon);
- ntotal = n / epsilon + n0 / epsilon0;
+ Int_t n = Int_t(fEvents->GetBinContent(1));
+ Int_t n0 = Int_t(fEvents->GetBinContent(2));
+ ntotal = n;
+ if (n0 > 0) {
+ ret->Reset();
+ DMSG(fDebug,1,
+ "Adding histograms %s(%d) and %s(%d) w/weights %f and %f resp.",
+ fSum0->GetName(), n, fSum->GetName(), n0, 1./epsilon,1./epsilon0);
+ ret->Add(fSum0, fSum, 1. / epsilon0, 1. / epsilon);
+ ntotal = n / epsilon + n0 / epsilon0;
+ }
TList* out = new TList;
out->SetOwner();
output->Add(out);
// Now make copies, normalize them, and store in output list
- TH2D* sumCopy = static_cast<TH2D*>(sum->Clone("sum"));
- TH2D* sum0Copy = static_cast<TH2D*>(sum0->Clone("sum0"));
+ // Note, these are the only ones normalized here
+ // These are mainly for diagnostics
+ TH2D* sumCopy = static_cast<TH2D*>(fSum->Clone("sum"));
+ TH2D* sum0Copy = static_cast<TH2D*>(fSum0->Clone("sum0"));
TH2D* retCopy = static_cast<TH2D*>(ret->Clone("sumAll"));
sumCopy->SetMarkerStyle(FlipHollowStyle(marker));
sumCopy->SetDirectory(0);
retCopy->SetMarkerStyle(marker);
retCopy->SetDirectory(0);
- TH1D* norm = ProjectX(sum, "norm", 0, 0, rootProj, corrEmpty, false);
- TH1D* norm0 = ProjectX(sum0, "norm0", 0, 0, rootProj, corrEmpty, false);
- TH1D* normAll = ProjectX(ret, "normAll", 0, 0, rootProj, corrEmpty, false);
+ Int_t nY = fSum->GetNbinsY();
+ Int_t o = 0; // nY+1;
+ TH1D* norm = ProjectX(fSum, "norm", o, o, rootProj, corrEmpty, false);
+ TH1D* norm0 = ProjectX(fSum0, "norm0", o, o, rootProj, corrEmpty, false);
+ TH1D* normAll = ProjectX(ret, "normAll", o, o, rootProj, corrEmpty, false);
+ norm->SetTitle("#eta coverage - >0-bin");
+ norm0->SetTitle("#eta coverage - 0-bin");
+ normAll->SetTitle("#eta coverage");
norm->SetDirectory(0);
norm0->SetDirectory(0);
normAll->SetDirectory(0);
- ScaleToCoverage(sumCopy, norm);
- ScaleToCoverage(sum0Copy, norm0);
- ScaleToCoverage(retCopy, normAll);
-
- Int_t nY = sum->GetNbinsY();
- TH1D* sumCopyPx = ProjectX(sumCopy, "average", 1, nY,rootProj,corrEmpty);
- TH1D* sum0CopyPx = ProjectX(sum0Copy, "average0", 1, nY,rootProj,corrEmpty);
- TH1D* retCopyPx = ProjectX(retCopy, "averageAll", 1, nY,rootProj,corrEmpty);
- sumCopyPx->SetDirectory(0);
+ TH1D* sumCopyPx = ProjectX(sumCopy, "average", 1,nY,rootProj,corrEmpty);
+ TH1D* sum0CopyPx = ProjectX(sum0Copy, "average0", 1,nY,rootProj,corrEmpty);
+ TH1D* retCopyPx = ProjectX(retCopy, "averageAll", 1,nY,rootProj,corrEmpty);
+ sumCopyPx-> SetTitle(Form("#sum_{i}^{N_{#phi}}%s", sumCopy->GetTitle()));
+ sum0CopyPx->SetTitle(Form("#sum_{i}^{N_{#phi}}%s", sum0Copy->GetTitle()));
+ retCopyPx-> SetTitle(Form("#sum_{i}^{N_{#phi}}%s", retCopy->GetTitle()));
+ sumCopyPx-> SetDirectory(0);
sum0CopyPx->SetDirectory(0);
- retCopyPx->SetDirectory(0);
+ retCopyPx-> SetDirectory(0);
+
+ TH1D* phi = ProjectX(fSum, "phi", nY+1,nY+1,rootProj,corrEmpty,false);
+ TH1D* phi0 = ProjectX(fSum0, "phi0", nY+1,nY+1,rootProj,corrEmpty,false);
+ TH1D* phiAll = ProjectX(ret, "phiAll", nY+1,nY+1,rootProj,corrEmpty,false);
+ phi ->SetTitle("#phi acceptance from dead strips - >0-bin");
+ phi0 ->SetTitle("#phi acceptance from dead strips - 0-bin");
+ phiAll->SetTitle("#phi acceptance from dead strips");
+ phi ->SetDirectory(0);
+ phi0 ->SetDirectory(0);
+ phiAll->SetDirectory(0);
+
+ const TH1D* cov = (corrEmpty ? norm : phi);
+ const TH1D* cov0 = (corrEmpty ? norm0 : phi0);
+ const TH1D* covAll = (corrEmpty ? normAll : phiAll);
+
+ // Here, we scale to the coverage (or phi acceptance)
+ ScaleToCoverage(sumCopy, cov);
+ ScaleToCoverage(sum0Copy, cov0);
+ ScaleToCoverage(retCopy, covAll);
// Scale our 1D histograms
- sumCopyPx->Scale(1., "width");
+ sumCopyPx ->Scale(1., "width");
sum0CopyPx->Scale(1., "width");
- retCopyPx->Scale(1., "width");
+ retCopyPx ->Scale(1., "width");
- AliInfo(Form("Maximum %f,%f changed to %f", sumCopyPx->GetMaximum(),
- sum0CopyPx->GetMaximum(), retCopyPx->GetMaximum()));
+ DMSG(fDebug,2,"Maximum %f,%f changed to %f", sumCopyPx->GetMaximum(),
+ sum0CopyPx->GetMaximum(), retCopyPx->GetMaximum());
// Scale the normalization - they should be 1 at the maximum
- norm->Scale(n > 0 ? 1. / n : 1);
- norm0->Scale(n0 > 0 ? 1. / n0 : 1);
+ norm ->Scale(n > 0 ? 1. / n : 1);
+ norm0 ->Scale(n0 > 0 ? 1. / n0 : 1);
normAll->Scale(ntotal > 0 ? 1. / ntotal : 1);
+ // Scale the normalization - they should be 1 at the maximum
+ phi ->Scale(n > 0 ? 1. / n : 1);
+ phi0 ->Scale(n0 > 0 ? 1. / n0 : 1);
+ phiAll->Scale(ntotal > 0 ? 1. / ntotal : 1);
+
out->Add(sumCopy);
out->Add(sum0Copy);
out->Add(retCopy);
out->Add(norm);
out->Add(norm0);
out->Add(normAll);
+ out->Add(phi);
+ out->Add(phi0);
+ out->Add(phiAll);
+
+ if (fDebug >= 1) {
+ if (n0 > 0)
+ DMSG(fDebug,1,"Returning (1/%f * %s + 1/%f * %s), "
+ "1/%f * %d + 1/%f * %d = %d",
+ epsilon0, fSum0->GetName(), epsilon, fSum->GetName(),
+ epsilon0, n0, epsilon, n, int(ntotal));
+ else
+ DMSG(fDebug,1, "Returning %s, %d", fSum->GetName(), int(ntotal));
+ }
- AliInfo(Form("Returning (1/%f * %s + 1/%f * %s), "
- "1/%f * %d + 1/%f * %d = %d",
- epsilon0, sum0->GetName(), epsilon, sum->GetName(),
- epsilon0, n0, epsilon, n, int(ntotal)));
#if 0
for (Int_t i = 1; i <= ret->GetNbinsX(); i++) {
Double_t nc = sum->GetBinContent(i, 0);
fSum(0),
fSumMC(0),
fTriggers(0),
+ fStatus(0),
fLow(0),
fHigh(0),
- fDoFinalMCCorrection(false),
+ fDoFinalMCCorrection(false),
+ fSatelliteVertices(false),
fDebug(0)
{
//
// Constructor
//
- DGUARD(0,0,"Centrality bin default construction");
+ DGUARD(fDebug,3,"Default CTOR of AliBasedNdeta::CentralityBin");
}
//____________________________________________________________________
AliBasedNdetaTask::CentralityBin::CentralityBin(const char* name,
fSum(0),
fSumMC(0),
fTriggers(0),
+ fStatus(0),
fLow(low),
fHigh(high),
fDoFinalMCCorrection(false),
+ fSatelliteVertices(false),
fDebug(0)
{
//
// low Lower centrality cut in percent
// high Upper centrality cut in percent
//
- DGUARD(0,0,"Named Centrality bin construction: %s [%d,%d]", name);
+ DGUARD(fDebug,3,"Named CTOR of AliBasedNdeta::CentralityBin: %s [%3d,%3d]",
+ name,low,high);
if (low <= 0 && high <= 0) {
fLow = 0;
fHigh = 0;
fSum(o.fSum),
fSumMC(o.fSumMC),
fTriggers(o.fTriggers),
+ fStatus(o.fStatus),
fLow(o.fLow),
fHigh(o.fHigh),
- fDoFinalMCCorrection(o.fDoFinalMCCorrection),
+ fDoFinalMCCorrection(o.fDoFinalMCCorrection),
+ fSatelliteVertices(o.fSatelliteVertices),
fDebug(o.fDebug)
{
//
// Parameters:
// other Object to copy from
//
- DGUARD(0,0,"Copy Centrality bin construction");
+ DGUARD(fDebug,3,"Copy CTOR of AliBasedNdeta::CentralityBin");
}
//____________________________________________________________________
AliBasedNdetaTask::CentralityBin::~CentralityBin()
//
// Destructor
//
- DGUARD(fDebug,3,"Centrality bin desctruction");
- if (fSums) fSums->Delete();
- if (fOutput) fOutput->Delete();
+ DGUARD(fDebug,3,"DTOR of AliBasedNdeta::CentralityBin");
+
+ // if (fSums) fSums->Delete();
+ // if (fOutput) fOutput->Delete();
}
//____________________________________________________________________
fSum = o.fSum;
fSumMC = o.fSumMC;
fTriggers = o.fTriggers;
+ fStatus = o.fStatus;
fLow = o.fLow;
fHigh = o.fHigh;
fDoFinalMCCorrection = o.fDoFinalMCCorrection;
+ fSatelliteVertices = o.fSatelliteVertices;
return *this;
}
}
//____________________________________________________________________
void
-AliBasedNdetaTask::CentralityBin::CreateOutputObjects(TList* dir)
+AliBasedNdetaTask::CentralityBin::CreateOutputObjects(TList* dir, Int_t mask)
{
//
// Create output objects
fSums->SetOwner();
dir->Add(fSums);
- fTriggers = AliAODForwardMult::MakeTriggerHistogram("triggers");
+ fTriggers = AliAODForwardMult::MakeTriggerHistogram("triggers", mask);
fTriggers->SetDirectory(0);
+
+ fStatus = AliAODForwardMult::MakeStatusHistogram("status");
+ fStatus->SetDirectory(0);
+
fSums->Add(fTriggers);
+ fSums->Add(fStatus);
}
//____________________________________________________________________
void
if (fSumMC) fSumMC->fDebug = lvl;
}
+//____________________________________________________________________
+Bool_t
+AliBasedNdetaTask::CentralityBin::ReadSum(TList* list, bool mc)
+{
+ const char* post = (mc ? "MC" : "");
+ TString sn = Sum::GetHistName(GetName(),0,post);
+ TString sn0 = Sum::GetHistName(GetName(),1,post);
+ TString ev = Sum::GetHistName(GetName(),2,post);
+ TH2D* sum = static_cast<TH2D*>(list->FindObject(sn));
+ TH2D* sum0 = static_cast<TH2D*>(list->FindObject(sn0));
+ TH1I* events = static_cast<TH1I*>(list->FindObject(ev));
+ if (!sum || !sum0 || !events) {
+ if (!mc)
+ AliWarningF("Failed to find one or more histograms: "
+ "%s (%p) %s (%p) %s (%p)",
+ sn.Data(), sum,
+ sn0.Data(), sum0,
+ ev.Data(), events);
+ return false;
+ }
+ Sum* ret = new Sum(GetName(), post);
+ ret->fSum = sum;
+ ret->fSum0 = sum0;
+ ret->fEvents = events;
+ ret->fDebug = fDebug;
+ if (mc) fSumMC = ret;
+ else fSum = ret;
+
+ return true;
+}
+
//____________________________________________________________________
void
AliBasedNdetaTask::CentralityBin::CreateSums(const TH2D* data, const TH2D* mc)
// data Data histogram to clone
// mc (optional) MC histogram to clone
//
- DGUARD(fDebug,1,"Create centrality bin sums from %s", data->GetName());
+ DGUARD(fDebug,1,"Create centrality bin sums from %s",
+ data ? data->GetName() : "(null)");
if (data) {
fSum = new Sum(GetName(),"");
fSum->Init(fSums, data, GetColor());
DGUARD(fDebug,2,"Check the event");
// We do not check for centrality here - it's already done
- return forward->CheckEvent(triggerMask, vzMin, vzMax, 0, 0, fTriggers);
+ return forward->CheckEvent(triggerMask, vzMin, vzMax, 0, 0,
+ fTriggers, fStatus);
}
//____________________________________________________________________
-void
+Bool_t
AliBasedNdetaTask::CentralityBin::ProcessEvent(const AliAODForwardMult* forward,
Int_t triggerMask, Bool_t isZero,
Double_t vzMin, Double_t vzMax,
// mc MC histogram
//
DGUARD(fDebug,1,"Process one event for %s a given centrality bin",
- data->GetName());
- if (!CheckEvent(forward, triggerMask, vzMin, vzMax)) return;
- if (!data) return;
+ data ? data->GetName() : "(null)");
+ if (!CheckEvent(forward, triggerMask, vzMin, vzMax)) return false;
+ if (!data) return false;
if (!fSum) CreateSums(data, mc);
fSum->Add(data, isZero);
if (mc) fSumMC->Add(mc, isZero);
+
+ return true;
}
//________________________________________________________________________
AliBasedNdetaTask::CentralityBin::Normalization(const TH1I& t,
UShort_t scheme,
Double_t trigEff,
- Double_t& ntotal) const
+ Double_t& ntotal,
+ TString* text) const
{
//
// Calculate normalization
// trigEff From MC
// ntotal On return, contains the number of events.
//
- DGUARD(fDebug,1,"Normalize centrality bin with %s", t.GetName());
+ DGUARD(fDebug,1,"Normalize centrality bin %s [%3d-%3d%%] with %s",
+ GetName(), fLow, fHigh, t.GetName());
Double_t nAll = t.GetBinContent(AliAODForwardMult::kBinAll);
Double_t nB = t.GetBinContent(AliAODForwardMult::kBinB);
Double_t nA = t.GetBinContent(AliAODForwardMult::kBinA);
Double_t nOffline = t.GetBinContent(AliAODForwardMult::kBinOffline);
Double_t nTriggered = t.GetBinContent(AliAODForwardMult::kWithTrigger);
Double_t nWithVertex = t.GetBinContent(AliAODForwardMult::kWithVertex);
- Double_t nAccepted = ntotal; // t.GetBinContent(AliAODForwardMult::kAccepted);
+ Double_t nAccepted = ntotal;
ntotal = 0;
if (nTriggered <= 0.1) {
Double_t scaler = 1;
Double_t beta = nA + nC - 2*nE;
- if (scheme & kEventLevel && !(scheme & kZeroBin)) {
- ntotal = nAccepted / vtxEff;
- scaler = vtxEff;
- AliInfo(Form("Calculating event normalisation as\n"
- " N = N_A * N_T / N_V = %d * %d / %d = %f (%f)",
- Int_t(nAccepted), Int_t(nTriggered), Int_t(nWithVertex),
- ntotal, scaler));
-
- if (scheme & kBackground) {
- // 1 E_V E_V
- // s = --------- = ------------- = ------------
- // 1 - beta 1 - beta E_V 1 - beta N_V
- // --- ---- -------- ---- ---
- // E_V N_V N_V N_V N_T
- //
- // E_V
- // = ------------
- // 1 - beta
- // ----
- // N_T
- //
- ntotal -= nAccepted * beta / nWithVertex;
- // This one is direct and correct.
- // scaler = 1. / (1. / vtxEff - beta / nWithVertex);
- // A simpler expresion
- scaler /= (1 - beta / nTriggered); // 0.831631 -> 0.780689
- AliInfo(Form("Calculating event normalisation as\n"
- " beta = N_a + N_c + 2 N_e = %d + %d - 2 * %d = %d\n"
- " N = N - N_A * beta / N_V = %f - %d * %d / %d = %f (%f)",
- Int_t(nA), Int_t(nC), Int_t(nE), Int_t(beta),
- nAccepted / vtxEff, Int_t(nAccepted), Int_t(beta),
- Int_t(nWithVertex), ntotal, scaler));
- }
- }
- if (scheme & kZeroBin) {
+
+ TString rhs("N = N_acc");
+ if (!(scheme & kZeroBin)) {
+ if (scheme & kEventLevel) {
+ ntotal = nAccepted / vtxEff;
+ scaler = vtxEff;
+ DMSG(fDebug,0,"Calculating event normalisation as\n"
+ " N = N_A * N_T / N_V = %d * %d / %d = %f (%f)",
+ Int_t(nAccepted), Int_t(nTriggered), Int_t(nWithVertex),
+ ntotal, scaler);
+ if (scheme & kBackground) {
+ // 1 E_V E_V
+ // s = --------- = ------------- = ------------
+ // 1 - beta 1 - beta E_V 1 - beta N_V
+ // --- ---- -------- ---- ---
+ // E_V N_V N_V N_V N_T
+ //
+ // E_V
+ // = ------------
+ // 1 - beta
+ // ----
+ // N_T
+ //
+ ntotal -= nAccepted * beta / nWithVertex;
+ // This one is direct and correct.
+ // scaler = 1. / (1. / vtxEff - beta / nWithVertex);
+ // A simpler expresion
+ scaler /= (1 - beta / nTriggered); // 0.831631 -> 0.780689
+ DMSG(fDebug,0,"Calculating event normalisation as\n"
+ " beta = N_a + N_c + 2 N_e = %d + %d - 2 * %d = %d\n"
+ " N = N - N_A * beta / N_V = %f - %d * %d / %d = %f (%f)",
+ Int_t(nA), Int_t(nC), Int_t(nE), Int_t(beta),
+ nAccepted / vtxEff, Int_t(nAccepted), Int_t(beta),
+ Int_t(nWithVertex), ntotal, scaler);
+ rhs.Append("(1/eps_V - beta/N_vtx)");
+ } // Background
+ else
+ rhs.Append("/eps_V");
+ } // Event-level
+ if (scheme & kTriggerEfficiency) {
+ Int_t old = Int_t(ntotal);
+ ntotal /= trigEff;
+ scaler *= trigEff;
+ DMSG(fDebug,0,"Correcting for trigger efficiency:\n"
+ " N = 1 / E_X * N = 1 / %f * %d = %f (%f)",
+ trigEff, old, ntotal, scaler);
+ rhs.Append("/eps_T");
+ } // Trigger efficiency
+ }
+ else {
// Calculate as
//
// N = N_A + 1/E_X * N_A / N_V (N_T - N_V - beta)
// = N_A (1 + 1/E_X (N_T/N_V - 1 - beta / N_V))
+ // = N_A (1 + 1/E_X (1/E_V - 1 - beta / N_V))
//
// s = N_A/N = 1 / (1 + 1/E_X (N_T/N_V - 1 - beta / N_V))
// = N_V / (N_V + 1/E_X (N_T - N_V - beta))
- beta / nWithVertex));
scaler = nWithVertex / (nWithVertex +
1/trigEff * (nTriggered-nWithVertex-beta));
- AliInfo(Form("Calculating event normalisation as\n"
- " beta = N_a + N_c + 2 N_e = %d + %d - 2 * %d = %d\n"
- " N = N_A (1 + 1/E_X (N_T/N_V - 1 - beta / N_V)) = "
- "%d (1 + 1 / %f (%d / %d - 1 - %d / %d)) = %f (%f)",
- Int_t(nA), Int_t(nC), Int_t(nE), Int_t(beta),
- Int_t(nAccepted), trigEff, Int_t(nTriggered),
- Int_t(nWithVertex), Int_t(beta), Int_t(nWithVertex),
- ntotal, scaler));
- }
- if (scheme & kTriggerEfficiency && !(scheme & kZeroBin)) {
- ntotal /= trigEff;
- scaler *= trigEff;
- AliInfo(Form("Correcting for trigger efficiency:\n"
- " N = 1 / E_X * N = 1 / %f * %d = %f (%f)",
- trigEff, Int_t(ntotal), ntotal / trigEff, scaler));
+ DMSG(fDebug,0,"Calculating event normalisation as\n"
+ " beta = N_a + N_c + 2 N_e = %d + %d - 2 * %d = %d\n"
+ " N = N_A (1 + 1/E_X (N_T/N_V - 1 - beta / N_V)) = "
+ "%d (1 + 1 / %f (%d / %d - 1 - %d / %d)) = %f (%f)",
+ Int_t(nA), Int_t(nC), Int_t(nE), Int_t(beta),
+ Int_t(nAccepted), trigEff, Int_t(nTriggered),
+ Int_t(nWithVertex), Int_t(beta), Int_t(nWithVertex),
+ ntotal, scaler);
+ rhs.Append("(1+1/eps_T(1/eps_V-1-beta/N_vtx))");
}
- AliInfo(Form("\n"
- " Total of %9d events for %s\n"
- " of these %9d have an offline trigger\n"
- " of these N_T = %9d has the selected trigger\n"
- " of these N_V = %9d has a vertex\n"
- " of these N_A = %9d were in the selected range\n"
- " Triggers by hardware type:\n"
- " N_b = %9d\n"
- " N_ac = %9d (%d+%d)\n"
+ if (text) {
+ text->Append(Form("%-40s = %d\n", "N_all", UInt_t(nAll)));
+ text->Append(Form("%-40s = %d\n", "N_acc", UInt_t(nAccepted)));
+ text->Append(Form("%-40s = %d\n", "N_trg", UInt_t(nTriggered)));
+ text->Append(Form("%-40s = %d\n", "N_vtx", UInt_t(nWithVertex)));
+ text->Append(Form("%-40s = %d\n", "N_B", UInt_t(nB)));
+ text->Append(Form("%-40s = %d\n", "N_A", UInt_t(nA)));
+ text->Append(Form("%-40s = %d\n", "N_C", UInt_t(nC)));
+ text->Append(Form("%-40s = %d\n", "N_E", UInt_t(nE)));
+ text->Append(Form("%-40s = %d\n", "beta = N_A + N_C - 2N_E",UInt_t(beta)));
+ text->Append(Form("%-40s = %f\n", "eps_V = N_vtx/N_trg",vtxEff));
+ text->Append(Form("%-40s = %f\n", "eps_T", trigEff));
+ text->Append(Form("%-40s = %f\n", rhs.Data(), ntotal));
+ }
+ TString tN = t.GetXaxis()->GetBinLabel(AliAODForwardMult::kWithTrigger);
+ tN.ReplaceAll("w/Selected trigger (","");
+ tN.ReplaceAll(")", "");
+ DMSG(fDebug,0,"\n"
+ " Total of %9d events for %s\n"
+ " of these %9d have an offline trigger\n"
+ " of these N_T = %9d has the selected trigger (%s)\n"
+ " of these N_V = %9d has a vertex\n"
+ " of these N_A = %9d were in the selected range\n"
+ " Triggers by hardware type:\n"
+ " N_b = %9d\n"
+ " N_ac = %9d (%d+%d)\n"
" N_e = %9d\n"
- " Vertex efficiency: %f\n"
- " Trigger efficiency: %f\n"
- " Total number of events: N = %f\n"
- " Scaler (N_A/N): %f",
- Int_t(nAll), GetTitle(), Int_t(nOffline),
- Int_t(nTriggered), Int_t(nWithVertex), Int_t(nAccepted),
- Int_t(nB), Int_t(nA+nC), Int_t(nA), Int_t(nC), Int_t(nE),
- vtxEff, trigEff, ntotal, scaler));
+ " Vertex efficiency: %f\n"
+ " Trigger efficiency: %f\n"
+ " Total number of events: N = %f\n"
+ " Scaler (N_A/N): %f\n"
+ " %25s = %f",
+ Int_t(nAll), GetTitle(), Int_t(nOffline),
+ Int_t(nTriggered), tN.Data(),
+ Int_t(nWithVertex), Int_t(nAccepted),
+ Int_t(nB), Int_t(nA+nC), Int_t(nA), Int_t(nC), Int_t(nE),
+ vtxEff, trigEff, ntotal, scaler, rhs.Data(), ntotal);
return scaler;
}
const char* postfix) const
{
static TString n;
- n = Form("dndeta%s%s",GetName(), postfix);
+ n = GetName();
+ n.ReplaceAll("dNdeta", "");
+ n.Prepend("dndeta");
+ n.Append(postfix);
+ // n = Form("dndeta%s%s",GetName(), postfix);
if (rebin > 1) n.Append(Form("_rebin%02d", rebin));
if (sym) n.Append("_mirror");
return n.Data();
}
//________________________________________________________________________
TH1*
-AliBasedNdetaTask::CentralityBin::GetResult(Int_t rebin,
- Bool_t sym,
- const char* postfix) const
+AliBasedNdetaTask::CentralityBin::GetResult(Int_t rebin,
+ Bool_t sym,
+ const char* postfix,
+ Bool_t verbose) const
{
if (!fOutput) {
- AliWarning(Form("No output list defined in %s [%3d,%3d]", GetName(),
- fLow, fHigh));
+ AliWarningF("No output list defined in %s [%3d,%3d]", GetName(),
+ fLow, fHigh);
return 0;
}
TString n = GetResultName(rebin, sym, postfix);
TObject* o = fOutput->FindObject(n.Data());
if (!o) {
- // AliWarning(Form("Object %s not found in output list", n.Data()));
+ if (verbose)
+ AliWarningF("Object %s not found in output list of %s",
+ n.Data(), GetName());
return 0;
}
return static_cast<TH1*>(o);
const char* postfix,
bool rootProj,
bool corrEmpty,
- const TH1* shapeCorr,
+ const TH2F* shapeCorr,
Double_t scaler,
bool symmetrice,
Int_t rebin,
// cutEdges Whether to cut edges when rebinning
//
DGUARD(fDebug,1,"Make centrality bin result from %s", sum->GetName());
- TH2D* copy = static_cast<TH2D*>(sum->Clone(Form("d2Ndetadphi%s%s",
- GetName(), postfix)));
- TH1D* accNorm = ProjectX(sum, Form("norm%s%s",GetName(), postfix), 0, 0,
- rootProj, corrEmpty, false);
+ TString base(GetName());
+ base.ReplaceAll("dNdeta", "");
+ base.Append(postfix);
+ TH2D* copy = static_cast<TH2D*>(sum->Clone(Form("d2Ndetadphi%s",
+ base.Data())));
+
+ TH1D* accNorm = 0;
+ Int_t nY = sum->GetNbinsY();
+ // Hack HHD Hans test code to manually remove FMD2 dead channel (but
+ // it is on outer?)
+ //
+ // cholm comment: The original hack has been moved to
+ // AliForwarddNdetaTask::CheckEventData - this simplifies things a
+ // great deal, and we could in principle use the new phi acceptance.
+ //
+ // However, since we may have filtered out the dead sectors in the
+ // AOD production already, we can't be sure we can recalculate the
+ // phi acceptance correctly, so for now, we rely on fCorrEmpty being set.
+ Int_t o = (corrEmpty ? 0 : nY+1);
+ accNorm = ProjectX(sum, Form("norm%s",base.Data()), o, o,
+ rootProj, corrEmpty, false);
accNorm->SetDirectory(0);
// ---- Scale by shape correction ----------------------------------
if (shapeCorr) copy->Divide(shapeCorr);
- else AliInfo("No shape correction specified, or disabled");
+ else DMSG(fDebug,1,"No shape correction specified, or disabled");
// --- Normalize to the coverage -----------------------------------
- ScaleToCoverage(copy, accNorm);
-
- // --- Event-level normalization -----------------------------------
- copy->Scale(scaler);
+ if (corrEmpty) {
+ ScaleToCoverage(copy, accNorm);
+ // --- Event-level normalization ---------------------------------
+ copy->Scale(scaler);
+ }
// --- Project on X axis -------------------------------------------
- TH1D* dndeta = ProjectX(copy, Form("dndeta%s%s",GetName(), postfix),
- 1, sum->GetNbinsY(), rootProj, corrEmpty);
+ TH1D* dndeta = ProjectX(copy, Form("dndeta%s",base.Data()),
+ 1, nY, rootProj, corrEmpty);
dndeta->SetDirectory(0);
// Event-level normalization
+ if (!corrEmpty) {
+ ScaleToCoverage(dndeta, accNorm);
+ dndeta->Scale(scaler);
+ }
dndeta->Scale(1., "width");
copy->Scale(1., "width");
- TH1D* dndetaMCCorrection = 0;
- TList* centlist = 0;
- TH1D* dndetaMCtruth = 0;
- TList* truthcentlist = 0;
+ TH1D* dndetaMCCorrection = 0;
+ TH1D* dndetaMCtruth = 0;
+ TList* centlist = 0;
+ TList* truthcentlist = 0;
- // Possible final correction to <MC analysis> / <MC truth>
- if(mclist)
+ // --- Possible final correction to <MC analysis> / <MC truth> -----
+ // we get the rebinned distribution for satellite to make the correction
+ TString rebinSuf(fSatelliteVertices ? "_rebin05" : "");
+ if(mclist) {
centlist = static_cast<TList*> (mclist->FindObject(GetListName()));
- if(centlist)
- dndetaMCCorrection = static_cast<TH1D*>(centlist->FindObject(Form("dndeta%s%s",GetName(), postfix)));
- if(truthlist)
- truthcentlist = static_cast<TList*> (truthlist->FindObject(GetListName()));
- if(truthcentlist)
- dndetaMCtruth = static_cast<TH1D*> (truthcentlist->FindObject("dndetaTruth"));
- //std::cout<<dndetaMCCorrection<<" "<<dndetaMCtruth<<std::endl;
- if(dndetaMCCorrection && dndetaMCtruth) {
+ if(centlist)
+ dndetaMCCorrection =
+ static_cast<TH1D*>(centlist->FindObject(Form("dndeta%s%s",
+ base.Data(),
+ rebinSuf.Data())));
+ }
+ if (truthlist) {
+ truthcentlist = static_cast<TList*>(truthlist->FindObject(GetListName()));
+ if (truthcentlist)
+ // TODO here new is "dndetaTruth"
+ dndetaMCtruth =
+ static_cast<TH1D*>(truthcentlist->FindObject(Form("dndetaMCTruth%s",
+ rebinSuf.Data())));
+ }
+
+ if (dndetaMCCorrection && dndetaMCtruth) {
AliInfo("Correcting with final MC correction");
+ TString testString(dndetaMCCorrection->GetName());
+
+ // We take the measured MC dN/deta and divide with truth
dndetaMCCorrection->Divide(dndetaMCtruth);
- dndeta->Divide(dndetaMCCorrection);
-
- //std::cout<<"histo "<<Form("dndeta%s%s",GetName(), postfix)<<" "<<GetListName()<<" "<<dndetaMCCorrection<<std::endl;
- //std::cout<<"truth "<<GetListName()<<" "<<dndetaMCtruth<<std::endl;
-
+ dndetaMCCorrection->SetTitle("Final MC correction");
+ dndetaMCCorrection->SetName("finalMCCorr");
+ for(Int_t m = 1; m <= dndetaMCCorrection->GetNbinsX(); m++) {
+ if(dndetaMCCorrection->GetBinContent(m) < 0.5 ||
+ dndetaMCCorrection->GetBinContent(m) > 1.75) {
+ dndetaMCCorrection->SetBinContent(m,1.);
+ dndetaMCCorrection->SetBinError(m,0.1);
+ }
+ }
+ // Applying the correction
+ if (!fSatelliteVertices)
+ // For non-satellites we took the same binning, so we do a straight
+ // division
+ dndeta->Divide(dndetaMCCorrection);
+ else {
+ // For satelitte events, we took coarser binned histograms, so
+ // we need to do a bit more
+ for(Int_t m = 1; m <= dndeta->GetNbinsX(); m++) {
+ if(dndeta->GetBinContent(m) <= 0.01 ) continue;
+
+ Double_t eta = dndeta->GetXaxis()->GetBinCenter(m);
+ Int_t bin = dndetaMCCorrection->GetXaxis()->FindBin(eta);
+ Double_t mccorr = dndetaMCCorrection->GetBinContent(bin);
+ Double_t mcerror = dndetaMCCorrection->GetBinError(bin);
+ if (mccorr < 1e-6) {
+ dndeta->SetBinContent(m, 0);
+ dndeta->SetBinError(m, 0);
+ }
+ Double_t value = dndeta->GetBinContent(m);
+ Double_t error = dndeta->GetBinError(m);
+ Double_t sumw2 = (error * error * mccorr * mccorr +
+ mcerror * mcerror * value * value);
+ dndeta->SetBinContent(m,value/mccorr) ;
+ dndeta->SetBinError(m,TMath::Sqrt(sumw2)/mccorr/mccorr);
+ }
+ }
}
- else AliInfo("No final MC correction applied");
+ else
+ DMSG(fDebug,1,"No final MC correction applied");
// --- Set some histogram attributes -------------------------------
TString post;
GetName(), post.Data()));
// --- Make symmetric extensions and rebinnings --------------------
- if (symmetrice) fOutput->Add(Symmetrice(dndeta));
+ if (symmetrice) fOutput->Add(Symmetrice(dndeta));
fOutput->Add(dndeta);
fOutput->Add(accNorm);
fOutput->Add(copy);
fOutput->Add(Rebin(dndeta, rebin, cutEdges));
- if (symmetrice) fOutput->Add(Symmetrice(Rebin(dndeta, rebin, cutEdges)));
+ if (symmetrice) fOutput->Add(Symmetrice(Rebin(dndeta, rebin, cutEdges)));
+ if (dndetaMCCorrection) fOutput->Add(dndetaMCCorrection);
+
+ // HHD Test of dN/deta in phi bins add flag later?
+ //
+ // cholm comment: We disable this for now
+#if 0
+ for (Int_t nn=1; nn <= sum->GetNbinsY(); nn++) {
+ TH1D* dndeta_phi = ProjectX(copy, Form("dndeta%s_phibin%d",
+ base.Data(), nn),
+ nn, nn, rootProj, corrEmpty);
+ dndeta_phi->SetDirectory(0);
+ // Event-level normalization
+ dndeta_phi->Scale(TMath::Pi()/10., "width");
+
+ if(centlist)
+ dndetaMCCorrection =
+ static_cast<TH1D*>(centlist->FindObject(Form("dndeta%s_phibin%d",
+ base.Data(),nn)));
+ if(truthcentlist)
+ dndetaMCtruth
+ = static_cast<TH1D*>(truthcentlist->FindObject("dndetaMCTruth"));
+
+ if (dndetaMCCorrection && dndetaMCtruth) {
+ AliInfo("Correcting with final MC correction");
+ TString testString(dndetaMCCorrection->GetName());
+ dndetaMCCorrection->Divide(dndetaMCtruth);
+ dndetaMCCorrection->SetTitle(Form("Final_MC_correction_phibin%d",nn));
+ dndetaMCCorrection->SetName(Form("Final_MC_correction_phibin%d",nn));
+ for(Int_t m = 1; m <= dndetaMCCorrection->GetNbinsX(); m++) {
+ if(dndetaMCCorrection->GetBinContent(m) < 0.25 ||
+ dndetaMCCorrection->GetBinContent(m) > 1.75) {
+ dndetaMCCorrection->SetBinContent(m,1.);
+ dndetaMCCorrection->SetBinError(m,0.1);
+ }
+ }
+ //Applying the correction
+ dndeta_phi->Divide(dndetaMCCorrection);
+ }
+ fOutput->Add(dndeta_phi);
+ fOutput->Add(Rebin(dndeta_phi, rebin, cutEdges));
+ if(dndetaMCCorrection) fOutput->Add(dndetaMCCorrection);
+ } // End of phi
+#endif
}
//________________________________________________________________________
AliBasedNdetaTask::CentralityBin::End(TList* sums,
TList* results,
UShort_t scheme,
- const TH1* shapeCorr,
+ const TH2F* shapeCorr,
Double_t trigEff,
+ Double_t trigEff0,
Bool_t symmetrice,
Int_t rebin,
Bool_t rootProj,
results->Add(fOutput);
if (!fSum) {
- AliInfo("This task did not produce any output");
- return;
+ if (!ReadSum(fSums, false)) {
+ AliInfo("This task did not produce any output");
+ return;
+ }
}
+ if (!fSumMC) ReadSum(fSums, true);
fTriggers = static_cast<TH1I*>(fSums->FindObject("triggers"));
AliError("Couldn't find histogram 'triggers' in list");
return;
}
- if (!fSum) {
- AliError(Form("No sum object for %s", GetName()));
- return;
- }
// --- Get normalization scaler ------------------------------------
Double_t epsilonT = trigEff;
- Double_t epsilonT0 = trigEff;
- // TEMPORARY FIX
- if (triggerMask == AliAODForwardMult::kNSD) {
- // This is a local change
- epsilonT = 0.96;
- AliWarning(Form("Using hard-coded NSD trigger efficiency of %f",epsilonT));
- }
- else if (triggerMask == AliAODForwardMult::kInel) {
- // This is a local change
- epsilonT = 0.934;
- AliWarning(Form("Using hard-coded Inel trigger efficiency of %f",epsilonT));
- }
- if (scheme & kZeroBin) {
- if (triggerMask==AliAODForwardMult::kInel)
- epsilonT0 = 0.785021; // 0.100240;
- else if (triggerMask==AliAODForwardMult::kInelGt0)
- epsilonT0 = 0;
- else if (triggerMask==AliAODForwardMult::kNSD)
- epsilonT0 = .706587;
- epsilonT = 1;
- AliWarning(Form("Using hard-coded NCluster>0 trigger efficiency of %f",
- epsilonT0));
- }
+ Double_t epsilonT0 = trigEff0;
+ DMSG(fDebug,2,"Using epsilonT=%f, epsilonT0=%f for 0x%x",
+ epsilonT, epsilonT0, triggerMask);
// Get our histograms
Double_t nSum = 0;
- TH2D* sum = fSum->GetSum(fSums, fOutput, nSum, epsilonT0, 1,
- marker, rootProj, corrEmpty);
+ TH2D* sum = fSum->CalcSum(fOutput, nSum, epsilonT0, epsilonT,
+ marker, rootProj, corrEmpty);
Double_t nSumMC = 0;
TH2D* sumMC = 0;
- if (fSumMC) sumMC = fSumMC->GetSum(fSums, fOutput, nSumMC,
- epsilonT0, 1, marker,
- rootProj, corrEmpty);
+ if (fSumMC) sumMC = fSumMC->CalcSum(fOutput, nSumMC,
+ epsilonT0, epsilonT, marker,
+ rootProj, corrEmpty);
if (!sum) {
AliError("Failed to get sum from summer - bailing out");
return;
}
+ TString text;
Double_t ntotal = nSum;
- Double_t scaler = Normalization(*fTriggers, scheme, epsilonT, ntotal);
+ Double_t scaler = Normalization(*fTriggers, scheme, epsilonT, ntotal, &text);
if (scaler < 0) {
AliError("Failed to calculate normalization - bailing out");
return;
}
fOutput->Add(fTriggers->Clone());
+ fOutput->Add(new TNamed("normCalc", text.Data()));
// --- Make result and store ---------------------------------------
MakeResult(sum, "", rootProj, corrEmpty, (scheme & kShape) ? shapeCorr : 0,
// if (!IsAllBin()) return;
}
+//____________________________________________________________________
+Bool_t
+AliBasedNdetaTask::ApplyEmpiricalCorrection(const AliAODForwardMult* aod,
+ TH2D* data)
+{
+ if (!fglobalempiricalcorrection || !data)
+ return true;
+ Float_t zvertex=aod->GetIpZ();
+ Int_t binzvertex=fglobalempiricalcorrection->GetXaxis()->FindBin(zvertex);
+ if(binzvertex<1||binzvertex>fglobalempiricalcorrection->GetNbinsX())
+ return false;
+ for (int i=1;i<=data->GetNbinsX();i++) {
+ Int_t bincorrection=fglobalempiricalcorrection->GetYaxis()
+ ->FindBin(data->GetXaxis()->GetBinCenter(i));
+ if(bincorrection<1||bincorrection>fglobalempiricalcorrection->GetNbinsY())
+ return false;
+ Float_t correction=fglobalempiricalcorrection
+ ->GetBinContent(binzvertex,bincorrection);
+ if(correction<0.001) {
+ data->SetBinContent(i,0,0);
+ data->SetBinContent(i,data->GetNbinsY()+1,0);
+ }
+ for(int j=1;j<=data->GetNbinsY();j++) {
+ if (data->GetBinContent(i,j)>0.0) {
+ data->SetBinContent(i,j,data->GetBinContent(i,j)*correction);
+ data->SetBinError(i,j,data->GetBinError(i,j)*correction);
+ }
+ }
+ }
+ return true;
+}
//
// EOF
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff