[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliForwardMCMultiplicityTask.cxx

// 
// Calculate the multiplicity in the forward regions event-by-event 
// 
// Inputs: 
//   - AliESDEvent 
//   - Kinematics
//   - Track references
//
// Outputs: 
//   - AliAODForwardMult 
// 
// Histograms 
//   
// Corrections used 
// 
#include "AliForwardMCMultiplicityTask.h"
#include "AliTriggerAnalysis.h"
#include "AliPhysicsSelection.h"
#include "AliLog.h"
#include "AliESDEvent.h"
#include "AliAODHandler.h"
#include "AliMultiplicity.h"
#include "AliInputEventHandler.h"
#include "AliForwardCorrectionManager.h"
#include "AliAnalysisManager.h"
#include <TH1.h>
#include <TDirectory.h>
#include <TTree.h>
#include <TROOT.h>
#include <iostream>
#include <iomanip>

//====================================================================
AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask()
  : AliForwardMultiplicityBase(),
    fHData(0),
    fESDFMD(),
    fHistos(),
    fAODFMD(),
    fMCESDFMD(),
    fMCHistos(),
    fMCAODFMD(),
    fEventInspector(),
    fEnergyFitter(),
    fSharingFilter(),
    fDensityCalculator(),
    fCorrections(),
    fHistCollector(),
    fList(0)
{
  // 
  // Constructor
  //
}

//____________________________________________________________________
AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const char* name)
  : AliForwardMultiplicityBase(name), 
    fHData(0),
    fESDFMD(),
    fHistos(),
    fAODFMD(kFALSE),
    fMCESDFMD(),
    fMCHistos(),
    fMCAODFMD(kTRUE),
    fEventInspector("event"),
    fEnergyFitter("energy"),
    fSharingFilter("sharing"), 
    fDensityCalculator("density"),
    fCorrections("corrections"),
    fHistCollector("collector"),
    fList(0)
{
  // 
  // Constructor 
  // 
  // Parameters:
  //    name Name of task 
  //
  DefineOutput(1, TList::Class());
}

//____________________________________________________________________
AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const AliForwardMCMultiplicityTask& o)
  : AliForwardMultiplicityBase(o),
    fHData(o.fHData),
    fESDFMD(o.fESDFMD),
    fHistos(o.fHistos),
    fAODFMD(o.fAODFMD),
    fMCESDFMD(o.fMCESDFMD),
    fMCHistos(o.fMCHistos),
    fMCAODFMD(o.fMCAODFMD),
    fEventInspector(o.fEventInspector),
    fEnergyFitter(o.fEnergyFitter),
   fSharingFilter(o.fSharingFilter),
    fDensityCalculator(o.fDensityCalculator),
    fCorrections(o.fCorrections),
    fHistCollector(o.fHistCollector),
    fList(o.fList) 
{
  // 
  // Copy constructor 
  // 
  // Parameters:
  //    o Object to copy from 
  //
  DefineOutput(1, TList::Class());
}

//____________________________________________________________________
AliForwardMCMultiplicityTask&
AliForwardMCMultiplicityTask::operator=(const AliForwardMCMultiplicityTask& o)
{
  // 
  // Assignment operator 
  // 
  // Parameters:
  //    o Object to assign from 
  // 
  // Return:
  //    Reference to this object 
  //
  AliForwardMultiplicityBase::operator=(o);

  fHData             = o.fHData;
  fEventInspector    = o.fEventInspector;
  fEnergyFitter      = o.fEnergyFitter;
  fSharingFilter     = o.fSharingFilter;
  fDensityCalculator = o.fDensityCalculator;
  fCorrections       = o.fCorrections;
  fHistCollector     = o.fHistCollector;
  fHistos            = o.fHistos;
  fAODFMD            = o.fAODFMD;
  fMCHistos          = o.fMCHistos;
  fMCAODFMD          = o.fMCAODFMD;
  fList              = o.fList;

  return *this;
}

//____________________________________________________________________
void
AliForwardMCMultiplicityTask::SetDebug(Int_t dbg)
{
  // 
  // Set debug level 
  // 
  // Parameters:
  //    dbg debug level
  //
  fEventInspector.SetDebug(dbg);
  fEnergyFitter.SetDebug(dbg);
  fSharingFilter.SetDebug(dbg);
  fDensityCalculator.SetDebug(dbg);
  fCorrections.SetDebug(dbg);
  fHistCollector.SetDebug(dbg);
}
//____________________________________________________________________
void
AliForwardMCMultiplicityTask::InitializeSubs()
{
  // 
  // Initialise the sub objects and stuff.  Called on first event 
  // 
  //
  AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  fcm.Init(fEventInspector.GetCollisionSystem(), 
	   fEventInspector.GetEnergy(),
	   fEventInspector.GetField(), 
	   true); // Last true is for MC flag 
  // Check that we have the energy loss fits, needed by 
  //   AliFMDSharingFilter 
  //   AliFMDDensityCalculator 
  if (!fcm.GetELossFit()) { 
    AliFatal(Form("No energy loss fits"));
    return;
  }
  // Check that we have the double hit correction - (optionally) used by 
  //  AliFMDDensityCalculator 
  if (!fcm.GetDoubleHit()) {
    AliWarning("No double hit corrections"); 
  }
  // Check that we have the secondary maps, needed by 
  //   AliFMDCorrections 
  //   AliFMDHistCollector
  if (!fcm.GetSecondaryMap()) {
    AliFatal("No secondary corrections");
    return;
  }
  // Check that we have the vertex bias correction, needed by 
  //   AliFMDCorrections 
  if (!fcm.GetVertexBias()) { 
    AliFatal("No event vertex bias corrections");
    return;
  }
  // Check that we have the merging efficiencies, optionally used by 
  //   AliFMDCorrections 
  if (!fcm.GetMergingEfficiency()) {
    AliWarning("No merging efficiencies");
  }


  const TAxis* pe = fcm.GetEtaAxis();
  const TAxis* pv = fcm.GetVertexAxis();
  if (!pe) AliFatal("No eta axis defined");
  if (!pv) AliFatal("No vertex axis defined");

  fHistos.Init(*pe);
  fAODFMD.Init(*pe);
  fMCHistos.Init(*pe);
  fMCAODFMD.Init(*pe);

  fHData = static_cast<TH2D*>(fAODFMD.GetHistogram().Clone("d2Ndetadphi"));
  fHData->SetStats(0);
  fHData->SetDirectory(0);
  fList->Add(fHData);

  fEnergyFitter.Init(*pe);
  fEventInspector.Init(*pv);
  fDensityCalculator.Init(*pe);
  fCorrections.Init(*pe);
  fHistCollector.Init(*pv);

  this->Print();
}

//____________________________________________________________________
void
AliForwardMCMultiplicityTask::UserCreateOutputObjects()
{
  // 
  // Create output objects 
  // 
  //
  fList = new TList;

  AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
  AliAODHandler*      ah = 
    dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
  if (!ah) AliFatal("No AOD output handler set in analysis manager");
    
    
  TObject* obj = &fAODFMD;
  ah->AddBranch("AliAODForwardMult", &obj);

  TObject* mcobj = &fMCAODFMD;
  ah->AddBranch("AliAODForwardMult", &mcobj);

  fEventInspector.DefineOutput(fList);
  fEnergyFitter.DefineOutput(fList);
  fSharingFilter.DefineOutput(fList);
  fDensityCalculator.DefineOutput(fList);
  fCorrections.DefineOutput(fList);

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

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

  // On the first event, initialize the parameters 
  if (fFirstEvent && esd->GetESDRun()) { 
    fEventInspector.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()));

    fFirstEvent = false;

    InitializeSubs();
  }
  // Clear stuff 
  fHistos.Clear();
  fESDFMD.Clear();
  fAODFMD.Clear();
  fMCHistos.Clear();
  fMCESDFMD.Clear();
  fMCAODFMD.Clear();

  Bool_t   lowFlux  = kFALSE;
  UInt_t   triggers = 0;
  UShort_t ivz      = 0;
  Double_t vz       = 0;
  UInt_t   found    = fEventInspector.Process(esd, triggers, lowFlux, ivz, vz);
  if (found & AliFMDEventInspector::kNoEvent)    return;
  if (found & AliFMDEventInspector::kNoTriggers) return;
 
  // Set trigger bits, and mark this event for storage 
  fAODFMD.SetTriggerBits(triggers);
  fMCAODFMD.SetTriggerBits(triggers);
  MarkEventForStore();

  if (found & AliFMDEventInspector::kNoSPD)     return;
  if (found & AliFMDEventInspector::kNoFMD)     return;
  if (found & AliFMDEventInspector::kNoVertex)  return;
  fAODFMD.SetIpZ(vz);
  fMCAODFMD.SetIpZ(vz);

  if (found & AliFMDEventInspector::kBadVertex) return;

  // We we do not want to use low flux specific code, we disable it here. 
  if (!fEnableLowFlux) lowFlux = false;

  // Get FMD data 
  AliESDFMD*  esdFMD  = esd->GetFMDData();
  AliMCEvent* mcEvent = MCEvent();
  // Apply the sharing filter (or hit merging or clustering if you like)
  if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD)) { 
    AliWarning("Sharing filter failed!");
    return;
  }
  if (!fSharingFilter.FilterMC(*esdFMD, *mcEvent, vz, fMCESDFMD)) { 
    AliWarning("MC Sharing filter failed!");
    return;
  }
  fSharingFilter.CompareResults(fESDFMD, fMCESDFMD);

  // Do the energy stuff 
  if (!fEnergyFitter.Accumulate(*esdFMD, triggers & AliAODForwardMult::kEmpty)){
    AliWarning("Energy fitter failed");
    return;
  }

  // Calculate the inclusive charged particle density 
  if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) { 
    AliWarning("Density calculator failed!");
    return;
  }
  if (!fDensityCalculator.CalculateMC(fMCESDFMD, fMCHistos)) { 
    AliWarning("MC Density calculator failed!");
    return;
  }
  fDensityCalculator.CompareResults(fHistos, fMCHistos);
  
  // Do the secondary and other corrections. 
  if (!fCorrections.Correct(fHistos, ivz)) { 
    AliWarning("Corrections failed");
    return;
  }
  if (!fCorrections.CorrectMC(fMCHistos, ivz)) { 
    AliWarning("MC Corrections failed");
    return;
  }
  fCorrections.CompareResults(fHistos, fMCHistos);
    
  if (!fHistCollector.Collect(fHistos, ivz, fAODFMD.GetHistogram())) {
    AliWarning("Histogram collector failed");
    return;
  }
  if (!fHistCollector.Collect(fMCHistos, ivz, fMCAODFMD.GetHistogram())) {
    AliWarning("MC Histogram collector failed");
    return;
  }

  if (fAODFMD.IsTriggerBits(AliAODForwardMult::kInel))
    fHData->Add(&(fAODFMD.GetHistogram()));

  PostData(1, fList);
}

//____________________________________________________________________
void
AliForwardMCMultiplicityTask::Terminate(Option_t*)
{
  // 
  // End of job
  // 
  // Parameters:
  //    option Not used 
  //
  TList* list = dynamic_cast<TList*>(GetOutputData(1));
  if (!list) {
    AliError(Form("No output list defined (%p)", GetOutputData(1)));
    if (GetOutputData(1)) GetOutputData(1)->Print();
    return;
  }
  
  // Get our histograms from the container 
  TH1I* hEventsTr    = 0;
  TH1I* hEventsTrVtx = 0;
  TH1I* hTriggers    = 0;
  if (!fEventInspector.FetchHistograms(list, hEventsTr, 
				       hEventsTrVtx, hTriggers)) { 
    AliError(Form("Didn't get histograms from event selector "
		  "(hEventsTr=%p,hEventsTrVtx=%p)", 
		  hEventsTr, hEventsTrVtx));
    list->ls();
    return;
  }

  TH2D* hData        = static_cast<TH2D*>(list->FindObject("d2Ndetadphi"));
  if (!hData) { 
    AliError(Form("Couldn't get our summed histogram from output "
		  "list %s (d2Ndetadphi=%p)", list->GetName(), hData));
    list->ls();
    return;
  }
  
  // TH1D* dNdeta = fHData->ProjectionX("dNdeta", 0, -1, "e");
  TH1D* dNdeta = hData->ProjectionX("dNdeta", 1, -1, "e");
  TH1D* norm   = hData->ProjectionX("norm",   0,  1,  "");
  dNdeta->SetTitle("dN_{ch}/d#eta in the forward regions");
  dNdeta->SetYTitle("#frac{1}{N}#frac{dN_{ch}}{d#eta}");
  dNdeta->Divide(norm);
  dNdeta->SetStats(0);
  dNdeta->Scale(Double_t(hEventsTrVtx->GetEntries())/hEventsTr->GetEntries(),
		"width");
  list->Add(dNdeta);
  list->Add(norm);

  fEnergyFitter.Fit(list);
  fSharingFilter.ScaleHistograms(list,hEventsTr->Integral());
  fDensityCalculator.ScaleHistograms(list,hEventsTrVtx->Integral());
  fCorrections.ScaleHistograms(list,hEventsTrVtx->Integral());
}

//____________________________________________________________________
void
AliForwardMCMultiplicityTask::Print(Option_t* option) const
{
  // 
  // Print information 
  // 
  // Parameters:
  //    option Not used
  //
  AliForwardMultiplicityBase::Print(option);
  gROOT->IncreaseDirLevel();
  fEventInspector   .Print(option);
  fEnergyFitter     .Print(option);    
  fSharingFilter    .Print(option);
  fDensityCalculator.Print(option);
  fCorrections      .Print(option);
  fHistCollector    .Print(option);
  gROOT->DecreaseDirLevel();
}

//
// EOF
//
Commit	Line	Data
7984e5f7	1	//
	2	// Calculate the multiplicity in the forward regions event-by-event
	3	//
	4	// Inputs:
	5	// - AliESDEvent
	6	// - Kinematics
	7	// - Track references
	8	//
	9	// Outputs:
	10	// - AliAODForwardMult
	11	//
	12	// Histograms
	13	//
	14	// Corrections used
	15	//
285e7b27	16	#include "AliForwardMCMultiplicityTask.h"
	17	#include "AliTriggerAnalysis.h"
	18	#include "AliPhysicsSelection.h"
	19	#include "AliLog.h"
285e7b27	20	#include "AliESDEvent.h"
	21	#include "AliAODHandler.h"
	22	#include "AliMultiplicity.h"
	23	#include "AliInputEventHandler.h"
	24	#include "AliForwardCorrectionManager.h"
	25	#include "AliAnalysisManager.h"
	26	#include <TH1.h>
	27	#include <TDirectory.h>
	28	#include <TTree.h>
	29	#include <TROOT.h>
	30	#include <iostream>
	31	#include <iomanip>
	32
	33	//====================================================================
	34	AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask()
	35	: AliForwardMultiplicityBase(),
	36	fHData(0),
	37	fESDFMD(),
	38	fHistos(),
	39	fAODFMD(),
	40	fMCESDFMD(),
	41	fMCHistos(),
	42	fMCAODFMD(),
	43	fEventInspector(),
	44	fEnergyFitter(),
	45	fSharingFilter(),
	46	fDensityCalculator(),
	47	fCorrections(),
	48	fHistCollector(),
	49	fList(0)
	50	{
7984e5f7	51	//
	52	// Constructor
	53	//
285e7b27	54	}
	55
	56	//____________________________________________________________________
	57	AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const char* name)
	58	: AliForwardMultiplicityBase(name),
	59	fHData(0),
	60	fESDFMD(),
	61	fHistos(),
	62	fAODFMD(kFALSE),
	63	fMCESDFMD(),
	64	fMCHistos(),
	65	fMCAODFMD(kTRUE),
	66	fEventInspector("event"),
	67	fEnergyFitter("energy"),
	68	fSharingFilter("sharing"),
	69	fDensityCalculator("density"),
	70	fCorrections("corrections"),
	71	fHistCollector("collector"),
	72	fList(0)
	73	{
7984e5f7	74	//
	75	// Constructor
	76	//
	77	// Parameters:
	78	// name Name of task
	79	//
285e7b27	80	DefineOutput(1, TList::Class());
	81	}
	82
	83	//____________________________________________________________________
	84	AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const AliForwardMCMultiplicityTask& o)
	85	: AliForwardMultiplicityBase(o),
	86	fHData(o.fHData),
	87	fESDFMD(o.fESDFMD),
	88	fHistos(o.fHistos),
	89	fAODFMD(o.fAODFMD),
	90	fMCESDFMD(o.fMCESDFMD),
	91	fMCHistos(o.fMCHistos),
	92	fMCAODFMD(o.fMCAODFMD),
	93	fEventInspector(o.fEventInspector),
	94	fEnergyFitter(o.fEnergyFitter),
	95	fSharingFilter(o.fSharingFilter),
	96	fDensityCalculator(o.fDensityCalculator),
	97	fCorrections(o.fCorrections),
	98	fHistCollector(o.fHistCollector),
	99	fList(o.fList)
	100	{
7984e5f7	101	//
	102	// Copy constructor
	103	//
	104	// Parameters:
	105	// o Object to copy from
	106	//
285e7b27	107	DefineOutput(1, TList::Class());
	108	}
	109
	110	//____________________________________________________________________
	111	AliForwardMCMultiplicityTask&
	112	AliForwardMCMultiplicityTask::operator=(const AliForwardMCMultiplicityTask& o)
	113	{
7984e5f7	114	//
	115	// Assignment operator
	116	//
	117	// Parameters:
	118	// o Object to assign from
	119	//
	120	// Return:
	121	// Reference to this object
	122	//
285e7b27	123	AliForwardMultiplicityBase::operator=(o);
	124
	125	fHData = o.fHData;
	126	fEventInspector = o.fEventInspector;
	127	fEnergyFitter = o.fEnergyFitter;
	128	fSharingFilter = o.fSharingFilter;
	129	fDensityCalculator = o.fDensityCalculator;
	130	fCorrections = o.fCorrections;
	131	fHistCollector = o.fHistCollector;
	132	fHistos = o.fHistos;
	133	fAODFMD = o.fAODFMD;
	134	fMCHistos = o.fMCHistos;
	135	fMCAODFMD = o.fMCAODFMD;
	136	fList = o.fList;
	137
	138	return *this;
	139	}
	140
	141	//____________________________________________________________________
	142	void
	143	AliForwardMCMultiplicityTask::SetDebug(Int_t dbg)
	144	{
7984e5f7	145	//
	146	// Set debug level
	147	//
	148	// Parameters:
	149	// dbg debug level
	150	//
285e7b27	151	fEventInspector.SetDebug(dbg);
	152	fEnergyFitter.SetDebug(dbg);
	153	fSharingFilter.SetDebug(dbg);
	154	fDensityCalculator.SetDebug(dbg);
	155	fCorrections.SetDebug(dbg);
	156	fHistCollector.SetDebug(dbg);
	157	}
	158	//____________________________________________________________________
	159	void
	160	AliForwardMCMultiplicityTask::InitializeSubs()
	161	{
7984e5f7	162	//
	163	// Initialise the sub objects and stuff. Called on first event
	164	//
	165	//
285e7b27	166	AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
	167	fcm.Init(fEventInspector.GetCollisionSystem(),
	168	fEventInspector.GetEnergy(),
	169	fEventInspector.GetField(),
	170	true); // Last true is for MC flag
	171	// Check that we have the energy loss fits, needed by
	172	// AliFMDSharingFilter
	173	// AliFMDDensityCalculator
	174	if (!fcm.GetELossFit()) {
	175	AliFatal(Form("No energy loss fits"));
	176	return;
	177	}
	178	// Check that we have the double hit correction - (optionally) used by
	179	// AliFMDDensityCalculator
	180	if (!fcm.GetDoubleHit()) {
	181	AliWarning("No double hit corrections");
	182	}
	183	// Check that we have the secondary maps, needed by
	184	// AliFMDCorrections
	185	// AliFMDHistCollector
	186	if (!fcm.GetSecondaryMap()) {
	187	AliFatal("No secondary corrections");
	188	return;
	189	}
	190	// Check that we have the vertex bias correction, needed by
	191	// AliFMDCorrections
	192	if (!fcm.GetVertexBias()) {
	193	AliFatal("No event vertex bias corrections");
	194	return;
	195	}
	196	// Check that we have the merging efficiencies, optionally used by
	197	// AliFMDCorrections
	198	if (!fcm.GetMergingEfficiency()) {
	199	AliWarning("No merging efficiencies");
	200	}
	201
	202
	203	const TAxis* pe = fcm.GetEtaAxis();
	204	const TAxis* pv = fcm.GetVertexAxis();
	205	if (!pe) AliFatal("No eta axis defined");
	206	if (!pv) AliFatal("No vertex axis defined");
	207
	208	fHistos.Init(*pe);
	209	fAODFMD.Init(*pe);
	210	fMCHistos.Init(*pe);
	211	fMCAODFMD.Init(*pe);
	212
	213	fHData = static_cast<TH2D*>(fAODFMD.GetHistogram().Clone("d2Ndetadphi"));
	214	fHData->SetStats(0);
	215	fHData->SetDirectory(0);
	216	fList->Add(fHData);
	217
	218	fEnergyFitter.Init(*pe);
	219	fEventInspector.Init(*pv);
	220	fDensityCalculator.Init(*pe);
	221	fCorrections.Init(*pe);
	222	fHistCollector.Init(*pv);
	223
	224	this->Print();
	225	}
	226
	227	//____________________________________________________________________
	228	void
	229	AliForwardMCMultiplicityTask::UserCreateOutputObjects()
230	{
7984e5f7	231	//
	232	// Create output objects
	233	//
	234	//
285e7b27	235	fList = new TList;
	236
	237	AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
	238	AliAODHandler* ah =
	239	dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
	240	if (!ah) AliFatal("No AOD output handler set in analysis manager");
	241
	242
	243	TObject* obj = &fAODFMD;
	244	ah->AddBranch("AliAODForwardMult", &obj);
	245
	246	TObject* mcobj = &fMCAODFMD;
	247	ah->AddBranch("AliAODForwardMult", &mcobj);
	248
	249	fEventInspector.DefineOutput(fList);
	250	fEnergyFitter.DefineOutput(fList);
	251	fSharingFilter.DefineOutput(fList);
	252	fDensityCalculator.DefineOutput(fList);
	253	fCorrections.DefineOutput(fList);
	254
	255	PostData(1, fList);
	256	}
	257	//____________________________________________________________________
	258	void
	259	AliForwardMCMultiplicityTask::UserExec(Option_t*)
	260	{
7984e5f7	261	//
	262	// Process each event
	263	//
	264	// Parameters:
	265	// option Not used
	266	//
	267
285e7b27	268	// Get the input data
	269	AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
	270	if (!esd) {
	271	AliWarning("No ESD event found for input event");
	272	return;
	273	}
	274
	275	// On the first event, initialize the parameters
	276	if (fFirstEvent && esd->GetESDRun()) {
	277	fEventInspector.ReadRunDetails(esd);
	278
	279	AliInfo(Form("Initializing with parameters from the ESD:\n"
	280	" AliESDEvent::GetBeamEnergy() ->%f\n"
	281	" AliESDEvent::GetBeamType() ->%s\n"
	282	" AliESDEvent::GetCurrentL3() ->%f\n"
	283	" AliESDEvent::GetMagneticField()->%f\n"
	284	" AliESDEvent::GetRunNumber() ->%d\n",
	285	esd->GetBeamEnergy(),
	286	esd->GetBeamType(),
	287	esd->GetCurrentL3(),
	288	esd->GetMagneticField(),
	289	esd->GetRunNumber()));
	290
	291	fFirstEvent = false;
	292
	293	InitializeSubs();
	294	}
	295	// Clear stuff
	296	fHistos.Clear();
	297	fESDFMD.Clear();
	298	fAODFMD.Clear();
	299	fMCHistos.Clear();
	300	fMCESDFMD.Clear();
	301	fMCAODFMD.Clear();
	302
	303	Bool_t lowFlux = kFALSE;
	304	UInt_t triggers = 0;
	305	UShort_t ivz = 0;
	306	Double_t vz = 0;
	307	UInt_t found = fEventInspector.Process(esd, triggers, lowFlux, ivz, vz);
	308	if (found & AliFMDEventInspector::kNoEvent) return;
	309	if (found & AliFMDEventInspector::kNoTriggers) return;
	310
	311	// Set trigger bits, and mark this event for storage
	312	fAODFMD.SetTriggerBits(triggers);
	313	fMCAODFMD.SetTriggerBits(triggers);
	314	MarkEventForStore();
	315
	316	if (found & AliFMDEventInspector::kNoSPD) return;
	317	if (found & AliFMDEventInspector::kNoFMD) return;
	318	if (found & AliFMDEventInspector::kNoVertex) return;
	319	fAODFMD.SetIpZ(vz);
	320	fMCAODFMD.SetIpZ(vz);
	321
	322	if (found & AliFMDEventInspector::kBadVertex) return;
	323
	324	// We we do not want to use low flux specific code, we disable it here.
	325	if (!fEnableLowFlux) lowFlux = false;
	326
	327	// Get FMD data
	328	AliESDFMD* esdFMD = esd->GetFMDData();
	329	AliMCEvent* mcEvent = MCEvent();
	330	// Apply the sharing filter (or hit merging or clustering if you like)
	331	if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD)) {
332	AliWarning("Sharing filter failed!");
333	return;
334	}
335	if (!fSharingFilter.FilterMC(esdFMD, mcEvent, vz, fMCESDFMD)) {
336	AliWarning("MC Sharing filter failed!");
337	return;
338	}
339	fSharingFilter.CompareResults(fESDFMD, fMCESDFMD);
340
341	// Do the energy stuff
342	if (!fEnergyFitter.Accumulate(*esdFMD, triggers & AliAODForwardMult::kEmpty)){
343	AliWarning("Energy fitter failed");
344	return;
345	}
346
347	// Calculate the inclusive charged particle density
348	if (!fDensityCalculator.Calculate(fESDFMD, fHistos, ivz, lowFlux)) {
349	AliWarning("Density calculator failed!");
350	return;
351	}
352	if (!fDensityCalculator.CalculateMC(fMCESDFMD, fMCHistos)) {
353	AliWarning("MC Density calculator failed!");
354	return;
355	}
356	fDensityCalculator.CompareResults(fHistos, fMCHistos);
357
358	// Do the secondary and other corrections.
359	if (!fCorrections.Correct(fHistos, ivz)) {
360	AliWarning("Corrections failed");
361	return;
362	}
363	if (!fCorrections.CorrectMC(fMCHistos, ivz)) {
364	AliWarning("MC Corrections failed");
365	return;
366	}
367	fCorrections.CompareResults(fHistos, fMCHistos);
368
369	if (!fHistCollector.Collect(fHistos, ivz, fAODFMD.GetHistogram())) {
370	AliWarning("Histogram collector failed");
371	return;
372	}
373	if (!fHistCollector.Collect(fMCHistos, ivz, fMCAODFMD.GetHistogram())) {
374	AliWarning("MC Histogram collector failed");
375	return;
376	}
377
378	if (fAODFMD.IsTriggerBits(AliAODForwardMult::kInel))
379	fHData->Add(&(fAODFMD.GetHistogram()));
380
381	PostData(1, fList);
382	}
383
384	//____________________________________________________________________
385	void
386	AliForwardMCMultiplicityTask::Terminate(Option_t*)
387	{
7984e5f7	388	//
	389	// End of job
	390	//
	391	// Parameters:
	392	// option Not used
	393	//
285e7b27	394	TList* list = dynamic_cast<TList*>(GetOutputData(1));
	395	if (!list) {
	396	AliError(Form("No output list defined (%p)", GetOutputData(1)));
	397	if (GetOutputData(1)) GetOutputData(1)->Print();
	398	return;
	399	}
	400
	401	// Get our histograms from the container
	402	TH1I* hEventsTr = 0;
	403	TH1I* hEventsTrVtx = 0;
	404	TH1I* hTriggers = 0;
	405	if (!fEventInspector.FetchHistograms(list, hEventsTr,
	406	hEventsTrVtx, hTriggers)) {
	407	AliError(Form("Didn't get histograms from event selector "
	408	"(hEventsTr=%p,hEventsTrVtx=%p)",
	409	hEventsTr, hEventsTrVtx));
	410	list->ls();
	411	return;
	412	}
	413
	414	TH2D* hData = static_cast<TH2D*>(list->FindObject("d2Ndetadphi"));
	415	if (!hData) {
	416	AliError(Form("Couldn't get our summed histogram from output "
	417	"list %s (d2Ndetadphi=%p)", list->GetName(), hData));
	418	list->ls();
	419	return;
	420	}
	421
	422	// TH1D* dNdeta = fHData->ProjectionX("dNdeta", 0, -1, "e");
	423	TH1D* dNdeta = hData->ProjectionX("dNdeta", 1, -1, "e");
	424	TH1D* norm = hData->ProjectionX("norm", 0, 1, "");
	425	dNdeta->SetTitle("dN_{ch}/d#eta in the forward regions");
	426	dNdeta->SetYTitle("#frac{1}{N}#frac{dN_{ch}}{d#eta}");
	427	dNdeta->Divide(norm);
	428	dNdeta->SetStats(0);
	429	dNdeta->Scale(Double_t(hEventsTrVtx->GetEntries())/hEventsTr->GetEntries(),
	430	"width");
	431	list->Add(dNdeta);
	432	list->Add(norm);
	433
	434	fEnergyFitter.Fit(list);
	435	fSharingFilter.ScaleHistograms(list,hEventsTr->Integral());
	436	fDensityCalculator.ScaleHistograms(list,hEventsTrVtx->Integral());
	437	fCorrections.ScaleHistograms(list,hEventsTrVtx->Integral());
	438	}
	439
	440	//____________________________________________________________________
	441	void
	442	AliForwardMCMultiplicityTask::Print(Option_t* option) const
	443	{
7984e5f7	444	//
	445	// Print information
	446	//
	447	// Parameters:
	448	// option Not used
	449	//
285e7b27	450	AliForwardMultiplicityBase::Print(option);
	451	gROOT->IncreaseDirLevel();
	452	fEventInspector .Print(option);
	453	fEnergyFitter .Print(option);
	454	fSharingFilter .Print(option);
	455	fDensityCalculator.Print(option);
	456	fCorrections .Print(option);
	457	fHistCollector .Print(option);
	458	gROOT->DecreaseDirLevel();
	459	}
	460
	461	//
	462	// EOF
	463	//