fastjet::ClusterSequence1GhostPassiveArea Class Reference

Class that behaves essentially like ClusterSequence except that it also provides access to the area of a jet (which will be a random quantity. More...

#include <ClusterSequence1GhostPassiveArea.hh>

Inheritance diagram for fastjet::ClusterSequence1GhostPassiveArea:

Inheritance graph
fastjet::ClusterSequencePassiveAreafastjet::ClusterSequenceActiveAreafastjet::ClusterSequenceAreaBasefastjet::ClusterSequence
[legend]
Collaboration diagram for fastjet::ClusterSequence1GhostPassiveArea:

Collaboration graph
fastjet::ClusterSequenceActiveAreafastjet::ClusterSequenceAreaBasefastjet::ClusterSequencefastjet::JetDefinitionfastjet::JetDefinition::DefaultRecombinerfastjet::JetDefinition::Recombinerfastjet::JetDefinition::PluginLimitedWarning
[legend]
List of all members.

Public Member Functions

 ClusterSequence1GhostPassiveArea ()
template<class L>
 ClusterSequence1GhostPassiveArea (const std::vector< L > &pseudojets, const JetDefinition &jet_def, const GhostedAreaSpec &area_spec, const bool &writeout_combinations=false)
 constructor based on JetDefinition and 1GhostPassiveAreaSpec
virtual double n_empty_jets (const RangeDefinition &range) const
 return an estimate for the number of empty jets -- one uses the AreaBase one rather than the ActiveArea one (which for which we do not have the information).

Protected Member Functions

void _initialise_and_run_1GPA (const JetDefinition &jet_def, const GhostedAreaSpec &area_spec, const bool &writeout_combinations=false)
 global routine for initialising and running a general passive area

Private Member Functions

void _run_1GPA (const GhostedAreaSpec &area_spec)
 routine for running a passive area one ghost at a time.

Detailed Description

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...)

Definition at line 49 of file ClusterSequence1GhostPassiveArea.hh.


Constructor & Destructor Documentation

fastjet::ClusterSequence1GhostPassiveArea::ClusterSequence1GhostPassiveArea (  )  [inline]

Definition at line 52 of file ClusterSequence1GhostPassiveArea.hh.

00052 {}

template<class L>
fastjet::ClusterSequence1GhostPassiveArea::ClusterSequence1GhostPassiveArea ( const std::vector< L > &  pseudojets,
const JetDefinition jet_def,
const GhostedAreaSpec area_spec,
const bool &  writeout_combinations = false 
) [inline]

constructor based on JetDefinition and 1GhostPassiveAreaSpec

Definition at line 84 of file ClusterSequence1GhostPassiveArea.hh.

00087                                      {
00088 
00089   // transfer the initial jets (type L) into our own array
00090   _transfer_input_jets(pseudojets);
00091 
00092   // run the clustering for passive areas
00093   _initialise_and_run_1GPA(jet_def, area_spec, writeout_combinations);
00094 
00095 }


Member Function Documentation

virtual double fastjet::ClusterSequence1GhostPassiveArea::n_empty_jets ( const RangeDefinition range  )  const [inline, virtual]

return an estimate for the number of empty jets -- one uses the AreaBase one rather than the ActiveArea one (which for which we do not have the information).

Reimplemented from fastjet::ClusterSequenceActiveArea.

Definition at line 64 of file ClusterSequence1GhostPassiveArea.hh.

00064                                                                    {
00065     return ClusterSequenceAreaBase::n_empty_jets(range);
00066   }

void fastjet::ClusterSequence1GhostPassiveArea::_initialise_and_run_1GPA ( const JetDefinition jet_def,
const GhostedAreaSpec area_spec,
const bool &  writeout_combinations = false 
) [protected]

global routine for initialising and running a general passive area

Definition at line 40 of file ClusterSequence1GhostPassiveArea.cc.

References fastjet::ClusterSequenceActiveArea::_initialise_AA(), fastjet::ClusterSequenceActiveArea::_postprocess_AA(), _run_1GPA(), and fastjet::ClusterSequence::jet_def().

Referenced by fastjet::ClusterSequencePassiveArea::_initialise_and_run_PA().

00043                                                     {
00044 
00045   bool continue_running;
00046   _initialise_AA(jet_def,  area_spec, writeout_combinations, continue_running);
00047   if (continue_running) {
00048     _run_1GPA(area_spec);
00049     _postprocess_AA(area_spec);
00050   }
00051 }

void fastjet::ClusterSequence1GhostPassiveArea::_run_1GPA ( const GhostedAreaSpec area_spec  )  [private]

routine for running a passive area one ghost at a time.

Definition at line 56 of file ClusterSequence1GhostPassiveArea.cc.

References fastjet::ClusterSequenceActiveArea::_average_area, fastjet::ClusterSequenceActiveArea::_average_area2, fastjet::ClusterSequence::_jets, fastjet::ClusterSequenceActiveArea::_transfer_areas(), fastjet::ClusterSequenceActiveArea::_transfer_ghost_free_history(), fastjet::GhostedAreaSpec::actual_ghost_area(), fastjet::GhostedAreaSpec::add_ghosts(), fastjet::ClusterSequence::jet_def(), fastjet::GhostedAreaSpec::repeat(), and fastjet::ClusterSequence::unique_history_order().

Referenced by _initialise_and_run_1GPA().

00056                                                                                    {
00057     // record the input jets as they are currently
00058   vector<PseudoJet> input_jets(_jets);
00059 
00060   // code for testing the unique tree
00061   vector<int> unique_tree;
00062 
00063   // initialise our temporary average area^2
00064   valarray<double> lcl_average_area2(0.0, _average_area.size());
00065   valarray<double> last_average_area(0.0, _average_area.size());
00066 
00067   // run the clustering multiple times so as to get areas of all the jets
00068   for (int irepeat = 0; irepeat < area_spec.repeat(); irepeat++) {
00069 
00070     // first establish list of all ghosts
00071     vector<PseudoJet> all_ghosts;
00072     area_spec.add_ghosts(all_ghosts);
00073 
00074     // then run many subsets of ghosts (actually each subset contains just one ghost)
00075     for (unsigned ig = 0; ig < all_ghosts.size(); ig++) {
00076       vector<PseudoJet> some_ghosts;
00077       some_ghosts.push_back(all_ghosts[ig]);
00078       ClusterSequenceActiveAreaExplicitGhosts clust_seq(input_jets, jet_def(), 
00079                                                        some_ghosts, area_spec.actual_ghost_area());
00080 
00081       if (irepeat == 0 && ig == 0) {
00082         // take the non-ghost part of the history and put into our own
00083         // history.
00084         _transfer_ghost_free_history(clust_seq);
00085         // get the "unique" order that will be used for transferring all areas. 
00086         unique_tree = unique_history_order();
00087       }
00088       
00089       // transfer areas from clust_seq into our object
00090       _transfer_areas(unique_tree, clust_seq);
00091     }
00092     // keep track of fluctuations in area
00093     lcl_average_area2 += (_average_area - last_average_area)*
00094                          (_average_area - last_average_area);
00095     last_average_area = _average_area;
00096   }
00097   _average_area2 = lcl_average_area2;
00098 }


The documentation for this class was generated from the following files:
Generated on Tue Dec 18 17:05:50 2007 for fastjet by  doxygen 1.5.2