FMD/flow/AliFMDFlowBin.cxx

   1 /* Copyright (C) 2007 Christian Holm Christensen <cholm@nbi.dk>
   2  *
   3  * This library is free software; you can redistribute it and/or
   4  * modify it under the terms of the GNU Lesser General Public License
   5  * as published by the Free Software Foundation; either version 2.1 of
   6  * the License, or (at your option) any later version.
   7  *
   8  * This library is distributed in the hope that it will be useful, but
   9  * WITHOUT ANY WARRANTY; without even the implied warranty of
  10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11  * Lesser General Public License for more details.
  12  *
  13  * You should have received a copy of the GNU Lesser General Public
  14  * License along with this library; if not, write to the Free Software
  15  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
  16  * USA
  17  */
  18 /** @file
  19     @brief implementation of a Bin in a Flow histogram */
  20 //____________________________________________________________________
  21 //
  22 // This contains an of class AliFMDFlowHarmonic and an object of
  23 // class AliFMDFlowEventPlane to calculate v_n and \Psi_k.  It contain
  24 // two objects of class AliFMDFlowEventPlane to calculate the
  25 // sub-event event planes Psi_A, \Psi_B.  It also contain 3 objects of
  26 // class AliFMDFlowResolution to calculate the event plane angle
  27 // resolution.
  28 //
  29 #include "flow/AliFMDFlowBin.h"
  30 #include <cmath>
  31 #include <iostream>
  32 #include <iomanip>
  33 #include <TBrowser.h>
  34
  35 //====================================================================
  36 AliFMDFlowBin::AliFMDFlowBin(const AliFMDFlowBin& o)
  37   : TObject(o),
  38     fPsi(o.fPsi),
  39     fPsiA(o.fPsiA),
  40     fPsiB(o.fPsiB),
  41     fRes(o.fRes),
  42     fResStar(o.fResStar),
  43     fResTdr(o.fResTdr),
  44     fHarmonic(o.fHarmonic)
  45 {
  46   // Copy constructor
  47   // Parameters:
  48   //    o Object to copy from
  49 }
  50
  51 //____________________________________________________________________
  52 AliFMDFlowBin&
  53 AliFMDFlowBin::operator=(const AliFMDFlowBin& o)
  54 {
  55   // Assignment operator
  56   // Parameters:
  57   //    o Object to assign from
  58   fPsi      = o.fPsi;
  59   fPsiA     = o.fPsiA;
  60   fPsiB     = o.fPsiB;
  61   fRes      = o.fRes;
  62   fResStar  = o.fResStar;
  63   fResTdr   = o.fResTdr;
  64   fHarmonic = o.fHarmonic;
  65   return *this;
  66 }
  67
  68 //____________________________________________________________________
  69 void
  70 AliFMDFlowBin::Begin()
  71 {
  72   // Clear event plane calculators
  73   fPsi.Clear();
  74   fPsiA.Clear();
  75   fPsiB.Clear();
  76 }
  77
  78 //____________________________________________________________________
  79 void
  80 AliFMDFlowBin::AddToEventPlane(Double_t phi, Double_t w, Bool_t a)
  81 {
  82   // Called to add a contribution to the event plane
  83   // Parameters
  84   //    phi   The angle phi in [0,2pi]
  85   //    w     Weight
  86   //    a    If true, add to sub-event A, otherwise to sub-event B.
  87   fPsi.Add(phi, w);
  88   if (a) fPsiA.Add(phi, w);
  89   else   fPsiB.Add(phi, w);
  90 }
  91
  92 //____________________________________________________________________
  93 void
  94 AliFMDFlowBin::AddToHarmonic(Double_t phi, Double_t w)
  95 {
  96   // Called to add a contribution to the harmonic.
  97   // Parameters:
  98   //   phi   The angle phi in [0,2pi]
  99   //   w     Weight of phi (only used in the calculation of
 100   //         the event plane).
 101   //
 102   // Disregard the obervation of phi from the event plane angle.
 103   Double_t psi   = fPsi.Psi(phi, w);
 104   fHarmonic.Add(phi, psi);
 105 }
 106
 107 //____________________________________________________________________
 108 void
 109 AliFMDFlowBin::End()
 110 {
 111   // Should be called at the end of an event
 112   Double_t psiA = fPsiA.Psi();
 113   Double_t psiB = fPsiB.Psi();
 114
 115   // Update the resolutions
 116   fRes.Add(psiA, psiB);
 117   fResStar.Add(psiA, psiB);
 118   fResTdr.Add(psiA, psiB);
 119 }
 120
 121 //____________________________________________________________________
 122 void
 123 AliFMDFlowBin::Event(Double_t* phis, Double_t* ws, UInt_t n)
 124 {
 125   // Analyse events
 126   // Parameters
 127   //  phis   Array of phi, (phi_1, ..., phi_n)
 128   //  ws     Weights (optional)
 129   //  n      Size of phis and possibly ws
 130   Begin();
 131
 132   // Calculate split.
 133   UInt_t split = n / 2;
 134   // First sub-event.
 135   for (UInt_t i = 0; i < split; i++)
 136     AddToEventPlane(phis[i], (ws ? ws[i] : 1), kTRUE);
 137   // Second sub-event.
 138   for (UInt_t i = split; i < n; i++)
 139     AddToEventPlane(phis[i], (ws ? ws[i] : 1), kFALSE);
 140   // Add contributions to the harmonic.
 141   for (UInt_t i = 0; i < n; i++)
 142     AddToHarmonic(phis[i], (ws ? ws[i] : 1));
 143
 144   End();
 145 }
 146
 147 //____________________________________________________________________
 148 Double_t
 149 AliFMDFlowBin::Value(CorType t) const
 150 {
 151   // Get the value in this bin
 152   // Parameters:
 153   //   t  Which type of correction
 154   //
 155   // return the value of the harmonic
 156   Double_t e;
 157   return Value(e, t);
 158 }
 159
 160 //____________________________________________________________________
 161 Double_t
 162 AliFMDFlowBin::EValue(CorType t) const
 163 {
 164   // Get the value in this bin
 165   // Parameters:
 166   //    t    Which type of correction
 167   //
 168   // return the error on the value of the harmonic
 169   Double_t e2;
 170   Value(e2, t);
 171   return sqrt(e2);
 172 }
 173
 174 //____________________________________________________________________
 175 Double_t
 176 AliFMDFlowBin::Value(Double_t& e2, CorType t) const
 177 {
 178   // Get the value in this bin
 179   // Parameters:
 180   //    e2   On return, the square error.
 181   //    t    Which type  of correction
 182   //
 183   // return the value of the harmonic
 184   Double_t r, er2;
 185   r = Correction(er2, t);
 186   return fHarmonic.Value(r, er2, e2);
 187 }
 188
 189 //____________________________________________________________________
 190 Double_t
 191 AliFMDFlowBin::Correction(Double_t& er2, CorType t) const
 192 {
 193   // Get the value in this bin
 194   // Parameters:
 195   //    e2   On return, the square error.
 196   //    t    Which type  of correction
 197   //
 198   // return the value of the Correction
 199   Double_t r = 1;
 200   UShort_t k = fHarmonic.Order()/fRes.Order();
 201   switch (t) {
 202   case kNaive: r = fRes.Correction(k, er2);     break;
 203   case kStar:  r = fResStar.Correction(k, er2); break;
 204   case kTdr:   r = fResTdr.Correction(k, er2);  break;
 205   default:     r = 1; er2 = 0;                  break;
 206   }
 207   return r;
 208 }
 209
 210 //____________________________________________________________________
 211 void
 212 AliFMDFlowBin::Finish()
 213 {
 214   // Called at the end of the event
 215 }
 216
 217 //____________________________________________________________________
 218 void
 219 AliFMDFlowBin::Browse(TBrowser* b)
 220 {
 221   // Browse this item
 222   b->Add(&fPsi,      "Full event plane");
 223   b->Add(&fPsiA,     "Sub-event A event plane");
 224   b->Add(&fPsiB,     "Sub-event A event plane");
 225   b->Add(&fRes,      "Naive resolution");
 226   b->Add(&fResStar,  "STAR resolution");
 227   b->Add(&fResTdr,   "TDR resolution");
 228   b->Add(&fHarmonic, "Harmonic");
 229 }
 230
 231 //____________________________________________________________________
 232 void
 233 AliFMDFlowBin::Print(Option_t*) const
 234 {
 235   // Print information
 236   Double_t e2v[4], v[4], r[4], e2r[4];
 237   const char* names[] = { "Bare", "Naive", "STAR", "TDR" };
 238   v[0] = 100 * Value(e2v[0], AliFMDFlowBin::kNone);
 239   v[1] = 100 * Value(e2v[1], AliFMDFlowBin::kNaive);
 240   v[2] = 100 * Value(e2v[2], AliFMDFlowBin::kStar);
 241   v[3] = 100 * Value(e2v[3], AliFMDFlowBin::kTdr);
 242   r[0] = 100 * Correction(e2r[0], AliFMDFlowBin::kNone);
 243   r[1] = 100 * Correction(e2r[1], AliFMDFlowBin::kNaive);
 244   r[2] = 100 * Correction(e2r[2], AliFMDFlowBin::kStar);
 245   r[3] = 100 * Correction(e2r[3], AliFMDFlowBin::kTdr);
 246
 247   std::streamsize         oldP = std::cout.precision(3);
 248   std::ios_base::fmtflags oldF = std::cout.setf(std::ios_base::fixed,
 249                                                 std::ios_base::floatfield);
 250   std::cout << "  v" << std::setw(1) << fHarmonic.Order() << ": ";
 251   for (UInt_t i = 0; i < 4; i++)
 252     std::cout << std::setw(6+(i == 0 ? 0 : 6)) << names[i] << ": "
 253               << std::setw(6) << v[i] << " +/- "
 254               << std::setw(6) << 100*sqrt(e2v[i]) << " ["
 255               << std::setw(7) << r[i] << " +/- "
 256               << std::setw(7) << 100*sqrt(e2r[i]) << "]\n";
 257   std::cout << std::flush;
 258   std::cout.precision(oldP);
 259   std::cout.setf(oldF, std::ios_base::floatfield);
 260 }
 261
 262
 263 //____________________________________________________________________
 264 //
 265 // EOF
 266 //