PWGLF/FORWARD/analysis2/AliForwardMultiplicityBase.cxx

   1 //====================================================================
   2 //
   3 // Base class for classes that calculate the multiplicity in the
   4 // forward regions event-by-event
   5 //
   6 // Inputs:
   7 //   - AliESDEvent
   8 //
   9 // Outputs:
  10 //   - AliAODForwardMult
  11 //
  12 // Histograms
  13 //
  14 // Corrections used
  15 #include "AliForwardMultiplicityBase.h"
  16 #include "AliForwardCorrectionManager.h"
  17 #include "AliForwardUtil.h"
  18 #include "AliFMDCorrELossFit.h"
  19 #include "AliLog.h"
  20 #include "AliAODHandler.h"
  21 #include "AliInputEventHandler.h"
  22 #include "AliAnalysisManager.h"
  23 #include "AliFMDEventInspector.h"
  24 #include "AliFMDSharingFilter.h"
  25 #include "AliFMDDensityCalculator.h"
  26 #include "AliFMDCorrector.h"
  27 #include "AliFMDHistCollector.h"
  28 #include "AliFMDEventPlaneFinder.h"
  29 #include "AliESDEvent.h"
  30 #include <TROOT.h>
  31 #include <TSystem.h>
  32 #include <TAxis.h>
  33 #include <TProfile.h>
  34 #include <THStack.h>
  35 #include <iostream>
  36 #include <iomanip>
  37
  38 //====================================================================
  39 AliForwardMultiplicityBase::AliForwardMultiplicityBase(const char* name)
  40   : AliAnalysisTaskSE(name),
  41     fEnableLowFlux(false),
  42     fFirstEvent(true),
  43     fStorePerRing(false),
  44     fList(0),
  45     fHData(0),
  46     fHistos(),
  47     fAODFMD(false),
  48     fAODEP(false),
  49     fRingSums(),
  50     fDoTiming(false),
  51     fHTiming(0),
  52     fCorrManager(0)
  53 {
  54   DGUARD(fDebug, 3,"Named CTOR of AliForwardMultiplicityBase %s",name);
  55   // Set our persistent pointer
  56   fCorrManager = &AliForwardCorrectionManager::Instance();
  57   fBranchNames =
  58     "ESD:AliESDRun.,AliESDHeader.,AliMultiplicity.,"
  59     "AliESDFMD.,SPDVertex.,PrimaryVertex.";
  60
  61   DefineOutput(1, TList::Class());
  62   DefineOutput(2, TList::Class());
  63 }
  64
  65 //____________________________________________________________________
  66 AliForwardMultiplicityBase&
  67 AliForwardMultiplicityBase::operator=(const AliForwardMultiplicityBase& o)
  68 {
  69   DGUARD(fDebug,2,"Assignment to AliForwardMultiplicityBase");
  70   if (&o == this) return *this;
  71   fEnableLowFlux = o.fEnableLowFlux;
  72   fFirstEvent    = o.fFirstEvent;
  73   fCorrManager   = o.fCorrManager;
  74   fHData             = o.fHData;
  75   fHistos            = o.fHistos;
  76   fAODFMD            = o.fAODFMD;
  77   fAODEP             = o.fAODEP;
  78   fRingSums          = o.fRingSums;
  79   fList              = o.fList;
  80   fStorePerRing      = o.fStorePerRing;
  81   fDoTiming          = o.fDoTiming;
  82   fHTiming           = o.fHTiming;
  83   return *this;
  84 }
  85
  86 //____________________________________________________________________
  87 void
  88 AliForwardMultiplicityBase::SetDebug(Int_t dbg)
  89 {
  90   //
  91   // Set debug level
  92   //
  93   // Parameters:
  94   //    dbg debug level
  95   //
  96   GetEventInspector()   .SetDebug(dbg);
  97   GetSharingFilter()    .SetDebug(dbg);
  98   GetDensityCalculator().SetDebug(dbg);
  99   GetCorrections()      .SetDebug(dbg);
 100   GetHistCollector()    .SetDebug(dbg);
 101   GetEventPlaneFinder() .SetDebug(dbg);
 102 }
 103 //____________________________________________________________________
 104 Bool_t
 105 AliForwardMultiplicityBase::Configure(const char* macro)
 106 {
 107   // --- Configure the task ------------------------------------------
 108   TString macroPath(gROOT->GetMacroPath());
 109   if (!macroPath.Contains("$(ALICE_ROOT)/PWGLF/FORWARD/analysis2")) {
 110     macroPath.Append(":$(ALICE_ROOT)/PWGLF/FORWARD/analysis2");
 111     gROOT->SetMacroPath(macroPath);
 112   }
 113   TString mac(macro);
 114   if (mac.EqualTo("-default-"))
 115     mac = "$(ALICE_ROOT)/PWGLF/FORWARD/analysis2/ForwardAODConfig.C";
 116   const char* config = gSystem->Which(gROOT->GetMacroPath(), mac.Data());
 117   if (!config) {
 118     AliWarningF("%s not found in %s", macro, gROOT->GetMacroPath());
 119     return false;
 120   }
 121
 122   AliInfoF("Loading configuration of '%s' from %s",  ClassName(), config);
 123   gROOT->Macro(Form("%s((AliForwardMultiplicityBase*)%p)", config, this));
 124   delete config;
 125
 126  return true;
 127 }
 128
 129 //____________________________________________________________________
 130 void
 131 AliForwardMultiplicityBase::UserCreateOutputObjects()
 132 {
 133   //
 134   // Create output objects
 135   //
 136   //
 137   DGUARD(fDebug,1,"Create user ouput");
 138   fList = new TList;
 139   fList->SetOwner();
 140
 141   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 142   AliAODHandler*      ah =
 143     dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
 144   //if (!ah) AliFatal("No AOD output handler set in analysis manager");
 145   if (ah)  CreateBranches(ah);
 146
 147   GetEventInspector()   .CreateOutputObjects(fList);
 148   GetSharingFilter()    .CreateOutputObjects(fList);
 149   GetDensityCalculator().CreateOutputObjects(fList);
 150   GetCorrections()      .CreateOutputObjects(fList);
 151   GetHistCollector()    .CreateOutputObjects(fList);
 152   GetEventPlaneFinder() .CreateOutputObjects(fList);
 153
 154   if (fDebug > 1) fDoTiming = true;
 155   if (fDoTiming) {
 156     fHTiming = new TProfile("timing", "Timing of task",
 157                             kTimingTotal, 0.5, kTimingTotal+.5);
 158     fHTiming->SetDirectory(0);
 159     fHTiming->SetFillColor(kRed+1);
 160     fHTiming->SetFillStyle(3001);
 161     fHTiming->SetMarkerStyle(20);
 162     fHTiming->SetMarkerColor(kBlack);
 163     fHTiming->SetLineColor(kBlack);
 164     fHTiming->SetXTitle("Part");
 165     fHTiming->SetYTitle("#LTt_{part}#GT [CPU]");
 166     fHTiming->SetStats(0);
 167     TAxis* xaxis = fHTiming->GetXaxis();
 168     xaxis->SetBinLabel(kTimingEventInspector,
 169                        GetEventInspector()   .GetName());
 170     xaxis->SetBinLabel(kTimingSharingFilter,
 171                        GetSharingFilter()    .GetName());
 172     xaxis->SetBinLabel(kTimingDensityCalculator,
 173                        GetDensityCalculator().GetName());
 174     xaxis->SetBinLabel(kTimingCorrections,
 175                        GetCorrections()      .GetName());
 176     xaxis->SetBinLabel(kTimingHistCollector,
 177                        GetHistCollector()    .GetName());
 178     xaxis->SetBinLabel(kTimingEventPlaneFinder,
 179                        GetEventPlaneFinder() .GetName());
 180     xaxis->SetBinLabel(kTimingTotal, "Total");
 181     fList->Add(fHTiming);
 182   }
 183   PostData(1, fList);
 184 }
 185 //____________________________________________________________________
 186 void
 187 AliForwardMultiplicityBase::CreateBranches(AliAODHandler* ah)
 188 {
 189   TObject* obj = &fAODFMD;
 190   ah->AddBranch("AliAODForwardMult", &obj);
 191   TObject* epobj = &fAODEP;
 192   ah->AddBranch("AliAODForwardEP", &epobj);
 193
 194   TAxis tmp(1, 0, 1);
 195   fHistos.Init(tmp);
 196
 197   if (!fStorePerRing) return;
 198
 199   AliWarning("Per-ring histograms in AOD\n"
 200              "*********************************************************\n"
 201              "* For each event 5 additional 2D histogram are stored   *\n"
 202              "* in separate branches of the AODs.  This will increase *\n"
 203              "* the size of the AODs - proceed with caution           *\n"
 204              "*********************************************************");
 205   TObject* hists[] = { fHistos.fFMD1i,
 206                        fHistos.fFMD2i, fHistos.fFMD2o,
 207                        fHistos.fFMD3i, fHistos.fFMD3o };
 208   for (Int_t i = 0; i < 5; i++) {
 209     ah->AddBranch("TH2D", &(hists[i]));
 210   }
 211 }
 212
 213 //____________________________________________________________________
 214 Bool_t
 215 AliForwardMultiplicityBase::SetupForData()
 216 {
 217   //
 218   // Initialise the sub objects and stuff.  Called on first event
 219   //
 220   //
 221   DGUARD(fDebug,1,"Initialize sub-algorithms");
 222   const TAxis* pe = 0;
 223   const TAxis* pv = 0;
 224
 225   Bool_t  mc  = this->IsA()->InheritsFrom("AliForwardMCMultiplicityTask");
 226   Bool_t  sat = false; // GetEventInspector().IsUseDisplacedVertices();
 227   if (!ReadCorrections(pe,pv,mc,sat)) return false;
 228
 229   InitMembers(pe,pv);
 230
 231   GetEventInspector()   .SetupForData(*pv);
 232   GetSharingFilter()    .SetupForData(*pe);
 233   GetDensityCalculator().SetupForData(*pe);
 234   GetCorrections()      .SetupForData(*pe);
 235   GetHistCollector()    .SetupForData(*pv,*pe);
 236   GetEventPlaneFinder() .SetupForData(*pe);
 237
 238   fAODFMD.SetBit(AliAODForwardMult::kSecondary,
 239                  GetCorrections().IsUseSecondaryMap());
 240   fAODFMD.SetBit(AliAODForwardMult::kVertexBias,
 241                  GetCorrections().IsUseVertexBias());
 242   fAODFMD.SetBit(AliAODForwardMult::kAcceptance,
 243                  GetCorrections().IsUseAcceptance());
 244   fAODFMD.SetBit(AliAODForwardMult::kMergingEfficiency,
 245                  GetCorrections().IsUseMergingEfficiency());
 246   fAODFMD.SetBit(AliAODForwardMult::kSum,
 247                  GetHistCollector().GetMergeMethod() ==
 248                  AliFMDHistCollector::kSum);
 249
 250   this->Print("R");
 251   return true;
 252 }
 253
 254 //____________________________________________________________________
 255 void
 256 AliForwardMultiplicityBase::InitMembers(const TAxis* pe, const TAxis* /*pv*/)
 257 {
 258   fHistos.ReInit(*pe);
 259   fAODFMD.Init(*pe);
 260   fAODEP.Init(*pe);
 261   fRingSums.Init(*pe);
 262
 263   fHData = static_cast<TH2D*>(fAODFMD.GetHistogram().Clone("d2Ndetadphi"));
 264   fHData->SetStats(0);
 265   fHData->SetDirectory(0);
 266   fList->Add(fHData);
 267
 268   TList* rings = new TList;
 269   rings->SetName("ringSums");
 270   rings->SetOwner();
 271   fList->Add(rings);
 272
 273   rings->Add(fRingSums.Get(1, 'I'));
 274   rings->Add(fRingSums.Get(2, 'I'));
 275   rings->Add(fRingSums.Get(2, 'O'));
 276   rings->Add(fRingSums.Get(3, 'I'));
 277   rings->Add(fRingSums.Get(3, 'O'));
 278   fRingSums.Get(1, 'I')->SetMarkerColor(AliForwardUtil::RingColor(1, 'I'));
 279   fRingSums.Get(2, 'I')->SetMarkerColor(AliForwardUtil::RingColor(2, 'I'));
 280   fRingSums.Get(2, 'O')->SetMarkerColor(AliForwardUtil::RingColor(2, 'O'));
 281   fRingSums.Get(3, 'I')->SetMarkerColor(AliForwardUtil::RingColor(3, 'I'));
 282   fRingSums.Get(3, 'O')->SetMarkerColor(AliForwardUtil::RingColor(3, 'O'));
 283 }
 284
 285 //____________________________________________________________________
 286 Bool_t
 287 AliForwardMultiplicityBase::CheckCorrections(UInt_t what) const
 288 {
 289   //
 290   // Check if all needed corrections are there and accounted for.  If not,
 291   // do a Fatal exit
 292   //
 293   // Parameters:
 294   //    what Which corrections is needed
 295   //
 296   // Return:
 297   //    true if all present, false otherwise
 298   //
 299   DGUARD(fDebug,1,"Checking corrections 0x%x", what);
 300
 301   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
 302   // Check that we have the energy loss fits, needed by
 303   //   AliFMDSharingFilter
 304   //   AliFMDDensityCalculator
 305   if (what & AliForwardCorrectionManager::kELossFits) {
 306     if (!fcm.GetELossFit()) {
 307       AliFatal(Form("No energy loss fits"));
 308       return false;
 309     }
 310     // Force this here so we select the proper quality
 311     const AliFMDCorrELossFit* fits = fcm.GetELossFit();
 312     fits->CacheBins(GetDensityCalculator().GetMinQuality());
 313   }
 314
 315   // Check that we have the double hit correction - (optionally) used by
 316   //  AliFMDDensityCalculator
 317   if (what & AliForwardCorrectionManager::kDoubleHit && !fcm.GetDoubleHit()) {
 318     AliFatal("No double hit corrections");
 319     return false;
 320   }
 321   // Check that we have the secondary maps, needed by
 322   //   AliFMDCorrector
 323   //   AliFMDHistCollector
 324   if (what & AliForwardCorrectionManager::kSecondaryMap &&
 325       !fcm.GetSecondaryMap()) {
 326     AliFatal("No secondary corrections");
 327     return false;
 328   }
 329   // Check that we have the vertex bias correction, needed by
 330   //   AliFMDCorrector
 331   if (what & AliForwardCorrectionManager::kVertexBias &&
 332       !fcm.GetVertexBias()) {
 333     AliFatal("No event vertex bias corrections");
 334     return false;
 335   }
 336   // Check that we have the merging efficiencies, optionally used by
 337   //   AliFMDCorrector
 338   if (what & AliForwardCorrectionManager::kMergingEfficiency &&
 339       !fcm.GetMergingEfficiency()) {
 340     AliFatal("No merging efficiencies");
 341     return false;
 342   }
 343   // Check that we have the acceptance correction, needed by
 344   //   AliFMDCorrector
 345   if (what & AliForwardCorrectionManager::kAcceptance &&
 346       !fcm.GetAcceptance()) {
 347     AliFatal("No acceptance corrections");
 348     return false;
 349   }
 350   return true;
 351 }
 352 //____________________________________________________________________
 353 Bool_t
 354 AliForwardMultiplicityBase::ReadCorrections(const TAxis*& pe,
 355                                             const TAxis*& pv,
 356                                             Bool_t        mc,
 357                                             Bool_t        sat)
 358 {
 359   //
 360   // Read corrections
 361   //
 362   //
 363
 364   UInt_t what = AliForwardCorrectionManager::kAll;
 365   if (!fEnableLowFlux)
 366     what ^= AliForwardCorrectionManager::kDoubleHit;
 367   if (!GetCorrections().IsUseVertexBias())
 368     what ^= AliForwardCorrectionManager::kVertexBias;
 369   if (!GetCorrections().IsUseAcceptance())
 370     what ^= AliForwardCorrectionManager::kAcceptance;
 371   if (!GetCorrections().IsUseMergingEfficiency())
 372     what ^= AliForwardCorrectionManager::kMergingEfficiency;
 373   DGUARD(fDebug,1,"Read corrections 0x%x", what);
 374
 375   AliForwardCorrectionManager& fcm = AliForwardCorrectionManager::Instance();
 376   if (!fcm.Init(GetEventInspector().GetRunNumber(),
 377                 GetEventInspector().GetCollisionSystem(),
 378                 GetEventInspector().GetEnergy(),
 379                 GetEventInspector().GetField(),
 380                 mc,
 381                 sat,
 382                 what,
 383                 false)) {
 384     AliWarning("Failed to read in some corrections, making task zombie");
 385     return false;
 386   }
 387   if (!CheckCorrections(what)) return false;
 388
 389   // Sett our persistency pointer
 390   // fCorrManager = &fcm;
 391
 392   // Get the eta axis from the secondary maps - if read in
 393   if (!pe) {
 394     pe = fcm.GetEtaAxis();
 395     if (!pe) AliFatal("No eta axis defined");
 396   }
 397   // Get the vertex axis from the secondary maps - if read in
 398   if (!pv) {
 399     pv = fcm.GetVertexAxis();
 400     if (!pv) AliFatal("No vertex axis defined");
 401   }
 402
 403   return true;
 404 }
 405 //____________________________________________________________________
 406 AliESDEvent*
 407 AliForwardMultiplicityBase::GetESDEvent()
 408 {
 409   //
 410   // Get the ESD event. IF this is the first event, initialise
 411   //
 412   DGUARD(fDebug,1,"Get the ESD event");
 413   if (IsZombie()) return 0;
 414   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
 415   if (!esd) {
 416     AliWarning("No ESD event found for input event");
 417     return 0;
 418   }
 419
 420   // --- Load the data -----------------------------------------------
 421   LoadBranches();
 422
 423   // On the first event, initialize the parameters
 424   if (fFirstEvent && esd->GetESDRun()) {
 425     GetEventInspector().ReadRunDetails(esd);
 426
 427     AliInfo(Form("Initializing with parameters from the ESD:\n"
 428                  "         AliESDEvent::GetBeamEnergy()   ->%f\n"
 429                  "         AliESDEvent::GetBeamType()     ->%s\n"
 430                  "         AliESDEvent::GetCurrentL3()    ->%f\n"
 431                  "         AliESDEvent::GetMagneticField()->%f\n"
 432                  "         AliESDEvent::GetRunNumber()    ->%d",
 433                  esd->GetBeamEnergy(),
 434                  esd->GetBeamType(),
 435                  esd->GetCurrentL3(),
 436                  esd->GetMagneticField(),
 437                  esd->GetRunNumber()));
 438
 439     fFirstEvent = false;
 440
 441     GetEventPlaneFinder().SetRunNumber(esd->GetRunNumber());
 442     if (!SetupForData()) {
 443       AliError("Failed to initialize sub-algorithms, making this a zombie");
 444       esd = 0; // Make sure we do nothing on this event
 445       Info("GetESDEvent", "ESD event pointer %p", esd);
 446       SetZombie(true);
 447       // return 0;
 448     }
 449   }
 450   return esd;
 451 }
 452 //____________________________________________________________________
 453 void
 454 AliForwardMultiplicityBase::MarkEventForStore() const
 455 {
 456   // Make sure the AOD tree is filled
 457   DGUARD(fDebug,3,"Mark AOD event for storage");
 458   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 459   AliAODHandler*      ah =
 460     dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
 461   if (ah)
 462         ah->SetFillAOD(kTRUE);
 463
 464 }
 465
 466 //____________________________________________________________________
 467 void
 468 AliForwardMultiplicityBase::Terminate(Option_t*)
 469 {
 470   //
 471   // End of job
 472   //
 473   // Parameters:
 474   //    option Not used
 475   //
 476   DGUARD(fDebug,1,"Processing the merged results");
 477
 478   TList* list = dynamic_cast<TList*>(GetOutputData(1));
 479   if (!list) {
 480     AliError(Form("No output list defined (%p)", GetOutputData(1)));
 481     if (GetOutputData(1)) GetOutputData(1)->Print();
 482     return;
 483   }
 484
 485   TList* output = new TList;
 486   output->SetName(Form("%sResults", GetName()));
 487   output->SetOwner();
 488
 489   Double_t nTr = 0, nTrVtx = 0, nAcc = 0;
 490   MakeSimpledNdeta(list, output, nTr, nTrVtx, nAcc);
 491
 492   EstimatedNdeta(list, output);
 493
 494   GetSharingFilter()    .Terminate(list,output,Int_t(nTr));
 495   GetDensityCalculator().Terminate(list,output,Int_t(nTrVtx));
 496   GetCorrections()      .Terminate(list,output,Int_t(nTrVtx));
 497
 498   TProfile* timing = static_cast<TProfile*>(list->FindObject("timing"));
 499   if (timing) {
 500     TProfile* p = static_cast<TProfile*>(timing->Clone());
 501     p->SetDirectory(0);
 502     p->Scale(100. / p->GetBinContent(p->GetNbinsX()));
 503     p->SetYTitle("#LTt_{part}#GT/#LTt_{total}#GT [%]");
 504     p->SetTitle("Relative timing of task");
 505     output->Add(p);
 506   }
 507   PostData(2, output);
 508 }
 509
 510 //____________________________________________________________________
 511 void
 512 AliForwardMultiplicityBase::EstimatedNdeta(const TList* input,
 513                                            TList*       output) const
 514 {
 515   MakeRingdNdeta(input, "ringSums", output, "ringResults");
 516 }
 517
 518 //____________________________________________________________________
 519 Bool_t
 520 AliForwardMultiplicityBase::MakeSimpledNdeta(const TList* input,
 521                                              TList*       output,
 522                                              Double_t&    nTr,
 523                                              Double_t&    nTrVtx,
 524                                              Double_t&    nAcc)
 525 {
 526   // Get our histograms from the container
 527   TH1I* hEventsTr    = 0;
 528   TH1I* hEventsTrVtx = 0;
 529   TH1I* hEventsAcc   = 0;
 530   TH1I* hTriggers    = 0;
 531   if (!GetEventInspector().FetchHistograms(input,
 532                                            hEventsTr,
 533                                            hEventsTrVtx,
 534                                            hEventsAcc,
 535                                            hTriggers)) {
 536     AliError(Form("Didn't get histograms from event selector "
 537                   "(hEventsTr=%p,hEventsTrVtx=%p,hEventsAcc=%p,hTriggers=%p)",
 538                   hEventsTr, hEventsTrVtx, hEventsAcc, hTriggers));
 539     input->ls();
 540     return false;
 541   }
 542   nTr             = hEventsTr->Integral();
 543   nTrVtx          = hEventsTrVtx->Integral();
 544   nAcc            = hEventsAcc->Integral();
 545   Double_t vtxEff = nTrVtx / nTr;
 546   TH2D*   hData   = static_cast<TH2D*>(input->FindObject("d2Ndetadphi"));
 547   if (!hData) {
 548     AliError(Form("Couldn't get our summed histogram from output "
 549                   "list %s (d2Ndetadphi=%p)", input->GetName(), hData));
 550     input->ls();
 551     return false;
 552   }
 553
 554   Int_t nY      = hData->GetNbinsY();
 555   TH1D* dNdeta  = hData->ProjectionX("dNdeta",  1,     nY, "e");
 556   TH1D* dNdeta_ = hData->ProjectionX("dNdeta_", 1,     nY, "e");
 557   TH1D* norm    = hData->ProjectionX("norm",    0,     0,  "");
 558   TH1D* phi     = hData->ProjectionX("phi",     nY+1,  nY+1,  "");
 559   dNdeta->SetTitle("dN_{ch}/d#eta in the forward regions");
 560   dNdeta->SetYTitle("#frac{1}{N}#frac{dN_{ch}}{d#eta}");
 561   dNdeta->SetMarkerColor(kRed+1);
 562   dNdeta->SetMarkerStyle(20);
 563   dNdeta->SetDirectory(0);
 564
 565   dNdeta_->SetTitle("dN_{ch}/d#eta in the forward regions");
 566   dNdeta_->SetYTitle("#frac{1}{N}#frac{dN_{ch}}{d#eta}");
 567   dNdeta_->SetMarkerColor(kMagenta+1);
 568   dNdeta_->SetMarkerStyle(21);
 569   dNdeta_->SetDirectory(0);
 570
 571   norm->SetTitle("Normalization to  #eta coverage");
 572   norm->SetYTitle("#eta coverage");
 573   norm->SetLineColor(kBlue+1);
 574   norm->SetMarkerColor(kBlue+1);
 575   norm->SetMarkerStyle(21);
 576   norm->SetFillColor(kBlue+1);
 577   norm->SetFillStyle(3005);
 578   norm->SetDirectory(0);
 579
 580   phi->SetTitle("Normalization to  #phi acceptance");
 581   phi->SetYTitle("#phi acceptance");
 582   phi->SetLineColor(kGreen+1);
 583   phi->SetMarkerColor(kGreen+1);
 584   phi->SetMarkerStyle(20);
 585   phi->SetFillColor(kGreen+1);
 586   phi->SetFillStyle(3004);
 587   // phi->Scale(1. / nAcc);
 588   phi->SetDirectory(0);
 589
 590   // dNdeta->Divide(norm);
 591   dNdeta->Divide(phi);
 592   dNdeta->SetStats(0);
 593   dNdeta->Scale(vtxEff, "width");
 594
 595   dNdeta_->Divide(norm);
 596   dNdeta_->SetStats(0);
 597   dNdeta_->Scale(vtxEff, "width");
 598
 599   output->Add(dNdeta);
 600   output->Add(dNdeta_);
 601   output->Add(norm);
 602   output->Add(phi);
 603
 604   return true;
 605 }
 606
 607
 608 //____________________________________________________________________
 609 void
 610 AliForwardMultiplicityBase::MakeRingdNdeta(const TList* input,
 611                                            const char*  inName,
 612                                            TList*       output,
 613                                            const char*  outName,
 614                                            Int_t        style) const
 615 {
 616   // Make dN/deta for each ring found in the input list.
 617   //
 618   // A stack of all the dN/deta is also made for easy drawing.
 619   //
 620   // Note, that the distributions are normalised to the number of
 621   // observed events only - they should be corrected for
 622   DGUARD(fDebug,3,"Make first-shot ring dN/deta");
 623
 624   if (!input) return;
 625   TList* list = static_cast<TList*>(input->FindObject(inName));
 626   if (!list) {
 627     AliWarning(Form("No list %s found in %s", inName, input->GetName()));
 628     return;
 629   }
 630
 631   TList* out = new TList;
 632   out->SetName(outName);
 633   out->SetOwner();
 634   output->Add(out);
 635
 636   THStack*     dndetaRings = new THStack("all", "dN/d#eta per ring");
 637   const char*  names[]     = { "FMD1I", "FMD2I", "FMD2O", "FMD3I", "FMD3O", 0 };
 638   const char** ptr         = names;
 639
 640   while (*ptr) {
 641     TList* thisList = new TList;
 642     thisList->SetOwner();
 643     thisList->SetName(*ptr);
 644     out->Add(thisList);
 645
 646     TH2D* h = static_cast<TH2D*>(list->FindObject(Form("%s_cache", *ptr)));
 647     if (!h) {
 648       AliWarning(Form("Didn't find %s_cache in %s", *ptr, list->GetName()));
 649       ptr++;
 650       continue;
 651     }
 652     TH2D* sumPhi = static_cast<TH2D*>(h->Clone("sum_phi"));
 653     sumPhi->SetDirectory(0);
 654     thisList->Add(sumPhi);
 655
 656     TH2D* sumEta = static_cast<TH2D*>(h->Clone("sum_eta"));
 657     sumEta->SetDirectory(0);
 658     thisList->Add(sumEta);
 659
 660     Int_t nY   = sumEta->GetNbinsY();
 661     TH1D* etaCov =static_cast<TH1D*>(h->ProjectionX("etaCov", 0,    0,    ""));
 662     TH1D* phiAcc =static_cast<TH1D*>(h->ProjectionX("phiAcc", nY+1, nY+1, ""));
 663
 664     etaCov->SetTitle("Normalization to #eta coverage");
 665     etaCov->SetYTitle("#eta coverage");
 666     etaCov->SetMarkerColor(kBlue+1);
 667     etaCov->SetFillColor(kBlue+1);
 668     etaCov->SetFillStyle(3005);
 669     etaCov->SetDirectory(0);
 670
 671     phiAcc->SetTitle("Normalization to #phi acceptance");
 672     phiAcc->SetYTitle("#phi acceptance");
 673     phiAcc->SetMarkerColor(kGreen+1);
 674     phiAcc->SetFillColor(kGreen+1);
 675     phiAcc->SetFillStyle(3004);
 676     // phiAcc->Scale(1. / nAcc);
 677     phiAcc->SetDirectory(0);
 678
 679     // Double_t s = (etaCov->GetMaximum() > 0 ? 1. / etaCov->GetMaximum() : 1);
 680     for (Int_t i = 1; i <= sumEta->GetNbinsX(); i++) {
 681       for (Int_t j = 1; j <= nY; j++) {
 682         Double_t c = sumEta->GetBinContent(i, j);
 683         Double_t e = sumEta->GetBinError(i, j);
 684         Double_t a = etaCov->GetBinContent(i);
 685         Double_t p = phiAcc->GetBinContent(i);
 686         // Double_t t = p; // * a
 687         sumEta->SetBinContent(i, j, a <= 0 ? 0 : c / a);
 688         sumEta->SetBinError(  i, j, a <= 0 ? 0 : e / a);
 689         sumPhi->SetBinContent(i, j, p <= 0 ? 0 : c / p);
 690         sumPhi->SetBinError(  i, j, p <= 0 ? 0 : e / p);
 691       }
 692     }
 693     // etaCov->Scale(s);
 694     // phiAcc->Scale(s);
 695
 696     TH1D* resPhi  =static_cast<TH1D*>(sumPhi->ProjectionX("dndeta_phi",
 697                                                           1,nY,"e"));
 698     resPhi->SetMarkerStyle(style);
 699     resPhi->SetDirectory(0);
 700     resPhi->Scale(1, "width");
 701
 702     TH1D* resEta  =static_cast<TH1D*>(sumEta->ProjectionX("dndeta_eta",
 703                                                           1,nY,"e"));
 704     resEta->SetMarkerStyle(style);
 705     resEta->SetDirectory(0);
 706     resEta->Scale(1, "width");
 707
 708     thisList->Add(resEta);
 709     thisList->Add(etaCov);
 710     thisList->Add(resPhi);
 711     thisList->Add(phiAcc);
 712     dndetaRings->Add(resPhi);
 713     ptr++;
 714   }
 715   out->Add(dndetaRings);
 716 }
 717 //____________________________________________________________________
 718 void
 719 AliForwardMultiplicityBase::Print(Option_t* option) const
 720 {
 721   //
 722   // Print information
 723   //
 724   // Parameters:
 725   //    option Not used
 726   //
 727
 728   std::cout << ClassName() << ": " << GetName() << "\n"
 729             << "  Enable low flux code:   " << (fEnableLowFlux ? "yes" : "no")
 730             << "\n"
 731             << "  Store per-ring hists:   " << (fStorePerRing ? "yes" : "no")
 732             << "\n"
 733             << "  Off-line trigger mask:  0x"
 734             << std::hex     << std::setfill('0')
 735             << std::setw (8) << fOfflineTriggerMask
 736             << std::dec     << std::setfill (' ') << "\n"
 737             << "  Make timing histogram:  " << std::boolalpha
 738             << fDoTiming << std::noboolalpha << std::endl;
 739   gROOT->IncreaseDirLevel();
 740   if (fCorrManager) fCorrManager->Print(option);
 741   else
 742     std::cout << "  Correction manager not set yet" << std::endl;
 743   GetEventInspector()   .Print(option);
 744   GetSharingFilter()    .Print(option);
 745   GetDensityCalculator().Print(option);
 746   GetCorrections()      .Print(option);
 747   GetHistCollector()    .Print(option);
 748   GetEventPlaneFinder() .Print(option);
 749   gROOT->DecreaseDirLevel();
 750 }
 751
 752
 753 //
 754 // EOF
 755 //