PWGLF/FORWARD/analysis2/AliFMDEnergyFitterTask.cxx

   1 //
   2 // Histogram and fit the energy loss distributions for the FMD
   3 //
   4 // Inputs:
   5 //   - AliESDEvent
   6 //
   7 // Outputs:
   8 //   - None
   9 //
  10 // Histograms:
  11 //
  12 // Corrections used:
  13 //   - None
  14 //
  15 //
  16 //
  17 #include "AliFMDEnergyFitterTask.h"
  18 #include "AliLog.h"
  19 #include "AliESDEvent.h"
  20 #include "AliAODForwardMult.h"
  21 #include "AliForwardCorrectionManager.h"
  22 #include "AliAnalysisManager.h"
  23 #include "AliAnalysisDataSlot.h"
  24 #include "AliAnalysisDataContainer.h"
  25 #include <TH1.h>
  26 #include <TDirectory.h>
  27 #include <TTree.h>
  28 #include <TFile.h>
  29 #include "AliMCEvent.h"
  30 #include "AliGenHijingEventHeader.h"
  31 #include "AliHeader.h"
  32 #include <iostream>
  33
  34 //====================================================================
  35 AliFMDEnergyFitterTask::AliFMDEnergyFitterTask()
  36   : AliAnalysisTaskSE(),
  37     fFirstEvent(true),
  38     fEventInspector(),
  39     fEnergyFitter(),
  40     fList(0),
  41     fbLow(0),
  42     fbHigh(100)
  43 {
  44   //
  45   // Constructor
  46   //
  47   DGUARD(fDebug, 3,"Default CTOR of AliFMDEnergyFitterTask");
  48 }
  49
  50 //____________________________________________________________________
  51 AliFMDEnergyFitterTask::AliFMDEnergyFitterTask(const char* name)
  52   : AliAnalysisTaskSE(name),
  53     fFirstEvent(true),
  54     fEventInspector("event"),
  55     fEnergyFitter("energy"),
  56     fList(0),
  57     fbLow(0),
  58     fbHigh(100)
  59 {
  60   //
  61   // Constructor
  62   //
  63   // Parameters:
  64   //    name Name of task
  65   //
  66   DGUARD(fDebug, 3,"Named CTOR of AliFMDEnergyFitterTask: %s", name);
  67   DefineOutput(1, TList::Class());
  68   DefineOutput(2, TList::Class());
  69   fBranchNames =
  70     "ESD:AliESDRun.,AliESDHeader.,AliMultiplicity.,"
  71     "AliESDFMD.,SPDVertex.,PrimaryVertex.";
  72 }
  73
  74 //____________________________________________________________________
  75 AliFMDEnergyFitterTask::AliFMDEnergyFitterTask(const AliFMDEnergyFitterTask& o)
  76   : AliAnalysisTaskSE(o),
  77     fFirstEvent(o.fFirstEvent),
  78     fEventInspector(o.fEventInspector),
  79     fEnergyFitter(o.fEnergyFitter),
  80     fList(o.fList),
  81     fbLow(o.fbLow),
  82     fbHigh(o.fbHigh)
  83 {
  84   //
  85   // Copy constructor
  86   //
  87   // Parameters:
  88   //    o Object to copy from
  89   //
  90   DGUARD(fDebug, 3,"COPY CTOR of AliFMDEnergyFitterTask");
  91   DefineOutput(1, TList::Class());
  92   DefineOutput(2, TList::Class());
  93 }
  94
  95 //____________________________________________________________________
  96 AliFMDEnergyFitterTask&
  97 AliFMDEnergyFitterTask::operator=(const AliFMDEnergyFitterTask& o)
  98 {
  99   //
 100   // Assignment operator
 101   //
 102   // Parameters:
 103   //    o Object to assign from
 104   //
 105   // Return:
 106   //    Reference to this object
 107   //
 108   fDebug = o.fDebug;
 109   DGUARD(fDebug,3,"Assignment of AliFMDEnergyFitterTask");
 110   if (&o == this) return *this;
 111   AliAnalysisTaskSE::operator=(o);
 112
 113   fFirstEvent        = o.fFirstEvent;
 114   fEventInspector    = o.fEventInspector;
 115   fEnergyFitter      = o.fEnergyFitter;
 116   fList              = o.fList;
 117   fbLow              = o.fbLow;
 118   fbHigh             = o.fbHigh;
 119
 120   return *this;
 121 }
 122
 123 //____________________________________________________________________
 124 void
 125 AliFMDEnergyFitterTask::SetDebug(Int_t dbg)
 126 {
 127   //
 128   // Set the debug level
 129   //
 130   // Parameters:
 131   //    dbg Debug level
 132   //
 133   fDebug = dbg;
 134   fEventInspector.SetDebug(dbg);
 135   fEnergyFitter.SetDebug(dbg);
 136 }
 137
 138 //____________________________________________________________________
 139 void
 140 AliFMDEnergyFitterTask::Init()
 141 {
 142   //
 143   // Initialize the task
 144   //
 145   //
 146   DGUARD(fDebug,1,"Initialize of AliFMDEnergyFitterTask");
 147   fFirstEvent = true;
 148 }
 149
 150 //____________________________________________________________________
 151 void
 152 AliFMDEnergyFitterTask::SetupForData()
 153 {
 154   //
 155   // Initialise the sub objects and stuff.  Called on first event
 156   //
 157   //
 158   DGUARD(fDebug,1,"Initialize subs of AliFMDEnergyFitterTask");
 159   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
 160   UShort_t sys = fEventInspector.GetCollisionSystem();
 161   UShort_t sNN = fEventInspector.GetEnergy();
 162   Short_t  fld = fEventInspector.GetField();
 163   fcm.Init(sys, sNN, fld, 0);
 164   TAxis eAxis(0,0,0);
 165   TAxis vAxis(10,-10,10);
 166   fEnergyFitter.SetupForData(eAxis);
 167   fEventInspector.SetupForData(vAxis);
 168
 169 }
 170
 171 //____________________________________________________________________
 172 void
 173 AliFMDEnergyFitterTask::UserCreateOutputObjects()
 174 {
 175   //
 176   // Create output objects
 177   //
 178   //
 179   DGUARD(fDebug,1,"Create output objects of AliFMDEnergyFitterTask");
 180   fList = new TList;
 181   fList->SetOwner();
 182
 183   fEventInspector.CreateOutputObjects(fList);
 184   fEnergyFitter.CreateOutputObjects(fList);
 185
 186   PostData(1, fList);
 187 }
 188 //____________________________________________________________________
 189 void
 190 AliFMDEnergyFitterTask::UserExec(Option_t*)
 191 {
 192   //
 193   // Process each event
 194   //
 195   // Parameters:
 196   //    option Not used
 197   //
 198
 199   // static Int_t cnt = 0;
 200   // cnt++;
 201   // Get the input data
 202   DGUARD(fDebug,3,"Analyse event of AliFMDEnergyFitterTask");
 203
 204   AliMCEvent* mcevent = MCEvent();
 205   if(mcevent) {
 206     AliHeader* header            = mcevent->Header();
 207     AliGenHijingEventHeader* hijingHeader =
 208       dynamic_cast<AliGenHijingEventHeader*>(header->GenEventHeader());
 209     if(hijingHeader) {
 210       Float_t b = hijingHeader->ImpactParameter();
 211       if(b<fbLow || b>fbHigh) return;
 212       else
 213         std::cout<<"Selecting event with impact parameter "<<b<<std::endl;
 214     }
 215
 216   }
 217
 218   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
 219   // AliInfo(Form("Event # %6d (esd=%p)", cnt, esd));
 220   if (!esd) {
 221     AliWarning("No ESD event found for input event");
 222     return;
 223   }
 224
 225   // --- Read in the data --------------------------------------------
 226   LoadBranches();
 227
 228   // On the first event, initialize the parameters
 229   if (fFirstEvent && esd->GetESDRun()) {
 230     fEventInspector.ReadRunDetails(esd);
 231
 232     AliInfo(Form("Initializing with parameters from the ESD:\n"
 233                  "         AliESDEvent::GetBeamEnergy()   ->%f\n"
 234                  "         AliESDEvent::GetBeamType()     ->%s\n"
 235                  "         AliESDEvent::GetCurrentL3()    ->%f\n"
 236                  "         AliESDEvent::GetMagneticField()->%f\n"
 237                  "         AliESDEvent::GetRunNumber()    ->%d\n",
 238                  esd->GetBeamEnergy(),
 239                  esd->GetBeamType(),
 240                  esd->GetCurrentL3(),
 241                  esd->GetMagneticField(),
 242                  esd->GetRunNumber()));
 243
 244
 245
 246     // AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
 247     // pars->SetParametersFromESD(esd);
 248     // pars->PrintStatus();
 249     fFirstEvent = false;
 250
 251     SetupForData();
 252   }
 253   Bool_t   lowFlux   = kFALSE;
 254   UInt_t   triggers  = 0;
 255   UShort_t ivz       = 0;
 256   TVector3 ip;
 257   Double_t cent      = 0;
 258   UShort_t nClusters = 0;
 259   UInt_t   found     = fEventInspector.Process(esd, triggers, lowFlux,
 260                                                ivz, ip, cent, nClusters);
 261   if (found & AliFMDEventInspector::kNoEvent)    return;
 262   if (found & AliFMDEventInspector::kNoTriggers) return;
 263   if (found & AliFMDEventInspector::kNoSPD)     return;
 264   if (found & AliFMDEventInspector::kNoFMD)     return;
 265   if (found & AliFMDEventInspector::kNoVertex)  return;
 266   if (found & AliFMDEventInspector::kBadVertex) return;
 267
 268   //  if(cent > 0) {
 269   //  if( cent < 40 || cent >50 ) return;
 270   //  else std::cout<<"selecting event with cent "<<cent<<std::endl;
 271   // }
 272
 273   // Get FMD data
 274   AliESDFMD* esdFMD = esd->GetFMDData();
 275   // Do the energy stuff
 276   if (!fEnergyFitter.Accumulate(*esdFMD, cent,
 277                                 triggers & AliAODForwardMult::kEmpty)){
 278     AliWarning("Energy fitter failed");
 279     return;
 280   }
 281   PostData(1, fList);
 282 }
 283
 284 //____________________________________________________________________
 285 void
 286 AliFMDEnergyFitterTask::FinishTaskOutput()
 287 {
 288   if (!fList) {
 289     Warning("FinishTaskOutput", "No list defined");
 290   }
 291   // else {
 292   //   fList->ls();
 293   // }
 294   // if (fDebug)
 295   // gDebug = 1;
 296   AliAnalysisTaskSE::FinishTaskOutput();
 297 }
 298
 299 //____________________________________________________________________
 300 void
 301 AliFMDEnergyFitterTask::Terminate(Option_t*)
 302 {
 303   //
 304   // End of job
 305   //
 306   // Parameters:
 307   //    option Not used
 308   //
 309   DGUARD(fDebug,1,"Processing merged output of AliFMDEnergyFitterTask");
 310
 311   TList* list = dynamic_cast<TList*>(GetOutputData(1));
 312   if (!list) {
 313     AliError(Form("No output list defined (%p)", GetOutputData(1)));
 314     if (GetOutputData(1)) GetOutputData(1)->Print();
 315     return;
 316   }
 317
 318   AliInfo("Fitting energy loss spectra");
 319   fEnergyFitter.Fit(list);
 320
 321   // Make a deep copy and post that as output 2
 322   TList* list2 = static_cast<TList*>(list->Clone(Form("%sResults",
 323                                                       list->GetName())));
 324   list2->SetOwner();
 325   PostData(2, list2);
 326 }
 327
 328 //____________________________________________________________________
 329 void
 330 AliFMDEnergyFitterTask::Print(Option_t*) const
 331 {
 332   //
 333   // Print information
 334   //
 335   // Parameters:
 336   //    option Not used
 337   //
 338 }
 339
 340 //
 341 // EOF
 342 //