PWG2/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 "AliInputEventHandler.h"
  23 #include "AliTriggerAnalysis.h"
  24 #include "AliPhysicsSelection.h"
  25 #include "AliAODForwardMult.h"
  26 #include "AliForwardUtil.h"
  27 #include <TH1.h>
  28 #include <TList.h>
  29 #include <TDirectory.h>
  30 #include <TROOT.h>
  31 #include <iostream>
  32 #include <iomanip>
  33
  34 //====================================================================
  35 AliFMDEventInspector::AliFMDEventInspector()
  36   : TNamed(),
  37     fHEventsTr(0),
  38     fHEventsTrVtx(0),
  39     fHTriggers(0),
  40     fHType(0),
  41     fLowFluxCut(1000),
  42     fMaxVzErr(0.1),
  43     fList(0),
  44     fEnergy(0),
  45     fField(999),
  46     fCollisionSystem(kUnknown),
  47     fDebug(0)
  48 {
  49   //
  50   // Constructor
  51   //
  52 }
  53
  54 //____________________________________________________________________
  55 AliFMDEventInspector::AliFMDEventInspector(const char* name)
  56   : TNamed("fmdEventInspector", name),
  57     fHEventsTr(0),
  58     fHEventsTrVtx(0),
  59     fHTriggers(0),
  60     fHType(0),
  61     fLowFluxCut(1000),
  62     fMaxVzErr(0.1),
  63     fList(0),
  64     fEnergy(0),
  65     fField(999),
  66     fCollisionSystem(kUnknown),
  67     fDebug(0)
  68 {
  69   //
  70   // Constructor
  71   //
  72   // Parameters:
  73   //   name Name of object
  74   //
  75 }
  76
  77 //____________________________________________________________________
  78 AliFMDEventInspector::AliFMDEventInspector(const AliFMDEventInspector& o)
  79   : TNamed(o),
  80     fHEventsTr(o.fHEventsTr),
  81     fHEventsTrVtx(o.fHEventsTrVtx),
  82     fHTriggers(o.fHTriggers),
  83     fHType(o.fHType),
  84     fLowFluxCut(o.fMaxVzErr),
  85     fMaxVzErr(o.fMaxVzErr),
  86     fList(o.fList),
  87     fEnergy(o.fEnergy),
  88     fField(o.fField),
  89     fCollisionSystem(o.fCollisionSystem),
  90     fDebug(0)
  91 {
  92   //
  93   // Copy constructor
  94   //
  95   // Parameters:
  96   //   o Object to copy from
  97   //
  98 }
  99
 100 //____________________________________________________________________
 101 AliFMDEventInspector::~AliFMDEventInspector()
 102 {
 103   //
 104   // Destructor
 105   //
 106   if (fHEventsTr)    delete fHEventsTr;
 107   if (fHEventsTrVtx) delete fHEventsTrVtx;
 108   if (fHTriggers)    delete fHTriggers;
 109   if (fHType)        delete fHType;
 110   if (fList)         delete fList;
 111 }
 112 //____________________________________________________________________
 113 AliFMDEventInspector&
 114 AliFMDEventInspector::operator=(const AliFMDEventInspector& o)
 115 {
 116   //
 117   // Assignement operator
 118   //
 119   // Parameters:
 120   //   o Object to assign from
 121   //
 122   // Return:
 123   //    Reference to this object
 124   //
 125   TNamed::operator=(o);
 126   fHEventsTr         = o.fHEventsTr;
 127   fHEventsTrVtx      = o.fHEventsTrVtx;
 128   fHTriggers         = o.fHTriggers;
 129   fHType             = o.fHType;
 130   fLowFluxCut        = o.fLowFluxCut;
 131   fMaxVzErr          = o.fMaxVzErr;
 132   fDebug             = o.fDebug;
 133   fList              = (o.fList ? new TList : 0);
 134   fEnergy            = o.fEnergy;
 135   fField             = o.fField;
 136   fCollisionSystem   = o.fCollisionSystem;
 137   if (fList) {
 138     fList->SetName(GetName());
 139     if (fHEventsTr)    fList->Add(fHEventsTr);
 140     if (fHEventsTrVtx) fList->Add(fHEventsTrVtx);
 141     if (fHTriggers)    fList->Add(fHTriggers);
 142     if (fHType)        fList->Add(fHType);
 143   }
 144   return *this;
 145 }
 146
 147 //____________________________________________________________________
 148 Bool_t
 149 AliFMDEventInspector::FetchHistograms(TList* d,
 150                                       TH1I*& hEventsTr,
 151                                       TH1I*& hEventsTrVtx,
 152                                       TH1I*& hTriggers) const
 153 {
 154   //
 155   // Fetch our histograms from the passed list
 156   //
 157   // Parameters:
 158   //   d             Input
 159   //   hEventsTr     On return, pointer to histogram, or null
 160   //   hEventsTrVtx  On return, pointer to histogram, or null
 161   //   hTriggers     On return, pointer to histogram, or null
 162   //
 163   // Return:
 164   //    true on success, false otherwise
 165   //
 166   hEventsTr    = 0;
 167   hEventsTrVtx = 0;
 168   hTriggers    = 0;
 169   TList* dd    = dynamic_cast<TList*>(d->FindObject(GetName()));
 170   if (!dd) return kFALSE;
 171
 172   hEventsTr    = dynamic_cast<TH1I*>(dd->FindObject("nEventsTr"));
 173   hEventsTrVtx = dynamic_cast<TH1I*>(dd->FindObject("nEventsTrVtx"));
 174   hTriggers    = dynamic_cast<TH1I*>(dd->FindObject("triggers"));
 175
 176   if (!hEventsTr || !hEventsTrVtx || !hTriggers) return kFALSE;
 177   return kTRUE;
 178 }
 179 //____________________________________________________________________
 180 void
 181 AliFMDEventInspector::Init(const TAxis& vtxAxis)
 182 {
 183   //
 184   // Initialize the object
 185   //
 186   // Parameters:
 187   //   vtxAxis Vertex axis in use
 188   //
 189   fHEventsTr = new TH1I("nEventsTr", "Number of events w/trigger",
 190                         vtxAxis.GetNbins(),
 191                         vtxAxis.GetXmin(),
 192                         vtxAxis.GetXmax());
 193   fHEventsTr->SetXTitle("v_{z} [cm]");
 194   fHEventsTr->SetYTitle("# of events");
 195   fHEventsTr->SetFillColor(kRed+1);
 196   fHEventsTr->SetFillStyle(3001);
 197   fHEventsTr->SetDirectory(0);
 198   // fHEventsTr->Sumw2();
 199   fList->Add(fHEventsTr);
 200
 201   fHEventsTrVtx = new TH1I("nEventsTrVtx",
 202                            "Number of events w/trigger and vertex",
 203                            vtxAxis.GetNbins(),
 204                            vtxAxis.GetXmin(),
 205                            vtxAxis.GetXmax());
 206   fHEventsTrVtx->SetXTitle("v_{z} [cm]");
 207   fHEventsTrVtx->SetYTitle("# of events");
 208   fHEventsTrVtx->SetFillColor(kBlue+1);
 209   fHEventsTrVtx->SetFillStyle(3001);
 210   fHEventsTrVtx->SetDirectory(0);
 211   // fHEventsTrVtx->Sumw2();
 212   fList->Add(fHEventsTrVtx);
 213
 214
 215   fHTriggers = new TH1I("triggers", "Triggers", 10, 0, 10);
 216   fHTriggers->SetFillColor(kRed+1);
 217   fHTriggers->SetFillStyle(3001);
 218   fHTriggers->SetStats(0);
 219   fHTriggers->SetDirectory(0);
 220   fHTriggers->GetXaxis()->SetBinLabel(kInel   +1,"INEL");
 221   fHTriggers->GetXaxis()->SetBinLabel(kInelGt0+1,"INEL>0");
 222   fHTriggers->GetXaxis()->SetBinLabel(kNSD    +1,"NSD");
 223   fHTriggers->GetXaxis()->SetBinLabel(kEmpty  +1,"Empty");
 224   fHTriggers->GetXaxis()->SetBinLabel(kA      +1,"A");
 225   fHTriggers->GetXaxis()->SetBinLabel(kB      +1,"B");
 226   fHTriggers->GetXaxis()->SetBinLabel(kC      +1,"C");
 227   fHTriggers->GetXaxis()->SetBinLabel(kE      +1,"E");
 228   fHTriggers->GetXaxis()->SetBinLabel(9,         "spare1");
 229   fHTriggers->GetXaxis()->SetBinLabel(10,        "spare2");
 230   fList->Add(fHTriggers);
 231
 232   fHType = new TH1I("type", Form("Event type (cut: SPD mult>%d)",
 233                                  fLowFluxCut), 2, -.5, 1.5);
 234   fHType->SetFillColor(kRed+1);
 235   fHType->SetFillStyle(3001);
 236   fHType->SetStats(0);
 237   fHType->SetDirectory(0);
 238   fHType->GetXaxis()->SetBinLabel(1,"Low-flux");
 239   fHType->GetXaxis()->SetBinLabel(2,"High-flux");
 240   fList->Add(fHType);
 241
 242 }
 243
 244 //____________________________________________________________________
 245 void
 246 AliFMDEventInspector::DefineOutput(TList* dir)
 247 {
 248   //
 249   // Define the output histograms.  These are put in a sub list of the
 250   // passed list.   The histograms are merged before the parent task calls
 251   // AliAnalysisTaskSE::Terminate
 252   //
 253   //   dir Directory to add to
 254   //
 255   fList = new TList;
 256   fList->SetName(GetName());
 257   dir->Add(fList);
 258 }
 259
 260 //____________________________________________________________________
 261 UInt_t
 262 AliFMDEventInspector::Process(const AliESDEvent* event,
 263                               UInt_t&            triggers,
 264                               Bool_t&            lowFlux,
 265                               UShort_t&          ivz,
 266                               Double_t&          vz)
 267 {
 268   //
 269   // Process the event
 270   //
 271   // Parameters:
 272   //   event     Input event
 273   //   triggers  On return, the triggers fired
 274   //   lowFlux   On return, true if the event is considered a low-flux
 275   //                  event (according to the setting of fLowFluxCut)
 276   //   ivz       On return, the found vertex bin (1-based).  A zero
 277   //                  means outside of the defined vertex range
 278   //   vz        On return, the z position of the interaction
 279   //
 280   // Return:
 281   //    0 (or kOk) on success, otherwise a bit mask of error codes
 282   //
 283
 284   // Check that we have an event
 285   if (!event) {
 286     AliWarning("No ESD event found for input event");
 287     return kNoEvent;
 288   }
 289
 290   // Read trigger information from the ESD and store in AOD object
 291   if (!ReadTriggers(event, triggers)) {
 292     if (fDebug > 2) {
 293       AliWarning("Failed to read triggers from ESD"); }
 294     return kNoTriggers;
 295   }
 296
 297   // Check if this is a high-flux event
 298   const AliMultiplicity* testmult = event->GetMultiplicity();
 299   if (!testmult) {
 300     if (fDebug > 3) {
 301       AliWarning("No central multiplicity object found"); }
 302     return kNoSPD;
 303   }
 304   lowFlux = testmult->GetNumberOfTracklets() < fLowFluxCut;
 305
 306   fHType->Fill(lowFlux ? 0 : 1);
 307
 308   // Check the FMD ESD data
 309   if (!event->GetFMDData()) {
 310     if (fDebug > 3) {
 311       AliWarning("No FMD data found in ESD"); }
 312     return kNoFMD;
 313   }
 314
 315   // Get the vertex information
 316   vz          = 0;
 317   Bool_t vzOk = ReadVertex(event, vz);
 318
 319   fHEventsTr->Fill(vz);
 320   if (!vzOk) {
 321     if (fDebug > 3) {
 322       AliWarning("Failed to read vertex from ESD"); }
 323     return kNoVertex;
 324   }
 325   fHEventsTrVtx->Fill(vz);
 326
 327   // Get the vertex bin
 328   ivz = fHEventsTr->GetXaxis()->FindBin(vz);
 329   if (ivz <= 0 || ivz > fHEventsTr->GetXaxis()->GetNbins()) {
 330     if (fDebug > 3) {
 331       AliWarning(Form("Vertex @ %f outside of range [%f,%f]",
 332                       vz, fHEventsTr->GetXaxis()->GetXmin(),
 333                       fHEventsTr->GetXaxis()->GetXmax())); }
 334     ivz = 0;
 335     return kBadVertex;
 336   }
 337   return kOk;
 338 }
 339
 340 //____________________________________________________________________
 341 Bool_t
 342 AliFMDEventInspector::ReadTriggers(const AliESDEvent* esd, UInt_t& triggers)
 343 {
 344   //
 345   // Read the trigger information from the ESD event
 346   //
 347   // Parameters:
 348   //   esd        ESD event
 349   //   triggers   On return, contains the trigger bits
 350   //
 351   // Return:
 352   //    @c true on success, @c false otherwise
 353   //
 354   triggers = 0;
 355
 356   // Get the analysis manager - should always be there
 357   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 358   if (!am) {
 359     AliWarning("No analysis manager defined!");
 360     return kFALSE;
 361   }
 362
 363   // Get the input handler - should always be there
 364   AliInputEventHandler* ih =
 365     static_cast<AliInputEventHandler*>(am->GetInputEventHandler());
 366   if (!ih) {
 367     AliWarning("No input handler");
 368     return kFALSE;
 369   }
 370
 371   // Get the physics selection - add that by using the macro
 372   // AddTaskPhysicsSelection.C
 373   AliPhysicsSelection* ps =
 374     static_cast<AliPhysicsSelection*>(ih->GetEventSelection());
 375   if (!ps) {
 376     AliWarning("No physics selection");
 377     return kFALSE;
 378   }
 379
 380   // Check if this is a collision candidate (INEL)
 381   Bool_t inel = ps->IsCollisionCandidate(esd);
 382   if (inel) {
 383     triggers |= AliAODForwardMult::kInel;
 384     fHTriggers->Fill(kInel+0.5);
 385   }
 386
 387   // IF this is inel, see if we have a tracklet
 388   if (inel) {
 389     const AliMultiplicity* spdmult = esd->GetMultiplicity();
 390     if (!spdmult) {
 391       AliWarning("No SPD multiplicity");
 392     }
 393     else {
 394       Int_t n = spdmult->GetNumberOfTracklets();
 395       for (Int_t j = 0; j < n; j++) {
 396         if(TMath::Abs(spdmult->GetEta(j)) < 1) {
 397           triggers |= AliAODForwardMult::kInelGt0;
 398           fHTriggers->Fill(kInelGt0+.5);
 399           break;
 400         }
 401       }
 402     }
 403   }
 404
 405   // Analyse some trigger stuff
 406   AliTriggerAnalysis ta;
 407   if (ta.IsOfflineTriggerFired(esd, AliTriggerAnalysis::kNSD1)) {
 408     triggers |= AliAODForwardMult::kNSD;
 409     fHTriggers->Fill(kNSD+.5);
 410   }
 411
 412   // Get trigger stuff
 413   TString trigStr = esd->GetFiredTriggerClasses();
 414   if (trigStr.Contains("CBEAMB-ABCE-NOPF-ALL")) {
 415     triggers |= AliAODForwardMult::kEmpty;
 416     fHTriggers->Fill(kEmpty+.5);
 417   }
 418
 419   if (trigStr.Contains("CINT1A-ABCE-NOPF-ALL")) {
 420     triggers |= AliAODForwardMult::kA;
 421     fHTriggers->Fill(kA+.5);
 422   }
 423
 424   if (trigStr.Contains("CINT1B-ABCE-NOPF-ALL")) {
 425     triggers |= AliAODForwardMult::kB;
 426     fHTriggers->Fill(kB+.5);
 427   }
 428
 429
 430   if (trigStr.Contains("CINT1C-ABCE-NOPF-ALL")) {
 431     triggers |= AliAODForwardMult::kC;
 432     fHTriggers->Fill(kC+.5);
 433   }
 434
 435   if (trigStr.Contains("CINT1-E-NOPF-ALL")) {
 436     triggers |= AliAODForwardMult::kE;
 437     fHTriggers->Fill(kE+.5);
 438   }
 439
 440   return kTRUE;
 441 }
 442 //____________________________________________________________________
 443 Bool_t
 444 AliFMDEventInspector::ReadVertex(const AliESDEvent* esd, Double_t& vz)
 445 {
 446   //
 447   // Read the vertex information from the ESD event
 448   //
 449   // Parameters:
 450   //   esd  ESD event
 451   //   vz   On return, the vertex Z position
 452   //
 453   // Return:
 454   //    @c true on success, @c false otherwise
 455   //
 456   vz = 0;
 457   // Get the vertex
 458   const AliESDVertex* vertex = esd->GetPrimaryVertexSPD();
 459   if (!vertex) {
 460     if (fDebug > 2) {
 461       AliWarning("No SPD vertex found in ESD"); }
 462     return kFALSE;
 463   }
 464
 465   // Check that enough tracklets contributed
 466   if(vertex->GetNContributors() <= 0) {
 467     if (fDebug > 2) {
 468       AliWarning(Form("Number of contributors to vertex is %d<=0",
 469                       vertex->GetNContributors())); }
 470     vz = 0;
 471     return kFALSE;
 472   }
 473
 474   // Check that the uncertainty isn't too large
 475   if (vertex->GetZRes() > fMaxVzErr) {
 476     if (fDebug > 2) {
 477       AliWarning(Form("Uncertaintity in Z of vertex is too large %f > %f",
 478                       vertex->GetZRes(), fMaxVzErr)); }
 479     return kFALSE;
 480   }
 481
 482   // Get the z coordiante
 483   vz = vertex->GetZ();
 484   return kTRUE;
 485 }
 486
 487 //____________________________________________________________________
 488 Bool_t
 489 AliFMDEventInspector::ReadRunDetails(const AliESDEvent* esd)
 490 {
 491   //
 492   // Read the collision system, collision energy, and L3 field setting
 493   // from the ESD
 494   //
 495   // Parameters:
 496   //   esd ESD to get information from
 497   //
 498   // Return:
 499   //    true on success, false
 500   //
 501   fCollisionSystem =
 502     AliForwardUtil::ParseCollisionSystem(esd->GetBeamType());
 503   fEnergy          =
 504     AliForwardUtil::ParseCenterOfMassEnergy(fCollisionSystem,
 505                                               2 * esd->GetBeamEnergy());
 506   fField           =
 507     AliForwardUtil::ParseMagneticField(esd->GetMagneticField());
 508
 509   if (fCollisionSystem   == AliForwardUtil::kUnknown ||
 510       fEnergy            <= 0                        ||
 511       TMath::Abs(fField) >  10)
 512     return kFALSE;
 513
 514   return kTRUE;
 515 }
 516
 517 //____________________________________________________________________
 518 void
 519 AliFMDEventInspector::Print(Option_t*) const
 520 {
 521   //
 522   // Print information
 523   //
 524   //   option Not used
 525   //
 526   char ind[gROOT->GetDirLevel()+1];
 527   for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
 528   ind[gROOT->GetDirLevel()] = '\0';
 529   TString sNN(AliForwardUtil::CenterOfMassEnergyString(fEnergy));
 530   sNN.Strip(TString::kBoth, '0');
 531   sNN.ReplaceAll("GeV", " GeV");
 532   TString field(AliForwardUtil::MagneticFieldString(fField));
 533   field.ReplaceAll("p",  "+");
 534   field.ReplaceAll("m",  "-");
 535   field.ReplaceAll("kG", " kG");
 536
 537   std::cout << ind << "AliFMDEventInspector: " << GetName() << '\n'
 538             << ind << " Low flux cut:           " << fLowFluxCut << '\n'
 539             << ind << " Max(delta v_z):         " << fMaxVzErr << " cm\n"
 540             << ind << " System:                 "
 541             << AliForwardUtil::CollisionSystemString(fCollisionSystem) << '\n'
 542             << ind << " CMS energy per nucleon: " << sNN << '\n'
 543             << ind << " Field:                  " <<  field << std::endl;
 544 }
 545
 546
 547 //
 548 // EOF
 549 //
 550