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