FastJet 3.0.0
BackgroundEstimatorBase.hh
00001 #ifndef __FASTJET_BACKGROUND_ESTIMATOR_BASE_HH__
00002 #define __FASTJET_BACKGROUND_ESTIMATOR_BASE_HH__
00003 
00004 //STARTHEADER
00005 // $Id: BackgroundEstimatorBase.hh 2616 2011-09-30 18:03:40Z salam $
00006 //
00007 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
00008 //
00009 //----------------------------------------------------------------------
00010 // This file is part of FastJet.
00011 //
00012 //  FastJet is free software; you can redistribute it and/or modify
00013 //  it under the terms of the GNU General Public License as published by
00014 //  the Free Software Foundation; either version 2 of the License, or
00015 //  (at your option) any later version.
00016 //
00017 //  The algorithms that underlie FastJet have required considerable
00018 //  development and are described in hep-ph/0512210. If you use
00019 //  FastJet as part of work towards a scientific publication, please
00020 //  include a citation to the FastJet paper.
00021 //
00022 //  FastJet is distributed in the hope that it will be useful,
00023 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00024 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00025 //  GNU General Public License for more details.
00026 //
00027 //  You should have received a copy of the GNU General Public License
00028 //  along with FastJet. If not, see <http://www.gnu.org/licenses/>.
00029 //----------------------------------------------------------------------
00030 //ENDHEADER
00031 
00032 #include <fastjet/ClusterSequenceAreaBase.hh>
00033 #include <fastjet/FunctionOfPseudoJet.hh>
00034 #include <fastjet/Selector.hh>
00035 #include <fastjet/Error.hh>
00036 #include <iostream>
00037 
00038 FASTJET_BEGIN_NAMESPACE     // defined in fastjet/internal/base.hh
00039 
00040 
00041 /// @ingroup tools_background
00042 /// \class BackgroundEstimatorBase
00043 ///
00044 /// Abstract base class that provides the basic interface for classes
00045 /// that estimate levels of background radiation in hadrion and
00046 /// heavy-ion collider events.
00047 ///
00048 ///
00049 class BackgroundEstimatorBase {
00050 public:
00051   /// @name  constructors and destructors
00052   //\{
00053   //----------------------------------------------------------------
00054   BackgroundEstimatorBase() : _rescaling_class(0) {}
00055   //\}
00056 
00057   /// a default virtual destructor that does nothing
00058   virtual ~BackgroundEstimatorBase() {}
00059 
00060 
00061   /// @name setting a new event
00062   //\{
00063   //----------------------------------------------------------------
00064 
00065   /// tell the background estimator that it has a new event, composed
00066   /// of the specified particles.
00067   virtual void set_particles(const std::vector<PseudoJet> & particles) = 0;
00068 
00069   //\}
00070 
00071   /// @name  retrieving fundamental information
00072   //\{
00073   //----------------------------------------------------------------
00074 
00075   /// get rho, the background density per unit area
00076   virtual double rho() const = 0;
00077 
00078   /// get sigma, the background fluctuations per unit area; must be
00079   /// multipled by sqrt(area) to get fluctuations for a region of a
00080   /// given area.
00081   virtual double sigma() const { 
00082     throw Error("sigma() not supported for this Background Estimator");
00083   }
00084 
00085   /// get rho, the background density per unit area, locally at the
00086   /// position of a given jet. Note that this is not const, because a
00087   /// user may then wish to query other aspects of the background that
00088   /// could depend on the position of the jet last used for a rho(jet)
00089   /// determination.
00090   virtual double rho(const PseudoJet & jet) = 0;
00091 
00092   /// get sigma, the background fluctuations per unit area, locally at
00093   /// the position of a given jet. As for rho(jet), it is non-const.
00094   virtual double sigma(const PseudoJet & jet) { 
00095     throw Error("sigma(jet) not supported for this Background Estimator");
00096   }
00097 
00098   /// returns true if this background estimator has support for
00099   /// determination of sigma
00100   virtual bool has_sigma() {return false;}
00101   //\}
00102   
00103 
00104   /// @name configuring the behaviour
00105   //\{
00106   //----------------------------------------------------------------
00107 
00108   /// Set a pointer to a class that calculates the rescaling factor as
00109   /// a function of the jet (position). Note that the rescaling factor
00110   /// is used both in the determination of the "global" rho (the pt/A
00111   /// of each jet is divided by this factor) and when asking for a
00112   /// local rho (the result is multiplied by this factor).
00113   ///
00114   /// The BackgroundRescalingYPolynomial class can be used to get a
00115   /// rescaling that depends just on rapidity.
00116   virtual void set_rescaling_class(const FunctionOfPseudoJet<double> * rescaling_class) { _rescaling_class = rescaling_class; }
00117 
00118   /// return the pointer to the jet density class
00119   const FunctionOfPseudoJet<double> *  rescaling_class() const{
00120     return _rescaling_class;
00121   }
00122 
00123   //\}
00124 
00125   /// @name description
00126   //\{
00127   //----------------------------------------------------------------
00128 
00129   /// returns a textual description of the background estimator
00130   virtual std::string description() const = 0;
00131 
00132   //\}
00133 
00134 protected:
00135   /// @name helpers for derived classes
00136   ///
00137   /// Note that these helpers are related to median-based estimation
00138   /// of the background, so there is no guarantee that they will
00139   /// remain in this base class in the long term
00140   //\{
00141   //----------------------------------------------------------------
00142 
00143   /// given a quantity in a vector (e.g. pt_over_area) and knowledge
00144   /// about the number of empty jets, calculate the median and
00145   /// stand_dev_if_gaussian (roughly from the 16th percentile)
00146   ///
00147   /// If do_fj2_calculation is set to true then this performs FastJet
00148   /// 2.X estimation of the standard deviation, which has a spurious
00149   /// offset in the limit of a small number of jets.
00150   void _median_and_stddev(const std::vector<double> & quantity_vector, 
00151                           double n_empty_jets, 
00152                           double & median, 
00153                           double & stand_dev_if_gaussian,
00154                           bool do_fj2_calculation = false
00155                           ) const;
00156 
00157   /// computes a percentile of a given _sorted_ vector
00158   ///  \param sorted_quantity_vector   the vector contains the data sample
00159   ///  \param percentile               the percentile (defined between 0 and 1) to compute
00160   ///  \param nempty                   an additional number of 0's
00161   ///                                  (considered at the beginning of 
00162   ///                                  the quantity vector)
00163   ///  \param do_fj2_calculation       carry out the calculation as it
00164   ///                                  was done in fj2 (suffers from "edge effects")
00165   double _percentile(const std::vector<double> & sorted_quantity_vector, 
00166                      const double percentile, 
00167                      const double nempty=0.0,
00168                      const bool do_fj2_calculation = false) const;
00169 
00170   //\}
00171 
00172   const FunctionOfPseudoJet<double> * _rescaling_class;
00173   static LimitedWarning _warnings_empty_area;
00174 };
00175 
00176 
00177 
00178 //----------------------------------------------------------------------
00179 /// @ingroup tools_background
00180 /// A background rescaling that is a simple polynomial in y
00181 class BackgroundRescalingYPolynomial : public FunctionOfPseudoJet<double> {
00182 public:
00183   /// construct a background rescaling polynomial of the form
00184   /// a0 + a1*y + a2*y^2 + a3*y^3 + a4*y^4
00185   ///
00186   /// The following values give a reasonable reproduction of the
00187   /// Pythia8 tune 4C background shape for pp collisions at
00188   /// sqrt(s)=7TeV:
00189   ///
00190   /// - a0 =  1.157
00191   /// - a1 =  0
00192   /// - a2 = -0.0266
00193   /// - a3 =  0
00194   /// - a4 =  0.000048
00195   ///
00196   BackgroundRescalingYPolynomial(double a0=1, 
00197                                  double a1=0, 
00198                                  double a2=0, 
00199                                  double a3=0, 
00200                                  double a4=0) : _a0(a0), _a1(a1), _a2(a2), _a3(a3), _a4(a4) {}
00201 
00202   /// return the rescaling factor associated with this jet
00203   virtual double result(const PseudoJet & jet) const;
00204 private:
00205   double _a0, _a1, _a2, _a3, _a4;
00206 };
00207 
00208 
00209 
00210 
00211 
00212 FASTJET_END_NAMESPACE
00213 
00214 #endif  // __BACKGROUND_ESTIMATOR_BASE_HH__
00215 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends