CMS 3D CMS Logo

HLTHemiDPhiFilter.cc
Go to the documentation of this file.
3 
5 
9 
12 
13 
16 
20 
21 //
22 // constructors and destructor
23 //
25  inputTag_ (iConfig.getParameter<edm::InputTag>("inputTag")),
26  min_dphi_ (iConfig.getParameter<double> ("minDPhi" )),
27  accept_NJ_ (iConfig.getParameter<bool> ("acceptNJ" ))
28 
29 {
30  m_theHemiToken = consumes<std::vector<math::XYZTLorentzVector>>(inputTag_);
31  LogDebug("") << "Inputs/minDphi/acceptNJ : "
32  << inputTag_.encode() << " "
33  << min_dphi_ << " "
34  << accept_NJ_ << ".";
35 }
36 
38 
39 void
43  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltRHemisphere"));
44  desc.add<double>("minDPhi",2.9415);
45  desc.add<bool>("acceptNJ",true);
46  descriptions.add("hltHemiDPhiFilter",desc);
47 }
48 
49 //
50 // member functions
51 //
52 
53 // ------------ method called to produce the data ------------
54 bool
56 {
57  using namespace std;
58  using namespace edm;
59  using namespace reco;
60  using namespace trigger;
61 
62  // get hold of collection of objects
64  iEvent.getByToken (m_theHemiToken,hemispheres);
65 
66  // check the the input collections are available
67  if (not hemispheres.isValid())
68  return false;
69 
70  if(hemispheres->empty()){ // the Hemisphere Maker will produce an empty collection of hemispheres if the number of jets in the
71  return accept_NJ_; // event is greater than the maximum number of jets
72  }
73 
74 
75 
76 
77  //***********************************
78  // Get the set of hemisphere axes
79 
80  TLorentzVector j1R(hemispheres->at(0).x(),hemispheres->at(0).y(),hemispheres->at(0).z(),hemispheres->at(0).t());
81  TLorentzVector j2R(hemispheres->at(1).x(),hemispheres->at(1).y(),hemispheres->at(1).z(),hemispheres->at(1).t());
82 
83  // compute the dPhi between them
84  double dphi = 50.;
85  dphi = deltaPhi(j1R.Phi(),j2R.Phi());
86 
87  // Dphi requirement
88 
89  if(dphi<=min_dphi_) return true;
90 
91  // filter decision
92  return false;
93 }
94 
95 
96 double HLTHemiDPhiFilter::deltaPhi(double v1, double v2)
97 {
98  // Computes the correctly normalized phi difference
99  // v1, v2 = phi of object 1 and 2
100  double diff = std::abs(v2 - v1);
101  return (diff < M_PI) ? diff : 2*M_PI - diff;
102 }
103 
#define LogDebug(id)
edm::InputTag inputTag_
HLTHemiDPhiFilter(const edm::ParameterSet &)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
std::string encode() const
Definition: InputTag.cc:159
int iEvent
Definition: GenABIO.cc:224
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::EDGetTokenT< std::vector< math::XYZTLorentzVector > > m_theHemiToken
static double deltaPhi(double, double)
ParameterDescriptionBase * add(U const &iLabel, T const &value)
bool isValid() const
Definition: HandleBase.h:74
~HLTHemiDPhiFilter() override
#define M_PI
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
void add(std::string const &label, ParameterSetDescription const &psetDescription)
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
fixed size matrix
HLT enums.