]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Important bug fixes for the PWGLF/FORWARD directory
authorcholm <Christian.Holm.Christensen@cern.ch>
Wed, 7 May 2014 11:47:22 +0000 (13:47 +0200)
committercholm <Christian.Holm.Christensen@cern.ch>
Wed, 7 May 2014 11:47:22 +0000 (13:47 +0200)
modified:   ../../../OADB/PWGLF/FORWARD/CORRECTIONS/data/fmd_corrections.root
modified:   ../../../OADB/PWGLF/FORWARD/CORRECTIONS/data/spd_corrections.root

Updated corrections and include new noise-gain correction for all possible runs

modified:   ../../CMakelibPWGLFforward2.pkg
modified:   ../../PWGLFforward2LinkDef.h

New classes added

modified:   AddTaskFMDELoss.C

Minor fixes for new AliFMDESDFixer

modified:   AliAODForwardMult.cxx

Allow trigger mask "CENT"

modified:   AliBaseESDTask.cxx
modified:   AliBaseESDTask.h

New method PreCorrections - meant to do stuff when we know the run number etc.,
but before reading the corrections

new file:   AliFMDCorrNoiseGain.h

A new type of correction.  This is used to add in the additional noise that
the reconstruction failed to add because of a bug in the AliFMDRawReader.  With
this, we can process old data with out fear of corruption.

modified:   AliFMDDensityCalculator.cxx
modified:   AliFMDDensityCalculator.h

Moved re-calculation of eta to new AliFMDESDFixer class

new file:   AliFMDESDFixer.cxx
new file:   AliFMDESDFixer.h

A class that fixes up the FMD ESD object for various flaws, including
bad values for 0's, missing noise from reconstruction, and possibly
re-calculation of eta.  This class takes an AliESDFMD object and
modifies the content in-place.  Unfortunately we cannot re-calculate
the phis in this, because we never store the phis in the ESD object.

modified:   AliFMDEnergyFitter.cxx
modified:   AliFMDEnergyFitter.h

Minor stuff

modified:   AliFMDEnergyFitterTask.cxx
modified:   AliFMDEnergyFitterTask.h

Use new class AliESDFMDFixer to fix up ESDs before processing

modified:   AliFMDEventInspector.cxx
modified:   AliFMDEventInspector.h
modified:   AliFMDMCEventInspector.cxx

Added new method for 2013 pA data.  Updated the method for
selecting which vertex method to use.

modified:   AliFMDSharingFilter.cxx
modified:   AliFMDSharingFilter.h

Moved ESD corrections into new class AliFMDESDFixer.

modified:   AliForwardCorrectionManager.cxx
modified:   AliForwardCorrectionManager.h

Possiblity for new correction object AliFMDCorrNoiseGain.

modified:   AliForwardMCMultiplicityTask.cxx
modified:   AliForwardMCMultiplicityTask.h
modified:   AliForwardMultiplicityBase.cxx
modified:   AliForwardMultiplicityBase.h
modified:   AliForwardMultiplicityTask.cxx
modified:   AliForwardMultiplicityTask.h

Updates for new correction ALiFMDCorrNoiseGain

modified:   DrawdNdeta.C

Minor things

modified:   ForwardAODConfig.C

Updates for changes in the code

new file:   corrs/Browse.C
modified:   corrs/ForwardOADBGui.C

New file and fixes for browsing corrections database

modified:   corrs/CorrDrawer.C

Updates for new noise-gain correction

modified:   corrs/RunTestF.C

Run test of Landau*Gaus shift of Delta_p

new file:   corrs/Upload.C

Generic script to upload corrections

modified:   qa/QAPlotter.C
modified:   qa/QATrender.C
modified:   qa/RunQAMT.sh
modified:   qa/runQA.sh

Updates for new centralized QA framework

modified:   scan/Trend.C
modified:   trains/MakeFMDELossTrain.C

Minor things.

43 files changed:
OADB/PWGLF/FORWARD/CORRECTIONS/data/fmd_corrections.root
OADB/PWGLF/FORWARD/CORRECTIONS/data/spd_corrections.root
PWGLF/CMakelibPWGLFforward2.pkg
PWGLF/FORWARD/analysis2/AddTaskFMDELoss.C
PWGLF/FORWARD/analysis2/AliAODForwardMult.cxx
PWGLF/FORWARD/analysis2/AliBaseESDTask.cxx
PWGLF/FORWARD/analysis2/AliBaseESDTask.h
PWGLF/FORWARD/analysis2/AliFMDCorrNoiseGain.h [new file with mode: 0644]
PWGLF/FORWARD/analysis2/AliFMDDensityCalculator.cxx
PWGLF/FORWARD/analysis2/AliFMDDensityCalculator.h
PWGLF/FORWARD/analysis2/AliFMDESDFixer.cxx [new file with mode: 0644]
PWGLF/FORWARD/analysis2/AliFMDESDFixer.h [new file with mode: 0644]
PWGLF/FORWARD/analysis2/AliFMDEnergyFitter.cxx
PWGLF/FORWARD/analysis2/AliFMDEnergyFitter.h
PWGLF/FORWARD/analysis2/AliFMDEnergyFitterTask.cxx
PWGLF/FORWARD/analysis2/AliFMDEnergyFitterTask.h
PWGLF/FORWARD/analysis2/AliFMDEventInspector.cxx
PWGLF/FORWARD/analysis2/AliFMDEventInspector.h
PWGLF/FORWARD/analysis2/AliFMDMCEventInspector.cxx
PWGLF/FORWARD/analysis2/AliFMDSharingFilter.cxx
PWGLF/FORWARD/analysis2/AliFMDSharingFilter.h
PWGLF/FORWARD/analysis2/AliForwardCorrectionManager.cxx
PWGLF/FORWARD/analysis2/AliForwardCorrectionManager.h
PWGLF/FORWARD/analysis2/AliForwardMCMultiplicityTask.cxx
PWGLF/FORWARD/analysis2/AliForwardMCMultiplicityTask.h
PWGLF/FORWARD/analysis2/AliForwardMultiplicityBase.cxx
PWGLF/FORWARD/analysis2/AliForwardMultiplicityBase.h
PWGLF/FORWARD/analysis2/AliForwardMultiplicityTask.cxx
PWGLF/FORWARD/analysis2/AliForwardMultiplicityTask.h
PWGLF/FORWARD/analysis2/DrawdNdeta.C
PWGLF/FORWARD/analysis2/ForwardAODConfig.C
PWGLF/FORWARD/analysis2/corrs/Browse.C [new file with mode: 0644]
PWGLF/FORWARD/analysis2/corrs/CorrDrawer.C
PWGLF/FORWARD/analysis2/corrs/ForwardOADBGui.C
PWGLF/FORWARD/analysis2/corrs/RunTestF.C
PWGLF/FORWARD/analysis2/corrs/Upload.C [new file with mode: 0644]
PWGLF/FORWARD/analysis2/qa/QAPlotter.C
PWGLF/FORWARD/analysis2/qa/QATrender.C
PWGLF/FORWARD/analysis2/qa/RunQAMT.sh
PWGLF/FORWARD/analysis2/qa/runQA.sh
PWGLF/FORWARD/analysis2/scan/Trend.C
PWGLF/FORWARD/analysis2/trains/MakeFMDELossTrain.C
PWGLF/PWGLFforward2LinkDef.h

index 40a2636b0a0e2f392d0205810ae0f46be33421c2..061c6221cdd62c9bc9b536e8ff27f42b88bfc31e 100644 (file)
Binary files a/OADB/PWGLF/FORWARD/CORRECTIONS/data/fmd_corrections.root and b/OADB/PWGLF/FORWARD/CORRECTIONS/data/fmd_corrections.root differ
index 8657b9b4f0ccb8b46db81d4f57369b967f9ae55a..bb7cda79b3aac104d806647893149bcc4edefca6 100644 (file)
Binary files a/OADB/PWGLF/FORWARD/CORRECTIONS/data/spd_corrections.root and b/OADB/PWGLF/FORWARD/CORRECTIONS/data/spd_corrections.root differ
index 30e848900eadbbfad9983842b83efc95e86408ea..06c1d9f73327610c34e79f4297755e9e3dc9ea01 100644 (file)
@@ -86,6 +86,7 @@ set ( SRCS
   FORWARD/analysis2/AliFMDEventInspector.cxx
   FORWARD/analysis2/AliFMDEventPlaneFinder.cxx
   FORWARD/analysis2/AliFMDHistCollector.cxx
+  FORWARD/analysis2/AliFMDESDFixer.cxx
   FORWARD/analysis2/AliFMDSharingFilter.cxx
   # FMD MC algorithms
   FORWARD/analysis2/AliFMDMCCorrector.cxx
@@ -119,6 +120,7 @@ set ( HDRS ${HDRS}
   FORWARD/analysis2/AliFMDStripIndex.h 
   FORWARD/analysis2/AliLandauGaus.h 
   FORWARD/analysis2/AliLandauGausFitter.h 
+  FORWARD/analysis2/AliFMDCorrNoiseGain.h
   )
 
 set ( EINCLUDE  
@@ -145,6 +147,7 @@ set ( EXPORT FORWARD/analysis2/AliAODForwardMult.h
             FORWARD/analysis2/AliFMDCorrAcceptance.h
             FORWARD/analysis2/AliFMDCorrSecondaryMap.h
             FORWARD/analysis2/AliFMDCorrELossFit.h
+            FORWARD/analysis2/AliFMDCorrNoiseGain.h
             FORWARD/analysis2/AliBaseESDTask.h
             FORWARD/analysis2/AliBaseAODTask.h
             FORWARD/analysis2/AliFMDEnergyFitter.h
index c49763d5ed6173ac628726c2eab69e7e9d6ce72d..84cd497c0eb836fd05b1781df019737f7e4f8f35 100644 (file)
@@ -34,7 +34,8 @@ AddTaskFMDELoss(Bool_t        mc,
                Bool_t        useCent,
                Bool_t        onlyMB=false, 
                Int_t         debug=0,
-               const Char_t* residuals="")
+               const Char_t* residuals="",
+               const Char_t* corrs="")
 {
   // --- Load libraries ----------------------------------------------
   gROOT->LoadClass("AliAODForwardMult", "libPWGLFforward2");
@@ -46,14 +47,30 @@ AddTaskFMDELoss(Bool_t        mc,
     return NULL;
   }   
   
+  // --- Set alternative corrections path ----------------------------
+  AliForwardCorrectionManager& cm = AliForwardCorrectionManager::Instance();
+  if (corrs && corrs[0] != '\0') cm.SetPrefix(corrs);
+
   // --- Make the task and add it to the manager ---------------------
   AliFMDEnergyFitterTask* task = new AliFMDEnergyFitterTask("ForwardELoss");
   // --- Set parameters on the algorithms ----------------------------
+
+  // --- Event inspector ---------------------------------------------
   // Set the number of SPD tracklets for which we consider the event a
   // low flux event
   task->GetEventInspector().SetLowFluxCut(1000); 
   // Set the maximum error on v_z [cm]
   task->GetEventInspector().SetMaxVzErr(0.2);
+
+  // --- ESD Fixer ---------------------------------------------------
+  // For MC input we explicitly disable the noise correction 
+  if (mc) task->GetESDFixer().SetRecoNoiseFactor(4);
+  // IF the noise correction is bigger than this, flag strip as dead 
+  task->GetESDFixer().SetMaxNoiseCorrection(0.05);
+  // Dead region in FMD2i
+  task->GetESDFixer().AddDeadRegion(2, 'I', 16, 17, 256, 511);  
+
+  // --- Energy loss fitter ------------------------------------------
   // Set the eta axis to use - note, this overrides whatever is used
   // by the rest of the algorithms - but only for the energy fitter
   // algorithm. 
@@ -78,12 +95,7 @@ AddTaskFMDELoss(Bool_t        mc,
   // Set the minimum number of entries in the distribution before
   // trying to fit to the data - 10K seems the absolute minimum
   task->GetEnergyFitter().SetMinEntries(10000);
-  // If set, only collect statistics for MB.  This is to prevent a
-  // bias when looping over data where the MB trigger is downscaled.
-  task->SetOnlyMB(onlyMB);
-  // Debug
-  task->SetDebug(debug);
-
+  // Check if we're to store the residuals 
   TString resi(residuals);
   resi.ToUpper();
   AliFMDEnergyFitter::EResidualMethod rm = AliFMDEnergyFitter::kNoResiduals;
@@ -97,6 +109,15 @@ AddTaskFMDELoss(Bool_t        mc,
   }
   task->GetEnergyFitter().SetStoreResiduals(rm);
 
+
+  // --- General -----------------------------------------------------
+  // If set, only collect statistics for MB.  This is to prevent a
+  // bias when looping over data where the MB trigger is downscaled.
+  task->SetOnlyMB(onlyMB);
+  // Debug
+  task->SetDebug(debug);
+
+
   // --- Set limits on fits the energy -------------------------------
   // DO NOT CHANGE THESE UNLESS YOU KNOW WHAT YOU ARE DOING
   // Maximum relative error on parameters 
index f177915ba5fed584ee12c20f21d0feafc5dcf339..896b41430b3b00fecea26821c1b184ebd50aad05 100644 (file)
@@ -319,6 +319,8 @@ AliAODForwardMult::MakeTriggerMask(const char* what)
     else if (s.CompareTo("SAT")        == 0) trgMask |= kSatellite;
     else if (s.CompareTo("E")          == 0) trgMask |= kE;
     else if (s.CompareTo("NCLUSTER>0") == 0) trgMask |= kNClusterGt0;
+    else if (s.CompareTo("CENT")       == 0) trgMask |= kInel;
+    // trgMask &= ~(kInel|kInelGt0|kNSD|kV0AND|kMCNSD);
     else 
       AliWarningGeneral("MakeTriggerMask", 
                        Form("Unknown trigger %s", s.Data()));
index bc96a5addf2fe7cb523014c09120e922c6b7948f..b996bc0e4276c17689749debaa619afc17675b68 100644 (file)
@@ -376,7 +376,9 @@ AliBaseESDTask::GetESDEvent()
           esd->GetCurrentL3(),
           esd->GetMagneticField(),
           esd->GetRunNumber());
-  
+
+  PreCorrections(esd);
+
   fFirstEvent = false;
   
   const   TAxis* pe = 0;
index 96d34501bb68d0d3e3cf932a2ad4e5d508a75e55..b97ae71ddfa739f7413ea3d99e6451fdbebd4363 100644 (file)
@@ -120,6 +120,14 @@ public:
    * @return true on success. 
    */
   virtual Bool_t Book() = 0;
+  /** 
+   * Called on first event _before_ reading corrections.  Here, the
+   * user class can do additional checking to see if the some (more or
+   * less) corrections are needed.
+   * 
+   * @param esd Event 
+   */
+  virtual void PreCorrections(const AliESDEvent* esd);
   /** 
    * Called after reading in the first event. Here we can setup stuff
    * depending on the conditions we're running under.
@@ -354,6 +362,7 @@ private:
 
   ClassDef(AliBaseESDTask,1);
 };
+inline void AliBaseESDTask::PreCorrections(const AliESDEvent*) {}
 #endif
 // Local Variables:
 //   mode: C++
diff --git a/PWGLF/FORWARD/analysis2/AliFMDCorrNoiseGain.h b/PWGLF/FORWARD/analysis2/AliFMDCorrNoiseGain.h
new file mode 100644 (file)
index 0000000..24e3ae4
--- /dev/null
@@ -0,0 +1,119 @@
+#ifndef ALIFMDCORRNOISEGAIN_H
+#define ALIFMDCORRNOISEGAIN_H
+#include <AliFMDFloatMap.h>
+
+/**
+ * Get the noise calibration.  That is, the ratio 
+ *
+ * @f[ 
+ *   \frac{\sigma_{i}}{g_{i}k}
+ * @f] 
+ *
+ * where @f$ k@f$ is a constant determined by the electronics of
+ * units DAC/MIP, and @f$ \sigma_i, g_i@f$ are the noise and gain of
+ * the @f$ i @f$ strip respectively.
+ *
+ * This correction is needed because some of the reconstructed data
+ * (what which have an AliESDFMD class version less than or equal to
+ * 3) used the wrong zero-suppression factor.  The zero suppression
+ * factor used by the on-line electronics was 4, but due to a coding
+ * error in the AliFMDRawReader a zero suppression factor of 1 was
+ * assumed during the reconstruction.  This shifts the zero of the
+ * energy loss distribution artificially towards the left (lover
+ * valued signals).
+ *
+ * So let's assume the real zero-suppression factor is @f$ f@f$ while
+ * the zero suppression factor @f$ f'@f$ assumed in the reconstruction
+ * was (wrongly) lower.  The number of ADC counts @f$ c_i'@f$ used in
+ * the reconstruction can be calculated from the reconstructed signal
+ * @f$ m_i'@f$ by
+ *
+ * @f[
+ *    c_i' = m_i \times g_i \times k / \cos\theta_i
+ * @f] 
+ *
+ * where @f$\theta_i@f$ is the incident angle of the @f$ i@f$ strip. 
+ * 
+ * This number of counts used the wrong noise factor @f$ f'@f$ so to
+ * correct to the on-line value, we need to do
+ *
+ * @f[ 
+ *   c_i = c_i' - \floor{f'\times n_i} + \floor{f\times n_i}
+ * @f] 
+ * 
+ * which gives the correct number of ADC counts over the pedestal. To
+ * convert back to the scaled energy loss signal we then need to
+ * calculate (noting that @f$ f,f'@f$ are integers)
+ *
+ * @f{eqnarray}
+ *    m_i &=& \frac{c_i \times \cos\theta_i}{g_i \times k}\\ 
+ *    &=& \left(c_i' - \lfloor f'\times n_i\rfloor + 
+ *           \lfloor f\times n_i\rfloor\right)\frac{\cos\theta}{g_i \times k}\\
+ *    &=& \left(\frac{m_i'\times g_i\times k}{\cos\theta} -
+ *            \lfloor f'\times n_i\rfloor + \lfloor f\times n_i\rfloor\right)
+ *         \frac{\cos\theta}{g_i \times k}\\
+ *    &=& m_i' + \frac{1}{g_i \times k}
+ *         \left(\lfloor f\times n_i}\rfloor-
+ *             \lfloor f'\times n_i\rfloor\right)\cos\theta\\
+ *    &=& m_i' + \frac{\lfloor n_i\rfloor}{g_i \times k}
+ *        \left(f-f'\right)\cos\theta
+ * @f{eqnarray}
+ * 
+ */
+class AliFMDCorrNoiseGain : public TObject 
+{
+public:
+  /**
+   * Default constructor 
+   */
+  AliFMDCorrNoiseGain() : fValues(0) { fValues.Reset(-1); }
+  /** 
+   * Constructor from a float map 
+   * 
+   * @param map Construct from this map 
+   */
+  AliFMDCorrNoiseGain(const AliFMDFloatMap& map) : fValues(map) {}
+  /** 
+   * Get the noise value for a particular strip 
+   * 
+   * @param d  Detector
+   * @param r  Ring 
+   * @param s  Sector 
+   * @param t  Strip 
+   * 
+   * @return Noise value for strip 
+   */
+  Float_t Get(UShort_t d, Char_t r, UShort_t s, UShort_t t) const 
+  { 
+    return fValues(d, r, s, t);
+  }
+  /** 
+   * Set the value for a strip. 
+   * 
+   * @param d Detector 
+   * @param r Ring 
+   * @param s Sector
+   * @param t Strip
+   * @param x Value 
+   */
+  void Set(UShort_t d, Char_t r, UShort_t s, UShort_t t, Float_t x) 
+  { 
+    fValues(d, r, s, t) = x;
+  }
+  /** 
+   * Get a reference to the noise map 
+   * 
+   * @return Noise map 
+   */
+  const AliFMDFloatMap& Values() { return fValues; }
+protected: 
+  AliFMDFloatMap fValues; // The noise-gain map 
+  ClassDef(AliFMDCorrNoiseGain,1); // Clone of AliFMDCalibPedestal
+};
+
+#endif
+// Local Variables:
+//   mode: C++
+// End:
+
index ba597be5a777758953d541089d12946ddaadfeb1..2ebf01b304cda27e31b41e7940260a8705838de3 100644 (file)
@@ -53,7 +53,6 @@ AliFMDDensityCalculator::AliFMDDensityCalculator()
     fPhiLumping(4),    
     fDebug(0),
     fCuts(),
-    fRecalculateEta(false),
     fRecalculatePhi(false),
     fMinQuality(10),
     fCache(),
@@ -91,7 +90,6 @@ AliFMDDensityCalculator::AliFMDDensityCalculator(const char* title)
     fPhiLumping(4),
     fDebug(0),
     fCuts(),
-    fRecalculateEta(false),
     fRecalculatePhi(false),
     fMinQuality(10),
     fCache(),
@@ -165,7 +163,6 @@ AliFMDDensityCalculator::AliFMDDensityCalculator(const
     fPhiLumping(o.fPhiLumping),
     fDebug(o.fDebug),
     fCuts(o.fCuts),
-    fRecalculateEta(o.fRecalculateEta),
     fRecalculatePhi(o.fRecalculatePhi),
     fMinQuality(o.fMinQuality),
     fCache(o.fCache),
@@ -229,7 +226,6 @@ AliFMDDensityCalculator::operator=(const AliFMDDensityCalculator& o)
   fEtaLumping         = o.fEtaLumping;
   fPhiLumping         = o.fPhiLumping;
   fCuts               = o.fCuts;
-  fRecalculateEta     = o.fRecalculateEta;
   fRecalculatePhi     = o.fRecalculatePhi;
   fMinQuality         = o.fMinQuality;
   fCache              = o.fCache;
@@ -347,7 +343,17 @@ AliFMDDensityCalculator::GetMultCut(UShort_t d, Char_t r, Double_t eta,
   return Rng2Cut(d, r, ieta, fLowCuts);
   // return fCuts.GetMultCut(d,r,eta,errors);
 }
-  
+
+#ifndef NO_TIMING
+# define START_TIMER(T) if (fDoTiming) T.Start(true)
+# define GET_TIMER(T,V) if (fDoTiming) V = T.CpuTime()
+# define ADD_TIMER(T,V) if (fDoTiming) V += T.CpuTime()
+#else
+# define START_TIMER(T) do {} while (false)
+# define GET_TIMER(T,V) do {} while (false)
+# define ADD_TIMER(T,V) do {} while (false)
+#endif
+
 //____________________________________________________________________
 Bool_t
 AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
@@ -380,14 +386,13 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
   //  Copy to cache       : fraction of sum  3.9%   of total  2.2%
   //  Poisson calculation : fraction of sum 18.7%   of total 10.6%
   //  Diagnostics         : fraction of sum  3.7%   of total  2.1%
-  Double_t reEtaTime  = 0;  
   Double_t nPartTime  = 0;
   Double_t corrTime   = 0;
   Double_t rePhiTime  = 0;
   Double_t copyTime   = 0;
   Double_t poissonTime= 0;
   Double_t diagTime   = 0;
-  if (fDoTiming) totalT.Start(true);
+  START_TIMER(totalT);
   
   Double_t etaCache[20*512]; // Same number of strips per ring 
   Double_t phiCache[20*512]; // whether it is inner our outer. 
@@ -433,10 +438,7 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
          Double_t oldPhi = phi;
 
          // --- Re-calculate eta - needed for satelittes ------------
-         if (fDoTiming) timer.Start(true);
-         if (eta == AliESDFMD::kInvalidEta || fRecalculateEta) 
-           eta = AliForwardUtil::GetEtaFromStrip(d,r,s,t,ip.Z());
-         if (fDoTiming) reEtaTime += timer.CpuTime();
+         START_TIMER(timer);
          etaCache[s*nt+t] = eta;
 
          // --- Check this strip ------------------------------------
@@ -455,13 +457,13 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
          rh->fGood->Fill(eta);
 
          // --- If we asked to re-calculate phi for (x,y) IP --------
-         if (fDoTiming) timer.Start(true);
+         START_TIMER(timer);
          if (fRecalculatePhi) {
            oldPhi = phi;
            phi = AliForwardUtil::GetPhiFromStrip(r, t, phi, ip.X(), ip.Y());
          }
          phiCache[s*nt+t] = phi;
-         if (fDoTiming) rePhiTime += timer.CpuTime();
+         ADD_TIMER(timer,rePhiTime);
 
          // --- Apply phi corner correction to eloss ----------------
          if (fUsePhiAcceptance == kPhiCorrectELoss) 
@@ -474,22 +476,22 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
                           d, r, s, t, eta);
 
          // --- Now caluculate Nch for this strip using fits --------
-         if (fDoTiming) timer.Start(true);
+         START_TIMER(timer);
          Double_t n   = 0;
          if (cut > 0 && mult > cut) n = NParticles(mult,d,r,eta,lowFlux);
          rh->fELoss->Fill(mult);
          // rh->fEvsN->Fill(mult,n);
          // rh->fEtaVsN->Fill(eta, n);
-         if (fDoTiming) nPartTime += timer.CpuTime();
+         ADD_TIMER(timer,nPartTime);
          
          // --- Calculate correction if needed ----------------------
-         if (fDoTiming) timer.Start(true);
+         START_TIMER(timer);
          // Temporary stuff - remove Correction call 
          Double_t c = 1;
          if (fUsePhiAcceptance == kPhiCorrectNch) 
            c = AcceptanceCorrection(r,t);
          // Double_t c = Correction(d,r,t,eta,lowFlux);
-         if (fDoTiming) corrTime += timer.CpuTime();
+         ADD_TIMER(timer,corrTime);
          fCorrections->Fill(c);
          if (c > 0) n /= c;
          // rh->fEvsM->Fill(mult,n);
@@ -518,7 +520,7 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
       rh->fGood->Divide(rh->fGood, rh->fTotal, 1, 1, "B");
 
       // --- Make a copy and reset as needed -------------------------
-      if (fDoTiming) timer.Start(true);
+      START_TIMER(timer);
       TH2D* hclone = fCache.Get(d,r);
       // hclone->Reset();
       // TH2D* hclone = static_cast<TH2D*>(h->Clone("hclone"));
@@ -533,10 +535,10 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
        // hclone->Add(h); 
        h->Reset(); 
       }
-      if (fDoTiming) copyTime += timer.CpuTime();
+      ADD_TIMER(timer,copyTime);
       
       // --- Store Poisson result ------------------------------------
-      if (fDoTiming) timer.Start(true);
+      START_TIMER(timer);
       TH2D* poisson = rh->fPoisson.Result();
       for (Int_t t=0; t < poisson->GetNbinsX(); t++) { 
        for (Int_t s=0; s < poisson->GetNbinsY(); s++) { 
@@ -548,10 +550,6 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
          Double_t  eta  = etaCache[s*nt+t]; 
          // Double_t  phi  = fmd.Phi(d,r,s,t) * TMath::DegToRad();
          // Double_t  eta  = fmd.Eta(d,r,s,t);
-         // if (fRecalculateEta)  
-         //   eta = AliForwardUtil::GetEtaFromStrip(d,r,s,t,ip.Z());
-         // if (fRecalculatePhi)
-         //   phi = AliForwardUtil::GetPhiFromStrip(r, t, phi, ip.X(), ip.Y());
          if (fUsePoisson) {
            h->Fill(eta,phi,poissonV);
            rh->fDensity->Fill(eta, phi, poissonV);
@@ -560,10 +558,10 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
            hclone->Fill(eta,phi,poissonV);
        }
       }
-      if (fDoTiming) poissonTime += timer.CpuTime();
+      ADD_TIMER(timer,poissonTime);
       
       // --- Make diagnostics - eloss vs poisson ---------------------
-      if (fDoTiming) timer.Start(true);
+      START_TIMER(timer);
       Int_t nY = h->GetNbinsY();
       Int_t nIn  = 0; // Count non-outliers
       Int_t nOut = 0; // Count outliers
@@ -611,14 +609,14 @@ AliFMDDensityCalculator::Calculate(const AliESDFMD&        fmd,
       rh->fOutliers->Fill(outRatio);
       if (outRatio < fMaxOutliers) rh->fPoisson.FillDiagnostics();
       else                         h->SetBit(AliForwardUtil::kSkipRing);
-      if (fDoTiming) diagTime += timer.CpuTime();
+      ADD_TIMER(timer,diagTime);
       // delete hclone;
       
     } // for q
   } // for d
 
   if (fDoTiming) {
-    fHTiming->Fill(1,reEtaTime);
+    // fHTiming->Fill(1,reEtaTime);
     fHTiming->Fill(2,nPartTime);
     fHTiming->Fill(3,corrTime);
     fHTiming->Fill(4,rePhiTime);
@@ -1093,7 +1091,6 @@ AliFMDDensityCalculator::CreateOutputObjects(TList* dir)
   d->Add(AliForwardUtil::MakeParameter("phiAcceptance",fUsePhiAcceptance));
   d->Add(AliForwardUtil::MakeParameter("etaLumping",   fEtaLumping));
   d->Add(AliForwardUtil::MakeParameter("phiLumping",   fPhiLumping));
-  d->Add(AliForwardUtil::MakeParameter("recalcEta",    fRecalculateEta));
   d->Add(AliForwardUtil::MakeParameter("recalcPhi",    fRecalculatePhi));
   d->Add(AliForwardUtil::MakeParameter("maxOutliers",  fMaxOutliers));
   d->Add(AliForwardUtil::MakeParameter("outlierCut",   fOutlierCut));
@@ -1165,12 +1162,11 @@ AliFMDDensityCalculator::Print(Option_t* option) const
 
   PFV("Max(particles)",                fMaxParticles );
   PFV("Min(quality)",           fMinQuality );
-  PFV("Poisson method",                fUsePoisson );
+  PFB("Poisson method",                fUsePoisson );
   PFV("Eta lumping",           fEtaLumping );
   PFV("Phi lumping",           fPhiLumping );
-  PFV("Recalculate eta",       fRecalculateEta );
-  PFV("Recalculate phi",       fRecalculatePhi );
-  PFV("Use phi acceptance",     phiM);
+  PFB("Recalculate phi",       fRecalculatePhi );
+  PFB("Use phi acceptance",     phiM);
   PFV("Lower cut", "");
   fCuts.Print();
 
index 454d0e0172e1b2530297e8176a4cef308a3b8d15..e9c6ada00d988e09962acb7d92a953f9de17f006 100644 (file)
@@ -157,13 +157,6 @@ public:
    * number of particles that has hit within a region.
    */
   void SetUsePoisson(Bool_t u) { fUsePoisson = u; }
-  /** 
-   * In case of a displaced vertices recalculate eta and angle correction
-   * 
-   * @param use recalculate or not
-   * 
-   */
-  void SetRecalculateEta(Bool_t use) { fRecalculateEta = use; }
   /** 
    * In case of a displaced vertices recalculate eta and angle correction
    * 
@@ -494,7 +487,6 @@ protected:
   Int_t    fPhiLumping;    //  How to lump phi bins for Poisson 
   Int_t    fDebug;         //  Debug level 
   AliFMDMultCuts fCuts;    // Cuts
-  Bool_t   fRecalculateEta;  // Whether to recalc eta and angle correction (disp vtx)
   Bool_t   fRecalculatePhi;  // Whether to correct for (X,Y) offset
   UShort_t fMinQuality;      // Least quality for fits
   AliForwardUtil::Histos fCache;
diff --git a/PWGLF/FORWARD/analysis2/AliFMDESDFixer.cxx b/PWGLF/FORWARD/analysis2/AliFMDESDFixer.cxx
new file mode 100644 (file)
index 0000000..85bb66f
--- /dev/null
@@ -0,0 +1,349 @@
+#include "AliFMDESDFixer.h"
+#include "AliESDFMD.h"
+#include "AliFMDStripIndex.h"
+#include "AliForwardUtil.h"
+#include "AliForwardCorrectionManager.h"
+#include "AliFMDCorrNoiseGain.h"
+#include <TH1.h>
+#include <TList.h>
+#include <TObjArray.h>
+#include <TROOT.h>
+#include <TMath.h>
+#include <iostream>
+#include <iomanip>
+
+
+//____________________________________________________________________
+AliFMDESDFixer::AliFMDESDFixer() 
+  : TObject(), 
+    fRecoFactor(1),
+    fMaxNoiseCorr(0.05),
+    fRecalculateEta(true),
+    fXtraDead(0),
+    fHasXtraDead(false),
+    fInvalidIsEmpty(false),
+    fNoiseChange(0),
+    fEtaChange(0),
+    fDeadChange(0)
+{}
+
+//____________________________________________________________________
+AliFMDESDFixer::AliFMDESDFixer(const char*) 
+  : TObject(), 
+    fRecoFactor(1),
+    fMaxNoiseCorr(0.05),
+    fRecalculateEta(false),
+    fXtraDead(AliFMDStripIndex::Pack(3,'O',19,511)+1),
+    fHasXtraDead(false),
+    fInvalidIsEmpty(false),
+    fNoiseChange(0),
+    fEtaChange(0),
+    fDeadChange(0)
+{}
+
+
+//____________________________________________________________________
+AliFMDESDFixer::AliFMDESDFixer(const AliFMDESDFixer& o) 
+  : TObject(), 
+    fRecoFactor(o.fRecoFactor),
+    fMaxNoiseCorr(o.fMaxNoiseCorr),
+    fRecalculateEta(o.fRecalculateEta),
+    fXtraDead(o.fXtraDead),
+    fHasXtraDead(o.fHasXtraDead),
+    fInvalidIsEmpty(o.fInvalidIsEmpty),
+    fNoiseChange(0),
+    fEtaChange(0),
+    fDeadChange(0)
+{}
+
+//____________________________________________________________________
+AliFMDESDFixer&
+AliFMDESDFixer::operator=(const AliFMDESDFixer& o) 
+{ 
+  if (&o == this) return *this;
+
+  fRecoFactor = o.fRecoFactor;
+  fMaxNoiseCorr = o.fMaxNoiseCorr;
+  fRecalculateEta = o.fRecalculateEta;
+  fXtraDead = o.fXtraDead;
+  fHasXtraDead = o.fHasXtraDead;
+  fInvalidIsEmpty = o.fInvalidIsEmpty;
+  fNoiseChange = 0;
+  fEtaChange  = 0;
+  fDeadChange = 0;
+
+  return *this; 
+};
+
+//____________________________________________________________________
+void
+AliFMDESDFixer::CreateOutputObjects(TList* l)
+{
+  TList* d = new TList;
+  d->SetName(GetName());
+  l->Add(d);
+
+  d->Add(AliForwardUtil::MakeParameter("recoFactor", fRecoFactor));
+  d->Add(AliForwardUtil::MakeParameter("recalcEta",  fRecalculateEta));
+  d->Add(AliForwardUtil::MakeParameter("invalidIsEmpty",  fInvalidIsEmpty));
+  
+  fNoiseChange = new TH1D("noiseChange", "#delta#Delta due to noise",30,0,.3);
+  fNoiseChange->SetDirectory(0);
+  fNoiseChange->SetFillColor(kYellow+1);
+  fNoiseChange->SetFillStyle(3001);
+  d->Add(fNoiseChange);
+
+  fEtaChange = new TH1D("etaChange", "#delta#eta",40,-1,1);
+  fEtaChange->SetDirectory(0);
+  fEtaChange->SetFillColor(kCyan+1);
+  fEtaChange->SetFillStyle(3001);
+  d->Add(fEtaChange);
+
+  fDeadChange = new TH1D("deadChange", "#deltaN_{dead}",100,0,51200);
+  fDeadChange->SetDirectory(0);
+  fDeadChange->SetFillColor(kMagenta+1);
+  fDeadChange->SetFillStyle(3001);
+  d->Add(fDeadChange);
+
+  TObjArray* extraDead = new TObjArray;
+  extraDead->SetOwner();
+  extraDead->SetName("extraDead");
+  fXtraDead.Compact();
+  UInt_t firstBit = fXtraDead.FirstSetBit();
+  UInt_t nBits    = fXtraDead.GetNbits();
+  for (UInt_t i = firstBit; i < nBits; i++) {
+    if (!fXtraDead.TestBitNumber(i)) continue;
+    UShort_t dd, s, t;
+    Char_t   r;
+    AliFMDStripIndex::Unpack(i, dd, r, s, t);
+    extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
+                                                     dd, r, s, t), 
+                                                UShort_t(i)));
+  }
+  d->Add(extraDead);
+  fHasXtraDead = nBits > 0;
+
+}
+
+//____________________________________________________________________
+void
+AliFMDESDFixer::AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t)
+{
+  if (d < 1 || d > 3) {
+    Warning("AddDead", "Invalid detector FMD%d", d);
+    return;
+  }
+  Bool_t inner = (r == 'I' || r == 'i');
+  if (d == 1 && !inner) { 
+    Warning("AddDead", "Invalid ring FMD%d%c", d, r);
+    return;
+  }
+  if ((inner && s >= 20) || (!inner && s >= 40)) { 
+    Warning("AddDead", "Invalid sector FMD%d%c[%02d]", d, r, s);
+    return;
+  }
+  if ((inner && t >= 512) || (!inner && t >= 256)) { 
+    Warning("AddDead", "Invalid strip FMD%d%c[%02d,%03d]", d, r, s, t);
+    return;
+  }
+    
+  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
+  // Int_t i  = 0;
+  fXtraDead.SetBitNumber(id, true);
+}
+//____________________________________________________________________
+void
+AliFMDESDFixer::AddDeadRegion(UShort_t d,  Char_t r, 
+                                  UShort_t s1, UShort_t s2, 
+                                  UShort_t t1, UShort_t t2)
+{
+  // Add a dead region spanning from FMD<d><r>[<s1>,<t1>] to 
+  // FMD<d><r>[<s2>,<t2>] (both inclusive)
+  for (Int_t s = s1; s <= s2; s++) 
+    for (Int_t t = t1; t <= t2; t++) 
+      AddDead(d, r, s, t);
+}
+//____________________________________________________________________
+void
+AliFMDESDFixer::AddDead(const Char_t* script)
+{
+  if (!script || script[0] == '\0') return;
+  
+  gROOT->Macro(Form("%s((AliFMDESDFixer*)%p);", script, this));
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDESDFixer::IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const
+{
+  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
+  return fXtraDead.TestBitNumber(id); 
+}
+
+//____________________________________________________________________
+Int_t
+AliFMDESDFixer::FindTargetNoiseFactor(const AliESDFMD& esd, Bool_t check) const
+{
+  if (!IsUseNoiseCorrection()) 
+    // If the reconstruction factor was high (4 or more), do nothing 
+    return 0;
+
+  if (!esd.NeedNoiseFix()) { 
+    // If the bit isn't set, doo nothing
+    return 0;
+  }
+  
+  // Get the target factor - even thought the method below returns a
+  // floating point value, we know that the noise factor is always
+  // integer, so we coerce it to be the same here. 
+  Int_t target = esd.GetNoiseFactor() - fRecoFactor;
+
+  // If the target factor is the same or smaller than the assumed
+  // factor, we have nothing to do here, and we return immediately
+  if (target <= 0) return 0;
+
+  // Get the scaled noise from the correction mananger 
+  if (check && !AliForwardCorrectionManager::Instance().GetNoiseGain()) 
+    return 0;
+
+  return target;
+}
+
+//____________________________________________________________________
+#define ETA2COS(ETA)                                           \
+  TMath::Cos(2*TMath::ATan(TMath::Exp(-TMath::Abs(ETA))))
+
+//____________________________________________________________________
+void
+AliFMDESDFixer::Fix(AliESDFMD& esd, Double_t zvtx)
+{
+
+  const AliFMDCorrNoiseGain* ng  = 0;
+  Int_t tgtFactor = FindTargetNoiseFactor(esd, false);
+  if (tgtFactor > 0) 
+    ng  = AliForwardCorrectionManager::Instance().GetNoiseGain();
+
+  if (!ng && !fHasXtraDead && !fRecalculateEta && !fInvalidIsEmpty) 
+    // We have nothing to do!
+    return;
+
+  UShort_t nDead = 0;
+  for (UShort_t d = 1; d <= 3; d++) { 
+    UShort_t nQ = d == 1 ? 1 : 2;
+
+    for (UShort_t q = 0; q < nQ; q++) { 
+      Char_t   r  = (q == 0 ? 'I' : 'O');
+      UShort_t nS = (q == 0 ?  20 :  40);
+      UShort_t nT = (q == 0 ? 512 : 256);
+      
+      for (UShort_t s = 0; s < nS; s++) { 
+       for (UShort_t t = 0; t < nT; t++) { 
+         Double_t mult     = esd.Multiplicity(d,r,s,t);
+         Double_t eta      = esd.Eta(d,r,s,t);
+         Double_t cosTheta = 0;
+
+         if (CheckDead(d,r,s,t,mult)) nDead++;
+         
+         // Possibly re-calculate eta 
+         if (fRecalculateEta) RecalculateEta(d,r,s,t,zvtx,eta,mult,cosTheta);
+
+         // Possibly correct for poor treatment of ZS in reconstruction. 
+         if (ng && mult != AliESDFMD::kInvalidMult) {
+           if (cosTheta <= 0) cosTheta = ETA2COS(eta);
+           if (!NoiseCorrect(tgtFactor,ng->Get(d,r,s,t), cosTheta, mult))
+             nDead++;
+         } 
+
+         // Write out final values to object 
+         if (mult >= AliESDFMD::kInvalidMult) mult = AliESDFMD::kInvalidMult;
+         esd.SetMultiplicity(d,r,s,t,mult);
+         esd.SetEta(d,r,s,t,eta);
+       } // for t
+      } // for s
+    } // for q
+  } // for d
+  fDeadChange->Fill(nDead);
+}
+
+//____________________________________________________________________
+Bool_t
+AliFMDESDFixer::CheckDead(UShort_t d, Char_t r, UShort_t s, UShort_t t,
+                         Double_t& mult)
+{
+  // Correct for zero's being flagged as invalid 
+  if (mult == AliESDFMD::kInvalidMult && fInvalidIsEmpty) mult = 0;
+
+  // Take into account what we're defined as dead 
+  if (IsDead(d,r,s,t)) {
+    mult = AliESDFMD::kInvalidMult;
+    return true;
+  }
+  return false;
+}
+
+//____________________________________________________________________
+void
+AliFMDESDFixer::RecalculateEta(UShort_t d, Char_t r, UShort_t s, UShort_t t,
+                              Double_t zvtx, Double_t& eta, Double_t& mult, 
+                              Double_t& cosTheta)
+{
+  Double_t oldEta = eta;
+  Double_t newEta = AliForwardUtil::GetEtaFromStrip(d,r,s,t, zvtx);
+  eta             = newEta;
+
+  fEtaChange->Fill(newEta-oldEta);
+
+  if (mult == AliESDFMD::kInvalidMult) return;
+
+  Double_t newCos = ETA2COS(newEta);
+  Double_t oldCos = ETA2COS(oldEta);
+  Double_t corr   = newCos / oldCos;
+  cosTheta        = newCos;
+  mult            *= corr;
+}
+//____________________________________________________________________
+Bool_t
+AliFMDESDFixer::NoiseCorrect(Int_t target, Double_t corr, Double_t cosTheta, 
+                            Double_t& mult)
+{
+  if (corr > fMaxNoiseCorr || corr <= 0) { 
+    mult = AliESDFMD::kInvalidMult;
+    return false;
+  }
+  Double_t add = corr * target * cosTheta;
+  fNoiseChange->Fill(add);
+
+  mult += add;
+  return true;
+}
+
+#define PF(N,V,...)                                    \
+  AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
+#define PFB(N,FLAG)                            \
+  do {                                                                 \
+    AliForwardUtil::PrintName(N);                                      \
+    std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
+  } while(false)
+#define PFV(N,VALUE)                                   \
+  do {                                                 \
+    AliForwardUtil::PrintName(N);                      \
+    std::cout << (VALUE) << std::endl; } while(false)
+
+//____________________________________________________________________
+void
+AliFMDESDFixer::Print(Option_t*) const
+{
+  AliForwardUtil::PrintTask(*this);
+  gROOT->IncreaseDirLevel();
+  PFB("Consider invalid null", fInvalidIsEmpty);
+  PFB("Has extra dead", fHasXtraDead);
+  PFV("Reco noise factor", fRecoFactor);
+  PFV("Max noise corr", fMaxNoiseCorr);
+  PFB("Recalc. eta", fRecalculateEta);
+  gROOT->DecreaseDirLevel();
+}
+
+//____________________________________________________________________
+// 
+// EOF
+// 
diff --git a/PWGLF/FORWARD/analysis2/AliFMDESDFixer.h b/PWGLF/FORWARD/analysis2/AliFMDESDFixer.h
new file mode 100644 (file)
index 0000000..b14c11a
--- /dev/null
@@ -0,0 +1,245 @@
+#ifndef ALIFMDESDFIXER_H
+#define ALIFMDESDFIXER_H 
+#include <TObject.h>
+#include <TBits.h>
+class TH1;
+class TList;
+class AliESDFMD;
+
+/**
+ * Class to fix up an ESD object for various small issues. 
+ *
+ * @par Input:
+ *    - AliESDFMD object - from reconstruction 
+ *    - Ip z coordinate - from reconstruction 
+ *    - Reco noise factor - Assumed noise factor used in reconstruction
+ *
+ * @par Output:
+ *    - The same AliESDFMD object but modified 
+ *
+ * @par Corrections used
+ *    - AliFMDCorrNoiseGain if the noise correction is enabled 
+ * 
+ * @par Histograms
+ *    - Change in @f$\Delta@f$ due to noise correction
+ *    - Number of channels declared dead 
+ *    - Change in @f$\eta@f$  
+ *
+ * @ingroup pwglf_forward_algo 
+ * @ingroup pwglf_forward_aod
+ */
+class AliFMDESDFixer : public TObject 
+{
+public:
+  /**
+   * Default CTOR - for ROOT I/O - do not use 
+   */
+  AliFMDESDFixer();
+
+  /** 
+   * User CTOR 
+   *
+   * @param name Dummy argument 
+   */
+  AliFMDESDFixer(const char* name);
+  /** 
+   * Get name of object 
+   * 
+   * @return Always the same 
+   */
+  const char* GetName() const { return "fmdESDFixer"; }
+  /**
+   * Create output objects 
+   */
+  void CreateOutputObjects(TList* l);
+  /** 
+   * Fix the ESD object 
+   * 
+   * @param esd ESD object 
+   * @param ipZ IP Z--coordinate
+   */
+  void Fix(AliESDFMD& esd, Double_t ipZ);
+
+  /** 
+   * @{ 
+   * @name Noise suppression fix 
+   */
+  /** 
+   * Set the factor assumed to be used in the reconstruction.  If this
+   * is set way high (>=4) then this corrector will effectively be
+   * disabled.
+   * 
+   * @param f 
+   */
+  void SetRecoNoiseFactor(Int_t f) { fRecoFactor = f; }
+  void SetMaxNoiseCorrection(Double_t x) { fMaxNoiseCorr = x; }
+  /** 
+   * Get the factor assumed in reconstruction pass
+   * 
+   * @return factor assumed in reconstruction pass
+   */
+  Int_t GetRecoNoiseFactor() const { return fRecoFactor; }
+  /** 
+   * Check if we're using the noise correction. 
+   * 
+   * @return true if fRecoFactor < 4
+   */
+  Bool_t IsUseNoiseCorrection() const { return fRecoFactor < 4; }
+  /** 
+   * Find the target noise factor 
+   * 
+   * @param esd   ESD object. 
+   * @param check If true, also check for correction object 
+   * 
+   * @return Needed noise factor, or 0 or less if no correction is needed
+   */
+  Int_t FindTargetNoiseFactor(const AliESDFMD& esd, Bool_t check=true) const;
+  /* @} */
+
+  /** 
+   * @{ 
+   * @name Recalculations 
+   */
+  /** 
+   * In case of a displaced vertices recalculate @f$\eta@f$ and angle
+   * correction
+   * 
+   * @param use recalculate or not
+   * 
+   */
+  void SetRecalculateEta(Bool_t use) { fRecalculateEta = use; }
+  /* @} */
+
+  /**
+   * @{ 
+   * @name Special treatmeant of invalid signals 
+   */
+  /** 
+   * Set whether to consider invalid multiplicities as null (or empty)
+   * signal. 
+   * 
+   * @param flag If true, count invalids as empty
+   */
+  void SetInvalidIsEmpty(Bool_t flag) { fInvalidIsEmpty = flag; }
+  /** @} */
+
+  /** 
+   * @{ 
+   * @name Dead strip handling 
+   */
+  /** 
+   * Add a dead strip
+   * 
+   * @param d  Detector
+   * @param r  Ring 
+   * @param s  Sector 
+   * @param t  Strip
+   */
+  void AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t);
+  /** 
+   * Add a dead region in a detector ring
+   * 
+   * @param d   Detector
+   * @param r   Ring
+   * @param s1  First sector (inclusive)
+   * @param s2  Last sector (inclusive)
+   * @param t1  First strip (inclusive)
+   * @param t2  Last strip (inclusive)
+   */
+  void AddDeadRegion(UShort_t d, Char_t r, UShort_t s1, UShort_t s2, 
+                    UShort_t t1, UShort_t t2);
+  /** 
+   * Add dead strips from a script.  The script is supposed to accept
+   * a pointer to this object (AliFMDSharingFilter) and then call
+   * AddDead or AddDeadRegion as needed.
+   * 
+   * @code 
+   * void deadstrips(AliFMDSharingFilter* filter)
+   * { 
+   *   filter->AddDead(...);
+   *   // ... and so on 
+   * }
+   * @endcode 
+   *
+   * @param script The script to read dead strips from. 
+   */
+  void AddDead(const Char_t* script);
+  /* @} */
+
+  void Print(Option_t* option="") const;
+protected:
+  AliFMDESDFixer(const AliFMDESDFixer&);
+  AliFMDESDFixer& operator=(const AliFMDESDFixer&);
+  
+  /** 
+   * @{ 
+   * @name Worker functions 
+   */
+  /** 
+   * Check if a strip is marked as dead 
+   * 
+   * @param d   Detector 
+   * @param r   Ring
+   * @param s   Sector 
+   * @param t   Strip 
+   * 
+   * @return true if dead 
+   */
+  virtual Bool_t IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const;
+  /** 
+   * Possibly raise a strip from the dead or kill it 
+   * 
+   * @param d   Detector 
+   * @param r   Ring
+   * @param s   Sector 
+   * @param t   Strip 
+   * @param m   Multiplicity
+   * 
+   * @return true if this was killed
+   */
+  Bool_t CheckDead(UShort_t d, Char_t r, UShort_t s, UShort_t t, Double_t& m);
+  /** 
+   * Re-calculate eta and correct the multiplicity accordingly 
+   * 
+   * @param d          Detector 
+   * @param r          Ring
+   * @param s          Sector 
+   * @param t          Strip 
+   * @param vz         Ip Z-coordinate
+   * @param mult       In/out multiplicity
+   * @param eta        In/out eta 
+   * @param cosTheta   On return, the cosine of theta or null
+   */
+  void RecalculateEta(UShort_t d, Char_t r, UShort_t s, UShort_t t, 
+                     Double_t vz, Double_t& mult, Double_t& eta, 
+                     Double_t& cosTheta);
+  /** 
+   * Correct for noise suppression
+   * 
+   * @param f            Factor to apply
+   * @param c            Correction 
+   * @param cosTheta     Cosine to theta 
+   * @param mult         In/Out multiplity 
+   *
+   * @return true if signal is good, otherwise false 
+   */
+  Bool_t NoiseCorrect(Int_t f, Double_t c, Double_t cosTheta, Double_t& mult);
+
+  Int_t    fRecoFactor;      // Noise factor used in Reco
+  Double_t fMaxNoiseCorr;    // If noise corr above this, flag as dead 
+  Bool_t   fRecalculateEta;  // Whether to recalc eta and angle cor (disp vtx)
+  TBits    fXtraDead;        // List of extra dead channels
+  Bool_t   fHasXtraDead;     // Whether we have xtra dead channels
+  Bool_t   fInvalidIsEmpty;  // Consider kInvalidMult as zero 
+  TH1*     fNoiseChange;     // Diagnostics
+  TH1*     fEtaChange;       // Diagnostics
+  TH1*     fDeadChange;      // Diagnostics
+  
+  ClassDef(AliFMDESDFixer,1);  // Fix FMD ESD object for issues 
+};
+
+#endif
+// 
+// EOF
+// 
+  
index 6e889f14dd9c4c95c86d994b64c32054ccd74c99..516866931da7fb4e56fca1082d873b32cec7421f 100644 (file)
@@ -202,6 +202,7 @@ AliFMDEnergyFitter::CreateOutputObjects(TList* dir)
   d->Add(AliForwardUtil::MakeParameter("regCut",        fRegularizationCut));
   d->Add(AliForwardUtil::MakeParameter("deltaShift", 
                                       AliLandauGaus::EnableSigmaShift()));
+
   if (fRingHistos.GetEntries() <= 0) { 
     AliFatal("No ring histograms where defined - giving up!");
     return;
@@ -325,11 +326,12 @@ AliFMDEnergyFitter::Accumulate(const AliESDFMD& input,
          Float_t mult = input.Multiplicity(d,r,s,t);
          
          // Keep dead-channel information. 
-         if(mult == AliESDFMD::kInvalidMult || mult > 10 || mult <= 0) 
+         if (mult == AliESDFMD::kInvalidMult || mult > 10 || mult <= 0) 
            continue;
 
          // Get the pseudo-rapidity 
          Double_t eta1 = input.Eta(d,r,s,t);
+
          // Int_t ieta = fEtaAxis.FindBin(eta1);
          // if (ieta < 1 || ieta >  fEtaAxis.GetNbins()) continue; 
 
@@ -524,6 +526,18 @@ AliFMDEnergyFitter::CheckSkip(UShort_t d, Char_t r, UShort_t skips)
   return (t & skips) == t;
 }
 
+#define PF(N,V,...)                                    \
+  AliForwardUtil::PrintField(N,V, ## __VA_ARGS__)
+#define PFB(N,FLAG)                            \
+  do {                                                                 \
+    AliForwardUtil::PrintName(N);                                      \
+    std::cout << std::boolalpha << (FLAG) << std::noboolalpha << std::endl; \
+  } while(false)
+#define PFV(N,VALUE)                                   \
+  do {                                                 \
+    AliForwardUtil::PrintName(N);                      \
+    std::cout << (VALUE) << std::endl; } while(false)
+
 //____________________________________________________________________
 void
 AliFMDEnergyFitter::Print(Option_t*) const
@@ -534,34 +548,31 @@ AliFMDEnergyFitter::Print(Option_t*) const
   // Parameters:
   //    option Not used 
   //
-  char ind[gROOT->GetDirLevel()+1];
-  for (Int_t i = 0; i < gROOT->GetDirLevel(); i++) ind[i] = ' ';
-  ind[gROOT->GetDirLevel()] = '\0';
-  std::cout << ind << ClassName() << ": " << GetName() << '\n'
-           << ind << " Low cut:                " << fLowCut << " E/E_mip\n"
-           << ind << " Max(particles):         " << fNParticles << '\n'
-           << ind << " Min(entries):           " << fMinEntries << '\n'
-           << ind << " Fit range:              " 
-           << fFitRangeBinWidth << " bins\n"
-           << ind << " Make fits:              " 
-           << (fDoFits ? "yes\n" : "no\n")
-           << ind << " Make object:            " 
-           << (fDoMakeObject ? "yes\n" : "no\n")
-           << ind << " Max E/E_mip:            " << fMaxE << '\n'
-           << ind << " N bins:                 " << fNEbins << '\n'
-           << ind << " Increasing bins:        " 
-           << (fUseIncreasingBins ?"yes\n" : "no\n")
-           << ind << " max(delta p/p):         " << fMaxRelParError << '\n'
-           << ind << " max(chi^2/nu):          " << fMaxChi2PerNDF << '\n'
-           << ind << " min(a_i):               " << fMinWeight << '\n'
-           << ind << " Residuals:              "; 
+  AliForwardUtil::PrintTask(*this);
+
+  gROOT->IncreaseDirLevel();
+  PFV("Low cut [E/E_mip]",     fLowCut);
+  PFV("Max(particles)",                fNParticles);
+  PFV("Min(entries)",          fMinEntries);
+  PFV("Fit range [bins]",       fFitRangeBinWidth);
+  PFB("Make fits",              fDoFits);
+  PFB("Make object",            fDoMakeObject);
+  PFV("Max E/E_mip",           fMaxE);
+  PFV("N bins",                        fNEbins);
+  PFB("Increasing bins",        fUseIncreasingBins);
+  PFV("max(delta p/p)",        fMaxRelParError);
+  PFV("max(chi^2/nu)",         fMaxChi2PerNDF);
+  PFV("min(a_i)",              fMinWeight);
+
+  TString r = "";
   switch (fResidualMethod) { 
-  case kNoResiduals: std::cout << "None"; break;
-  case kResidualDifference: std::cout << "Difference"; break;
-  case kResidualScaledDifference: std::cout << "Scaled difference"; break;
-  case kResidualSquareDifference: std::cout << "Square difference"; break;
+  case kNoResiduals:              r = "None";       break;
+  case kResidualDifference:       r = "Difference"; break;
+  case kResidualScaledDifference: r = "Scaled difference"; break;
+  case kResidualSquareDifference: r = "Square difference"; break;
   }
-  std::cout << std::endl;
+  PFV("Residuals", r);
+  gROOT->DecreaseDirLevel();
 }
   
 //====================================================================
index 00ae2bc314ee900b34b54baa7299f8ea3d56e6f4..c0c488027d4c0e5e90f0cbb2eb2f2ea4f9247588 100644 (file)
@@ -735,7 +735,7 @@ protected:
   UShort_t        fSkips;             // Rings to skip when fitting 
   Double_t        fRegularizationCut; // When to regularize the chi^2
 
-  ClassDef(AliFMDEnergyFitter,6); //
+  ClassDef(AliFMDEnergyFitter,8); //
 };
 
 #endif
index 263f3960ccca3d89a9e26fac60d3dde45cd9c889..deeca0e5d0d5e2e3063104fb79bd574aa0b9f8b0 100644 (file)
@@ -32,7 +32,8 @@
 //====================================================================
 AliFMDEnergyFitterTask::AliFMDEnergyFitterTask()
   : AliBaseESDTask(),
-    fEventInspector(),
+    fEventInspector(),                         
+    fESDFixer(),
     fEnergyFitter(),
     fOnlyMB(false)
 {
@@ -47,6 +48,7 @@ AliFMDEnergyFitterTask::AliFMDEnergyFitterTask()
 AliFMDEnergyFitterTask::AliFMDEnergyFitterTask(const char* name)
   : AliBaseESDTask(name, "", &(AliForwardCorrectionManager::Instance())), 
     fEventInspector("event"),
+    fESDFixer("esdFizer"),
     fEnergyFitter("energy"),
     fOnlyMB(false)
 {
@@ -108,14 +110,38 @@ AliFMDEnergyFitterTask::Book()
   DGUARD(fDebug,1,"Create output objects of AliFMDEnergyFitterTask");
 
   // We don't need any corrections for this task 
-  fNeededCorrections = 0;
+  fNeededCorrections = 0; 
   fExtraCorrections  = 0;
+  if (fESDFixer.IsUseNoiseCorrection()) 
+    fNeededCorrections = AliForwardCorrectionManager::kNoiseGain;
 
+  fESDFixer    .CreateOutputObjects(fList);
   fEnergyFitter.CreateOutputObjects(fList);
 
   fList->Add(AliForwardUtil::MakeParameter("onlyMB", fOnlyMB));
   return true;
 }
+//____________________________________________________________________
+void
+AliFMDEnergyFitterTask::PreCorrections(const AliESDEvent* esd)
+{
+  if (!esd) return; 
+  
+  AliESDFMD* esdFMD = esd->GetFMDData();  
+  if (!esdFMD) return;
+
+  // TODO: We should always disable this on MC!
+  Int_t tgt = fESDFixer.FindTargetNoiseFactor(*esdFMD, false);
+  if (tgt <= 0) {
+    // If the target noise factor is 0 or less, disable the noise/gain
+    // correction.
+    fESDFixer.SetRecoNoiseFactor(4);
+    fNeededCorrections ^= AliForwardCorrectionManager::kNoiseGain;
+  }
+  else 
+    AliWarning("The noise corrector has been enabled!");
+}
+
 //____________________________________________________________________
 Bool_t
 AliFMDEnergyFitterTask::PreData(const TAxis& /*vertex*/, const TAxis& eta)
@@ -128,6 +154,7 @@ AliFMDEnergyFitterTask::PreData(const TAxis& /*vertex*/, const TAxis& eta)
 
   fEnergyFitter.SetupForData(eta);
 
+  Print();
   return true;
 }
 
@@ -183,6 +210,10 @@ AliFMDEnergyFitterTask::Event(AliESDEvent& esd)
   
   // Get FMD data 
   AliESDFMD* esdFMD = esd.GetFMDData();
+
+  // Fix up ESD 
+  fESDFixer.Fix(*esdFMD, ip.Z());
+
   // Do the energy stuff 
   if (!fEnergyFitter.Accumulate(*esdFMD, cent, 
                                triggers & AliAODForwardMult::kEmpty)){
@@ -229,6 +260,7 @@ AliFMDEnergyFitterTask::Print(Option_t* option) const
   AliBaseESDTask::Print(option);
   gROOT->IncreaseDirLevel();
   PFB("Only MB", fOnlyMB);
+  fESDFixer    .Print(option);
   fEnergyFitter.Print(option);
   gROOT->DecreaseDirLevel();
 }
index a04aef52fceb6dbb6308cc06bab9c3c0c68ea5d2..6a4e0b326cc58832bf4a39484af34ea0aa553dde 100644 (file)
@@ -16,6 +16,7 @@
 #include "AliBaseESDTask.h"
 #include "AliFMDEventInspector.h"
 #include "AliFMDEnergyFitter.h"
+#include "AliFMDESDFixer.h"
 class AliESDEvent;
 class TH2D;
 class TList;
@@ -75,6 +76,14 @@ public:
    * @return true on success. 
    */
   virtual Bool_t Book();
+  /** 
+   * Called on first event _before_ reading corrections.  Here, the
+   * user class can do additional checking to see if the some (more or
+   * less) corrections are needed.
+   * 
+   * @param esd Event 
+   */
+  virtual void PreCorrections(const AliESDEvent* esd);
   /** 
    * Called after reading in the first event. Here we can setup stuff
    * depending on the conditions we're running under.
@@ -121,6 +130,12 @@ public:
    * @return Reference to AliFMDEventInspector object 
    */
   const AliFMDEventInspector& GetEventInspector() const{return fEventInspector;}
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  AliFMDESDFixer& GetESDFixer() { return fESDFixer; }
   /**
    * Get reference to the EnergyFitter algorithm 
    * 
@@ -188,6 +203,7 @@ protected:
   AliFMDEnergyFitterTask& operator=(const AliFMDEnergyFitterTask& o);
 
   AliFMDEventInspector fEventInspector; // Algorithm
+  AliFMDESDFixer       fESDFixer;       // Algorithm
   AliFMDEnergyFitter   fEnergyFitter;   // Algorithm
   Bool_t               fOnlyMB;         // Only MB flag
 
index 9d0cf32731f6de8b582b1001c4df0356f1551b20..4a0619575d05b5933fc05d4d73d99658fef7d443 100644 (file)
@@ -95,18 +95,19 @@ AliFMDEventInspector::AliFMDEventInspector()
     fDebug(0),
     fCentAxis(0),
     fVtxAxis(10,-10,10),
-    fUseFirstPhysicsVertex(false),
+    fVtxMethod(kNormal),
+    // fUseFirstPhysicsVertex(false),
     fUseV0AND(false),
     fMinPileupContrib(3), 
     fMinPileupDistance(0.8),
-    fUseDisplacedVertices(false),
+// fUseDisplacedVertices(false),
     fDisplacedVertex(),
     fCollWords(),
     fBgWords(),
     fCentMethod("V0M"),
     fMinCent(-1.0),
     fMaxCent(-1.0),
-    fUsepA2012Vertex(false),
+// fUsepA2012Vertex(false),
     fRunNumber(0),
     fMC(false),
   fProdYear(-1),
@@ -146,18 +147,19 @@ AliFMDEventInspector::AliFMDEventInspector(const char* name)
     fDebug(0),
     fCentAxis(0),
     fVtxAxis(10,-10,10),
-    fUseFirstPhysicsVertex(false),
+    fVtxMethod(kNormal),
+// fUseFirstPhysicsVertex(false),
     fUseV0AND(false),
     fMinPileupContrib(3), 
     fMinPileupDistance(0.8),
-    fUseDisplacedVertices(false),
+// fUseDisplacedVertices(false),
   fDisplacedVertex(),
   fCollWords(),
   fBgWords(),
   fCentMethod("V0M"),
   fMinCent(-1.0),
   fMaxCent(-1.0),
-  fUsepA2012Vertex(false),
+// fUsepA2012Vertex(false),
   fRunNumber(0),
   fMC(false),
   fProdYear(-1),
@@ -574,7 +576,7 @@ AliFMDEventInspector::SetupForData(const TAxis& vtxAxis)
   xAxis->SetBinLabel(kOther,           "Other");
   fList->Add(fHTrgStatus);
 
-  if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", false);
+  if (AllowDisplaced()) fDisplacedVertex.SetupForData(fList, "", false);
 }
 
 //____________________________________________________________________
@@ -592,11 +594,12 @@ AliFMDEventInspector::StoreInformation()
   fList->Add(AliForwardUtil::MakeParameter("field", fField));
   fList->Add(AliForwardUtil::MakeParameter("runNo", fRunNumber));
   fList->Add(AliForwardUtil::MakeParameter("lowFlux", fLowFluxCut));
-  fList->Add(AliForwardUtil::MakeParameter("fpVtx",fUseFirstPhysicsVertex));
+  fList->Add(AliForwardUtil::MakeParameter("ipMethod", fVtxMethod));
+  //fList->Add(AliForwardUtil::MakeParameter("fpVtx",fUseFirstPhysicsVertex));
   fList->Add(AliForwardUtil::MakeParameter("v0and",fUseV0AND));
   fList->Add(AliForwardUtil::MakeParameter("nPileUp", fMinPileupContrib));
   fList->Add(AliForwardUtil::MakeParameter("dPileup", fMinPileupDistance));
-  fList->Add(AliForwardUtil::MakeParameter("satellite", fUseDisplacedVertices));
+  fList->Add(AliForwardUtil::MakeParameter("satellite", AllowDisplaced()));
   fList->Add(AliForwardUtil::MakeParameter("alirootRev", 
                                           AliForwardUtil::AliROOTRevision()));
   fList->Add(AliForwardUtil::MakeParameter("alirootBranch", 
@@ -713,7 +716,7 @@ AliFMDEventInspector::Process(const AliESDEvent* event,
   }
 
   // --- Process satellite event information is requested ------------
-  if (fUseDisplacedVertices) { 
+  if (AllowDisplaced()) { 
     if (!fDisplacedVertex.Process(event)) 
       AliWarning("Failed to process satellite event");
   }
@@ -823,7 +826,7 @@ AliFMDEventInspector::ReadCentrality(const AliESDEvent& esd,
 
   // We overwrite with satellite events, so we can be sure to get the
   // centrality determination from the satellite vertex selection
-  if (fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
+  if (AllowDisplaced() && fDisplacedVertex.IsSatellite()) {
     cent = fDisplacedVertex.GetCentralityPercentile();
     qual = 0;
   }
@@ -893,7 +896,7 @@ AliFMDEventInspector::ReadTriggers(const AliESDEvent& esd, UInt_t& triggers,
   if (trigStr.IsNull()) fHTrgStatus->Fill(kNoTrgWords);
   if (fHWords) fHWords->Fill(trigStr.Data(), 1);
   
-  if(fUseDisplacedVertices) {
+  if(AllowDisplaced()) {
     DMSG(fDebug,3,"Using displaced vertex stuff");
     // if (TMath::Abs(fDisplacedVertex.GetVertexZ()) >= 999) offline = false;
     if (fDisplacedVertex.IsSatellite()) 
@@ -1227,14 +1230,16 @@ AliFMDEventInspector::ReadVertex(const AliESDEvent& esd, TVector3& ip)
   ip.SetXYZ(1024, 1024, 0);
   
   EVtxStatus s = kNoVtx;
-  if (fUseDisplacedVertices && fDisplacedVertex.IsSatellite()) {
+  if (AllowDisplaced() && fDisplacedVertex.IsSatellite()) {
     s = kVtxOK;
     ip.SetZ(fDisplacedVertex.GetVertexZ());
   }
-  else if (fUseFirstPhysicsVertex
+  else if (fVtxMethod == kPWGUD
     s = CheckPWGUDVertex(esd, ip);
-  else if (fUsepA2012Vertex
+  else if (fVtxMethod == kpA2012
     s = CheckpA2012Vertex(esd,ip);     
+  else if (fVtxMethod == kpA2013) 
+    s = CheckpA2013Vertex(esd,ip);     
   else 
     s = CheckVertex(esd, ip);
   
@@ -1285,9 +1290,10 @@ AliFMDEventInspector::EVtxStatus
 AliFMDEventInspector::CheckpA2012Vertex(const AliESDEvent& esd, 
                                        TVector3& ip)  const
 {      
+  const Int_t nMinContrib = 0;
   const AliESDVertex *vertex = esd.GetPrimaryVertexSPD();
   if (!vertex) return kNoSPDVtx;
-  if (vertex->GetNContributors() <= 0) return kFewContrib;
+  if (vertex->GetNContributors() <= nMinContrib) return kFewContrib;
   
   TString vtxTyp = vertex->GetTitle();
   if (vtxTyp.Contains("vertexer: Z")) return kNotVtxZ;
@@ -1301,6 +1307,39 @@ AliFMDEventInspector::CheckpA2012Vertex(const AliESDEvent& esd,
 
   return kVtxOK;
 }
+//____________________________________________________________________
+AliFMDEventInspector::EVtxStatus
+AliFMDEventInspector::CheckpA2013Vertex(const AliESDEvent& esd, 
+                                       TVector3& ip)  const
+{      
+  // This code is adopted from 
+  // 
+  //   AliAnalysisUtils::IsVertexSelected2013pA(AliVEvent *event)
+  // 
+  const Int_t nMinContrib = 0;
+  const AliESDVertex* primVtx = esd.GetPrimaryVertex();
+  if (!primVtx) return kNoVtx;
+  if (primVtx->GetNContributors() <= nMinContrib) return kFewContrib;
+  
+  const AliESDVertex* spdVtx = esd.GetPrimaryVertex();
+  if (!spdVtx) return kNoSPDVtx;
+  if (spdVtx->GetNContributors() <= nMinContrib) return kFewContrib;
+  
+  TString vtxTyp = spdVtx->GetTitle();
+  if (vtxTyp.Contains("vertexer: Z")) {
+    if (spdVtx->GetZRes()>0.25) return kUncertain;
+  }
+  Bool_t correlateVtx = true;
+  if (correlateVtx) { 
+    if (TMath::Abs(spdVtx->GetZ() - primVtx->GetZ()) > 0.5) return kUncertain;
+  }
+
+  ip.SetX(primVtx->GetX());
+  ip.SetY(primVtx->GetY());
+  ip.SetZ(primVtx->GetZ());            
+
+  return kVtxOK;
+}
 
 //____________________________________________________________________
 AliFMDEventInspector::EVtxStatus
@@ -1450,9 +1489,18 @@ AliFMDEventInspector::Print(Option_t*) const
   PFV("System", AliForwardUtil::CollisionSystemString(fCollisionSystem));
   PFV("CMS energy per nucleon",        sNN);
   PFV("Field",                 field);
-  PFB("Satellite events",      fUseDisplacedVertices);
-  PFB("Use 2012 pA vertex",    fUsepA2012Vertex );
-  PFB("Use PWG-UD vertex",     fUseFirstPhysicsVertex);
+  PFB("Satellite events",      AllowDisplaced());
+  // PFB("Use 2012 pA vertex", fUsepA2012Vertex );
+  // PFB("Use PWG-UD vertex",  fUseFirstPhysicsVertex);
+  TString vtxMethod("normal");
+  switch(fVtxMethod) { 
+  case kNormal:    vtxMethod = "normal" ;    break;
+  case kpA2012:    vtxMethod = "pA 2012";    break;
+  case kpA2013:    vtxMethod = "pA 2013";    break;
+  case kPWGUD:     vtxMethod = "PWG-UD";     break;
+  case kDisplaced: vtxMethod = "Satellite";  break;
+  }
+  PFV("IP method",              vtxMethod);
   PFB("Simulation input",      fMC );
   PFV("Centrality method",     fCentMethod);
   PFV("Centrality axis",        (!fCentAxis ? "none" : ""));
index 9d89019eb73e7e9fae57569ddf6840161812d076..7167342d49f35e8b42c2415ef6f2e1575b951129 100644 (file)
@@ -88,6 +88,13 @@ public:
     kSatellite,
     kOffline
   };
+  enum EVtxType { 
+    kNormal, 
+    kpA2012, 
+    kpA2013, 
+    kPWGUD, 
+    kDisplaced
+  };
   /** 
    * Centrality methods 
    */
@@ -208,18 +215,28 @@ public:
    * @param c Maximum error (in centimeters)
    */
   void SetMaxVzErr(Double_t c=0.1) { fMaxVzErr = c; }
+  /** 
+   * Set the vertex method to use 
+   * 
+   * @param t Method 
+   */
+  void SetVertexMethod(EVtxType t) { fVtxMethod = t; }
   /** 
    * Use the first physics vtx code.   
    * 
    * @param use Use it or not 
+   * @deprecated Use SetVertexMethod 
    */
-  void SetUseFirstPhysicsVtx(Bool_t use) {fUseFirstPhysicsVertex = use; }
+  void SetUseFirstPhysicsVtx(Bool_t use) {
+    SetVertexMethod(use ? kPWGUD : kNormal); }
  /** 
    * Use the first physics vtx code.   
    * 
    * @param use Use it or not 
+   * @deprecated Use SetVertexMethod 
    */
-  void SetpA2012Vtx(Bool_t use) {fUsepA2012Vertex= use; }
+  void SetpA2012Vtx(Bool_t use) {
+    SetVertexMethod(use ? kpA2012 : kNormal); }
  /** 
    * Use the 2012 pA vtx code.   
    * 
@@ -249,12 +266,11 @@ public:
    * Enable selection of displaced vertices. 
    * 
    * @param use whether to use
+   * @deprecated Use SetVertexMethod 
    */
-  void SetUseDisplacedVertices(Bool_t use=true)
-  {
-    fUseDisplacedVertices = use;
-  }  
-  Bool_t IsUseDisplacedVertices() const { return fUseDisplacedVertices; }
+  void SetUseDisplacedVertices(Bool_t use=true) { 
+    SetVertexMethod(use ? kDisplaced : kNormal); }  
+  Bool_t IsUseDisplacedVertices() const { return AllowDisplaced(); }
   /** 
    * Set the lower centrality cut - if negative, do not use 
    *
@@ -447,6 +463,7 @@ protected:
    * @return Reference to this object
    */
   AliFMDEventInspector& operator=(const AliFMDEventInspector& o);
+  Bool_t AllowDisplaced() const { return fVtxMethod == kDisplaced; }
   /** 
    * Cache the configure trigger classes from the physis selection.  
    * 
@@ -588,7 +605,24 @@ protected:
    * @return status
    */
   virtual EVtxStatus CheckpA2012Vertex(const AliESDEvent& esd, 
-                                    TVector3& ip) const;
+                                      TVector3& ip) const;
+  /** 
+   * Check the vertex. That is
+   *
+   * - Check if we have an normal vertex and that it's status is OK 
+   * - Check that we have enough contributors 
+   * - Check if we have an SPD vertex
+   * - Check that we have enough contributors 
+   * - If from Z-vertexer, check the resolution
+   * - Check that the two found vertices are within 0.5cm of each other
+   * 
+   * @param esd Data 
+   * @param ip  On return, the coordinates of the IP
+   * 
+   * @return status
+   */
+  virtual EVtxStatus CheckpA2013Vertex(const AliESDEvent& esd, 
+                                      TVector3& ip) const;
   /** 
    * Check the vertex for pA 2012 settings. That is
    *
@@ -598,10 +632,6 @@ protected:
    * 
    * @return true if the vertex was found and met the requirements
    */
-
-
-
-
   virtual EVtxStatus CheckVertex(const AliESDEvent& esd, TVector3& ip) const;
   /** 
    * Read centrality from event 
@@ -637,18 +667,19 @@ protected:
   Int_t    fDebug;                //  Debug level 
   TAxis*   fCentAxis;             // Centrality axis used in histograms
   TAxis    fVtxAxis;              // IP_z Axis 
-  Bool_t   fUseFirstPhysicsVertex;//Use the vtx code from p+p first physics
+  EVtxType fVtxMethod;            // Vertex method to use 
+  // Bool_t   fUseFirstPhysicsVertex;//Use the vtx code from p+p first physics
   Bool_t   fUseV0AND;             // Use the vtx code from p+p first physics
   UShort_t fMinPileupContrib;     // Min contributors to 2nd pile-up IP
   Double_t fMinPileupDistance;    // Min distance of 2nd pile-up IP
-  Bool_t   fUseDisplacedVertices; // Analyze displaced vertices?
+  // Bool_t   fUseDisplacedVertices; // Analyze displaced vertices?
   AliDisplacedVertexSelection fDisplacedVertex; //Displaced vertex selector
   TList    fCollWords;            //! Configured collision words 
   TList    fBgWords;              //! Configured background words 
   TString  fCentMethod;           // Centrality method
   Double_t fMinCent;              // min centrality
   Double_t fMaxCent;              // max centrailty
-  Bool_t   fUsepA2012Vertex;      // flag to use pA2012 Veretx selection
+  // Bool_t   fUsepA2012Vertex;      // flag to use pA2012 Veretx selection
   ULong_t  fRunNumber;            // Current run number 
   Bool_t   fMC;                   // Is this MC input
   Short_t  fProdYear;             // Production year 
@@ -657,7 +688,7 @@ protected:
   Int_t    fProdSVN;              // AliROOT revision used in production
   Bool_t   fProdMC;               // True if anchor production
 
-  ClassDef(AliFMDEventInspector,12); // Inspect the event 
+  ClassDef(AliFMDEventInspector,13); // Inspect the event 
 };
 
 #endif
index e1367e6dee6037eb004dd05876cd7cd8ad705b7c..0d8dca4e3717cf020786869a657c89d4c56d5c96 100644 (file)
@@ -145,10 +145,10 @@ AliFMDMCEventInspector::SetupForData(const TAxis& vtxAxis)
 
   // We temporary disable the displaced vertices so we can initialize
   // the routine ourselves.
-  Bool_t saveDisplaced  = fUseDisplacedVertices;
-  fUseDisplacedVertices = false;
+  Bool_t saveDisplaced  = AllowDisplaced();
+  if (saveDisplaced) SetVertexMethod(kNormal);
   AliFMDEventInspector::SetupForData(vtxAxis);
-  fUseDisplacedVertices = saveDisplaced;
+  if (saveDisplaced) SetVertexMethod(kDisplaced);
 
   Int_t    maxPart = 450;
   Int_t    maxBin  = 225;
@@ -253,7 +253,7 @@ AliFMDMCEventInspector::SetupForData(const TAxis& vtxAxis)
   fHCentVsMcC->SetZTitle("Events");
   fList->Add(fHCentVsMcC);
 
-  if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", true);
+  if (AllowDisplaced()) fDisplacedVertex.SetupForData(fList, "", true);
 }
 
 namespace
@@ -459,7 +459,7 @@ AliFMDMCEventInspector::ProcessMC(AliMCEvent*       event,
   fHBvsPart->Fill(b, npart);
   fHBvsBin->Fill(b, nbin);
 
-  if(fUseDisplacedVertices) {
+  if(AllowDisplaced()) {
 #if 0
     // Put the vertex at fixed locations 
     Double_t zvtx  = vz;
index fda63dd65c01a342c743c6c98a80c339c7d5f966..f3dae49a3a4f10bac03e12d0d3dce26d5a0e4811 100644 (file)
@@ -64,10 +64,6 @@ AliFMDSharingFilter::AliFMDSharingFilter()
     fHCuts(),
     fUseSimpleMerging(false),
     fThreeStripSharing(true),
-    fRecalculateEta(false),
-    // fExtraDead(0),
-    fXtraDead(0),
-    fInvalidIsEmpty(false),
     fMergingDisabled(false),
     fIgnoreESDForAngleCorrection(false)
 {
@@ -92,10 +88,6 @@ AliFMDSharingFilter::AliFMDSharingFilter(const char* title)
     fHCuts(),
     fUseSimpleMerging(false),
     fThreeStripSharing(true),
-    fRecalculateEta(false), 
-    // fExtraDead(51200),
-    fXtraDead(AliFMDStripIndex::Pack(3,'O',19,511)+1),
-    fInvalidIsEmpty(false),
     fMergingDisabled(false),
     fIgnoreESDForAngleCorrection(false)
 {
@@ -156,84 +148,6 @@ AliFMDSharingFilter::GetRingHistos(UShort_t d, Char_t r) const
   return static_cast<RingHistos*>(fRingHistos.At(idx));
 }
 
-//____________________________________________________________________
-void
-AliFMDSharingFilter::AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t)
-{
-  if (d < 1 || d > 3) {
-    Warning("AddDead", "Invalid detector FMD%d", d);
-    return;
-  }
-  Bool_t inner = (r == 'I' || r == 'i');
-  if (d == 1 && !inner) { 
-    Warning("AddDead", "Invalid ring FMD%d%c", d, r);
-    return;
-  }
-  if ((inner && s >= 20) || (!inner && s >= 40)) { 
-    Warning("AddDead", "Invalid sector FMD%d%c[%02d]", d, r, s);
-    return;
-  }
-  if ((inner && t >= 512) || (!inner && t >= 256)) { 
-    Warning("AddDead", "Invalid strip FMD%d%c[%02d,%03d]", d, r, s, t);
-    return;
-  }
-    
-  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
-  // Int_t i  = 0;
-  fXtraDead.SetBitNumber(id, true);
-#if 0
-  for (i = 0; i < fExtraDead.GetSize(); i++) {
-    Int_t j = fExtraDead.At(i);
-    if (j == id) return; // Already there 
-    if (j <  0) break; // Free slot 
-  }
-  if (i >= fExtraDead.GetSize()) { 
-    Warning("AddDead", "No free slot to add FMD%d%c[%02d,%03d] at", 
-           d, r, s, t);
-    return;
-  }
-  fExtraDead[i] = id;
-#endif
-}
-//____________________________________________________________________
-void
-AliFMDSharingFilter::AddDeadRegion(UShort_t d,  Char_t r, 
-                                  UShort_t s1, UShort_t s2, 
-                                  UShort_t t1, UShort_t t2)
-{
-  // Add a dead region spanning from FMD<d><r>[<s1>,<t1>] to 
-  // FMD<d><r>[<s2>,<t2>] (both inclusive)
-  for (Int_t s = s1; s <= s2; s++) 
-    for (Int_t t = t1; t <= t2; t++) 
-      AddDead(d, r, s, t);
-}
-//____________________________________________________________________
-void
-AliFMDSharingFilter::AddDead(const Char_t* script)
-{
-  if (!script || script[0] == '\0') return;
-  
-  gROOT->Macro(Form("%s((AliFMDSharingFilter*)%p);", script, this));
-}
-
-//____________________________________________________________________
-Bool_t
-AliFMDSharingFilter::IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const
-{
-  Int_t id = AliFMDStripIndex::Pack(d, r, s, t);
-  return fXtraDead.TestBitNumber(id); 
-#if 0
-  for (Int_t i = 0; i < fExtraDead.GetSize(); i++) {
-    Int_t j = fExtraDead.At(i);
-    if (j == id) {
-      //Info("IsDead", "FMD%d%c[%02d,%03d] marked as dead here", d, r, s, t);
-      return true;
-    }
-    if (j < 0) break; // High water mark 
-  }
-  return false;
-#endif
-}
 //____________________________________________________________________
 void
 AliFMDSharingFilter::SetupForData(const TAxis& axis)
@@ -243,9 +157,6 @@ AliFMDSharingFilter::SetupForData(const TAxis& axis)
   AliForwardCorrectionManager& fcm  = AliForwardCorrectionManager::Instance();
   const AliFMDCorrELossFit*    fits = fcm.GetELossFit();
 
-  // Compactify the xtra dead bits 
-  fXtraDead.Compact();
-
   // Get the high cut.  The high cut is defined as the 
   // most-probably-value peak found from the energy distributions, minus 
   // 2 times the width of the corresponding Landau.
@@ -287,7 +198,7 @@ Bool_t
 AliFMDSharingFilter::Filter(const AliESDFMD& input, 
                            Bool_t           /*lowFlux*/,
                            AliESDFMD&       output, 
-                           Double_t         zvtx)
+                           Double_t         /*zvtx*/)
 {
   // 
   // Filter the input AliESDFMD object
@@ -348,42 +259,19 @@ AliFMDSharingFilter::Filter(const AliESDFMD& input,
          Double_t phi = input.Phi(d,r,s,t) * TMath::Pi() / 180.;
          if (s == 0) output.SetEta(d,r,s,t,eta);
          
-         if(fRecalculateEta) { 
-           Double_t etaOld  = eta;
-           Double_t etaCalc = AliForwardUtil::GetEtaFromStrip(d,r,s,t,zvtx);
-           eta              = etaCalc;
-           
-           if (mult > 0 && mult != AliESDFMD::kInvalidMult ) {
-             Double_t cosOld =  ETA2COS(etaOld);
-             Double_t cosNew =  ETA2COS(etaCalc);
-             Double_t corr   =  cosNew / cosOld;
-             mult            *= corr;
-             multNext        *= corr;
-             multNextNext    *= corr;
-           }
-         } // Recalculate eta 
-
-         // Special case for pre revision 43611 AliFMDReconstructor.
-         // If fInvalidIsEmpty and we get an invalid signal from the
-         // ESD, then we need to set this signal to zero.  Note, dead
-         // strips added in the ForwardAODConfig.C file are not
-         // effected by this, and we can use that to set the proper
-         // dead strips.
-         if (mult == AliESDFMD::kInvalidMult && fInvalidIsEmpty) 
-           mult = 0;
-
          // Keep dead-channel information - either from the ESD (but
          // see above for older data) or from the settings in the
          // ForwardAODConfig.C file.
-         if (mult == AliESDFMD::kInvalidMult || IsDead(d,r,s,t)) {
+         if (mult == AliESDFMD::kInvalidMult) {
            output.SetMultiplicity(d,r,s,t,AliESDFMD::kInvalidMult);
            histos->fBefore->Fill(-1);
            mult = AliESDFMD::kInvalidMult;
          }
          
-         if (mult != AliESDFMD::kInvalidMult) 
+         if (mult != AliESDFMD::kInvalidMult) {
            // Always fill the ESD sum histogram 
            histos->fSumESD->Fill(eta, phi, mult);
+         }
 
          // If no signal or dead strip, go on. 
          if (mult == AliESDFMD::kInvalidMult || mult == 0) {
@@ -789,35 +677,7 @@ AliFMDSharingFilter::CreateOutputObjects(TList* dir)
   TParameter<int>* nFiles = new TParameter<int>("nFiles", 1);
   nFiles->SetMergeMode('+');
   d->Add(nFiles);
-  
-  TObjArray* extraDead = new TObjArray;
-  extraDead->SetOwner();
-  extraDead->SetName("extraDead");
-#if 0
-  for (Int_t i = 0; i < fExtraDead.GetSize(); i++) { 
-    if (fExtraDead.At(i) < 0) break;
-    UShort_t dd, s, t;
-    Char_t  r;
-    Int_t   id = fExtraDead.At(i);
-    AliFMDStripIndex::Unpack(id, dd, r, s, t);
-    extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
-                                                     dd, r, s, t), id));
-  }
-#endif
-  fXtraDead.Compact();
-  UInt_t firstBit = fXtraDead.FirstSetBit();
-  UInt_t nBits    = fXtraDead.GetNbits();
-  for (UInt_t i = firstBit; i < nBits; i++) {
-    if (!fXtraDead.TestBitNumber(i)) continue;
-    UShort_t dd, s, t;
-    Char_t  r;
-    AliFMDStripIndex::Unpack(i, dd, r, s, t);
-    extraDead->Add(AliForwardUtil::MakeParameter(Form("FMD%d%c[%02d,%03d]",
-                                                     dd, r, s, t), 
-                                                UShort_t(i)));
-  }
-  // d->Add(&fXtraDead);
-  d->Add(extraDead);
+
   fLCuts.Output(d,"lCuts");
   fHCuts.Output(d,"hCuts");
 
@@ -855,7 +715,6 @@ AliFMDSharingFilter::Print(Option_t* /*option*/) const
   PFB("Use corrected angles",  fCorrectAngles);
   PFB("Zero below threshold",  fZeroSharedHitsBelowThreshold);
   PFB("Use simple sharing",    fUseSimpleMerging);
-  PFB("Consider invalid null", fInvalidIsEmpty);
   PFB("Allow 3 strip merging", fThreeStripSharing);
   PF("Low cuts",       "");
   fLCuts.Print();
index 1e3c8723a2fd158cc721dd5f29eaa29d9436fa75..37ea26178b29789bbc1f9aabacefcf353c24ffc1 100644 (file)
@@ -21,7 +21,6 @@
 #include <TList.h>
 #include "AliForwardUtil.h"
 #include "AliFMDMultCuts.h"
-#include <TBits.h>
 class AliESDFMD;
 class TAxis;
 class TList;
@@ -133,20 +132,6 @@ public:
    * 
    */
   void SetAllow3Strips(Bool_t use) { fThreeStripSharing = use; }  
-  /** 
-   * In case of a displaced vertices recalculate eta and angle correction
-   * 
-   * @param use recalculate or not
-   * 
-   */
-  void SetRecalculateEta(Bool_t use) { fRecalculateEta = use; }
-  /** 
-   * Set whether to consider invalid multiplicities as null (or empty)
-   * signal. 
-   * 
-   * @param flag If true, count invalids as empty
-   */
-  void SetInvalidIsEmpty(Bool_t flag) { fInvalidIsEmpty = flag; }
   /**
    * Set whether to ignore the ESD info when angle correcting, this
    * is to counter a known issue where the info in the ESD is incorrect
@@ -245,48 +230,6 @@ public:
    */  
   void SetHCuts(const AliFMDMultCuts& c) { fHCuts = c; }
   /* @} */
-  /** 
-   * @{ 
-   * @name Dead strip handling 
-   */
-  /** 
-   * Add a dead strip
-   * 
-   * @param d  Detector
-   * @param r  Ring 
-   * @param s  Sector 
-   * @param t  Strip
-   */
-  void AddDead(UShort_t d, Char_t r, UShort_t s, UShort_t t);
-  /** 
-   * Add a dead region in a detector ring
-   * 
-   * @param d   Detector
-   * @param r   Ring
-   * @param s1  First sector (inclusive)
-   * @param s2  Last sector (inclusive)
-   * @param t1  First strip (inclusive)
-   * @param t2  Last strip (inclusive)
-   */
-  void AddDeadRegion(UShort_t d, Char_t r, UShort_t s1, UShort_t s2, 
-                    UShort_t t1, UShort_t t2);
-  /** 
-   * Add dead strips from a script.  The script is supposed to accept
-   * a pointer to this object (AliFMDSharingFilter) and then call
-   * AddDead or AddDeadRegion as needed.
-   * 
-   * @code 
-   * void deadstrips(AliFMDSharingFilter* filter)
-   * { 
-   *   filter->AddDead(...);
-   *   // ... and so on 
-   * }
-   * @endcode 
-   *
-   * @param script The script to read dead strips from. 
-   */
-  void AddDead(const Char_t* script);
-  /* @} */
 protected:
   /** 
    * Copy constructor - not implemented
@@ -443,17 +386,6 @@ protected:
    * @return 
    */
   virtual Double_t GetLowCut(UShort_t d, Char_t r, Double_t eta) const;
-  /** 
-   * Check if a strip is marked as dead 
-   * 
-   * @param d   Detector 
-   * @param r   Ring
-   * @param s   Sector 
-   * @param t   Strip 
-   * 
-   * @return true if dead 
-   */
-  virtual Bool_t IsDead(UShort_t d, Char_t r, UShort_t s, UShort_t t) const;
   TList    fRingHistos;      // List of histogram containers
   Bool_t   fCorrectAngles;   // Whether to work on angle corrected signals
   TH2*     fHighCuts;        // High cuts used
@@ -464,12 +396,9 @@ protected:
   AliFMDMultCuts fHCuts;     // Cuts object for high cuts
   Bool_t   fUseSimpleMerging;// enable simple sharing by HHD
   Bool_t   fThreeStripSharing; //In case of simple sharing allow 3 strips
-  Bool_t   fRecalculateEta;  // Whether to recalc eta and angle cor (disp vtx)
-  TBits    fXtraDead;
-  Bool_t   fInvalidIsEmpty;  // Consider kInvalidMult as zero 
   Bool_t   fMergingDisabled; // If true, do not merge
   Bool_t   fIgnoreESDForAngleCorrection; // Ignore ESD information when angle correcting
-  ClassDef(AliFMDSharingFilter,10); //
+  ClassDef(AliFMDSharingFilter,11); //
 };
 
 #endif
index e2360ba59a11d16bb466dcb9ee86353ac943d213..dad4cbf4eaeb6f246a513df9af9738cc763acfc4 100644 (file)
@@ -8,6 +8,7 @@
 #include "AliFMDCorrVertexBias.h"
 #include "AliFMDCorrMergingEfficiency.h"
 #include "AliFMDCorrAcceptance.h"
+#include "AliFMDCorrNoiseGain.h"
 #include "AliForwardUtil.h"
 #include "AliOADBForward.h"
 #include <TString.h>
@@ -28,6 +29,7 @@ const char* AliForwardCorrectionManager::fgkELossFitsSkel    = "elossfits";
 const char* AliForwardCorrectionManager::fgkVertexBiasSkel   = "vertexbias";
 const char* AliForwardCorrectionManager::fgkMergingEffSkel   = "merging";
 const char* AliForwardCorrectionManager::fgkAcceptanceSkel   = "acceptance";
+const char* AliForwardCorrectionManager::fgkNoiseGainSkel    = "noisegain";
 
 #define PREFIX  "$(ALICE_ROOT)/OADB/PWGLF/FORWARD/CORRECTIONS/data/"
 #define DB_NAME "fmd_corrections.root"
@@ -79,17 +81,19 @@ AliForwardCorrectionManager::AliForwardCorrectionManager(Bool_t d)
   RegisterCorrection(kIdAcceptance, fgkAcceptanceSkel, 
                     PREFIX DB_NAME, AliFMDCorrAcceptance::Class(),
                     kRun|kSys|kSNN|kSatellite);
+  RegisterCorrection(kIdNoiseGain, fgkNoiseGainSkel,
+                    PREFIX DB_NAME, AliFMDCorrNoiseGain::Class(), kRun);
 }
 //____________________________________________________________________
 Bool_t
 AliForwardCorrectionManager::Init(ULong_t     runNo, 
-                                     const char* sys, 
-                                     Float_t     sNN, 
-                                     Float_t     field,
-                                     Bool_t      mc,
-                                     Bool_t      sat,
-                                     UInt_t      what,
-                                     Bool_t      force)
+                                 const char* sys, 
+                                 Float_t     sNN, 
+                                 Float_t     field,
+                                 Bool_t      mc,
+                                 Bool_t      sat,
+                                 UInt_t      what,
+                                 Bool_t      force)
 {
   // 
   // Read in correction based on passed parameters
@@ -145,6 +149,7 @@ AliForwardCorrectionManager::Init(ULong_t  runNo,
   EnableCorrection(kIdAcceptance,      what & kAcceptance);
   EnableCorrection(kIdVertexBias,      what & kVertexBias);
   EnableCorrection(kIdMergingEfficiency,what & kMergingEfficiency);
+  EnableCorrection(kIdNoiseGain,       what & kNoiseGain);
   
   return InitCorrections(runNo, sys, sNN, field, mc, sat, force);
 }
@@ -176,6 +181,8 @@ AliForwardCorrectionManager::ParseFields(const TString& fields)
       ret |= kMergingEfficiency;
     else if (str.Contains(fgkAcceptanceSkel, TString::kIgnoreCase))
       ret |= kAcceptance;
+    else if (str.Contains(fgkNoiseGainSkel, TString::kIgnoreCase))
+      ret |= kNoiseGain;
     else 
       AliWarningClassF("Unknown correction: %s", str.Data());
   }
@@ -265,6 +272,17 @@ AliForwardCorrectionManager::GetAcceptance() const
    */
   return static_cast<const AliFMDCorrAcceptance*>(Get(kIdAcceptance)); 
 }
+//____________________________________________________________________
+const AliFMDCorrNoiseGain*
+AliForwardCorrectionManager::GetNoiseGain() const 
+{
+  /** 
+   * Get the noisegain calibration
+   * 
+   * @return NoiseGain calibration
+   */
+  return static_cast<const AliFMDCorrNoiseGain*>(Get(kIdNoiseGain)); 
+}
 
 //____________________________________________________________________
 const TAxis* 
index b840c87e17924e60cfe8a0fb911c3a557b411abd..87444175abedbbfe38057be3576c3fd37779937c 100644 (file)
@@ -22,6 +22,7 @@ class AliFMDCorrVertexBias;
 class AliFMDCorrMergingEfficiency;
 class AliFMDCorrAcceptance;
 class AliFMDCorrSecondaryMap;
+class AliFMDCorrNoiseGain;
 class TAxis;
 
 /**
@@ -46,7 +47,8 @@ private:
     kIdVertexBias,
     kIdMergingEfficiency,
     kIdDoubleHit,
-    kIdAcceptance
+    kIdAcceptance,
+    kIdNoiseGain
   };
 public:
   /**
@@ -59,13 +61,15 @@ public:
     kMergingEfficiency         = 0x08,
     kDoubleHit                 = 0x10,
     kAcceptance                = 0x20,
+    kNoiseGain                 = 0x40,
     kDefault                   = (kSecondaryMap|kELossFits|kAcceptance),
     kAll                       = (kSecondaryMap| 
                                  kELossFits|
                                  kVertexBias|
                                  kMergingEfficiency|
                                  kDoubleHit|
-                                 kAcceptance)
+                                 kAcceptance|
+                                 kNoiseGain)
   };
   /** 
    * Default constructor.  This is public for the sake of the ROOT I/O
@@ -247,17 +251,29 @@ public:
   /** 
    * Get the merging efficiency 
    * 
-   * 
    * @return Get the vertex efficiency correction 
    */
   const AliFMDCorrMergingEfficiency* GetMergingEfficiency() const;
   /** 
    * Get the acceptance correction due to dead channels 
    * 
-   * 
    * @return Acceptance correction due to dead channels 
    */
   const AliFMDCorrAcceptance* GetAcceptance() const;
+  /** 
+   * Get the noise calibration.  That is, the ratio 
+   *
+   * @f[ 
+   *   \frac{\sigma_{i}}{g_{i}k}
+   * @f] 
+   *
+   * where @f$ k@f$ is a constant determined by the electronics of
+   * units DAC/MIP, and @f$ \sigma_i, g_i@f$ are the noise and gain of
+   * the @f$ i @f$ strip respectively
+   * 
+   * @return Noise-gain calibration
+   */
+  const AliFMDCorrNoiseGain* GetNoiseGain() const;
 private:
   /** 
    * Non-default constructor - initializes corrections - used by
@@ -280,10 +296,11 @@ private:
   static const Char_t* fgkVertexBiasSkel;    // Name of correction object 
   static const Char_t* fgkMergingEffSkel;    // Name of correction object 
   static const Char_t* fgkAcceptanceSkel;    // Name of correction object 
+  static const Char_t* fgkNoiseGainSkel;     // Name of correction object 
   /* 
    * @} 
    */
-  ClassDef(AliForwardCorrectionManager,4) // Manager of corrections 
+  ClassDef(AliForwardCorrectionManager,5) // Manager of corrections 
 };
 
 #endif
index 670d3ec3dc3c350c8d4dcaac9b15284938f324c6..4901268d721559d7e94f3407812b5d245d70a6e4 100644 (file)
@@ -47,6 +47,7 @@ AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask()
     fMCRingSums(),
     fPrimary(0),
     fEventInspector(),
+    fESDFixer(),
     fSharingFilter(),
     fDensityCalculator(),
     fCorrections(),
@@ -68,6 +69,7 @@ AliForwardMCMultiplicityTask::AliForwardMCMultiplicityTask(const char* name)
     fMCRingSums(),
     fPrimary(0),
     fEventInspector("event"),
+    fESDFixer("esdFizer"),
     fSharingFilter("sharing"), 
     fDensityCalculator("density"),
     fCorrections("corrections"),
@@ -118,6 +120,9 @@ AliForwardMCMultiplicityTask::CreateBranches(AliAODHandler* ah)
 Bool_t
 AliForwardMCMultiplicityTask::Book()
 {
+  // We do this to explicitly disable the noise corrector for MC
+  GetESDFixer().SetRecoNoiseFactor(5);
+
   Bool_t ret = AliForwardMultiplicityBase::Book();
   POST(MCAOD_SLOT, &fMCAODFMD);
   POST(PRIM_SLOT,  fPrimary);
@@ -259,6 +264,9 @@ AliForwardMCMultiplicityTask::Event(AliESDEvent& esd)
 
   // Get FMD data 
   AliESDFMD*  esdFMD  = esd.GetFMDData();
+  
+  // Fix up the the ESD 
+  GetESDFixer().Fix(*esdFMD, ip.Z());
 
   // Apply the sharing filter (or hit merging or clustering if you like)
   if (isAccepted && !fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD,ip.Z())){
index a17b2e8b1b41a9e4e79a5defea60dc79ad10a53f..89a0e61b06fa211e372abb98aa11268c82ec4ad4 100644 (file)
@@ -15,6 +15,7 @@
  */
 #include "AliForwardMultiplicityBase.h"
 #include "AliFMDMCEventInspector.h"
+#include "AliFMDESDFixer.h"
 #include "AliFMDMCSharingFilter.h"
 #include "AliFMDMCDensityCalculator.h"
 #include "AliFMDMCCorrector.h"
@@ -114,6 +115,12 @@ public:
    * @return Reference to AliFMDEventInspector object 
    */
   AliFMDEventInspector& GetEventInspector() { return fEventInspector; }
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  AliFMDESDFixer& GetESDFixer() { return fESDFixer; }
   /**
    * Get reference to the SharingFilter algorithm 
    * 
@@ -144,6 +151,12 @@ public:
    * @return Reference to AliFMDEventInspector object 
    */
   const AliFMDEventInspector& GetEventInspector() const { return fEventInspector; }
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  const AliFMDESDFixer& GetESDFixer() const { return fESDFixer; }
   /**
    * Get reference to the SharingFilter algorithm 
    * 
@@ -227,6 +240,7 @@ protected:
   TH2D*                  fPrimary;      // Per event primary particles 
 
   AliFMDMCEventInspector    fEventInspector;    // Algorithm
+  AliFMDESDFixer            fESDFixer;          // Algorithm
   AliFMDMCSharingFilter     fSharingFilter;     // Algorithm
   AliFMDMCDensityCalculator fDensityCalculator; // Algorithm
   AliFMDMCCorrector         fCorrections;       // Algorithm
index ddf7cd5e713b128706e91ae3cbc80ff94149ba50..1dc038ada8ec0ddb88f44b1e3a40412b112a46cf 100644 (file)
@@ -21,6 +21,7 @@
 #include "AliInputEventHandler.h"
 #include "AliAnalysisManager.h"
 #include "AliFMDEventInspector.h"
+#include "AliFMDESDFixer.h"
 #include "AliFMDSharingFilter.h"
 #include "AliFMDDensityCalculator.h"
 #include "AliFMDCorrector.h"
@@ -93,6 +94,8 @@ AliForwardMultiplicityBase::Book()
   UInt_t what = AliForwardCorrectionManager::kAll;
   if (!fEnableLowFlux)
     what ^= AliForwardCorrectionManager::kDoubleHit;
+  if (!GetESDFixer().IsUseNoiseCorrection()) 
+    what ^= AliForwardCorrectionManager::kNoiseGain;
   if (!GetCorrections().IsUseVertexBias())
     what ^= AliForwardCorrectionManager::kVertexBias;
   if (!GetCorrections().IsUseAcceptance())
@@ -101,6 +104,7 @@ AliForwardMultiplicityBase::Book()
     what ^= AliForwardCorrectionManager::kMergingEfficiency;
   fNeededCorrections = what;
 
+  GetESDFixer()         .CreateOutputObjects(fList);
   GetSharingFilter()   .CreateOutputObjects(fList);
   GetDensityCalculator().CreateOutputObjects(fList);
   GetCorrections()     .CreateOutputObjects(fList);
@@ -541,6 +545,7 @@ AliForwardMultiplicityBase::Print(Option_t* option) const
   PFB("Store per-ring hists", fStorePerRing);
   PFB("Make timing histogram", fDoTiming);
   // gROOT->IncreaseDirLevel();
+  GetESDFixer()         .Print(option);        
   GetSharingFilter()    .Print(option);
   GetDensityCalculator().Print(option);
   GetCorrections()      .Print(option);
index 86417cb258f54ba426093d29268ec4566babd666..2fa4f54cbfab7e42683deca28f3247702ddba863 100644 (file)
@@ -19,6 +19,7 @@
 #include "AliAODForwardMult.h"
 #include "AliAODForwardEP.h"
 // class AliFMDEnergyFitter;
+class AliFMDESDFixer;
 class AliFMDSharingFilter;
 class AliFMDDensityCalculator;
 class AliFMDCorrector;
@@ -146,6 +147,12 @@ public:
    * @{ 
    * @name Access to sub-algorithms 
    */
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  virtual AliFMDESDFixer& GetESDFixer() = 0;
   /**
    * Get reference to the SharingFilter algorithm 
    * 
@@ -170,6 +177,12 @@ public:
    * @return Reference to AliFMDHistCollector object 
    */
   virtual AliFMDHistCollector& GetHistCollector() = 0;
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  virtual const AliFMDESDFixer& GetESDFixer() const = 0;
   /**
    * Get reference to the SharingFilter algorithm 
    * 
index 15209fb55ecb3687c27629e7585a7cff1cb9c00e..092aff7462dfbf87636ff27be58da3c38fd69412 100644 (file)
@@ -34,6 +34,7 @@ AliForwardMultiplicityTask::AliForwardMultiplicityTask()
   : AliForwardMultiplicityBase(),
     fESDFMD(),
     fEventInspector(),
+    fESDFixer(),
     fSharingFilter(),
     fDensityCalculator(),
     fCorrections(),
@@ -51,6 +52,7 @@ AliForwardMultiplicityTask::AliForwardMultiplicityTask(const char* name)
   : AliForwardMultiplicityBase(name),
     fESDFMD(),
     fEventInspector("event"),
+    fESDFixer("esdFizer"),
     fSharingFilter("sharing"), 
     fDensityCalculator("density"),
     fCorrections("corrections"),
@@ -67,6 +69,25 @@ AliForwardMultiplicityTask::AliForwardMultiplicityTask(const char* name)
 }
 
 
+//____________________________________________________________________
+void
+AliForwardMultiplicityTask::PreCorrections(const AliESDEvent* esd)
+{
+  if (!esd) return; 
+  
+  AliESDFMD* esdFMD = esd->GetFMDData();  
+  if (!esdFMD) return;
+
+  Int_t tgt = GetESDFixer().FindTargetNoiseFactor(*esdFMD, false);
+  if (tgt <= 0) {
+    // If the target noise factor is 0 or less, disable the noise/gain
+    // correction.
+    GetESDFixer().SetRecoNoiseFactor(4);
+    fNeededCorrections ^= AliForwardCorrectionManager::kNoiseGain;
+  }
+  else 
+    AliWarning("The noise corrector has been enabled!");
+}
 //____________________________________________________________________
 Bool_t
 AliForwardMultiplicityTask::PreEvent()
@@ -131,6 +152,9 @@ AliForwardMultiplicityTask::Event(AliESDEvent& esd)
   // Get FMD data 
   AliESDFMD* esdFMD = esd.GetFMDData();  
 
+  // Fix up the the ESD 
+  GetESDFixer().Fix(*esdFMD, ip.Z());
+
   // Apply the sharing filter (or hit merging or clustering if you like)
   if (fDoTiming) individual.Start(true);
   if (!fSharingFilter.Filter(*esdFMD, lowFlux, fESDFMD, ip.Z())) { 
index a52cca5804bfa96a32ac8b0cf4d6027a3df65dd5..338efbc43ab445209bd648de9077ce288a988773 100644 (file)
@@ -16,6 +16,7 @@
 #include "AliForwardMultiplicityBase.h"
 #include "AliForwardUtil.h"
 #include "AliFMDEventInspector.h"
+#include "AliFMDESDFixer.h"
 #include "AliFMDSharingFilter.h"
 #include "AliFMDDensityCalculator.h"
 #include "AliFMDCorrector.h"
@@ -62,6 +63,14 @@ public:
    * @{ 
    * @name Interface methods 
    */
+  /** 
+   * Called on first event _before_ reading corrections.  Here, the
+   * user class can do additional checking to see if the some (more or
+   * less) corrections are needed.
+   * 
+   * @param esd Event 
+   */
+  virtual void PreCorrections(const AliESDEvent* esd);
   /** 
    * Called before processing a single event - should not do anything
    * but clear data, etc.
@@ -90,6 +99,12 @@ public:
    * @return Reference to AliFMDEventInspector object 
    */
   AliFMDEventInspector& GetEventInspector() { return fEventInspector; }
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  AliFMDESDFixer& GetESDFixer() { return fESDFixer; }
   /**
    * Get reference to the SharingFilter algorithm 
    * 
@@ -120,6 +135,12 @@ public:
    * @return Reference to AliFMDEventInspector object 
    */
   const AliFMDEventInspector& GetEventInspector() const { return fEventInspector; }
+  /**
+   * Get reference to the ESDFixer algorithm 
+   * 
+   * @return Reference to AliFMDESDFixer object 
+   */
+  const AliFMDESDFixer& GetESDFixer() const { return fESDFixer; }
   /**
    * Get reference to the SharingFilter algorithm 
    * 
@@ -177,6 +198,7 @@ protected:
 
   AliESDFMD               fESDFMD;            // Sharing corrected ESD object
   AliFMDEventInspector    fEventInspector;    // Algorithm
+  AliFMDESDFixer          fESDFixer;          // Algorithm
   AliFMDSharingFilter     fSharingFilter;     // Algorithm
   AliFMDDensityCalculator fDensityCalculator; // Algorithm
   AliFMDCorrector         fCorrections;       // Algorithm
index 9f007e35e3257a0959c07d699e659a37ee986937..890005c7b824ddfdec55e35e71d5a86367953e3d 100644 (file)
@@ -1,6 +1,6 @@
 /**
  * @file   DrawdNdeta.C
- * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk>
+ * @author Christian Holm Christensen <cholm@nbi.dk>
  * @date   Wed Mar 23 14:07:10 2011
  * 
  * @brief  Script to visualise the dN/deta for pp and PbPb
@@ -780,42 +780,51 @@ struct dNdetaDrawer
                  TObjArray&    truths)
   {
     if (!list) return 0;
-    UShort_t   n = HasCent() ? fCentAxis->GetNbins() : 0;
-    // Info("FetchTopResults","got %d centrality bins", n);
-    if (n == 0) {
-      TH1*  h = FetchOne(list, mcList, empCorr, name, "all",
-                        FetchOthers(0,0), -1000, 0, 
-                        max, rmax, amax, fTruth);
-      if (!h) return 0;
-      TObjArray* a = new TObjArray;
-      // Info("FetchTopResults", "Adding %s to result stack", h->GetName());
-      a->AddAt(h, 0);
-      return a;
-    }
+    UShort_t   n = HasCent() ? fCentAxis->GetNbins() : 1;
+    // // Info("FetchTopResults","got %d centrality bins", n);
+    // if (n == 0) {
+    //   TH1* h  = FetchOne(list, mcList, empCorr, name, "all",
+    //                          FetchOthers(0,0), -1000, 0, 
+    //                          max, rmax, amax, truths);
+    //   if (!h) return 0;
+    //   TObjArray* a = new TObjArray;
+    //   // Info("FetchTopResults", "Adding %s to result stack", h->GetName());
+    //   a->AddAt(h, 0);
+    //   return a;
+    // }
     
     TObjArray* a = new TObjArray;
     truths.Expand(n);
     for (UShort_t i = 0; i < n; i++) { 
-      UShort_t centLow  = fCentAxis->GetBinLowEdge(i+1);
-      UShort_t centHigh = fCentAxis->GetBinUpEdge(i+1);
-      TString  lname    = Form("cent%03d_%03d", centLow, centHigh);
-      Int_t    col      = GetCentralityColor(i+1);
-      TString  centTxt  = Form("%3d%%-%3d%% central", centLow, centHigh);
-
+      UShort_t centLow  = 0;
+      UShort_t centHigh = 0;
+      TString  lname    = "all";
+      Int_t    col      = -1000;
+      TString  centTxt  = "";
+      if (HasCent()) {
+       centLow  = fCentAxis->GetBinLowEdge(i+1);
+       centHigh = fCentAxis->GetBinUpEdge(i+1);
+       lname    = Form("cent%03d_%03d", centLow, centHigh);
+       col      = GetCentralityColor(i+1);
+        centTxt  = Form("%3d%%-%3d%% central", centLow, centHigh);
+      }
       TH1* tt = static_cast<TH1*>(truths.At(i));
       TH1* ot = tt;
       TH1* h  = FetchOne(list, mcList, empCorr, name, lname,
                         FetchOthers(centLow,centHigh), col, 
-                        centTxt.Data(), max, rmax, amax, fTruth);
+                        centTxt.Data(), max, rmax, amax, tt);
       if (!h) continue;
-      if (ot != tt) { 
-       //Info("FetchTopResults", "old truth=%p new truth=%p (%s)", 
-       //     ot, tt, name);
+
+      if (tt != ot) { 
        truths.AddAt(tt, i);
       }
       // Info("FetchTopResults", "Adding %p to result stack", h);
       a->AddAt(h, i);
     }
+    if (a->GetEntries() <= 0) {
+      delete a;
+      a = 0;
+    }
     return a;
   } 
   //__________________________________________________________________
@@ -908,6 +917,8 @@ struct dNdetaDrawer
     TH1* dndeta      = FetchHistogram(list, Form("dndeta%s", name));
     TH1* dndetaMC    = FetchHistogram(list, Form("dndeta%sMC", name));
     TH1* dndetaTruth = FetchHistogram(list, "dndetaTruth");
+    Info("", "dN/deta truth from %s: %p", list->GetName(), dndetaTruth);
+    Info("", "dN/deta truth from external: %p", truth);
 
     if (mcList && FetchHistogram(mcList, "finalMCCorr")) 
       Warning("FetchCentResults", "dNdeta already corrected for final MC");
@@ -962,7 +973,7 @@ struct dNdetaDrawer
          while ((hist = static_cast<TH1*>(next()))) 
            max = TMath::Max(max, AddHistogram(fResults, hist));
        }
-      }
+      } // If show rings
       // Info("FetchCentResults", "Got %p, %p, %p from %s with name %s, max=%f",
       //      dndeta, dndetaMC, dndetaTruth, list->GetName(), name, max);
       
@@ -986,13 +997,16 @@ struct dNdetaDrawer
          }
        }
        // fOthers->Add(thisOther);
-      }
-    }
+      } // if others for this 
+    } // if not truth 
     if (dndetaMC) { 
       fRatios->Add(Ratio(dndeta,    dndetaMC,    rmax));
       fRatios->Add(Ratio(dndetaSym, dndetaMCSym, rmax));
     }
     if (truth) {
+      Info("", "Forming ratio to truth:\n\t%s\n\t%s",
+          dndeta->GetName(), 
+          truth->GetName());
       fRatios->Add(Ratio(dndeta,      truth, rmax));
       fRatios->Add(Ratio(dndetaSym,   truth, rmax));
     }
index 0588d0fb0095dc3efed166ef898f840d7c5a83f4..f75193d9dfa837da41d9ce608baa5fda74fc7711 100644 (file)
@@ -43,8 +43,8 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
   // pedestal value in them, so the absolute number of high-value
   // pedestal signals is large - despite the low probablity).
   AliFMDMultCuts cSharingLow(AliFMDMultCuts::kFixed,0.15);
-  AliFMDMultCuts cSharingHigh(AliFMDMultCuts::kMPVFraction,0.8,0.8,0.8,0.8,0.8);
-  AliFMDMultCuts cDensity(AliFMDMultCuts::kMPVFraction,0.7);
+  AliFMDMultCuts cSharingHigh(AliFMDMultCuts::kLandauSigmaWidth,1);
+  AliFMDMultCuts cDensity(AliFMDMultCuts::kLandauSigmaWidth,1);
   
   // --- Event inspector ---------------------------------------------
   // Set the number of SPD tracklets for which we consider the event a
@@ -58,39 +58,26 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
   task->GetEventInspector().SetMinPileupDistance(.8);
   // V0-AND triggered events flagged as NSD 
   task->GetEventInspector().SetUseV0AndForNSD(false);
-  // Use primary vertex selection from 1st physics WG
-  task->GetEventInspector().SetUseFirstPhysicsVtx(false);
-  // Use satellite collisions
-  task->GetEventInspector().SetUseDisplacedVertices(false);
+  // Set the kind of vertex to look for.  Can be one of 
+  //  
+  //   - kNormal:    SPD vertex 
+  //   - kpA2012:    Selection tuned for 2012 pA data 
+  //   - kpA2013:    Selection tuned for 2013 pA data 
+  //   - kPWGUD:     Selection used by 'first physics' 
+  //   - kDisplaced: Satellite collisions, with kNormal fall-back 
+  // 
+  task->GetEventInspector().SetVertexMethod(AliFMDEventInspector::kNormal);
   // Which centrality estimator to use 
   task->GetEventInspector().SetCentralityMethod("V0M");
 
-  // --- Sharing filter ----------------------------------------------
-  // If the following is set to true, then the merging of shared
-  // signals is disabled completely
-  task->GetSharingFilter().SetDisableMerging(false);
-  // Enable use of angle corrected signals in the algorithm 
-  task->GetSharingFilter().SetUseAngleCorrectedSignals(true);
-  // Ignore the ESD information when angle correcting.
-  // 
-  // *IMPORTANT* 
-  // 
-  // This is to counter a known issue with AliESDFMD with ClassDef 
-  // version < 4, where the angle correction flag is incorrectly set.
-  // A fix is coming to AliESDFMD to handle it directly in the class. 
-  // Only set the flag below to true if you know it to be necessary for
-  // your data set.
-  task->GetSharingFilter().SetIgnoreESDWhenAngleCorrecting(false);
-  // Disable use of angle corrected signals in the algorithm 
-  task->GetSharingFilter().SetZeroSharedHitsBelowThreshold(false);
-  // Whether to use simple merging algorithm
-  task->GetSharingFilter().SetUseSimpleSharing(true);
-  // Whether to allow for 3 strip hits - deprecated
-  task->GetSharingFilter().SetAllow3Strips(false);
-  // Set upper sharing cut 
-  task->GetSharingFilter().SetHCuts(cSharingHigh);
-  // Enable use of angle corrected signals in the algorithm 
-  task->GetSharingFilter().SetLCuts(cSharingLow);
+  // --- ESD fixer ---------------------------------------------------
+  // Sets the noise factor that was used during reconstruction.  If
+  // this is set to 4 or more, then this correction will be disabled.
+  task->GetESDFixer().SetRecoNoiseFactor(1);
+  // IF the noise correction is bigger than this, flag strip as dead 
+  task->GetESDFixer().SetMaxNoiseCorrection(0.05);
+  // Sets whether to recalculate eta 
+  task->GetESDFixer().SetRecalculateEta(false);
   // If true, consider AliESDFMD::kInvalidMult as a zero signal.  This
   // has the unfortunate side effect, that we cannot use the
   // on-the-fly calculation of the phi acceptance.  
@@ -117,9 +104,9 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
   // LHC10c-7TeV is effected up-to and including pass2
   // LHC10c-CPass0 should be OK, but has limited statistics 
   // LHC10c_11a_FMD should be OK, but has few runs  
-  task->GetSharingFilter().SetInvalidIsEmpty(false);
+  task->GetESDFixer().SetInvalidIsEmpty(false);
   // Dead region in FMD2i
-  task->GetSharingFilter().AddDeadRegion(2, 'I', 16, 17, 256, 511);  
+  task->GetESDFixer().AddDeadRegion(2, 'I', 16, 17, 256, 511);  
   // One can add extra dead strips from a script like 
   // 
   //   void deadstrips(AliFMDSharingFilter* filter)
@@ -130,7 +117,34 @@ ForwardAODConfig(AliForwardMultiplicityBase* task)
   //
   // and then do here 
   // 
-  // task->GetSharingFilter().AddDead("deadstrips.C");
+  // task->GetESDFixer().AddDead("deadstrips.C");
+
+  // --- Sharing filter ----------------------------------------------
+  // If the following is set to true, then the merging of shared
+  // signals is disabled completely
+  // task->GetSharingFilter().SetMergingDisabled(false);
+  // Enable use of angle corrected signals in the algorithm 
+  task->GetSharingFilter().SetUseAngleCorrectedSignals(true);
+  // Ignore the ESD information when angle correcting.
+  // 
+  // *IMPORTANT* 
+  // 
+  // This is to counter a known issue with AliESDFMD with ClassDef 
+  // version < 4, where the angle correction flag is incorrectly set.
+  // A fix is coming to AliESDFMD to handle it directly in the class. 
+  // Only set the flag below to true if you know it to be necessary for
+  // your data set.
+  task->GetSharingFilter().SetIgnoreESDWhenAngleCorrecting(false);
+  // Disable use of angle corrected signals in the algorithm 
+  task->GetSharingFilter().SetZeroSharedHitsBelowThreshold(false);
+  // Whether to use simple merging algorithm
+  task->GetSharingFilter().SetUseSimpleSharing(true);
+  // Whether to allow for 3 strip hits - deprecated
+  task->GetSharingFilter().SetAllow3Strips(false);
+  // Set upper sharing cut 
+  task->GetSharingFilter().SetHCuts(cSharingHigh);
+  // Enable use of angle corrected signals in the algorithm 
+  task->GetSharingFilter().SetLCuts(cSharingLow);
    
   // --- Density calculator ------------------------------------------
   // Set the maximum number of particle to try to reconstruct 
diff --git a/PWGLF/FORWARD/analysis2/corrs/Browse.C b/PWGLF/FORWARD/analysis2/corrs/Browse.C
new file mode 100644 (file)
index 0000000..6bf58d8
--- /dev/null
@@ -0,0 +1,15 @@
+       TObject* Browse()
+{
+  const char* fwd = "/opt/alice/aliroot/trunk-inst/PWGLF/FORWARD/analysis2";
+  if (!gROOT->GetClass("AliOADBForward"))
+    gROOT->Macro(Form("%s/scripts/LoadLibs.C", fwd));
+  gSystem->AddIncludePath(Form("-I%s -I$ALICE_ROOT/include", fwd));
+  gROOT->LoadMacro(Form("%s/corrs/ForwardOADBGui.C", fwd));
+  
+  AliOADBForward* db = new AliOADBForward;
+  db->Open("fmd_corrections.root", "*");
+  
+  ForwardOADBGui(db);
+
+  return db;
+}
index 7d529a80b9c463bb7d818102f1b7a42d33e13654..90ae76d6c93a1a0504ab0597f9f9d982c7296408 100644 (file)
@@ -3,6 +3,7 @@
 # include "AliFMDCorrAcceptance.h"
 # include "AliFMDCorrSecondaryMap.h"
 # include "AliFMDCorrELossFit.h"
+# include "AliFMDCorrNoiseGain.h"
 # include "AliForwardUtil.h"
 # include "AliForwardCorrectionManager.h"
 # include "AliLog.h"
@@ -15,6 +16,7 @@ class TObject;
 class AliFMDCorrAcceptance;
 class AliFMDCorrSecondaryMap;
 class AliFMDCorrELossFit;
+class AliFMDCorrNoiseGain;
 #include <TString.h>
 #endif
 
@@ -28,7 +30,7 @@ public:
    */
   CorrDrawer() 
   {
-    fELossExtra = "forward_eloss.root";
+    fELossExtra = ""; // forward_eloss.root";
     fMinQuality = 8;
   }
   /** 
@@ -154,6 +156,8 @@ public:
       AppendName(name, AliForwardCorrectionManager::kAcceptance);
     if (what & AliForwardCorrectionManager::kELossFits) 
       AppendName(name, AliForwardCorrectionManager::kELossFits);
+    if (what & AliForwardCorrectionManager::kNoiseGain) 
+      AppendName(name, AliForwardCorrectionManager::kNoiseGain);
     if (what & AliForwardCorrectionManager::kVertexBias) 
       Warning("CorrDrawer","Vertex bias not implemented yet");
     if (what & AliForwardCorrectionManager::kDoubleHit) 
@@ -164,7 +168,8 @@ public:
     // Filter the ones we can handle 
     UShort_t flags = what & (AliForwardCorrectionManager::kELossFits|
                             AliForwardCorrectionManager::kAcceptance|
-                            AliForwardCorrectionManager::kSecondaryMap);
+                            AliForwardCorrectionManager::kSecondaryMap|
+                            AliForwardCorrectionManager::kNoiseGain);
     if (!mgr.Init(runNo, sys, sNN, field, mc, sat, flags, true)) {
       Error("CorrDrawer", "Failed to initialize for flags=0x%02x"
                "run=%lu, sys=%hu, sNN=%hu, field=%hd, mc=%d, sat=%d",
@@ -200,6 +205,8 @@ public:
       DrawIt(mgr.GetAcceptance(), details);
     if (what & AliForwardCorrectionManager::kELossFits)
       DrawIt(mgr.GetELossFit(), details, few);
+    if (what & AliForwardCorrectionManager::kNoiseGain)
+      DrawIt(mgr.GetNoiseGain(), details);
 
     // Done
     CloseCanvas();
@@ -326,6 +333,19 @@ public:
     DrawIt(sec, pdf); 
     if (pdf) CloseCanvas();
   }
+  /** 
+   * Draw a single summary plot multiple plots of the energy loss
+   * fits.  A new canvas is created for this.
+   * 
+   * @param fits Energy loss fits
+   * @param pdf  If true, do multiple plots. Otherwise a single summary plot
+   */
+  virtual void Summarize(const AliFMDCorrNoiseGain* corr, Bool_t pdf=true) 
+  { 
+    CreateCanvas(CanvasName("noisegain.pdf"), true, pdf);
+    DrawIt(corr, pdf); 
+    if (pdf) CloseCanvas();
+  }
   /** 
    * Draw a single summary plot multiple plots of the energy loss
    * fits.  A new canvas is created for this.
@@ -351,10 +371,10 @@ public:
    *
    * @deprecated Use Run instead 
    */
-  static void Summarize(const TString& what   = TString("")
-                       Bool_t         /*mc*/ = false,
-                       const TString& output = TString("forward_eloss.root"), 
-                       const TString& local  = TString("fmd_corrections.root"),
+  static void Summarize(const TString& what   = ""
+                       Bool_t         mc     = false,
+                       const TString& output = "",
+                       const TString& local  = "fmd_corrections.root",
                        Option_t*      options= "")
   {
     CorrDrawer* drawer = new CorrDrawer;
@@ -374,8 +394,8 @@ public:
    * @deprecated Use Run instead 
    */
   static void Summarize(UShort_t       what, 
-                       Bool_t         /*mc*/ = false,
-                       const TString& output = "forward_eloss.root", 
+                       Bool_t         mc = false,
+                       const TString& output = "",
                        const TString& local  = "fmd_corrections.root",
                        Option_t*      options= "")
   {
@@ -399,6 +419,8 @@ protected:
       what.Append("acceptance"); break;
     case AliForwardCorrectionManager::kELossFits:                
       what.Append("elossfits"); break;
+    case AliForwardCorrectionManager::kNoiseGain:                
+      what.Append("noisegain"); break;
     default:
       what.Append("unknown"); break;
     }
@@ -701,6 +723,46 @@ protected:
       PrintCanvas("#LTsecondary map#GT");
     }
   }
+  virtual void DrawIt(const AliFMDCorrNoiseGain* corr, bool details) 
+  {
+    if (!corr) {
+      Warning("CorrDrawer","No noise-gain correction available");
+      return;
+    }
+
+    if (!fCanvas) {
+      Warning("CorrDrawer", "No canvas");
+      return;
+    }
+
+    DivideForRings(false,false);
+    
+    for (UShort_t d = 1; d <= 3; d++) { 
+      UShort_t nQ = d == 1 ? 1 : 2;
+      for (UShort_t q = 0; q < nQ; q++) { 
+       Char_t   r  = q == 0 ? 'I' : 'O';
+       UShort_t nS = q == 0 ?  20 :  40;
+       UShort_t nT = q == 0 ? 512 : 256;
+       
+       TH2* h = new TH2D(Form("fmd%d%c", d, r),
+                         Form("FMD%d%c", d, r), 
+                         nT, -.5, nT-.5,  nS, -.5, nS-.5);
+       h->SetDirectory(0);
+       h->SetXTitle("Strip");
+       h->SetYTitle("Sector");
+
+       for (UShort_t s = 0; s < nS; s++) { 
+         for (UShort_t t = 0; t < nT; t++) { 
+           Float_t c = corr->Get(d,r,s,t);
+           h->Fill(t,s,c);
+         }
+       }
+       h->GetZaxis()->SetRangeUser(0,0.05);
+       DrawInRingPad(d,r,h,"COLZ");
+      }
+    }
+    PrintCanvas("Noise correction");
+  }
   /** 
    * Draw the energy loss fits correction 
    * 
index d1045f8c378293b8d36d6e302a0d2f0ad682d2a2..f040909f48b23bdd1f57ca73166052bc0f95b3c1 100644 (file)
@@ -434,7 +434,8 @@ struct ForwardOADBGUI
     TDatime dt(e->fTimestamp);
     lve->SetSubnames(Form("%lu", e->fRunNo), 
                     (e->fSys == 1 ? "p-p" : 
-                     e->fSys == 2 ? "Pb-Pb" : "p-Pb"),
+                     e->fSys == 2 ? "Pb-Pb" : 
+                     e->fSys == 3 ? "p-Pb" : "?"),
                     Form("%4huGeV",e->fSNN), 
                     Form("%+2hdkG", e->fField), 
                     (e->fMC ? "MC" : "Real"),
index 20a7073efebf3dca7e3dac0852135fcfb0dbf02e..8ea2835cfa77273882b557d945222fee76d2cd47 100644 (file)
@@ -1,5 +1,5 @@
 void
-testF(UShort_t single=0, Bool_t old=false)
+RunTestF(UShort_t single=0, Bool_t old=false)
 {
   gSystem->AddIncludePath(Form("%s-DTEST_SHIFT -DTEST_FITTER -I$ANA_SRC -I.", 
                               (old ? "-DNO_SIGMA_SHIFT " : "")));
diff --git a/PWGLF/FORWARD/analysis2/corrs/Upload.C b/PWGLF/FORWARD/analysis2/corrs/Upload.C
new file mode 100644 (file)
index 0000000..ba8c27f
--- /dev/null
@@ -0,0 +1,22 @@
+// Generated by ExtractELoss.C
+TString MakeDest(const TString& dest, const TString& fname)
+{
+  TString tmp(dest);
+  if (!tmp.IsNull()) {
+    if (!tmp.EndsWith("/")) tmp.Append("/");
+    tmp.Append(fname);
+  }
+  return tmp;
+}
+
+void Upload(const TString& dest="")
+{
+  gROOT->Macro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/LoadLibs.C");
+  
+  const char* fmdFile = "fmd_corrections.root";
+  TString fdest = MakeDest(dest, fmdFile);
+  
+  AliForwardCorrectionManager::Instance().Append(fmdFile, fdest);
+}
+// EOF
+
index f7c3a863e1bacde1001f0768de6047063868a951..883ea805b4107b4d089b77c04a6af372b798834b 100644 (file)
@@ -204,7 +204,7 @@ struct QAPlotter : public QABase
    */
   QAPlotter(Long_t prodYear, Char_t prodLetter, Bool_t useVar) 
     : QABase("", (prodYear < 2000 ? 2000 : 0) + prodYear, 
-            Form("LHC%02d%c", prodYear % 100, prodLetter, "pass0")), 
+            Form("LHC%02d%c", int(prodYear % 100), prodLetter), "pass0"), 
       fNAccepted(0),
       fVz(0), 
       fUseVar(useVar)
index 583e9fb61c7e2c986fc4454c082c0889b40d2e26..3fce9f2e6cd20855045e3624a3689363b2211b7c 100644 (file)
@@ -558,7 +558,7 @@ public:
            Int_t  prodYear, 
            char   prodLetter) 
     : QABase("data", (prodYear < 2000 ? 2000 : 0) + prodYear,
-            Form("LHC%02d%c", (prodYear%2000), prodLetter), "pass0")
+            Form("LHC%02d%c", (prodYear%2000), prodLetter), "pass0"),
       fRunNo(-1),
       fCurrentFile(0),
       fSharingFilter(0),
index 038e2c6e79223d95e58448cfc6f559b518214cdf..5395c777390475228230ba64e7c262f662d4b474 100755 (executable)
@@ -222,6 +222,8 @@ files=
 path=
 numf=0
 from_local=0
+type=data
+passid=
 get_filelist()
 {
     mess 3 "Getting file list" 
@@ -240,6 +242,8 @@ get_filelist()
        x) ;; 
        *)  mess 3 "Assuming simulation output"
            datd=sim 
+           type=sim
+           passfull=passMC
            esdd= 
            ;; 
     esac
@@ -254,6 +258,7 @@ get_filelist()
        passpre=
        post=
     fi
+    passid=${paid}
     local post=${passpost}
     case x$post in 
        x_*) ;; 
@@ -332,7 +337,9 @@ analyse_file()
 {
     local dir=`dirname $1` 
     local inp=`basename $1` 
-    local out=`echo $inp | sed 's/trending_/tree_/'` 
+    local r=$2
+    local out=trending.root 
+    # `echo $inp | sed 's/trending_/tree_/'` 
     local ret=0
     mess 2 -n "Analysing $inp -> $out ... "
 
@@ -344,16 +351,24 @@ analyse_file()
        rm -f $dir/$out
     fi
 
+    
+    mess 1 "Running runQA.sh '$inp' '$type' '$prodyear' '$prodfull' '$passid' '$r' in '$dir'"
     (cd $dir 
-       root -l -b  <<EOF 
-.L RunFileQA.C
-RunFileQA("$inp", "$out", $prodyear, "$prodletter");
-.q
-EOF
+       $ALICE_ROOT/PWGLF/FORWARD/analysis2/qa/runQA.sh \
+           "$inp" "$type" $prodyear "$prodfull" "$passid" "$r" 
        ret=$? 
-       mess 2 " -> $ret"
-       # rm -f ${scr}.C 
-    ) 2>> $redir
+       mess 2 " -> $ret")
+#     (cd $dir 
+#      root -l -b  <<EOF 
+# .L RunFileQA.C
+# RunFileQA("$inp", "$out", $prodyear, "$prodletter");
+# .q
+# EOF
+#      ret=$? 
+#      mess 2 " -> $ret"
+#      # rm -f ${scr}.C 
+#     ) 2>> $redir
+
     return $ret
 }
 
@@ -364,7 +379,9 @@ analyse_run()
     local source=$1 ; shift 
     local store=$1 ; shift 
     local r=$1 ; shift 
-    local o=${store}/`basename $file .root`_${r}.root 
+    local rr=`printf %09d $r`
+    local o=${store}/${rr}/input.root
+    mkdir -p ${store}/${rr}
 
     mess 2 -n "$source -> $o ... "
     if test -f $o && test $force -lt 2; then 
@@ -399,11 +416,11 @@ analyse_run()
     check_file ${o} 
     local ret=$? 
     case $ret in 
-       0|2) ;; 
+       0|2) ;; 
        1|3|4|5|6) return 2 ;; 
     esac
 
-    analyse_file ${o}
+    analyse_file ${o} ${r}
 
     return 0
 }
@@ -510,44 +527,12 @@ make_trend()
     local ret=0
     mess 1 -n "Analysing for trend $dir ... "
     (cd $dir 
-       rm -f trend_*_*.html 
-       rm -f trend_*_*.pdf
-       rm -f trend_*_*.root
-
-       root -l -b <<EOF 
-.L RunFinalQA.C
-RunFinalQA(".", $prodyear, "$prodletter", $variance);
-.q
-EOF
-       mess 1 -n " ... "
-       # mess 3 -n "root -l -b -q ${scr}.C "
-       # root -l -b -q ${scr}.C  > /dev/null 2>&1 
-       # local ret=$? 
-       # mess 1 " -> $ret"
-       # rm -f ${scr}.C 
-
-       # do the index file 
-       local idx=`ls trend_*_*.html 2> /dev/null` 
-       for i in $idx ; do 
-           mess 1 "Making index.html point to $i" 
-           sed -e 's,index.html,../index.html,' \
-               -e "s,!--JOBID--,a target='_blank' href='${jobu}${jobid}'>Job</a," \
-               < $i > index.html 
-           cp index.html $i
-       done
-       
-       if test ! -f index.html ; then 
-           echo "No index file found" 
-           ret=1
-       else 
-           fix_perm index.html 
-           fix_perm . > /dev/null 2>&1 
-       fi
-
-       copy_style
-       copy_aliroot_file $favicon
-       copy_aliroot_file $logo
-    ) 2>> $redir
+       echo "hadd trending.root 000*/trending.root"
+       rm -f trending.root 
+       hadd -k trending.root 000*/trending.root 
+       $ALICE_ROOT/PWGLF/FORWARD/analysis2/qa/periodQA.sh trending.root 
+       ret=$?
+    )
     return $ret
 }
 
index 19377ce430e9a17891608440efcaf7b104db93b8..6c45a89aa5236d7bc61a6929074fb1892c0dd615 100755 (executable)
@@ -10,7 +10,7 @@ fwd=$ANA_SRC
 # 6: run
 aliroot -l -b -q ${fwd}/qa/RunQA.C\(\"$1\",\"$2\",$3,\"$4\",\"$5\",$6\)
 cp ${fwd}/qa/style.css .
-cp ${fwd}/qa/script.css .
+cp ${fwd}/qa/script.js .
 cp ${fwd}/qa/fmd_favicon.png . 
 cp ${fwd}/qa/fmd_logo.png .
 
index ec5fd999405121465d7ac72d5a34380d1f35d41f..94fcb71c68401f9758f6133d17c91383cb37b649 100644 (file)
@@ -338,6 +338,8 @@ struct Trend : public SummaryDrawer
     TCollection* results = GetCollection(file, "ForwarddNdetaResults");
     if (!results) return false;
 
+    TCollection* mcResults = GetCollection(file, "MCTruthdNdetaResults");
+
     THStack* all = new THStack("all",    "All");
     THStack* rat = new THStack("ratios", "Ratios");
     
@@ -345,18 +347,24 @@ struct Trend : public SummaryDrawer
     TObject* pCent = 0;
     Int_t    jCent = 0;
     while ((pCent = iCent())) {
-      TGraph*  graph = static_cast<TGraph*>(fOthers.At(jCent));
-      if (!graph) break;
-
       TString folderName(pCent->GetName());
       TCollection* centFolder = GetCollection(results, folderName);
-      if (!centFolder) break;
+      if (!centFolder) {
+       Warning("", "Didn't get the centrality %s folder from %s",
+               folderName.Data(), results->GetName());
+       // results->ls();
+       break;
+      }
 
+      TCollection* mcCentFolder = 0;
+      if (mcResults) mcCentFolder = GetCollection(mcResults, folderName);
+      
       TH1* dNdeta = GetH1(centFolder, Form("dndetaForward%s", 
                                           fRebinned ? "_rebin05" : ""));
       if (!dNdeta) {
        Warning("", "Didn't get histogram for jCent=%d in %s", 
                jCent, path.Data());
+       // results->ls();
        break;
       }
       dNdeta->SetDirectory(out);
@@ -364,7 +372,22 @@ struct Trend : public SummaryDrawer
       dNdeta->SetMarkerColor(jCent+1);
       dNdeta->SetLineColor(jCent+1);
 
-      TH1* other = G2H(graph, *(dNdeta->GetXaxis()));
+      TH1* other = 0; 
+      if (mcCentFolder) 
+       other = GetH1(mcCentFolder, Form("dndetaMCTruth%s", 
+                                        fRebinned ? "_rebin05" : ""));
+      if (!other) {
+       TGraph*  graph = static_cast<TGraph*>(fOthers.At(jCent));
+       if (!graph) break;
+
+       other = G2H(graph, *(dNdeta->GetXaxis()));
+      }
+
+      if (!other) {
+       Warning("", "No other data found for %s", path.Data());
+       break;
+      }
+
       other->SetMarkerColor(dNdeta->GetMarkerColor());
       other->SetMarkerSize(dNdeta->GetMarkerSize());
       other->SetLineColor(dNdeta->GetLineColor());
@@ -600,6 +623,10 @@ struct Trend : public SummaryDrawer
       Warning("DrawStacks", "Stack is missing!");
       return;
     }
+    if (!stack->GetHists() || stack->GetHists()->GetEntries() <= 0) { 
+      Warning("DrawStacks", "Stack is empty");
+      return;
+    }
     TH1*    first = static_cast<TH1*>(stack->GetHists()->At(0));
     TString xT    = first->GetXaxis()->GetTitle();
     
index 1988169a031dc4daeca5f76d84edf2796c6b2eeb..d644d0a832f0e0ebe80aed6a34bb229aa1230749 100644 (file)
@@ -31,6 +31,7 @@ public:
     fOptions.Add("cent", "Use centrality");
     fOptions.Add("only-mb", "Only collect statistics from MB events");
     fOptions.Add("residuals", "MODE", "Optional calculation of residuals", "");
+    fOptions.Add("corr", "DIR", "Corrections dir", "");
     fOptions.Set("type", "ESD");
   }
 protected:
@@ -71,12 +72,21 @@ protected:
     Bool_t   cent   = fOptions.Has("cent");
     Bool_t   onlyMB = fOptions.AsBool("only-mb");
     Int_t    verb   = fOptions.AsInt("verbose");
+    TString  corrs  = "";
+    if (fOptions.Has("corr")) {
+      corrs = fOptions.Get("corr"); 
+    }
     TString  resi   = "";
     if (fOptions.Has("residuals")) resi = fOptions.Get("residuals"); 
 
     // --- Add the task ----------------------------------------------
-    AddTask("AddTaskFMDELoss.C", Form("%d,%d,%d,%d,\"%s\"", 
-                                     mc, cent, onlyMB, verb, resi.Data()));
+    AddTask("AddTaskFMDELoss.C", Form("%d,%d,%d,%d,\"%s\",\"%s\"", 
+                                     mc, cent, onlyMB, verb, 
+                                     resi.Data(), corrs.Data()));
+
+    if (!corrs.IsNull()) {
+      fHelper->LoadAux(Form("%s/fmd_corrections.root",corrs.Data()), true);
+    }
   }
   /** 
    * Create entrality selection if enabled 
@@ -167,7 +177,7 @@ protected:
     f << std::boolalpha 
       << "// Generated by " << ClassName() << "\n"
       << "// If force=true, then force set parameters\n"
-      << "// If shift=true, enable extra shift in Delta from \sigma/xi
+      << "// If shift=true, enable extra shift in Delta from sigma/xi\n"
       << "//\n"
       << "void ReFit(Bool_t      force=false,\n"
       << "           Bool_t      shift=true,\n"
index 6e261e70fa48ed294cf96fc8327c13ec2961bf5b..bbd9d9d7cb6cf90cb4175171db1518f1e5e7840a 100644 (file)
 #pragma link C++ class AliFMDCorrMergingEfficiency+;
 #pragma link C++ class AliFMDCorrSecondaryMap+;
 #pragma link C++ class AliFMDCorrVertexBias+;
+#pragma link C++ class AliFMDCorrNoiseGain+;
 
 // FMD algorithms 
 #pragma link C++ class AliFMDDensityCalculator+;
 #pragma link C++ class AliFMDEventInspector+;
 #pragma link C++ class AliFMDEventPlaneFinder+;
 #pragma link C++ class AliFMDHistCollector+;
+#pragma link C++ class AliFMDESDFixer+;
 #pragma link C++ class AliFMDSharingFilter+;
 #if ROOT_VERSION_CODE < 0x56300 // ROOT_VERSION(5,99,0)
 #pragma link C++ class AliFMDSharingFilter::RingHistos+;