PWGLF/FORWARD/analysis2/AliFMDSharingFilter.cxx

   1 //
   2 // Class to do the sharing correction.  That is, a filter that merges
   3 // adjacent strip signals presumably originating from a single particle
   4 // that impinges on the detector in such a way that it deposite energy
   5 // into two or more strips.
   6 //
   7 // Input:
   8 //    - AliESDFMD object  - from reconstruction
   9 //
  10 // Output:
  11 //    - AliESDFMD object  - copy of input, but with signals merged
  12 //
  13 // Corrections used:
  14 //    - AliFMDCorrELossFit
  15 //
  16 // Histograms:
  17 //    - For each ring (FMD1i, FMD2i, FMD2o, FMD3i, FMD3o) the distribution of
  18 //      signals before and after the filter.
  19 //    - For each ring (see above), an array of distributions of number of
  20 //      hit strips for each vertex bin (if enabled - see Init method)
  21 //
  22 //
  23 //
  24 #include "AliFMDSharingFilter.h"
  25 #include "AliFMDStripIndex.h"
  26 #include <AliESDFMD.h>
  27 #include <TAxis.h>
  28 #include <TList.h>
  29 #include <TH1.h>
  30 #include <TMath.h>
  31 #include "AliForwardCorrectionManager.h"
  32 #include "AliFMDCorrELossFit.h"
  33 #include <AliLog.h>
  34 #include <TROOT.h>
  35 #include <THStack.h>
  36 #include <TParameter.h>
  37 #include <iostream>
  38 #include <iomanip>
  39
  40 ClassImp(AliFMDSharingFilter)
  41 #if 0
  42 ; // This is for Emacs
  43 #endif
  44
  45 #define DBG(L,M)                                                        \
  46   do { if (L>fDebug)break; std::cout << (M) << std::flush;} while(false)
  47 #define DBGL(L,M)                                                       \
  48   do { if (L>fDebug)break; std::cout << (M) << std::endl;} while(false)
  49
  50
  51
  52 //____________________________________________________________________
  53 AliFMDSharingFilter::AliFMDSharingFilter()
  54   : TNamed(),
  55     fRingHistos(),
  56     fCorrectAngles(kFALSE),
  57     // fSummed(0),
  58     fHighCuts(0),
  59     fLowCuts(0),
  60     // fOper(0),
  61     fDebug(0),
  62     fZeroSharedHitsBelowThreshold(false),
  63     fLCuts(),
  64     fHCuts(),
  65     fUseSimpleMerging(false),
  66     fThreeStripSharing(true),
  67     fRecalculateEta(false),
  68     // fExtraDead(0),
  69     fXtraDead(0),
  70     fInvalidIsEmpty(false),
  71     fMergingDisabled(false)
  72 {
  73   //
  74   // Default Constructor - do not use
  75   //
  76   DGUARD(fDebug,1, "Default CTOR for AliFMDSharingFilter");
  77 }
  78
  79 //____________________________________________________________________
  80 AliFMDSharingFilter::AliFMDSharingFilter(const char* title)
  81   : TNamed("fmdSharingFilter", title),
  82     fRingHistos(),
  83     fCorrectAngles(kFALSE),
  84     // fSummed(0),
  85     fHighCuts(0),
  86     fLowCuts(0),
  87     // fOper(0),
  88     fDebug(0),
  89     fZeroSharedHitsBelowThreshold(false),
  90     fLCuts(),
  91     fHCuts(),
  92     fUseSimpleMerging(false),
  93     fThreeStripSharing(true),
  94     fRecalculateEta(false),
  95     // fExtraDead(51200),
  96     fXtraDead(AliFMDStripIndex::Pack(3,'O',19,511)+1),
  97     fInvalidIsEmpty(false),
  98     fMergingDisabled(false)
  99 {
 100   //
 101   // Constructor
 102   //
 103   // Parameters:
 104   //    title Title of object  - not significant
 105   //
 106   DGUARD(fDebug,1, "Named CTOR for AliFMDSharingFilter: %s", title);
 107   fRingHistos.SetName(GetName());
 108   fRingHistos.SetOwner();
 109
 110   fRingHistos.Add(new RingHistos(1, 'I'));
 111   fRingHistos.Add(new RingHistos(2, 'I'));
 112   fRingHistos.Add(new RingHistos(2, 'O'));
 113   fRingHistos.Add(new RingHistos(3, 'I'));
 114   fRingHistos.Add(new RingHistos(3, 'O'));
 115
 116   fHCuts.SetNXi(1);
 117   fHCuts.SetIncludeSigma(1);
 118   fLCuts.SetMultCuts(.15);
 119
 120   // fExtraDead.Reset(-1);
 121 }
 122
 123 //____________________________________________________________________
 124 AliFMDSharingFilter::~AliFMDSharingFilter()
 125 {
 126   //
 127   // Destructor
 128   //
 129   DGUARD(fDebug,3, "DTOR for AliFMDSharingFilter");
 130   // fRingHistos.Delete();
 131 }
 132
 133 //____________________________________________________________________
 134 AliFMDSharingFilter::RingHistos*
 135 AliFMDSharingFilter::GetRingHistos(UShort_t d, Char_t r) const
 136 {
 137   //
 138   // Get the ring histogram container
 139   //
 140   // Parameters:
 141   //    d Detector
 142   //    r Ring
 143   //
 144   // Return:
 145   //    Ring histogram container
 146   //
 147   Int_t idx = -1;
 148   switch (d) {
 149   case 1: idx = 0; break;
 150   case 2: idx = 1 + (r == 'I' || r == 'i' ? 0 : 1); break;
 151   case 3: idx = 3 + (r == 'I' || r == 'i' ? 0 : 1); break;
 152   }
 153   if (idx < 0 || idx >= fRingHistos.GetEntries()) return 0;
 154
 155   return static_cast<RingHistos*>(fRingHistos.At(idx));
 156 }
 157
 158 //____________________________________________________________________
 159 void
 160 AliFMDSharingFilter::AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t)
 161 {
 162   if (d < 1 || d > 3) {
 163     Warning("AddDead", "Invalid detector FMD%d", d);
 164     return;
 165   }
 166   Bool_t inner = (r == 'I' || r == 'i');
 167   if (d == 1 && !inner) {
 168     Warning("AddDead", "Invalid ring FMD%d%c", d, r);
 169     return;
 170   }
 171   if ((inner && s >= 20) || (!inner && s >= 40)) {
 172     Warning("AddDead", "Invalid sector FMD%d%c[%02d]", d, r, s);
 173     return;
 174   }
 175   if ((inner && t >= 512) || (!inner && t >= 256)) {
 176     Warning("AddDead", "Invalid strip FMD%d%c[%02d,%03d]", d, r, s, t);
 177     return;
 178   }
 179
 180   Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
 181   // Int_t i  = 0;
 182   fXtraDead.SetBitNumber(id, true);
 183 #if 0
 184   for (i = 0; i < fExtraDead.GetSize(); i++) {
 185     Int_t j = fExtraDead.At(i);
 186     if (j == id) return; // Already there
 187     if (j <  0) break; // Free slot
 188   }
 189   if (i >= fExtraDead.GetSize()) {
 190     Warning("AddDead", "No free slot to add FMD%d%c[%02d,%03d] at",
 191             d, r, s, t);
 192     return;
 193   }
 194   fExtraDead[i] = id;
 195 #endif
 196 }
 197 //____________________________________________________________________
 198 void
 199 AliFMDSharingFilter::AddDeadRegion(UShort_t d,  Char_t r,
 200                                    UShort_t s1, UShort_t s2,
 201                                    UShort_t t1, UShort_t t2)
 202 {
 203   // Add a dead region spanning from FMD<d><r>[<s1>,<t1>] to
 204   // FMD<d><r>[<s2>,<t2>] (both inclusive)
 205   for (Int_t s = s1; s <= s2; s++)
 206     for (Int_t t = t1; t <= t2; t++)
 207       AddDead(d, r, s, t);
 208 }
 209 //____________________________________________________________________
 210 void
 211 AliFMDSharingFilter::AddDead(const Char_t* script)
 212 {
 213   if (!script || script[0] == '\0') return;
 214
 215   gROOT->Macro(Form("%s((AliFMDSharingFilter*)%p);", script, this));
 216 }
 217
 218 //____________________________________________________________________
 219 Bool_t
 220 AliFMDSharingFilter::IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const
 221 {
 222   Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
 223   return fXtraDead.TestBitNumber(id);
 224 #if 0
 225   for (Int_t i = 0; i < fExtraDead.GetSize(); i++) {
 226     Int_t j = fExtraDead.At(i);
 227     if (j == id) {
 228       //Info("IsDead", "FMD%d%c[%02d,%03d] marked as dead here", d, r, s, t);
 229       return true;
 230     }
 231     if (j < 0) break; // High water mark
 232   }
 233   return false;
 234 #endif
 235 }
 236 //____________________________________________________________________
 237 void
 238 AliFMDSharingFilter::SetupForData(const TAxis& axis)
 239 {
 240   // Initialise - called on first event
 241   DGUARD(fDebug,1, "Initialize for AliFMDSharingFilter");
 242   AliForwardCorrectionManager& fcm  = AliForwardCorrectionManager::Instance();
 243   const AliFMDCorrELossFit*    fits = fcm.GetELossFit();
 244
 245   // Compactify the xtra dead bits
 246   fXtraDead.Compact();
 247
 248   // Get the high cut.  The high cut is defined as the
 249   // most-probably-value peak found from the energy distributions, minus
 250   // 2 times the width of the corresponding Landau.
 251
 252   TAxis eAxis(axis.GetNbins(),
 253               axis.GetXmin(),
 254               axis.GetXmax());
 255   if(fits)
 256     eAxis.Set(fits->GetEtaAxis().GetNbins(),
 257               fits->GetEtaAxis().GetXmin(),
 258               fits->GetEtaAxis().GetXmax());
 259
 260   UShort_t nEta = eAxis.GetNbins();
 261
 262   fHighCuts->SetBins(nEta, eAxis.GetXmin(), eAxis.GetXmax(), 5, .5, 5.5);
 263   fHighCuts->GetYaxis()->SetBinLabel(1, "FMD1i");
 264   fHighCuts->GetYaxis()->SetBinLabel(2, "FMD2i");
 265   fHighCuts->GetYaxis()->SetBinLabel(3, "FMD2o");
 266   fHighCuts->GetYaxis()->SetBinLabel(4, "FMD3i");
 267   fHighCuts->GetYaxis()->SetBinLabel(5, "FMD3o");
 268
 269   fLowCuts->SetBins(nEta, eAxis.GetXmin(), eAxis.GetXmax(), 5, .5, 5.5);
 270   fLowCuts->GetYaxis()->SetBinLabel(1, "FMD1i");
 271   fLowCuts->GetYaxis()->SetBinLabel(2, "FMD2i");
 272   fLowCuts->GetYaxis()->SetBinLabel(3, "FMD2o");
 273   fLowCuts->GetYaxis()->SetBinLabel(4, "FMD3i");
 274   fLowCuts->GetYaxis()->SetBinLabel(5, "FMD3o");
 275
 276   // Cache our cuts in histograms
 277   fLCuts.FillHistogram(fLowCuts);
 278   fHCuts.FillHistogram(fHighCuts);
 279 }
 280
 281 //____________________________________________________________________
 282 #define ETA2COS(ETA)                                            \
 283   TMath::Cos(2*TMath::ATan(TMath::Exp(-TMath::Abs(ETA))))
 284
 285 Bool_t
 286 AliFMDSharingFilter::Filter(const AliESDFMD& input,
 287                             Bool_t           /*lowFlux*/,
 288                             AliESDFMD&       output,
 289                             Double_t         zvtx)
 290 {
 291   //
 292   // Filter the input AliESDFMD object
 293   //
 294   // Parameters:
 295   //    input     Input
 296   //    lowFlux   If this is a low-flux event
 297   //    output    Output AliESDFMD object
 298   //
 299   // Return:
 300   //    True on success, false otherwise
 301   //
 302   DGUARD(fDebug,1, "Filter event in AliFMDSharingFilter");
 303   output.Clear();
 304   TIter    next(&fRingHistos);
 305   RingHistos* o      = 0;
 306   while ((o = static_cast<RingHistos*>(next()))) o->Clear();
 307
 308   Int_t nSingle    = 0;
 309   Int_t nDouble    = 0;
 310   Int_t nTriple    = 0;
 311
 312   for(UShort_t d = 1; d <= 3; d++) {
 313     Int_t nRings = (d == 1 ? 1 : 2);
 314     for (UShort_t q = 0; q < nRings; q++) {
 315       Char_t      r      = (q == 0 ? 'I' : 'O');
 316       UShort_t    nsec   = (q == 0 ?  20 :  40);
 317       UShort_t    nstr   = (q == 0 ? 512 : 256);
 318       RingHistos* histos = GetRingHistos(d, r);
 319
 320       for(UShort_t s = 0; s < nsec;  s++) {
 321         // `used' flags if the _current_ strip was used by _previous_
 322         // iteration.
 323         Bool_t   used            = kFALSE;
 324         // `eTotal' contains the current sum of merged signals so far
 325         Double_t eTotal          = -1;
 326         // Int_t    nDistanceBefore = -1;
 327         // Int_t    nDistanceAfter  = -1;
 328         // `twoLow' flags if we saw two consequtive strips with a
 329         // signal between the two cuts.
 330         Bool_t   twoLow          = kFALSE;
 331         Int_t    nStripsAboveCut = 0;
 332
 333         for(UShort_t t = 0; t < nstr; t++) {
 334           // nDistanceBefore++;
 335           // nDistanceAfter++;
 336
 337           output.SetMultiplicity(d,r,s,t,0.);
 338           Float_t mult         = SignalInStrip(input,d,r,s,t);
 339           Float_t multNext     = (t<nstr-1) ? SignalInStrip(input,d,r,s,t+1) :0;
 340           Float_t multNextNext = (t<nstr-2) ? SignalInStrip(input,d,r,s,t+2) :0;
 341           if (multNext     ==  AliESDFMD::kInvalidMult) multNext     = 0;
 342           if (multNextNext ==  AliESDFMD::kInvalidMult) multNextNext = 0;
 343           if(!fThreeStripSharing) multNextNext = 0;
 344
 345           // Get the pseudo-rapidity
 346           Double_t eta = input.Eta(d,r,s,t);
 347           Double_t phi = input.Phi(d,r,s,t) * TMath::Pi() / 180.;
 348           if (s == 0) output.SetEta(d,r,s,t,eta);
 349
 350           if(fRecalculateEta) {
 351             Double_t etaOld  = eta;
 352             Double_t etaCalc = AliForwardUtil::GetEtaFromStrip(d,r,s,t,zvtx);
 353             eta              = etaCalc;
 354
 355             if (mult > 0 && mult != AliESDFMD::kInvalidMult ) {
 356               Double_t cosOld =  ETA2COS(etaOld);
 357               Double_t cosNew =  ETA2COS(etaCalc);
 358               Double_t corr   =  cosNew / cosOld;
 359               mult            *= corr;
 360               multNext        *= corr;
 361               multNextNext    *= corr;
 362             }
 363           } // Recalculate eta
 364
 365           // Special case for pre revision 43611 AliFMDReconstructor.
 366           // If fInvalidIsEmpty and we get an invalid signal from the
 367           // ESD, then we need to set this signal to zero.  Note, dead
 368           // strips added in the ForwardAODConfig.C file are not
 369           // effected by this, and we can use that to set the proper
 370           // dead strips.
 371           if (mult == AliESDFMD::kInvalidMult && fInvalidIsEmpty)
 372             mult = 0;
 373
 374           // Keep dead-channel information - either from the ESD (but
 375           // see above for older data) or from the settings in the
 376           // ForwardAODConfig.C file.
 377           if (mult == AliESDFMD::kInvalidMult || IsDead(d,r,s,t)) {
 378             output.SetMultiplicity(d,r,s,t,AliESDFMD::kInvalidMult);
 379             histos->fBefore->Fill(-1);
 380             mult = AliESDFMD::kInvalidMult;
 381           }
 382
 383           if (mult != AliESDFMD::kInvalidMult)
 384             // Always fill the ESD sum histogram
 385             histos->fSumESD->Fill(eta, phi, mult);
 386
 387           // If no signal or dead strip, go on.
 388           if (mult == AliESDFMD::kInvalidMult || mult == 0) {
 389             if (mult == 0) histos->fSum->Fill(eta,phi,mult);
 390             // Flush a possible signal
 391             if (eTotal > 0 && t > 0)
 392               output.SetMultiplicity(d,r,s,t-1,eTotal);
 393             // Reset states so we do not try to merge over a dead strip.
 394             eTotal = -1;
 395             used   = false;
 396             twoLow = false;
 397             if (t > 0)
 398               histos->fNConsecutive->Fill(nStripsAboveCut);
 399             if (mult == AliESDFMD::kInvalidMult)
 400               // Why not fill immidiately here?
 401               nStripsAboveCut = -1;
 402             else
 403               // Why not fill immidiately here?
 404               nStripsAboveCut = 0;
 405             continue;
 406           }
 407
 408           // Fill the diagnostics histogram
 409           histos->fBefore->Fill(mult);
 410
 411           Double_t mergedEnergy = mult;
 412           // it seems to me that this logic could be condensed a bit
 413           if(mult > GetLowCut(d, r, eta)) {
 414             if(nStripsAboveCut < 1) {
 415               if(t > 0)
 416                 histos->fNConsecutive->Fill(nStripsAboveCut);
 417               nStripsAboveCut=0;
 418             }
 419             nStripsAboveCut++;
 420           }
 421           else {
 422             if (t > 0)
 423               histos->fNConsecutive->Fill(nStripsAboveCut);
 424             nStripsAboveCut=0;
 425           }
 426
 427           if (!fMergingDisabled) {
 428             mergedEnergy = 0;
 429
 430             // The current sum
 431             Float_t etot = 0;
 432
 433             // Fill in neighbor information
 434             if (t < nstr-1) histos->fNeighborsBefore->Fill(mult,multNext);
 435
 436             Bool_t thisValid = mult     > GetLowCut(d, r, eta);
 437             Bool_t nextValid = multNext > GetLowCut(d, r, eta);
 438             Bool_t thisSmall = mult     < GetHighCut(d, r, eta ,false);
 439             Bool_t nextSmall = multNext < GetHighCut(d, r, eta ,false);
 440
 441             // If this strips signal is above the high cut, reset distance
 442             // if (!thisSmall) {
 443             //    histos->fDistanceBefore->Fill(nDistanceBefore);
 444             //    nDistanceBefore = -1;
 445             // }
 446
 447             // If the total signal in the past 1 or 2 strips are non-zero
 448             // we need to check
 449             if (eTotal > 0) {
 450               // Here, we have already flagged one strip as a candidate
 451
 452               // If 3-strip merging is enabled, then check the next
 453               // strip to see that it falls within cut, or if we have
 454               // two low signals
 455               if (fThreeStripSharing && nextValid && (nextSmall || twoLow)) {
 456                 eTotal = eTotal + multNext;
 457                 used = kTRUE;
 458                 histos->fTriple->Fill(eTotal);
 459                 nTriple++;
 460                 twoLow = kFALSE;
 461               }
 462               // Otherwise, we got a double hit before, and that
 463               // should be stored.
 464               else {
 465                 used = kFALSE;
 466                 histos->fDouble->Fill(eTotal);
 467                 nDouble++;
 468               }
 469               // Store energy loss and reset sum
 470               etot   = eTotal;
 471               eTotal = -1;
 472             } // if (eTotal>0)
 473             else {
 474               // If we have no current sum
 475
 476               // Check if this is marked as used, and if so, continue
 477               if (used) {used = kFALSE; continue; }
 478
 479               // If the signal is abvoe the cut, set current
 480               if (thisValid) etot = mult;
 481
 482               // If the signal is abiove the cut, and so is the next
 483               // signal and either of them are below the high cut,
 484               if (thisValid  && nextValid  && (thisSmall || nextSmall)) {
 485
 486                 // If this is below the high cut, and the next is too, then
 487                 // we have two low signals
 488                 if (thisSmall && nextSmall) twoLow = kTRUE;
 489
 490                 // If this signal is bigger than the next, and the
 491                 // one after that is below the low-cut, then update
 492                 // the sum
 493                 if (mult>multNext && multNextNext < GetLowCut(d, r, eta)) {
 494                   etot = mult + multNext;
 495                   used = kTRUE;
 496                   histos->fDouble->Fill(etot);
 497                   nDouble++;
 498                 }
 499                 // Otherwise, we may need to merge with a third strip
 500                 else {
 501                   etot   = 0;
 502                   eTotal = mult + multNext;
 503                 }
 504               }
 505               // This is a signle hit
 506               else if(etot > 0) {
 507                 histos->fSingle->Fill(etot);
 508                 histos->fSinglePerStrip->Fill(etot,t);
 509                 nSingle++;
 510               }
 511             } // else if (etotal >= 0)
 512
 513             mergedEnergy = etot;
 514             // if (mergedEnergy > GetHighCut(d, r, eta ,false)) {
 515             //   histos->fDistanceAfter->Fill(nDistanceAfter);
 516             //   nDistanceAfter    = -1;
 517             // }
 518             //if(mult>0 && multNext >0)
 519             //  std::cout<<mult<<"  "<<multNext<<"  "<<mergedEnergy<<std::endl;
 520           } // if (!fMergingDisabled)
 521
 522           if (!fCorrectAngles)
 523             mergedEnergy = AngleCorrect(mergedEnergy, eta);
 524           // if (mergedEnergy > 0) histos->Incr();
 525
 526           if (t != 0)
 527             histos->fNeighborsAfter->Fill(output.Multiplicity(d,r,s,t-1),
 528                                           mergedEnergy);
 529           histos->fBeforeAfter->Fill(mult, mergedEnergy);
 530           if(mergedEnergy > 0)
 531             histos->fAfter->Fill(mergedEnergy);
 532           histos->fSum->Fill(eta,phi,mergedEnergy);
 533
 534           output.SetMultiplicity(d,r,s,t,mergedEnergy);
 535         } // for strip
 536         histos->fNConsecutive->Fill(nStripsAboveCut); // fill the last sector
 537       } // for sector
 538     } // for ring
 539   } // for detector
 540   DMSG(fDebug, 3,"single=%9d, double=%9d, triple=%9d",
 541        nSingle, nDouble, nTriple);
 542   next.Reset();
 543   // while ((o = static_cast<RingHistos*>(next()))) o->Finish();
 544
 545   return kTRUE;
 546 }
 547
 548 //_____________________________________________________________________
 549 Double_t
 550 AliFMDSharingFilter::SignalInStrip(const AliESDFMD& input,
 551                                    UShort_t         d,
 552                                    Char_t           r,
 553                                    UShort_t         s,
 554                                    UShort_t         t) const
 555 {
 556   //
 557   // Get the signal in a strip
 558   //
 559   // Parameters:
 560   //    fmd   ESD object
 561   //    d     Detector
 562   //    r     Ring
 563   //    s     Sector
 564   //    t     Strip
 565   //
 566   // Return:
 567   //    The energy signal
 568   //
 569   Double_t mult = input.Multiplicity(d,r,s,t);
 570   // In case of
 571   //  - bad value (invalid or 0)
 572   //  - we want angle corrected and data is
 573   //  - we don't want angle corrected and data isn't
 574   // just return read value
 575   if (mult == AliESDFMD::kInvalidMult               ||
 576       mult == 0                                     ||
 577       (fCorrectAngles && input.IsAngleCorrected()) ||
 578       (!fCorrectAngles && !input.IsAngleCorrected()))
 579     return mult;
 580
 581   // If we want angle corrected data, correct it,
 582   // otherwise de-correct it
 583   if (fCorrectAngles) mult = AngleCorrect(mult, input.Eta(d,r,s,t));
 584   else                mult = DeAngleCorrect(mult, input.Eta(d,r,s,t));
 585   return mult;
 586 }
 587
 588 namespace {
 589   Double_t Rng2Cut(UShort_t d, Char_t r, Double_t eta, TH2* h) {
 590     Double_t ret = 1024;
 591     Int_t ybin = 0;
 592     switch(d) {
 593     case 1: ybin = 1; break;
 594     case 2: ybin = (r=='i' || r=='I') ? 2 : 3; break;
 595     case 3: ybin = (r=='i' || r=='I') ? 4 : 5; break;
 596     default: return ret;
 597     }
 598     Int_t xbin = h->GetXaxis()->FindBin(eta);
 599     if (xbin < 1 && xbin > h->GetXaxis()->GetNbins()) return ret;
 600     ret = h->GetBinContent(xbin,ybin);
 601     return ret;
 602   }
 603 }
 604
 605   //_____________________________________________________________________
 606 Double_t
 607 AliFMDSharingFilter::GetLowCut(UShort_t d, Char_t r, Double_t eta) const
 608 {
 609   //
 610   // Get the low cut.  Normally, the low cut is taken to be the lower
 611   // value of the fit range used when generating the energy loss fits.
 612   // However, if fLowCut is set (using SetLowCit) to a value greater
 613   // than 0, then that value is used.
 614   //
 615   return Rng2Cut(d, r, eta, fLowCuts);
 616   // return fLCuts.GetMultCut(d,r,eta,false);
 617 }
 618
 619 //_____________________________________________________________________
 620 Double_t
 621 AliFMDSharingFilter::GetHighCut(UShort_t d, Char_t r,
 622                                 Double_t eta, Bool_t /*errors*/) const
 623 {
 624   //
 625   // Get the high cut.  The high cut is defined as the
 626   // most-probably-value peak found from the energy distributions, minus
 627   // 2 times the width of the corresponding Landau.
 628   //
 629   return Rng2Cut(d, r, eta, fHighCuts);
 630   // return fHCuts.GetMultCut(d,r,eta,errors);
 631 }
 632
 633 //____________________________________________________________________
 634 Double_t
 635 AliFMDSharingFilter::AngleCorrect(Double_t mult, Double_t eta) const
 636 {
 637   //
 638   // Angle correct the signal
 639   //
 640   // Parameters:
 641   //    mult Angle Un-corrected Signal
 642   //    eta  Pseudo-rapidity
 643   //
 644   // Return:
 645   //    Angle corrected signal
 646   //
 647   Double_t theta =  2 * TMath::ATan(TMath::Exp(-eta));
 648   if (eta < 0) theta -= TMath::Pi();
 649   return mult * TMath::Cos(theta);
 650 }
 651 //____________________________________________________________________
 652 Double_t
 653 AliFMDSharingFilter::DeAngleCorrect(Double_t mult, Double_t eta) const
 654 {
 655   //
 656   // Angle de-correct the signal
 657   //
 658   // Parameters:
 659   //    mult Angle corrected Signal
 660   //    eta  Pseudo-rapidity
 661   //
 662   // Return:
 663   //    Angle un-corrected signal
 664   //
 665   Double_t theta =  2 * TMath::ATan(TMath::Exp(-eta));
 666   if (eta < 0) theta -= TMath::Pi();
 667   return mult / TMath::Cos(theta);
 668 }
 669
 670 //____________________________________________________________________
 671 void
 672 AliFMDSharingFilter::Terminate(const TList* dir, TList* output, Int_t nEvents)
 673 {
 674   //
 675   // Scale the histograms to the total number of events
 676   //
 677   // Parameters:
 678   //    dir     Where the output is
 679   //    nEvents Number of events
 680   //
 681   DGUARD(fDebug,1, "Scale histograms in AliFMDSharingFilter");
 682   if (nEvents <= 0) return;
 683   TList* d = static_cast<TList*>(dir->FindObject(GetName()));
 684   if (!d) return;
 685
 686   TList* out = new TList;
 687   out->SetName(d->GetName());
 688   out->SetOwner();
 689
 690   TParameter<int>* nFiles =
 691     static_cast<TParameter<int>*>(d->FindObject("nFiles"));
 692
 693   TH2* lowCuts  = static_cast<TH2*>(d->FindObject("lowCuts"));
 694   TH2* highCuts = static_cast<TH2*>(d->FindObject("highCuts"));
 695   if (lowCuts && nFiles) {
 696     lowCuts->Scale(1. / nFiles->GetVal());
 697     out->Add(lowCuts->Clone());
 698   }
 699   else
 700     AliWarning("low cuts histogram not found in input list");
 701   if (highCuts && nFiles) {
 702     highCuts->Scale(1. / nFiles->GetVal());
 703     out->Add(highCuts->Clone());
 704   }
 705   else
 706     AliWarning("high cuts histogram not found in input list");
 707
 708   TIter    next(&fRingHistos);
 709   RingHistos* o = 0;
 710   THStack* sums = new THStack("sums", "Sum of ring signals");
 711   THStack* sumsESD = new THStack("sumsESD", "Sum of ring ESD signals");
 712   while ((o = static_cast<RingHistos*>(next()))) {
 713     o->Terminate(d, nEvents);
 714     if (!o->fSum) {
 715       Warning("Terminate", "No sum histogram found for ring %s", o->GetName());
 716       continue;
 717     }
 718     TH1D* sum = o->fSum->ProjectionX(o->GetName(), 1, o->fSum->GetNbinsY(),"e");
 719     sum->Scale(1., "width");
 720     sum->SetTitle(o->GetName());
 721     sum->SetDirectory(0);
 722     sum->SetYTitle("#sum_{c} #Delta/#Delta_{mip}");
 723     sums->Add(sum);
 724
 725     sum = o->fSumESD->ProjectionX(o->GetName(), 1, o->fSumESD->GetNbinsY(),"e");
 726     sum->Scale(1., "width");
 727     sum->SetTitle(o->GetName());
 728     sum->SetDirectory(0);
 729     sum->SetYTitle("#sum_{s} #Delta/#Delta_{mip}");
 730     sumsESD->Add(sum);
 731   }
 732   out->Add(sums);
 733   out->Add(sumsESD);
 734   output->Add(out);
 735 }
 736
 737 //____________________________________________________________________
 738 void
 739 AliFMDSharingFilter::CreateOutputObjects(TList* dir)
 740 {
 741   //
 742   // Define the output histograms.  These are put in a sub list of the
 743   // passed list.   The histograms are merged before the parent task calls
 744   // AliAnalysisTaskSE::Terminate
 745   //
 746   // Parameters:
 747   //    dir Directory to add to
 748   //
 749   DGUARD(fDebug,1, "Define output in AliFMDSharingFilter");
 750   TList* d = new TList;
 751   d->SetName(GetName());
 752   dir->Add(d);
 753
 754 #if 0
 755   fSummed = new TH2I("operations", "Strip operations",
 756                      kMergedInto, kNone-.5, kMergedInto+.5,
 757                      51201, -.5, 51200.5);
 758   fSummed->SetXTitle("Operation");
 759   fSummed->SetYTitle("# of strips");
 760   fSummed->SetZTitle("Events");
 761   fSummed->GetXaxis()->SetBinLabel(kNone,            "None");
 762   fSummed->GetXaxis()->SetBinLabel(kCandidate,       "Candidate");
 763   fSummed->GetXaxis()->SetBinLabel(kMergedWithOther, "Merged w/other");
 764   fSummed->GetXaxis()->SetBinLabel(kMergedInto,      "Merged into");
 765   fSummed->SetDirectory(0);
 766   d->Add(fSummed);
 767 #endif
 768
 769   fHighCuts = new TH2D("highCuts", "High cuts used", 1,0,1, 1,0,1);
 770   fHighCuts->SetXTitle("#eta");
 771   fHighCuts->SetDirectory(0);
 772   d->Add(fHighCuts);
 773
 774   fLowCuts = new TH2D("lowCuts", "Low cuts used", 1,0,1, 1,0,1);
 775   fLowCuts->SetXTitle("#eta");
 776   fLowCuts->SetDirectory(0);
 777   d->Add(fLowCuts);
 778
 779   // d->Add(lowCut);
 780   // d->Add(nXi);
 781   // d->Add(sigma);
 782   d->Add(AliForwardUtil::MakeParameter("angle", fCorrectAngles));
 783   d->Add(AliForwardUtil::MakeParameter("lowSignal",
 784                                        fZeroSharedHitsBelowThreshold));
 785   d->Add(AliForwardUtil::MakeParameter("simple", fUseSimpleMerging));
 786   d->Add(AliForwardUtil::MakeParameter("sumThree", fThreeStripSharing));
 787   d->Add(AliForwardUtil::MakeParameter("disabled", fMergingDisabled));
 788   TParameter<int>* nFiles = new TParameter<int>("nFiles", 1);
 789   nFiles->SetMergeMode('+');
 790   d->Add(nFiles);
 791
 792   TObjArray* extraDead = new TObjArray;
 793   extraDead->SetOwner();
 794   extraDead->SetName("extraDead");
 795 #if 0
 796   for (Int_t i = 0; i < fExtraDead.GetSize(); i++) {
 797     if (fExtraDead.At(i) < 0) break;
 798     UShort_t dd, s, t;
 799     Char_t  r;
 800     Int_t   id = fExtraDead.At(i);
 801     AliFMDStripIndex::Unpack(id, dd, r, s, t);
 802     extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
 803                                                       dd, r, s, t), id));
 804   }
 805 #endif
 806   fXtraDead.Compact();
 807   UInt_t firstBit = fXtraDead.FirstSetBit();
 808   UInt_t nBits    = fXtraDead.GetNbits();
 809   for (UInt_t i = firstBit; i < nBits; i++) {
 810     if (!fXtraDead.TestBitNumber(i)) continue;
 811     UShort_t dd, s, t;
 812     Char_t  r;
 813     AliFMDStripIndex::Unpack(i, dd, r, s, t);
 814     extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
 815                                                       dd, r, s, t),
 816                                                  UShort_t(i)));
 817   }
 818   // d->Add(&fXtraDead);
 819   d->Add(extraDead);
 820   fLCuts.Output(d,"lCuts");
 821   fHCuts.Output(d,"hCuts");
 822
 823   TIter    next(&fRingHistos);
 824   RingHistos* o = 0;
 825   while ((o = static_cast<RingHistos*>(next()))) {
 826     o->CreateOutputObjects(d);
 827   }
 828 }
 829 #define PF(N,V,...)                                     \
 830   AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
 831 #define PFB(N,FLAG)                             \
 832   do {                                                                  \
 833     AliForwardUtil::PrintName(N);                                       \
 834     std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
 835   } while(false)
 836 #define PFV(N,VALUE)                                    \
 837   do {                                                  \
 838     AliForwardUtil::PrintName(N);                       \
 839     std::cout << (VALUE) << std::endl; } while(false)
 840
 841 //____________________________________________________________________
 842 void
 843 AliFMDSharingFilter::Print(Option_t* /*option*/) const
 844 {
 845   //
 846   // Print information
 847   //
 848   // Parameters:
 849   //    option Not used
 850   //
 851   AliForwardUtil::PrintTask(*this);
 852   gROOT->IncreaseDirLevel();
 853
 854   PFB("Use corrected angles",  fCorrectAngles);
 855   PFB("Zero below threshold",  fZeroSharedHitsBelowThreshold);
 856   PFB("Use simple sharing",    fUseSimpleMerging);
 857   PFB("Consider invalid null", fInvalidIsEmpty);
 858   PFB("Allow 3 strip merging", fThreeStripSharing);
 859   PF("Low cuts",        "");
 860   fLCuts.Print();
 861   PF("High cuts",       "");
 862   fHCuts.Print();
 863   gROOT->DecreaseDirLevel();
 864 }
 865
 866 //====================================================================
 867 AliFMDSharingFilter::RingHistos::RingHistos()
 868   : AliForwardUtil::RingHistos(),
 869     fBefore(0),
 870     fAfter(0),
 871     fSingle(0),
 872     fDouble(0),
 873     fTriple(0),
 874     fSinglePerStrip(0),
 875     // fDistanceBefore(0),
 876     // fDistanceAfter(0),
 877     fBeforeAfter(0),
 878     fNeighborsBefore(0),
 879     fNeighborsAfter(0),
 880     fSumESD(0),
 881     fSum(0),
 882     fNConsecutive(0)
 883      // ,
 884     // fHits(0),
 885     // fNHits(0)
 886 {
 887   //
 888   // Default CTOR
 889   //
 890
 891 }
 892
 893 //____________________________________________________________________
 894 AliFMDSharingFilter::RingHistos::RingHistos(UShort_t d, Char_t r)
 895   : AliForwardUtil::RingHistos(d,r),
 896     fBefore(0),
 897     fAfter(0),
 898     fSingle(0),
 899     fDouble(0),
 900     fTriple(0),
 901     fSinglePerStrip(0),
 902     // fDistanceBefore(0),
 903     // fDistanceAfter(0),
 904     fBeforeAfter(0),
 905     fNeighborsBefore(0),
 906     fNeighborsAfter(0),
 907     fSumESD(0),
 908     fSum(0),
 909     fNConsecutive(0)
 910      //,
 911     // fHits(0),
 912     // fNHits(0)
 913 {
 914   //
 915   // Constructor
 916   //
 917   // Parameters:
 918   //    d detector
 919   //    r ring
 920   //
 921   fBefore = new TH1D("esdEloss", Form("Energy loss in %s (reconstruction)",
 922                                       GetName()), 640, -1, 15);
 923   fBefore->SetXTitle("#Delta E/#Delta E_{mip}");
 924   fBefore->SetYTitle("P(#Delta E/#Delta E_{mip})");
 925   fBefore->SetFillColor(Color());
 926   fBefore->SetFillStyle(3001);
 927   fBefore->SetLineColor(kBlack);
 928   fBefore->SetLineStyle(2);
 929   fBefore->SetDirectory(0);
 930
 931   fAfter  = static_cast<TH1D*>(fBefore->Clone("anaEloss"));
 932   fAfter->SetTitle(Form("Energy loss in %s (sharing corrected)", GetName()));
 933   fAfter->SetFillColor(Color()+2);
 934   fAfter->SetLineStyle(1);
 935   fAfter->SetDirectory(0);
 936
 937   fSingle = new TH1D("singleEloss", "Energy loss (single strips)",
 938                      600, 0, 15);
 939   fSingle->SetXTitle("#Delta/#Delta_{mip}");
 940   fSingle->SetYTitle("P(#Delta/#Delta_{mip})");
 941   fSingle->SetFillColor(Color());
 942   fSingle->SetFillStyle(3001);
 943   fSingle->SetLineColor(kBlack);
 944   fSingle->SetLineStyle(2);
 945   fSingle->SetDirectory(0);
 946
 947   fDouble = static_cast<TH1D*>(fSingle->Clone("doubleEloss"));
 948   fDouble->SetTitle("Energy loss (two strips)");
 949   fDouble->SetFillColor(Color()+1);
 950   fDouble->SetDirectory(0);
 951
 952   fTriple = static_cast<TH1D*>(fSingle->Clone("tripleEloss"));
 953   fTriple->SetTitle("Energy loss (three strips)");
 954   fTriple->SetFillColor(Color()+2);
 955   fTriple->SetDirectory(0);
 956
 957   //Int_t nBinsForInner = (r == 'I' ? 32 : 16);
 958   Int_t nBinsForInner = (r == 'I' ? 512 : 256);
 959   Int_t nStrips       = (r == 'I' ? 512 : 256);
 960
 961   fSinglePerStrip = new TH2D("singlePerStrip", "SinglePerStrip",
 962                              600,0,15, nBinsForInner,0,nStrips);
 963   fSinglePerStrip->SetXTitle("#Delta/#Delta_{mip}");
 964   fSinglePerStrip->SetYTitle("Strip #");
 965   fSinglePerStrip->SetZTitle("Counts");
 966   fSinglePerStrip->SetDirectory(0);
 967
 968 #if 0
 969   fDistanceBefore = new TH1D("distanceBefore", "Distance before sharing",
 970                              nStrips , 0,nStrips );
 971   fDistanceBefore->SetXTitle("Distance");
 972   fDistanceBefore->SetYTitle("Counts");
 973   fDistanceBefore->SetFillColor(kGreen+2);
 974   fDistanceBefore->SetFillStyle(3001);
 975   fDistanceBefore->SetLineColor(kBlack);
 976   fDistanceBefore->SetLineStyle(2);
 977   fDistanceBefore->SetDirectory(0);
 978
 979   fDistanceAfter = static_cast<TH1D*>(fDistanceBefore->Clone("distanceAfter"));
 980   fDistanceAfter->SetTitle("Distance after sharing");
 981   fDistanceAfter->SetFillColor(kGreen+1);
 982   fDistanceAfter->SetDirectory(0);
 983 #endif
 984
 985   Double_t max = 15;
 986   Double_t min = -1;
 987   Int_t    n   = int((max-min) / (max / 300));
 988   fBeforeAfter = new TH2D("beforeAfter", "Before and after correlation",
 989                           n, min, max, n, min, max);
 990   fBeforeAfter->SetXTitle("#Delta E/#Delta E_{mip} before");
 991   fBeforeAfter->SetYTitle("#Delta E/#Delta E_{mip} after");
 992   fBeforeAfter->SetZTitle("Correlation");
 993   fBeforeAfter->SetDirectory(0);
 994
 995   fNeighborsBefore = static_cast<TH2D*>(fBeforeAfter->Clone("neighborsBefore"));
 996   fNeighborsBefore->SetTitle("Correlation of neighbors before");
 997   fNeighborsBefore->SetXTitle("#Delta E_{i}/#Delta E_{mip}");
 998   fNeighborsBefore->SetYTitle("#Delta E_{i+1}/#Delta E_{mip}");
 999   fNeighborsBefore->SetDirectory(0);
1000
1001   fNeighborsAfter =
1002     static_cast<TH2D*>(fNeighborsBefore->Clone("neighborsAfter"));
1003   fNeighborsAfter->SetTitle("Correlation of neighbors after");
1004   fNeighborsAfter->SetDirectory(0);
1005
1006   fSumESD = new TH2D("summedESD", "Summed ESD signal", 200, -4, 6,
1007                   NSector(), 0, 2*TMath::Pi());
1008   fSumESD->SetDirectory(0);
1009   fSumESD->Sumw2();
1010   fSumESD->SetMarkerColor(Color());
1011   // fSum->SetFillColor(Color());
1012   fSumESD->SetXTitle("#eta");
1013   fSumESD->SetYTitle("#varphi [radians]");
1014   fSumESD->SetZTitle("#sum_{strip} #Delta/#Delta_{mip}(#eta,#varphi) ");
1015
1016   fSum = static_cast<TH2D*>(fSumESD->Clone("summed"));
1017   fSum->SetTitle("Summed cluster signal");
1018   fSum->SetZTitle("#sum_{cluster} #Delta/#Delta_{mip}(#eta,#varphi) ");
1019   fSum->SetDirectory(0);
1020
1021   // Perhaps we need to ensure that this histogram has enough range to
1022   // accommondate all possible ranges - that is, from -1 to the number
1023   // of strips in this ring(-type) - i.e., NStrips().  Perhaps the
1024   // axis should be defined with increasin bin size - e.g.,
1025   //
1026   //   -1.5,-.5,.5,1.5,...,100.5,128.5,192.5,...,NStrips()
1027   //
1028   fNConsecutive = new TH1D("nConsecutive","# consecutive strips above low cut",
1029                            201,-1.5,199.5);
1030   fNConsecutive->SetXTitle("N_{strips}");
1031   fNConsecutive->SetYTitle("N_{entries}");
1032   fNConsecutive->SetFillColor(kYellow+2);
1033   fNConsecutive->SetFillStyle(3001);
1034   fNConsecutive->SetDirectory(0);
1035
1036
1037 #if 0
1038   fHits = new TH1D("hits", "Number of hits", 200, 0, 200000);
1039   fHits->SetDirectory(0);
1040   fHits->SetXTitle("# of hits");
1041   fHits->SetFillColor(kGreen+1);
1042 #endif
1043 }
1044 //____________________________________________________________________
1045 AliFMDSharingFilter::RingHistos::RingHistos(const RingHistos& o)
1046   : AliForwardUtil::RingHistos(o),
1047     fBefore(o.fBefore),
1048     fAfter(o.fAfter),
1049     fSingle(o.fSingle),
1050     fDouble(o.fDouble),
1051     fTriple(o.fTriple),
1052     fSinglePerStrip(o.fSinglePerStrip),
1053     // fDistanceBefore(o.fDistanceBefore),
1054     // fDistanceAfter(o.fDistanceAfter),
1055     fBeforeAfter(o.fBeforeAfter),
1056     fNeighborsBefore(o.fNeighborsBefore),
1057     fNeighborsAfter(o.fNeighborsAfter),
1058     fSumESD(o.fSumESD), //,
1059     fSum(o.fSum),
1060     fNConsecutive(o.fNConsecutive)
1061      //,
1062     // fHits(o.fHits),
1063     // fNHits(o.fNHits)
1064 {
1065   //
1066   // Copy constructor
1067   //
1068   // Parameters:
1069   //    o Object to copy from
1070   //
1071 }
1072
1073 //____________________________________________________________________
1074 AliFMDSharingFilter::RingHistos&
1075 AliFMDSharingFilter::RingHistos::operator=(const RingHistos& o)
1076 {
1077   //
1078   // Assignment operator
1079   //
1080   // Parameters:
1081   //    o Object to assign from
1082   //
1083   // Return:
1084   //    Reference to this
1085   //
1086   if (&o == this) return *this;
1087   AliForwardUtil::RingHistos::operator=(o);
1088   fDet = o.fDet;
1089   fRing = o.fRing;
1090
1091   if (fBefore)         delete  fBefore;
1092   if (fAfter)          delete  fAfter;
1093   if (fSingle)         delete  fSingle;
1094   if (fDouble)         delete  fDouble;
1095   if (fTriple)         delete  fTriple;
1096   if (fSinglePerStrip) delete  fSinglePerStrip;
1097   if (fNConsecutive)   delete  fNConsecutive;
1098   // if (fDistanceBefore) delete fDistanceBefore;
1099   // if (fDistanceAfter)  delete fDistanceAfter;
1100   // if (fHits)                delete fHits;
1101
1102
1103   fBefore          = static_cast<TH1D*>(o.fBefore->Clone());
1104   fAfter           = static_cast<TH1D*>(o.fAfter->Clone());
1105   fSingle          = static_cast<TH1D*>(o.fSingle->Clone());
1106   fDouble          = static_cast<TH1D*>(o.fDouble->Clone());
1107   fTriple          = static_cast<TH1D*>(o.fTriple->Clone());
1108   fSinglePerStrip  = static_cast<TH2D*>(o.fSinglePerStrip->Clone());
1109   // fDistanceBefore  = static_cast<TH1D*>(o.fDistanceBefore->Clone());
1110   // fDistanceAfter   = static_cast<TH1D*>(o.fDistanceAfter->Clone());
1111   fBeforeAfter     = static_cast<TH2D*>(o.fBeforeAfter->Clone());
1112   fNeighborsBefore = static_cast<TH2D*>(o.fNeighborsBefore->Clone());
1113   fNeighborsAfter  = static_cast<TH2D*>(o.fNeighborsAfter->Clone());
1114   // fHits            = static_cast<TH1D*>(o.fHits->Clone());
1115   fSumESD          = static_cast<TH2D*>(o.fSumESD->Clone());
1116   fSum             = static_cast<TH2D*>(o.fSum->Clone());
1117   fNConsecutive    = static_cast<TH1D*>(o.fNConsecutive->Clone());
1118
1119   return *this;
1120 }
1121 //____________________________________________________________________
1122 AliFMDSharingFilter::RingHistos::~RingHistos()
1123 {
1124   //
1125   // Destructor
1126   //
1127 }
1128 #if 0
1129 //____________________________________________________________________
1130 void
1131 AliFMDSharingFilter::RingHistos::Finish()
1132 {
1133   //
1134   // Finish off
1135   //
1136   //
1137   // fHits->Fill(fNHits);
1138 }
1139 #endif
1140 //____________________________________________________________________
1141 void
1142 AliFMDSharingFilter::RingHistos::Terminate(const TList* dir, Int_t nEvents)
1143 {
1144   //
1145   // Scale the histograms to the total number of events
1146   //
1147   // Parameters:
1148   //    nEvents Number of events
1149   //    dir     Where the output is
1150   //
1151   TList* l = GetOutputList(dir);
1152   if (!l) return;
1153
1154 #if 0
1155   TH1D* before = static_cast<TH1D*>(l->FindObject("esdEloss"));
1156   TH1D* after  = static_cast<TH1D*>(l->FindObject("anaEloss"));
1157   if (before) before->Scale(1./nEvents);
1158   if (after)  after->Scale(1./nEvents);
1159 #endif
1160
1161   TH2D* summed = static_cast<TH2D*>(l->FindObject("summed"));
1162   if (summed) summed->Scale(1./nEvents);
1163   fSum = summed;
1164
1165   TH2D* summedESD = static_cast<TH2D*>(l->FindObject("summedESD"));
1166   if (summedESD) summedESD->Scale(1./nEvents);
1167   fSumESD = summedESD;
1168
1169   TH1D* consecutive = static_cast<TH1D*>(l->FindObject("nConsecutive"));
1170   if (consecutive) consecutive->Scale(1./nEvents);
1171   fNConsecutive= consecutive;
1172 }
1173
1174 //____________________________________________________________________
1175 void
1176 AliFMDSharingFilter::RingHistos::CreateOutputObjects(TList* dir)
1177 {
1178   //
1179   // Make output
1180   //
1181   // Parameters:
1182   //    dir where to store
1183   //
1184   TList* d = DefineOutputList(dir);
1185
1186   d->Add(fBefore);
1187   d->Add(fAfter);
1188   d->Add(fSingle);
1189   d->Add(fDouble);
1190   d->Add(fTriple);
1191   d->Add(fSinglePerStrip);
1192   // d->Add(fDistanceBefore);
1193   // d->Add(fDistanceAfter);
1194   d->Add(fBeforeAfter);
1195   d->Add(fNeighborsBefore);
1196   d->Add(fNeighborsAfter);
1197   // d->Add(fHits);
1198   d->Add(fSumESD);
1199   d->Add(fSum);
1200   d->Add(fNConsecutive);
1201
1202   // Removed to avoid doubly adding the list which destroys
1203   // the merging
1204   //dir->Add(d);
1205 }
1206
1207 //____________________________________________________________________
1208 //
1209 // EOF
1210 //