PWG2/FORWARD/analysis2/AliAODForwardMult.h

   1 //
   2 // See implementation or Doxygen comments for more information
   3 //
   4 #ifndef ALIAODFORWARDMULT_H
   5 #define ALIAODFORWARDMULT_H
   6 /**
   7  * @file   AliAODForwardMult.h
   8  * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk>
   9  * @date   Wed Mar 23 13:58:00 2011
  10  *
  11  * @brief
  12  *
  13  * @ingroup pwg2_forward_aod
  14  *
  15  */
  16 #include <TObject.h>
  17 #include <TH2D.h>
  18 class TBrowser;
  19 class TH1I;
  20 /**
  21  * Class that contains the forward multiplicity data per event
  22  *
  23  * This class contains a histogram of
  24  * @f[
  25  *   \frac{d^2N_{ch}}{d\eta d\phi}\quad,
  26  * @f]
  27  * as well as a trigger mask for each analysed event.
  28  *
  29  * The eta acceptance of the event is stored in the underflow bins of
  30  * the histogram.  So to build the final histogram, one needs to
  31  * correct for this acceptance (properly weighted by the events), and
  32  * the vertex efficiency.  This simply boils down to defining a 2D
  33  * histogram and summing the event histograms in that histogram.  One
  34  * should of course also do proper book-keeping of the accepted event.
  35  *
  36  * @code
  37  * TTree* GetAODTree()
  38  * {
  39  *    TFile* file = TFile::Open("AliAODs.root","READ");
  40  *    TTree* tree = static_cast<TTree*>(file->Get("aodTree"));
  41  *    return tree;
  42  * }
  43  *
  44  * void Analyse()
  45  * {
  46  *   TH2D*              sum        = 0;                  // Summed hist
  47  *   TTree*             tree       = GetAODTree();       // AOD tree
  48  *   AliAODForwardMult* mult       = 0;                  // AOD object
  49  *   Int_t              nTriggered = 0;                  // # of triggered ev.
  50  *   Int_t              nWithVertex= 0;                  // # of ev. w/vertex
  51  *   Int_t              nAccepted  = 0;                  // # of ev. used
  52  *   Int_t              nAvailable = tree->GetEntries(); // How many entries
  53  *   Float_t            vzLow      = -10;                // Lower ip cut
  54  *   Float_t            vzHigh     =  10;                // Upper ip cut
  55  *   Int_t              mask       = AliAODForwardMult::kInel;// Trigger mask
  56  *   tree->SetBranchAddress("forward", &forward);        // Set the address
  57  *
  58  *   for (int i = 0; i < nAvailable; i++) {
  59  *     // Create sum histogram on first event - to match binning to input
  60  *     if (!sum) sum = static_cast<TH2D*>(mult->Clone("d2ndetadphi"));
  61  *
  62  *     tree->GetEntry(i);
  63  *
  64  *     // Other trigger/event requirements could be defined
  65  *     if (!mult->IsTriggerBits(mask)) continue;
  66  *     nTriggered++;
  67  *
  68  *     // Check if we have vertex
  69  *     if (!mult->HasIpZ()) continue;
  70  *     nWithVertex++;
  71  *
  72  *     // Select vertex range (in centimeters)
  73  *     if (!mult->InRange(vzLow, vzHigh) continue;
  74  *     nAccepted++;
  75  *
  76  *     // Add contribution from this event
  77  *     sum->Add(&(mult->GetHistogram()));
  78  *   }
  79  *
  80  *   // Get acceptance normalisation from underflow bins
  81  *   TH1D* norm   = sum->Projection("norm", 0, 1, "");
  82  *   // Project onto eta axis - _ignoring_underflow_bins_!
  83  *   TH1D* dndeta = sum->Projection("dndeta", 1, -1, "e");
  84  *   // Normalize to the acceptance
  85  *   dndeta->Divide(norm);
  86  *   // Scale by the vertex efficiency
  87  *   dndeta->Scale(Double_t(nWithVertex)/nTriggered, "width");
  88  *   // And draw the result
  89  *   dndeta->Draw();
  90  * }
  91  * @endcode
  92  *
  93  * The above code will draw the final @f$ dN_{ch}/d\eta@f$ for the
  94  * selected event class and vertex range
  95  *
  96  * The histogram can be used as input for other kinds of analysis too,
  97  * like flow, event-plane, centrality, and so on.
  98  *
  99  * @ingroup pwg2_forward
 100  * @ingroup pwg2_forward_aod
 101  */
 102 class AliAODForwardMult : public TObject
 103 {
 104 public:
 105   /**
 106    * Bits of the trigger pattern
 107    */
 108   enum {
 109     /** In-elastic collision */
 110     kInel     = 0x001,
 111     /** In-elastic collision with at least one SPD tracklet */
 112     kInelGt0  = 0x002,
 113     /** Non-single diffractive collision */
 114     kNSD      = 0x004,
 115     /** Empty bunch crossing */
 116     kEmpty    = 0x008,
 117     /** A-side trigger */
 118     kA        = 0x010,
 119     /** B(arrel) trigger */
 120     kB        = 0x020,
 121     /** C-side trigger */
 122     kC        = 0x080,
 123     /** Empty trigger */
 124     kE        = 0x100,
 125     /** pileup from SPD */
 126     kPileUp   = 0x200,
 127     /** true NSD from MC */
 128     kMCNSD    = 0x400,
 129     /** Offline MB triggered */
 130     kOffline  = 0x800
 131   };
 132   /**
 133    * Bin numbers in trigger histograms
 134    */
 135   enum {
 136     kBinAll=1,
 137     kBinInel,
 138     kBinInelGt0,
 139     kBinNSD,
 140     kBinA,
 141     kBinB,
 142     kBinC,
 143     kBinE,
 144     kBinPileUp,
 145     kBinMCNSD,
 146     kBinOffline,
 147     kWithTrigger,
 148     kWithVertex,
 149     kAccepted
 150   };
 151   /**
 152    * Default constructor
 153    *
 154    * Used by ROOT I/O sub-system - do not use
 155    */
 156   AliAODForwardMult();
 157   /**
 158    * Constructor
 159    *
 160    * @param isMC Whether this was from MC or not
 161    */
 162   AliAODForwardMult(Bool_t isMC);
 163   /**
 164    * Destructor
 165    */
 166   virtual ~AliAODForwardMult() {} // Destructor
 167   /**
 168    * Initialize
 169    *
 170    * @param etaAxis  Pseudo-rapidity axis
 171    */
 172   void Init(const TAxis& etaAxis);
 173   /**
 174    * Get the @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
 175    *
 176    * @return @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
 177    */
 178   const TH2D& GetHistogram() const { return fHist; } // Get histogram
 179   /**
 180    * Get the @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
 181    *
 182    * @return @f$ d^2N_{ch}/d\eta d\phi@f$ histogram,
 183    */
 184   TH2D& GetHistogram() { return fHist; } // Get histogram
 185   /**
 186    * Get the trigger mask
 187    *
 188    * @return Trigger mask
 189    */
 190   UInt_t GetTriggerMask() const { return fTriggers; } // Get triggers
 191   /**
 192    * Set the trigger mask
 193    *
 194    * @param trg Trigger mask
 195    */
 196   void SetTriggerMask(UInt_t trg) { fTriggers = trg; } // Set triggers
 197   /**
 198    * Set bit(s) in trigger mask
 199    *
 200    * @param bits bit(s) to set
 201    */
 202   void SetTriggerBits(UInt_t bits) { fTriggers |= bits; } // Set trigger bits
 203   /**
 204    * Check if bit(s) are set in the trigger mask
 205    *
 206    * @param bits Bits to test for
 207    *
 208    * @return
 209    */
 210   Bool_t IsTriggerBits(UInt_t bits) const;
 211   /**
 212    * Whether we have any trigger bits
 213    */
 214   Bool_t HasTrigger() const { return fTriggers != 0; } // Check for triggers
 215   /**
 216    * Clear all data
 217    *
 218    * @param option  Passed on to TH2::Reset verbatim
 219    */
 220   void Clear(Option_t* option="");
 221   /**
 222    * browse this object
 223    *
 224    * @param b Browser
 225    */
 226   void Browse(TBrowser* b);
 227   /**
 228    * This is a folder
 229    *
 230    * @return Always true
 231    */
 232   Bool_t IsFolder() const { return kTRUE; } // Always true
 233   /**
 234    * Print content
 235    *
 236    * @param option Passed verbatim to TH2::Print
 237    */
 238   void Print(Option_t* option="") const;
 239   /**
 240    * Set the z coordinate of the interaction point
 241    *
 242    * @param ipZ Interaction point z coordinate
 243    */
 244   void SetIpZ(Float_t ipZ) { fIpZ = ipZ; } // Set Ip's Z coordinate
 245   /**
 246    * Set the center of mass energy per nucleon-pair.  This is stored
 247    * in the (0,0) of the histogram
 248    *
 249    * @param sNN Center of mass energy per nucleon pair (GeV)
 250    */
 251   void SetSNN(UShort_t sNN);
 252   /**
 253    * Get the collision system number
 254    * - 0: Unknown
 255    * - 1: pp
 256    * - 2: PbPb
 257    *
 258    * @param sys Collision system number
 259    */
 260   void SetSystem(UShort_t sys);
 261   /**
 262    * Set the event centrality
 263    *
 264    * @param c Centrality
 265    */
 266   void SetCentrality(Float_t c) { fCentrality = c; }
 267   /**
 268    * Set the z coordinate of the interaction point
 269    *
 270    * @return Interaction point z coordinate
 271    */
 272   Float_t GetIpZ() const { return fIpZ; } // Get Ip's Z coordinate
 273   /**
 274    * Check if we have a valid z coordinate of the interaction point
 275    *
 276    * @return True if we have a valid interaction point z coordinate
 277    */
 278   Bool_t HasIpZ() const;
 279   /**
 280    * Get the center of mass energy per nucleon pair (GeV)
 281    *
 282    * @return Center of mass energy per nucleon pair (GeV)
 283    */
 284   UShort_t GetSNN() const;
 285   /**
 286    * Get the collision system number
 287    * - 0: Unknown
 288    * - 1: pp
 289    * - 2: PbPb
 290    *
 291    * @return Collision system number
 292    */
 293   UShort_t GetSystem() const;
 294   /**
 295    * Check if the z coordinate of the interaction point is within the
 296    * given limits.  Note that the convention used corresponds to the
 297    * convention used in ROOTs TAxis.
 298    *
 299    * @param low  Lower cut (inclusive)
 300    * @param high Upper cut (exclusive)
 301    *
 302    * @return true if @f$ low \ge ipz < high@f$
 303    */
 304   Bool_t InRange(Float_t low, Float_t high) const;
 305   /**
 306    * Get the event centrality
 307    *
 308    *
 309    * @return
 310    */
 311   Float_t GetCentrality() const { return fCentrality; }
 312   /**
 313    * Check if we have a valid centrality
 314    *
 315    *
 316    * @return
 317    */
 318   Bool_t  HasCentrality() const { return !(fCentrality  < 0); }
 319
 320   /**
 321    * Get the name of the object
 322    *
 323    * @return Name of object
 324    */
 325   const Char_t* GetName() const { return (fIsMC ? "ForwardMC" : "Forward"); }
 326   /**
 327    * Check if event meets the passses requirements.
 328    *
 329    * It returns true if @e all of the following is true
 330    *
 331    * - The trigger is within the bit mask passed.
 332    * - The vertex is within the specified limits.
 333    * - The centrality is within the specified limits, or if lower
 334    *   limit is equal to or larger than the upper limit.
 335    *
 336    * If a histogram is passed in the last parameter, then that
 337    * histogram is filled with the trigger bits.
 338    *
 339    * @param triggerMask  Trigger mask
 340    * @param vzMin        Minimum @f$ v_z@f$ (in centimeters)
 341    * @param vzMax        Maximum @f$ v_z@f$ (in centimeters)
 342    * @param cMin         Minimum centrality (in percent)
 343    * @param cMax         Maximum centrality (in percent)
 344    * @param hist         Histogram to fill
 345    *
 346    * @return @c true if the event meets the requirements
 347    */
 348   Bool_t CheckEvent(Int_t    triggerMask=kInel,
 349                     Double_t vzMin=-10, Double_t vzMax=10,
 350                     UShort_t cMin=0,    UShort_t cMax=100,
 351                     TH1*     hist=0) const;
 352   /**
 353    * Get a string correspondig to the trigger mask
 354    *
 355    * @param mask Trigger mask
 356    *
 357    * @return Static string (copy before use)
 358    */
 359   static const Char_t* GetTriggerString(UInt_t mask);
 360   /**
 361    * Make a histogram to record triggers in.
 362    *
 363    * The bins defined by the trigger enumeration in this class.  One
 364    * can use this enumeration to retrieve the number of triggers for
 365    * each class.
 366    *
 367    * @param name Name of the histogram
 368    *
 369    * @return Newly allocated histogram
 370    */
 371   static TH1I* MakeTriggerHistogram(const char* name="triggers");
 372 protected:
 373   Bool_t  fIsMC;     // Whether this is from MC
 374   TH2D    fHist;     // Histogram of d^2N_{ch}/(deta dphi) for this event
 375   UInt_t  fTriggers; // Trigger bit mask
 376   Float_t fIpZ;      // Z coordinate of the interaction point
 377   Float_t fCentrality; // Event centrality
 378
 379   static const Float_t fgkInvalidIpZ; // Invalid IpZ value
 380   ClassDef(AliAODForwardMult,2); // AOD forward multiplicity
 381 };
 382
 383 //____________________________________________________________________
 384 inline Bool_t
 385 AliAODForwardMult::InRange(Float_t low, Float_t high) const
 386 {
 387   return HasIpZ() && fIpZ >= low && fIpZ < high;
 388 }
 389
 390 //____________________________________________________________________
 391 inline Bool_t
 392 AliAODForwardMult::IsTriggerBits(UInt_t bits) const
 393 {
 394   return HasTrigger() && ((fTriggers & bits) == bits);
 395 }
 396
 397 #endif
 398 // Local Variables:
 399 //  mode: C++
 400 // End:
 401