FastJet 3.0.2
ClusterSequenceAreaBase.hh
00001 //STARTHEADER
00002 // $Id: ClusterSequenceAreaBase.hh 2687 2011-11-14 11:17:51Z soyez $
00003 //
00004 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
00005 //
00006 //----------------------------------------------------------------------
00007 // This file is part of FastJet.
00008 //
00009 //  FastJet is free software; you can redistribute it and/or modify
00010 //  it under the terms of the GNU General Public License as published by
00011 //  the Free Software Foundation; either version 2 of the License, or
00012 //  (at your option) any later version.
00013 //
00014 //  The algorithms that underlie FastJet have required considerable
00015 //  development and are described in hep-ph/0512210. If you use
00016 //  FastJet as part of work towards a scientific publication, please
00017 //  include a citation to the FastJet paper.
00018 //
00019 //  FastJet is distributed in the hope that it will be useful,
00020 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00021 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00022 //  GNU General Public License for more details.
00023 //
00024 //  You should have received a copy of the GNU General Public License
00025 //  along with FastJet. If not, see <http://www.gnu.org/licenses/>.
00026 //----------------------------------------------------------------------
00027 //ENDHEADER
00028 
00029 #ifndef __FASTJET_CLUSTERSEQUENCEAREABASE_HH__
00030 #define __FASTJET_CLUSTERSEQUENCEAREABASE_HH__
00031 
00032 #include "fastjet/ClusterSequence.hh"
00033 #include "fastjet/LimitedWarning.hh"
00034 #include "fastjet/Selector.hh"
00035 
00036 FASTJET_BEGIN_NAMESPACE
00037 
00038 /// @ingroup area_classes
00039 /// \class ClusterSequenceAreaBase
00040 /// base class that sets interface for extensions of ClusterSequence
00041 /// that provide information about the area of each jet
00042 ///
00043 /// the virtual functions here all return 0, since no area determination
00044 /// is implemented.
00045 class ClusterSequenceAreaBase : public ClusterSequence {
00046 public:
00047   
00048   /// a constructor which just carries out the construction of the
00049   /// parent class
00050   template<class L> ClusterSequenceAreaBase
00051          (const std::vector<L> & pseudojets, 
00052           const JetDefinition & jet_def_in,
00053           const bool & writeout_combinations = false) :
00054            ClusterSequence(pseudojets, jet_def_in, writeout_combinations) {}
00055 
00056 
00057   /// default constructor
00058   ClusterSequenceAreaBase() {}
00059 
00060 
00061   /// destructor
00062   virtual ~ClusterSequenceAreaBase() {}
00063 
00064 
00065   /// return the area associated with the given jet; this base class
00066   /// returns 0.
00067   virtual double area       (const PseudoJet & ) const {return 0.0;}
00068 
00069   /// return the error (uncertainty) associated with the determination
00070   /// of the area of this jet; this base class returns 0.
00071   virtual double area_error (const PseudoJet & ) const {return 0.0;}
00072 
00073   /// return a PseudoJet whose 4-vector is defined by the following integral
00074   ///
00075   ///       \int drap d\phi PseudoJet("rap,phi,pt=one") *
00076   ///                           * Theta("rap,phi inside jet boundary")
00077   ///
00078   /// where PseudoJet("rap,phi,pt=one") is a 4-vector with the given
00079   /// rapidity (rap), azimuth (phi) and pt=1, while Theta("rap,phi
00080   /// inside jet boundary") is a function that is 1 when rap,phi
00081   /// define a direction inside the jet boundary and 0 otherwise.
00082   ///
00083   /// This base class returns a null 4-vector.
00084   virtual PseudoJet area_4vector(const PseudoJet & ) const {
00085     return PseudoJet(0.0,0.0,0.0,0.0);}
00086 
00087   /// true if a jet is made exclusively of ghosts
00088   ///
00089   /// NB: most area classes do not give any explicit ghost jets, but
00090   /// some do, and they should replace this function with their own
00091   /// version.
00092   virtual bool is_pure_ghost(const PseudoJet & ) const {
00093     return false;
00094   }
00095 
00096   /// returns true if ghosts are explicitly included within 
00097   /// jets for this ClusterSequence; 
00098   ///
00099   /// Derived classes that do include explicit ghosts should provide
00100   /// an alternative version of this routine and set it properly.
00101   virtual bool has_explicit_ghosts() const {
00102     return false;
00103   }
00104 
00105   /// return the total area, corresponding to the given Selector, that
00106   /// is free of jets, in general based on the inclusive jets.
00107   /// 
00108   /// The selector passed as an argument has to have a finite area and
00109   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00110   /// tools for more generic usages)
00111   virtual double empty_area(const Selector & selector) const;
00112 
00113   /// return the total area, corresponding to the given Selector, that
00114   /// is free of jets, based on the supplied all_jets
00115   /// 
00116   /// The selector passed as an argument has to have a finite area and
00117   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00118   /// tools for more generic usages)
00119   double empty_area_from_jets(const std::vector<PseudoJet> & all_jets,
00120                               const Selector & selector) const;
00121 
00122   /// return something similar to the number of pure ghost jets
00123   /// in the given selector's range in an active area case.
00124   /// For the local implementation we return empty_area/(0.55 pi R^2),
00125   /// based on measured properties of ghost jets with kt and cam
00126   /// (cf arXiv:0802.1188).
00127   ///
00128   /// Note that the number returned is a double.
00129   /// 
00130   /// The selector passed as an argument has to have a finite area and
00131   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00132   /// tools for more generic usages)
00133   virtual double n_empty_jets(const Selector & selector) const {
00134     double R = jet_def().R();
00135     return empty_area(selector)/(0.55*pi*R*R);
00136   }
00137 
00138   /// the median of (pt/area) for jets contained within the selector
00139   /// range, making use also of the info on n_empty_jets
00140   /// 
00141   /// The selector passed as an argument has to have a finite area and
00142   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00143   /// tools for more generic usages)
00144   double median_pt_per_unit_area(const Selector & selector) const;
00145 
00146   /// the median of (pt/area_4vector) for jets contained within the
00147   /// selector range, making use also of the info on n_empty_jets
00148   /// 
00149   /// The selector passed as an argument has to have a finite area and
00150   /// apply jet-by-jet
00151   double median_pt_per_unit_area_4vector(const Selector & selector) const;
00152   
00153   /// the function that does the work for median_pt_per_unit_area and 
00154   /// median_pt_per_unit_area_4vector: 
00155   /// - something_is_area_4vect = false -> use plain area
00156   /// - something_is_area_4vect = true  -> use 4-vector area
00157   double median_pt_per_unit_something(
00158                     const Selector & selector, bool use_area_4vector) const;
00159 
00160   /// using jets withing the selector range (and with 4-vector areas if
00161   /// use_area_4vector), calculate the median pt/area, as well as an
00162   /// "error" (uncertainty), which is defined as the 1-sigma
00163   /// half-width of the distribution of pt/A, obtained by looking for
00164   /// the point below which we have (1-0.6827)/2 of the jets
00165   /// (including empty jets).
00166   ///
00167   /// The subtraction for a jet with uncorrected pt pt^U and area A is
00168   ///
00169   ///   pt^S = pt^U - median*A +- sigma*sqrt(A)
00170   ///
00171   /// where the error is only that associated with the fluctuations
00172   /// in the noise and not that associated with the noise having 
00173   /// caused changes in the hard-particle content of the jet.
00174   /// 
00175   /// The selector passed as an argument has to have a finite area and
00176   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00177   /// tools for more generic usages)
00178   ///
00179   /// NB: subtraction may also be done with 4-vector area of course,
00180   /// and this is recommended for jets with larger values of R, as
00181   /// long as rho has also been determined with a 4-vector area;
00182   /// using a scalar area causes one to neglect terms of relative
00183   /// order $R^2/8$ in the jet $p_t$.
00184   virtual void get_median_rho_and_sigma(const Selector & selector, 
00185                                         bool use_area_4vector,
00186                                         double & median, double & sigma,
00187                                         double & mean_area) const;
00188 
00189   /// a more advanced version of get_median_rho_and_sigma, which allows
00190   /// one to use any "view" of the event containing all jets (so that, 
00191   /// e.g. one might use Cam on a different resolution scale without
00192   /// have to rerun the algorithm).
00193   ///
00194   /// By default it will assume that "all" are not inclusive jets, 
00195   /// so that in dealing with empty area it has to calculate
00196   /// the number of empty jets based on the empty area and the
00197   /// the observed <area> of jets rather than a surmised area
00198   ///
00199   /// Note that for small effective radii, this can cause problems
00200   /// because the harder jets get an area >> <ghost-jet-area>
00201   /// and so the estimate comes out all wrong. In these situations
00202   /// it is highly advisable to use an area with explicit ghosts, since
00203   /// then the "empty" jets are actually visible.
00204   /// 
00205   /// The selector passed as an argument has to have a finite area and
00206   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00207   /// tools for more generic usages)
00208   virtual void get_median_rho_and_sigma(const std::vector<PseudoJet> & all_jets,
00209                                         const Selector & selector, 
00210                                         bool use_area_4vector,
00211                                         double & median, double & sigma,
00212                                         double & mean_area,
00213                                         bool all_are_inclusive = false) const;
00214 
00215   /// same as the full version of get_median_rho_and_error, but without
00216   /// access to the mean_area
00217   /// 
00218   /// The selector passed as an argument has to have a finite area and
00219   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00220   /// tools for more generic usages)
00221   virtual void get_median_rho_and_sigma(const Selector & selector, 
00222                                 bool use_area_4vector,
00223                                 double & median, double & sigma) const {
00224     double mean_area;
00225     get_median_rho_and_sigma(selector,  use_area_4vector,
00226                              median,  sigma, mean_area);
00227   }
00228   
00229 
00230   /// fits a form pt_per_unit_area(y) = a + b*y^2 in the selector range. 
00231   /// exclude_above allows one to exclude large values of pt/area from fit. 
00232   ///               (if negative, the cut is discarded)
00233   /// use_area_4vector = true uses the 4vector areas.
00234   /// 
00235   /// The selector passed as an argument has to have a finite area and
00236   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00237   /// tools for more generic usages)
00238   virtual void parabolic_pt_per_unit_area(double & a, double & b, 
00239                                           const Selector & selector, 
00240                                           double exclude_above=-1.0, 
00241                                           bool use_area_4vector=false) const;
00242 
00243   /// return a vector of all subtracted jets, using area_4vector, given rho.
00244   /// Only inclusive_jets above ptmin are subtracted and returned.
00245   /// the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()),
00246   /// i.e. not necessarily ordered in pt once subtracted
00247   std::vector<PseudoJet> subtracted_jets(const double rho,
00248                                          const double ptmin=0.0) const;
00249 
00250   /// return a vector of subtracted jets, using area_4vector.
00251   /// Only inclusive_jets above ptmin are subtracted and returned.
00252   /// the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()),
00253   /// i.e. not necessarily ordered in pt once subtracted
00254   /// 
00255   /// The selector passed as an argument has to have a finite area and
00256   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00257   /// tools for more generic usages)
00258   std::vector<PseudoJet> subtracted_jets(const Selector & selector, 
00259                                          const double ptmin=0.0) const;
00260 
00261   /// return a subtracted jet, using area_4vector, given rho
00262   PseudoJet subtracted_jet(const PseudoJet & jet,
00263                            const double rho) const;
00264 
00265   /// return a subtracted jet, using area_4vector; note
00266   /// that this is potentially inefficient if repeatedly used for many
00267   /// different jets, because rho will be recalculated each time
00268   /// around.
00269   /// 
00270   /// The selector passed as an argument has to have a finite area and
00271   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00272   /// tools for more generic usages)
00273   PseudoJet subtracted_jet(const PseudoJet & jet,
00274                            const Selector & selector) const;
00275 
00276   /// return the subtracted pt, given rho
00277   double subtracted_pt(const PseudoJet & jet,
00278                        const double rho,
00279                        bool use_area_4vector=false) const;
00280 
00281   /// return the subtracted pt; note that this is
00282   /// potentially inefficient if repeatedly used for many different
00283   /// jets, because rho will be recalculated each time around.
00284   /// 
00285   /// The selector passed as an argument has to have a finite area and
00286   /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
00287   /// tools for more generic usages)
00288   double subtracted_pt(const PseudoJet & jet,
00289                        const Selector & selector,
00290                        bool use_area_4vector=false) const;
00291 
00292 protected:
00293   /// check the selector is suited for the computations i.e. applies jet by jet and has a finite area
00294   void _check_selector_good_for_median(const Selector &selector) const;
00295 
00296 
00297 private:
00298   /// handle warning messages
00299   static LimitedWarning _warnings;
00300   static LimitedWarning _warnings_zero_area;
00301   static LimitedWarning _warnings_empty_area;
00302 
00303   /// check the jet algorithm is suitable (and if not issue a warning)
00304   void _check_jet_alg_good_for_median() const;
00305   
00306 };
00307 
00308 
00309 
00310 FASTJET_END_NAMESPACE
00311 
00312 #endif // __FASTJET_CLUSTERSEQUENCEAREABASE_HH__
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends