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

#include "AliForwardMCMultiplicityTask.h"
#include "AliTriggerAnalysis.h"
#include "AliPhysicsSelection.h"
#include "AliLog.h"
// #include "AliFMDAnaParameters.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)
{
}

//____________________________________________________________________
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)
{
  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) 
{
  DefineOutput(1, TList::Class());
}

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