gen::JetMatchingMGFastJet Class Reference

#include <JetMatchingMGFastJet.h>

Inheritance diagram for gen::JetMatchingMGFastJet:

Public Member Functions
const std::vector< int > *	getPartonList ()
	JetMatchingMGFastJet (const edm::ParameterSet &params)
	~JetMatchingMGFastJet ()
Public Member Functions inherited from gen::JetMatchingMadgraph
	JetMatchingMadgraph (const edm::ParameterSet &params)
	~JetMatchingMadgraph ()
Public Member Functions inherited from gen::JetMatching
bool	isMatchingDone ()
	JetMatching (const edm::ParameterSet &params)
void	resetMatchingStatus ()
virtual	~JetMatching ()

Protected Member Functions
void	beforeHadronisation (const lhef::LHEEvent *)
void	beforeHadronisationExec ()
virtual void	init (const lhef::LHERunInfo *runInfo)
bool	initAfterBeams ()
int	match (const lhef::LHEEvent *partonLevel, const std::vector< fastjet::PseudoJet > *jetInput)
Protected Member Functions inherited from gen::JetMatchingMadgraph
std::set< std::string >	capabilities () const
virtual double	getJetEtaMax () const
template<typename T >
T	getParameter (const std::string &var, const T &defValue=T()) const

Private Types
enum	partonTypes { ID_TOP =6, ID_GLUON =21, ID_PHOTON =22 }
enum	vetoStatus {   NONE, LESS_JETS, MORE_JETS, HARD_JET,   UNMATCHED_PARTON }

Private Attributes
double	clFact
double	coneRadius
bool	doMerge
double	etaJetMax
std::vector< fastjet::PseudoJet >	fClusJets
int	fDJROutFlag
std::ofstream	fDJROutput
bool	fExcLocal
bool	fIsInit
fastjet::JetDefinition *	fJetFinder
std::vector< fastjet::PseudoJet >	fPtSortedJets
int	jetAlgoPower
int	nJetMax
int	nJetMin
int	nQmatch
double	qCut
double	qCutSq
std::vector< int >	typeIdx [3]

Additional Inherited Members
Static Public Member Functions inherited from gen::JetMatching
static std::auto_ptr< JetMatching >	create (const edm::ParameterSet &params)
Static Protected Member Functions inherited from gen::JetMatchingMadgraph
template<typename T >
static T	getParameter (const std::map< std::string, std::string > &params, const std::string &var, const T &defValue=T())
template<typename T >
static T	parseParameter (const std::string &value)
template<typename T >
static void	updateOrDie (const std::map< std::string, std::string > &params, T &param, const std::string &name)
Protected Attributes inherited from gen::JetMatchingMadgraph
bool	eventInitialized
bool	exclusive
std::map< std::string, std::string >	mgParams
bool	runInitialized
bool	soup
Protected Attributes inherited from gen::JetMatching
bool	fMatchingStatus

Detailed Description

Definition at line 36 of file JetMatchingMGFastJet.h.

Member Enumeration Documentation

enum gen::JetMatchingMGFastJet::partonTypes

private

Enumerator
ID_TOP
ID_GLUON
ID_PHOTON

Definition at line 61 of file JetMatchingMGFastJet.h.

{ ID_TOP=6, ID_GLUON=21, ID_PHOTON=22 };

enum gen::JetMatchingMGFastJet::vetoStatus

private

Enumerator
NONE
LESS_JETS
MORE_JETS
HARD_JET
UNMATCHED_PARTON

Definition at line 60 of file JetMatchingMGFastJet.h.

{ NONE, LESS_JETS, MORE_JETS, HARD_JET, UNMATCHED_PARTON };

Constructor & Destructor Documentation

gen::JetMatchingMGFastJet::JetMatchingMGFastJet ( const edm::ParameterSet &  params )

inline

Definition at line 41 of file JetMatchingMGFastJet.h.

References fDJROutFlag, and edm::ParameterSet::getParameter().

                                                          : JetMatchingMadgraph(params), 
                                                              fJetFinder(0),  
                                  fIsInit(false) 
                                  { fDJROutFlag = params.getParameter<int>("outTree_flag"); }

gen::JetMatchingMGFastJet::~JetMatchingMGFastJet ( )

inline

Definition at line 45 of file JetMatchingMGFastJet.h.

References fJetFinder.

{ if (fJetFinder) delete fJetFinder; }

Member Function Documentation

void gen::JetMatchingMGFastJet::beforeHadronisation ( const lhef::LHEEvent *  lhee )

protectedvirtual

Reimplemented from gen::JetMatchingMadgraph.

Definition at line 76 of file JetMatchingMGFastJet.cc.

References edm::hlt::Exception, gen::JetMatchingMadgraph::exclusive, fExcLocal, lhef::LHEEvent::getHEPEUP(), i, ID_GLUON, ID_TOP, lhef::HEPEUP::IDUP, customizeTrackingMonitorSeedNumber::idx, lhef::HEPEUP::ISTUP, lhef::HEPEUP::MOTHUP, nJetMax, nQmatch, lhef::HEPEUP::NUP, gen::JetMatchingMadgraph::runInitialized, gen::JetMatchingMadgraph::soup, and typeIdx.

{

   if (!runInitialized)
        throw cms::Exception("Generator|PartonShowerVeto")
            << "Run not initialized in JetMatchingMadgraph"
            << std::endl;
   
   for (int i = 0; i < 3; i++) 
   {
    typeIdx[i].clear();
   }

   // Sort original process final state into light/heavy jets and 'other'.
   // Criteria:
   //   1 <= ID <= 5 and massless, or ID == 21 --> light jet (typeIdx[0])
   //   4 <= ID <= 6 and massive               --> heavy jet (typeIdx[1])
   //   All else                               --> other     (typeIdx[2])
   //
   const lhef::HEPEUP& hepeup = *lhee->getHEPEUP(); 
   int idx = 2; 
   for ( int i=0; i < hepeup.NUP; i++ )
   {
      if ( hepeup.ISTUP[i] < 0 ) continue;
      if ( hepeup.MOTHUP[i].first != 1 && hepeup.MOTHUP[i].second !=2 ) continue; // this way we skip entries that come
                                                                                  // from resonance decays;
                                          // we only take those that descent
                                          // directly from "incoming partons"
      idx = 2;
      if ( hepeup.IDUP[i] == ID_GLUON || (fabs(hepeup.IDUP[i]) <= nQmatch) ) // light jet
         // light jet
     idx = 0;
      else if ( fabs(hepeup.IDUP[i]) > nQmatch && fabs(hepeup.IDUP[i]) <= ID_TOP) // heavy jet
         idx = 1;
      // Store
      typeIdx[idx].push_back(i); 
   } 
   
   // NOTE: In principle, I should use exclusive, inclusive, or soup !!!   
   // should be like this:
   if ( soup )
   {
      int NPartons = typeIdx[0].size(); 
      fExcLocal = ( NPartons < nJetMax );
   }
   else
      fExcLocal = exclusive;
        
   return;

}

void gen::JetMatchingMGFastJet::beforeHadronisationExec ( )

inlineprotectedvirtual

Reimplemented from gen::JetMatchingMadgraph.

Definition at line 54 of file JetMatchingMGFastJet.h.

{ return; }

const std::vector<int>* gen::JetMatchingMGFastJet::getPartonList ( )

inlinevirtual

Reimplemented from gen::JetMatching.

Definition at line 47 of file JetMatchingMGFastJet.h.

References typeIdx.

{ return typeIdx; }

virtual void gen::JetMatchingMGFastJet::init ( const lhef::LHERunInfo *  runInfo )

inlineprotectedvirtual

Reimplemented from gen::JetMatchingMadgraph.

Definition at line 50 of file JetMatchingMGFastJet.h.

References fIsInit, gen::JetMatchingMadgraph::init(), and initAfterBeams().

                                                         { if (fIsInit) return; JetMatchingMadgraph::init(runInfo); 
                                                             initAfterBeams(); fIsInit=true; return; }

bool gen::JetMatchingMGFastJet::initAfterBeams ( )

protectedvirtual

Reimplemented from gen::JetMatching.

Definition at line 32 of file JetMatchingMGFastJet.cc.

References clFact, coneRadius, doMerge, gen::MEMAIN::etaclmax, etaJetMax, fClusJets, fExcLocal, fIsInit, fJetFinder, fPtSortedJets, gen::UPPRIV::ickkw, jetAlgoPower, gen::MEMAIN::maxjets, gen::memain_, gen::MEMAIN::minjets, nJetMax, nJetMin, gen::MEMAIN::nqmatch, nQmatch, funct::pow(), gen::MEMAIN::qcut, qCut, qCutSq, and gen::uppriv_.

Referenced by init().

{
  
  if ( fIsInit ) return true;
  
  //
  doMerge = uppriv_.ickkw; 
  // doMerge = true;
  qCut    = memain_.qcut; // 
  nQmatch        = memain_.nqmatch;
  clFact         = 1.; // Default value 
                      // NOTE: ME2pythia seems to default to 1.5 - need to check !!!
              // in general, needs to read key ALPSFACT from LHE file - fix CMSSW code !!!
  nJetMin        = memain_.minjets;
  nJetMax        = memain_.maxjets;

  etaJetMax      = memain_.etaclmax; 

  coneRadius     = 1.0; 
  jetAlgoPower   = 1; //   this is the kT algorithm !!! 

  // Matching procedure
  //
  qCutSq         = pow(qCut,2);  
  // this should be something like memaev_.iexc
  fExcLocal = true;


   // If not merging, then done (?????)
   //
   // if (!doMerge) return true;
  
   // Initialise chosen jet algorithm. 
   //
   fJetFinder = new fastjet::JetDefinition( fastjet::kt_algorithm, coneRadius );
   fClusJets.clear();
   fPtSortedJets.clear();  
   
   fIsInit = true;
   
   return true;

}

int gen::JetMatchingMGFastJet::match ( const lhef::LHEEvent *  partonLevel, const std::vector< fastjet::PseudoJet > *  jetInput  )

protectedvirtual

Reimplemented from gen::JetMatchingMadgraph.

Definition at line 128 of file JetMatchingMGFastJet.cc.

References clFact, funct::false, fClusJets, fDJROutFlag, fDJROutput, fExcLocal, fJetFinder, fPtSortedJets, lhef::LHEEvent::getHEPEUP(), i, customizeTrackingMonitorSeedNumber::idx, LESS_JETS, max(), MORE_JETS, NONE, convertSQLitetoXML_cfg::output, funct::pow(), lhef::HEPEUP::PUP, qCut, qCutSq, mathSSE::sqrt(), typeIdx, and UNMATCHED_PARTON.

{
         
   // Number of hard partons
   //
   int NPartons = typeIdx[0].size();
   
   fClusJets.clear();
   fPtSortedJets.clear();
   
   int ClusSeqNJets = 0;

   fastjet::ClusterSequence ClusSequence( *jetInput, *fJetFinder );
   
   if ( fExcLocal )
   {
      fClusJets = ClusSequence.exclusive_jets( qCutSq );
   }
   else
   {
      fClusJets = ClusSequence.inclusive_jets( qCut ); 
   }
   
   ClusSeqNJets = fClusJets.size(); 
   
   
   if ( ClusSeqNJets < NPartons ) return LESS_JETS;
   
   double localQcutSq = qCutSq;
   
   if ( fExcLocal ) // exclusive
   {
      if( ClusSeqNJets > NPartons ) return MORE_JETS;
   }
   else // inclusive 
   {
      fPtSortedJets = fastjet::sorted_by_pt( *jetInput );
      localQcutSq = std::max( qCutSq, fPtSortedJets[0].pt2() );
      fClusJets = ClusSequence.exclusive_jets( NPartons ); // override
      ClusSeqNJets = NPartons;
   }
   
  if( clFact != 0 ) localQcutSq *= pow(clFact,2);
  
  std::vector<fastjet::PseudoJet> MatchingInput;
  
  std::vector<bool> jetAssigned;
  jetAssigned.assign( fClusJets.size(), false );

  int counter = 0;
  
  const lhef::HEPEUP& hepeup = *partonLevel->getHEPEUP();
    
  while ( counter < NPartons )
  {

     MatchingInput.clear();

     for ( int i=0; i<ClusSeqNJets; i++ )
     {
        if ( jetAssigned[i] ) continue;
        if ( i == NPartons )  break;
    //
    // this looks "awkward" but this way we do NOT pass in cluster_hist_index
    // which would actually indicate position in "history" of the "master" ClusSeq
    //
    fastjet::PseudoJet exjet = fClusJets[i];
        MatchingInput.push_back( fastjet::PseudoJet( exjet.px(), exjet.py(), exjet.pz(), exjet.e() ) );
    MatchingInput.back().set_user_index(i);
     }

     int idx = typeIdx[0][counter];
     MatchingInput.push_back( fastjet::PseudoJet( hepeup.PUP[idx][0],
                                                  hepeup.PUP[idx][1],
                          hepeup.PUP[idx][2],
                          hepeup.PUP[idx][3]) );

     //
     // in principle, one can use ClusterSequence::n_particles() 
     //
     int NBefore = MatchingInput.size();

     // Create new clustering object - which includes the 1st clustering run !!!
     // NOTE-1: it better be a new object for each try, or history will be "too crowded".
     // NOTE-2: when created, the it ALWAYS makes at least 1 clustering step, so in most
     //         cases at least 1 new jet object will be added; although in some cases 
     //         no jet is added, if the system doesn't cluster to anything (for example,
     //         input jet(s) and parton(s) are going opposite directions)
     //
     fastjet::ClusterSequence ClusSeq( MatchingInput, *fJetFinder );
     
     const std::vector<fastjet::PseudoJet>&        output = ClusSeq.jets();
     int NClusJets = output.size() - NBefore;
     
     // 
     // JVY - I think this is the right idea:
     // at least 1 (one) new jet needs to be added
     // however, need to double check details and refine, especially for inclusive mode
     //
     // 
     if ( NClusJets < 1 )
     {
        return UNMATCHED_PARTON;
     }
     //
     // very unlikely case but let's do it just to be safe
     //
     if ( NClusJets >= NBefore )
     {
        return MORE_JETS;
     }
      
     // Now browse history and see how close the clustering distance
     //
     // NOTE: Remember, there maybe more than one new jet in the list (for example,
     //       for process=2,3,...);
     //       in this case we take the ** first ** one that satisfies the distance/cut,
     //       which is ** typically ** also the best one 
     //
     bool matched = false;
     const std::vector<fastjet::ClusterSequence::history_element>&  history = ClusSeq.history();
              
     // for ( unsigned int i=nBefore; i<history.size(); i++ )
     for ( unsigned int i=NBefore; i<output.size(); i++ )
     {

        int hidx = output[i].cluster_sequence_history_index();
    double dNext = history[hidx].dij;
    if ( dNext < localQcutSq )
    {
       //
       // the way we form input, parent1 is always jet, and parent2 can be any,
       // but if it's a good match/cluster, it should point at the parton
       //
       int parent1 = history[hidx].parent1;
       int parent2 = history[hidx].parent2;
       if ( parent1 < 0 || parent2 < 0 ) break; // bad cluster, no match
       //
       // pull up jet's "global" index 
       //
       int pidx = MatchingInput[parent1].user_index();
       jetAssigned[pidx] = true; 
       matched = true;
       break;
    }
     }
     if ( !matched )
     {
        return UNMATCHED_PARTON;
     }
     
     counter++;   
  }

   // Now save some info for DJR analysis (if requested).
   // This mimics what is done in ME2pythia.f
   // Basically, NJets and matching scale for these 4 cases:
   // 1->0, 2->1, 3->2, and 4->3
   //
   if ( fDJROutFlag > 0 )
   {   
      std::vector<double> MergingScale;
      MergingScale.clear();
      for ( int nj=0; nj<4; nj++ )
      {
         double dmscale2 = ClusSequence.exclusive_dmerge( nj );
         double dmscale = sqrt( dmscale2 );
         MergingScale.push_back( dmscale );
      }
      fDJROutput.open( "events.tree", std::ios_base::app );  
      double dNJets = (double)NPartons;
      fDJROutput << " " << dNJets << " " << MergingScale[0] << " " 
                                         << MergingScale[1] << " "
                                         << MergingScale[2] << " "
                         << MergingScale[3] << std::endl;  
      fDJROutput.close();
   }
 
  return NONE;

}

Member Data Documentation

double gen::JetMatchingMGFastJet::clFact

private

Definition at line 64 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), and match().

double gen::JetMatchingMGFastJet::coneRadius

private

Definition at line 81 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams().

bool gen::JetMatchingMGFastJet::doMerge

private

Definition at line 74 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams().

double gen::JetMatchingMGFastJet::etaJetMax

private

Definition at line 81 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams().

std::vector<fastjet::PseudoJet> gen::JetMatchingMGFastJet::fClusJets

private

Definition at line 95 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), and match().

int gen::JetMatchingMGFastJet::fDJROutFlag

private

Definition at line 100 of file JetMatchingMGFastJet.h.

Referenced by JetMatchingMGFastJet(), and match().

std::ofstream gen::JetMatchingMGFastJet::fDJROutput

private

Definition at line 99 of file JetMatchingMGFastJet.h.

Referenced by match().

bool gen::JetMatchingMGFastJet::fExcLocal

private

Definition at line 86 of file JetMatchingMGFastJet.h.

Referenced by beforeHadronisation(), initAfterBeams(), and match().

bool gen::JetMatchingMGFastJet::fIsInit

private

Definition at line 102 of file JetMatchingMGFastJet.h.

Referenced by init(), and initAfterBeams().

fastjet::JetDefinition* gen::JetMatchingMGFastJet::fJetFinder

private

Definition at line 94 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), match(), and ~JetMatchingMGFastJet().

std::vector<fastjet::PseudoJet> gen::JetMatchingMGFastJet::fPtSortedJets

private

Definition at line 95 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), and match().

int gen::JetMatchingMGFastJet::jetAlgoPower

private

Definition at line 80 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams().

int gen::JetMatchingMGFastJet::nJetMax

private

Definition at line 77 of file JetMatchingMGFastJet.h.

Referenced by beforeHadronisation(), and initAfterBeams().

int gen::JetMatchingMGFastJet::nJetMin

private

Definition at line 77 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams().

int gen::JetMatchingMGFastJet::nQmatch

private

Definition at line 65 of file JetMatchingMGFastJet.h.

Referenced by beforeHadronisation(), and initAfterBeams().

double gen::JetMatchingMGFastJet::qCut

private

Definition at line 63 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), and match().

double gen::JetMatchingMGFastJet::qCutSq

private

Definition at line 63 of file JetMatchingMGFastJet.h.

Referenced by initAfterBeams(), and match().

std::vector< int > gen::JetMatchingMGFastJet::typeIdx[3]

private

Definition at line 89 of file JetMatchingMGFastJet.h.

Referenced by beforeHadronisation(), getPartonList(), and match().