[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / AliForwardMCCorrectionsTask.cxx

// 
// Calculate the corrections in the forward regions
// 
// Inputs: 
//   - AliESDEvent 
//
// Outputs: 
//   - AliAODForwardMult 
// 
// Histograms 
//   
// Corrections used 
// 
#include "AliForwardMCCorrectionsTask.h"
#include "AliTriggerAnalysis.h"
#include "AliPhysicsSelection.h"
#include "AliLog.h"
#include "AliHeader.h"
#include "AliGenEventHeader.h"
#include "AliESDEvent.h"
#include "AliAODHandler.h"
#include "AliMultiplicity.h"
#include "AliInputEventHandler.h"
#include "AliStack.h"
#include "AliMCEvent.h"
#include "AliAODForwardMult.h"
#include "AliFMDStripIndex.h"
#include "AliFMDCorrSecondaryMap.h"
#include <TH1.h>
#include <TH2D.h>
#include <TDirectory.h>
#include <TTree.h>
#include <TList.h>
#include <TROOT.h>
#include <iostream>

//====================================================================
namespace {
  const char* GetEventName(Bool_t tr, Bool_t vtx) 
  {
    return Form("nEvents%s%s", (tr ? "Tr" : ""), (vtx ? "Vtx" : ""));
  }
#if 0
  const char* GetHitsName(UShort_t d, Char_t r) 
  {
    return Form("hitsFMD%d%c", d, r);
  }
  const char* GetStripsName(UShort_t d, Char_t r)
  {
    return Form("stripsFMD%d%c", d, r);
  }
  const char* GetPrimaryName(Char_t r, Bool_t trVtx)
  {
    return Form("primaries%s%s", 
		((r == 'I' || r == 'i') ? "Inner" : "Outer"), 
		(trVtx ? "TrVtx" : "All"));
  }
#endif
}

//====================================================================
AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask()
  : AliAnalysisTaskSE(),
    fInspector(),
    fTrackDensity(),
    fESDFMD(),
    fVtxBins(0),
    fFirstEvent(true),
    fHEvents(0), 
    fHEventsTr(0), 
    fHEventsTrVtx(0),
    fVtxAxis(),
    fEtaAxis(),
    fList()
{
  // 
  // Constructor 
  // 
  // Parameters:
  //    name Name of task 
  //
}

//____________________________________________________________________
AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const char* name)
  : AliAnalysisTaskSE(name),
    fInspector("eventInspector"), 
    fTrackDensity("trackDensity"),
    fESDFMD(),
    fVtxBins(0),
    fFirstEvent(true),
    fHEvents(0), 
    fHEventsTr(0), 
    fHEventsTrVtx(0),
    fVtxAxis(10,-10,10), 
    fEtaAxis(200,-4,6),
    fList()
{
  // 
  // Constructor 
  // 
  // Parameters:
  //    name Name of task 
  //
  DefineOutput(1, TList::Class());
  DefineOutput(2, TList::Class());
}

//____________________________________________________________________
AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const AliForwardMCCorrectionsTask& o)
  : AliAnalysisTaskSE(o),
    fInspector(o.fInspector),
    fTrackDensity(),
    fESDFMD(o.fESDFMD),
    fVtxBins(0),
    fFirstEvent(o.fFirstEvent),
    fHEvents(o.fHEvents), 
    fHEventsTr(o.fHEventsTr), 
    fHEventsTrVtx(o.fHEventsTrVtx),
    fVtxAxis(10,-10,10), 
    fEtaAxis(200,-4,6),
    fList(o.fList)
{
  // 
  // Copy constructor 
  // 
  // Parameters:
  //    o Object to copy from 
  //
}

//____________________________________________________________________
AliForwardMCCorrectionsTask&
AliForwardMCCorrectionsTask::operator=(const AliForwardMCCorrectionsTask& o)
{
  // 
  // Assignment operator 
  // 
  // Parameters:
  //    o Object to assign from 
  // 
  // Return:
  //    Reference to this object 
  //
  if (&o == this) return *this;
  fInspector         = o.fInspector;
  fTrackDensity      = o.fTrackDensity;
  fESDFMD            = o.fESDFMD;
  fVtxBins           = o.fVtxBins;
  fFirstEvent        = o.fFirstEvent;
  fHEvents           = o.fHEvents;
  fHEventsTr         = o.fHEventsTr;
  fHEventsTrVtx      = o.fHEventsTrVtx;
  SetVertexAxis(o.fVtxAxis);
  SetEtaAxis(o.fEtaAxis);

  return *this;
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::Init()
{
  // 
  // Initialize the task 
  // 
  //
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::SetVertexAxis(Int_t nBin, Double_t min, 
					   Double_t max)
{
  // 
  // Set the vertex axis to use
  // 
  // Parameters:
  //    nBins Number of bins
  //    vzMin Least @f$z@f$ coordinate of interation point
  //    vzMax Largest @f$z@f$ coordinate of interation point
  //
  if (max < min) { 
    Double_t tmp = min;
    min          = max;
    max          = tmp;
  }
/*
  if (min < -10) 
    AliWarning(Form("Minimum vertex %f < -10, make sure you want this",min));
  if (max > +10) 
    AliWarning(Form("Minimum vertex %f > +10, make sure you want this",max));
*/
  fVtxAxis.Set(nBin, min, max);
}
//____________________________________________________________________
void
AliForwardMCCorrectionsTask::SetVertexAxis(const TAxis& axis)
{
  // 
  // Set the vertex axis to use
  // 
  // Parameters:
  //    axis Axis
  //
  SetVertexAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::SetEtaAxis(Int_t nBin, Double_t min, Double_t max)
{
  // 
  // Set the eta axis to use
  // 
  // Parameters:
  //    nBins Number of bins
  //    vzMin Least @f$\eta@f$ 
  //    vzMax Largest @f$\eta@f$ 
  //
  if (max < min) { 
    Double_t tmp = min;
    min          = max;
    max          = tmp;
  }
  if (min < -4) 
    AliWarning(Form("Minimum eta %f < -4, make sure you want this",min));
  if (max > +6) 
    AliWarning(Form("Minimum eta %f > +6, make sure you want this",max));
  fEtaAxis.Set(nBin, min, max);
}
//____________________________________________________________________
void
AliForwardMCCorrectionsTask::SetEtaAxis(const TAxis& axis)
{
  // 
  // Set the eta axis to use
  // 
  // Parameters:
  //    axis Axis
  //
  SetEtaAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::DefineBins(TList* l)
{
  if (!fVtxBins) fVtxBins = new TObjArray(fVtxAxis.GetNbins(), 1);
  if (fVtxBins->GetEntries() > 0) return;

  fVtxBins->SetOwner();
  for (Int_t  i    = 1; i <= fVtxAxis.GetNbins(); i++) { 
    Double_t  low  = fVtxAxis.GetBinLowEdge(i);
    Double_t  high = fVtxAxis.GetBinUpEdge(i);
    VtxBin*   bin  = new VtxBin(low, high, fEtaAxis);
    fVtxBins->AddAt(bin, i);
    bin->DefineOutput(l);
  }
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::UserCreateOutputObjects()
{
  // 
  // Create output objects 
  // 
  //
  fList = new TList;
  fList->SetOwner();
  fList->SetName(Form("%sSums", GetName()));

  DefineBins(fList);

  fHEvents = new TH1I(GetEventName(false,false),
		      "Number of all events", 
		      fVtxAxis.GetNbins(), 
		      fVtxAxis.GetXmin(), 
		      fVtxAxis.GetXmax());
  fHEvents->SetXTitle("v_{z} [cm]");
  fHEvents->SetYTitle("# of events");
  fHEvents->SetFillColor(kBlue+1);
  fHEvents->SetFillStyle(3001);
  fHEvents->SetDirectory(0);
  fList->Add(fHEvents);

  fHEventsTr = new TH1I(GetEventName(true, false), 
			"Number of triggered events",
			fVtxAxis.GetNbins(), 
			fVtxAxis.GetXmin(), 
			fVtxAxis.GetXmax());
  fHEventsTr->SetXTitle("v_{z} [cm]");
  fHEventsTr->SetYTitle("# of events");
  fHEventsTr->SetFillColor(kRed+1);
  fHEventsTr->SetFillStyle(3001);
  fHEventsTr->SetDirectory(0);
  fList->Add(fHEventsTr);

  fHEventsTrVtx = new TH1I(GetEventName(true,true),
			   "Number of events w/trigger and vertex", 
			   fVtxAxis.GetNbins(), 
			   fVtxAxis.GetXmin(), 
			   fVtxAxis.GetXmax());
  fHEventsTrVtx->SetXTitle("v_{z} [cm]");
  fHEventsTrVtx->SetYTitle("# of events");
  fHEventsTrVtx->SetFillColor(kBlue+1);
  fHEventsTrVtx->SetFillStyle(3001);
  fHEventsTrVtx->SetDirectory(0);
  fList->Add(fHEventsTrVtx);

  // Copy axis objects to output 
  TH1* vtxAxis = new TH1D("vtxAxis", "Vertex axis", 
			  fVtxAxis.GetNbins(), 
			  fVtxAxis.GetXmin(), 
			  fVtxAxis.GetXmax());
  TH1* etaAxis = new TH1D("etaAxis", "Eta axis", 
			  fEtaAxis.GetNbins(), 
			  fEtaAxis.GetXmin(), 
			  fEtaAxis.GetXmax());
  fList->Add(vtxAxis);
  fList->Add(etaAxis);

  AliInfo(Form("Initialising sub-routines: %p, %p", 
	       &fInspector, &fTrackDensity));
  fInspector.DefineOutput(fList);
  fInspector.Init(fVtxAxis);
  fTrackDensity.DefineOutput(fList);

  PostData(1, fList);
}
//____________________________________________________________________
void
AliForwardMCCorrectionsTask::UserExec(Option_t*)
{
  // 
  // Process each event 
  // 
  // Parameters:
  //    option Not used
  //  

  // Get the input data - MC event
  AliMCEvent*  mcEvent = MCEvent();
  if (!mcEvent) { 
    AliWarning("No MC event found");
    return;
  }

  // Get the input data - ESD event
  AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
  if (!esd) { 
    AliWarning("No ESD event found for input event");
    return;
  }

  //--- Read run information -----------------------------------------
  if (fFirstEvent && esd->GetESDRun()) {
    fInspector.ReadRunDetails(esd);
    
    AliInfo(Form("Initializing with parameters from the ESD:\n"
		 "         AliESDEvent::GetBeamEnergy()   ->%f\n"
		 "         AliESDEvent::GetBeamType()     ->%s\n"
		 "         AliESDEvent::GetCurrentL3()    ->%f\n"
		 "         AliESDEvent::GetMagneticField()->%f\n"
		 "         AliESDEvent::GetRunNumber()    ->%d\n",
		 esd->GetBeamEnergy(),
		 esd->GetBeamType(),
		 esd->GetCurrentL3(),
		 esd->GetMagneticField(),
		 esd->GetRunNumber()));

    Print();
    fFirstEvent = false;
  }

  // Some variables 
  UInt_t   triggers; // Trigger bits
  Bool_t   lowFlux;  // Low flux flag
  UShort_t iVz;      // Vertex bin from ESD
  Double_t vZ;       // Z coordinate from ESD
  Double_t cent;     // Centrality 
  UShort_t iVzMc;    // Vertex bin from MC
  Double_t vZMc;     // Z coordinate of IP vertex from MC
  Double_t b;        // Impact parameter
  Int_t    nPart;    // Number of participants 
  Int_t    nBin;     // Number of binary collisions 
  Double_t phiR;     // Reaction plane from MC
  UShort_t nClusters;// Number of SPD clusters 
  // Process the data 
  UInt_t retESD = fInspector.Process(esd, triggers, lowFlux, iVz, vZ, 
				     cent, nClusters);
  fInspector.ProcessMC(mcEvent, triggers, iVzMc, vZMc, 
		       b, nPart, nBin, phiR);

  Bool_t isInel   = triggers & AliAODForwardMult::kInel;
  Bool_t hasVtx   = retESD == AliFMDMCEventInspector::kOk;

  // Fill the event count histograms 
  if (isInel)           fHEventsTr->Fill(vZMc);
  if (isInel && hasVtx) fHEventsTrVtx->Fill(vZMc);
  fHEvents->Fill(vZMc);

  // Now find our vertex bin object 
  VtxBin* bin = 0;
  if (iVzMc > 0 && iVzMc <= fVtxAxis.GetNbins()) 
    bin = static_cast<VtxBin*>(fVtxBins->At(iVzMc));
  if (!bin) { 
    // AliError(Form("No vertex bin object @ %d (%f)", iVzMc, vZMc));
    return;
  }

  // Clear our ESD object 
  fESDFMD.Clear();

  // Get FMD data 
  AliESDFMD* esdFMD = esd->GetFMDData();
  
  // Now process our input data and store in own ESD object 
  fTrackDensity.Calculate(*esdFMD, *mcEvent, vZMc, fESDFMD, bin->fPrimary);
  bin->fCounts->Fill(0.5);

  // And then bin the data in our vtxbin 
  for (UShort_t d=1; d<=3; d++) { 
    UShort_t nr = (d == 1 ? 1 : 2);
    for (UShort_t q=0; q<nr; q++) { 
      Char_t      r = (q == 0 ? 'I' : 'O');
      UShort_t    ns= (q == 0 ?  20 :  40);
      UShort_t    nt= (q == 0 ? 512 : 256);
      TH2D*       h = bin->fHists.Get(d,r);

      for (UShort_t s=0; s<ns; s++) { 
	for (UShort_t t=0; t<nt; t++) {
	  Float_t mult = fESDFMD.Multiplicity(d,r,s,t);
	  
	  if (mult == 0 || mult > 20) continue;

	  Float_t phi = fESDFMD.Phi(d,r,s,t) / 180 * TMath::Pi();
	  Float_t eta = fESDFMD.Eta(d,r,s,t);
	  h->Fill(eta,phi,mult);
	} // for t
      } // for s 
    } // for q 
  } // for d

  PostData(1, fList);
}
//____________________________________________________________________
void
AliForwardMCCorrectionsTask::Terminate(Option_t*)
{
  // 
  // End of job
  // 
  // Parameters:
  //    option Not used 
  //

  fList = dynamic_cast<TList*>(GetOutputData(1));
  if (!fList) {
    AliError("No output list defined");
    return;
  }

  DefineBins(fList);

  // Output list 
  TList* output = new TList;
  output->SetOwner();
  output->SetName(Form("%sResults", GetName()));

  // --- Fill correction object --------------------------------------
  AliFMDCorrSecondaryMap* corr = new AliFMDCorrSecondaryMap;
  corr->SetVertexAxis(fVtxAxis);
  corr->SetEtaAxis(fEtaAxis);
  
  TIter     next(fVtxBins);
  VtxBin*   bin = 0;
  UShort_t  iVz = 1;
  while ((bin = static_cast<VtxBin*>(next()))) 
    bin->Finish(fList, output, iVz++, corr);

  output->Add(corr);

  PostData(2, output);
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::Print(Option_t* option) const
{
  std::cout << ClassName() << "\n"
	    << "  Vertex bins:      " << fVtxAxis.GetNbins() << '\n'
	    << "  Vertex range:     [" << fVtxAxis.GetXmin() 
	    << "," << fVtxAxis.GetXmax() << "]\n"
	    << "  Eta bins:         " << fEtaAxis.GetNbins() << '\n'
	    << "  Eta range:        [" << fEtaAxis.GetXmin() 
	    << "," << fEtaAxis.GetXmax() << "]"
	    << std::endl;
  gROOT->IncreaseDirLevel();
  fInspector.Print(option);
  fTrackDensity.Print(option);
  gROOT->DecreaseDirLevel();
}

//====================================================================
const char*
AliForwardMCCorrectionsTask::VtxBin::BinName(Double_t low, 
					     Double_t high) 
{
  static TString buf;
  buf = Form("vtx%+05.1f_%+05.1f", low, high);
  buf.ReplaceAll("+", "p");
  buf.ReplaceAll("-", "m");
  buf.ReplaceAll(".", "d");
  return buf.Data();
}


//____________________________________________________________________
AliForwardMCCorrectionsTask::VtxBin::VtxBin()
  : fHists(), 
    fPrimary(0),
    fCounts(0)
{
}
//____________________________________________________________________
AliForwardMCCorrectionsTask::VtxBin::VtxBin(Double_t low, 
					    Double_t high, 
					    const TAxis& axis)
  : TNamed(BinName(low, high), 
	   Form("%+5.1fcm<v_{z}<%+5.1fcm", low, high)),
    fHists(), 
    fPrimary(0),
    fCounts(0)
{
  fHists.Init(axis);

  fPrimary = new TH2D("primary", "Primaries", 
		      axis.GetNbins(), axis.GetXmin(), axis.GetXmax(), 
		      40, 0, 2*TMath::Pi());
  fPrimary->SetXTitle("#eta");
  fPrimary->SetYTitle("#varphi [radians]");
  fPrimary->Sumw2();
  fPrimary->SetDirectory(0);

  fCounts = new TH1D("counts", "Counts", 1, 0, 1);
  fCounts->SetXTitle("Events");
  fCounts->SetYTitle("# of Events");
  fCounts->SetDirectory(0);
}

//____________________________________________________________________
AliForwardMCCorrectionsTask::VtxBin::VtxBin(const VtxBin& o)
  : TNamed(o),
    fHists(o.fHists),
    fPrimary(0), 
    fCounts(0)
{
  if (o.fPrimary) {
    fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
    fPrimary->SetDirectory(0);
  }
  if (o.fCounts) {
    fCounts = static_cast<TH1D*>(o.fCounts->Clone());
    fCounts->SetDirectory(0);
  }
}

//____________________________________________________________________
AliForwardMCCorrectionsTask::VtxBin&
AliForwardMCCorrectionsTask::VtxBin::operator=(const VtxBin& o)
{
  if (&o == this) return *this;
  TNamed::operator=(o);
  fHists   = o.fHists;
  fPrimary = 0;
  fCounts  = 0;
  if (o.fPrimary) {
    fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
    fPrimary->SetDirectory(0);
  }
  if (o.fCounts) {
    fCounts = static_cast<TH1D*>(o.fCounts->Clone());
    fCounts->SetDirectory(0);
  }
  return *this;
}

//____________________________________________________________________
void
AliForwardMCCorrectionsTask::VtxBin::DefineOutput(TList* l)
{
  TList* d = new TList;
  d->SetName(GetName());
  d->SetOwner();
  l->Add(d);

  d->Add(fHists.fFMD1i);
  d->Add(fHists.fFMD2i);
  d->Add(fHists.fFMD2o);
  d->Add(fHists.fFMD3i);
  d->Add(fHists.fFMD3o);

  d->Add(fPrimary);
  d->Add(fCounts);
}

//____________________________________________________________________
TH2D*
AliForwardMCCorrectionsTask::VtxBin::MakeBg(const TH2D* hits, 
					    const TH2D* primary) const
{
  TH2D* h = static_cast<TH2D*>(hits->Clone());
  h->SetDirectory(0);
  TString n(h->GetName());
  n.ReplaceAll("_cache", "");
  h->SetName(n);
  h->Divide(primary);
  
  return h;
}
  
//____________________________________________________________________
void
AliForwardMCCorrectionsTask::VtxBin::Finish(const TList* input, 
					    TList* output, 
					    UShort_t iVz,
					    AliFMDCorrSecondaryMap* map)
{
  TList* out = new TList;
  out->SetName(GetName());
  out->SetOwner();
  output->Add(out);

  TList* l = static_cast<TList*>(input->FindObject(GetName()));
  if (!l) { 
    AliError(Form("List %s not found in %s", GetName(), input->GetName()));
    return;
  }

  TH2D*   fmd1i = static_cast<TH2D*>(l->FindObject("FMD1I_cache"));
  TH2D*   fmd2i = static_cast<TH2D*>(l->FindObject("FMD2I_cache"));
  TH2D*   fmd2o = static_cast<TH2D*>(l->FindObject("FMD2O_cache"));
  TH2D*   fmd3i = static_cast<TH2D*>(l->FindObject("FMD3I_cache"));
  TH2D*   fmd3o = static_cast<TH2D*>(l->FindObject("FMD3O_cache"));
  TH2D*   primO = static_cast<TH2D*>(l->FindObject("primary"));
  if (!fmd1i || !fmd2i || !fmd2o || !fmd3i || !fmd3o || !primO) {
    AliError(Form("Missing histogram(s): %p,%p,%p,%p,%p,%p",
		  fmd1i, fmd2i, fmd2o, fmd3i, fmd3o, primO));
    return;
  }

  // Half coverage in phi for inners
  TH2D*   primI = static_cast<TH2D*>(primO->Clone());
  primI->SetDirectory(0);
  primI->RebinY(2); 

  TH2D* bg1i = MakeBg(fmd1i, primI);
  TH2D* bg2i = MakeBg(fmd2i, primI);
  TH2D* bg2o = MakeBg(fmd2o, primO);
  TH2D* bg3i = MakeBg(fmd3i, primI);
  TH2D* bg3o = MakeBg(fmd3o, primO);
  map->SetCorrection(1, 'I', iVz, bg1i);
  map->SetCorrection(2, 'I', iVz, bg2i);
  map->SetCorrection(2, 'O', iVz, bg2o);
  map->SetCorrection(3, 'I', iVz, bg3i);
  map->SetCorrection(3, 'O', iVz, bg3o);
  out->Add(bg1i);
  out->Add(bg2i);
  out->Add(bg2o);
  out->Add(bg3i);
  out->Add(bg3o);
 
}

//
// EOF
//
Commit	Line	Data
7984e5f7	1	//
	2	// Calculate the corrections in the forward regions
	3	//
	4	// Inputs:
	5	// - AliESDEvent
	6	//
	7	// Outputs:
	8	// - AliAODForwardMult
	9	//
	10	// Histograms
	11	//
	12	// Corrections used
	13	//
0bd4b00f	14	#include "AliForwardMCCorrectionsTask.h"
563a673f	15	#include "AliTriggerAnalysis.h"
	16	#include "AliPhysicsSelection.h"
	17	#include "AliLog.h"
0bd4b00f	18	#include "AliHeader.h"
0bd4b00f	19	#include "AliGenEventHeader.h"
563a673f	20	#include "AliESDEvent.h"
	21	#include "AliAODHandler.h"
	22	#include "AliMultiplicity.h"
	23	#include "AliInputEventHandler.h"
0bd4b00f	24	#include "AliStack.h"
0bd4b00f	25	#include "AliMCEvent.h"
5bb5d1f6	26	#include "AliAODForwardMult.h"
5bb5d1f6	27	#include "AliFMDStripIndex.h"
f91358e1	28	#include "AliFMDCorrSecondaryMap.h"
563a673f	29	#include <TH1.h>
fb3430ac	30	#include <TH2D.h>
563a673f	31	#include <TDirectory.h>
563a673f	32	#include <TTree.h>
fb3430ac	33	#include <TList.h>
f91358e1	34	#include <TROOT.h>
f91358e1	35	#include <iostream>
563a673f	36
	37	//====================================================================
	38	namespace {
	39	const char* GetEventName(Bool_t tr, Bool_t vtx)
	40	{
0bd4b00f	41	return Form("nEvents%s%s", (tr ? "Tr" : ""), (vtx ? "Vtx" : ""));
563a673f	42	}
badf1074	43	#if 0
563a673f	44	const char* GetHitsName(UShort_t d, Char_t r)
	45	{
	46	return Form("hitsFMD%d%c", d, r);
	47	}
	48	const char* GetStripsName(UShort_t d, Char_t r)
	49	{
	50	return Form("stripsFMD%d%c", d, r);
	51	}
	52	const char* GetPrimaryName(Char_t r, Bool_t trVtx)
	53	{
	54	return Form("primaries%s%s",
	55	((r == 'I' \|\| r == 'i') ? "Inner" : "Outer"),
	56	(trVtx ? "TrVtx" : "All"));
	57	}
badf1074	58	#endif
563a673f	59	}
	60
	61	//====================================================================
0bd4b00f	62	AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask()
563a673f	63	: AliAnalysisTaskSE(),
5bb5d1f6	64	fInspector(),
f91358e1	65	fTrackDensity(),
	66	fESDFMD(),
	67	fVtxBins(0),
5bb5d1f6	68	fFirstEvent(true),
0bd4b00f	69	fHEvents(0),
563a673f	70	fHEventsTr(0),
563a673f	71	fHEventsTrVtx(0),
563a673f	72	fVtxAxis(),
0bd4b00f	73	fEtaAxis(),
0bd4b00f	74	fList()
563a673f	75	{
7984e5f7	76	//
	77	// Constructor
	78	//
	79	// Parameters:
	80	// name Name of task
	81	//
563a673f	82	}
	83
	84	//____________________________________________________________________
0bd4b00f	85	AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const char* name)
5bb5d1f6	86	: AliAnalysisTaskSE(name),
5bb5d1f6	87	fInspector("eventInspector"),
f91358e1	88	fTrackDensity("trackDensity"),
	89	fESDFMD(),
	90	fVtxBins(0),
5bb5d1f6	91	fFirstEvent(true),
0bd4b00f	92	fHEvents(0),
563a673f	93	fHEventsTr(0),
0bd4b00f	94	fHEventsTrVtx(0),
563a673f	95	fVtxAxis(10,-10,10),
0bd4b00f	96	fEtaAxis(200,-4,6),
0bd4b00f	97	fList()
563a673f	98	{
7984e5f7	99	//
	100	// Constructor
	101	//
	102	// Parameters:
	103	// name Name of task
	104	//
563a673f	105	DefineOutput(1, TList::Class());
5bb5d1f6	106	DefineOutput(2, TList::Class());
563a673f	107	}
	108
	109	//____________________________________________________________________
0bd4b00f	110	AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const AliForwardMCCorrectionsTask& o)
563a673f	111	: AliAnalysisTaskSE(o),
5bb5d1f6	112	fInspector(o.fInspector),
f91358e1	113	fTrackDensity(),
	114	fESDFMD(o.fESDFMD),
	115	fVtxBins(0),
5bb5d1f6	116	fFirstEvent(o.fFirstEvent),
0bd4b00f	117	fHEvents(o.fHEvents),
563a673f	118	fHEventsTr(o.fHEventsTr),
0bd4b00f	119	fHEventsTrVtx(o.fHEventsTrVtx),
f91358e1	120	fVtxAxis(10,-10,10),
f91358e1	121	fEtaAxis(200,-4,6),
0bd4b00f	122	fList(o.fList)
563a673f	123	{
7984e5f7	124	//
	125	// Copy constructor
	126	//
	127	// Parameters:
	128	// o Object to copy from
	129	//
563a673f	130	}
	131
	132	//____________________________________________________________________
0bd4b00f	133	AliForwardMCCorrectionsTask&
0bd4b00f	134	AliForwardMCCorrectionsTask::operator=(const AliForwardMCCorrectionsTask& o)
563a673f	135	{
7984e5f7	136	//
	137	// Assignment operator
	138	//
	139	// Parameters:
	140	// o Object to assign from
	141	//
	142	// Return:
	143	// Reference to this object
	144	//
d015ecfe	145	if (&o == this) return *this;
5bb5d1f6	146	fInspector = o.fInspector;
f91358e1	147	fTrackDensity = o.fTrackDensity;
	148	fESDFMD = o.fESDFMD;
	149	fVtxBins = o.fVtxBins;
5bb5d1f6	150	fFirstEvent = o.fFirstEvent;
f91358e1	151	fHEvents = o.fHEvents;
563a673f	152	fHEventsTr = o.fHEventsTr;
563a673f	153	fHEventsTrVtx = o.fHEventsTrVtx;
0bd4b00f	154	SetVertexAxis(o.fVtxAxis);
0bd4b00f	155	SetEtaAxis(o.fEtaAxis);
563a673f	156
	157	return *this;
	158	}
	159
	160	//____________________________________________________________________
	161	void
0bd4b00f	162	AliForwardMCCorrectionsTask::Init()
563a673f	163	{
7984e5f7	164	//
	165	// Initialize the task
	166	//
	167	//
563a673f	168	}
	169
	170	//____________________________________________________________________
	171	void
0bd4b00f	172	AliForwardMCCorrectionsTask::SetVertexAxis(Int_t nBin, Double_t min,
0bd4b00f	173	Double_t max)
563a673f	174	{
7984e5f7	175	//
	176	// Set the vertex axis to use
	177	//
	178	// Parameters:
	179	// nBins Number of bins
	180	// vzMin Least @f$z@f$ coordinate of interation point
	181	// vzMax Largest @f$z@f$ coordinate of interation point
	182	//
bcd522ff	183	if (max < min) {
563a673f	184	Double_t tmp = min;
	185	min = max;
	186	max = tmp;
	187	}
bcd522ff	188	/*
563a673f	189	if (min < -10)
	190	AliWarning(Form("Minimum vertex %f < -10, make sure you want this",min));
	191	if (max > +10)
	192	AliWarning(Form("Minimum vertex %f > +10, make sure you want this",max));
bcd522ff	193	*/
563a673f	194	fVtxAxis.Set(nBin, min, max);
	195	}
	196	//____________________________________________________________________
	197	void
0bd4b00f	198	AliForwardMCCorrectionsTask::SetVertexAxis(const TAxis& axis)
563a673f	199	{
7984e5f7	200	//
	201	// Set the vertex axis to use
	202	//
	203	// Parameters:
	204	// axis Axis
	205	//
563a673f	206	SetVertexAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
	207	}
	208
	209	//____________________________________________________________________
	210	void
0bd4b00f	211	AliForwardMCCorrectionsTask::SetEtaAxis(Int_t nBin, Double_t min, Double_t max)
563a673f	212	{
7984e5f7	213	//
	214	// Set the eta axis to use
	215	//
	216	// Parameters:
	217	// nBins Number of bins
	218	// vzMin Least @f$\eta@f$
	219	// vzMax Largest @f$\eta@f$
	220	//
bcd522ff	221	if (max < min) {
563a673f	222	Double_t tmp = min;
	223	min = max;
	224	max = tmp;
	225	}
	226	if (min < -4)
	227	AliWarning(Form("Minimum eta %f < -4, make sure you want this",min));
	228	if (max > +6)
	229	AliWarning(Form("Minimum eta %f > +6, make sure you want this",max));
bcd522ff	230	fEtaAxis.Set(nBin, min, max);
563a673f	231	}
	232	//____________________________________________________________________
	233	void
0bd4b00f	234	AliForwardMCCorrectionsTask::SetEtaAxis(const TAxis& axis)
563a673f	235	{
7984e5f7	236	//
	237	// Set the eta axis to use
	238	//
	239	// Parameters:
	240	// axis Axis
	241	//
563a673f	242	SetEtaAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
563a673f	243	}
f91358e1	244
563a673f	245	//____________________________________________________________________
f91358e1	246	void
f91358e1	247	AliForwardMCCorrectionsTask::DefineBins(TList* l)
563a673f	248	{
f91358e1	249	if (!fVtxBins) fVtxBins = new TObjArray(fVtxAxis.GetNbins(), 1);
	250	if (fVtxBins->GetEntries() > 0) return;
	251
	252	fVtxBins->SetOwner();
	253	for (Int_t i = 1; i <= fVtxAxis.GetNbins(); i++) {
	254	Double_t low = fVtxAxis.GetBinLowEdge(i);
	255	Double_t high = fVtxAxis.GetBinUpEdge(i);
	256	VtxBin* bin = new VtxBin(low, high, fEtaAxis);
	257	fVtxBins->AddAt(bin, i);
	258	bin->DefineOutput(l);
	259	}
563a673f	260	}
	261
	262	//____________________________________________________________________
	263	void
0bd4b00f	264	AliForwardMCCorrectionsTask::UserCreateOutputObjects()
563a673f	265	{
7984e5f7	266	//
	267	// Create output objects
	268	//
	269	//
563a673f	270	fList = new TList;
5bb5d1f6	271	fList->SetOwner();
5bb5d1f6	272	fList->SetName(Form("%sSums", GetName()));
563a673f	273
f91358e1	274	DefineBins(fList);
f91358e1	275
0bd4b00f	276	fHEvents = new TH1I(GetEventName(false,false),
563a673f	277	"Number of all events",
0bd4b00f	278	fVtxAxis.GetNbins(),
	279	fVtxAxis.GetXmin(),
	280	fVtxAxis.GetXmax());
563a673f	281	fHEvents->SetXTitle("v_{z} [cm]");
	282	fHEvents->SetYTitle("# of events");
	283	fHEvents->SetFillColor(kBlue+1);
	284	fHEvents->SetFillStyle(3001);
	285	fHEvents->SetDirectory(0);
	286	fList->Add(fHEvents);
	287
0bd4b00f	288	fHEventsTr = new TH1I(GetEventName(true, false),
0bd4b00f	289	"Number of triggered events",
563a673f	290	fVtxAxis.GetNbins(),
	291	fVtxAxis.GetXmin(),
	292	fVtxAxis.GetXmax());
	293	fHEventsTr->SetXTitle("v_{z} [cm]");
	294	fHEventsTr->SetYTitle("# of events");
	295	fHEventsTr->SetFillColor(kRed+1);
	296	fHEventsTr->SetFillStyle(3001);
	297	fHEventsTr->SetDirectory(0);
	298	fList->Add(fHEventsTr);
	299
0bd4b00f	300	fHEventsTrVtx = new TH1I(GetEventName(true,true),
563a673f	301	"Number of events w/trigger and vertex",
	302	fVtxAxis.GetNbins(),
	303	fVtxAxis.GetXmin(),
	304	fVtxAxis.GetXmax());
	305	fHEventsTrVtx->SetXTitle("v_{z} [cm]");
	306	fHEventsTrVtx->SetYTitle("# of events");
	307	fHEventsTrVtx->SetFillColor(kBlue+1);
	308	fHEventsTrVtx->SetFillStyle(3001);
	309	fHEventsTrVtx->SetDirectory(0);
	310	fList->Add(fHEventsTrVtx);
563a673f	311
f91358e1	312	// Copy axis objects to output
	313	TH1* vtxAxis = new TH1D("vtxAxis", "Vertex axis",
	314	fVtxAxis.GetNbins(),
	315	fVtxAxis.GetXmin(),
	316	fVtxAxis.GetXmax());
	317	TH1* etaAxis = new TH1D("etaAxis", "Eta axis",
	318	fEtaAxis.GetNbins(),
	319	fEtaAxis.GetXmin(),
	320	fEtaAxis.GetXmax());
	321	fList->Add(vtxAxis);
	322	fList->Add(etaAxis);
	323
	324	AliInfo(Form("Initialising sub-routines: %p, %p",
	325	&fInspector, &fTrackDensity));
5bb5d1f6	326	fInspector.DefineOutput(fList);
f91358e1	327	fInspector.Init(fVtxAxis);
f91358e1	328	fTrackDensity.DefineOutput(fList);
5bb5d1f6	329
563a673f	330	PostData(1, fList);
	331	}
	332	//____________________________________________________________________
	333	void
0bd4b00f	334	AliForwardMCCorrectionsTask::UserExec(Option_t*)
563a673f	335	{
7984e5f7	336	//
	337	// Process each event
	338	//
	339	// Parameters:
	340	// option Not used
	341	//
	342
563a673f	343	// Get the input data - MC event
563a673f	344	AliMCEvent* mcEvent = MCEvent();
cc83fca2	345	if (!mcEvent) {
563a673f	346	AliWarning("No MC event found");
	347	return;
	348	}
	349
	350	// Get the input data - ESD event
	351	AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
	352	if (!esd) {
	353	AliWarning("No ESD event found for input event");
	354	return;
	355	}
563a673f	356
5bb5d1f6	357	//--- Read run information -----------------------------------------
	358	if (fFirstEvent && esd->GetESDRun()) {
	359	fInspector.ReadRunDetails(esd);
	360
	361	AliInfo(Form("Initializing with parameters from the ESD:\n"
	362	" AliESDEvent::GetBeamEnergy() ->%f\n"
	363	" AliESDEvent::GetBeamType() ->%s\n"
	364	" AliESDEvent::GetCurrentL3() ->%f\n"
	365	" AliESDEvent::GetMagneticField()->%f\n"
	366	" AliESDEvent::GetRunNumber() ->%d\n",
	367	esd->GetBeamEnergy(),
	368	esd->GetBeamType(),
	369	esd->GetCurrentL3(),
	370	esd->GetMagneticField(),
	371	esd->GetRunNumber()));
f91358e1	372
f91358e1	373	Print();
5bb5d1f6	374	fFirstEvent = false;
563a673f	375	}
563a673f	376
5bb5d1f6	377	// Some variables
	378	UInt_t triggers; // Trigger bits
	379	Bool_t lowFlux; // Low flux flag
	380	UShort_t iVz; // Vertex bin from ESD
	381	Double_t vZ; // Z coordinate from ESD
	382	Double_t cent; // Centrality
	383	UShort_t iVzMc; // Vertex bin from MC
	384	Double_t vZMc; // Z coordinate of IP vertex from MC
	385	Double_t b; // Impact parameter
	386	Int_t nPart; // Number of participants
	387	Int_t nBin; // Number of binary collisions
	388	Double_t phiR; // Reaction plane from MC
	389	UShort_t nClusters;// Number of SPD clusters
	390	// Process the data
	391	UInt_t retESD = fInspector.Process(esd, triggers, lowFlux, iVz, vZ,
	392	cent, nClusters);
f91358e1	393	fInspector.ProcessMC(mcEvent, triggers, iVzMc, vZMc,
f91358e1	394	b, nPart, nBin, phiR);
5bb5d1f6	395
5bb5d1f6	396	Bool_t isInel = triggers & AliAODForwardMult::kInel;
5bb5d1f6	397	Bool_t hasVtx = retESD == AliFMDMCEventInspector::kOk;
f91358e1	398
563a673f	399	// Fill the event count histograms
5bb5d1f6	400	if (isInel) fHEventsTr->Fill(vZMc);
	401	if (isInel && hasVtx) fHEventsTrVtx->Fill(vZMc);
	402	fHEvents->Fill(vZMc);
563a673f	403
f91358e1	404	// Now find our vertex bin object
	405	VtxBin* bin = 0;
	406	if (iVzMc > 0 && iVzMc <= fVtxAxis.GetNbins())
	407	bin = static_cast<VtxBin*>(fVtxBins->At(iVzMc));
	408	if (!bin) {
	409	// AliError(Form("No vertex bin object @ %d (%f)", iVzMc, vZMc));
	410	return;
563a673f	411	}
563a673f	412
f91358e1	413	// Clear our ESD object
	414	fESDFMD.Clear();
	415
	416	// Get FMD data
	417	AliESDFMD* esdFMD = esd->GetFMDData();
	418
	419	// Now process our input data and store in own ESD object
	420	fTrackDensity.Calculate(esdFMD, mcEvent, vZMc, fESDFMD, bin->fPrimary);
	421	bin->fCounts->Fill(0.5);
	422
	423	// And then bin the data in our vtxbin
	424	for (UShort_t d=1; d<=3; d++) {
	425	UShort_t nr = (d == 1 ? 1 : 2);
	426	for (UShort_t q=0; q<nr; q++) {
	427	Char_t r = (q == 0 ? 'I' : 'O');
	428	UShort_t ns= (q == 0 ? 20 : 40);
	429	UShort_t nt= (q == 0 ? 512 : 256);
	430	TH2D* h = bin->fHists.Get(d,r);
	431
	432	for (UShort_t s=0; s<ns; s++) {
	433	for (UShort_t t=0; t<nt; t++) {
	434	Float_t mult = fESDFMD.Multiplicity(d,r,s,t);
	435
	436	if (mult == 0 \|\| mult > 20) continue;
	437
	438	Float_t phi = fESDFMD.Phi(d,r,s,t) / 180 * TMath::Pi();
	439	Float_t eta = fESDFMD.Eta(d,r,s,t);
	440	h->Fill(eta,phi,mult);
	441	} // for t
	442	} // for s
	443	} // for q
	444	} // for d
563a673f	445
bcd522ff	446	PostData(1, fList);
bcd522ff	447	}
563a673f	448	//____________________________________________________________________
563a673f	449	void
f91358e1	450	AliForwardMCCorrectionsTask::Terminate(Option_t*)
563a673f	451	{
7984e5f7	452	//
f91358e1	453	// End of job
7984e5f7	454	//
7984e5f7	455	// Parameters:
f91358e1	456	// option Not used
7984e5f7	457	//
7984e5f7	458
f91358e1	459	fList = dynamic_cast<TList*>(GetOutputData(1));
	460	if (!fList) {
	461	AliError("No output list defined");
	462	return;
563a673f	463	}
563a673f	464
bcd522ff	465	DefineBins(fList);
bcd522ff	466
f91358e1	467	// Output list
	468	TList* output = new TList;
	469	output->SetOwner();
	470	output->SetName(Form("%sResults", GetName()));
563a673f	471
f91358e1	472	// --- Fill correction object --------------------------------------
	473	AliFMDCorrSecondaryMap* corr = new AliFMDCorrSecondaryMap;
	474	corr->SetVertexAxis(fVtxAxis);
	475	corr->SetEtaAxis(fEtaAxis);
563a673f	476
f91358e1	477	TIter next(fVtxBins);
	478	VtxBin* bin = 0;
	479	UShort_t iVz = 1;
	480	while ((bin = static_cast<VtxBin*>(next())))
	481	bin->Finish(fList, output, iVz++, corr);
563a673f	482
f91358e1	483	output->Add(corr);
563a673f	484
f91358e1	485	PostData(2, output);
563a673f	486	}
563a673f	487
563a673f	488	//____________________________________________________________________
563a673f	489	void
f91358e1	490	AliForwardMCCorrectionsTask::Print(Option_t* option) const
563a673f	491	{
f91358e1	492	std::cout << ClassName() << "\n"
	493	<< " Vertex bins: " << fVtxAxis.GetNbins() << '\n'
	494	<< " Vertex range: [" << fVtxAxis.GetXmin()
	495	<< "," << fVtxAxis.GetXmax() << "]\n"
	496	<< " Eta bins: " << fEtaAxis.GetNbins() << '\n'
	497	<< " Eta range: [" << fEtaAxis.GetXmin()
	498	<< "," << fEtaAxis.GetXmax() << "]"
	499	<< std::endl;
	500	gROOT->IncreaseDirLevel();
	501	fInspector.Print(option);
	502	fTrackDensity.Print(option);
	503	gROOT->DecreaseDirLevel();
	504	}
7984e5f7	505
f91358e1	506	//====================================================================
	507	const char*
	508	AliForwardMCCorrectionsTask::VtxBin::BinName(Double_t low,
	509	Double_t high)
	510	{
	511	static TString buf;
	512	buf = Form("vtx%+05.1f_%+05.1f", low, high);
	513	buf.ReplaceAll("+", "p");
	514	buf.ReplaceAll("-", "m");
	515	buf.ReplaceAll(".", "d");
	516	return buf.Data();
	517	}
563a673f	518
563a673f	519
f91358e1	520	//____________________________________________________________________
	521	AliForwardMCCorrectionsTask::VtxBin::VtxBin()
	522	: fHists(),
	523	fPrimary(0),
	524	fCounts(0)
	525	{
	526	}
	527	//____________________________________________________________________
	528	AliForwardMCCorrectionsTask::VtxBin::VtxBin(Double_t low,
	529	Double_t high,
	530	const TAxis& axis)
	531	: TNamed(BinName(low, high),
	532	Form("%+5.1fcm<v_{z}<%+5.1fcm", low, high)),
	533	fHists(),
	534	fPrimary(0),
	535	fCounts(0)
	536	{
	537	fHists.Init(axis);
	538
	539	fPrimary = new TH2D("primary", "Primaries",
	540	axis.GetNbins(), axis.GetXmin(), axis.GetXmax(),
	541	40, 0, 2*TMath::Pi());
	542	fPrimary->SetXTitle("#eta");
	543	fPrimary->SetYTitle("#varphi [radians]");
	544	fPrimary->Sumw2();
	545	fPrimary->SetDirectory(0);
	546
	547	fCounts = new TH1D("counts", "Counts", 1, 0, 1);
	548	fCounts->SetXTitle("Events");
	549	fCounts->SetYTitle("# of Events");
	550	fCounts->SetDirectory(0);
	551	}
563a673f	552
f91358e1	553	//____________________________________________________________________
	554	AliForwardMCCorrectionsTask::VtxBin::VtxBin(const VtxBin& o)
	555	: TNamed(o),
	556	fHists(o.fHists),
	557	fPrimary(0),
	558	fCounts(0)
	559	{
	560	if (o.fPrimary) {
	561	fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
	562	fPrimary->SetDirectory(0);
563a673f	563	}
f91358e1	564	if (o.fCounts) {
	565	fCounts = static_cast<TH1D*>(o.fCounts->Clone());
	566	fCounts->SetDirectory(0);
563a673f	567	}
f91358e1	568	}
563a673f	569
f91358e1	570	//____________________________________________________________________
	571	AliForwardMCCorrectionsTask::VtxBin&
	572	AliForwardMCCorrectionsTask::VtxBin::operator=(const VtxBin& o)
	573	{
d015ecfe	574	if (&o == this) return *this;
f91358e1	575	TNamed::operator=(o);
	576	fHists = o.fHists;
	577	fPrimary = 0;
	578	fCounts = 0;
	579	if (o.fPrimary) {
	580	fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
	581	fPrimary->SetDirectory(0);
563a673f	582	}
f91358e1	583	if (o.fCounts) {
	584	fCounts = static_cast<TH1D*>(o.fCounts->Clone());
	585	fCounts->SetDirectory(0);
	586	}
	587	return *this;
	588	}
563a673f	589
f91358e1	590	//____________________________________________________________________
	591	void
	592	AliForwardMCCorrectionsTask::VtxBin::DefineOutput(TList* l)
	593	{
	594	TList* d = new TList;
	595	d->SetName(GetName());
	596	d->SetOwner();
	597	l->Add(d);
	598
	599	d->Add(fHists.fFMD1i);
	600	d->Add(fHists.fFMD2i);
	601	d->Add(fHists.fFMD2o);
	602	d->Add(fHists.fFMD3i);
	603	d->Add(fHists.fFMD3o);
	604
	605	d->Add(fPrimary);
	606	d->Add(fCounts);
563a673f	607	}
563a673f	608
f91358e1	609	//____________________________________________________________________
	610	TH2D*
	611	AliForwardMCCorrectionsTask::VtxBin::MakeBg(const TH2D* hits,
	612	const TH2D* primary) const
	613	{
	614	TH2D* h = static_cast<TH2D*>(hits->Clone());
	615	h->SetDirectory(0);
	616	TString n(h->GetName());
	617	n.ReplaceAll("_cache", "");
	618	h->SetName(n);
	619	h->Divide(primary);
	620
	621	return h;
	622	}
	623
563a673f	624	//____________________________________________________________________
563a673f	625	void
f91358e1	626	AliForwardMCCorrectionsTask::VtxBin::Finish(const TList* input,
	627	TList* output,
	628	UShort_t iVz,
	629	AliFMDCorrSecondaryMap* map)
563a673f	630	{
f91358e1	631	TList* out = new TList;
	632	out->SetName(GetName());
	633	out->SetOwner();
	634	output->Add(out);
	635
	636	TList* l = static_cast<TList*>(input->FindObject(GetName()));
	637	if (!l) {
	638	AliError(Form("List %s not found in %s", GetName(), input->GetName()));
	639	return;
	640	}
	641
	642	TH2D* fmd1i = static_cast<TH2D*>(l->FindObject("FMD1I_cache"));
	643	TH2D* fmd2i = static_cast<TH2D*>(l->FindObject("FMD2I_cache"));
	644	TH2D* fmd2o = static_cast<TH2D*>(l->FindObject("FMD2O_cache"));
	645	TH2D* fmd3i = static_cast<TH2D*>(l->FindObject("FMD3I_cache"));
	646	TH2D* fmd3o = static_cast<TH2D*>(l->FindObject("FMD3O_cache"));
	647	TH2D* primO = static_cast<TH2D*>(l->FindObject("primary"));
	648	if (!fmd1i \|\| !fmd2i \|\| !fmd2o \|\| !fmd3i \|\| !fmd3o \|\| !primO) {
	649	AliError(Form("Missing histogram(s): %p,%p,%p,%p,%p,%p",
	650	fmd1i, fmd2i, fmd2o, fmd3i, fmd3o, primO));
	651	return;
	652	}
	653
	654	// Half coverage in phi for inners
	655	TH2D* primI = static_cast<TH2D*>(primO->Clone());
	656	primI->SetDirectory(0);
	657	primI->RebinY(2);
	658
	659	TH2D* bg1i = MakeBg(fmd1i, primI);
	660	TH2D* bg2i = MakeBg(fmd2i, primI);
	661	TH2D* bg2o = MakeBg(fmd2o, primO);
	662	TH2D* bg3i = MakeBg(fmd3i, primI);
	663	TH2D* bg3o = MakeBg(fmd3o, primO);
	664	map->SetCorrection(1, 'I', iVz, bg1i);
	665	map->SetCorrection(2, 'I', iVz, bg2i);
	666	map->SetCorrection(2, 'O', iVz, bg2o);
	667	map->SetCorrection(3, 'I', iVz, bg3i);
	668	map->SetCorrection(3, 'O', iVz, bg3o);
	669	out->Add(bg1i);
	670	out->Add(bg2i);
	671	out->Add(bg2o);
	672	out->Add(bg3i);
	673	out->Add(bg3o);
	674
563a673f	675	}
	676
	677	//
	678	// EOF
	679	//