PWG2/FORWARD/analysis2/AliCentralMultiplicityTask.cxx

   1 //====================================================================
   2 //
   3 // Base class for classes that calculate the multiplicity in the
   4 // central region event-by-event
   5 //
   6 // Inputs:
   7 //   - AliESDEvent
   8 //
   9 // Outputs:
  10 //   - AliAODCentralMult
  11 //
  12 // Histograms
  13 //
  14 // Corrections used
  15 #include "AliCentralMultiplicityTask.h"
  16 #include "AliAODForwardMult.h"
  17 #include "AliForwardUtil.h"
  18 #include "AliLog.h"
  19 #include "AliAODHandler.h"
  20 #include "AliAnalysisManager.h"
  21 #include "AliESDEvent.h"
  22 #include "AliMultiplicity.h"
  23 #include <TROOT.h>
  24 #include <TFile.h>
  25 #include <TError.h>
  26 #include <TSystem.h>
  27 #include <iostream>
  28 #include <iomanip>
  29
  30 //====================================================================
  31 AliCentralMultiplicityTask::AliCentralMultiplicityTask(const char* name)
  32   : AliAnalysisTaskSE(name),
  33     fInspector("centralEventInspector"),
  34     fData(0),
  35     fList(0),
  36     fAODCentral(kFALSE),
  37     fManager(),
  38     fUseSecondary(true),
  39     fUseAcceptance(true),
  40     fFirstEventSeen(false),
  41     fIvz(0),
  42     fNClusterTracklet(0),
  43     fClusterPerTracklet(0),
  44     fNCluster(0),
  45     fNTracklet(0),
  46     fEtaMin(0),
  47     fEtaMax(0)
  48 {
  49   //
  50   // Constructor
  51   //
  52   DefineOutput(1, TList::Class());
  53   fBranchNames =
  54     "ESD:AliESDRun.,AliESDHeader.,AliMultiplicity.,"
  55     "SPDVertex.,PrimaryVertex.";
  56 }
  57 //____________________________________________________________________
  58 AliCentralMultiplicityTask::AliCentralMultiplicityTask()
  59   : AliAnalysisTaskSE(),
  60     fInspector(),
  61     fData(0),
  62     fList(0),
  63     fAODCentral(),
  64     fManager(),
  65     fUseSecondary(true),
  66     fUseAcceptance(true),
  67     fFirstEventSeen(false),
  68     fIvz(0),
  69     fNClusterTracklet(0),
  70     fClusterPerTracklet(0),
  71     fNCluster(0),
  72     fNTracklet(0),
  73     fEtaMin(0),
  74     fEtaMax(0)
  75 {
  76   //
  77   // Constructor
  78   //
  79 }
  80 //____________________________________________________________________
  81 AliCentralMultiplicityTask::AliCentralMultiplicityTask(const AliCentralMultiplicityTask& o)
  82   : AliAnalysisTaskSE(o),
  83     fInspector(o.fInspector),
  84     fData(o.fData),
  85     fList(o.fList),
  86     fAODCentral(o.fAODCentral),
  87     fManager(o.fManager),
  88     fUseSecondary(o.fUseSecondary),
  89     fUseAcceptance(o.fUseAcceptance),
  90     fFirstEventSeen(o.fFirstEventSeen),
  91     fIvz(0),
  92     fNClusterTracklet(o.fNClusterTracklet),
  93     fClusterPerTracklet(o.fClusterPerTracklet),
  94     fNCluster(o.fNCluster),
  95     fNTracklet(o.fNTracklet),
  96     fEtaMin(o.fEtaMin),
  97     fEtaMax(o.fEtaMax)
  98 {
  99   //
 100   // Copy constructor
 101   //
 102 }
 103 //____________________________________________________________________
 104 AliCentralMultiplicityTask&
 105 AliCentralMultiplicityTask::operator=(const AliCentralMultiplicityTask& o)
 106 {
 107   //
 108   // Assignment operator
 109   //
 110   fInspector         = o.fInspector;
 111   fData              = o.fData;
 112   fList              = o.fList;
 113   fAODCentral        = o.fAODCentral;
 114   fManager           = o.fManager;
 115   fUseSecondary      = o.fUseSecondary;
 116   fUseAcceptance     = o.fUseAcceptance;
 117   fFirstEventSeen    = o.fFirstEventSeen;
 118   fIvz               = 0;
 119   fNClusterTracklet  = o.fNClusterTracklet;
 120   fClusterPerTracklet= o.fClusterPerTracklet;
 121   fNCluster          = o.fNCluster;
 122   fNTracklet         = o.fNTracklet;
 123   fEtaMin            = o.fEtaMin;
 124   fEtaMax            = o.fEtaMax;
 125   return *this;
 126 }
 127 //____________________________________________________________________
 128 void AliCentralMultiplicityTask::UserCreateOutputObjects()
 129 {
 130   //
 131   // Create output objects
 132   //
 133   //
 134
 135   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 136   AliAODHandler*      ah =
 137     dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
 138   if (!ah) AliFatal("No AOD output handler set in analysis manager");
 139
 140
 141   TObject* obj = &fAODCentral;
 142   ah->AddBranch("AliAODCentralMult", &obj);
 143
 144
 145
 146   fList = new TList();
 147   fList->SetOwner();
 148
 149   fInspector.DefineOutput(fList);
 150
 151   PostData(1,fList);
 152 }
 153
 154 //____________________________________________________________________
 155 AliESDEvent*
 156 AliCentralMultiplicityTask::GetESDEvent()
 157 {
 158   //
 159   // Get the ESD event. IF this is the first event, initialise
 160   //
 161   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
 162   if (!esd) {
 163     AliWarning("No ESD event found for input event");
 164     return 0;
 165   }
 166
 167   // IF we've read the first event already, just return the event
 168   if (fFirstEventSeen) return esd;
 169
 170   // Read the details of the rung
 171   fInspector.ReadRunDetails(esd);
 172
 173   // If we weren't initialised before (i.e., in the setup), do so now.
 174   if (!GetManager().IsInit()) {
 175     GetManager().Init(fInspector.GetCollisionSystem(),
 176                       fInspector.GetEnergy(),
 177                       fInspector.GetField());
 178     AliInfo("Manager of corrections in AliCentralMultiplicityTask init");
 179   }
 180
 181   // Check for existence and get secondary map
 182   AliCentralCorrSecondaryMap* secMap = GetManager().GetSecMap();
 183   if (!secMap) AliFatal("No secondary map defined!");
 184   const TAxis& vaxis = secMap->GetVertexAxis();
 185
 186   FindEtaLimits();
 187
 188   fNClusterTracklet = new TH2D("nClusterVsnTracklet",
 189                                "Total number of cluster vs number of tracklets",
 190                                100, 0, 100, 100, 0, 100);
 191   fNClusterTracklet->SetDirectory(0);
 192   fNClusterTracklet->SetXTitle("# of free clusters");
 193   fNClusterTracklet->SetYTitle("# of tracklets");
 194   fNClusterTracklet->SetStats(0);
 195   fList->Add(fNClusterTracklet);
 196
 197   Int_t    nEta = 80;
 198   Double_t lEta = 2;
 199   fClusterPerTracklet = new TH2D("clusterPerTracklet",
 200                                  "N_{free cluster}/N_{tracklet} vs. #eta",
 201                                  nEta,-lEta,lEta, 101, -.05, 10.05);
 202   fClusterPerTracklet->SetDirectory(0);
 203   fClusterPerTracklet->SetXTitle("#eta");
 204   fClusterPerTracklet->SetYTitle("N_{free cluster}/N_{tracklet}");
 205   fClusterPerTracklet->SetStats(0);
 206   fList->Add(fClusterPerTracklet);
 207
 208   // Cache histograms
 209   fNCluster = new TH1D("cacheCluster", "", nEta,-lEta,lEta);
 210   fNCluster->SetDirectory(0);
 211   fNCluster->Sumw2();
 212
 213   fNTracklet = new TH1D("cacheTracklet", "", nEta,-lEta,lEta);
 214   fNTracklet->SetDirectory(0);
 215   fNTracklet->Sumw2();
 216
 217   // Initialize the inspecto
 218   fInspector.Init(vaxis);
 219   fFirstEventSeen = kTRUE;
 220
 221   // Print some information
 222   Print();
 223
 224   return esd;
 225 }
 226 //____________________________________________________________________
 227 void
 228 AliCentralMultiplicityTask::MarkEventForStore() const
 229 {
 230   // Make sure the AOD tree is filled
 231   AliAnalysisManager* am = AliAnalysisManager::GetAnalysisManager();
 232   AliAODHandler*      ah =
 233     dynamic_cast<AliAODHandler*>(am->GetOutputEventHandler());
 234   if (!ah)
 235     AliFatal("No AOD output handler set in analysis manager");
 236
 237   ah->SetFillAOD(kTRUE);
 238 }
 239
 240 //____________________________________________________________________
 241 void AliCentralMultiplicityTask::FindEtaLimits()
 242 {
 243   AliCentralCorrSecondaryMap* secMap = GetManager().GetSecMap();
 244
 245   const TAxis& vaxis = secMap->GetVertexAxis();
 246
 247   fEtaMin.Set(vaxis.GetNbins());
 248   fEtaMax.Set(vaxis.GetNbins());
 249
 250   TList* hits = new TList;
 251   hits->SetOwner();
 252   hits->SetName("hitMaps");
 253   fList->Add(hits);
 254
 255   TList* secs = new TList;
 256   secs->SetOwner();
 257   secs->SetName("secondaryMaps");
 258   fList->Add(secs);
 259   TH2D* hCoverage = new TH2D("coverage", "#eta coverage per v_{z}",
 260                              secMap->GetCorrection(UShort_t(1))->GetXaxis()->GetNbins(),
 261                              secMap->GetCorrection(UShort_t(1))->GetXaxis()->GetXmin(),
 262                              secMap->GetCorrection(UShort_t(1))->GetXaxis()->GetXmax(),
 263                              vaxis.GetNbins(),vaxis.GetXmin(),vaxis.GetXmax());
 264   hCoverage->SetDirectory(0);
 265   hCoverage->SetXTitle("#eta");
 266   hCoverage->SetYTitle("v_{z} [cm]");
 267   hCoverage->SetZTitle("n_{bins}");
 268   fList->Add(hCoverage);
 269
 270   for (Int_t v = 1; v <= vaxis.GetNbins(); v++) {
 271     TH2D* corr = secMap->GetCorrection(UShort_t(v));
 272     TH1D* proj = corr->ProjectionX(Form("secCor%02d", v));
 273     proj->Scale(1. / corr->GetNbinsY());
 274     proj->SetTitle(Form("Projection of secondary correction "
 275                         "for %+5.1f<v_{z}<%+5.1f",
 276                         vaxis.GetBinLowEdge(v), vaxis.GetBinUpEdge(v)));
 277     proj->SetYTitle("#LT 2^{nd} correction#GT");
 278     proj->SetDirectory(0);
 279     proj->SetMarkerStyle(20);
 280     proj->SetMarkerColor(kBlue+1);
 281     secs->Add(proj);
 282
 283     TH2D* obg = static_cast<TH2D*>(corr->Clone(Form("secCor2DFiducial%02d",v)));
 284     obg->SetDirectory(0);
 285     secs->Add(obg);
 286
 287     TH1D* after = static_cast<TH1D*>(proj->Clone(Form("secCorFiducial%02d",v)));
 288     after->SetDirectory(0);
 289     after->SetMarkerColor(kRed+1);
 290     secs->Add(after);
 291
 292     TH2D* data = static_cast<TH2D*>(corr->Clone(Form("hitMap%02d",v)));
 293     //d->SetTitle(Form("hitMap%02d",v));
 294     data->SetTitle(Form("d^{2}N/d#eta d#phi "
 295                         "for %+5.1f<v_{z}<%+5.1f",
 296                         vaxis.GetBinLowEdge(v), vaxis.GetBinUpEdge(v)));
 297     data->GetZaxis()->SetTitle("");
 298     data->SetMarkerColor(kBlack);
 299     data->SetMarkerStyle(1);
 300     hits->Add(data);
 301
 302     TH1D* hAcceptance = fManager.GetAcceptanceCorrection(v);
 303     TH1D* accClone   = static_cast<TH1D*>(hAcceptance->Clone(Form("acceptance%02d",v)));
 304     secs->Add(accClone);
 305
 306     Double_t prev = 0;
 307     for (Int_t e = 1; e <= proj->GetNbinsX(); e++) {
 308       Double_t c = proj->GetBinContent(e);
 309       if (c > .5 /*&& TMath::Abs(c - prev) < .1*c*/) {
 310         fEtaMin[v-1] = e;
 311         break;
 312       }
 313       prev = c;
 314       after->SetBinContent(e, 0);
 315       after->SetBinError(e, 0);
 316       for(Int_t nn =1; nn <=obg->GetNbinsY();nn++)
 317         obg->SetBinContent(e,nn,0);
 318
 319     }
 320     for (Int_t e = proj->GetNbinsX(); e >= 1; e--) {
 321       Double_t c = proj->GetBinContent(e);
 322       if (c > .5 /*&& TMath::Abs(c - prev) < .1*c*/) {
 323         fEtaMax[v-1] = e;
 324         break;
 325       }
 326       prev = c;
 327       after->SetBinContent(e, 0);
 328       after->SetBinError(e, 0);
 329       for(Int_t nn =1; nn <=obg->GetNbinsY();nn++)
 330         obg->SetBinContent(e,nn,0);
 331
 332     }
 333
 334     for (Int_t nn = fEtaMin[v-1]; nn<=fEtaMax[v-1]; nn++) {
 335       hCoverage->SetBinContent(nn,v,1);
 336     }
 337
 338   }
 339 }
 340
 341 //____________________________________________________________________
 342 void AliCentralMultiplicityTask::UserExec(Option_t* /*option*/)
 343 {
 344   //
 345   // Process each event
 346   //
 347   // Parameters:
 348   //    option Not used
 349   //
 350   fAODCentral.Clear("");
 351   fIvz = 0;
 352
 353   AliESDEvent* esd = GetESDEvent();
 354
 355   Bool_t   lowFlux   = kFALSE;
 356   UInt_t   triggers  = 0;
 357   UShort_t ivz       = 0;
 358   Double_t vz        = 0;
 359   Double_t cent      = -1;
 360   UShort_t nClusters = 0;
 361   UInt_t   found     = fInspector.Process(esd, triggers, lowFlux,
 362                                           ivz, vz, cent, nClusters);
 363
 364   // No event or no trigger
 365   if (found & AliFMDEventInspector::kNoEvent)    return;
 366   if (found & AliFMDEventInspector::kNoTriggers) return;
 367
 368   // Make sure AOD is filled
 369   MarkEventForStore();
 370
 371   if (found == AliFMDEventInspector::kNoSPD)      return;
 372   if (found == AliFMDEventInspector::kNoVertex)   return;
 373   if (triggers & AliAODForwardMult::kPileUp)      return;
 374   if (found == AliFMDEventInspector::kBadVertex)  return; // Out of range
 375
 376   //Doing analysis
 377   fIvz = ivz;
 378   const AliMultiplicity* spdmult = esd->GetMultiplicity();
 379
 380   TH2D& aodHist = fAODCentral.GetHistogram();
 381
 382   ProcessESD(aodHist, spdmult);
 383   CorrectData(aodHist, ivz);
 384   //Producing hit maps
 385   TList* hitList = static_cast<TList*>(fList->FindObject("hitMaps"));
 386   TH2D* data  = 0;
 387   if(hitList)
 388     data = static_cast<TH2D*>(hitList->At(ivz-1));
 389   if(data)
 390     data->Add(&aodHist);
 391
 392   PostData(1,fList);
 393 }
 394 //____________________________________________________________________
 395 void
 396 AliCentralMultiplicityTask::ProcessESD(TH2D& aodHist,
 397                                        const AliMultiplicity* spdmult) const
 398 {
 399   fNTracklet->Reset();
 400   fNCluster->Reset();
 401
 402   //Filling clusters in layer 1 used for tracklets...
 403   for(Int_t j = 0; j< spdmult->GetNumberOfTracklets();j++) {
 404     Double_t eta = spdmult->GetEta(j);
 405     fNTracklet->Fill(eta);
 406     aodHist.Fill(eta,spdmult->GetPhi(j));
 407   }
 408
 409   //...and then the unused ones in layer 1
 410   for(Int_t j = 0; j< spdmult->GetNumberOfSingleClusters();j++) {
 411     Double_t eta = -TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.));
 412     fNCluster->Fill(eta);
 413     aodHist.Fill(eta, spdmult->GetPhiSingle(j));
 414   }
 415   fNClusterTracklet->Fill(fNCluster->GetEntries(),
 416                           fNTracklet->GetEntries());
 417
 418   fNCluster->Divide(fNTracklet);
 419   for (Int_t j = 1; j <= fNCluster->GetNbinsX(); j++)
 420     fClusterPerTracklet->Fill(fNCluster->GetXaxis()->GetBinCenter(j),
 421                               fNCluster->GetBinContent(j));
 422
 423 }
 424
 425 //____________________________________________________________________
 426 void
 427 AliCentralMultiplicityTask::CorrectData(TH2D& aodHist, UShort_t vtxbin) const
 428 {
 429   // Corrections
 430   TH1D* hAcceptance = fManager.GetAcceptanceCorrection(vtxbin);
 431   TH2D* hSecMap     = fManager.GetSecMapCorrection(vtxbin);
 432
 433   if (!hSecMap)     AliFatal("No secondary map!");
 434   if (!hAcceptance) AliFatal("No acceptance!");
 435
 436   if (fUseSecondary && hSecMap) aodHist.Divide(hSecMap);
 437
 438   for(Int_t nx = 1; nx <= aodHist.GetNbinsX(); nx++) {
 439     Float_t accCor = hAcceptance->GetBinContent(nx);
 440     Float_t accErr = hAcceptance->GetBinError(nx);
 441
 442     Bool_t fiducial = true;
 443     if (nx < fEtaMin[vtxbin-1] || nx > fEtaMax[vtxbin-1])
 444       fiducial = false;
 445     //  Bool_t etabinSeen = kFALSE;
 446     for(Int_t ny = 1; ny <= aodHist.GetNbinsY(); ny++) {
 447 #if 1
 448       if (!fiducial) {
 449         aodHist.SetBinContent(nx, ny, 0);
 450         aodHist.SetBinError(nx, ny, 0);
 451         continue;
 452       }
 453 #endif
 454       // Get currrent value
 455       Float_t aodValue = aodHist.GetBinContent(nx,ny);
 456       Float_t aodErr   = aodHist.GetBinError(nx,ny);
 457
 458 #if 0 // This is done once in the FindEtaBins function
 459       // Set underflow bin
 460       Float_t secCor   = 0;
 461       if(hSecMap)       secCor     = hSecMap->GetBinContent(nx,ny);
 462       if (secCor > 0.5) etabinSeen = kTRUE;
 463 #endif
 464       if (aodValue < 0.000001) {
 465         aodHist.SetBinContent(nx,ny, 0);
 466         aodHist.SetBinError(nx,ny, 0);
 467         continue;
 468       }
 469       if (!fUseAcceptance) continue;
 470
 471       // Acceptance correction
 472       if (accCor   < 0.000001) accCor = 1;
 473       Float_t aodNew   = aodValue / accCor ;
 474       Float_t error    = aodNew*TMath::Sqrt(TMath::Power(aodErr/aodValue,2) +
 475                                             TMath::Power(accErr/accCor,2) );
 476       aodHist.SetBinContent(nx,ny, aodNew);
 477       //test
 478       aodHist.SetBinError(nx,ny,error);
 479       aodHist.SetBinError(nx,ny,aodErr);
 480     }
 481     //Filling underflow bin if we eta bin is in range
 482     if (fiducial) aodHist.SetBinContent(nx,0, 1.);
 483     // if (etabinSeen) aodHist.SetBinContent(nx,0, 1.);
 484   }
 485 }
 486
 487 //____________________________________________________________________
 488 void AliCentralMultiplicityTask::Terminate(Option_t* /*option*/)
 489 {
 490   //
 491   // End of job
 492   //
 493   // Parameters:
 494   //    option Not used
 495   //
 496 }
 497 //____________________________________________________________________
 498 void
 499 AliCentralMultiplicityTask::Print(Option_t* option) const
 500 {
 501   //
 502   // Print information
 503   //
 504   // Parameters:
 505   //    option Not used
 506   //
 507   std::cout << ClassName() << ": " << GetName() << "\n"
 508             << std::boolalpha
 509             << "  Use secondary correction:  " << fUseSecondary << '\n'
 510             << "  Use acceptance correction: " << fUseAcceptance << '\n'
 511             << "  Off-line trigger mask:  0x"
 512             << std::hex     << std::setfill('0')
 513             << std::setw (8) << fOfflineTriggerMask
 514             << std::dec     << std::setfill (' ')
 515             << std::noboolalpha << std::endl;
 516   AliCentralCorrSecondaryMap* secMap = GetManager().GetSecMap();
 517   if (secMap) {
 518     const TAxis& vaxis = secMap->GetVertexAxis();
 519     std::cout << "  Eta ranges:\n"
 520             << "     Vertex        | Eta bins\n"
 521               << "   bin     range   | \n"
 522               << "   ----------------+-----------" << std::endl;
 523     for (Int_t v = 1; v <= vaxis.GetNbins(); v++) {
 524       std::cout << "   " << std::setw(2) << v << "  "
 525                 << std::setw(5) << vaxis.GetBinLowEdge(v) << "-"
 526                 << std::setw(5) << vaxis.GetBinUpEdge(v) << " | "
 527                 << std::setw(3) << fEtaMin[v-1] << "-"
 528                 << std::setw(3) << fEtaMax[v-1] << std::endl;
 529     }
 530   }
 531
 532   gROOT->IncreaseDirLevel();
 533   fManager.Print(option);
 534   fInspector.Print(option);
 535   gROOT->DecreaseDirLevel();
 536
 537 }
 538 //====================================================================
 539 AliCentralMultiplicityTask::Manager::Manager() :
 540   fAcceptancePath("$ALICE_ROOT/PWG2/FORWARD/corrections/CentralAcceptance"),
 541   fSecMapPath("$ALICE_ROOT/PWG2/FORWARD/corrections/CentralSecMap"),
 542   fAcceptance(),
 543   fSecmap(),
 544   fAcceptanceName("centralacceptance"),
 545   fSecMapName("centralsecmap"),
 546   fIsInit(kFALSE)
 547 {
 548   //
 549   // Constructor
 550   //
 551 }
 552 //____________________________________________________________________
 553 AliCentralMultiplicityTask::Manager::Manager(const Manager& o)
 554   :fAcceptancePath(o.fAcceptancePath),
 555    fSecMapPath(o.fSecMapPath),
 556    fAcceptance(o.fAcceptance),
 557    fSecmap(o.fSecmap),
 558    fAcceptanceName(o.fAcceptanceName),
 559    fSecMapName(o.fSecMapName),
 560    fIsInit(o.fIsInit)
 561 {
 562   //
 563   // Copy Constructor
 564   //
 565 }
 566 //____________________________________________________________________
 567 AliCentralMultiplicityTask::Manager&
 568 AliCentralMultiplicityTask::Manager::operator=(const Manager& o)
 569 {
 570   //
 571   // Assignment operator
 572   //
 573   fAcceptancePath = o.fAcceptancePath;
 574   fSecMapPath     = o.fSecMapPath;
 575   fAcceptance     = o.fAcceptance;
 576   fSecmap         = o.fSecmap;
 577   fAcceptanceName = o.fAcceptanceName;
 578   fSecMapName     = o.fSecMapName;
 579   fIsInit         = o.fIsInit;
 580   return *this;
 581 }
 582
 583 //____________________________________________________________________
 584 const char*
 585 AliCentralMultiplicityTask::Manager::GetFullFileName(UShort_t what,
 586                                                      UShort_t sys,
 587                                                      UShort_t sNN,
 588                                                      Short_t  field) const
 589 {
 590   //
 591   // Get full path name to object file
 592   //
 593   // Parameters:
 594   //    what   What to get
 595   //    sys    Collision system
 596   //    sNN    Center of mass energy
 597   //    field  Magnetic field
 598   //
 599   // Return:
 600   //
 601   //
 602   return Form("%s/%s",
 603               what == 0 ? GetSecMapPath() : GetAcceptancePath(),
 604               GetFileName(what, sys, sNN, field));
 605 }
 606
 607 //____________________________________________________________________
 608 const char*
 609 AliCentralMultiplicityTask::Manager::GetFileName(UShort_t  what ,
 610                                                  UShort_t  sys,
 611                                                  UShort_t  sNN,
 612                                                  Short_t   field) const
 613 {
 614   //
 615   // Get the full path name
 616   //
 617   // Parameters:
 618   //    what   What to get
 619   //    sys    Collision system
 620   //    sNN    Center of mass energy
 621   //    field  Magnetic field
 622   //
 623   // Return:
 624   //
 625   //
 626   // Must be static - otherwise the data may disappear on return from
 627   // this member function
 628   static TString fname = "";
 629
 630   switch(what) {
 631   case 0:  fname = fSecMapName;     break;
 632   case 1:  fname = fAcceptanceName; break;
 633   default:
 634     ::Error("GetFileName",
 635             "Invalid indentifier %d for central object, must be 0 or 1!", what);
 636     break;
 637   }
 638   fname.Append(Form("_%s_%04dGeV_%c%1dkG.root",
 639                     AliForwardUtil::CollisionSystemString(sys),
 640                     sNN, (field < 0 ? 'm' : 'p'), TMath::Abs(field)));
 641
 642   return fname.Data();
 643 }
 644
 645 //____________________________________________________________________
 646 TH2D*
 647 AliCentralMultiplicityTask::Manager::GetSecMapCorrection(UShort_t vtxbin) const
 648 {
 649   //
 650   // Get the secondary map
 651   //
 652   // Parameters:
 653   //    vtxbin
 654   //
 655   // Return:
 656   //
 657   //
 658   if (!fSecmap) {
 659     ::Warning("GetSecMapCorrection","No secondary map defined");
 660     return 0;
 661   }
 662   return fSecmap->GetCorrection(vtxbin);
 663 }
 664 //____________________________________________________________________
 665 TH1D*
 666 AliCentralMultiplicityTask::Manager::GetAcceptanceCorrection(UShort_t vtxbin)
 667   const
 668 {
 669   //
 670   // Get the acceptance correction
 671   //
 672   // Parameters:
 673   //    vtxbin
 674   //
 675   // Return:
 676   //
 677   //
 678   if (!fAcceptance) {
 679     ::Warning("GetAcceptanceCorrection","No acceptance map defined");
 680     return 0;
 681   }
 682   return fAcceptance->GetCorrection(vtxbin);
 683 }
 684
 685 //____________________________________________________________________
 686 void
 687 AliCentralMultiplicityTask::Manager::Init(UShort_t  sys,
 688                                           UShort_t  sNN,
 689                                           Short_t   field)
 690 {
 691   //
 692   // Initialize
 693   //
 694   // Parameters:
 695   //    sys    Collision system (1: pp, 2: PbPb)
 696   //    sNN    Center of mass energy per nucleon pair [GeV]
 697   //    field  Magnetic field [kG]
 698   //
 699   if(fIsInit) ::Warning("Init","Already initialised - overriding...");
 700
 701   TFile fsec(GetFullFileName(0,sys,sNN,field));
 702   fSecmap =
 703     dynamic_cast<AliCentralCorrSecondaryMap*>(fsec.Get(fSecMapName.Data()));
 704   if(!fSecmap) {
 705     ::Error("Init", "no central Secondary Map found!") ;
 706     return;
 707   }
 708   TFile facc(GetFullFileName(1,sys,sNN,field));
 709   fAcceptance =
 710     dynamic_cast<AliCentralCorrAcceptance*>(facc.Get(fAcceptanceName.Data()));
 711   if(!fAcceptance) {
 712     ::Error("Init", "no central Acceptance found!") ;
 713     return;
 714   }
 715
 716   if(fSecmap && fAcceptance) {
 717     fIsInit = kTRUE;
 718     ::Info("Init",
 719            "Central Manager initialised for %s, energy %dGeV, field %dkG",
 720            sys == 1 ? "pp" : sys == 2 ? "PbPb" : "unknown", sNN,field);
 721   }
 722 }
 723 //____________________________________________________________________
 724 Bool_t
 725 AliCentralMultiplicityTask::Manager::WriteFile(UShort_t what,
 726                                                UShort_t sys,
 727                                                UShort_t sNN,
 728                                                Short_t  fld,
 729                                                TObject* obj,
 730                                                Bool_t   full) const
 731 {
 732   //
 733   // Write correction output to (a temporary) file
 734   //
 735   // Parameters:
 736   //   What     What to write
 737   //   sys      Collision system (1: pp, 2: PbPb)
 738   //   sNN      Center of mass energy per nucleon (GeV)
 739   //   fld      Field (kG)
 740   //   obj      Object to write
 741   //   full     if true, write to full path, otherwise locally
 742   //
 743   // Return:
 744   //   true on success.
 745   TString ofName;
 746   if (!full)
 747     ofName = GetFileName(what, sys, sNN, fld);
 748   else
 749     ofName = GetFullFileName(what, sys, sNN, fld);
 750   if (ofName.IsNull()) {
 751     AliErrorGeneral("Manager",Form("Unknown object type %d", what));
 752     return false;
 753   }
 754   TFile* output = TFile::Open(ofName, "RECREATE");
 755   if (!output) {
 756     AliErrorGeneral("Manager",Form("Failed to open file %s", ofName.Data()));
 757     return false;
 758   }
 759
 760   TString oName(GetObjectName(what));
 761   Int_t ret = obj->Write(oName);
 762   if (ret <= 0) {
 763     AliErrorGeneral("Manager",Form("Failed to write %p to %s/%s (%d)",
 764                                    obj, ofName.Data(), oName.Data(), ret));
 765     return false;
 766   }
 767
 768   ret = output->Write();
 769   if (ret < 0) {
 770     AliErrorGeneral("Manager",
 771                     Form("Failed to write %s to disk (%d)", ofName.Data(),ret));
 772     return false;
 773   }
 774   output->ls();
 775   output->Close();
 776
 777   TString cName(obj->IsA()->GetName());
 778   AliInfoGeneral("Manager",
 779                  Form("Wrote %s object %s to %s\n",
 780                       cName.Data(),oName.Data(), ofName.Data()));
 781   if (!full) {
 782     TString dName(GetFileDir(what));
 783     AliInfoGeneral("Manager",
 784                    Form("%s should be copied to %s\n"
 785                         "Do for example\n\t"
 786                         "aliroot $ALICE_ROOT/PWG2/FORWARD/analysis2/scripts/"
 787                         "MoveCorrections.C\\(%d\\)\nor\n\t"
 788                         "cp %s %s/",
 789                         ofName.Data(),dName.Data(),
 790                         what, ofName.Data(),
 791                         gSystem->ExpandPathName(dName.Data())));
 792
 793
 794   }
 795   return true;
 796 }
 797
 798 //____________________________________________________________________
 799 void
 800 AliCentralMultiplicityTask::Manager::Print(Option_t* option) const
 801 {
 802   //
 803   // Print information to standard output
 804   //
 805   std::cout << " AliCentralMultiplicityTask::Manager\n"
 806             << std::boolalpha
 807             << "  Initialized:     " << fIsInit << '\n'
 808             << "  Acceptance path: " << fAcceptancePath << '\n'
 809             << "  Acceptance name: " << fAcceptanceName << '\n'
 810             << "  Acceptance:      " << fAcceptance << '\n'
 811             << "  Secondary path:  " << fSecMapPath << '\n'
 812             << "  Secondary name:  " << fSecMapName << '\n'
 813             << "  Secondary map:   " << fSecmap
 814             << std::noboolalpha << std::endl;
 815   if (fAcceptance) fAcceptance->Print(option);
 816   if (fSecmap)     fSecmap->Print(option);
 817 }
 818
 819 //
 820 // EOF
 821 //