PWGLF/FORWARD/analysis2/AliFMDEventInspector.cxx

   1 //
   2 // This class inspects the event
   3 //
   4 // Input:
   5 //   - AliESDFMD object possibly corrected for sharing
   6 //
   7 // Output:
   8 //   - A histogram of v_z of events with triggers.
   9 //   - A histogram of v_z of events with vertex and triggers
  10 //   - A histogram of trigger counters
  11 //
  12 // Note, that these are added to the master output list
  13 //
  14 // Corrections used:
  15 //   - None
  16 //
  17 #include "AliFMDEventInspector.h"
  18 #include "AliLog.h"
  19 #include "AliESDEvent.h"
  20 #include "AliMultiplicity.h"
  21 #include "AliAnalysisManager.h"
  22 #include "AliMCEventHandler.h"
  23 #include "AliInputEventHandler.h"
  24 #include "AliTriggerAnalysis.h"
  25 #include "AliPhysicsSelection.h"
  26 #include "AliOADBPhysicsSelection.h"
  27 #include "AliAODForwardMult.h"
  28 #include "AliForwardUtil.h"
  29 #include "AliCentrality.h"
  30 #include <TH1.h>
  31 #include <TList.h>
  32 #include <TDirectory.h>
  33 #include <TROOT.h>
  34 #include <TParameter.h>
  35 #include <iostream>
  36 #include <iomanip>
  37 #include "AliMCEvent.h"
  38 #include "AliHeader.h"
  39 #include "AliGenEventHeader.h"
  40 #include "AliCollisionGeometry.h"
  41 #include "AliVVZERO.h"
  42
  43 //====================================================================
  44 const char* AliFMDEventInspector::fgkFolderName = "fmdEventInspector";
  45
  46 //____________________________________________________________________
  47 AliFMDEventInspector::AliFMDEventInspector()
  48   : TNamed(),
  49     fHEventsTr(0),
  50     fHEventsTrVtx(0),
  51     fHEventsAccepted(0),
  52     fHEventsAcceptedXY(0),
  53     fHTriggers(0),
  54     fHTriggerCorr(0),
  55     fHType(0),
  56     fHWords(0),
  57     fHCent(0),
  58     fHCentVsQual(0),
  59     fHStatus(0),
  60     fHVtxStatus(0),
  61     fHTrgStatus(0),
  62     fLowFluxCut(1000),
  63     fMaxVzErr(0.2),
  64     fList(0),
  65     fEnergy(0),
  66     fField(999),
  67     fCollisionSystem(kUnknown),
  68     fDebug(0),
  69     fCentAxis(0),
  70     fVtxAxis(10,-10,10),
  71     fUseFirstPhysicsVertex(false),
  72     fUseV0AND(false),
  73     fMinPileupContrib(3),
  74   fMinPileupDistance(0.8),
  75   fUseDisplacedVertices(false),
  76   fDisplacedVertex(),
  77   fCollWords(),
  78   fBgWords(),
  79   fCentMethod("V0M"),
  80   fMinCent(-1.0),
  81   fMaxCent(-1.0),
  82   fUsepA2012Vertex(false),
  83   fRunNumber(0)
  84 {
  85   //
  86   // Constructor
  87   //
  88   DGUARD(fDebug,1,"Default CTOR of AliFMDEventInspector");
  89 }
  90
  91 //____________________________________________________________________
  92 AliFMDEventInspector::AliFMDEventInspector(const char* name)
  93   : TNamed(fgkFolderName, name),
  94     fHEventsTr(0),
  95     fHEventsTrVtx(0),
  96     fHEventsAccepted(0),
  97     fHEventsAcceptedXY(0),
  98     fHTriggers(0),
  99     fHTriggerCorr(0),
 100     fHType(0),
 101     fHWords(0),
 102     fHCent(0),
 103     fHCentVsQual(0),
 104     fHStatus(0),
 105     fHVtxStatus(0),
 106     fHTrgStatus(0),
 107     fLowFluxCut(1000),
 108     fMaxVzErr(0.2),
 109     fList(0),
 110     fEnergy(0),
 111     fField(999),
 112     fCollisionSystem(kUnknown),
 113     fDebug(0),
 114     fCentAxis(0),
 115     fVtxAxis(10,-10,10),
 116     fUseFirstPhysicsVertex(false),
 117     fUseV0AND(false),
 118     fMinPileupContrib(3),
 119     fMinPileupDistance(0.8),
 120     fUseDisplacedVertices(false),
 121   fDisplacedVertex(),
 122   fCollWords(),
 123   fBgWords(),
 124   fCentMethod("V0M"),
 125   fMinCent(-1.0),
 126   fMaxCent(-1.0),
 127   fUsepA2012Vertex(false),
 128   fRunNumber(0)
 129 {
 130   //
 131   // Constructor
 132   //
 133   // Parameters:
 134   //   name Name of object
 135   //
 136   DGUARD(fDebug,1,"Named CTOR of AliFMDEventInspector: %s", name);
 137 }
 138
 139 //____________________________________________________________________
 140 AliFMDEventInspector::AliFMDEventInspector(const AliFMDEventInspector& o)
 141   : TNamed(o),
 142     fHEventsTr(o.fHEventsTr),
 143     fHEventsTrVtx(o.fHEventsTrVtx),
 144     fHEventsAccepted(o.fHEventsAccepted),
 145     fHEventsAcceptedXY(o.fHEventsAcceptedXY),
 146     fHTriggers(o.fHTriggers),
 147     fHTriggerCorr(o.fHTriggerCorr),
 148     fHType(o.fHType),
 149     fHWords(o.fHWords),
 150     fHCent(o.fHCent),
 151     fHCentVsQual(o.fHCentVsQual),
 152     fHStatus(o.fHStatus),
 153     fHVtxStatus(o.fHVtxStatus),
 154     fHTrgStatus(o.fHTrgStatus),
 155     fLowFluxCut(o.fLowFluxCut),
 156     fMaxVzErr(o.fMaxVzErr),
 157     fList(o.fList),
 158     fEnergy(o.fEnergy),
 159     fField(o.fField),
 160     fCollisionSystem(o.fCollisionSystem),
 161     fDebug(0),
 162     fCentAxis(0),
 163     fVtxAxis(o.fVtxAxis),
 164     fUseFirstPhysicsVertex(o.fUseFirstPhysicsVertex),
 165     fUseV0AND(o.fUseV0AND),
 166     fMinPileupContrib(o.fMinPileupContrib),
 167     fMinPileupDistance(o.fMinPileupDistance),
 168     fUseDisplacedVertices(o.fUseDisplacedVertices),
 169     fDisplacedVertex(o.fDisplacedVertex),
 170   fCollWords(),
 171   fBgWords(),
 172   fCentMethod(o.fCentMethod),
 173   fMinCent(o.fMinCent),
 174   fMaxCent(o.fMaxCent),
 175   fUsepA2012Vertex(o.fUsepA2012Vertex),
 176   fRunNumber(o.fRunNumber)
 177
 178 {
 179   //
 180   // Copy constructor
 181   //
 182   // Parameters:
 183   //   o Object to copy from
 184   //
 185   DGUARD(fDebug,1,"Copy CTOR of AliFMDEventInspector");
 186 }
 187
 188 //____________________________________________________________________
 189 AliFMDEventInspector::~AliFMDEventInspector()
 190 {
 191   //
 192   // Destructor
 193   //
 194   DGUARD(fDebug,1,"DTOR of AliFMDEventInspector");
 195   // if (fList)         delete fList;
 196 }
 197 //____________________________________________________________________
 198 AliFMDEventInspector&
 199 AliFMDEventInspector::operator=(const AliFMDEventInspector& o)
 200 {
 201   //
 202   // Assignement operator
 203   //
 204   // Parameters:
 205   //   o Object to assign from
 206   //
 207   // Return:
 208   //    Reference to this object
 209   //
 210   DGUARD(fDebug,3,"Assignment of AliFMDEventInspector");
 211   if (&o == this) return *this;
 212   TNamed::operator=(o);
 213   fHEventsTr         = o.fHEventsTr;
 214   fHEventsTrVtx      = o.fHEventsTrVtx;
 215   fHEventsAccepted   = o.fHEventsAccepted;
 216   fHEventsAcceptedXY = o.fHEventsAcceptedXY;
 217   fHTriggers         = o.fHTriggers;
 218   fHTriggerCorr      = o.fHTriggerCorr;
 219   fHType             = o.fHType;
 220   fHWords            = o.fHWords;
 221   fHCent             = o.fHCent;
 222   fHCentVsQual       = o.fHCentVsQual;
 223   fHStatus           = o.fHStatus;
 224   fHVtxStatus        = o.fHVtxStatus;
 225   fHTrgStatus        = o.fHTrgStatus;
 226   fLowFluxCut        = o.fLowFluxCut;
 227   fMaxVzErr          = o.fMaxVzErr;
 228   fDebug             = o.fDebug;
 229   fList              = (o.fList ? new TList : 0);
 230   fEnergy            = o.fEnergy;
 231   fField             = o.fField;
 232   fCollisionSystem   = o.fCollisionSystem;
 233   fVtxAxis.Set(o.fVtxAxis.GetNbins(), o.fVtxAxis.GetXmin(),
 234                o.fVtxAxis.GetXmax());
 235
 236   fUseFirstPhysicsVertex = o.fUseFirstPhysicsVertex;
 237   fUseV0AND              = o.fUseV0AND;
 238   fMinPileupContrib      = o.fMinPileupContrib;
 239   fMinPileupDistance     = o.fMinPileupDistance;
 240   fUseDisplacedVertices  = o.fUseDisplacedVertices;
 241   fDisplacedVertex       = o.fDisplacedVertex;
 242   fCentMethod            = o.fCentMethod;
 243   fMinCent               = o.fMinCent;
 244   fMaxCent               = o.fMaxCent;
 245   fUsepA2012Vertex       = o.fUsepA2012Vertex;
 246   fRunNumber             = o.fRunNumber;
 247
 248   if (fList) {
 249     fList->SetName(GetName());
 250     if (fHEventsTr)    fList->Add(fHEventsTr);
 251     if (fHEventsTrVtx) fList->Add(fHEventsTrVtx);
 252     if (fHTriggers)    fList->Add(fHTriggers);
 253     if (fHTriggerCorr) fList->Add(fHTriggerCorr);
 254     if (fHType)        fList->Add(fHType);
 255     if (fHWords)       fList->Add(fHWords);
 256     if (fHCent)        fList->Add(fHCent);
 257     if (fHCentVsQual)  fList->Add(fHCentVsQual);
 258     if (fHStatus)      fList->Add(fHStatus);
 259     if (fHVtxStatus)   fList->Add(fHVtxStatus);
 260     if (fHTrgStatus)   fList->Add(fHTrgStatus);
 261   }
 262   return *this;
 263 }
 264
 265 //____________________________________________________________________
 266 void
 267 AliFMDEventInspector::SetCentralityMethod(ECentMethod m)
 268 {
 269   switch (m) {
 270   case kV0Multiplicity: fCentMethod = "VOM"; break; // VZERO multiplicity
 271   case kV0Amplitude:    fCentMethod = "V0A"; break; // VZERO amplitude
 272   case kV0Charge:       fCentMethod = "V0C"; break; // VZERO charge
 273   case kFMDRough:       fCentMethod = "FMD"; break; // FMD scaled energy l
 274   case kNTracks:        fCentMethod = "TRK"; break; // Number of tracks
 275   case kLTracks:        fCentMethod = "TKL"; break; // Number of tracks
 276   case kCL0:            fCentMethod = "CL0"; break; //
 277   case kCL1:            fCentMethod = "CL1"; break; //
 278   case kCND:            fCentMethod = "CND"; break; //
 279   case kNParticles:     fCentMethod = "NPA"; break; // Neutral particles
 280   case kNeutrons:       fCentMethod = "ZNA"; break; // ZDC neutron amplitu
 281   case kV0vsFMD:        fCentMethod = "V0MvsFMD"; break; // VZERO versus FMD
 282   case kV0vsNTracks:    fCentMethod = "TKLvsVOM"; break; // Tracks versus VZERO
 283   case kZEMvsZDC:       fCentMethod = "ZEMvsZDC"; break; // ZDC
 284   default:              fCentMethod = "V0M"; break;
 285   }
 286 }
 287
 288 //____________________________________________________________________
 289 void
 290 AliFMDEventInspector::SetMinCentrality(Double_t minCent)
 291 {
 292   AliWarning("\n"
 293              "*******************************************************\n"
 294              "* Setting centrality cuts in this stage is deprecated *\n"
 295              "*******************************************************");
 296   fMinCent = minCent;
 297 }
 298 //____________________________________________________________________
 299 void
 300 AliFMDEventInspector::SetMaxCentrality(Double_t maxCent)
 301 {
 302   AliWarning("\n"
 303              "*******************************************************\n"
 304              "* Setting centrality cuts in this stage is deprecated *\n"
 305              "*******************************************************");
 306   fMaxCent = maxCent;
 307 }
 308
 309 //____________________________________________________________________
 310 Bool_t
 311 AliFMDEventInspector::FetchHistograms(const TList* d,
 312                                       TH1I*& hEventsTr,
 313                                       TH1I*& hEventsTrVtx,
 314                                       TH1I*& hEventsAcc,
 315                                       TH1I*& hTriggers) const
 316 {
 317   //
 318   // Fetch our histograms from the passed list
 319   //
 320   // Parameters:
 321   //   d             Input
 322   //   hEventsTr     On return, pointer to histogram, or null
 323   //   hEventsTrVtx  On return, pointer to histogram, or null
 324   //   hTriggers     On return, pointer to histogram, or null
 325   //
 326   // Return:
 327   //    true on success, false otherwise
 328   //
 329   DGUARD(fDebug,3,"Fetch histograms in AliFMDEventInspector");
 330   hEventsTr    = 0;
 331   hEventsTrVtx = 0;
 332   hEventsAcc   = 0;
 333   hTriggers    = 0;
 334   TList* dd    = dynamic_cast<TList*>(d->FindObject(GetName()));
 335   if (!dd) return kFALSE;
 336
 337   hEventsTr    = dynamic_cast<TH1I*>(dd->FindObject("nEventsTr"));
 338   hEventsTrVtx = dynamic_cast<TH1I*>(dd->FindObject("nEventsTrVtx"));
 339   hEventsAcc   = dynamic_cast<TH1I*>(dd->FindObject("nEventsAccepted"));
 340   hTriggers    = dynamic_cast<TH1I*>(dd->FindObject("triggers"));
 341
 342   if (!hEventsTr    ||
 343       !hEventsTrVtx ||
 344       !hEventsAcc   ||
 345       !hTriggers) return kFALSE;
 346   return kTRUE;
 347 }
 348 //____________________________________________________________________
 349 void
 350 AliFMDEventInspector::CacheConfiguredTriggerClasses(TList& cache,
 351                                                     const TList* classes,
 352                                                     AliOADBPhysicsSelection* o)
 353 {
 354   TIter nextClass(classes);
 355   TObjString* trigClass = 0;
 356   // Loop over all trigger classes.  Trigger classes have the format
 357   //
 358   //   class          := positive_words SPACE(s) negative_words
 359   //   positive_words :=
 360   //                  |  '+' words
 361   //   negative_words :=
 362   //                  |  '-' words
 363   //   words          := word
 364   //                  |  word ',' words
 365   //
 366   while ((trigClass = static_cast<TObjString*>(nextClass()))) {
 367     // Tokenize on space to get positive and negative parts
 368     TString     side   = o->GetBeamSide(trigClass->String());
 369     TObjArray*  parts  = trigClass->String().Tokenize(" ");
 370     TObjString* part   = 0;
 371     TIter       nextPart(parts);
 372     while ((part = static_cast<TObjString*>(nextPart()))) {
 373       // We only care about the positive ones
 374       if (part->GetName()[0] != '+') continue;
 375       part->String().Remove(0,1);
 376
 377       // Tokenize on a comma to get the words
 378       TObjArray*  words = part->String().Tokenize(",");
 379       TObjString* word  = 0;
 380       TIter       nextWord(words);
 381       while ((word = static_cast<TObjString*>(nextWord()))) {
 382         TNamed* store = new TNamed(word->String(), side);
 383         cache.Add(store);
 384         DMSG(fDebug,3,"Caching %s trigger word %s",
 385              store->GetTitle(), store->GetName());
 386       } // while (word)
 387       delete words;
 388     }
 389     delete parts;
 390   }
 391 }
 392
 393 //____________________________________________________________________
 394 void
 395 AliFMDEventInspector::SetupForData(const TAxis& vtxAxis)
 396 {
 397   //
 398   // Initialize the object - this is called on the first seen event.
 399   //
 400   // Parameters:
 401   //   vtxAxis Vertex axis in use
 402   //
 403   DGUARD(fDebug,1,"Initialize in AliFMDEventInspector");
 404
 405   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 406
 407   // Get the input handler - should always be there
 408   AliInputEventHandler* ih =
 409     static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
 410   if (!ih) {
 411     AliWarning("No input handler");
 412     return;
 413   }
 414   // Get the physics selection - should always be there
 415   AliPhysicsSelection* ps =
 416     static_cast<AliPhysicsSelection*>(ih->GetEventSelection());
 417   if (!ps) {
 418     AliWarning("No physics selection");
 419     return;
 420   }
 421   // Get the configured triggers
 422   AliOADBPhysicsSelection* oadb =
 423     const_cast<AliOADBPhysicsSelection*>(ps->GetOADBPhysicsSelection());
 424   if (!oadb) {
 425     AliWarning("No OADB physics selection object");
 426     return;
 427   }
 428   // Get the configured trigger words from the physics selection
 429   const TList* collTriggClasses = ps->GetCollisionTriggerClasses();
 430   const TList* bgTriggClasses   = ps->GetBGTriggerClasses();
 431   if (!collTriggClasses) {
 432     AliWarning("No configured collision trigger classes");
 433     return;
 434   }
 435   if (!bgTriggClasses) {
 436     AliWarning("No configured background trigger classes");
 437     return;
 438   }
 439   CacheConfiguredTriggerClasses(fCollWords, collTriggClasses, oadb);
 440   CacheConfiguredTriggerClasses(fBgWords,   bgTriggClasses,   oadb);
 441   // fCollWords.ls();
 442   // fBgWords.ls();
 443
 444
 445   TArrayD limits;
 446   if ((fMinCent < 0 && fMaxCent < 0) || fMaxCent <= fMinCent) {
 447     // -1.5 -0.5 0.5 1.5 ... 89.5 ... 100.5
 448     // ----- 92 number --------- ---- 1 ---
 449     limits.Set(93);
 450     for (Int_t i = 0; i < 92; i++) limits[i] = -1.5 + i;
 451     limits[92] = 100.5;
 452   }
 453   else {
 454     Int_t n = fMaxCent-fMinCent+2;
 455     limits.Set(n);
 456     for (Int_t i = 0; i < n; i++) {
 457       limits[i] = fMinCent + i - .5;
 458     }
 459   }
 460
 461   fVtxAxis.Set(vtxAxis.GetNbins(), vtxAxis.GetXmin(), vtxAxis.GetXmax());
 462
 463   fCentAxis  = new TAxis(limits.GetSize()-1, limits.GetArray());
 464   fHEventsTr = new TH1I("nEventsTr", "Number of events w/trigger",
 465                         4*vtxAxis.GetNbins(),
 466                         2*vtxAxis.GetXmin(),
 467                         2*vtxAxis.GetXmax());
 468   fHEventsTr->SetXTitle("v_{z} [cm]");
 469   fHEventsTr->SetYTitle("# of events");
 470   fHEventsTr->SetFillColor(kRed+1);
 471   fHEventsTr->SetFillStyle(3001);
 472   fHEventsTr->SetDirectory(0);
 473   // fHEventsTr->Sumw2();
 474   fList->Add(fHEventsTr);
 475
 476   fHEventsTrVtx = static_cast<TH1I*>(fHEventsTr->Clone("nEventsTrVtx"));
 477   fHEventsTrVtx->SetTitle("Number of events w/trigger and vertex");
 478   fHEventsTrVtx->SetFillColor(kBlue+1);
 479   fHEventsTrVtx->SetDirectory(0);
 480   // fHEventsTrVtx->Sumw2();
 481   fList->Add(fHEventsTrVtx);
 482
 483   fHEventsAccepted = new TH1I("nEventsAccepted",
 484                               "Number of events  w/trigger and vertex in range",
 485                               2*vtxAxis.GetNbins(),
 486                               2*vtxAxis.GetXmin(),
 487                               2*vtxAxis.GetXmax());
 488   fHEventsAccepted->SetXTitle("v_{z} [cm]");
 489   fHEventsAccepted->SetYTitle("# of events");
 490   fHEventsAccepted->SetFillColor(kGreen+1);
 491   fHEventsAccepted->SetFillStyle(3001);
 492   fHEventsAccepted->SetDirectory(0);
 493   // fHEventsAccepted->Sumw2();
 494   fList->Add(fHEventsAccepted);
 495
 496   fHEventsAcceptedXY = new TH2D("nEventsAcceptedXY",
 497                               "XY vertex w/trigger and Z vertex in range",
 498                                 1000,-1,1,1000,-1,1);
 499
 500   fHEventsAcceptedXY->SetXTitle("v_{x} [cm]");
 501   fHEventsAcceptedXY->SetYTitle("v_{y} [cm]");
 502   fHEventsAcceptedXY->SetDirectory(0);
 503   // fHEventsAccepted->Sumw2();
 504   fList->Add(fHEventsAcceptedXY);
 505
 506
 507   fHTriggers = new TH1I("triggers", "Triggers", kOffline+1, 0, kOffline+1);
 508   fHTriggers->SetFillColor(kRed+1);
 509   fHTriggers->SetFillStyle(3001);
 510   fHTriggers->SetStats(0);
 511   fHTriggers->SetDirectory(0);
 512
 513   fHTriggerCorr = new TH2I("triggerCorr", "Trigger correlation",
 514                            kOffline+1, 0, kOffline+1,
 515                            kOffline+1, 0, kOffline+1);
 516   fHTriggerCorr->SetStats(0);
 517   fHTriggerCorr->SetDirectory(0);
 518   fHTriggerCorr->SetXTitle("Requirement");
 519   fHTriggerCorr->SetYTitle("Companion");
 520
 521   Int_t binNum[] = { kInel   +1,
 522                      kInelGt0+1,
 523                      kNSD    +1,
 524                      kV0AND  +1,
 525                      kEmpty  +1,
 526                      kA      +1,
 527                      kB      +1,
 528                      kC      +1,
 529                      kE      +1,
 530                      kPileUp +1,
 531                      kMCNSD  +1,
 532                      kSatellite+1,
 533                      kOffline+1 };
 534   const char* binLbl[] = { "INEL",
 535                            "INEL>0",
 536                            "NSD",
 537                            "VOAND",
 538                            "Empty",
 539                            "A",
 540                            "B",
 541                            "C",
 542                            "E",
 543                            "Pileup",
 544                            "NSD_{MC}",
 545                            "Satellite",
 546                            "Offline" };
 547   for (Int_t i = 0; i < kOffline+1; i++) {
 548     fHTriggers->GetXaxis()->SetBinLabel(binNum[i], binLbl[i]);
 549     fHTriggerCorr->GetXaxis()->SetBinLabel(binNum[i], binLbl[i]);
 550     fHTriggerCorr->GetYaxis()->SetBinLabel(binNum[i], binLbl[i]);
 551   }
 552   fList->Add(fHTriggers);
 553   fList->Add(fHTriggerCorr);
 554
 555
 556   fHType = new TH1I("type", Form("Event type (cut: SPD mult>%d)",
 557                                  fLowFluxCut), 2, -.5, 1.5);
 558   fHType->SetFillColor(kRed+1);
 559   fHType->SetFillStyle(3001);
 560   fHType->SetStats(0);
 561   fHType->SetDirectory(0);
 562   fHType->GetXaxis()->SetBinLabel(1,"Low-flux");
 563   fHType->GetXaxis()->SetBinLabel(2,"High-flux");
 564   fList->Add(fHType);
 565
 566 #if 0
 567   // This histogram disabled as it causes problems in the merge
 568   fHWords = new TH1I("words", "Trigger words seen", 1, 0, 0);
 569   fHWords->SetFillColor(kBlue+1);
 570   fHWords->SetFillStyle(3001);
 571   fHWords->SetStats(0);
 572   fHWords->SetDirectory(0);
 573   fHWords->SetBit(TH1::kCanRebin);
 574   fList->Add(fHWords);
 575 #endif
 576
 577   fHCent = new TH1F("cent", "Centrality", limits.GetSize()-1,limits.GetArray());
 578   fHCent->SetFillColor(kBlue+1);
 579   fHCent->SetFillStyle(3001);
 580   fHCent->SetStats(0);
 581   fHCent->SetDirectory(0);
 582   fHCent->SetXTitle("Centrality [%]");
 583   fHCent->SetYTitle("Events");
 584   fList->Add(fHCent);
 585
 586   fHCentVsQual = new TH2F("centVsQuality", "Quality vs Centrality",
 587                           5, 0, 5, limits.GetSize()-1, limits.GetArray());
 588   fHCentVsQual->SetXTitle("Quality");
 589   fHCentVsQual->SetYTitle("Centrality [%]");
 590   fHCentVsQual->SetZTitle("Events");
 591   fHCentVsQual->GetXaxis()->SetBinLabel(1, "OK");
 592   fHCentVsQual->GetXaxis()->SetBinLabel(2, "Outside v_{z} cut");
 593   fHCentVsQual->GetXaxis()->SetBinLabel(3, "V0 vs SPD outlier");
 594   fHCentVsQual->GetXaxis()->SetBinLabel(4, "V0 vs TPC outlier");
 595   fHCentVsQual->GetXaxis()->SetBinLabel(5, "V0 vs ZDC outlier");
 596   fHCentVsQual->SetDirectory(0);
 597   fList->Add(fHCentVsQual);
 598
 599   fHStatus = new TH1I("status", "Status", 7, 1, 8);
 600   fHStatus->SetFillColor(kBlue+1);
 601   fHStatus->SetFillStyle(3001);
 602   fHStatus->SetStats(0);
 603   fHStatus->SetDirectory(0);
 604   TAxis* xAxis = fHStatus->GetXaxis();
 605   xAxis->SetBinLabel(1, "OK");
 606   xAxis->SetBinLabel(2, "No event");
 607   xAxis->SetBinLabel(3, "No triggers");
 608   xAxis->SetBinLabel(4, "No SPD");
 609   xAxis->SetBinLabel(5, "No FMD");
 610   xAxis->SetBinLabel(6, "No vertex");
 611   xAxis->SetBinLabel(7, "Bad vertex");
 612   fList->Add(fHStatus);
 613
 614   fHVtxStatus = new TH1I("vtxStatus","Vertex Status",
 615                          kNotVtxZ,kVtxOK,kNotVtxZ+1);
 616   fHVtxStatus->SetFillColor(kGreen+1);
 617   fHVtxStatus->SetFillStyle(3001);
 618   fHVtxStatus->SetStats(0);
 619   fHVtxStatus->SetDirectory(0);
 620   xAxis = fHVtxStatus->GetXaxis();
 621   xAxis->SetBinLabel(kVtxOK,      "OK");
 622   xAxis->SetBinLabel(kNoVtx,      "None/bad status");
 623   xAxis->SetBinLabel(kNoSPDVtx,   "No SPD/bad status");
 624   xAxis->SetBinLabel(kFewContrib, "N_{contrib} <= 0");
 625   xAxis->SetBinLabel(kUncertain,  Form("#delta z > %4.2f", fMaxVzErr));
 626   xAxis->SetBinLabel(kNotVtxZ,    "Not Z vertexer");
 627   fList->Add(fHVtxStatus);
 628
 629   fHTrgStatus = new TH1I("trgStatus", "Trigger Status",
 630                          kOther, kNoTrgWords, kOther+1);
 631   fHTrgStatus->SetFillColor(kMagenta+1);
 632   fHTrgStatus->SetFillStyle(3001);
 633   fHTrgStatus->SetStats(0);
 634   fHTrgStatus->SetDirectory(0);
 635   xAxis = fHTrgStatus->GetXaxis();
 636   xAxis->SetBinLabel(kNoTrgWords,       "No words");
 637   xAxis->SetBinLabel(kPP2760Fast,       "FAST in pp@#sqrt{s}=2.76TeV");
 638   xAxis->SetBinLabel(kMUON,             "Muon trigger");
 639   xAxis->SetBinLabel(kTriggered,        "Triggered");
 640   xAxis->SetBinLabel(kMinBias,          "CINT1 (V0A||V0C||FASTOR)");
 641   xAxis->SetBinLabel(kMinBiasNoSPD,     "CINT5 (V0A||V0C)");
 642   xAxis->SetBinLabel(kV0AndTrg,         "CINT7 (V0A&&V0C)");
 643   xAxis->SetBinLabel(kHighMult,         "N>>0");
 644   xAxis->SetBinLabel(kCentral,          "Central");
 645   xAxis->SetBinLabel(kSemiCentral,      "Semi-central");
 646   xAxis->SetBinLabel(kDiffractive,      "Diffractive");
 647   xAxis->SetBinLabel(kUser,             "User");
 648   xAxis->SetBinLabel(kOther,            "Other");
 649   fList->Add(fHTrgStatus);
 650
 651   if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", false);
 652 }
 653
 654 //____________________________________________________________________
 655 void
 656 AliFMDEventInspector::StoreInformation()
 657 {
 658   // Write TNamed objects to output list containing information about
 659   // the running conditions
 660   DGUARD(fDebug,2,"Store information from AliFMDEventInspector");
 661   if (!fList) return;
 662
 663
 664   fList->Add(AliForwardUtil::MakeParameter("sys", fCollisionSystem));
 665   fList->Add(AliForwardUtil::MakeParameter("sNN", fEnergy));
 666   fList->Add(AliForwardUtil::MakeParameter("field", fField));
 667   fList->Add(AliForwardUtil::MakeParameter("runNo", fRunNumber));
 668   fList->Add(AliForwardUtil::MakeParameter("lowFlux", fLowFluxCut));
 669   fList->Add(AliForwardUtil::MakeParameter("fpVtx",fUseFirstPhysicsVertex));
 670   fList->Add(AliForwardUtil::MakeParameter("v0and",fUseV0AND));
 671   fList->Add(AliForwardUtil::MakeParameter("nPileUp", fMinPileupContrib));
 672   fList->Add(AliForwardUtil::MakeParameter("dPileup", fMinPileupDistance));
 673   fList->Add(AliForwardUtil::MakeParameter("satellite", fUseDisplacedVertices));
 674   fList->Add(AliForwardUtil::MakeParameter("alirootRev",
 675                                            AliForwardUtil::AliROOTRevision()));
 676   fList->Add(AliForwardUtil::MakeParameter("alirootBranch",
 677                                            AliForwardUtil::AliROOTBranch()));
 678
 679 }
 680
 681 //____________________________________________________________________
 682 void
 683 AliFMDEventInspector::CreateOutputObjects(TList* dir)
 684 {
 685   //
 686   // Define the output histograms.  These are put in a sub list of the
 687   // passed list.   The histograms are merged before the parent task calls
 688   // AliAnalysisTaskSE::Terminate
 689   //
 690   //   dir Directory to add to
 691   //
 692   DGUARD(fDebug,1,"Define output from AliFMDEventInspector");
 693   fList = new TList;
 694   fList->SetName(GetName());
 695   fList->SetOwner();
 696   dir->Add(fList);
 697 }
 698
 699 //____________________________________________________________________
 700 UInt_t
 701 AliFMDEventInspector::Process(const AliESDEvent* event,
 702                               UInt_t&            triggers,
 703                               Bool_t&            lowFlux,
 704                               UShort_t&          ivz,
 705                               TVector3&          ip,
 706                               Double_t&          cent,
 707                               UShort_t&          nClusters)
 708 {
 709   //
 710   // Process the event
 711   //
 712   // Parameters:
 713   //   event     Input event
 714   //   triggers  On return, the triggers fired
 715   //   lowFlux   On return, true if the event is considered a low-flux
 716   //                  event (according to the setting of fLowFluxCut)
 717   //   ivz       On return, the found vertex bin (1-based).  A zero
 718   //                  means outside of the defined vertex range
 719   //   vz        On return, the z position of the interaction
 720   //   cent      On return, the centrality - if not available < 0
 721   //
 722   // Return:
 723   //    0 (or kOk) on success, otherwise a bit mask of error codes
 724   //
 725   DGUARD(fDebug,1,"Process event in AliFMDEventInspector");
 726   // --- Check that we have an event ---------------------------------
 727   if (!event) {
 728     AliWarning("No ESD event found for input event");
 729     fHStatus->Fill(2);
 730     return kNoEvent;
 731   }
 732
 733   // --- Process satellite event information is requested ------------
 734   if (fUseDisplacedVertices) {
 735     if (!fDisplacedVertex.Process(event))
 736       AliWarning("Failed to process satellite event");
 737   }
 738
 739   // --- Read trigger information from the ESD and store in AOD object
 740   if (!ReadTriggers(*event, triggers, nClusters)) {
 741     if (fDebug > 2) {
 742       AliWarning("Failed to read triggers from ESD"); }
 743     fHStatus->Fill(3);
 744     return kNoTriggers;
 745   }
 746
 747   // --- Check if this is a high-flux event --------------------------
 748   const AliMultiplicity* testmult = event->GetMultiplicity();
 749   if (!testmult) {
 750     if (fDebug > 3) {
 751       AliWarning("No central multiplicity object found"); }
 752   }
 753   else
 754     lowFlux = testmult->GetNumberOfTracklets() < fLowFluxCut;
 755
 756   fHType->Fill(lowFlux ? 0 : 1);
 757
 758   // --- Get the interaction point -----------------------------------
 759   Bool_t vzOk = ReadVertex(*event, ip);
 760   fHEventsTr->Fill(ip.Z());
 761   if (!vzOk) {
 762     if (fDebug > 3) {
 763       AliWarning("Failed to read vertex from ESD"); }
 764     fHStatus->Fill(6);
 765     return kNoVertex;
 766   }
 767
 768    // --- Read centrality information
 769   cent          = -10;
 770   UShort_t qual = 0;
 771   if (!ReadCentrality(*event, cent, qual)) {
 772     if (fDebug > 3)
 773       AliWarning("Failed to get centrality");
 774   }
 775   // --- check centrality cut
 776
 777   if(fMinCent > -0.0001 && cent < fMinCent) return  kNoEvent;
 778   if(fMaxCent > -0.0001 && cent > fMaxCent) return  kNoEvent;
 779   fHCent->Fill(cent);
 780   if (qual == 0) fHCentVsQual->Fill(0., cent);
 781   else {
 782     for (UShort_t i = 0; i < 4; i++)
 783       if (qual & (1 << i)) fHCentVsQual->Fill(Double_t(i+1), cent);
 784   }
 785
 786
 787   fHEventsTrVtx->Fill(ip.Z());
 788
 789   // --- Get the vertex bin ------------------------------------------
 790   ivz = fVtxAxis.FindBin(ip.Z());
 791   if (ivz <= 0 || ivz > fVtxAxis.GetNbins()) {
 792     if (fDebug > 3) {
 793       AliWarning(Form("Vertex @ %f outside of range [%f,%f]",
 794                       ip.Z(), fVtxAxis.GetXmin(), fVtxAxis.GetXmax()));
 795     }
 796     ivz = 0;
 797     fHStatus->Fill(7);
 798     return kBadVertex;
 799   }
 800   fHEventsAccepted->Fill(ip.Z());
 801   fHEventsAcceptedXY->Fill(ip.X(),ip.Y());
 802
 803   // --- Check the FMD ESD data --------------------------------------
 804   if (!event->GetFMDData()) {
 805     if (fDebug > 3) {
 806       AliWarning("No FMD data found in ESD"); }
 807     fHStatus->Fill(5);
 808     return kNoFMD;
 809   }
 810
 811   fHStatus->Fill(1);
 812   return kOk;
 813 }
 814
 815 //____________________________________________________________________
 816 Bool_t
 817 AliFMDEventInspector::ReadCentrality(const AliESDEvent& esd,
 818                                      Double_t& cent,
 819                                      UShort_t& qual) const
 820 {
 821   //
 822   // Read centrality from event
 823   //
 824   // Parameters:
 825   //    esd  Event
 826   //    cent On return, the centrality or negative if not found
 827   //
 828   // Return:
 829   //    False on error, true otherwise
 830   //
 831   DGUARD(fDebug,2,"Read the centrality in AliFMDEventInspector");
 832
 833   cent = -1;
 834   qual = 1;
 835   AliCentrality* centObj = const_cast<AliESDEvent&>(esd).GetCentrality();
 836   if (centObj) {
 837     cent = centObj->GetCentralityPercentile(fCentMethod);
 838     qual = centObj->GetQuality();
 839   }
 840
 841   if (qual > 0 && fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
 842     cent = fDisplacedVertex.GetCentralityPercentile();
 843     qual = 0;
 844   }
 845
 846   return true;
 847 }
 848
 849 //____________________________________________________________________
 850 Bool_t
 851 AliFMDEventInspector::CheckpAExtraV0(const AliESDEvent& esd) const
 852 {
 853   if (fCollisionSystem != AliForwardUtil::kPPb) return true;
 854
 855    AliVVZERO* esdV0 = esd.GetVZEROData();
 856    if ((esdV0->GetV0ADecision()!=1) || (esdV0->GetV0CDecision()!=1))
 857      return false;
 858    return true;
 859 }
 860
 861 //____________________________________________________________________
 862 Bool_t
 863 AliFMDEventInspector::ReadTriggers(const AliESDEvent& esd, UInt_t& triggers,
 864                                    UShort_t& nClusters)
 865 {
 866   //
 867   // Read the trigger information from the ESD event
 868   //
 869   // Parameters:
 870   //   esd        ESD event
 871   //   triggers   On return, contains the trigger bits
 872   //
 873   // Return:
 874   //    @c true on success, @c false otherwise
 875   //
 876   DGUARD(fDebug,2,"Read the triggers in AliFMDEventInspector");
 877   triggers = 0;
 878
 879   // Get the analysis manager - should always be there
 880   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 881   DMSG(fDebug,10,"Got analysis manager %p", am);
 882   if (!am) {
 883     AliWarning("No analysis manager defined!");
 884     return kFALSE;
 885   }
 886
 887   // Get the input handler - should always be there
 888   AliInputEventHandler* ih =
 889     static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
 890   DMSG(fDebug,10,"Got input handler %p", ih);
 891   if (!ih) {
 892     AliWarning("No input handler");
 893     return kFALSE;
 894   }
 895
 896   // Check if this is a collision candidate (MB)
 897   ///
 898   // Historic remark: Note, that we should use the value cached in the
 899   //   input handler rather than calling IsCollisionCandiate directly
 900   //   on the AliPhysicsSelection obejct.  If we called the latter
 901   //   then the AliPhysicsSelection object would overcount by a factor
 902   //   of 2! :-(
 903   UInt_t  trgMask  = ih->IsEventSelected();
 904   Bool_t  offline  = trgMask;
 905   Bool_t  fastonly = (trgMask & AliVEvent::kFastOnly);
 906   TString trigStr  = esd.GetFiredTriggerClasses();
 907
 908   if (trigStr.IsNull()) fHTrgStatus->Fill(kNoTrgWords);
 909   if (fHWords) fHWords->Fill(trigStr.Data(), 1);
 910
 911   if(fUseDisplacedVertices) {
 912     DMSG(fDebug,3,"Using displaced vertex stuff");
 913     // if (TMath::Abs(fDisplacedVertex.GetVertexZ()) >= 999) offline = false;
 914     if (fDisplacedVertex.IsSatellite())
 915       triggers |= AliAODForwardMult::kSatellite;
 916   }
 917
 918   if (CheckFastPartition(fastonly)) {
 919     fHTrgStatus->Fill(kPP2760Fast);
 920     offline = false;
 921   }
 922
 923   if (offline && CheckCosmics(trigStr)) {
 924     fHTrgStatus->Fill(kMUON);
 925     offline = false;
 926   }
 927   if (offline) fHTrgStatus->Fill(kTriggered);
 928   Int_t f = 0;
 929   if (trgMask & AliVEvent::kMB)          f += fHTrgStatus->Fill(kMinBias);
 930   if (trgMask & AliVEvent::kCINT5)       f += fHTrgStatus->Fill(kMinBiasNoSPD);
 931   if (trgMask & AliVEvent::kINT7)        f += fHTrgStatus->Fill(kV0AndTrg);
 932   if (trgMask & AliVEvent::kHighMult)    f += fHTrgStatus->Fill(kHighMult);
 933   if (trgMask & AliVEvent::kCentral)     f += fHTrgStatus->Fill(kCentral);
 934   if (trgMask & AliVEvent::kSemiCentral) f += fHTrgStatus->Fill(kSemiCentral);
 935   if (trgMask & AliVEvent::kDG5)         f += fHTrgStatus->Fill(kDiffractive);
 936   if (trgMask & AliVEvent::kUserDefined) f += fHTrgStatus->Fill(kUser);
 937   if (f <= 0) fHTrgStatus->Fill(kOther);
 938
 939   // if (!CheckpAExtraV0(esd))             offline = false;
 940
 941   DMSG(fDebug,2,"Event is %striggered by off-line", offline ? "" : "NOT ");
 942
 943   if (offline) {
 944     triggers |= AliAODForwardMult::kOffline;
 945     triggers |= AliAODForwardMult::kInel;
 946     CheckINELGT0(esd, nClusters, triggers);
 947   }
 948
 949   CheckNSD(esd,triggers);
 950   CheckPileup(esd, triggers);
 951   CheckEmpty(trigStr, triggers);
 952   // if (CheckPileup(esd, triggers)) fHTriggers->Fill(kPileUp+.5);
 953   // if (CheckEmpty(trigStr, triggers)) fHTriggers->Fill(kEmpty+.5);
 954
 955   CheckWords(esd, triggers);
 956
 957 #define TEST_TRIG_BIN(RET,BIN,TRIGGERS) \
 958   do { switch (BIN) { \
 959     case kInel:     RET = triggers & AliAODForwardMult::kInel;      break; \
 960     case kInelGt0:  RET = triggers & AliAODForwardMult::kInelGt0;   break; \
 961     case kNSD:      RET = triggers & AliAODForwardMult::kNSD;       break; \
 962     case kV0AND:    RET = triggers & AliAODForwardMult::kV0AND;     break; \
 963     case kEmpty:    RET = triggers & AliAODForwardMult::kEmpty;     break; \
 964     case kA:        RET = triggers & AliAODForwardMult::kA;         break; \
 965     case kB:        RET = triggers & AliAODForwardMult::kB;         break; \
 966     case kC:        RET = triggers & AliAODForwardMult::kC;         break; \
 967     case kE:        RET = triggers & AliAODForwardMult::kE;         break; \
 968     case kPileUp:   RET = triggers & AliAODForwardMult::kPileUp;    break; \
 969     case kMCNSD:    RET = triggers & AliAODForwardMult::kMCNSD;     break; \
 970     case kSatellite:RET = triggers & AliAODForwardMult::kSatellite; break; \
 971     case kOffline:  RET = triggers & AliAODForwardMult::kOffline;   break; \
 972     default:        RET = false; } } while(false)
 973
 974   if (!fHTriggers) {
 975     AliWarning("Histogram of triggers not defined - has init been called");
 976     return false;
 977   }
 978
 979   for (Int_t i = 0; i < kOffline+1; i++) {
 980     Bool_t hasX = false;
 981     TEST_TRIG_BIN(hasX, i, triggers);
 982     if (!hasX) continue;
 983     fHTriggers->Fill(i+.5);
 984     for (Int_t j = 0; j < kOffline+1; j++) {
 985       Bool_t hasY = false;
 986       TEST_TRIG_BIN(hasY, j, triggers);
 987       if (!hasY) continue;
 988
 989       fHTriggerCorr->Fill(i+.5, j+.5);
 990     }
 991   }
 992   return kTRUE;
 993 }
 994
 995 //____________________________________________________________________
 996 Bool_t
 997 AliFMDEventInspector::CheckFastPartition(bool fastonly) const
 998 {
 999   // For the 2.76 TeV p+p run, the FMD ran in the slow partition
1000   // so it received no triggers from the fast partition. Therefore
1001   // the fast triggers are removed here but not for MC where all
1002   // triggers are fast.
1003   if (TMath::Abs(fEnergy - 2750.) > 20) return false;
1004   if (fCollisionSystem != AliForwardUtil::kPP) return false;
1005   if (fastonly)
1006     DMSG(fDebug,1,"Fast trigger in pp @ sqrt(s)=2.76TeV removed");
1007
1008   return fastonly;
1009 }
1010
1011 //____________________________________________________________________
1012 Bool_t
1013 AliFMDEventInspector::CheckCosmics(const TString& trigStr) const
1014 {
1015   // MUON triggers are not strictly minimum bias (MB) so they are
1016   // removed (HHD)
1017   if(trigStr.Contains("CMUS1")) {
1018     DMSG(fDebug,1,"Cosmic trigger ins't min-bias, removed");
1019     return true;
1020   }
1021   return false;
1022 }
1023
1024 //____________________________________________________________________
1025 Bool_t
1026 AliFMDEventInspector::CheckINELGT0(const AliESDEvent& esd,
1027                                    UShort_t& nClusters,
1028                                    UInt_t& triggers) const
1029 {
1030   nClusters = 0;
1031
1032   // If this is inel, see if we have a tracklet
1033   const AliMultiplicity* spdmult = esd.GetMultiplicity();
1034   if (!spdmult) {
1035     AliWarning("No SPD multiplicity");
1036     return false;
1037   }
1038
1039   // Check if we have one or more tracklets
1040   // in the range -1 < eta < 1 to set the INEL>0
1041   // trigger flag.
1042   //
1043   // Also count tracklets as a single cluster
1044   Int_t n = spdmult->GetNumberOfTracklets();
1045   for (Int_t j = 0; j < n; j++) {
1046     if(TMath::Abs(spdmult->GetEta(j)) < 1) {
1047       triggers |= AliAODForwardMult::kInelGt0;
1048       nClusters++;
1049     }
1050   }
1051   n = spdmult->GetNumberOfSingleClusters();
1052   for (Int_t j = 0; j < n; j++) {
1053     Double_t eta = -TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.));
1054     if (TMath::Abs(eta) < 1) nClusters++;
1055   }
1056   if (nClusters > 0) triggers |= AliAODForwardMult::kNClusterGt0;
1057
1058   return triggers & AliAODForwardMult::kNClusterGt0;
1059 }
1060
1061 //____________________________________________________________________
1062 Bool_t
1063 AliFMDEventInspector::CheckNSD(const AliESDEvent& esd, UInt_t& triggers) const
1064 {
1065   // Analyse some trigger stuff
1066   AliTriggerAnalysis ta;
1067   if (ta.IsOfflineTriggerFired(&esd, AliTriggerAnalysis::kV0AND)) {
1068     triggers |= AliAODForwardMult::kV0AND;
1069     if (fUseV0AND)
1070       triggers |= AliAODForwardMult::kNSD;
1071   }
1072   if (ta.IsOfflineTriggerFired(&esd, AliTriggerAnalysis::kNSD1))
1073     triggers |= AliAODForwardMult::kNSD;
1074   return triggers & AliAODForwardMult::kNSD;
1075 }
1076 //____________________________________________________________________
1077 Bool_t
1078 AliFMDEventInspector::CheckPileup(const AliESDEvent& esd,
1079                                   UInt_t& triggers) const
1080 {
1081   // Check for multiple vertices (pile-up) with at least 3
1082   // contributors and at least 0.8cm from the primary vertex
1083   if(fCollisionSystem != AliForwardUtil::kPP) return false;
1084
1085   Bool_t pileup =  esd.IsPileupFromSPD(fMinPileupContrib,fMinPileupDistance);
1086   if (pileup) triggers |= AliAODForwardMult::kPileUp;
1087   return pileup;
1088 }
1089
1090 //____________________________________________________________________
1091 Bool_t
1092 AliFMDEventInspector::CheckEmpty(const TString& trigStr, UInt_t& triggers) const
1093 {
1094   if (trigStr.Contains("CBEAMB-ABCE-NOPF-ALL")) {
1095     triggers |= AliAODForwardMult::kEmpty;
1096     return true;
1097   }
1098   return false;
1099 }
1100 //____________________________________________________________________
1101 Bool_t
1102 AliFMDEventInspector::CheckWords(const AliESDEvent& esd, UInt_t& triggers) const
1103 {
1104   TObject* word = 0;
1105   TIter nextColl(&fCollWords);
1106   while ((word = nextColl())) {
1107     DMSG(fDebug,10,"Checking if %s trigger %s is fired",
1108          word->GetTitle(), word->GetName());
1109     if (!esd.IsTriggerClassFired(word->GetName())) continue;
1110
1111     TString beamSide = word->GetTitle();
1112     DMSG(fDebug,10,"Found it - this is a %s trigger", beamSide.Data());
1113
1114     if (!beamSide.EqualTo("B")) continue;
1115     triggers |= AliAODForwardMult::kB;
1116     break; // No more to do here
1117   }
1118   TIter nextBg(&fBgWords);
1119   UInt_t all = (AliAODForwardMult::kA |
1120                 AliAODForwardMult::kC |
1121                 AliAODForwardMult::kE);
1122   while ((word = nextBg())) {
1123     DMSG(fDebug,10,"Checking if %s trigger %s is fired",
1124          word->GetTitle(), word->GetName());
1125     if (!esd.IsTriggerClassFired(word->GetName())) continue;
1126
1127     TString beamSide = word->GetTitle();
1128     DMSG(fDebug,10,"Found it - this is a %s trigger", beamSide.Data());
1129
1130     if (beamSide.Contains("A")) triggers |= AliAODForwardMult::kA;
1131     if (beamSide.Contains("C")) triggers |= AliAODForwardMult::kC;
1132     if (beamSide.Contains("E")) triggers |= AliAODForwardMult::kE;
1133
1134     if ((triggers & all) == all) break; // No more to do
1135   }
1136   return true;
1137 }
1138
1139
1140 //____________________________________________________________________
1141 Bool_t
1142 AliFMDEventInspector::ReadVertex(const AliESDEvent& esd, TVector3& ip)
1143 {
1144   //
1145   // Read the vertex information from the ESD event
1146   //
1147   // Parameters:
1148   //   esd  ESD event
1149   //   vz   On return, the vertex Z position
1150   //
1151   // Return:
1152   //    @c true on success, @c false otherwise
1153   //
1154   DGUARD(fDebug,2,"Read the vertex in AliFMDEventInspector");
1155   ip.SetXYZ(1024, 1024, 0);
1156
1157   EVtxStatus s = kNoVtx;
1158   if (fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
1159     s = kVtxOK;
1160     ip.SetZ(fDisplacedVertex.GetVertexZ());
1161   }
1162   else if (fUseFirstPhysicsVertex)
1163     s = CheckPWGUDVertex(esd, ip);
1164   else if (fUsepA2012Vertex)
1165     s = CheckpA2012Vertex(esd,ip);
1166   else
1167     s = CheckVertex(esd, ip);
1168
1169   fHVtxStatus->Fill(s);
1170
1171   return s == kVtxOK;
1172 }
1173
1174 //____________________________________________________________________
1175 AliFMDEventInspector::EVtxStatus
1176 AliFMDEventInspector::CheckPWGUDVertex(const AliESDEvent& esd,
1177                                        TVector3& ip)  const
1178 {
1179   // This is the code used by the 1st physics people
1180   const AliESDVertex* vertex    = esd.GetPrimaryVertex();
1181   if (!vertex  || !vertex->GetStatus()) {
1182     DMSG(fDebug,2,"No primary vertex (%p) or bad status %d",
1183          vertex, (vertex ? vertex->GetStatus() : -1));
1184     return kNoVtx;
1185   }
1186   const AliESDVertex* vertexSPD = esd.GetPrimaryVertexSPD();
1187   if (!vertexSPD || !vertexSPD->GetStatus()) {
1188     DMSG(fDebug,2,"No primary SPD vertex (%p) or bad status %d",
1189          vertexSPD, (vertexSPD ? vertexSPD->GetStatus() : -1));
1190     return kNoSPDVtx;
1191   }
1192
1193   // if vertex is from SPD vertexZ, require more stringent cuts
1194   if (vertex->IsFromVertexerZ()) {
1195     if (vertex->GetDispersion() > fMaxVzErr ||
1196         vertex->GetZRes() > 1.25 * fMaxVzErr) {
1197       DMSG(fDebug,2,"Dispersion %f > %f or resolution %f > %f",
1198            vertex->GetDispersion(), fMaxVzErr,
1199            vertex->GetZRes(), 1.25 * fMaxVzErr);
1200       return kUncertain;
1201     }
1202   }
1203   ip.SetZ(vertex->GetZ());
1204
1205   if(!vertex->IsFromVertexerZ()) {
1206     ip.SetX(vertex->GetX());
1207     ip.SetY(vertex->GetY());
1208   }
1209   return kVtxOK;
1210 }
1211 //____________________________________________________________________
1212 AliFMDEventInspector::EVtxStatus
1213 AliFMDEventInspector::CheckpA2012Vertex(const AliESDEvent& esd,
1214                                         TVector3& ip)  const
1215 {
1216   const AliESDVertex *vertex = esd.GetPrimaryVertexSPD();
1217   if (!vertex) return kNoSPDVtx;
1218   if (vertex->GetNContributors() <= 0) return kFewContrib;
1219
1220   TString vtxTyp = vertex->GetTitle();
1221   if (vtxTyp.Contains("vertexer: Z")) return kNotVtxZ;
1222
1223   if (vertex->GetDispersion() >= 0.04 || vertex->GetZRes()>=0.25)
1224     return kUncertain;
1225
1226   ip.SetX(vertex->GetX());
1227   ip.SetY(vertex->GetY());
1228   ip.SetZ(vertex->GetZ());
1229
1230   return kVtxOK;
1231 }
1232
1233 //____________________________________________________________________
1234 AliFMDEventInspector::EVtxStatus
1235 AliFMDEventInspector::CheckVertex(const AliESDEvent& esd,
1236                                   TVector3& ip) const
1237 {
1238   // Use standard SPD vertex (perhaps preferable for Pb+Pb)
1239   // Get the vertex
1240   const AliESDVertex* vertex = esd.GetPrimaryVertexSPD();
1241   if (!vertex) {
1242     if (fDebug > 2) {
1243       AliWarning("No SPD vertex found in ESD"); }
1244     return kNoSPDVtx;
1245   }
1246
1247   // #if 0 // Check disabled - seem to kill a lot of PbPb events
1248   // Check that enough tracklets contributed
1249   if(vertex->GetNContributors() <= 0) {
1250     DMSG(fDebug,2,"Number of contributors to vertex is %d<=0",
1251          vertex->GetNContributors());
1252     ip.SetZ(0);
1253     return kFewContrib;
1254   }
1255   // #endif
1256
1257   // Check that the uncertainty isn't too large
1258   if (vertex->GetZRes() > fMaxVzErr) {
1259     DMSG(fDebug,2,"Uncertaintity in Z of vertex is too large %f > %f",
1260          vertex->GetZRes(), fMaxVzErr);
1261     return kUncertain;
1262   }
1263
1264   // Get the z coordiante
1265   ip.SetZ(vertex->GetZ());
1266   const AliESDVertex* vertexXY = esd.GetPrimaryVertex();
1267
1268
1269   if(!vertexXY->IsFromVertexerZ()) {
1270     ip.SetX(vertexXY->GetX());
1271     ip.SetY(vertexXY->GetY());
1272   }
1273   return kVtxOK;
1274 }
1275
1276 //____________________________________________________________________
1277 Bool_t
1278 AliFMDEventInspector::ReadRunDetails(const AliESDEvent* esd)
1279 {
1280   //
1281   // Read the collision system, collision energy, and L3 field setting
1282   // from the ESD
1283   //
1284   // Parameters:
1285   //   esd ESD to get information from
1286   //
1287   // Return:
1288   //    true on success, false
1289   //
1290   // AliInfo(Form("Parameters from 1st ESD event: cms=%s, sNN=%f, field=%f",
1291   //           esd->GetBeamType(), 2*esd->GetBeamEnergy(),
1292   //           esd->GetMagneticField()));
1293   DGUARD(fDebug,2,"Read the run details in AliFMDEventInspector");
1294   const char* sys  = esd->GetBeamType();
1295   Float_t     cms  = 2 * esd->GetBeamEnergy();
1296   Float_t     fld  = esd->GetMagneticField();
1297   fCollisionSystem = AliForwardUtil::ParseCollisionSystem(sys);
1298   fEnergy          = AliForwardUtil::ParseCenterOfMassEnergy(fCollisionSystem,
1299                                                              cms);
1300   fField           = AliForwardUtil::ParseMagneticField(fld);
1301   fRunNumber       = esd->GetRunNumber();
1302   StoreInformation();
1303   if (fCollisionSystem   == AliForwardUtil::kUnknown) {
1304     AliWarningF("Unknown collision system: %s - please check", sys);
1305     return false;
1306   }
1307   if (fEnergy            <= 0) {
1308     AliWarningF("Unknown CMS energy: %f (%d) - please check", cms, fEnergy);
1309     return false;
1310   }
1311   if (TMath::Abs(fField) >  10) {
1312     AliWarningF("Unknown L3 field setting: %f (%d) - please check", fld,fField);
1313     return false;
1314   }
1315
1316   return true;
1317 }
1318
1319
1320 //____________________________________________________________________
1321 const Char_t*
1322 AliFMDEventInspector::CodeString(UInt_t code)
1323 {
1324   static TString s;
1325   s = "";
1326   if (code & kNoEvent)    s.Append("NOEVENT ");
1327   if (code & kNoTriggers) s.Append("NOTRIGGERS ");
1328   if (code & kNoSPD)      s.Append("NOSPD ");
1329   if (code & kNoFMD)      s.Append("NOFMD ");
1330   if (code & kNoVertex)   s.Append("NOVERTEX ");
1331   if (code & kBadVertex)  s.Append("BADVERTEX ");
1332   return s.Data();
1333 }
1334 //____________________________________________________________________
1335 void
1336 AliFMDEventInspector::Print(Option_t*) const
1337 {
1338   //
1339   // Print information
1340   //
1341   //   option Not used
1342   //
1343   char ind[gROOT->GetDirLevel()+1];
1344   for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
1345   ind[gROOT->GetDirLevel()] = '\0';
1346   TString sNN(AliForwardUtil::CenterOfMassEnergyString(fEnergy));
1347   sNN.Strip(TString::kBoth, '0');
1348   sNN.ReplaceAll("GeV", " GeV");
1349   TString field(AliForwardUtil::MagneticFieldString(fField));
1350   field.ReplaceAll("p",  "+");
1351   field.ReplaceAll("m",  "-");
1352   field.ReplaceAll("kG", " kG");
1353
1354   std::cout << std::boolalpha
1355             << ind << ClassName() << ": " << GetName() << '\n'
1356             << ind << " Vertex bins:            " << fVtxAxis.GetNbins() << '\n'
1357             << ind << " Vertex range:           [" << fVtxAxis.GetXmin()
1358             << "," << fVtxAxis.GetXmax() << "]\n"
1359             << ind << " Low flux cut:           " << fLowFluxCut << '\n'
1360             << ind << " Max(delta v_z):         " << fMaxVzErr << " cm\n"
1361             << ind << " Min(nContrib_pileup):   " << fMinPileupContrib << '\n'
1362             << ind << " Min(v-pileup):          " << fMinPileupDistance << '\n'
1363             << ind << " System:                 "
1364             << AliForwardUtil::CollisionSystemString(fCollisionSystem) << '\n'
1365             << ind << " CMS energy per nucleon: " << sNN << '\n'
1366             << ind << " Field:                  " <<  field << '\n'
1367             << ind << " Satellite events:       " << fUseDisplacedVertices<<'\n'
1368             << ind << " Centrality method:      " << fCentMethod << '\n'
1369             << std::noboolalpha;
1370   if (!fCentAxis) { std::cout << std::flush; return; }
1371   Int_t nBin = fCentAxis->GetNbins();
1372   std::cout << ind << " Centrality axis:        " << nBin << " bins"
1373             << std::flush;
1374   for (Int_t i = 0; i < nBin; i++) {
1375     if ((i % 10) == 0) std::cout << '\n' << ind << "  ";
1376     std::cout << std::setw(5) << fCentAxis->GetBinLowEdge(i+1) << '-';
1377   }
1378   std::cout << std::setw(5) << fCentAxis->GetBinUpEdge(nBin) << std::endl;
1379 }
1380
1381
1382 //
1383 // EOF
1384 //
1385