//STARTHEADER // $Id: ClusterSequenceActiveArea.hh,v 1.1 2008-11-06 14:32:07 ovyn Exp $ // // Copyright (c) 2005-2006, Matteo Cacciari and Gavin Salam // //---------------------------------------------------------------------- // This file is part of FastJet. // // FastJet is free software; you can redistribute it and/or modify // it under the terms of the GNU General Public License as published by // the Free Software Foundation; either version 2 of the License, or // (at your option) any later version. // // The algorithms that underlie FastJet have required considerable // development and are described in hep-ph/0512210. If you use // FastJet as part of work towards a scientific publication, please // include a citation to the FastJet paper. // // FastJet is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the // GNU General Public License for more details. // // You should have received a copy of the GNU General Public License // along with FastJet; if not, write to the Free Software // Foundation, Inc.: // 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA //---------------------------------------------------------------------- //ENDHEADER #ifndef __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__ #define __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__ #include "Utilities/Fastjet/include/fastjet/PseudoJet.hh" #include "Utilities/Fastjet/include/fastjet/ClusterSequenceAreaBase.hh" #include "Utilities/Fastjet/include/fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh" #include #include //------------ backwards compatibility with version 2.1 ------------- // for backwards compatibility make ActiveAreaSpec name available #include "Utilities/Fastjet/include/fastjet/ActiveAreaSpec.hh" #include "Utilities/Fastjet/include/fastjet/ClusterSequenceWithArea.hh" //-------------------------------------------------------------------- FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh //using namespace std; /// Class that behaves essentially like ClusterSequence except /// that it also provides access to the area of a jet (which /// will be a random quantity... Figure out what to do about seeds /// later...) class ClusterSequenceActiveArea : public ClusterSequenceAreaBase { public: /// default constructor ClusterSequenceActiveArea() {} /// constructor based on JetDefinition and GhostedAreaSpec template ClusterSequenceActiveArea (const std::vector & pseudojets, const JetDefinition & jet_def, const GhostedAreaSpec & area_spec, const bool & writeout_combinations = false) ; virtual double area (const PseudoJet & jet) const { return _average_area[jet.cluster_hist_index()];}; virtual double area_error (const PseudoJet & jet) const { return _average_area2[jet.cluster_hist_index()];}; virtual PseudoJet area_4vector (const PseudoJet & jet) const { return _average_area_4vector[jet.cluster_hist_index()];}; /// enum providing a variety of tentative strategies for estimating /// the background (e.g. non-jet) activity in a highly populated event; the /// one that has been most extensively tested is median. enum mean_pt_strategies{median=0, non_ghost_median, pttot_over_areatot, pttot_over_areatot_cut, mean_ratio_cut, play, median_4vector}; /// return the transverse momentum per unit area according to one /// of the above strategies; for some strategies (those with "cut" /// in their name) the parameter "range" allows one to exclude a /// subset of the jets for the background estimation, those that /// have pt/area > median(pt/area)*range. /// /// NB: This call is OBSOLETE; use media_pt_per_unit_area from the // ClusterSequenceAreaBase class instead double pt_per_unit_area(mean_pt_strategies strat=median, double range=2.0 ) const; // following code removed -- now dealt with by AreaBase class (and // this definition here conflicts with it). // /// fits a form pt_per_unit_area(y) = a + b*y^2 in the range // /// abs(y) & unique_hist_order, const ClusterSequenceActiveAreaExplicitGhosts & ); /// child classes benefit from having these at their disposal std::valarray _average_area, _average_area2; std::valarray _average_area_4vector; /// returns true if there are any particles whose transverse momentum /// if so low that there's a risk of the ghosts having modified the /// clustering sequence bool has_dangerous_particles() const {return _has_dangerous_particles;} private: double _non_jet_area, _non_jet_area2, _non_jet_number; double _maxrap_for_area; // max rap where we put ghosts double _safe_rap_for_area; // max rap where we trust jet areas bool _has_dangerous_particles; /// routine for extracting the tree in an order that will be independent /// of any degeneracies in the recombination sequence that don't /// affect the composition of the final jets void _extract_tree(std::vector &) const; /// do the part of the extraction associated with pos, working /// through its children and their parents void _extract_tree_children(int pos, std::valarray &, const std::valarray &, std::vector &) const; /// do the part of the extraction associated with the parents of pos. void _extract_tree_parents (int pos, std::valarray &, const std::valarray &, std::vector &) const; /// check if two jets have the same momentum to within the /// tolerance (and if pt's are not the same we're forgiving and /// look to see if the energy is the same) void _throw_unless_jets_have_same_perp_or_E(const PseudoJet & jet, const PseudoJet & refjet, double tolerance, const ClusterSequenceActiveAreaExplicitGhosts & jets_ghosted_seq ) const; /// since we are playing nasty games with seeds, we should warn /// the user a few times //static int _n_seed_warnings; //const static int _max_seed_warnings = 10; // record the number of repeats int _area_spec_repeat; /// a class for our internal storage of ghost jets class GhostJet : public PseudoJet { public: GhostJet(const PseudoJet & j, double a) : PseudoJet(j), area(a){} double area; }; std::vector _ghost_jets; std::vector _unclustered_ghosts; }; template ClusterSequenceActiveArea::ClusterSequenceActiveArea (const std::vector & pseudojets, const JetDefinition & jet_def, const GhostedAreaSpec & area_spec, const bool & writeout_combinations) { // transfer the initial jets (type L) into our own array _transfer_input_jets(pseudojets); // run the clustering for active areas _initialise_and_run_AA(jet_def, area_spec, writeout_combinations); } FASTJET_END_NAMESPACE #endif // __FASTJET_CLUSTERSEQUENCEACTIVEAREA_HH__