PWG2/FORWARD/analysis2/AliFMDCorrELossFit.h

   1 // Object holding the Energy loss fit 'correction'
   2 //
   3 // These are generated from Monte-Carlo or real ESDs.
   4 #ifndef ALIFMDCORRELOSSFIT_H
   5 #define ALIFMDCORRELOSSFIT_H
   6 #include <TObject.h>
   7 #include <TAxis.h>
   8 #include <TObjArray.h>
   9 class TF1;
  10 class TBrowser;
  11
  12 /**
  13  * @defgroup pwg2_forward_corr Corrections
  14  *
  15  * @ingroup pwg2_forward
  16  */
  17 /**
  18  * Object holding the Energy loss fit 'correction'
  19  *
  20  * These are generated from Monte-Carlo or real ESDs.
  21  *
  22  * @ingroup pwg2_forward_corr
  23  */
  24 class AliFMDCorrELossFit : public TObject
  25 {
  26 public:
  27   /**
  28    * POD structure to hold data from fits
  29    *
  30    * @ingroup pwg2_forward_corr
  31    */
  32   struct ELossFit : public TObject
  33   {
  34     Int_t     fN;      // Number of peaks fitted
  35     UShort_t  fNu;     // Number degrees of freedom
  36     Double_t  fChi2;   // Chi square from fit
  37     Double_t  fC;      // Scaling constant
  38     Double_t  fDelta;  // Most probable value
  39     Double_t  fXi;     // Width parameter of Landau
  40     Double_t  fSigma;  // Sigma on folded gaussian
  41     Double_t  fSigmaN; // Sigma of detector noise
  42     Double_t* fA;      // [fN] Weights
  43     Double_t  fEC;     // Error on C
  44     Double_t  fEDelta; // Error on Delta
  45     Double_t  fEXi;    // Error on Xi
  46     Double_t  fESigma; // Error on sigma
  47     Double_t  fESigmaN;// Error on sigma (noise)
  48     Double_t* fEA;     // [fN] Error on weights
  49     Int_t     fQuality;// Assigned quality
  50     UShort_t  fDet;    // Detector
  51     Char_t    fRing;   // Ring
  52     UShort_t  fBin;    // Eta bin
  53
  54     static Double_t fgMaxRelError;  // Global default max relative error
  55     static Double_t fgLeastWeight;  // Global default least weight
  56     static Double_t fgMaxChi2nu;    // Global default maximum reduced chi^2
  57     /**
  58      * Default constructor
  59      *
  60      */
  61     ELossFit();
  62     /**
  63      * Construct from a function
  64      *
  65      * @param quality Quality flag
  66      * @param f       Function
  67      */
  68     ELossFit(Int_t quality,const TF1& f);
  69     /**
  70      * Constructor with full parameter set
  71      *
  72      * @param quality   Quality flag
  73      * @param n         @f$ N@f$ - Number of fitted peaks
  74      * @param chi2      @f$ \chi^2 @f$
  75      * @param nu        @f$ \nu @f$ - number degrees of freedom
  76      * @param c         @f$ C@f$ - scale constant
  77      * @param ec        @f$ \delta C@f$ - error on @f$ C@f$
  78      * @param delta     @f$ \Delta@f$ - Most probable value
  79      * @param edelta    @f$ \delta\Delta@f$ - error on @f$\Delta@f$
  80      * @param xi        @f$ \xi@f$ - width
  81      * @param exi       @f$ \delta\xi@f$ - error on @f$\xi@f$
  82      * @param sigma     @f$ \sigma@f$ - Width of Gaussian
  83      * @param esigma    @f$ \delta\sigma@f$ - error on @f$\sigma@f$
  84      * @param sigman    @f$ \sigma_n@f$ - Noise width
  85      * @param esigman   @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$
  86      * @param a         Array of @f$ N-1@f$ weights @f$ a_i@f$ for
  87      *                  @f$ i=2,\ldots@f$
  88      * @param ea        Array of @f$ N-1@f$ error on the weights @f$ a_i@f$ for
  89      *                  @f$ i=2,\ldots@f$
  90      */
  91     ELossFit(Int_t     quality,UShort_t  n,
  92              Double_t  chi2,   UShort_t  nu,
  93              Double_t  c,      Double_t  ec,
  94              Double_t  delta,  Double_t  edelta,
  95              Double_t  xi,     Double_t  exi,
  96              Double_t  sigma,  Double_t  esigma,
  97              Double_t  sigman, Double_t  esigman,
  98              Double_t* a,      Double_t* ea);
  99     /**
 100      * Copy constructor
 101      *
 102      * @param o Object to copy from
 103      */
 104     ELossFit(const ELossFit& o);
 105     /**
 106      * Assignment operator
 107      *
 108      * @param o Object to assign from
 109      *
 110      * @return Reference to this object
 111      */
 112     ELossFit& operator=(const ELossFit& o);
 113     /**
 114      * Destructor
 115      */
 116     ~ELossFit();
 117     /**
 118      * @{
 119      * @name Access to parameters
 120      */
 121     /**
 122      * @return Number of peaks fitted
 123      */
 124     Int_t GetN() const { return fN; }
 125     /**
 126      * @return Number degrees of freedom
 127      */
 128     UShort_t GetNu() const { return fNu; }
 129     /**
 130      * @return Chi square from fit
 131      */
 132     Double_t GetChi2() const { return fChi2; }
 133     /**
 134      * @return Scaling constant
 135      */
 136     Double_t GetC() const { return fC; }
 137     /**
 138      * @return Most probable value
 139      */
 140     Double_t GetDelta() const { return fDelta; }
 141     /**
 142      * @return Width parameter of Landau
 143      */
 144     Double_t GetXi() const { return fXi; }
 145     /**
 146      * @return Sigma on folded gaussian
 147      */
 148     Double_t GetSigma() const { return fSigma; }
 149     /**
 150      * @return Sigma of detector noise
 151      */
 152     Double_t GetSigmaN() const { return fSigmaN; }
 153     /**
 154      * @return Weights
 155      */
 156     Double_t* GetAs() const { return fA; }
 157     /**
 158      * @return Weights
 159      */
 160     Double_t GetA(UShort_t i) const;
 161     /**
 162      * @return Error on C
 163      */
 164     Double_t GetEC() const { return fEC; }
 165     /**
 166      * @return Error on Delta
 167      */
 168     Double_t GetEDelta() const { return fEDelta; }
 169     /**
 170      * @return Error on Xi
 171      */
 172     Double_t GetEXi() const { return fEXi; }
 173     /**
 174      * @return Error on sigma
 175      */
 176     Double_t GetESigma() const { return fESigma; }
 177     /**
 178      * @return Error on sigma (noise)
 179      */
 180     Double_t GetESigmaN() const { return fESigmaN; }
 181     /**
 182      * @return Error on weights
 183      */
 184     Double_t* GetEAs() const { return fEA; }
 185     /**
 186      * @return Error on weights
 187      */
 188     Double_t GetEA(UShort_t i) const;
 189     /**
 190      * @return Assigned quality
 191      */
 192     Int_t GetQuality() const { return fQuality; }
 193     /**
 194      * @return Detector
 195      */
 196     UShort_t GetDet() const { return fDet; }
 197     /**
 198      * @return Ring
 199      */
 200     Char_t GetRing() const { return fRing; }
 201     /**
 202      * @return Eta bin
 203      */
 204     UShort_t GetBin() const { return fBin; }
 205     /* @} */
 206
 207     /**
 208      * @{
 209      * @name Evaluation
 210      */
 211     /**
 212      * Evaluate
 213      * @f[
 214      *  f_N(x;\Delta,\xi,\sigma') =
 215      *     \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i')
 216      * @f]
 217      *
 218      * (see AliForwardUtil::NLandauGaus) for the maximum @f$ N @f$
 219      * that fulfills the requirements
 220      *
 221      * @param x           Where to evaluate
 222      * @param maxN        @f$ \max{N}@f$
 223      *
 224      * @return @f$ f_N(x;\Delta,\xi,\sigma')@f$
 225      */
 226     Double_t Evaluate(Double_t x,
 227                       UShort_t maxN=999) const;
 228     /**
 229      * Evaluate
 230      * @f[
 231      *   f_W(x;\Delta,\xi,\sigma') =
 232      *   \frac{\sum_{i=1}^{n} i a_i f_i(x;\Delta,\xi,\sigma')}{
 233      *     f_N(x;\Delta,\xi,\sigma')} =
 234      *   \frac{\sum_{i=1}^{n} i a_i f(x;\Delta_i,\xi_i,\sigma_i')}{
 235      *     \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i')}
 236      * @f]
 237      * where @f$ n@f$ fulfills the requirements (see FindMaxWeight).
 238      *
 239      * If the denominator is zero, then 1 is returned.
 240      *
 241      * See also AliForwardUtil::ILandauGaus and AliForwardUtil::NLandauGaus
 242      * for more information on the evaluated functions.
 243      *
 244      * @param x           Where to evaluate
 245      * @param maxN        @f$ \max{N}@f$
 246      *
 247      * @return @f$ f_W(x;\Delta,\xi,\sigma')@f$.
 248      */
 249     Double_t EvaluateWeighted(Double_t x,
 250                               UShort_t maxN=9999) const;
 251     /**
 252      * Find the maximum weight to use.  The maximum weight is the
 253      * largest i for which
 254      *
 255      * - @f$ i \leq \max{N}@f$
 256      * - @f$ a_i > \min{a}@f$
 257      * - @f$ \delta a_i/a_i > \delta_{max}@f$
 258      *
 259      * @param maxRelError @f$ \min{a}@f$
 260      * @param leastWeight @f$ \delta_{max}@f$
 261      * @param maxN        @f$ \max{N}@f$
 262      *
 263      * @return The largest index @f$ i@f$ for which the above
 264      * conditions hold.  Will never return less than 1.
 265      */
 266     Int_t FindMaxWeight(Double_t maxRelError=2*fgMaxRelError,
 267                         Double_t leastWeight=fgLeastWeight,
 268                         UShort_t maxN=999) const;
 269     /* @} */
 270     /**
 271      * @{
 272      * @name TObject Sortable interface
 273      */
 274     /**
 275      * Declare this object as sortable
 276      *
 277      * @return Always true
 278      */
 279     Bool_t IsSortable() const { return kTRUE; }
 280     /**
 281      * Compare to another ELossFit object.
 282      *
 283      * - +1, if this quality is better (larger) than other objects quality
 284      * - -1, if this quality is worse (smaller) than other objects quality
 285      * - +1, if this @f$|\chi^2/\nu-1|@f$ is smaller than the same for other
 286      * - -1, if this @f$|\chi^2/\nu-1|@f$ is larger than the same for other
 287      * - 0 otherwise
 288      *
 289      * @param o Other object to compare to
 290      */
 291     Int_t Compare(const TObject* o) const;
 292     /* @} */
 293     /**
 294      * @{
 295      * name Auxiliary member functions
 296      */
 297     /**
 298      * Information to standard output
 299      *
 300      * @param option Not used
 301      */
 302     void Print(Option_t* option) const; // *MENU*
 303     /**
 304      * Draw this fit
 305      *
 306      * @param option Options
 307      *  - COMP  Draw components too
 308      */
 309     void Draw(Option_t* option="comp"); // *MENU*
 310     /**
 311      * Browse this object
 312      *
 313      * @param b Browser
 314      */
 315     void Browse(TBrowser* b);
 316     /**
 317      * Get the name of this object
 318      *
 319      *
 320      * @return
 321      */
 322     const Char_t* GetName() const;
 323     /**
 324      * Calculate the quality
 325      */
 326     void CalculateQuality(Double_t maxChi2nu=fgMaxChi2nu,
 327                           Double_t maxRelError=fgMaxRelError,
 328                           Double_t leastWeight=fgLeastWeight);
 329     /* @} */
 330     ClassDef(ELossFit,1); // Result of fit
 331   };
 332
 333   /**
 334    * Default constructor
 335    */
 336   AliFMDCorrELossFit();
 337   /**
 338    * Copy constructor
 339    *
 340    * @param o Object to copy from
 341    */
 342   AliFMDCorrELossFit(const AliFMDCorrELossFit& o);
 343   /**
 344    * Destructor
 345    */
 346   virtual ~AliFMDCorrELossFit();
 347   /**
 348    * Assignment operator
 349    *
 350    * @param o Object to assign from
 351    *
 352    * @return Reference to this object
 353    */
 354   AliFMDCorrELossFit& operator=(const AliFMDCorrELossFit& o);
 355
 356   /**
 357    * @{
 358    * @name Set fits
 359    */
 360   /**
 361    * Set the fit parameters from a function
 362    *
 363    * @param d        Detector
 364    * @param r        Ring
 365    * @param eta      Eta
 366    * @param quality  Quality flag
 367    * @param f        Function from fit
 368    */
 369   Bool_t SetFit(UShort_t d, Char_t r, Double_t eta, Int_t quality,
 370                 const TF1& f);
 371   /**
 372    * Set the fit parameters from a function
 373    *
 374    * @param d    Detector
 375    * @param r    Ring
 376    * @param eta  Eta
 377    * @param f    ELoss fit result - note, the object will take ownership
 378    */
 379   Bool_t SetFit(UShort_t d, Char_t r, Double_t eta, ELossFit* f);
 380   /**
 381    * Set the fit parameters from a function
 382    *
 383    * @param d       Detector
 384    * @param r       Ring
 385    * @param etaBin  Eta (bin number, 1->nBins)
 386    * @param f       ELoss fit result - note, the object will take ownership
 387    */
 388   Bool_t SetFit(UShort_t d, Char_t r, Int_t etaBin, ELossFit* f);
 389   /**
 390    * Set the fit parameters from a function
 391    *
 392    * @param d         Detector number
 393    * @param r         Ring identifier
 394    * @param eta       Eta value
 395    * @param quality   Quality flag
 396    * @param n         @f$ N@f$ - Number of fitted peaks
 397    * @param chi2      @f$ \chi^2 @f$
 398    * @param nu        @f$ \nu @f$ - number degrees of freedom
 399    * @param c         @f$ C@f$ - scale constant
 400    * @param ec        @f$ \delta C@f$ - error on @f$ C@f$
 401    * @param delta     @f$ \Delta@f$ - most probable value
 402    * @param edelta    @f$ \delta\Delta@f$ - error on @f$\Delta@f$
 403    * @param xi        @f$ \xi@f$ - Landau width
 404    * @param exi       @f$ \delta\xi@f$ - error on @f$\xi@f$
 405    * @param sigma     @f$ \sigma@f$ - Gaussian width
 406    * @param esigma    @f$ \delta\sigma@f$ - error on @f$\sigma@f$
 407    * @param sigman    @f$ \sigma_n@f$ - Noise width
 408    * @param esigman   @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$
 409    * @param a         Array of @f$ N-1@f$ weights @f$ a_i@f$ for
 410    *                  @f$ i=2,\ldots@f$
 411    * @param ea        Array of @f$ N-1@f$ errors on weights @f$ a_i@f$ for
 412    *                  @f$ i=2,\ldots@f$
 413    */
 414   Bool_t SetFit(UShort_t  d,      Char_t    r, Double_t eta,
 415                 Int_t     quality,UShort_t  n,
 416                 Double_t  chi2,   UShort_t  nu,
 417                 Double_t  c,      Double_t  ec,
 418                 Double_t  delta,  Double_t  edelta,
 419                 Double_t  xi,     Double_t  exi,
 420                 Double_t  sigma,  Double_t  esigma,
 421                 Double_t  sigman, Double_t  esigman,
 422                 Double_t* a,      Double_t* ea);
 423   /* @} */
 424
 425   /**
 426    * @{
 427    * @name Set and get eta axis
 428    */
 429   /**
 430    * Set the eta axis to use
 431    *
 432    * @param axis Eta axis
 433    */
 434   void SetEtaAxis(const TAxis& axis);
 435   /**
 436    * Set the eta axis to use
 437    *
 438    * @param nBins Number of bins
 439    * @param min   Minimum @f$ \eta@f$
 440    * @param max   maximum @f$ \eta@f$
 441    */
 442   void SetEtaAxis(Int_t nBins, Double_t min, Double_t max);
 443   /**
 444    * Get the eta axis used
 445    *
 446    * @return
 447    */
 448   const TAxis& GetEtaAxis() const { return fEtaAxis; }
 449   /**
 450    * Set the low cut used when fitting
 451    *
 452    * @param cut Cut value
 453    */
 454   void SetLowCut(Double_t cut) { fLowCut = cut; }
 455   /**
 456    * Get the low cut used when fitting
 457    *
 458    * @return Low cut used for fitting
 459    */
 460   Double_t GetLowCut() const { return fLowCut; }
 461   /**
 462    * Find the eta bin corresponding to the given eta
 463    *
 464    * @param eta  Eta value
 465    *
 466    * @return Bin (in @f$[1,N_{bins}]@f$) corresponding to the given
 467    * eta, or 0 if out of range.
 468    */
 469   Int_t FindEtaBin(Double_t eta) const;
 470   /* @} */
 471
 472   /**
 473    * @{
 474    * @name Find fits
 475    */
 476   /**
 477    * Find the fit corresponding to the specified parameters
 478    *
 479    * @param d   Detector
 480    * @param r   Ring
 481    * @param eta Eta value
 482    *
 483    * @return Fit parameters or null in case of problems
 484    */
 485   ELossFit* FindFit(UShort_t d, Char_t r, Double_t eta) const;
 486   /**
 487    * Find the fit corresponding to the specified parameters
 488    *
 489    * @param d      Detector
 490    * @param r      Ring
 491    * @param etabin Eta bin (1 based)
 492    *
 493    * @return Fit parameters or null in case of problems
 494    */
 495   ELossFit* FindFit(UShort_t d, Char_t r, Int_t etabin) const;
 496   /* @} */
 497
 498   /**
 499    * @{
 500    * @name Miscellaneous
 501    */
 502   /**
 503    * Get the ring array corresponding to the specified ring
 504    *
 505    * @param d Detector
 506    * @param r Ring
 507    *
 508    * @return Pointer to ring array, or null in case of problems
 509    */
 510   TObjArray* GetRingArray(UShort_t d, Char_t r) const;
 511   /**
 512    * Signal that this is a folder
 513    *
 514    * @return Always true
 515    */
 516   Bool_t IsFolder() const { return true; }
 517   /**
 518    * Browse this object
 519    *
 520    * @param b
 521    */
 522   void Browse(TBrowser* b);
 523   /**
 524    * Draw this object
 525    *
 526    * @param option Options.  Possible values are
 527    *  - err Plot error bars
 528    */
 529   void Draw(Option_t* option=""); //*MENU*
 530   /**
 531    * Print this object.
 532    *
 533    * @param option Options
 534    *   - R   Print recursive
 535    *
 536    */
 537   void Print(Option_t* option="R") const; //*MENU*
 538   /* @} */
 539 protected:
 540   /**
 541    * Get the ring array corresponding to the specified ring
 542    *
 543    * @param d Detector
 544    * @param r Ring
 545    *
 546    * @return Pointer to ring array, or newly created container
 547    */
 548   TObjArray* GetOrMakeRingArray(UShort_t d, Char_t r);
 549
 550   TObjArray  fRings;    // Array of rings
 551   TAxis      fEtaAxis;  // Eta axis used
 552   Double_t   fLowCut;   // Low cut used when fitting
 553
 554   ClassDef(AliFMDCorrELossFit,2);
 555 };
 556
 557 //____________________________________________________________________
 558 inline void
 559 AliFMDCorrELossFit::SetEtaAxis(Int_t nBins, Double_t min, Double_t max)
 560 {
 561   fEtaAxis.Set(nBins, min, max);
 562 }
 563 //____________________________________________________________________
 564 inline void
 565 AliFMDCorrELossFit::SetEtaAxis(const TAxis& e)
 566 {
 567   fEtaAxis.Set(e.GetNbins(), e.GetXmin(), e.GetXmax());
 568 }
 569 //____________________________________________________________________
 570 inline Double_t
 571 AliFMDCorrELossFit::ELossFit::GetA(UShort_t i) const
 572 {
 573   if (i <  1)   return 0;
 574   if (i >  fN)  return 0;
 575   if (i == 1)   return 1;
 576   return fA[i-2];
 577 }
 578 //____________________________________________________________________
 579 inline Double_t
 580 AliFMDCorrELossFit::ELossFit::GetEA(UShort_t i) const
 581 {
 582   if (i <  1)   return 0;
 583   if (i >  fN)  return 0;
 584   if (i == 1)   return 1;
 585   return fEA[i-2];
 586 }
 587
 588
 589 #endif
 590 // Local Variables:
 591 //   mode: C++
 592 // End:
 593