FastJet 3.0.2
ClusterSequenceVoronoiArea.hh
00001 //STARTHEADER
00002 // $Id: ClusterSequenceVoronoiArea.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_CLUSTERSEQUENCEVORONOIAREA_HH__
00030 #define __FASTJET_CLUSTERSEQUENCEVORONOIAREA_HH__
00031 
00032 #include "fastjet/PseudoJet.hh"
00033 #include "fastjet/AreaDefinition.hh"
00034 #include "fastjet/ClusterSequenceAreaBase.hh"
00035 #include <memory>
00036 #include <vector>
00037 
00038 FASTJET_BEGIN_NAMESPACE      // defined in fastjet/internal/base.hh
00039 
00040 /// @ingroup sec_area_classes
00041 /// \class ClusterSequenceVoronoiArea
00042 /// Like ClusterSequence with computation of the Voronoi jet area
00043 ///
00044 /// Handle the computation of Voronoi jet area.
00045 ///
00046 /// This class should not be used directly. Rather use
00047 /// ClusterSequenceArea with the appropriate AreaDefinition
00048 class ClusterSequenceVoronoiArea : public ClusterSequenceAreaBase {
00049 public:
00050   /// template ctor
00051   /// \param pseudojet              list of jets (template type)
00052   /// \param jet_def                jet definition
00053   /// \param effective_Rfact        effective radius
00054   /// \param writeout_combinations  ??????
00055   template<class L> ClusterSequenceVoronoiArea
00056          (const std::vector<L> & pseudojets, 
00057           const JetDefinition & jet_def,
00058           const VoronoiAreaSpec & spec = VoronoiAreaSpec(),
00059           const bool & writeout_combinations = false);
00060   
00061   /// default dtor
00062   ~ClusterSequenceVoronoiArea();
00063 
00064   /// return the area associated with the given jet
00065   virtual inline double area(const PseudoJet & jet) const {
00066     return _voronoi_area[jet.cluster_hist_index()];}
00067 
00068   /// return a 4-vector area associated with the given jet -- stricly
00069   /// this is not the exact 4-vector area, but rather an approximation
00070   /// made of sums of centres of all Voronoi cells in jet, each
00071   /// contributing with a normalisation equal to the area of the cell
00072   virtual inline PseudoJet area_4vector(const PseudoJet & jet) const {
00073     return _voronoi_area_4vector[jet.cluster_hist_index()];}
00074 
00075   /// return the error of the area associated with the given jet
00076   /// (0 by definition for a voronoi area)
00077   virtual inline double area_error(const PseudoJet & /*jet*/) const {
00078     return 0.0;}
00079 
00080   /// passive area calculator -- to be defined in the .cc file (it will do
00081   /// the true hard work)
00082   class VoronoiAreaCalc; 
00083   
00084 
00085 private:  
00086   /// initialisation of the Voronoi Area
00087   void _initializeVA();
00088 
00089   std::vector<double> _voronoi_area;  ///< vector containing the result
00090   std::vector<PseudoJet> _voronoi_area_4vector; ///< vector containing approx 4-vect areas
00091   VoronoiAreaCalc *_pa_calc;          ///< area calculator
00092   double _effective_Rfact;            ///< effective radius
00093 };
00094 
00095 
00096 
00097 
00098 /// template constructor need to be specified in the header!
00099 //----------------------------------------------------------------------
00100 template<class L> ClusterSequenceVoronoiArea::ClusterSequenceVoronoiArea
00101 (const std::vector<L> &pseudojets, 
00102  const JetDefinition &jet_def_in,
00103  const VoronoiAreaSpec & spec,
00104  const bool & writeout_combinations) :
00105   _effective_Rfact(spec.effective_Rfact()) {
00106 
00107   // transfer the initial jets (type L) into our own array
00108   _transfer_input_jets(pseudojets);
00109 
00110   // run the clustering
00111   _initialise_and_run(jet_def_in,writeout_combinations);
00112 
00113   // the jet clustering's already been done, now worry about areas...
00114   _initializeVA();
00115 }
00116 
00117 FASTJET_END_NAMESPACE
00118 
00119 #endif // __FASTJET_CLUSTERSEQUENCEVORONOIAREA_HH__
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends