PWGLF/FORWARD/analysis2/AliForwardUtil.h

   1 //
   2 // Utilities used in the forward multiplcity analysis
   3 //
   4 //
   5 #ifndef ALIFORWARDUTIL_H
   6 #define ALIFORWARDUTIL_H
   7 /**
   8  * @file   AliForwardUtil.h
   9  * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk>
  10  * @date   Wed Mar 23 14:06:54 2011
  11  *
  12  * @brief
  13  *
  14  *
  15  * @ingroup pwglf_forward
  16  */
  17 #include <TObject.h>
  18 #include <TString.h>
  19 #include <TObjArray.h>
  20 class TH2D;
  21 class TH1I;
  22 class TH1;
  23 class TF1;
  24 class TAxis;
  25 class AliESDEvent;
  26 class AliAODEvent;
  27 class AliAnalysisTaskSE;
  28
  29 /**
  30  * Utilities used in the forward multiplcity analysis
  31  *
  32  * @ingroup pwglf_forward
  33  */
  34 class AliForwardUtil : public TObject
  35 {
  36 public:
  37   /**
  38    * Get the standard color for a ring
  39    *
  40    * @param d Detector
  41    * @param r Ring
  42    *
  43    * @return
  44    */
  45   static Color_t RingColor(UShort_t d, Char_t r)
  46   {
  47     return ((d == 1 ? kRed : (d == 2 ? kGreen : kBlue))
  48             + ((r == 'I' || r == 'i') ? 2 : -3));
  49   }
  50   //==================================================================
  51   /**
  52    * @{
  53    * @name AliROOT version
  54    */
  55   /**
  56    * Get the revision number of AliROOT
  57    *
  58    * @return Subversion revision number of AliROOT used
  59    */
  60   static ULong_t AliROOTRevision();
  61   /**
  62    * Get the branch identifier of AliROOT.  In case of trunk, return
  63    * 0xFFFFFFFF, while for @b vM-N-R{-S}, we get
  64    *
  65    * @code
  66    *    ((M & 0xFF) << 12 | (N & 0xFF) << 8 | (R & 0xFF) << 3 | (X))
  67    * @endcode
  68    * where @c X is 0xAA if @b S is specified (e.g., analysis tag).
  69    *
  70    * @return branch identifer encoded in bits
  71    */
  72   static ULong_t AliROOTBranch();
  73   //==================================================================
  74   /**
  75    * @{
  76    * @name Collision/run parameters
  77    */
  78   /**
  79    * Defined collision types
  80    */
  81   enum ECollisionSystem {
  82     kUnknown,
  83     kPP,
  84     kPbPb,
  85     kPPb
  86   };
  87   //__________________________________________________________________
  88   /**
  89    * Calculate the beam rapidity.
  90    *
  91    * @b Note: The beam energy is given in GeV/charge
  92    *
  93    * @param beam Beam energy in GeV/charge
  94    * @param z    Charge number of projectile
  95    * @param a    Mass number of projectile
  96    *
  97    * @return The rapidity of the beam
  98    */
  99   static Float_t BeamRapidity(Float_t beam, UShort_t z, UShort_t a);
 100   /**
 101    * Calculate the center of mass energy from the beam energy per
 102    * charge and the nucleus numbers.
 103    *
 104    * @param beam Beam energy in GeV/charge
 105    * @param z1   Charge number of projectile
 106    * @param a1   Mass number of projectile
 107    * @param z2   Charge number of projectile
 108    * @param a2   Mass number of projectile
 109    *
 110    * @return The center of mass energy in GeV/nucleon
 111    */
 112   static Float_t CenterOfMassEnergy(Float_t beam, UShort_t z1, UShort_t a1,
 113                                     Short_t z2=-1, Short_t a2=-1);
 114   /**
 115    * Calculate the center of mass rapidity (shift)
 116    *
 117    * @param z1 Charge number of projectile
 118    * @param a1 Mass number of projectile
 119    * @param z2 Charge number of projectile
 120    * @param a2 Mass number of projectile
 121    * x
 122    * @return Rapidity of the center of mass
 123    */
 124   static Float_t CenterOfMassRapidity(UShort_t z1, UShort_t a1,
 125                                       Short_t z2=-1, Short_t a2=-1);
 126   /**
 127    * Parse a collision system spec given in a string.   Known values are
 128    *
 129    *  - "pp", "p-p" which returns kPP
 130    *  - "PbPb", "Pb-Pb", "A-A", which returns kPbPb
 131    *  - "pPb", "p-Pb", "pA", p-A" which returns kPPb
 132    *  - Everything else gives kUnknown
 133    *
 134    * @param sys Collision system spec
 135    *
 136    * @return Collision system id
 137    */
 138   static UShort_t ParseCollisionSystem(const char* sys);
 139   /**
 140    * Get a string representation of the collision system
 141    *
 142    * @param sys  Collision system
 143    * - kPP -> "pp"
 144    * - kPbPb -> "PbPb"
 145    * - kPPb -> "pPb"
 146    * - anything else gives "unknown"
 147    *
 148    * @return String representation of the collision system
 149    */
 150   static const char* CollisionSystemString(UShort_t sys);
 151   //__________________________________________________________________
 152   /**
 153    * Parse the center of mass energy given as a float and return known
 154    * values as a unsigned integer
 155    *
 156    * @param sys  Collision system (needed for AA)
 157    * @param cms  Center of mass energy * total charge
 158    *
 159    * @return Center of mass energy per nucleon
 160    */
 161   static UShort_t ParseCenterOfMassEnergy(UShort_t sys, Float_t cms);
 162   /**
 163    * Get a string representation of the center of mass energy per nuclean
 164    *
 165    * @param cms  Center of mass energy per nucleon
 166    *
 167    * @return String representation of the center of mass energy per nuclean
 168    */
 169   static const char* CenterOfMassEnergyString(UShort_t cms);
 170   //__________________________________________________________________
 171   /**
 172    * Parse the magnetic field (in kG) as given by a floating point number
 173    *
 174    * @param field  Magnetic field in kG
 175    *
 176    * @return Short integer value of magnetic field in kG
 177    */
 178   static Short_t ParseMagneticField(Float_t field);
 179   /**
 180    * Get a string representation of the magnetic field
 181    *
 182    * @param field Magnetic field in kG
 183    *
 184    * @return String representation of the magnetic field
 185    */
 186   static const char* MagneticFieldString(Short_t field);
 187   /* @} */
 188
 189   //==================================================================
 190   /**
 191    * @{
 192    * @name Recalculate @f$\eta@f$, @f$\phi@f$, etc.
 193    */
 194   /**
 195    * Get the radius of a strip.
 196    *
 197    * @param ring  Ring identifier 'I' or 'O'
 198    * @param strip Strip number
 199    *
 200    * @return Radial distance from beam of the strip
 201    */
 202   static Double_t GetStripR(Char_t ring, UShort_t strip);
 203   /**
 204    * Get eta from strip
 205    *
 206    * @param det, ring, sec, strip, zvtx
 207    *
 208    * @return eta
 209    */
 210   static Double_t GetEtaFromStrip(UShort_t det, Char_t ring,
 211                                   UShort_t sec, UShort_t strip, Double_t zvtx);
 212   /**
 213    * Get the azimuthal angle of a strip
 214    *
 215    * @param ring  Ring identifier 'I' or 'O'
 216    * @param strip Strip number
 217    * @param phi   Straight forward strip phi
 218    * @param xvtx  Ip X coordinate
 219    * @param yvtx  Ip Y coordinate
 220    *
 221    * @return The phi angle correctef for (X,Y) off set.
 222    */
 223   static Double_t GetPhiFromStrip(Char_t ring, UShort_t strip,
 224                                   Double_t phi, Double_t xvtx, Double_t yvtx);
 225   /* @} */
 226
 227   //==================================================================
 228   /**
 229    * @{
 230    * @name Manager related tasks
 231    */
 232   /**
 233    * Get the AOD event - either from the input (AOD analysis) or the
 234    * output (ESD analysis)
 235    *
 236    * @param task Task to do the investigation for
 237    *
 238    * @return Found AOD event or null
 239    */
 240   static AliAODEvent* GetAODEvent(AliAnalysisTaskSE* task);
 241   /**
 242    * Check if we have something that will provide and AOD event
 243    *
 244    * @return 0 if there's nothing that provide an AOD event, 1 if it
 245    * is provided on the input (AOD analysis) or 2 if it is provided on
 246    * the output (ESD analysis)
 247    */
 248   static UShort_t CheckForAOD();
 249   /**
 250    * Check if we have a particular (kind) of task in our train
 251    *
 252    * @param clsOrName  Class name or name of task
 253    * @param cls If true, look for a task of a particular class -
 254    * otherwise search for a speficially name task
 255    *
 256    * @return true if the needed task was found
 257    */
 258   static Bool_t CheckForTask(const char* clsOrName, Bool_t cls=true);
 259   /* @} */
 260
 261   //==================================================================
 262   /**
 263    * @{
 264    * @name Member functions to store and retrieve analysis parameters
 265    */
 266   static TObject* MakeParameter(const char* name, UShort_t value);
 267   static TObject* MakeParameter(const char* name, Int_t value);
 268   static TObject* MakeParameter(const char* name, Double_t value);
 269   static TObject* MakeParameter(const char* name, Bool_t value);
 270   static TObject* MakeParameter(const char* name, ULong_t value);
 271   static void GetParameter(TObject* o, UShort_t& value);
 272   static void GetParameter(TObject* o, Int_t& value);
 273   static void GetParameter(TObject* o, Double_t& value);
 274   static void GetParameter(TObject* o, Bool_t& value);
 275   static void GetParameter(TObject* o, ULong_t& value);
 276   /* @} */
 277
 278   //==================================================================
 279   /**
 280    * @{
 281    * @name Energy stragling functions
 282    */
 283   //__________________________________________________________________
 284   /**
 285    * Number of steps to do in the Landau, Gaussiam convolution
 286    */
 287   static Int_t fgConvolutionSteps; // Number of convolution steps
 288   //------------------------------------------------------------------
 289   /**
 290    * How many sigma's of the Gaussian in the Landau, Gaussian
 291    * convolution to integrate over
 292    */
 293   static Double_t fgConvolutionNSigma; // Number of convolution sigmas
 294   //------------------------------------------------------------------
 295   /**
 296    * Calculate the shifted Landau
 297    * @f[
 298    *    f'_{L}(x;\Delta,\xi) = f_L(x;\Delta+0.22278298\xi)
 299    * @f]
 300    *
 301    * where @f$ f_{L}@f$ is the ROOT implementation of the Landau
 302    * distribution (known to have @f$ \Delta_{p}=-0.22278298@f$ for
 303    * @f$\Delta=0,\xi=1@f$.
 304    *
 305    * @param x      Where to evaluate @f$ f'_{L}@f$
 306    * @param delta  Most probable value
 307    * @param xi     The 'width' of the distribution
 308    *
 309    * @return @f$ f'_{L}(x;\Delta,\xi) @f$
 310    */
 311   static Double_t Landau(Double_t x, Double_t delta, Double_t xi);
 312
 313   //------------------------------------------------------------------
 314   /**
 315    * Calculate the value of a Landau convolved with a Gaussian
 316    *
 317    * @f[
 318    * f(x;\Delta,\xi,\sigma') = \frac{1}{\sigma' \sqrt{2 \pi}}
 319    *    \int_{-\infty}^{+\infty} d\Delta' f'_{L}(x;\Delta',\xi)
 320    *    \exp{-\frac{(\Delta-\Delta')^2}{2\sigma'^2}}
 321    * @f]
 322    *
 323    * where @f$ f'_{L}@f$ is the Landau distribution, @f$ \Delta@f$ the
 324    * energy loss, @f$ \xi@f$ the width of the Landau, and
 325    * @f$ \sigma'^2=\sigma^2-\sigma_n^2 @f$.  Here, @f$\sigma@f$ is the
 326    * variance of the Gaussian, and @f$\sigma_n@f$ is a parameter modelling
 327    * noise in the detector.
 328    *
 329    * Note that this function uses the constants fgConvolutionSteps and
 330    * fgConvolutionNSigma
 331    *
 332    * References:
 333    *  - <a href="http://dx.doi.org/10.1016/0168-583X(84)90472-5">Nucl.Instrum.Meth.B1:16</a>
 334    *  - <a href="http://dx.doi.org/10.1103/PhysRevA.28.615">Phys.Rev.A28:615</a>
 335    *  - <a href="http://root.cern.ch/root/htmldoc/tutorials/fit/langaus.C.html">ROOT implementation</a>
 336    *
 337    * @param x         where to evaluate @f$ f@f$
 338    * @param delta     @f$ \Delta@f$ of @f$ f(x;\Delta,\xi,\sigma')@f$
 339    * @param xi        @f$ \xi@f$ of @f$ f(x;\Delta,\xi,\sigma')@f$
 340    * @param sigma     @f$ \sigma@f$ of @f$\sigma'^2=\sigma^2-\sigma_n^2 @f$
 341    * @param sigma_n   @f$ \sigma_n@f$ of @f$\sigma'^2=\sigma^2-\sigma_n^2 @f$
 342    *
 343    * @return @f$ f@f$ evaluated at @f$ x@f$.
 344    */
 345   static Double_t LandauGaus(Double_t x, Double_t delta, Double_t xi,
 346                              Double_t sigma, Double_t sigma_n);
 347
 348   //------------------------------------------------------------------
 349   /**
 350    * Evaluate
 351    * @f[
 352    *    f_i(x;\Delta,\xi,\sigma') = f(x;\Delta_i,\xi_i,\sigma_i')
 353    * @f]
 354    * corresponding to @f$ i@f$ particles i.e., with the substitutions
 355    * @f{eqnarray*}{
 356    *    \Delta    \rightarrow \Delta_i    &=& i(\Delta + \xi\log(i))\\
 357    *    \xi       \rightarrow \xi_i       &=& i \xi\\
 358    *    \sigma    \rightarrow \sigma_i    &=& \sqrt{i}\sigma\\
 359    *    \sigma'^2 \rightarrow \sigma_i'^2 &=& \sigma_n^2 + \sigma_i^2
 360    * @f}
 361    *
 362    * @param x        Where to evaluate
 363    * @param delta    @f$ \Delta@f$
 364    * @param xi       @f$ \xi@f$
 365    * @param sigma    @f$ \sigma@f$
 366    * @param sigma_n  @f$ \sigma_n@f$
 367    * @param i        @f$ i @f$
 368    *
 369    * @return @f$ f_i @f$ evaluated
 370    */
 371   static Double_t ILandauGaus(Double_t x, Double_t delta, Double_t xi,
 372                               Double_t sigma, Double_t sigma_n, Int_t i);
 373
 374   //------------------------------------------------------------------
 375   /**
 376    * Numerically evaluate
 377    * @f[
 378    *    \left.\frac{\partial f_i}{\partial p_i}\right|_{x}
 379    * @f]
 380    * where @f$ p_i@f$ is the @f$ i^{\mbox{th}}@f$ parameter.  The mapping
 381    * of the parameters is given by
 382    *
 383    * - 0: @f$\Delta@f$
 384    * - 1: @f$\xi@f$
 385    * - 2: @f$\sigma@f$
 386    * - 3: @f$\sigma_n@f$
 387    *
 388    * This is the partial derivative with respect to the parameter of
 389    * the response function corresponding to @f$ i@f$ particles i.e.,
 390    * with the substitutions
 391    * @f[
 392    *    \Delta    \rightarrow \Delta_i    = i(\Delta + \xi\log(i))\\
 393    *    \xi       \rightarrow \xi_i       = i \xi\\
 394    *    \sigma    \rightarrow \sigma_i    = \sqrt{i}\sigma\\
 395    *    \sigma'^2 \rightarrow \sigma_i'^2 = \sigma_n^2 + \sigma_i^2
 396    * @f]
 397    *
 398    * @param x        Where to evaluate
 399    * @param ipar     Parameter number
 400    * @param dp       @f$ \epsilon\delta p_i@f$ for some value of @f$\epsilon@f$
 401    * @param delta    @f$ \Delta@f$
 402    * @param xi       @f$ \xi@f$
 403    * @param sigma    @f$ \sigma@f$
 404    * @param sigma_n  @f$ \sigma_n@f$
 405    * @param i        @f$ i@f$
 406    *
 407    * @return @f$ f_i@f$ evaluated
 408    */
 409   static Double_t IdLandauGausdPar(Double_t x, UShort_t ipar, Double_t dp,
 410                                    Double_t delta, Double_t xi,
 411                                    Double_t sigma, Double_t sigma_n, Int_t i);
 412
 413   //------------------------------------------------------------------
 414   /**
 415    * Evaluate
 416    * @f[
 417    *   f_N(x;\Delta,\xi,\sigma') = \sum_{i=1}^N a_i f_i(x;\Delta,\xi,\sigma'a)
 418    * @f]
 419    *
 420    * where @f$ f(x;\Delta,\xi,\sigma')@f$ is the convolution of a
 421    * Landau with a Gaussian (see LandauGaus).  Note that
 422    * @f$ a_1 = 1@f$, @f$\Delta_i = i(\Delta_1 + \xi\log(i))@f$,
 423    * @f$\xi_i=i\xi_1@f$, and @f$\sigma_i'^2 = \sigma_n^2 + i\sigma_1^2@f$.
 424    *
 425    * References:
 426    *  - <a href="http://dx.doi.org/10.1016/0168-583X(84)90472-5">Nucl.Instrum.Meth.B1:16</a>
 427    *  - <a href="http://dx.doi.org/10.1103/PhysRevA.28.615">Phys.Rev.A28:615</a>
 428    *  - <a href="http://root.cern.ch/root/htmldoc/tutorials/fit/langaus.C.html">ROOT implementation</a>
 429    *
 430    * @param x        Where to evaluate @f$ f_N@f$
 431    * @param delta    @f$ \Delta_1@f$
 432    * @param xi       @f$ \xi_1@f$
 433    * @param sigma    @f$ \sigma_1@f$
 434    * @param sigma_n  @f$ \sigma_n@f$
 435    * @param n        @f$ N@f$ in the sum above.
 436    * @param a        Array of size @f$ N-1@f$ of the weights @f$ a_i@f$ for
 437    *                 @f$ i > 1@f$
 438    *
 439    * @return @f$ f_N(x;\Delta,\xi,\sigma')@f$
 440    */
 441   static Double_t NLandauGaus(Double_t x, Double_t delta, Double_t xi,
 442                               Double_t sigma, Double_t sigma_n, Int_t n,
 443                               const Double_t* a);
 444   /**
 445    * Generate a TF1 object of @f$ f_I@f$
 446    *
 447    * @param c        Constant
 448    * @param delta    @f$ \Delta@f$
 449    * @param xi       @f$ \xi_1@f$
 450    * @param sigma    @f$ \sigma_1@f$
 451    * @param sigma_n  @f$ \sigma_n@f$
 452    * @param i        @f$ i@f$ - the number of particles
 453    * @param xmin     Least value of range
 454    * @param xmax     Largest value of range
 455    *
 456    * @return Newly allocated TF1 object
 457    */
 458   static TF1* MakeILandauGaus(Double_t c,
 459                               Double_t delta, Double_t xi,
 460                               Double_t sigma, Double_t sigma_n,
 461                               Int_t    i,
 462                               Double_t xmin,  Double_t  xmax);
 463   /**
 464    * Generate a TF1 object of @f$ f_N@f$
 465    *
 466    * @param c         Constant
 467    * @param delta     @f$ \Delta@f$
 468    * @param xi        @f$ \xi_1@f$
 469    * @param sigma     @f$ \sigma_1@f$
 470    * @param sigma_n   @f$ \sigma_n@f$
 471    * @param n         @f$ N@f$ - how many particles to sum to
 472    * @param a         Array of size @f$ N-1@f$ of the weights @f$ a_i@f$ for
 473    *                  @f$ i > 1@f$
 474    * @param xmin      Least value of range
 475    * @param xmax      Largest value of range
 476    *
 477    * @return Newly allocated TF1 object
 478    */
 479   static TF1* MakeNLandauGaus(Double_t c,
 480                               Double_t delta, Double_t  xi,
 481                               Double_t sigma, Double_t  sigma_n,
 482                               Int_t    n,     const Double_t* a,
 483                               Double_t xmin,  Double_t  xmax);
 484
 485   //__________________________________________________________________
 486   /**
 487    * Structure to do fits to the energy loss spectrum
 488    *
 489    * @ingroup pwglf_forward
 490    */
 491   struct ELossFitter
 492   {
 493     enum {
 494       kC     = 0,
 495       kDelta,
 496       kXi,
 497       kSigma,
 498       kSigmaN,
 499       kN,
 500       kA
 501     };
 502     /**
 503      * Constructor
 504      *
 505      * @param lowCut     Lower cut of spectrum - data below this cuts is ignored
 506      * @param maxRange   Maximum range to fit to
 507      * @param minusBins  The number of bins below maximum to use
 508      */
 509     ELossFitter(Double_t lowCut, Double_t maxRange, UShort_t minusBins);
 510     /**
 511      * Destructor
 512      *
 513      */
 514     virtual ~ELossFitter();
 515     void SetDebug(Bool_t debug=true) { fDebug = debug; }
 516     /**
 517      * Clear internal arrays
 518      *
 519      */
 520     void Clear();
 521     /**
 522      * Fit a 1-particle signal to the passed energy loss distribution
 523      *
 524      * Note that this function clears the internal arrays first
 525      *
 526      * @param dist    Data to fit the function to
 527      * @param sigman If larger than zero, the initial guess of the
 528      *               detector induced noise. If zero or less, then this
 529      *               parameter is ignored in the fit (fixed at 0)
 530      *
 531      * @return The function fitted to the data
 532      */
 533     TF1* Fit1Particle(TH1* dist, Double_t sigman=-1);
 534     /**
 535      * Fit a N-particle signal to the passed energy loss distribution
 536      *
 537      * If there's no 1-particle fit present, it does that first
 538      *
 539      * @param dist   Data to fit the function to
 540      * @param n      Number of particle signals to fit
 541      * @param sigman If larger than zero, the initial guess of the
 542      *               detector induced noise. If zero or less, then this
 543      *               parameter is ignored in the fit (fixed at 0)
 544      *
 545      * @return The function fitted to the data
 546      */
 547     TF1* FitNParticle(TH1* dist, UShort_t n, Double_t sigman=-1);
 548     /**
 549      * Fit a composite distribution of energy loss from both primaries
 550      * and secondaries
 551      *
 552      * @param dist   Distribution
 553      * @param sigman If larger than zero, the initial guess of the
 554      *                detector included noise.  If zero or less this
 555      *                parameter is fixed to 0.
 556      *
 557      * @return Function fitted to the data
 558      */
 559     TF1* FitComposite(TH1* dist, Double_t sigman);
 560     /**
 561      * Get Lower cut on data
 562      *
 563      * @return Lower cut on data
 564      */
 565     Double_t GetLowCut() const { return fLowCut; }
 566     /**
 567      * Get Maximum range to fit
 568      *
 569      * @return Maximum range to fit
 570      */
 571     Double_t GetMaxRange() const { return fMaxRange; }
 572     /**
 573      * Get Number of bins from maximum to fit 1st peak
 574      *
 575      * @return Number of bins from maximum to fit 1st peak
 576      */
 577     UShort_t GetMinusBins() const { return fMinusBins; }
 578     /**
 579      * Get Array of fit results
 580      *
 581      * @return Array of fit results
 582      */
 583     const TObjArray& GetFitResults() const { return fFitResults; }
 584     /**
 585      * Get Array of fit results
 586      *
 587      * @return Array of fit results
 588      */
 589     TObjArray& GetFitResults() { return fFitResults; }
 590     /**
 591      * Get Array of functions
 592      *
 593      * @return Array of functions
 594      */
 595     const TObjArray& GetFunctions() const { return fFunctions; }
 596     /**
 597      * Get Array of functions
 598      *
 599      * @return Array of functions
 600      */
 601     TObjArray& GetFunctions() { return fFunctions; }
 602   private:
 603     const Double_t fLowCut;     // Lower cut on data
 604     const Double_t fMaxRange;   // Maximum range to fit
 605     const UShort_t fMinusBins;  // Number of bins from maximum to fit 1st peak
 606     TObjArray fFitResults;      // Array of fit results
 607     TObjArray fFunctions;       // Array of functions
 608     Bool_t    fDebug;
 609   };
 610   /* @} */
 611
 612
 613   //==================================================================
 614   /**
 615    * @{
 616    * @name Convenience containers
 617    */
 618   /**
 619    * Structure to hold histograms
 620    *
 621    * @ingroup pwglf_forward
 622    */
 623   struct Histos : public TObject
 624   {
 625     /**
 626      * Constructor
 627      *
 628      *
 629      */
 630     Histos() : fFMD1i(0), fFMD2i(0), fFMD2o(0), fFMD3i(0), fFMD3o(0) {}
 631     /**
 632      * Copy constructor
 633      *
 634      * @param o Object to copy from
 635      */
 636     Histos(const Histos& o)
 637       : TObject(o),
 638         fFMD1i(o.fFMD1i),
 639         fFMD2i(o.fFMD2i),
 640         fFMD2o(o.fFMD2o),
 641         fFMD3i(o.fFMD3i),
 642         fFMD3o(o.fFMD3o)
 643     {}
 644     /**
 645      * Assignement operator
 646      *
 647      * @return Reference to this
 648      */
 649     Histos& operator=(const Histos&) { return *this;}
 650     /**
 651      * Destructor.  This does not delete the interally allocated
 652      * memory.  Use the member function Delete for that.
 653      */
 654     ~Histos();
 655     /**
 656      * Clear internal memory.  Note, if the internal histograms are
 657      * added to an output container, then we must not free this
 658      * memory.
 659      */
 660     void Delete(Option_t* opt="");
 661     /**
 662      * Initialize the object
 663      *
 664      * @param etaAxis Eta axis to use
 665      */
 666     void Init(const TAxis& etaAxis);
 667     /**
 668      * Re-initialize the object with new @f$\eta@f$ axis
 669      *
 670      * @param etaAxis Eta axis to use
 671      */
 672     void ReInit(const TAxis& etaAxis);
 673     /**
 674      * Make a histogram
 675      *
 676      * @param d        Detector
 677      * @param r        Ring
 678      * @param etaAxis  Eta axis to use
 679      *
 680      * @return Newly allocated histogram
 681      */
 682     static TH2D* Make(UShort_t d, Char_t r, const TAxis& etaAxis);
 683     /**
 684      * Set the @f$\eta@f$ axis
 685      *
 686      * @param hist    Histogram
 687      * @param etaAxis @f$\eta@f$ axis to use
 688      */
 689     static void RebinEta(TH2D* hist, const TAxis& etaAxis);
 690     /**
 691      * Clear data
 692      *
 693      * @param option Not used
 694      */
 695     void  Clear(Option_t* option="");
 696     // const TH2D* Get(UShort_t d, Char_t r) const;
 697     /**
 698      * Get the histogram for a particular detector,ring
 699      *
 700      * @param d Detector
 701      * @param r Ring
 702      *
 703      * @return Histogram for detector,ring or nul
 704      */
 705     TH2D* Get(UShort_t d, Char_t r) const;
 706     TH2D* fFMD1i; // Histogram for FMD1i
 707     TH2D* fFMD2i; // Histogram for FMD2i
 708     TH2D* fFMD2o; // Histogram for FMD2o
 709     TH2D* fFMD3i; // Histogram for FMD3i
 710     TH2D* fFMD3o; // Histogram for FMD3o
 711
 712     ClassDef(Histos,2)
 713   };
 714
 715   //__________________________________________________________________
 716   /**
 717    * Base class for structure holding ring specific histograms
 718    *
 719    * @ingroup pwglf_forward
 720    */
 721   struct RingHistos : public TObject
 722   {
 723     /**
 724      * Constructor
 725      *
 726      */
 727     RingHistos() : fDet(0), fRing('\0'), fName("") {}
 728     /**
 729      *
 730      *
 731      * @param d Detector
 732      * @param r Ring
 733      */
 734     RingHistos(UShort_t d, Char_t r)
 735       : fDet(d), fRing(r), fName(TString::Format("FMD%d%c", d, r))
 736     {}
 737     /**
 738      * Copy constructor
 739      *
 740      * @param o Object to copy from
 741      */
 742     RingHistos(const RingHistos& o)
 743       : TObject(o), fDet(o.fDet), fRing(o.fRing), fName(o.fName)
 744     {}
 745     /**
 746      *
 747      */
 748     virtual ~RingHistos() {}
 749     /**
 750      * Assignement operator
 751      *
 752      * @param o Object to assign from
 753      *
 754      * @return Reference to this
 755      */
 756     RingHistos& operator=(const RingHistos& o)
 757     {
 758       if (&o == this) return *this;
 759       TObject::operator=(o);
 760       fDet  = o.fDet;
 761       fRing = o.fRing;
 762       fName = o.fName;
 763       return *this;
 764     }
 765     /**
 766      * Define the outout list in @a d
 767      *
 768      * @param d Where to put the output list
 769      *
 770      * @return Newly allocated TList object or null
 771      */
 772     TList* DefineOutputList(TList* d) const;
 773     /**
 774      * Get our output list from the container @a d
 775      *
 776      * @param d where to get the output list from
 777      *
 778      * @return The found TList or null
 779      */
 780     TList* GetOutputList(const TList* d) const;
 781     /**
 782      * Find a specific histogram in the source list @a d
 783      *
 784      * @param d     (top)-container
 785      * @param name  Name of histogram
 786      *
 787      * @return Found histogram or null
 788      */
 789     TH1* GetOutputHist(const TList* d, const char* name) const;
 790     /**
 791      * Get the colour of this ring
 792      *
 793      *
 794      * @return
 795      */
 796     Color_t Color() const
 797     {
 798       return AliForwardUtil::RingColor(fDet, fRing);
 799     }
 800     /**
 801      * The name of this ring
 802      *
 803      * @return Name of this ring
 804      */
 805     const char* GetName() const { return fName.Data(); }
 806     UShort_t fDet;   // Detector
 807     Char_t   fRing;  // Ring
 808     TString  fName;  // Name
 809
 810     ClassDef(RingHistos,1)
 811   };
 812   /* @} */
 813
 814   //__________________________________________________________________
 815   /**
 816    * A guard idom for producing debug output
 817    *
 818    */
 819   struct DebugGuard
 820   {
 821     /**
 822      * Constructor
 823      *
 824      * @param lvl       Current level
 825      * @param msgLvl    Target level
 826      * @param format    @c printf -like format
 827      *
 828      * @return
 829      */
 830     DebugGuard(Int_t lvl, Int_t msgLvl, const char* format, ...);
 831     /**
 832      * Destructor
 833      */
 834     ~DebugGuard();
 835     /**
 836      * Make a message
 837      *
 838      * @param lvl    Current level
 839      * @param msgLvl Target level
 840      * @param format @c printf -like format
 841      */
 842     static void Message(Int_t lvl, Int_t msgLvl, const char* format, ...);
 843   private:
 844     /**
 845      * Output the message
 846      *
 847      * @param in    Direction
 848      * @param msg   Message
 849      */
 850     static void Output(int in, TString& msg);
 851     /**
 852      * Format a message
 853      *
 854      * @param out     Output is stored here
 855      * @param format  @c printf -like format
 856      * @param ap      List of arguments
 857      */
 858     static void Format(TString& out, const char* format, va_list ap);
 859     TString fMsg;
 860   };
 861 private:
 862   /**
 863    * Constructor
 864    */
 865   AliForwardUtil() {}
 866   /**
 867    * Copy constructor
 868    *
 869    * @param o Object to copy from
 870    */
 871   AliForwardUtil(const AliForwardUtil& o) : TObject(o) {}
 872   /**
 873    * Assingment operator
 874    *
 875    *
 876    * @return Reference to this object
 877    */
 878   AliForwardUtil& operator=(const AliForwardUtil&) { return *this; }
 879   /**
 880    * Destructor
 881    */
 882   ~AliForwardUtil() {}
 883
 884
 885   ClassDef(AliForwardUtil,1) // Utilities - do not make object
 886 };
 887
 888 // #ifdef LOG_NO_DEBUG
 889 // # define DGUARD(L,N,F,...) do {} while(false)
 890 // #else
 891 /**
 892  * Macro to declare a DebugGuard
 893  *
 894  * @param L Current debug level
 895  * @param N Target debug level
 896  * @param F @c printf -like Format
 897  */
 898 # define DGUARD(L,N,F,...)                                      \
 899   AliForwardUtil::DebugGuard _GUARD(L,N,F, ## __VA_ARGS__)
 900 /**
 901  * Macro to make a debug message, using DebugGuard::Message
 902  *
 903  * @param L Current debug level
 904  * @param N Target debug level
 905  * @param F @c printf -like Format
 906  */
 907 # define DMSG(L,N,F,...)                                        \
 908   AliForwardUtil::DebugGuard::Message(L,N,F, ## __VA_ARGS__)
 909 // #endif
 910 #endif
 911 // Local Variables:
 912 //  mode: C++
 913 // End:
 914