#include <PFJetFilter.h>

Inheritance diagram for PFJetFilter:

Public Member Functions
	PFJetFilter (const edm::ParameterSet &)
	~PFJetFilter ()
Private Member Functions
virtual void	beginJob ()
virtual void	endJob ()
virtual bool	filter (edm::Event &, const edm::EventSetup &)
double	resolution (double, bool)
double	response (double, bool)
Private Attributes
PFBenchmarkAlgo *	algo_
double	deltaEt_max
double	deltaEt_min
double	deltaR_max
double	deltaR_min
unsigned int	entry
double	eta_max
double	eta_min
double	genPt_cut
edm::InputTag	inputRecoLabel_
edm::InputTag	inputTruthLabel_
double	recPt_cut
bool	verbose

Detailed Description

Definition at line 12 of file PFJetFilter.h.

Constructor & Destructor Documentation

PFJetFilter::PFJetFilter ( const edm::ParameterSet & iConfig ) [explicit]

Definition at line 16 of file PFJetFilter.cc.

References edm::ParameterSet::getParameter().

{
  inputTruthLabel_             = iConfig.getParameter<edm::InputTag>("InputTruthLabel");
  inputRecoLabel_              = iConfig.getParameter<edm::InputTag>("InputRecoLabel");

  recPt_cut                    = iConfig.getParameter<double>("recPt");
  genPt_cut                    = iConfig.getParameter<double>("genPt");

  eta_min                      = iConfig.getParameter<double>("minEta");
  eta_max                      = iConfig.getParameter<double>("maxEta");

  deltaR_min                   = iConfig.getParameter<double>("deltaRMin");
  deltaR_max                   = iConfig.getParameter<double>("deltaRMax");

  deltaEt_min                  = iConfig.getParameter<double>("minDeltaEt");
  deltaEt_max                  = iConfig.getParameter<double>("maxDeltaEt");

  verbose                      = iConfig.getParameter<bool>("verbose");

  entry = 0;
}

PFJetFilter::~PFJetFilter ( )

Definition at line 38 of file PFJetFilter.cc.

{}

Member Function Documentation

void PFJetFilter::beginJob ( void ) [private, virtual]

Reimplemented from edm::EDFilter.

Definition at line 43 of file PFJetFilter.cc.

{}

void PFJetFilter::endJob ( void ) [private, virtual]

Reimplemented from edm::EDFilter.

Definition at line 46 of file PFJetFilter.cc.

{}

bool PFJetFilter::filter	(	edm::Event &	iEvent,
		const edm::EventSetup &	iESetup
	)		`[private, virtual]`

Implements edm::EDFilter.

Definition at line 49 of file PFJetFilter.cc.

References abs, sipixelsummary::Barrel, gather_cfg::cout, reco::deltaR(), reco::Candidate::eta(), edm::Event::getByLabel(), i, j, NULL, reco::Candidate::phi(), edm::Handle< T >::product(), reco::Candidate::pt(), ptmin, and dtT0WireCalibration_cfg::resolution.

{

// Typedefs to use views
  typedef edm::View<reco::Candidate> candidateCollection ;
  typedef edm::View<reco::Candidate> candidateCollection ;
  
  const candidateCollection *truth_candidates;
  const candidateCollection *reco_candidates;
 
  // ==========================================================
  // Retrieve!
  // ==========================================================

 // Get Truth Candidates (GenCandidates, GenJets, etc.)
    Handle<candidateCollection> truth_hnd;
    bool isGen = iEvent.getByLabel(inputTruthLabel_, truth_hnd);
    if ( !isGen ) { 
      std::cout << "Warning : no Gen jets in input !" << std::endl;
      return false;
    }

    truth_candidates = truth_hnd.product();

   // Get Reco Candidates (PFlow, CaloJet, etc.)
    Handle<candidateCollection> reco_hnd;
    bool isReco = iEvent.getByLabel(inputRecoLabel_, reco_hnd);
    if ( !isReco ) { 
      std::cout << "Warning : no Reco jets in input !" << std::endl;
      return false; 
    }
    reco_candidates = reco_hnd.product();
    if (!truth_candidates || !reco_candidates) return false;

    
    bool pass = false;
    
    for (unsigned int i = 0; i < reco_candidates->size(); i++) {
      
      const reco::Candidate *particle = &(*reco_candidates)[i];
      
      // This protection is meant at not being used !
      assert( particle!=NULL ); 

      double rec_pt = particle->pt();
      double rec_eta = particle->eta();
      double rec_phi = particle->phi();

      // skip PFjets with pt < recPt_cut GeV
      if (rec_pt < recPt_cut ) continue;

      // skip PFjets with eta > maxEta_cut or eta < minEta_cut
      if ( fabs(rec_eta) > eta_max ) continue;
      if ( fabs(rec_eta) < eta_min ) continue;

      bool Barrel = false;
      bool Endcap = false;
      bool Forward = false;
      if (std::abs(rec_eta) < 1.4 ) Barrel = true;
      if (std::abs (rec_eta) > 1.4 && std::abs (rec_eta) < 2.6 ) Endcap = true;
      if (std::abs (rec_eta) > 3.5 && std::abs (rec_eta) < 4.5 ) Forward = true;

      // Keep only jets in the barrel or the endcaps, within the tracker acceptance
      if ( !Barrel && !Endcap ) continue;

      // Find the closets recJet
      double deltaRmin = 999.;
      double ptmin = 0.;
      for (unsigned int j = 0; j < reco_candidates->size(); j++) {
      
        if ( i == j ) continue;
        const reco::Candidate *other = &(*reco_candidates)[j];
        double deltaR = algo_->deltaR(particle, other);
        if ( deltaR < deltaRmin && 
             other->pt() > 0.25*particle->pt() && 
             other->pt() > recPt_cut ) { 
          deltaRmin = deltaR;
          ptmin = other->pt();
        } 
        if ( deltaRmin < deltaR_min ) break;  
      }
      if ( deltaRmin < deltaR_min ) continue;  
      

      // Find the closest genJet.
      const reco::Candidate *gen_particle = 
        algo_->matchByDeltaR(particle,
                             truth_candidates);

      // Check there is a genJet associated to the recoJet
      if(gen_particle==NULL) continue; 

      // check deltaR is small enough
      double deltaR = algo_->deltaR(particle, gen_particle);
      if ( deltaR > deltaR_max ) continue;

      // double true_E = gen_particle->p();
      double true_pt = gen_particle->pt();
      double true_eta = gen_particle->eta();
      double true_phi = gen_particle->phi();

      // skip PFjets with pt < genPt_cut GeV
      if (true_pt < genPt_cut ) continue;
      
      // Find the pT residual
      double resPt = (rec_pt -true_pt)/true_pt;
      double sigma = resolution(true_pt,Barrel);
      double avera = response(true_pt,Barrel);
      double nSig  = (resPt-avera)/sigma;

      if ( nSig > deltaEt_max || nSig < deltaEt_min ) { 
        /* */
        if ( verbose ) 
          std::cout << "Entry " << entry 
                    << " resPt = " << resPt 
                    << " sigma/avera/nSig = " << sigma << "/" << avera <<"/" << nSig
                    << " pT (T/R) " << true_pt << "/" << rec_pt 
                    << " Eta (T/R) " << true_eta << "/" << rec_eta 
                    << " Phi (T/R) " << true_phi << "/" << rec_phi 
                    << " deltaRMin/ptmin " << deltaRmin << "/" << ptmin  
                    << std::endl;
        /* */
        pass = true;
      }
      
      if ( pass ) break;

    }

    entry++;
    return pass;

}

double PFJetFilter::resolution	(	double	pt,
		bool	barrel
	)		`[private]`

Definition at line 123 of file JetResolution.cc.

References newFWLiteAna::bin, i, MultiGaussianStateTransform::N, JetResolution::parameterFncs_, JetResolution::parameters_, JetResolution::resolutionFnc_, and findQualityFiles::size.

{
  unsigned N(y.size());
  for (unsigned iPar=0;iPar<parameters_.size();iPar++) {
    int bin = parameters_[iPar]->binIndex(x);
    assert(bin>=0);
    assert(bin<(int)parameters_[iPar]->size());
    const std::vector<float>& pars = parameters_[iPar]->record(bin).parameters();
    for (unsigned i=2*N;i<pars.size();i++)
      parameterFncs_[iPar]->SetParameter(i-2*N,pars[i]);
    float yy[4];
    for (unsigned i=0;i<N;i++)
      yy[i] = (y[i] < pars[2*i]) ? pars[2*i] : (y[i] > pars[2*i+1]) ? pars[2*i+1] : y[i];
    resolutionFnc_->SetParameter(iPar+1,
                                 parameterFncs_[iPar]->Eval(yy[0],yy[1],yy[2],yy[3]));
  }
  return resolutionFnc_;
}

double PFJetFilter::response	(	double	pt,
		bool	barrel
	)		`[private]`

Definition at line 196 of file PFJetFilter.cc.

References funct::exp(), p1, and p2.

                                            {

  double p0 = barrel ?  1.09906E-1 :  6.91939E+1;
  double p1 = barrel ? -1.61443E-1 : -6.92733E+1;
  double p2 = barrel ?  3.45489E+3 :  1.58207E+6;

  double resp = p0 + p1 * std::exp(-pt/p2);
  return resp;

}