Fork me on GitHub

Changeset de6d698 in git for modules/FastJetFinder.cc


Ignore:
Timestamp:
Mar 16, 2015, 4:39:51 PM (10 years ago)
Author:
Michele Selvaggi <michele.selvaggi@…>
Branches:
ImprovedOutputFile, Timing, dual_readout, llp, master
Children:
986d9d5
Parents:
666d795
Message:

added Trimming, Pruning and SoftDrop in FastJetFinder

File:
1 edited

Legend:

Unmodified
Added
Removed
  • modules/FastJetFinder.cc

    r666d795 rde6d698  
    6666#include "fastjet/contribs/Nsubjettiness/WinnerTakeAllRecombiner.hh"
    6767
     68#include "fastjet/tools/Filter.hh"
     69#include "fastjet/tools/Pruner.hh"
     70#include "fastjet/contribs/RecursiveTools/SoftDrop.hh"
     71
    6872using namespace std;
    6973using namespace fastjet;
     
    121125  fN = GetInt("N", 2);                  // used only if Njettiness is used as jet clustering algo (case 8)
    122126
     127  //-- Trimming parameters --
     128 
     129  fComputeTrimming = GetBool("ComputeTrimming", false);
     130  fRTrim = GetDouble("RTrim", 0.2);
     131  fPtFracTrim = GetDouble("PtFracTrim", 0.05);
     132 
     133
     134  //-- Pruning parameters --
     135 
     136  fComputePruning = GetBool("ComputePruning", false);
     137  fZcutPrun = GetDouble("ZcutPrun", 0.1);
     138  fRcutPrun = GetDouble("RcutPrun", 0.5);
     139  fRPrun = GetDouble("RPrun", 0.8);
     140 
     141  //-- SoftDrop parameters --
     142 
     143  fComputeSoftDrop     = GetBool("ComputeSoftDrop", false);
     144  fBetaSoftDrop        = GetDouble("BetaSoftDrop", 0.0);
     145  fSymmetryCutSoftDrop = GetDouble("SymmetryCutSoftDrop", 0.1);
     146  fR0SoftDrop= GetDouble("R0SoftDrop=", 0.8);
     147 
     148
    123149  // ---  Jet Area Parameters ---
    124150  fAreaAlgorithm = GetInt("AreaAlgorithm", 0);
     
    260286  PseudoJet jet, area;
    261287  ClusterSequence *sequence;
    262   vector< PseudoJet > inputList, outputList;
     288  vector< PseudoJet > inputList, outputList, subjets;
    263289  vector< PseudoJet >::iterator itInputList, itOutputList;
    264290  vector< TEstimatorStruct >::iterator itEstimators;
     
    352378    candidate->DeltaEta = detaMax;
    353379    candidate->DeltaPhi = dphiMax;
    354 
     380       
     381    //------------------------------------
     382    // Trimming
     383    //------------------------------------
     384
     385     if(fComputeTrimming)
     386    {
     387
     388      fastjet::Filter    trimmer(fastjet::JetDefinition(fastjet::kt_algorithm,fRTrim),fastjet::SelectorPtFractionMin(fPtFracTrim));
     389      fastjet::PseudoJet trimmed_jet = trimmer(*itOutputList);
     390     
     391      trimmed_jet = join(trimmed_jet.constituents());
     392     
     393      candidate->TrimmedP4[0].SetPtEtaPhiM(trimmed_jet.pt(), trimmed_jet.eta(), trimmed_jet.phi(), trimmed_jet.m());
     394       
     395      // four hardest subjets
     396      subjets.clear();
     397      subjets = trimmed_jet.pieces();
     398      subjets = sorted_by_pt(subjets);
     399     
     400      candidate->NSubJetsTrimmed = subjets.size();
     401
     402      for (size_t i = 0; i < subjets.size() and i < 4; i++){
     403        if(subjets.at(i).pt() < 0) continue ;
     404        candidate->TrimmedP4[i+1].SetPtEtaPhiM(subjets.at(i).pt(), subjets.at(i).eta(), subjets.at(i).phi(), subjets.at(i).m());
     405      }
     406    }
     407   
     408   
     409    //------------------------------------
     410    // Pruning
     411    //------------------------------------
     412   
     413   
     414    if(fComputePruning)
     415    {
     416
     417      fastjet::Pruner    pruner(fastjet::JetDefinition(fastjet::cambridge_algorithm,fRPrun),fZcutPrun,fRcutPrun);
     418      fastjet::PseudoJet pruned_jet = pruner(*itOutputList);
     419
     420      candidate->PrunedP4[0].SetPtEtaPhiM(pruned_jet.pt(), pruned_jet.eta(), pruned_jet.phi(), pruned_jet.m());
     421         
     422      // four hardest subjet
     423      subjets.clear();
     424      subjets = pruned_jet.pieces();
     425      subjets = sorted_by_pt(subjets);
     426     
     427      candidate->NSubJetsPruned = subjets.size();
     428
     429      for (size_t i = 0; i < subjets.size() and i < 4; i++){
     430        if(subjets.at(i).pt() < 0) continue ;
     431        candidate->PrunedP4[i+1].SetPtEtaPhiM(subjets.at(i).pt(), subjets.at(i).eta(), subjets.at(i).phi(), subjets.at(i).m());
     432      }
     433
     434    }
     435     
     436    //------------------------------------
     437    // SoftDrop
     438    //------------------------------------
     439   
     440    if(fComputeSoftDrop)
     441    {
     442   
     443      contrib::SoftDrop  softDrop(fBetaSoftDrop,fSymmetryCutSoftDrop,fR0SoftDrop);
     444      fastjet::PseudoJet softdrop_jet = softDrop(*itOutputList);
     445     
     446      candidate->SoftDroppedP4[0].SetPtEtaPhiM(softdrop_jet.pt(), softdrop_jet.eta(), softdrop_jet.phi(), softdrop_jet.m());
     447       
     448      // four hardest subjet
     449     
     450      subjets.clear();
     451      subjets    = softdrop_jet.pieces();
     452      subjets    = sorted_by_pt(subjets);
     453      candidate->NSubJetsSoftDropped = softdrop_jet.pieces().size();
     454
     455      for (size_t i = 0; i < subjets.size()  and i < 4; i++){
     456        if(subjets.at(i).pt() < 0) continue ;
     457        candidate->SoftDroppedP4[i+1].SetPtEtaPhiM(subjets.at(i).pt(), subjets.at(i).eta(), subjets.at(i).phi(), subjets.at(i).m());
     458      }
     459    }
     460 
    355461    // --- compute N-subjettiness with N = 1,2,3,4,5 ----
    356462
Note: See TracChangeset for help on using the changeset viewer.