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

// 
// Calculate the multiplicity in the forward regions event-by-event 
// 
// Inputs: 
//   - AliESDEvent 
//
// Outputs: 
//   - AliAODForwardMult 
// 
// Histograms 
//   
// Corrections used 
//
#include "AliForwardQATask.h"
#include "AliForwardUtil.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 "AliAODForwardMult.h"
#include <TH1.h>
#include <TDirectory.h>
#include <TTree.h>
#include <TROOT.h>
#include <TStopwatch.h>

//====================================================================
AliForwardQATask::AliForwardQATask()
  : AliAnalysisTaskSE(),
    fEnableLowFlux(false), 
    fFirstEvent(true),
    fCorrManager(0),
    fESDFMD(),
    fHistos(),
    fEventInspector(),
    fEnergyFitter(),
    fSharingFilter(),
    fDensityCalculator(),
    fList(0),
    fDebug(0)
{
  // 
  // Constructor
  //
}

//____________________________________________________________________
AliForwardQATask::AliForwardQATask(const char* name)
  : AliAnalysisTaskSE(name),
    fEnableLowFlux(false), 
    fFirstEvent(true),
    fCorrManager(0),
    fESDFMD(),
    fHistos(),
    fEventInspector("event"),
    fEnergyFitter("energy"),
    fSharingFilter("sharing"), 
    fDensityCalculator("density"),
    fList(0),
    fDebug(0)
{
  // 
  // Constructor 
  // 
  // Parameters:
  //    name Name of task 
  //
  DefineOutput(1, TList::Class());
  DefineOutput(2, TList::Class());
  fCorrManager = &AliForwardCorrectionManager::Instance(); 
  fEnergyFitter.SetNParticles(1); // Just find the 1st peak 
  fEnergyFitter.SetDoMakeObject(false); 
  fEnergyFitter.SetUseIncreasingBins(true);
  fEnergyFitter.SetDoFits(kTRUE);
  fEnergyFitter.SetLowCut(0.4);
  fEnergyFitter.SetFitRangeBinWidth(4);
  fEnergyFitter.SetMinEntries(1000);
}

//____________________________________________________________________
AliForwardQATask::AliForwardQATask(const AliForwardQATask& o)
  : AliAnalysisTaskSE(o),
    fEnableLowFlux(o.fEnableLowFlux), 
    fFirstEvent(o.fFirstEvent),
    fCorrManager(o.fCorrManager),
    fESDFMD(o.fESDFMD),
    fHistos(o.fHistos),
    fEventInspector(o.fEventInspector),
    fEnergyFitter(o.fEnergyFitter),
    fSharingFilter(o.fSharingFilter),
    fDensityCalculator(o.fDensityCalculator),
    fList(o.fList),
    fDebug(o.fDebug) 
{
  // 
  // Copy constructor 
  // 
  // Parameters:
  //    o Object to copy from 
  //
  DefineOutput(1, TList::Class());
  DefineOutput(2, TList::Class());
}

//____________________________________________________________________
AliForwardQATask&
AliForwardQATask::operator=(const AliForwardQATask& o)
{
  // 
  // Assignment operator 
  // 
  // Parameters:
  //    o Object to assign from 
  // 
  // Return:
  //    Reference to this object 
  //
  if (&o == this) return *this;
  AliAnalysisTaskSE::operator=(o);

  fEnableLowFlux     = o.fEnableLowFlux;
  fFirstEvent        = o.fFirstEvent;
  fCorrManager       = o.fCorrManager;
  fEventInspector    = o.fEventInspector;
  fEnergyFitter      = o.fEnergyFitter;
  fSharingFilter     = o.fSharingFilter;
  fDensityCalculator = o.fDensityCalculator;
  fHistos            = o.fHistos;
  fList              = o.fList;
  fDebug             = o.fDebug;

  return *this;
}

//____________________________________________________________________
void
AliForwardQATask::SetDebug(Int_t dbg)
{
  // 
  // Set debug level 
  // 
  // Parameters:
  //    dbg Debug level
  //
  fDebug = dbg;
  fEventInspector.SetDebug(dbg);
  fEnergyFitter.SetDebug(dbg);
  fSharingFilter.SetDebug(dbg);
  fDensityCalculator.SetDebug(dbg);
}

//____________________________________________________________________
Bool_t 
AliForwardQATask::CheckCorrections(UInt_t what)
{
  // 
  // Check if all needed corrections are there and accounted for.  If not,
  // do a Fatal exit 
  // 
  // Parameters:
  //    what Which corrections is needed
  // 
  // Return:
  //    true if all present, false otherwise
  //  

  AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  // Check that we have the energy loss fits, needed by 
  //   AliFMDSharingFilter 
  //   AliFMDDensityCalculator 
  if (what & AliForwardCorrectionManager::kELossFits && !fcm.GetELossFit()) { 
    AliWarning("No energy loss fits");

    // Fall-back values if we do not have the energy loss fits 
    AliFMDMultCuts& sfLCuts = GetSharingFilter().GetLCuts();
    if (sfLCuts.GetMethod() != AliFMDMultCuts::kFixed) { 
      Double_t cut = 0.3;
      AliWarningF("Using fixed cut @ %f for the lower bound "
		 "of the sharing filter", cut);
      sfLCuts.SetMultCuts(cut);
    }
    AliFMDMultCuts& sfHCuts = GetSharingFilter().GetHCuts();
    if (sfHCuts.GetMethod() != AliFMDMultCuts::kFixed) { 
      Double_t cut = 100;
      AliWarningF("Using fixed cut @ %f for the upper bound "
		 "of the sharing filter", cut);
      sfHCuts.SetMultCuts(cut);
    }
    AliFMDMultCuts& dcCuts  = GetDensityCalculator().GetCuts();
    if (dcCuts.GetMethod() != AliFMDMultCuts::kFixed) { 
      Double_t cut = 0.3;
      AliWarningF("Using fixed cut @ %f for the lower bound "
		 "of the density calculator", cut);
      dcCuts.SetMultCuts(cut);
    }
  }
  return true;
}

//____________________________________________________________________
Bool_t
AliForwardQATask::ReadCorrections(const TAxis*& pe, 
				  const TAxis*& pv, 
				  Bool_t        mc)
{
  //
  // Read corrections
  //
  //
  UInt_t what = AliForwardCorrectionManager::kAll;
  what ^= AliForwardCorrectionManager::kDoubleHit;
  what ^= AliForwardCorrectionManager::kVertexBias;
  what ^= AliForwardCorrectionManager::kAcceptance;
  what ^= AliForwardCorrectionManager::kMergingEfficiency;

  AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
  if (!fcm.Init(GetEventInspector().GetCollisionSystem(),
		GetEventInspector().GetEnergy(),
		GetEventInspector().GetField(),
		mc,
		what)) return false;
  if (!CheckCorrections(what)) {
    return false;
  }

  // Sett our persistency pointer 
  // fCorrManager = &fcm;

  // Get the eta axis from the secondary maps - if read in
  if (!pe) {
    pe = fcm.GetEtaAxis();
    if (!pe) AliFatal("No eta axis defined");
  }
  // Get the vertex axis from the secondary maps - if read in
  if (!pv) {
    pv = fcm.GetVertexAxis();
    if (!pv) AliFatal("No vertex axis defined");
  }

  return true;
}

//____________________________________________________________________
AliESDEvent*
AliForwardQATask::GetESDEvent()
{
  //
  // Get the ESD event. IF this is the first event, initialise
  //
  DGUARD(fDebug,2,"Get the ESD event");
  if (IsZombie()) {
    DMSG(fDebug,3,"We're a Zombie - bailing out");
    return 0;
  }
  AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
  if (!esd) {
    AliWarning("No ESD event found for input event");
    return 0;
  }

  // On the first event, initialize the parameters
  if (fFirstEvent && esd->GetESDRun()) {
    GetEventInspector().ReadRunDetails(esd);

    AliInfoF("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",
	     esd->GetBeamEnergy(),
	     esd->GetBeamType(),
	     esd->GetCurrentL3(),
	     esd->GetMagneticField(),
	     esd->GetRunNumber());


    if (!SetupForData()) {
      AliWarning("Initialisation of sub algorithms failed!");
      SetZombie(true);
      esd = 0;
      return 0;
    }
    AliInfoF("Clearing first event flag from %s to false", 
	     fFirstEvent ? "true" : "false");
    fFirstEvent = false;
  }
  return esd;
}
//____________________________________________________________________
Bool_t
AliForwardQATask::SetupForData()
{
  // 
  // Initialise the sub objects and stuff.  Called on first event 
  // 
  //
  const TAxis* pe = 0;
  const TAxis* pv = 0;

  if (!ReadCorrections(pe,pv))  { 
    AliWarning("Using default binning");
    pv = new TAxis(10,-10, 10);
    pe = new TAxis(240,-6,6);
  }

  fHistos.Init(*pe);

  fEventInspector.SetupForData(*pv);
  fEnergyFitter.SetupForData(*pe);
  fSharingFilter.SetupForData(*pe);
  fDensityCalculator.SetupForData(*pe);

  this->Print();

  return true;
}

//____________________________________________________________________
void
AliForwardQATask::UserCreateOutputObjects()
{
  // 
  // Create output objects 
  // 
  //
  fList = new TList;
  fList->SetOwner();
  
  fEventInspector.CreateOutputObjects(fList);
  fEnergyFitter.CreateOutputObjects(fList);
  fSharingFilter.CreateOutputObjects(fList);
  fDensityCalculator.CreateOutputObjects(fList);

  PostData(1, fList);
}
//____________________________________________________________________
void
AliForwardQATask::UserExec(Option_t*)
{
  // 
  // Process each event 
  // 
  // Parameters:
  //    option Not used
  //  
  DGUARD(fDebug,1,"Process the input event");

  // static Int_t cnt = 0;
  // cnt++;
  // Get the input data 
  AliESDEvent* esd = GetESDEvent();
  if (!esd) { 
    AliWarning("Got no ESD event");
    return;
  }
  if (fFirstEvent) { 
    // If the first event flag wasn't cleared in the above call to
    // GetESDEvent, we should not do anything, since nothing has been
    // initialised yet, so we opt out here (with a warning) 
    AliWarning("Nothing has been initialized yet, opt'ing out");
    return;
  }

  // Clear stuff 
  fHistos.Clear();
  fESDFMD.Clear();
  
  Bool_t   lowFlux   = kFALSE;
  UInt_t   triggers  = 0;
  UShort_t ivz       = 0;
  TVector3 ip;
  Double_t cent      = -1;
  UShort_t nClusters = 0;
  UInt_t   found     = fEventInspector.Process(esd, triggers, lowFlux, 
					       ivz, ip, cent, nClusters);
  
  Bool_t ok = true;
  if (found & AliFMDEventInspector::kNoEvent)    ok = false;
  if (found & AliFMDEventInspector::kNoTriggers) ok = false;
  if (found & AliFMDEventInspector::kNoSPD)      ok = false;
  if (found & AliFMDEventInspector::kNoFMD)      ok = false;
  if (found & AliFMDEventInspector::kNoVertex)   ok = false;
  if (triggers & AliAODForwardMult::kPileUp)     ok = false;
  if (found & AliFMDEventInspector::kBadVertex)  ok = false;
  if (!ok) { 
    DMSG(fDebug,2,"Event failed selection: %s", 
	 AliFMDEventInspector::CodeString(found));
    return;
  }
  DMSG(fDebug,2,"Event triggers: %s", 
       AliAODForwardMult::GetTriggerString(triggers));

  // 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();
  
  // Run the energy loss fitter 
  if (!fEnergyFitter.Accumulate(*esdFMD, cent, 
				triggers & AliAODForwardMult::kEmpty)) {
    AliWarning("Energy fitter failed");
    return;
  }
  
  //  // Apply the sharing filter (or hit merging or clustering if you like)
  if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD, ip.Z())) { 
    AliWarning("Sharing filter failed!");
    return;
  }
 
  // Calculate the inclusive charged particle density 
  if (!fDensityCalculator.Calculate(fESDFMD, fHistos, lowFlux, cent, ip)) { 
    // if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) { 
    AliWarning("Density calculator failed!");
    return;
  }
  
  PostData(1, fList);
}

//____________________________________________________________________
void
AliForwardQATask::Terminate(Option_t*)
{
  // 
  // End of job
  // 
  // Parameters:
  //    option Not used 
  //
  if (fDebug) AliInfo("In Forwards terminate");
  TStopwatch swt;
  swt.Start();

  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;
  }
  
  // Make a deep copy and post that as output 2 
  TList* list2 = static_cast<TList*>(list->Clone(Form("%sResults", 
						      list->GetName())));
  list2->SetOwner();

  // Get our histograms from the container 
  TH1I* hEventsTr    = 0;//static_cast<TH1I*>(list->FindObject("nEventsTr"));
  TH1I* hEventsTrVtx = 0;//static_cast<TH1I*>(list->FindObject("nEventsTrVtx"));
  TH1I* hEventsAcc   = 0;
  TH1I* hTriggers    = 0;
  if (!fEventInspector.FetchHistograms(list, 
				       hEventsTr, 
				       hEventsTrVtx, 
				       hEventsAcc,
				       hTriggers)) { 
    AliError(Form("Didn't get histograms from event selector "
		  "(hEventsTr=%p,hEventsTrVtx=%p,hEventsAcc=%p)", 
		  hEventsTr, hEventsTrVtx,hEventsAcc));
    return;
  }

  TStopwatch swf;
  swf.Start();
  fEnergyFitter.Fit(list2);
  swf.Stop();
  AliInfoF("Fitting took %d real-time seconds, and %f CPU seconds", 
	   Int_t(swf.RealTime()), swf.CpuTime());

  fSharingFilter.Terminate(list,list2,Int_t(hEventsTr->Integral()));
  fDensityCalculator.Terminate(list,list2,Int_t(hEventsTrVtx->Integral()));

  if (fDebug) AliInfoF("Posting post processing results to %s", 
		       list2->GetName());
  PostData(2, list2);

  swt.Stop();
  AliInfoF("Terminate took %d real-time seconds, and %f CPU seconds", 
	   Int_t(swt.RealTime()), swt.CpuTime());

  
}

//____________________________________________________________________
void
AliForwardQATask::Print(Option_t* option) const
{
  // 
  // Print information 
  // 
  // Parameters:
  //    option Not used
  //
  
  std::cout << ClassName() << ": " << GetName() << "\n" 
	    << "  Enable low flux code:   " << (fEnableLowFlux ? "yes" : "no") 
	    << "\n"
	    << "  Off-line trigger mask:  0x" 
	    << std::hex     << std::setfill('0') 
	    << std::setw (8) << fOfflineTriggerMask 
	    << std::dec     << std::setfill (' ') << std::endl;
  gROOT->IncreaseDirLevel();
  if (fCorrManager) fCorrManager->Print();
  else  
    std::cout << "  Correction manager not set yet" << std::endl;
  GetEventInspector()   .Print(option);
  GetEnergyFitter()     .Print(option);
  GetSharingFilter()    .Print(option);
  GetDensityCalculator().Print(option);
  gROOT->DecreaseDirLevel();
}

//
// EOF
//
Commit	Line	Data
96624385	1	//
	2	// Calculate the multiplicity in the forward regions event-by-event
	3	//
	4	// Inputs:
	5	// - AliESDEvent
	6	//
	7	// Outputs:
	8	// - AliAODForwardMult
	9	//
	10	// Histograms
	11	//
	12	// Corrections used
	13	//
	14	#include "AliForwardQATask.h"
	15	#include "AliForwardUtil.h"
	16	#include "AliTriggerAnalysis.h"
	17	#include "AliPhysicsSelection.h"
	18	#include "AliLog.h"
	19	#include "AliESDEvent.h"
	20	#include "AliAODHandler.h"
	21	#include "AliMultiplicity.h"
	22	#include "AliInputEventHandler.h"
	23	#include "AliForwardCorrectionManager.h"
	24	#include "AliAnalysisManager.h"
	25	#include "AliAODForwardMult.h"
	26	#include <TH1.h>
	27	#include <TDirectory.h>
	28	#include <TTree.h>
	29	#include <TROOT.h>
e2213ed5	30	#include <TStopwatch.h>
96624385	31
	32	//====================================================================
	33	AliForwardQATask::AliForwardQATask()
	34	: AliAnalysisTaskSE(),
	35	fEnableLowFlux(false),
	36	fFirstEvent(true),
	37	fCorrManager(0),
	38	fESDFMD(),
	39	fHistos(),
	40	fEventInspector(),
	41	fEnergyFitter(),
	42	fSharingFilter(),
	43	fDensityCalculator(),
9b2f2e39	44	fList(0),
9b2f2e39	45	fDebug(0)
96624385	46	{
	47	//
	48	// Constructor
	49	//
	50	}
	51
	52	//____________________________________________________________________
	53	AliForwardQATask::AliForwardQATask(const char* name)
	54	: AliAnalysisTaskSE(name),
	55	fEnableLowFlux(false),
	56	fFirstEvent(true),
	57	fCorrManager(0),
	58	fESDFMD(),
	59	fHistos(),
	60	fEventInspector("event"),
	61	fEnergyFitter("energy"),
	62	fSharingFilter("sharing"),
	63	fDensityCalculator("density"),
9b2f2e39	64	fList(0),
9b2f2e39	65	fDebug(0)
96624385	66	{
	67	//
	68	// Constructor
	69	//
	70	// Parameters:
	71	// name Name of task
	72	//
	73	DefineOutput(1, TList::Class());
	74	DefineOutput(2, TList::Class());
	75	fCorrManager = &AliForwardCorrectionManager::Instance();
	76	fEnergyFitter.SetNParticles(1); // Just find the 1st peak
	77	fEnergyFitter.SetDoMakeObject(false);
	78	fEnergyFitter.SetUseIncreasingBins(true);
	79	fEnergyFitter.SetDoFits(kTRUE);
	80	fEnergyFitter.SetLowCut(0.4);
	81	fEnergyFitter.SetFitRangeBinWidth(4);
	82	fEnergyFitter.SetMinEntries(1000);
	83	}
	84
	85	//____________________________________________________________________
	86	AliForwardQATask::AliForwardQATask(const AliForwardQATask& o)
	87	: AliAnalysisTaskSE(o),
	88	fEnableLowFlux(o.fEnableLowFlux),
	89	fFirstEvent(o.fFirstEvent),
	90	fCorrManager(o.fCorrManager),
	91	fESDFMD(o.fESDFMD),
	92	fHistos(o.fHistos),
	93	fEventInspector(o.fEventInspector),
	94	fEnergyFitter(o.fEnergyFitter),
	95	fSharingFilter(o.fSharingFilter),
	96	fDensityCalculator(o.fDensityCalculator),
9b2f2e39	97	fList(o.fList),
9b2f2e39	98	fDebug(o.fDebug)
96624385	99	{
	100	//
	101	// Copy constructor
	102	//
	103	// Parameters:
	104	// o Object to copy from
	105	//
	106	DefineOutput(1, TList::Class());
	107	DefineOutput(2, TList::Class());
	108	}
	109
	110	//____________________________________________________________________
	111	AliForwardQATask&
	112	AliForwardQATask::operator=(const AliForwardQATask& o)
	113	{
	114	//
	115	// Assignment operator
	116	//
	117	// Parameters:
	118	// o Object to assign from
	119	//
	120	// Return:
	121	// Reference to this object
	122	//
d015ecfe	123	if (&o == this) return *this;
96624385	124	AliAnalysisTaskSE::operator=(o);
	125
	126	fEnableLowFlux = o.fEnableLowFlux;
	127	fFirstEvent = o.fFirstEvent;
	128	fCorrManager = o.fCorrManager;
	129	fEventInspector = o.fEventInspector;
	130	fEnergyFitter = o.fEnergyFitter;
	131	fSharingFilter = o.fSharingFilter;
	132	fDensityCalculator = o.fDensityCalculator;
	133	fHistos = o.fHistos;
	134	fList = o.fList;
9b2f2e39	135	fDebug = o.fDebug;
96624385	136
	137	return *this;
	138	}
	139
	140	//____________________________________________________________________
	141	void
	142	AliForwardQATask::SetDebug(Int_t dbg)
	143	{
	144	//
	145	// Set debug level
	146	//
	147	// Parameters:
	148	// dbg Debug level
	149	//
9b2f2e39	150	fDebug = dbg;
96624385	151	fEventInspector.SetDebug(dbg);
	152	fEnergyFitter.SetDebug(dbg);
	153	fSharingFilter.SetDebug(dbg);
	154	fDensityCalculator.SetDebug(dbg);
	155	}
	156
	157	//____________________________________________________________________
	158	Bool_t
57522224	159	AliForwardQATask::CheckCorrections(UInt_t what)
96624385	160	{
	161	//
	162	// Check if all needed corrections are there and accounted for. If not,
	163	// do a Fatal exit
	164	//
	165	// Parameters:
	166	// what Which corrections is needed
	167	//
	168	// Return:
	169	// true if all present, false otherwise
	170	//
	171
	172	AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
	173	// Check that we have the energy loss fits, needed by
	174	// AliFMDSharingFilter
	175	// AliFMDDensityCalculator
	176	if (what & AliForwardCorrectionManager::kELossFits && !fcm.GetELossFit()) {
460a5c02	177	AliWarning("No energy loss fits");
57522224	178
	179	// Fall-back values if we do not have the energy loss fits
	180	AliFMDMultCuts& sfLCuts = GetSharingFilter().GetLCuts();
	181	if (sfLCuts.GetMethod() != AliFMDMultCuts::kFixed) {
	182	Double_t cut = 0.3;
	183	AliWarningF("Using fixed cut @ %f for the lower bound "
	184	"of the sharing filter", cut);
	185	sfLCuts.SetMultCuts(cut);
	186	}
	187	AliFMDMultCuts& sfHCuts = GetSharingFilter().GetHCuts();
	188	if (sfHCuts.GetMethod() != AliFMDMultCuts::kFixed) {
	189	Double_t cut = 100;
	190	AliWarningF("Using fixed cut @ %f for the upper bound "
	191	"of the sharing filter", cut);
	192	sfHCuts.SetMultCuts(cut);
	193	}
	194	AliFMDMultCuts& dcCuts = GetDensityCalculator().GetCuts();
	195	if (dcCuts.GetMethod() != AliFMDMultCuts::kFixed) {
	196	Double_t cut = 0.3;
	197	AliWarningF("Using fixed cut @ %f for the lower bound "
	198	"of the density calculator", cut);
	199	dcCuts.SetMultCuts(cut);
	200	}
96624385	201	}
	202	return true;
	203	}
	204
	205	//____________________________________________________________________
	206	Bool_t
	207	AliForwardQATask::ReadCorrections(const TAxis*& pe,
	208	const TAxis*& pv,
	209	Bool_t mc)
	210	{
	211	//
	212	// Read corrections
	213	//
	214	//
	215	UInt_t what = AliForwardCorrectionManager::kAll;
	216	what ^= AliForwardCorrectionManager::kDoubleHit;
	217	what ^= AliForwardCorrectionManager::kVertexBias;
	218	what ^= AliForwardCorrectionManager::kAcceptance;
	219	what ^= AliForwardCorrectionManager::kMergingEfficiency;
	220
	221	AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
	222	if (!fcm.Init(GetEventInspector().GetCollisionSystem(),
	223	GetEventInspector().GetEnergy(),
	224	GetEventInspector().GetField(),
	225	mc,
	226	what)) return false;
460a5c02	227	if (!CheckCorrections(what)) {
	228	return false;
	229	}
96624385	230
	231	// Sett our persistency pointer
	232	// fCorrManager = &fcm;
	233
	234	// Get the eta axis from the secondary maps - if read in
	235	if (!pe) {
	236	pe = fcm.GetEtaAxis();
	237	if (!pe) AliFatal("No eta axis defined");
	238	}
	239	// Get the vertex axis from the secondary maps - if read in
	240	if (!pv) {
	241	pv = fcm.GetVertexAxis();
	242	if (!pv) AliFatal("No vertex axis defined");
	243	}
	244
	245	return true;
	246	}
	247
	248	//____________________________________________________________________
	249	AliESDEvent*
	250	AliForwardQATask::GetESDEvent()
	251	{
	252	//
	253	// Get the ESD event. IF this is the first event, initialise
	254	//
241cca4d	255	DGUARD(fDebug,2,"Get the ESD event");
	256	if (IsZombie()) {
	257	DMSG(fDebug,3,"We're a Zombie - bailing out");
	258	return 0;
	259	}
96624385	260	AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
	261	if (!esd) {
	262	AliWarning("No ESD event found for input event");
	263	return 0;
	264	}
	265
	266	// On the first event, initialize the parameters
	267	if (fFirstEvent && esd->GetESDRun()) {
	268	GetEventInspector().ReadRunDetails(esd);
	269
422a78c8	270	AliInfoF("Initializing with parameters from the ESD:\n"
	271	" AliESDEvent::GetBeamEnergy() ->%f\n"
	272	" AliESDEvent::GetBeamType() ->%s\n"
	273	" AliESDEvent::GetCurrentL3() ->%f\n"
	274	" AliESDEvent::GetMagneticField()->%f\n"
6ff251d8	275	" AliESDEvent::GetRunNumber() ->%d",
422a78c8	276	esd->GetBeamEnergy(),
	277	esd->GetBeamType(),
	278	esd->GetCurrentL3(),
	279	esd->GetMagneticField(),
	280	esd->GetRunNumber());
96624385	281
96624385	282
5934a3e3	283	if (!SetupForData()) {
422a78c8	284	AliWarning("Initialisation of sub algorithms failed!");
6ff251d8	285	SetZombie(true);
6ff251d8	286	esd = 0;
422a78c8	287	return 0;
422a78c8	288	}
460a5c02	289	AliInfoF("Clearing first event flag from %s to false",
	290	fFirstEvent ? "true" : "false");
	291	fFirstEvent = false;
96624385	292	}
	293	return esd;
	294	}
	295	//____________________________________________________________________
422a78c8	296	Bool_t
5934a3e3	297	AliForwardQATask::SetupForData()
96624385	298	{
	299	//
	300	// Initialise the sub objects and stuff. Called on first event
	301	//
	302	//
	303	const TAxis* pe = 0;
	304	const TAxis* pv = 0;
	305
d70d2a44	306	if (!ReadCorrections(pe,pv)) {
	307	AliWarning("Using default binning");
	308	pv = new TAxis(10,-10, 10);
	309	pe = new TAxis(240,-6,6);
460a5c02	310	}
96624385	311
	312	fHistos.Init(*pe);
	313
5934a3e3	314	fEventInspector.SetupForData(*pv);
	315	fEnergyFitter.SetupForData(*pe);
	316	fSharingFilter.SetupForData(*pe);
	317	fDensityCalculator.SetupForData(*pe);
96624385	318
96624385	319	this->Print();
422a78c8	320
422a78c8	321	return true;
96624385	322	}
	323
	324	//____________________________________________________________________
	325	void
	326	AliForwardQATask::UserCreateOutputObjects()
	327	{
	328	//
	329	// Create output objects
	330	//
	331	//
	332	fList = new TList;
	333	fList->SetOwner();
	334
5934a3e3	335	fEventInspector.CreateOutputObjects(fList);
	336	fEnergyFitter.CreateOutputObjects(fList);
	337	fSharingFilter.CreateOutputObjects(fList);
	338	fDensityCalculator.CreateOutputObjects(fList);
96624385	339
	340	PostData(1, fList);
	341	}
	342	//____________________________________________________________________
	343	void
	344	AliForwardQATask::UserExec(Option_t*)
	345	{
	346	//
	347	// Process each event
	348	//
	349	// Parameters:
	350	// option Not used
	351	//
241cca4d	352	DGUARD(fDebug,1,"Process the input event");
96624385	353
	354	// static Int_t cnt = 0;
	355	// cnt++;
	356	// Get the input data
	357	AliESDEvent* esd = GetESDEvent();
422a78c8	358	if (!esd) {
	359	AliWarning("Got no ESD event");
	360	return;
	361	}
460a5c02	362	if (fFirstEvent) {
	363	// If the first event flag wasn't cleared in the above call to
	364	// GetESDEvent, we should not do anything, since nothing has been
	365	// initialised yet, so we opt out here (with a warning)
	366	AliWarning("Nothing has been initialized yet, opt'ing out");
	367	return;
	368	}
96624385	369
	370	// Clear stuff
	371	fHistos.Clear();
	372	fESDFMD.Clear();
	373
	374	Bool_t lowFlux = kFALSE;
	375	UInt_t triggers = 0;
	376	UShort_t ivz = 0;
5ca83fee	377	TVector3 ip;
96624385	378	Double_t cent = -1;
	379	UShort_t nClusters = 0;
	380	UInt_t found = fEventInspector.Process(esd, triggers, lowFlux,
5ca83fee	381	ivz, ip, cent, nClusters);
96624385	382
241cca4d	383	Bool_t ok = true;
	384	if (found & AliFMDEventInspector::kNoEvent) ok = false;
	385	if (found & AliFMDEventInspector::kNoTriggers) ok = false;
	386	if (found & AliFMDEventInspector::kNoSPD) ok = false;
	387	if (found & AliFMDEventInspector::kNoFMD) ok = false;
	388	if (found & AliFMDEventInspector::kNoVertex) ok = false;
	389	if (triggers & AliAODForwardMult::kPileUp) ok = false;
	390	if (found & AliFMDEventInspector::kBadVertex) ok = false;
	391	if (!ok) {
	392	DMSG(fDebug,2,"Event failed selection: %s",
	393	AliFMDEventInspector::CodeString(found));
	394	return;
	395	}
5ca83fee	396	DMSG(fDebug,2,"Event triggers: %s",
5ca83fee	397	AliAODForwardMult::GetTriggerString(triggers));
96624385	398
	399	// We we do not want to use low flux specific code, we disable it here.
	400	if (!fEnableLowFlux) lowFlux = false;
	401
	402	// Get FMD data
	403	AliESDFMD* esdFMD = esd->GetFMDData();
	404
	405	// Run the energy loss fitter
	406	if (!fEnergyFitter.Accumulate(*esdFMD, cent,
	407	triggers & AliAODForwardMult::kEmpty)) {
	408	AliWarning("Energy fitter failed");
	409	return;
	410	}
4077f3e8	411
96624385	412	// // Apply the sharing filter (or hit merging or clustering if you like)
5ca83fee	413	if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD, ip.Z())) {
96624385	414	AliWarning("Sharing filter failed!");
	415	return;
	416	}
4077f3e8	417
96624385	418	// Calculate the inclusive charged particle density
57522224	419	if (!fDensityCalculator.Calculate(fESDFMD, fHistos, lowFlux, cent, ip)) {
96624385	420	// if (!fDensityCalculator.Calculate(*esdFMD, fHistos, ivz, lowFlux)) {
	421	AliWarning("Density calculator failed!");
	422	return;
	423	}
4077f3e8	424
96624385	425	PostData(1, fList);
	426	}
	427
	428	//____________________________________________________________________
	429	void
	430	AliForwardQATask::Terminate(Option_t*)
	431	{
	432	//
	433	// End of job
	434	//
	435	// Parameters:
	436	// option Not used
	437	//
9b2f2e39	438	if (fDebug) AliInfo("In Forwards terminate");
e2213ed5	439	TStopwatch swt;
e2213ed5	440	swt.Start();
96624385	441
	442	TList* list = dynamic_cast<TList*>(GetOutputData(1));
	443	if (!list) {
	444	AliError(Form("No output list defined (%p)", GetOutputData(1)));
	445	if (GetOutputData(1)) GetOutputData(1)->Print();
	446	return;
	447	}
	448
5934a3e3	449	// Make a deep copy and post that as output 2
	450	TList* list2 = static_cast<TList*>(list->Clone(Form("%sResults",
	451	list->GetName())));
	452	list2->SetOwner();
	453
96624385	454	// Get our histograms from the container
	455	TH1I* hEventsTr = 0;//static_cast<TH1I*>(list->FindObject("nEventsTr"));
	456	TH1I* hEventsTrVtx = 0;//static_cast<TH1I*>(list->FindObject("nEventsTrVtx"));
5ca83fee	457	TH1I* hEventsAcc = 0;
96624385	458	TH1I* hTriggers = 0;
5ca83fee	459	if (!fEventInspector.FetchHistograms(list,
	460	hEventsTr,
	461	hEventsTrVtx,
	462	hEventsAcc,
	463	hTriggers)) {
96624385	464	AliError(Form("Didn't get histograms from event selector "
5ca83fee	465	"(hEventsTr=%p,hEventsTrVtx=%p,hEventsAcc=%p)",
5ca83fee	466	hEventsTr, hEventsTrVtx,hEventsAcc));
96624385	467	return;
	468	}
	469
e2213ed5	470	TStopwatch swf;
e2213ed5	471	swf.Start();
57522224	472	fEnergyFitter.Fit(list2);
e2213ed5	473	swf.Stop();
	474	AliInfoF("Fitting took %d real-time seconds, and %f CPU seconds",
	475	Int_t(swf.RealTime()), swf.CpuTime());
96624385	476
5934a3e3	477	fSharingFilter.Terminate(list,list2,Int_t(hEventsTr->Integral()));
5934a3e3	478	fDensityCalculator.Terminate(list,list2,Int_t(hEventsTrVtx->Integral()));
96624385	479
9b2f2e39	480	if (fDebug) AliInfoF("Posting post processing results to %s",
9b2f2e39	481	list2->GetName());
96624385	482	PostData(2, list2);
96624385	483
e2213ed5	484	swt.Stop();
	485	AliInfoF("Terminate took %d real-time seconds, and %f CPU seconds",
	486	Int_t(swt.RealTime()), swt.CpuTime());
	487
5934a3e3	488
96624385	489	}
	490
	491	//____________________________________________________________________
	492	void
	493	AliForwardQATask::Print(Option_t* option) const
	494	{
	495	//
	496	// Print information
	497	//
	498	// Parameters:
	499	// option Not used
	500	//
	501
	502	std::cout << ClassName() << ": " << GetName() << "\n"
	503	<< " Enable low flux code: " << (fEnableLowFlux ? "yes" : "no")
	504	<< "\n"
	505	<< " Off-line trigger mask: 0x"
	506	<< std::hex << std::setfill('0')
	507	<< std::setw (8) << fOfflineTriggerMask
	508	<< std::dec << std::setfill (' ') << std::endl;
	509	gROOT->IncreaseDirLevel();
	510	if (fCorrManager) fCorrManager->Print();
	511	else
	512	std::cout << " Correction manager not set yet" << std::endl;
	513	GetEventInspector() .Print(option);
	514	GetEnergyFitter() .Print(option);
	515	GetSharingFilter() .Print(option);
57522224	516	GetDensityCalculator().Print(option);
96624385	517	gROOT->DecreaseDirLevel();
	518	}
	519
	520	//
	521	// EOF
	522	//