Loading...
Searching...
No Matches
14-groomers.cc
Go to the documentation of this file.
//----------------------------------------------------------------------
/// \file
/// \page Example14 14 - unified use of transformers
///
/// fastjet example program to illustrate the use of the
/// fastjet::Filter and fastjet::Pruner classes in a unified way
/// through their derivation from fastjet::Transformer
///
/// The two hardest jets in a boosted top event, clustered with an
/// (abnormally) large R, are then groomed using different tools. One
/// notes the reduction in the mass of the jets after grooming.
///
/// run it with : ./14-groomers < data/boosted_top_event.dat
///
/// Source code: 14-groomers.cc
//----------------------------------------------------------------------
//STARTHEADER
// $Id$
//
// Copyright (c) 2005-2018, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
//
//----------------------------------------------------------------------
// 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, see <http://www.gnu.org/licenses/>.
//----------------------------------------------------------------------
//ENDHEADER
#include <fastjet/PseudoJet.hh>
#include <fastjet/ClusterSequence.hh>
#include <fastjet/Selector.hh>
#include <iostream>
#include "fastjet/tools/Filter.hh"
#include "fastjet/tools/Pruner.hh"
#include <cstdio> // needed for io
using namespace fastjet;
using namespace std;
/// an example program showing how to use Filter and Pruner in FastJet
int main(){
// read in input particles
//----------------------------------------------------------
vector<PseudoJet> input_particles;
double px, py , pz, E;
while (cin >> px >> py >> pz >> E) {
// create a PseudoJet with these components and put it onto
// back of the input_particles vector
input_particles.push_back(PseudoJet(px,py,pz,E));
}
// get the resulting jets ordered in pt
//----------------------------------------------------------
JetDefinition jet_def(cambridge_algorithm, 1.5);
ClusterSequence clust_seq(input_particles, jet_def);
vector<PseudoJet> inclusive_jets =
sorted_by_pt(clust_seq.inclusive_jets(5.0));
// label the columns
printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", "phi", "pt", "mass");
// print out the details for each jet
for (unsigned int i = 0; i < inclusive_jets.size(); i++) {
printf("%5u %15.8f %15.8f %15.8f %15.8f\n",
i, inclusive_jets[i].rap(), inclusive_jets[i].phi(),
inclusive_jets[i].perp(),inclusive_jets[i].m());
}
// simple test to avoid that the example below crashes:
// make sure there is at least 3 jets above our 5 GeV
if (inclusive_jets.size()<2){
cout << "Please provide an event with at least 2 jets above 5 GeV" << endl;
return 1;
}
// We will groom the two hardest jets of the event
//----------------------------------------------------------
vector<PseudoJet> candidates = SelectorNHardest(2)(inclusive_jets);
// create 3 groomers
//----------------------------------------------------------
vector<Transformer *> groomers;
// 1.
// the Cambridge/Aachen filter with Rfilt=0.3
// (simplified version of arXiv:0802.2470)
double Rfilt = 0.3;
unsigned int nfilt = 3;
groomers.push_back(new Filter(JetDefinition(cambridge_algorithm, Rfilt),
SelectorNHardest(nfilt) ) );
// 2.
// Filtering with a pt cut as for trimming (arXiv:0912.1342)
double Rtrim = 0.2;
double ptfrac = 0.03;
groomers.push_back(new Filter(JetDefinition(kt_algorithm, Rtrim),
SelectorPtFractionMin(ptfrac) ) );
// 3.
// Pruning (arXiv:0903.5081)
double zcut = 0.1;
double rcut_factor = 0.5;
groomers.push_back(new Pruner(cambridge_algorithm, zcut, rcut_factor));
// apply the various groomers to the test PseudoJet's
// and show the result
//----------------------------------------------------------
// print out original jet candidates
cout << "\nOriginal jets that will be grooomed: " << endl;
for (vector<PseudoJet>::iterator jit=candidates.begin(); jit!=candidates.end(); jit++){
const PseudoJet & c = *jit;
cout << " rap = " << c.rap() << ", phi = " << c.phi() << ", pt = " << c.perp()
<< ", mass = " << c.m()
<< " [" << c.description() << "]" << endl;
}
// loop on groomers
for (unsigned int i=0; i < groomers.size(); i++){
const Transformer & f = *groomers[i];
cout << "\nUsing groomer: " << f.description() << endl;
// loop on jet candidates
for (vector<PseudoJet>::iterator jit=candidates.begin(); jit!=candidates.end(); jit++){
const PseudoJet & c = *jit;
// apply groomer f to jet c
PseudoJet j = f(c);
// access properties specific to the given transformer
//
// We first make sure that the jet indeed has a structure
// compatible with the result of a Filter or Pruner (using
// has_structure_of()), and then retrieve the pieces rejected by the
// groomer (using structure_of())
int n_rejected;
if (j.has_structure_of<Filter>()) {
const Filter::StructureType & fj_struct = j.structure_of<Filter>();
n_rejected = fj_struct.rejected().size();
}else {
assert(j.has_structure_of<Pruner>()); // make sure
const Pruner::StructureType & fj_struct = j.structure_of<Pruner>();
n_rejected = fj_struct.rejected().size();
}
// write out result
cout << " rap = " << j.rap() << ", phi = " << j.phi() << ", pt = " << j.perp()
<< " mass = " << j.m()
<< " [kept: " << j.pieces().size()
<< ", rejected: " << n_rejected;
if (j.has_structure_of<Pruner>()) {
cout << ", Rcut: " << j.structure_of<Pruner>().Rcut();
}
cout << "]" << endl;
}
}
// a bit of memory cleaning
for (unsigned int i=0; i < groomers.size(); i++) delete groomers[i];
return 0;
}
main
int main()
an example program showing how to use Filter and Pruner in FastJet
Definition 14-groomers.cc:59
fastjet::ClusterSequence
deals with clustering
Definition ClusterSequence.hh:63
fastjet::ClusterSequence::inclusive_jets
std::vector< PseudoJet > inclusive_jets(const double ptmin=0.0) const
return a vector of all jets (in the sense of the inclusive algorithm) with pt >= ptmin.
Definition ClusterSequence.cc:916
fastjet::FilterStructure::rejected
const std::vector< PseudoJet > & rejected() const
returns the subjets that were not kept during the filtering procedure (subtracted if the filter reque...
Definition Filter.hh:228
fastjet::Filter
Class that helps perform filtering (Butterworth, Davison, Rubin and Salam, arXiv:0802....
Definition Filter.hh:97
fastjet::JetDefinition
class that is intended to hold a full definition of the jet clusterer
Definition JetDefinition.hh:250
fastjet::PrunerStructure::rejected
std::vector< PseudoJet > rejected() const
return the constituents that have been rejected
Definition Pruner.hh:200
fastjet::Pruner
Transformer that prunes a jet.
Definition Pruner.hh:107
fastjet::PseudoJet
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
Definition PseudoJet.hh:68
fastjet::PseudoJet::rap
double rap() const
returns the rapidity or some large value when the rapidity is infinite
Definition PseudoJet.hh:138
fastjet::PseudoJet::structure_of
const TransformerType::StructureType & structure_of() const
this is a helper to access any structure created by a Transformer (that is, of type Transformer::Stru...
Definition PseudoJet.hh:1122
fastjet::PseudoJet::phi
double phi() const
returns phi (in the range 0..2pi)
Definition PseudoJet.hh:123
fastjet::PseudoJet::perp
double perp() const
returns the scalar transverse momentum
Definition PseudoJet.hh:158
fastjet::PseudoJet::has_structure_of
bool has_structure_of() const
check if the PseudoJet has the structure resulting from a Transformer (that is, its structure is comp...
Definition PseudoJet.hh:1112
fastjet::PseudoJet::m
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP)
Definition PseudoJet.hh:1060
fastjet::PseudoJet::description
std::string description() const
return a string describing what kind of PseudoJet we are dealing with
Definition PseudoJet.cc:517
fastjet::PseudoJet::pieces
virtual std::vector< PseudoJet > pieces() const
retrieve the pieces that make up the jet.
Definition PseudoJet.cc:792
fastjet::Transformer
Base (abstract) class for a jet transformer.
Definition Transformer.hh:71
fastjet::Transformer::description
virtual std::string description() const =0
This should be overloaded to return a description of the Transformer.
fastjet::SelectorNHardest
Selector SelectorNHardest(unsigned int n)
select the n hardest objects
Definition Selector.cc:1074
fastjet::SelectorPtFractionMin
Selector SelectorPtFractionMin(double fraction)
select objects that carry at least a fraction "fraction" of the reference jet.
Definition Selector.cc:1365
fastjet
the FastJet namespace
Definition AreaDefinition.cc:37
fastjet::cambridge_algorithm
@ cambridge_algorithm
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
Definition JetDefinition.hh:142
fastjet::kt_algorithm
@ kt_algorithm
the longitudinally invariant kt algorithm
Definition JetDefinition.hh:139
fastjet::sorted_by_pt
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
Definition PseudoJet.cc:882