PWGLF/FORWARD/analysis2/AliFMDHistCollector.h

   1 //
   2 // This class collects the event histograms into single histograms,
   3 // one for each ring in each vertex bin.
   4 //
   5 #ifndef ALIFMDHISTCOLLECTOR_H
   6 #define ALIFMDHISTCOLLECTOR_H
   7 /**
   8  * @file   AliFMDHistCollector.h
   9  * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk>
  10  * @date   Wed Mar 23 14:03:01 2011
  11  *
  12  * @brief
  13  *
  14  *
  15  * @ingroup pwglf_forward_aod
  16  */
  17 #include <TNamed.h>
  18 #include <TList.h>
  19 #include <TArrayI.h>
  20 #include "AliForwardUtil.h"
  21 class AliESDFMD;
  22 class TH2D;
  23
  24 /**
  25  * This class collects the event histograms into single histograms,
  26  * one for each ring in each vertex bin.
  27  *
  28  * @par Input:
  29  *   - AliESDFMD object possibly corrected for sharing
  30  *
  31  * @par Output:
  32  *   - 5 RingHistos objects - each with a number of vertex dependent
  33  *     2D histograms of the inclusive charge particle density
  34  *
  35  * @par HistCollector used:
  36  *   - AliFMDCorrSecondaryMap
  37  *
  38  * @ingroup pwglf_forward_algo
  39  * @ingroup pwglf_forward_aod
  40  */
  41 class AliFMDHistCollector : public TNamed
  42 {
  43 public:
  44   /**
  45    * Methods to use when merging overlapping bins @f$b_1@f$, @f$b_2@f$
  46    * with content @f$c_1@f$, @f$c_2@f$, and errors @f$e_1@f$,
  47    * @f$e_2@f$ into bin @f$ b@f$ with content @f$c@f$ and error @f$e@f$
  48    */
  49   enum MergeMethod {
  50     /**
  51      * @f[
  52      *   c = \frac{1}{2}(c_1+c_2)
  53      * @f]
  54      * @f[
  55      *   e = \sqrt{e_1^2+e_2^2}
  56      * @f]
  57      */
  58     kStraightMean,
  59     /**
  60      * As above, exept zero's are ignored
  61      */
  62     kStraightMeanNoZero,
  63     /**
  64      * @f[
  65      *   c = \frac{\frac{c_1}{e_1^2}+\frac{c_2}{e_2^2}}{
  66      *             \frac{1}{e_1^2}+\frac{1}{e_2^2}}
  67      * @f]
  68      * @f[
  69      *   e = \sqrt{\frac{1}{\frac{1}{e_1^2}+\frac{1}{e_2^2}}}
  70      * @f]
  71      */
  72     kWeightedMean,
  73     /**
  74      * @f[
  75      *     c = \left\{\begin{array}{cl}
  76      *          c_1 & \mbox{if $e_1 < e_2$} \\
  77      *          c_2 & \mbox{otherwise}\end{array}\right.
  78      * @f]
  79      */
  80     kLeastError
  81   };
  82   /**
  83    * How to obtain the fiducial cuts
  84    */
  85   enum FiducialMethod {
  86     /**
  87      * Select bins by fixed cut.  Bins with a secondary correction
  88      * less than the cut is considered as non-valid
  89      */
  90     kByCut,
  91     /**
  92      * A bin is considered non-valid, if it is less then twice as
  93      * large as it's neighbors (in eta)
  94      */
  95     kDistance
  96   };
  97   /**
  98    * FMD ring bits for skipping
  99    */
 100    enum FMDRingBits {
 101      kFMD1I=0x01,
 102      kFMD1 =0x01,
 103      kFMD2I=0x02,
 104      kFMD2O=0x04,
 105      kFMD2 =0x06,
 106      kFMD3I=0x08,
 107      kFMD3O=0x10,
 108      kFMD3 =0x18
 109   };
 110   /**
 111    * Constructor
 112    */
 113   AliFMDHistCollector();
 114   /**
 115    * Constructor
 116    *
 117    * @param title Name of object
 118    */
 119   AliFMDHistCollector(const char* title);
 120   /**
 121    * Copy constructor
 122    *
 123    * @param o Object to copy from
 124    */
 125   AliFMDHistCollector(const AliFMDHistCollector& o);
 126
 127   /**
 128    * Destructor
 129    */
 130   virtual ~AliFMDHistCollector();
 131   /**
 132    * Assignement operator
 133    *
 134    * @param o Object to assign from
 135    *
 136    * @return Reference to this object
 137    */
 138   AliFMDHistCollector& operator=(const AliFMDHistCollector&);
 139   /**
 140    * Intialise
 141    *
 142    * @param vtxAxis  @f$ v_z@f$ axis
 143    * @param etaAxis  @f$ \eta@f$ axis
 144    */
 145   virtual void Init(const TAxis& vtxAxis,
 146                     const TAxis& etaAxis);
 147   /**
 148    * Do the calculations
 149    *
 150    * @param hists    Cache of histograms
 151    * @param sums     Cache to sum ring histograms in
 152    * @param vtxBin   Vertex bin (1 based)
 153    * @param out      Output histogram
 154    *
 155    * @return true on successs
 156    */
 157   virtual Bool_t Collect(const AliForwardUtil::Histos& hists,
 158                          AliForwardUtil::Histos&       sums,
 159                          UShort_t                      vtxBin,
 160                          TH2D&                         out);
 161   /**
 162    * Output diagnostic histograms to directory
 163    *
 164    * @param dir List to write in
 165    */
 166   virtual void DefineOutput(TList* dir);
 167   /**
 168    * Set the merge method
 169    *
 170    * @param m Method
 171    */
 172   void SetMergeMethod(MergeMethod m) { fMergeMethod = m; }
 173   /**
 174    * Set the method for finding the fidicual area of the secondary maps
 175    *
 176    * @param m Method
 177    */
 178   void SetFiducialMethod(FiducialMethod m) { fFiducialMethod = m; }
 179   /**
 180    * Set the number of extra bins (beyond the secondary map border)
 181    * to cut away.
 182    *
 183    * @param n Number of bins
 184    */
 185   void SetNCutBins(UInt_t n=2) { fNCutBins = n; }
 186   /**
 187    * Set the correction cut, that is, when bins in the secondary
 188    * correction maps have a value lower than this cut, they are
 189    * considered uncertain and not used
 190    *
 191    * @param cut Cut-off
 192    */
 193   void SetCorrectionCut(Float_t cut=0.5) { fCorrectionCut = cut; }
 194   /**
 195    * Set FMD rings to skip. Argument should be kFirstRingToSkip|kSecondRingToSkip...
 196    *
 197    * @param mask bit pattern
 198    */
 199   void SetFMDRingsToSkip(UShort_t mask) { fSkipFMDRings = mask; }
 200  /**
 201    * Set whether to make bg maps or not
 202    *
 203    * @param use make them
 204    */
 205   void SetMakeBGHitMaps(Bool_t use) { fBgAndHitMaps = use; }
 206
 207   /**
 208    * Set the debug level. The higher the value the more output
 209    *
 210    * @param dbg Debug level
 211    */
 212
 213   void SetDebug(Int_t dbg=1) { fDebug = dbg; }
 214   /**
 215    * Print information
 216    *
 217    * @param option Not used
 218    */
 219   void Print(Option_t* option="") const;
 220 protected:
 221   /**
 222    * Get the first and last eta bin to use for a given ring and vertex
 223    *
 224    * @param d        Detector
 225    * @param r        Ring
 226    * @param vtxBin   Vertex bin (1 based)
 227    * @param first    On return, the first eta bin to use
 228    * @param last     On return, the last eta bin to use
 229    */
 230   virtual void GetFirstAndLast(UShort_t d, Char_t r, UShort_t vtxBin,
 231                                Int_t& first, Int_t& last) const;
 232   /**
 233    * Get the first and last eta bin to use for a given ring and vertex
 234    *
 235    * @param idx      Ring index as given by GetIdx
 236    * @param vtxBin   Vertex bin (1 based)
 237    * @param first    On return, the first eta bin to use
 238    * @param last     On return, the last eta bin to use
 239    */
 240   virtual void GetFirstAndLast(Int_t idx, UShort_t vtxBin,
 241                                Int_t& first, Int_t& last) const;
 242   /**
 243    * Get the first eta bin to use for a given ring and vertex
 244    *
 245    * @param d Detector
 246    * @param r Ring
 247    * @param v vertex bin (1 based)
 248    *
 249    * @return First eta bin to use, or -1 in case of problems
 250    */
 251   Int_t GetFirst(UShort_t d, Char_t r, UShort_t v) const;
 252   /**
 253    * Get the first eta bin to use for a given ring and vertex
 254    *
 255    * @param idx Ring index as given by GetIdx
 256    * @param v vertex bin (1 based)
 257    *
 258    * @return First eta bin to use, or -1 in case of problems
 259    */
 260   Int_t GetFirst(Int_t idx, UShort_t v) const;
 261   /**
 262    * Get the last eta bin to use for a given ring and vertex
 263    *
 264    * @param d Detector
 265    * @param r Ring
 266    * @param v vertex bin (1 based)
 267    *
 268    * @return Last eta bin to use, or -1 in case of problems
 269    */
 270   Int_t GetLast(UShort_t d, Char_t r, UShort_t v) const;
 271   /**
 272    * Get the last eta bin to use for a given ring and vertex
 273    *
 274    * @param idx Ring index as given by GetIdx
 275    * @param v vertex bin (1 based)
 276    *
 277    * @return Last eta bin to use, or -1 in case of problems
 278    */
 279   Int_t GetLast(Int_t idx, UShort_t v) const;
 280   /**
 281    * Get the detector and ring from the ring index
 282    *
 283    * @param idx Ring index
 284    * @param d   On return, the detector or 0 in case of errors
 285    * @param r   On return, the ring id or '0' in case of errors
 286    */
 287   void GetDetRing(Int_t idx, UShort_t& d, Char_t& r) const;
 288   /**
 289    * Get the ring index from detector number and ring identifier
 290    *
 291    * @param d Detector
 292    * @param r Ring identifier
 293    *
 294    * @return ring index or -1 in case of problems
 295    */
 296   Int_t GetIdx(UShort_t d, Char_t r) const;
 297   /**
 298    * Get the possibly overlapping histogram of eta bin @a e in
 299    * detector and ring
 300    *
 301    * @param d Detector
 302    * @param r Ring
 303    * @param e Eta bin
 304    * @param v Vertex bin (1 based)
 305    *
 306    * @return Overlapping histogram index or -1
 307    */
 308   Int_t GetOverlap(UShort_t d, Char_t r, Int_t e, UShort_t v) const;
 309   /**
 310    * Get the possibly overlapping histogram of eta bin @a e in
 311    * detector and ring
 312    *
 313    * @param i Ring index
 314    * @param e Eta bin
 315    * @param v Vertex bin (1 based)
 316    *
 317    * @return Overlapping histogram index or -1
 318    */
 319   Int_t GetOverlap(Int_t i, Int_t e, UShort_t v) const;
 320   /**
 321    * Check if there's an overlapping histogram with this eta bin of
 322    * the detector and ring
 323    *
 324    * @param d Detector
 325    * @param r Ring
 326    * @param e eta bin
 327    * @param v Vertex bin (1 based)
 328    *
 329    * @return True if there's an overlapping histogram
 330    */
 331   Bool_t HasOverlap(UShort_t d, Char_t r, Int_t e, UShort_t v) const;
 332   /**
 333    * Check if there's an overlapping histogram with this eta bin of
 334    * ring
 335    *
 336    * @param i Ring index
 337    * @param e eta bin
 338    * @param v Vertex bin
 339    *
 340    * @return True if there's an overlapping histogram
 341    */
 342   Bool_t HasOverlap(Int_t i, Int_t e, UShort_t v) const;
 343   /**
 344    * Check if we should include the bin in the data range
 345    *
 346    * @param bg Secondary map histogram
 347    * @param ie Eta bin
 348    * @param ip Phi bin
 349    *
 350    * @return True if to be used
 351    */
 352   Bool_t CheckCorrection(const TH2D* bg, Int_t ie, Int_t ip) const;
 353   /**
 354    * Merge bins accoring to set method
 355    *
 356    * @param c   Current content
 357    * @param e   Current error
 358    * @param oc  Old content
 359    * @param oe  Old error
 360    * @param rc  On return, the new content
 361    * @param re  On return, tne new error
 362    */
 363   void MergeBins(Double_t c,   Double_t e,
 364                  Double_t oc,  Double_t oe,
 365                  Double_t& rc, Double_t& re) const;
 366
 367
 368   Int_t       fNCutBins;        // Number of additional bins to cut away
 369   Float_t     fCorrectionCut;   // Cut-off on secondary corrections
 370   TArrayI     fFirstBins;       // Array of first eta bins
 371   TArrayI     fLastBins;        // Array of last eta bins
 372   Int_t       fDebug;           // Debug level
 373   TList*      fList;            // Output list
 374   TH2D*       fSumRings;        // Sum per ring (on y-axis)
 375   TH2D*       fCoverage;        // Sum per ring (on y-axis)
 376   MergeMethod fMergeMethod;     // Merge methiod for overlapping bins
 377   FiducialMethod fFiducialMethod; // Fidicual method
 378   UShort_t    fSkipFMDRings;    // FMD rings to ignore
 379   Bool_t      fBgAndHitMaps;    // Make hit/bg maps or not
 380
 381   ClassDef(AliFMDHistCollector,2); // Calculate Nch density
 382 };
 383
 384 //____________________________________________________________________
 385 inline void
 386 AliFMDHistCollector::GetFirstAndLast(UShort_t d, Char_t r, UShort_t vtxbin,
 387                                      Int_t& first, Int_t& last) const
 388 {
 389   GetFirstAndLast(GetIdx(d,r), vtxbin, first, last);
 390 }
 391 //____________________________________________________________________
 392 inline Int_t
 393 AliFMDHistCollector::GetFirst(UShort_t d, Char_t r, UShort_t v) const
 394 {
 395   return GetFirst(GetIdx(d,r), v);
 396 }
 397 //____________________________________________________________________
 398 inline Int_t
 399 AliFMDHistCollector::GetLast(UShort_t d, Char_t r, UShort_t v) const
 400 {
 401   return GetLast(GetIdx(d, r), v);
 402 }
 403 //____________________________________________________________________
 404 inline Bool_t
 405 AliFMDHistCollector::HasOverlap(UShort_t d, Char_t r, Int_t e, UShort_t v) const
 406 {
 407   return GetOverlap(d,r,e,v) >= 0;
 408 }
 409 //____________________________________________________________________
 410 inline Bool_t
 411 AliFMDHistCollector::HasOverlap(Int_t i, Int_t e, UShort_t v) const
 412 {
 413   return GetOverlap(i,e,v) >= 0;
 414 }
 415
 416 #endif
 417 // Local Variables:
 418 //   mode: C++
 419 // End:
 420