CMS 3D CMS Logo

HLTRFilter.cc
Go to the documentation of this file.
2 
4 
8 
11 
13 
16 
20 
23 
24 //
25 // constructors and destructor
26 //
28  inputTag_ (iConfig.getParameter<edm::InputTag>("inputTag")),
29  inputMetTag_ (iConfig.getParameter<edm::InputTag>("inputMetTag")),
30  doMuonCorrection_(iConfig.getParameter<bool> ("doMuonCorrection" )),
31  min_R_ (iConfig.getParameter<double> ("minR" )),
32  min_MR_ (iConfig.getParameter<double> ("minMR" )),
33  DoRPrime_ (iConfig.getParameter<bool> ("doRPrime" )),
34  accept_NJ_ (iConfig.getParameter<bool> ("acceptNJ" )),
35  R_offset_ ((iConfig.existsAs<double>("R2Offset") ? iConfig.getParameter<double>("R2Offset"):0)),
36  MR_offset_ ((iConfig.existsAs<double>("MROffset") ? iConfig.getParameter<double>("MROffset"):0)),
37  R_MR_cut_ ((iConfig.existsAs<double>("RMRCut") ? iConfig.getParameter<double>("RMRCut"):-999999.))
38 
39 {
40  m_theInputToken = consumes<std::vector<math::XYZTLorentzVector>>(inputTag_);
41  m_theMETToken = consumes<edm::View<reco::MET> >(inputMetTag_);
42  LogDebug("") << "Inputs/minR/minMR/doRPrime/acceptNJ/R2Offset/MROffset/RMRCut : "
43  << inputTag_.encode() << " "
44  << inputMetTag_.encode() << " "
45  << min_R_ << " "
46  << min_MR_ << " "
47  << accept_NJ_
48  << R_offset_ << " "
49  << MR_offset_ << " "
50  << R_MR_cut_
51  << ".";
52 
53  //put a dummy METCollection into the event, holding values for MR and Rsq
54  produces<reco::METCollection>();
55 
56 }
57 
58 HLTRFilter::~HLTRFilter() = default;
59 
60 void
64  desc.add<edm::InputTag>("inputTag",edm::InputTag("hltRHemisphere"));
65  desc.add<edm::InputTag>("inputMetTag",edm::InputTag("hltMet"));
66  desc.add<bool>("doMuonCorrection",false);
67  desc.add<double>("minR",0.3);
68  desc.add<double>("minMR",100.0);
69  desc.add<bool>("doRPrime",false);
70  desc.add<bool>("acceptNJ",true);
71  desc.add<double>("R2Offset",0.0);
72  desc.add<double>("MROffset",0.0);
73  desc.add<double>("RMRCut",-999999.0);
74  descriptions.add("hltRFilter",desc);
75 }
76 
77 //
78 // member functions
79 //
80 
81 // ------------ method called to produce the data ------------
82 bool
84 {
85  using namespace std;
86  using namespace edm;
87  using namespace reco;
88 
89  // get hold of collection of objects
91  iEvent.getByToken (m_theInputToken,hemispheres);
92 
93  // get hold of the MET Collection
95  iEvent.getByToken(m_theMETToken,inputMet);
96 
97  // check the the input collections are available
98  if (not hemispheres.isValid() or not inputMet.isValid())
99  return false;
100 
101  if(hemispheres->empty()){ // the Hemisphere Maker will produce an empty collection of hemispheres if the number of jets in the
102  return accept_NJ_; // event is greater than the maximum number of jets
103  }
104 
105  //***********************************
106  //Calculate R
107 
108  int nMuons;
109  switch(hemispheres->size()){
110  case 2:
111  nMuons=0; break;
112  case 5:
113  nMuons=1; break;
114  case 10:
115  nMuons=2; break;
116  default:
117  return false; //invalid hemisphere collection
118  }
119 
120  //muons as MET
121  TLorentzVector ja(hemispheres->at(0).x(),hemispheres->at(0).y(),hemispheres->at(0).z(),hemispheres->at(0).t());
122  TLorentzVector jb(hemispheres->at(1).x(),hemispheres->at(1).y(),hemispheres->at(1).z(),hemispheres->at(1).t());
123 
124  std::vector<math::XYZTLorentzVector> muonVec;
125 
126  double MR = CalcMR(ja,jb);
127  double R = CalcR(MR,ja,jb,inputMet,muonVec);
128 
129  if(MR>=min_MR_ && R>=min_R_
130  && ( (R*R - R_offset_)*(MR-MR_offset_) )>=R_MR_cut_) {
131  addObjects(iEvent, filterproduct, MR, R*R);
132  return true;
133  }
134  if(nMuons==0) {
135  addObjects(iEvent, filterproduct, MR, R*R);
136  return false; // if no muons and we get here, reject event
137  }
138  //Lead Muon as Jet
139  ja.SetXYZT(hemispheres->at(3).x(),hemispheres->at(3).y(),hemispheres->at(3).z(),hemispheres->at(3).t());
140  jb.SetXYZT(hemispheres->at(4).x(),hemispheres->at(4).y(),hemispheres->at(4).z(),hemispheres->at(4).t());
141  muonVec.push_back(hemispheres->at(2)); // lead muon at position 2
142 
143  MR = CalcMR(ja,jb);
144  R = CalcR(MR,ja,jb,inputMet,muonVec);
145  if(MR>=min_MR_ && R>=min_R_
146  && ( (R*R - R_offset_)*(MR-MR_offset_) )>=R_MR_cut_){
147  addObjects(iEvent, filterproduct, MR, R*R);
148  return true;
149  }
150 
151  if(nMuons==1){
152  addObjects(iEvent, filterproduct, MR, R*R);
153  return false; // if no muons and we get here, reject event
154  }
155 
156  muonVec.pop_back();
157  //Second Muon as Jet
158  ja.SetXYZT(hemispheres->at(6).x(),hemispheres->at(6).y(),hemispheres->at(6).z(),hemispheres->at(6).t());
159  jb.SetXYZT(hemispheres->at(7).x(),hemispheres->at(7).y(),hemispheres->at(7).z(),hemispheres->at(7).t());
160  muonVec.push_back(hemispheres->at(5)); // sublead muon at position 5
161 
162  MR = CalcMR(ja,jb);
163  R = CalcR(MR,ja,jb,inputMet,muonVec);
164  if(MR>=min_MR_ && R>=min_R_
165  && ( (R*R - R_offset_)*(MR-MR_offset_) )>=R_MR_cut_){
166  addObjects(iEvent, filterproduct, MR, R*R);
167  return true;
168  }
169 
170  ja.SetXYZT(hemispheres->at(8).x(),hemispheres->at(8).y(),hemispheres->at(8).z(),hemispheres->at(8).t());
171  jb.SetXYZT(hemispheres->at(9).x(),hemispheres->at(9).y(),hemispheres->at(9).z(),hemispheres->at(9).t());
172  muonVec.push_back(hemispheres->at(2)); // lead muon at position 2
173 
174  MR = CalcMR(ja,jb);
175  R = CalcR(MR,ja,jb,inputMet,muonVec);
176  if(MR>=min_MR_ && R>=min_R_
177  && ( (R*R - R_offset_)*(MR-MR_offset_) )>=R_MR_cut_){
178  addObjects(iEvent, filterproduct, MR, R*R);
179  return true;
180  }
181 
182  // filter decision
183  addObjects(iEvent, filterproduct, MR, R*R);
184  return false;
185 }
186 
187 void HLTRFilter::addObjects(edm::Event& iEvent, trigger::TriggerFilterObjectWithRefs & filterproduct, double MR, double Rsq) const{
188 
189  //create METCollection for storing MR and Rsq results
190  std::unique_ptr<reco::METCollection> razorObject(new reco::METCollection());
191 
192  reco::MET::LorentzVector mrRsqP4(MR,Rsq,0,0);
193  reco::MET::Point vtx(0,0,0);
194 
195  reco::MET mrRsq(mrRsqP4, vtx);
196  razorObject->push_back(mrRsq);
197 
199 
200  //put the METCollection into the event (necessary because of how addCollectionTag works...)
201  iEvent.put(std::move(razorObject));
202  edm::Ref<reco::METCollection> mrRsqRef(ref_before_put, 0);
203 
204  //add it to the trigger object collection
205  if (saveTags()) filterproduct.addCollectionTag(edm::InputTag( *moduleLabel()));
206  filterproduct.addObject(trigger::TriggerMET, mrRsqRef); //give it the ID of a MET object
207 }
208 
209 double
210 HLTRFilter::CalcMR(TLorentzVector ja, TLorentzVector jb){
211  if(ja.Pt()<=0.1) return -1;
212 
213  ja.SetPtEtaPhiM(ja.Pt(),ja.Eta(),ja.Phi(),0.0);
214  jb.SetPtEtaPhiM(jb.Pt(),jb.Eta(),jb.Phi(),0.0);
215 
216  if(ja.Pt() > jb.Pt()){
217  TLorentzVector temp = ja;
218  ja = jb;
219  jb = temp;
220  }
221 
222  double A = ja.P();
223  double B = jb.P();
224  double az = ja.Pz();
225  double bz = jb.Pz();
226  TVector3 jaT, jbT;
227  jaT.SetXYZ(ja.Px(),ja.Py(),0.0);
228  jbT.SetXYZ(jb.Px(),jb.Py(),0.0);
229  double ATBT = (jaT+jbT).Mag2();
230 
231  double MR = sqrt((A+B)*(A+B)-(az+bz)*(az+bz)-
232  (jbT.Dot(jbT)-jaT.Dot(jaT))*(jbT.Dot(jbT)-jaT.Dot(jaT))/(jaT+jbT).Mag2());
233 
234  double mybeta = (jbT.Dot(jbT)-jaT.Dot(jaT))/
235  sqrt(ATBT*((A+B)*(A+B)-(az+bz)*(az+bz)));
236 
237  double mygamma = 1./sqrt(1.-mybeta*mybeta);
238 
239  //use gamma times MRstar
240  return MR*mygamma;
241 }
242 
243 double
244  HLTRFilter::CalcR(double MR, TLorentzVector ja, TLorentzVector jb, edm::Handle<edm::View<reco::MET> > inputMet, const std::vector<math::XYZTLorentzVector>& muons){
245  //now we can calculate MTR
246  TVector3 met;
247  met.SetPtEtaPhi((inputMet->front()).pt(),0.0,(inputMet->front()).phi());
248 
249  std::vector<math::XYZTLorentzVector>::const_iterator muonIt;
250  for(muonIt = muons.begin(); muonIt!=muons.end(); muonIt++){
251  TVector3 tmp;
252  tmp.SetPtEtaPhi(muonIt->pt(),0,muonIt->phi());
253  met-=tmp;
254  }
255 
256  double MTR = sqrt(0.5*(met.Mag()*(ja.Pt()+jb.Pt()) - met.Dot(ja.Vect()+jb.Vect())));
257 
258  //filter events
259  return float(MTR)/float(MR); //R
260 
261 }
263 
#define LogDebug(id)
bool accept_NJ_
Definition: HLTRFilter.h:50
double R_MR_cut_
Definition: HLTRFilter.h:53
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:125
double min_MR_
Definition: HLTRFilter.h:48
double min_R_
Definition: HLTRFilter.h:47
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
double MR_offset_
Definition: HLTRFilter.h:52
~HLTRFilter() override
bool hltFilter(edm::Event &, const edm::EventSetup &, trigger::TriggerFilterObjectWithRefs &filterproduct) const override
Definition: HLTRFilter.cc:83
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: HLTRFilter.cc:61
std::vector< reco::MET > METCollection
collection of MET objects
Definition: METCollection.h:23
std::string encode() const
Definition: InputTag.cc:159
void addObjects(edm::Event &, trigger::TriggerFilterObjectWithRefs &filterproduct, double MR, double Rsq) const
Definition: HLTRFilter.cc:187
void addObject(int id, const reco::RecoEcalCandidateRef &ref)
setters for L3 collections: (id=physics type, and Ref<C>)
double R_offset_
Definition: HLTRFilter.h:51
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
Definition: MET.h:42
T sqrt(T t)
Definition: SSEVec.h:18
static double CalcMR(TLorentzVector ja, TLorentzVector jb)
Definition: HLTRFilter.cc:210
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
static double CalcR(double MR, TLorentzVector ja, TLorentzVector jb, edm::Handle< edm::View< reco::MET > > met, const std::vector< math::XYZTLorentzVector > &muons)
Definition: HLTRFilter.cc:244
HLTRFilter(const edm::ParameterSet &)
Definition: HLTRFilter.cc:27
static const std::string B
ParameterDescriptionBase * add(U const &iLabel, T const &value)
bool isValid() const
Definition: HandleBase.h:74
const std::string * moduleLabel() const
Definition: HLTFilter.cc:66
RefProd< PROD > getRefBeforePut()
Definition: Event.h:150
static void makeHLTFilterDescription(edm::ParameterSetDescription &desc)
Definition: HLTFilter.cc:29
edm::InputTag inputMetTag_
Definition: HLTRFilter.h:45
edm::EDGetTokenT< std::vector< math::XYZTLorentzVector > > m_theInputToken
Definition: HLTRFilter.h:42
met
===> hadronic RAZOR
edm::InputTag inputTag_
Definition: HLTRFilter.h:44
void addCollectionTag(const edm::InputTag &collectionTag)
collectionTags
void add(std::string const &label, ParameterSetDescription const &psetDescription)
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:37
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
edm::EDGetTokenT< edm::View< reco::MET > > m_theMETToken
Definition: HLTRFilter.h:43
fixed size matrix
bool saveTags() const
Definition: HLTFilter.h:45
HLT enums.
math::XYZPoint Point
point in the space
Definition: Candidate.h:41
def move(src, dest)
Definition: eostools.py:511