PWGLF/FORWARD/analysis2/AliForwardMCCorrectionsTask.cxx

   1 //
   2 // Calculate the corrections in the forward regions
   3 //
   4 // Inputs:
   5 //   - AliESDEvent
   6 //
   7 // Outputs:
   8 //   - AliAODForwardMult
   9 //
  10 // Histograms
  11 //
  12 // Corrections used
  13 //
  14 #include "AliForwardMCCorrectionsTask.h"
  15 #include "AliTriggerAnalysis.h"
  16 #include "AliPhysicsSelection.h"
  17 #include "AliLog.h"
  18 #include "AliHeader.h"
  19 #include "AliGenEventHeader.h"
  20 #include "AliESDEvent.h"
  21 #include "AliAODHandler.h"
  22 #include "AliMultiplicity.h"
  23 #include "AliInputEventHandler.h"
  24 #include "AliStack.h"
  25 #include "AliMCEvent.h"
  26 #include "AliAODForwardMult.h"
  27 #include "AliFMDStripIndex.h"
  28 #include "AliFMDCorrSecondaryMap.h"
  29 #include <TH1.h>
  30 #include <TH2D.h>
  31 #include <TDirectory.h>
  32 #include <TTree.h>
  33 #include <TList.h>
  34 #include <TROOT.h>
  35 #include <iostream>
  36
  37 //====================================================================
  38 namespace {
  39   const char* GetEventName(Bool_t tr, Bool_t vtx)
  40   {
  41     return Form("nEvents%s%s", (tr ? "Tr" : ""), (vtx ? "Vtx" : ""));
  42   }
  43 #if 0
  44   const char* GetHitsName(UShort_t d, Char_t r)
  45   {
  46     return Form("hitsFMD%d%c", d, r);
  47   }
  48   const char* GetStripsName(UShort_t d, Char_t r)
  49   {
  50     return Form("stripsFMD%d%c", d, r);
  51   }
  52   const char* GetPrimaryName(Char_t r, Bool_t trVtx)
  53   {
  54     return Form("primaries%s%s",
  55                 ((r == 'I' || r == 'i') ? "Inner" : "Outer"),
  56                 (trVtx ? "TrVtx" : "All"));
  57   }
  58 #endif
  59 }
  60
  61 //====================================================================
  62 AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask()
  63   : AliAnalysisTaskSE(),
  64     fInspector(),
  65     fTrackDensity(),
  66     fESDFMD(),
  67     fVtxBins(0),
  68     fFirstEvent(true),
  69     fHEvents(0),
  70     fHEventsTr(0),
  71     fHEventsTrVtx(0),
  72     fVtxAxis(),
  73     fEtaAxis(),
  74     fList()
  75 {
  76   //
  77   // Constructor
  78   //
  79   // Parameters:
  80   //    name Name of task
  81   //
  82 }
  83
  84 //____________________________________________________________________
  85 AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const char* name)
  86   : AliAnalysisTaskSE(name),
  87     fInspector("eventInspector"),
  88     fTrackDensity("trackDensity"),
  89     fESDFMD(),
  90     fVtxBins(0),
  91     fFirstEvent(true),
  92     fHEvents(0),
  93     fHEventsTr(0),
  94     fHEventsTrVtx(0),
  95     fVtxAxis(10,-10,10),
  96     fEtaAxis(200,-4,6),
  97     fList()
  98 {
  99   //
 100   // Constructor
 101   //
 102   // Parameters:
 103   //    name Name of task
 104   //
 105   DefineOutput(1, TList::Class());
 106   DefineOutput(2, TList::Class());
 107 }
 108
 109 //____________________________________________________________________
 110 AliForwardMCCorrectionsTask::AliForwardMCCorrectionsTask(const AliForwardMCCorrectionsTask& o)
 111   : AliAnalysisTaskSE(o),
 112     fInspector(o.fInspector),
 113     fTrackDensity(),
 114     fESDFMD(o.fESDFMD),
 115     fVtxBins(0),
 116     fFirstEvent(o.fFirstEvent),
 117     fHEvents(o.fHEvents),
 118     fHEventsTr(o.fHEventsTr),
 119     fHEventsTrVtx(o.fHEventsTrVtx),
 120     fVtxAxis(10,-10,10),
 121     fEtaAxis(200,-4,6),
 122     fList(o.fList)
 123 {
 124   //
 125   // Copy constructor
 126   //
 127   // Parameters:
 128   //    o Object to copy from
 129   //
 130 }
 131
 132 //____________________________________________________________________
 133 AliForwardMCCorrectionsTask&
 134 AliForwardMCCorrectionsTask::operator=(const AliForwardMCCorrectionsTask& o)
 135 {
 136   //
 137   // Assignment operator
 138   //
 139   // Parameters:
 140   //    o Object to assign from
 141   //
 142   // Return:
 143   //    Reference to this object
 144   //
 145   if (&o == this) return *this;
 146   fInspector         = o.fInspector;
 147   fTrackDensity      = o.fTrackDensity;
 148   fESDFMD            = o.fESDFMD;
 149   fVtxBins           = o.fVtxBins;
 150   fFirstEvent        = o.fFirstEvent;
 151   fHEvents           = o.fHEvents;
 152   fHEventsTr         = o.fHEventsTr;
 153   fHEventsTrVtx      = o.fHEventsTrVtx;
 154   SetVertexAxis(o.fVtxAxis);
 155   SetEtaAxis(o.fEtaAxis);
 156
 157   return *this;
 158 }
 159
 160 //____________________________________________________________________
 161 void
 162 AliForwardMCCorrectionsTask::Init()
 163 {
 164   //
 165   // Initialize the task
 166   //
 167   //
 168 }
 169
 170 //____________________________________________________________________
 171 void
 172 AliForwardMCCorrectionsTask::SetVertexAxis(Int_t nBin, Double_t min,
 173                                            Double_t max)
 174 {
 175   //
 176   // Set the vertex axis to use
 177   //
 178   // Parameters:
 179   //    nBins Number of bins
 180   //    vzMin Least @f$z@f$ coordinate of interation point
 181   //    vzMax Largest @f$z@f$ coordinate of interation point
 182   //
 183   if (max < min) {
 184     Double_t tmp = min;
 185     min          = max;
 186     max          = tmp;
 187   }
 188 /*
 189   if (min < -10)
 190     AliWarning(Form("Minimum vertex %f < -10, make sure you want this",min));
 191   if (max > +10)
 192     AliWarning(Form("Minimum vertex %f > +10, make sure you want this",max));
 193 */
 194   fVtxAxis.Set(nBin, min, max);
 195 }
 196 //____________________________________________________________________
 197 void
 198 AliForwardMCCorrectionsTask::SetVertexAxis(const TAxis& axis)
 199 {
 200   //
 201   // Set the vertex axis to use
 202   //
 203   // Parameters:
 204   //    axis Axis
 205   //
 206   SetVertexAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
 207 }
 208
 209 //____________________________________________________________________
 210 void
 211 AliForwardMCCorrectionsTask::SetEtaAxis(Int_t nBin, Double_t min, Double_t max)
 212 {
 213   //
 214   // Set the eta axis to use
 215   //
 216   // Parameters:
 217   //    nBins Number of bins
 218   //    vzMin Least @f$\eta@f$
 219   //    vzMax Largest @f$\eta@f$
 220   //
 221   if (max < min) {
 222     Double_t tmp = min;
 223     min          = max;
 224     max          = tmp;
 225   }
 226   if (min < -4)
 227     AliWarning(Form("Minimum eta %f < -4, make sure you want this",min));
 228   if (max > +6)
 229     AliWarning(Form("Minimum eta %f > +6, make sure you want this",max));
 230   fEtaAxis.Set(nBin, min, max);
 231 }
 232 //____________________________________________________________________
 233 void
 234 AliForwardMCCorrectionsTask::SetEtaAxis(const TAxis& axis)
 235 {
 236   //
 237   // Set the eta axis to use
 238   //
 239   // Parameters:
 240   //    axis Axis
 241   //
 242   SetEtaAxis(axis.GetNbins(),axis.GetXmin(),axis.GetXmax());
 243 }
 244
 245 //____________________________________________________________________
 246 void
 247 AliForwardMCCorrectionsTask::DefineBins(TList* l)
 248 {
 249   if (!fVtxBins) fVtxBins = new TObjArray(fVtxAxis.GetNbins(), 1);
 250   if (fVtxBins->GetEntries() > 0) return;
 251
 252   fVtxBins->SetOwner();
 253   for (Int_t  i    = 1; i <= fVtxAxis.GetNbins(); i++) {
 254     Double_t  low  = fVtxAxis.GetBinLowEdge(i);
 255     Double_t  high = fVtxAxis.GetBinUpEdge(i);
 256     VtxBin*   bin  = new VtxBin(low, high, fEtaAxis);
 257     fVtxBins->AddAt(bin, i);
 258     bin->DefineOutput(l);
 259   }
 260 }
 261
 262 //____________________________________________________________________
 263 void
 264 AliForwardMCCorrectionsTask::UserCreateOutputObjects()
 265 {
 266   //
 267   // Create output objects
 268   //
 269   //
 270   fList = new TList;
 271   fList->SetOwner();
 272   fList->SetName(Form("%sSums", GetName()));
 273
 274   DefineBins(fList);
 275
 276   fHEvents = new TH1I(GetEventName(false,false),
 277                       "Number of all events",
 278                       fVtxAxis.GetNbins(),
 279                       fVtxAxis.GetXmin(),
 280                       fVtxAxis.GetXmax());
 281   fHEvents->SetXTitle("v_{z} [cm]");
 282   fHEvents->SetYTitle("# of events");
 283   fHEvents->SetFillColor(kBlue+1);
 284   fHEvents->SetFillStyle(3001);
 285   fHEvents->SetDirectory(0);
 286   fList->Add(fHEvents);
 287
 288   fHEventsTr = new TH1I(GetEventName(true, false),
 289                         "Number of triggered events",
 290                         fVtxAxis.GetNbins(),
 291                         fVtxAxis.GetXmin(),
 292                         fVtxAxis.GetXmax());
 293   fHEventsTr->SetXTitle("v_{z} [cm]");
 294   fHEventsTr->SetYTitle("# of events");
 295   fHEventsTr->SetFillColor(kRed+1);
 296   fHEventsTr->SetFillStyle(3001);
 297   fHEventsTr->SetDirectory(0);
 298   fList->Add(fHEventsTr);
 299
 300   fHEventsTrVtx = new TH1I(GetEventName(true,true),
 301                            "Number of events w/trigger and vertex",
 302                            fVtxAxis.GetNbins(),
 303                            fVtxAxis.GetXmin(),
 304                            fVtxAxis.GetXmax());
 305   fHEventsTrVtx->SetXTitle("v_{z} [cm]");
 306   fHEventsTrVtx->SetYTitle("# of events");
 307   fHEventsTrVtx->SetFillColor(kBlue+1);
 308   fHEventsTrVtx->SetFillStyle(3001);
 309   fHEventsTrVtx->SetDirectory(0);
 310   fList->Add(fHEventsTrVtx);
 311
 312   // Copy axis objects to output
 313   TH1* vtxAxis = new TH1D("vtxAxis", "Vertex axis",
 314                           fVtxAxis.GetNbins(),
 315                           fVtxAxis.GetXmin(),
 316                           fVtxAxis.GetXmax());
 317   TH1* etaAxis = new TH1D("etaAxis", "Eta axis",
 318                           fEtaAxis.GetNbins(),
 319                           fEtaAxis.GetXmin(),
 320                           fEtaAxis.GetXmax());
 321   fList->Add(vtxAxis);
 322   fList->Add(etaAxis);
 323
 324   AliInfo(Form("Initialising sub-routines: %p, %p",
 325                &fInspector, &fTrackDensity));
 326   fInspector.DefineOutput(fList);
 327   fInspector.Init(fVtxAxis);
 328   fTrackDensity.DefineOutput(fList);
 329
 330   PostData(1, fList);
 331 }
 332 //____________________________________________________________________
 333 void
 334 AliForwardMCCorrectionsTask::UserExec(Option_t*)
 335 {
 336   //
 337   // Process each event
 338   //
 339   // Parameters:
 340   //    option Not used
 341   //
 342
 343   // Get the input data - MC event
 344   AliMCEvent*  mcEvent = MCEvent();
 345   if (!mcEvent) {
 346     AliWarning("No MC event found");
 347     return;
 348   }
 349
 350   // Get the input data - ESD event
 351   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
 352   if (!esd) {
 353     AliWarning("No ESD event found for input event");
 354     return;
 355   }
 356
 357   //--- Read run information -----------------------------------------
 358   if (fFirstEvent && esd->GetESDRun()) {
 359     fInspector.ReadRunDetails(esd);
 360
 361     AliInfo(Form("Initializing with parameters from the ESD:\n"
 362                  "         AliESDEvent::GetBeamEnergy()   ->%f\n"
 363                  "         AliESDEvent::GetBeamType()     ->%s\n"
 364                  "         AliESDEvent::GetCurrentL3()    ->%f\n"
 365                  "         AliESDEvent::GetMagneticField()->%f\n"
 366                  "         AliESDEvent::GetRunNumber()    ->%d\n",
 367                  esd->GetBeamEnergy(),
 368                  esd->GetBeamType(),
 369                  esd->GetCurrentL3(),
 370                  esd->GetMagneticField(),
 371                  esd->GetRunNumber()));
 372
 373     Print();
 374     fFirstEvent = false;
 375   }
 376
 377   // Some variables
 378   UInt_t   triggers; // Trigger bits
 379   Bool_t   lowFlux;  // Low flux flag
 380   UShort_t iVz;      // Vertex bin from ESD
 381   Double_t vZ;       // Z coordinate from ESD
 382   Double_t cent;     // Centrality
 383   UShort_t iVzMc;    // Vertex bin from MC
 384   Double_t vZMc;     // Z coordinate of IP vertex from MC
 385   Double_t b;        // Impact parameter
 386   Int_t    nPart;    // Number of participants
 387   Int_t    nBin;     // Number of binary collisions
 388   Double_t phiR;     // Reaction plane from MC
 389   UShort_t nClusters;// Number of SPD clusters
 390   // Process the data
 391   UInt_t retESD = fInspector.Process(esd, triggers, lowFlux, iVz, vZ,
 392                                      cent, nClusters);
 393   fInspector.ProcessMC(mcEvent, triggers, iVzMc, vZMc,
 394                        b, nPart, nBin, phiR);
 395
 396   Bool_t isInel   = triggers & AliAODForwardMult::kInel;
 397   Bool_t hasVtx   = retESD == AliFMDMCEventInspector::kOk;
 398
 399   // Fill the event count histograms
 400   if (isInel)           fHEventsTr->Fill(vZMc);
 401   if (isInel && hasVtx) fHEventsTrVtx->Fill(vZMc);
 402   fHEvents->Fill(vZMc);
 403
 404   // Now find our vertex bin object
 405   VtxBin* bin = 0;
 406   if (iVzMc > 0 && iVzMc <= fVtxAxis.GetNbins())
 407     bin = static_cast<VtxBin*>(fVtxBins->At(iVzMc));
 408   if (!bin) {
 409     // AliError(Form("No vertex bin object @ %d (%f)", iVzMc, vZMc));
 410     return;
 411   }
 412
 413   // Clear our ESD object
 414   fESDFMD.Clear();
 415
 416   // Get FMD data
 417   AliESDFMD* esdFMD = esd->GetFMDData();
 418
 419   // Now process our input data and store in own ESD object
 420   fTrackDensity.Calculate(*esdFMD, *mcEvent, vZMc, fESDFMD, bin->fPrimary);
 421   bin->fCounts->Fill(0.5);
 422
 423   // And then bin the data in our vtxbin
 424   for (UShort_t d=1; d<=3; d++) {
 425     UShort_t nr = (d == 1 ? 1 : 2);
 426     for (UShort_t q=0; q<nr; q++) {
 427       Char_t      r = (q == 0 ? 'I' : 'O');
 428       UShort_t    ns= (q == 0 ?  20 :  40);
 429       UShort_t    nt= (q == 0 ? 512 : 256);
 430       TH2D*       h = bin->fHists.Get(d,r);
 431
 432       for (UShort_t s=0; s<ns; s++) {
 433         for (UShort_t t=0; t<nt; t++) {
 434           Float_t mult = fESDFMD.Multiplicity(d,r,s,t);
 435
 436           if (mult == 0 || mult > 20) continue;
 437
 438           Float_t phi = fESDFMD.Phi(d,r,s,t) / 180 * TMath::Pi();
 439           Float_t eta = fESDFMD.Eta(d,r,s,t);
 440           h->Fill(eta,phi,mult);
 441         } // for t
 442       } // for s
 443     } // for q
 444   } // for d
 445
 446   PostData(1, fList);
 447 }
 448 //____________________________________________________________________
 449 void
 450 AliForwardMCCorrectionsTask::Terminate(Option_t*)
 451 {
 452   //
 453   // End of job
 454   //
 455   // Parameters:
 456   //    option Not used
 457   //
 458
 459   fList = dynamic_cast<TList*>(GetOutputData(1));
 460   if (!fList) {
 461     AliError("No output list defined");
 462     return;
 463   }
 464
 465   DefineBins(fList);
 466
 467   // Output list
 468   TList* output = new TList;
 469   output->SetOwner();
 470   output->SetName(Form("%sResults", GetName()));
 471
 472   // --- Fill correction object --------------------------------------
 473   AliFMDCorrSecondaryMap* corr = new AliFMDCorrSecondaryMap;
 474   corr->SetVertexAxis(fVtxAxis);
 475   corr->SetEtaAxis(fEtaAxis);
 476
 477   TIter     next(fVtxBins);
 478   VtxBin*   bin = 0;
 479   UShort_t  iVz = 1;
 480   while ((bin = static_cast<VtxBin*>(next())))
 481     bin->Finish(fList, output, iVz++, corr);
 482
 483   output->Add(corr);
 484
 485   PostData(2, output);
 486 }
 487
 488 //____________________________________________________________________
 489 void
 490 AliForwardMCCorrectionsTask::Print(Option_t* option) const
 491 {
 492   std::cout << ClassName() << "\n"
 493             << "  Vertex bins:      " << fVtxAxis.GetNbins() << '\n'
 494             << "  Vertex range:     [" << fVtxAxis.GetXmin()
 495             << "," << fVtxAxis.GetXmax() << "]\n"
 496             << "  Eta bins:         " << fEtaAxis.GetNbins() << '\n'
 497             << "  Eta range:        [" << fEtaAxis.GetXmin()
 498             << "," << fEtaAxis.GetXmax() << "]"
 499             << std::endl;
 500   gROOT->IncreaseDirLevel();
 501   fInspector.Print(option);
 502   fTrackDensity.Print(option);
 503   gROOT->DecreaseDirLevel();
 504 }
 505
 506 //====================================================================
 507 const char*
 508 AliForwardMCCorrectionsTask::VtxBin::BinName(Double_t low,
 509                                              Double_t high)
 510 {
 511   static TString buf;
 512   buf = Form("vtx%+05.1f_%+05.1f", low, high);
 513   buf.ReplaceAll("+", "p");
 514   buf.ReplaceAll("-", "m");
 515   buf.ReplaceAll(".", "d");
 516   return buf.Data();
 517 }
 518
 519
 520 //____________________________________________________________________
 521 AliForwardMCCorrectionsTask::VtxBin::VtxBin()
 522   : fHists(),
 523     fPrimary(0),
 524     fCounts(0)
 525 {
 526 }
 527 //____________________________________________________________________
 528 AliForwardMCCorrectionsTask::VtxBin::VtxBin(Double_t low,
 529                                             Double_t high,
 530                                             const TAxis& axis)
 531   : TNamed(BinName(low, high),
 532            Form("%+5.1fcm<v_{z}<%+5.1fcm", low, high)),
 533     fHists(),
 534     fPrimary(0),
 535     fCounts(0)
 536 {
 537   fHists.Init(axis);
 538
 539   fPrimary = new TH2D("primary", "Primaries",
 540                       axis.GetNbins(), axis.GetXmin(), axis.GetXmax(),
 541                       40, 0, 2*TMath::Pi());
 542   fPrimary->SetXTitle("#eta");
 543   fPrimary->SetYTitle("#varphi [radians]");
 544   fPrimary->Sumw2();
 545   fPrimary->SetDirectory(0);
 546
 547   fCounts = new TH1D("counts", "Counts", 1, 0, 1);
 548   fCounts->SetXTitle("Events");
 549   fCounts->SetYTitle("# of Events");
 550   fCounts->SetDirectory(0);
 551 }
 552
 553 //____________________________________________________________________
 554 AliForwardMCCorrectionsTask::VtxBin::VtxBin(const VtxBin& o)
 555   : TNamed(o),
 556     fHists(o.fHists),
 557     fPrimary(0),
 558     fCounts(0)
 559 {
 560   if (o.fPrimary) {
 561     fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
 562     fPrimary->SetDirectory(0);
 563   }
 564   if (o.fCounts) {
 565     fCounts = static_cast<TH1D*>(o.fCounts->Clone());
 566     fCounts->SetDirectory(0);
 567   }
 568 }
 569
 570 //____________________________________________________________________
 571 AliForwardMCCorrectionsTask::VtxBin&
 572 AliForwardMCCorrectionsTask::VtxBin::operator=(const VtxBin& o)
 573 {
 574   if (&o == this) return *this;
 575   TNamed::operator=(o);
 576   fHists   = o.fHists;
 577   fPrimary = 0;
 578   fCounts  = 0;
 579   if (o.fPrimary) {
 580     fPrimary = static_cast<TH2D*>(o.fPrimary->Clone());
 581     fPrimary->SetDirectory(0);
 582   }
 583   if (o.fCounts) {
 584     fCounts = static_cast<TH1D*>(o.fCounts->Clone());
 585     fCounts->SetDirectory(0);
 586   }
 587   return *this;
 588 }
 589
 590 //____________________________________________________________________
 591 void
 592 AliForwardMCCorrectionsTask::VtxBin::DefineOutput(TList* l)
 593 {
 594   TList* d = new TList;
 595   d->SetName(GetName());
 596   d->SetOwner();
 597   l->Add(d);
 598
 599   d->Add(fHists.fFMD1i);
 600   d->Add(fHists.fFMD2i);
 601   d->Add(fHists.fFMD2o);
 602   d->Add(fHists.fFMD3i);
 603   d->Add(fHists.fFMD3o);
 604
 605   d->Add(fPrimary);
 606   d->Add(fCounts);
 607 }
 608
 609 //____________________________________________________________________
 610 TH2D*
 611 AliForwardMCCorrectionsTask::VtxBin::MakeBg(const TH2D* hits,
 612                                             const TH2D* primary) const
 613 {
 614   TH2D* h = static_cast<TH2D*>(hits->Clone());
 615   h->SetDirectory(0);
 616   TString n(h->GetName());
 617   n.ReplaceAll("_cache", "");
 618   h->SetName(n);
 619   h->Divide(primary);
 620
 621   return h;
 622 }
 623
 624 //____________________________________________________________________
 625 void
 626 AliForwardMCCorrectionsTask::VtxBin::Finish(const TList* input,
 627                                             TList* output,
 628                                             UShort_t iVz,
 629                                             AliFMDCorrSecondaryMap* map)
 630 {
 631   TList* out = new TList;
 632   out->SetName(GetName());
 633   out->SetOwner();
 634   output->Add(out);
 635
 636   TList* l = static_cast<TList*>(input->FindObject(GetName()));
 637   if (!l) {
 638     AliError(Form("List %s not found in %s", GetName(), input->GetName()));
 639     return;
 640   }
 641
 642   TH2D*   fmd1i = static_cast<TH2D*>(l->FindObject("FMD1I_cache"));
 643   TH2D*   fmd2i = static_cast<TH2D*>(l->FindObject("FMD2I_cache"));
 644   TH2D*   fmd2o = static_cast<TH2D*>(l->FindObject("FMD2O_cache"));
 645   TH2D*   fmd3i = static_cast<TH2D*>(l->FindObject("FMD3I_cache"));
 646   TH2D*   fmd3o = static_cast<TH2D*>(l->FindObject("FMD3O_cache"));
 647   TH2D*   primO = static_cast<TH2D*>(l->FindObject("primary"));
 648   if (!fmd1i || !fmd2i || !fmd2o || !fmd3i || !fmd3o || !primO) {
 649     AliError(Form("Missing histogram(s): %p,%p,%p,%p,%p,%p",
 650                   fmd1i, fmd2i, fmd2o, fmd3i, fmd3o, primO));
 651     return;
 652   }
 653
 654   // Half coverage in phi for inners
 655   TH2D*   primI = static_cast<TH2D*>(primO->Clone());
 656   primI->SetDirectory(0);
 657   primI->RebinY(2);
 658
 659   TH2D* bg1i = MakeBg(fmd1i, primI);
 660   TH2D* bg2i = MakeBg(fmd2i, primI);
 661   TH2D* bg2o = MakeBg(fmd2o, primO);
 662   TH2D* bg3i = MakeBg(fmd3i, primI);
 663   TH2D* bg3o = MakeBg(fmd3o, primO);
 664   map->SetCorrection(1, 'I', iVz, bg1i);
 665   map->SetCorrection(2, 'I', iVz, bg2i);
 666   map->SetCorrection(2, 'O', iVz, bg2o);
 667   map->SetCorrection(3, 'I', iVz, bg3i);
 668   map->SetCorrection(3, 'O', iVz, bg3o);
 669   out->Add(bg1i);
 670   out->Add(bg2i);
 671   out->Add(bg2o);
 672   out->Add(bg3i);
 673   out->Add(bg3o);
 674
 675 }
 676
 677 //
 678 // EOF
 679 //