ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MuJetAnalysis/AnalysisTools/FitNtuple/src/FitNtuple.cc
Revision: 1.32
Committed: Thu Jun 27 08:16:37 2013 UTC (11 years, 10 months ago) by pakhotin
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.31: +15 -11 lines
Error occurred while calculating annotation data.
Log Message:
remove pT cut on PF isolation

File Contents

# Content
1 // -*- C++ -*-
2 //
3 // Package: FitNtuple
4 // Class: FitNtuple
5 //
6 /**\class FitNtuple FitNtuple.cc AnalysisTools/FitNtuple/src/FitNtuple.cc
7
8 Description: [one line class summary]
9
10 Implementation:
11 [Notes on implementation]
12 */
13 //
14 // Original Author: James Pivarski
15 // Created: Fri Nov 19 17:04:33 CST 2010
16 // $Id: FitNtuple.cc,v 1.31 2013/06/27 08:05:11 pakhotin Exp $
17 //
18 //
19
20 // system include files
21 #include <memory>
22
23 // user include files
24 #include "TTree.h"
25 #include "TRandom3.h"
26 #include "TMath.h"
27
28 #include "CommonTools/UtilAlgos/interface/TFileService.h"
29 #include "DataFormats/Candidate/interface/Candidate.h"
30 #include "DataFormats/DetId/interface/DetId.h"
31 #include "DataFormats/GeometrySurface/interface/Cylinder.h"
32 #include "DataFormats/GeometrySurface/interface/Plane.h"
33 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
34 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
35 #include "DataFormats/MuonReco/interface/MuonChamberMatch.h"
36 #include "DataFormats/MuonReco/interface/Muon.h"
37 #include "DataFormats/MuonReco/interface/MuonSegmentMatch.h"
38 #include "DataFormats/PatCandidates/interface/Muon.h"
39 #include "DataFormats/PatCandidates/interface/TriggerEvent.h"
40 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
41 #include "DataFormats/TrackReco/interface/HitPattern.h"
42 #include "DataFormats/TrackReco/interface/TrackFwd.h"
43 #include "DataFormats/TrackReco/interface/Track.h"
44 #include "DataFormats/VertexReco/interface/VertexFwd.h"
45 #include "DataFormats/VertexReco/interface/Vertex.h"
46 #include "FWCore/Framework/interface/EDAnalyzer.h"
47 #include "FWCore/Framework/interface/ESHandle.h"
48 #include "FWCore/Framework/interface/Event.h"
49 #include "FWCore/Framework/interface/Frameworkfwd.h"
50 #include "FWCore/Framework/interface/MakerMacros.h"
51 #include "FWCore/ParameterSet/interface/ParameterSet.h"
52 #include "FWCore/ServiceRegistry/interface/Service.h"
53 #include "FWCore/Utilities/interface/InputTag.h"
54 #include "Geometry/CommonDetUnit/interface/GlobalTrackingGeometry.h"
55 #include "Geometry/CommonTopologies/interface/PixelTopology.h"
56 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
57 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
58 #include "Geometry/TrackerGeometryBuilder/interface/PixelGeomDetUnit.h"
59 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
60 #include "MagneticField/Engine/interface/MagneticField.h"
61 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
62 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
63 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
64 #include "TrackingTools/Records/interface/TrackingComponentsRecord.h"
65 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
66 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
67 #include "TrackPropagation/SteppingHelixPropagator/interface/SteppingHelixPropagator.h"
68
69 #include "AnalysisDataFormats/MuJetAnalysis/interface/MultiMuon.h"
70
71 // class declaration
72 class FitNtuple : public edm::EDAnalyzer {
73
74 public:
75 explicit FitNtuple(const edm::ParameterSet&);
76 ~FitNtuple();
77
78 private:
79 virtual void beginJob() ;
80 virtual void analyze(const edm::Event&, const edm::EventSetup&);
81 virtual void endJob() ;
82
83 // ----------member data ---------------------------
84
85 // input parameters
86 edm::InputTag m_muons;
87 edm::InputTag m_muonsTriggerMatch;
88 edm::InputTag m_muJets;
89 edm::InputTag m_muJetOrphans;
90 double m_trigpt;
91 std::string m_dataset;
92
93 TRandom3 m_trandom3;
94
95 // all ntuples need these variables
96 Int_t m_run; // run lumi, and event number so that we can look at
97 Int_t m_lumi; // interesting cases in the event display
98 Int_t m_event; //
99
100 Int_t m_trigger; // satisfied trigger? 0 = failed all
101 // 2 = passed HLT_Mu9
102 // 4 = passed HLT_Mu11
103 // 8 = passed HLT_Mu15_v1
104 // 16 = passed HLT_Mu30_v*
105 // it is negative if there are no pT > 15 GeV/c, |eta| < 0.9 muons matched to L3
106 Int_t m_trigger15;
107 Int_t m_trigger20;
108 Int_t m_trigger24;
109 Int_t m_trigger30;
110 Int_t m_trigger40;
111 Int_t m_trigger40eta;
112
113
114 Int_t m_primaryverticessize;
115 Int_t m_generaltrackssize;
116
117 Int_t m_maxtrackspervtx;
118
119 Float_t m_muon1pt;
120 Float_t m_muon1eta;
121 Float_t m_muon2pt;
122 Float_t m_muon2eta;
123 Float_t m_muon3pt;
124 Float_t m_muon3eta;
125 Float_t m_muon4pt;
126 Float_t m_muon4eta;
127 Float_t m_muontrigpt;
128 Float_t m_muontrigeta;
129
130 // control sample "lowdimuon" (single mu-jet with two muons and low pT)
131 TTree *m_lowdimuon;
132 Int_t m_lowdimuon_containstrig;
133 Int_t m_lowdimuon_containstrig2;
134 Float_t m_lowdimuon_genmass;
135 Float_t m_lowdimuon_mass;
136 Float_t m_lowdimuon_corr_mass;
137 Float_t m_lowdimuon_recomass;
138 Float_t m_lowdimuon_pt;
139 Float_t m_lowdimuon_eta;
140 Float_t m_lowdimuon_phi;
141 Float_t m_lowdimuon_dr;
142 Float_t m_lowdimuon_pluspx;
143 Float_t m_lowdimuon_pluspy;
144 Float_t m_lowdimuon_pluspz;
145 Float_t m_lowdimuon_minuspx;
146 Float_t m_lowdimuon_minuspy;
147 Float_t m_lowdimuon_minuspz;
148 Float_t m_lowdimuon_vprob;
149 Float_t m_lowdimuon_vnchi2;
150 Float_t m_lowdimuon_vx;
151 Float_t m_lowdimuon_vy;
152 Float_t m_lowdimuon_vz;
153 Float_t m_lowdimuon_isoTk;
154 Float_t m_lowdimuon_isoTk_3mm;
155 Float_t m_lowdimuon_isoTk_2mm;
156 Float_t m_lowdimuon_isoTk_1mm;
157 Int_t m_lowdimuon_plusmatches;
158 Int_t m_lowdimuon_plushits;
159 Float_t m_lowdimuon_plusnormchi2;
160 Int_t m_lowdimuon_minusmatches;
161 Int_t m_lowdimuon_minushits;
162 Float_t m_lowdimuon_minusnormchi2;
163 Float_t m_lowdimuon_lxy;
164 Int_t m_lowdimuon_bbbarlike;
165
166 Float_t m_lowdimuon_trigpt1;
167 Float_t m_lowdimuon_trigeta1;
168 Float_t m_lowdimuon_trigpt2;
169 Float_t m_lowdimuon_trigeta2;
170
171 /*
172 Float_t m_lowdimuon_dphi100;
173 Float_t m_lowdimuon_dphi200;
174 Float_t m_lowdimuon_dphi300;
175 Float_t m_lowdimuon_dphi425;
176 Float_t m_lowdimuon_dphi470;
177 Float_t m_lowdimuon_dphi510;
178 Float_t m_lowdimuon_dphi565;
179 Float_t m_lowdimuon_dphi620;
180 Float_t m_lowdimuon_dphi670;
181 Float_t m_lowdimuon_dphi720;
182 Float_t m_lowdimuon_dphi800;
183 Float_t m_lowdimuon_dphi900;
184 */
185
186 Int_t m_lowdimuon_plus_thits;
187 Int_t m_lowdimuon_plus_mhits;
188 Int_t m_lowdimuon_plus_phits;
189
190 Int_t m_lowdimuon_minus_thits;
191 Int_t m_lowdimuon_minus_mhits;
192 Int_t m_lowdimuon_minus_phits;
193
194 Float_t m_lowdimuon_plus_dxy;
195 Float_t m_lowdimuon_minus_dxy;
196
197 Float_t m_lowdimuon_plus_pt;
198 Float_t m_lowdimuon_minus_pt;
199 Float_t m_lowdimuon_plus_eta;
200 Float_t m_lowdimuon_minus_eta;
201 Float_t m_lowdimuon_plus_phi;
202 Float_t m_lowdimuon_minus_phi;
203
204 Float_t m_lowdimuon_Deltaphi;
205
206 Float_t m_lowdimuon_dz;
207
208 // control sample "dimuorphan" (single, two-muon mu-jet and a single muon)
209 TTree *m_dimuorphan;
210 Float_t m_dimuorphan_deltaphi;
211 Float_t m_dimuorphan_orphanpt;
212 Float_t m_dimuorphan_orphaneta;
213 Float_t m_dimuorphan_orphanphi;
214 Int_t m_dimuorphan_orphanisglobal;
215 Int_t m_dimuorphan_orphanmatches;
216 Int_t m_dimuorphan_stationmask;
217 Int_t m_dimuorphan_orphanhits;
218 Float_t m_dimuorphan_orphanchi2;
219 Int_t m_dimuorphan_containstrig;
220 Int_t m_dimuorphan_containstrig2;
221 Float_t m_dimuorphan_mass;
222 Float_t m_dimuorphan_corr_mass;
223 Float_t m_dimuorphan_recomass;
224 Float_t m_dimuorphan_pt;
225 Float_t m_dimuorphan_eta;
226 Float_t m_dimuorphan_phi;
227 Float_t m_dimuorphan_dphi;
228 Float_t m_dimuorphan_dr;
229 Float_t m_dimuorphan_vprob;
230 Float_t m_dimuorphan_vnchi2;
231 Float_t m_dimuorphan_lxy;
232 Float_t m_dimuorphan_lxyz;
233 Float_t m_dimuorphan_caloiso;
234 Float_t m_dimuorphan_isoTk;
235 Float_t m_dimuorphan_isoTk_3mm;
236 Float_t m_dimuorphan_isoTk_2mm;
237 Float_t m_dimuorphan_isoTk_1mm;
238 Float_t m_dimuorphan_isoPF_1mm;
239 Float_t m_dimuorphan_orphan_isoTk_3mm;
240 Float_t m_dimuorphan_orphan_isoTk_2mm;
241 Float_t m_dimuorphan_orphan_isoTk_1mm;
242 Float_t m_dimuorphan_orphan_isoPF_1mm;
243 Int_t m_dimuorphan_dimuoncontainstrig;
244 Int_t m_dimuorphan_dimuonbbbarlike;
245
246 Float_t m_dimuorphan_trigpt1;
247 Float_t m_dimuorphan_trigeta1;
248 Float_t m_dimuorphan_trigpt2;
249 Float_t m_dimuorphan_trigeta2;
250
251 Int_t m_dimuorphan_trackdensity;
252
253 Int_t m_dimuorphan_plusmatches;
254 Int_t m_dimuorphan_plushits;
255 Float_t m_dimuorphan_plusnormchi2;
256 Int_t m_dimuorphan_minusmatches;
257 Int_t m_dimuorphan_minushits;
258 Float_t m_dimuorphan_minusnormchi2;
259
260 Int_t m_dimuorphan_plusstationmask;
261 Int_t m_dimuorphan_minusstationmask;
262
263 Float_t m_dimuorphan_plus_qoverpError;
264 Float_t m_dimuorphan_plus_ptError;
265 Float_t m_dimuorphan_plus_phiError;
266 Float_t m_dimuorphan_plus_etaError;
267 Int_t m_dimuorphan_plus_isGlobal;
268 Int_t m_dimuorphan_plus_isStandAlone;
269 Int_t m_dimuorphan_plus_GlobalHits;
270 Float_t m_dimuorphan_plus_GlobalChi2;
271 Float_t m_dimuorphan_plus_StandAloneHits;
272 Float_t m_dimuorphan_plus_StandAloneChi2;
273
274 Float_t m_dimuorphan_minus_qoverpError;
275 Float_t m_dimuorphan_minus_ptError;
276 Float_t m_dimuorphan_minus_phiError;
277 Float_t m_dimuorphan_minus_etaError;
278 Int_t m_dimuorphan_minus_isGlobal;
279 Int_t m_dimuorphan_minus_isStandAlone;
280 Int_t m_dimuorphan_minus_GlobalHits;
281 Float_t m_dimuorphan_minus_GlobalChi2;
282 Float_t m_dimuorphan_minus_StandAloneHits;
283 Float_t m_dimuorphan_minus_StandAloneChi2;
284
285 /*
286 Float_t m_dimuorphan_dphi100;
287 Float_t m_dimuorphan_dphi200;
288 Float_t m_dimuorphan_dphi300;
289 Float_t m_dimuorphan_dphi425;
290 Float_t m_dimuorphan_dphi470;
291 Float_t m_dimuorphan_dphi510;
292 Float_t m_dimuorphan_dphi565;
293 Float_t m_dimuorphan_dphi620;
294 Float_t m_dimuorphan_dphi670;
295 Float_t m_dimuorphan_dphi720;
296 Float_t m_dimuorphan_dphi800;
297 Float_t m_dimuorphan_dphi900;
298 */
299
300 Int_t m_dimuorphan_plus_thits;
301 Int_t m_dimuorphan_plus_mhits;
302 Int_t m_dimuorphan_plus_phits;
303
304 Int_t m_dimuorphan_minus_thits;
305 Int_t m_dimuorphan_minus_mhits;
306 Int_t m_dimuorphan_minus_phits;
307
308 Float_t m_dimuorphan_plus_dxy;
309 Float_t m_dimuorphan_minus_dxy;
310
311 Float_t m_dimuorphan_plus_pt;
312 Float_t m_dimuorphan_minus_pt;
313 Float_t m_dimuorphan_plus_eta;
314 Float_t m_dimuorphan_minus_eta;
315 Float_t m_dimuorphan_plus_phi;
316 Float_t m_dimuorphan_minus_phi;
317
318 Float_t m_dimuorphan_Deltaphi;
319 Float_t m_dimuorphan_Deltaphi_orphan;
320
321 Float_t m_dimuorphan_dz1;
322 Float_t m_dimuorphan_dz2;
323 Float_t m_dimuorphan_deltaz;
324
325 Int_t m_dimuorphan_muonssize;
326
327 Float_t m_dimuorphan_dr1;
328 Float_t m_dimuorphan_dr2;
329
330 // signal region "dimudimu" (two mu-jets with two muons each)
331 TTree *m_dimudimu;
332 Int_t m_dimudimu_orphans;
333 Float_t m_dimudimu_wholemass;
334 Float_t m_dimudimu_wholept;
335 Float_t m_dimudimu_deltaphi;
336 Int_t m_dimudimu_containstrig;
337 Int_t m_dimudimu_containstrig2;
338 Float_t m_dimudimu_massC;
339 Float_t m_dimudimu_ptC;
340 Float_t m_dimudimu_etaC;
341 Float_t m_dimudimu_phiC;
342 Float_t m_dimudimu_drC;
343 Float_t m_dimudimu_vprobC;
344 Float_t m_dimudimu_lxyC;
345 Float_t m_dimudimu_lxyzC;
346 Float_t m_dimudimu_C_isoTk;
347 Float_t m_dimudimu_C_isoTk_3mm;
348 Float_t m_dimudimu_C_isoTk_2mm;
349 Float_t m_dimudimu_C_isoTk_1mm;
350 Float_t m_dimudimu_C_isoPF_1mm;
351 Float_t m_dimudimu_massF;
352 Float_t m_dimudimu_ptF;
353 Float_t m_dimudimu_etaF;
354 Float_t m_dimudimu_phiF;
355 Float_t m_dimudimu_drF;
356 Float_t m_dimudimu_vprobF;
357 Float_t m_dimudimu_lxyF;
358 Float_t m_dimudimu_lxyzF;
359 Float_t m_dimudimu_F_isoTk;
360 Float_t m_dimudimu_F_isoTk_3mm;
361 Float_t m_dimudimu_F_isoTk_2mm;
362 Float_t m_dimudimu_F_isoTk_1mm;
363 Float_t m_dimudimu_F_isoPF_1mm;
364
365 Float_t m_dimudimu_dz1;
366 Float_t m_dimudimu_dz2;
367 Float_t m_dimudimu_deltaz;
368
369 };
370
371 //
372 // constants, enums and typedefs
373 //
374
375 //
376 // static data member definitions
377 //
378
379 //
380 // constructors and destructor
381 //
382 FitNtuple::FitNtuple(const edm::ParameterSet& iConfig)
383 : m_muons(iConfig.getParameter<edm::InputTag>("muons")),
384 m_muonsTriggerMatch(iConfig.getParameter<edm::InputTag>("muonsTriggerMatch")),
385 m_muJets(iConfig.getParameter<edm::InputTag>("muJets")),
386 m_muJetOrphans(iConfig.getParameter<edm::InputTag>("muJetOrphans")),
387 m_trigpt(iConfig.getParameter<double>("trigpt")),
388 m_dataset(iConfig.getParameter<std::string>("dataset"))
389 {
390 //now do what ever initialization is needed
391 edm::Service<TFileService> tFile;
392
393 m_lowdimuon = tFile->make<TTree>("lowdimuon", "lowdimuon");
394 m_lowdimuon->Branch("run", &m_run, "run/I");
395 m_lowdimuon->Branch("lumi", &m_lumi, "lumi/I");
396 m_lowdimuon->Branch("event", &m_event, "event/I");
397 m_lowdimuon->Branch("trigger", &m_trigger, "trigger/I");
398 m_lowdimuon->Branch("trigger15", &m_trigger15, "trigger15/I");
399 m_lowdimuon->Branch("trigger20", &m_trigger20, "trigger20/I");
400 m_lowdimuon->Branch("trigger24", &m_trigger24, "trigger24/I");
401 m_lowdimuon->Branch("trigger30", &m_trigger30, "trigger30/I");
402 m_lowdimuon->Branch("trigger40", &m_trigger40, "trigger40/I");
403 m_lowdimuon->Branch("trigger40eta", &m_trigger40eta, "trigger40eta/I");
404 m_lowdimuon->Branch("muon1pt", &m_muon1pt, "muon1pt/F");
405 m_lowdimuon->Branch("muon1eta", &m_muon1eta, "muon1eta/F");
406 m_lowdimuon->Branch("muon2pt", &m_muon2pt, "muon2pt/F");
407 m_lowdimuon->Branch("muon2eta", &m_muon2eta, "muon2eta/F");
408 m_lowdimuon->Branch("muontrigpt", &m_muontrigpt, "muontrigpt/F");
409 m_lowdimuon->Branch("muontrigeta", &m_muontrigeta, "muontrigeta/F");
410 m_lowdimuon->Branch("containstrig", &m_lowdimuon_containstrig, "containstrig/I");
411 m_lowdimuon->Branch("containstrig2", &m_lowdimuon_containstrig2, "containstrig2/I");
412 m_lowdimuon->Branch("genmass", &m_lowdimuon_genmass, "genmass/F");
413 m_lowdimuon->Branch("mass", &m_lowdimuon_mass, "mass/F");
414 m_lowdimuon->Branch("corr_mass", &m_lowdimuon_corr_mass, "corr_mass/F");
415 m_lowdimuon->Branch("recomass", &m_lowdimuon_recomass, "recomass/F");
416 m_lowdimuon->Branch("pt", &m_lowdimuon_pt, "pt/F");
417 m_lowdimuon->Branch("eta", &m_lowdimuon_eta, "eta/F");
418 m_lowdimuon->Branch("phi", &m_lowdimuon_phi, "phi/F");
419 m_lowdimuon->Branch("dr", &m_lowdimuon_dr, "dr/F");
420 m_lowdimuon->Branch("pluspx", &m_lowdimuon_pluspx, "pluspx/F");
421 m_lowdimuon->Branch("pluspy", &m_lowdimuon_pluspy, "pluspy/F");
422 m_lowdimuon->Branch("pluspz", &m_lowdimuon_pluspz, "pluspz/F");
423 m_lowdimuon->Branch("minuspx", &m_lowdimuon_minuspx, "minuspx/F");
424 m_lowdimuon->Branch("minuspy", &m_lowdimuon_minuspy, "minuspy/F");
425 m_lowdimuon->Branch("minuspz", &m_lowdimuon_minuspz, "minuspz/F");
426 m_lowdimuon->Branch("vprob", &m_lowdimuon_vprob, "vprob/F");
427 m_lowdimuon->Branch("vnchi2", &m_lowdimuon_vnchi2, "vnchi2/F");
428 m_lowdimuon->Branch("vx", &m_lowdimuon_vx, "vx/F");
429 m_lowdimuon->Branch("vy", &m_lowdimuon_vy, "vy/F");
430 m_lowdimuon->Branch("vz", &m_lowdimuon_vz, "vz/F");
431 m_lowdimuon->Branch("isoTk", &m_lowdimuon_isoTk, "isoTk/F");
432 m_lowdimuon->Branch("isoTk_3mm", &m_lowdimuon_isoTk_3mm, "isoTk_3mm/F");
433 m_lowdimuon->Branch("isoTk_2mm", &m_lowdimuon_isoTk_2mm, "isoTk_2mm/F");
434 m_lowdimuon->Branch("isoTk_1mm", &m_lowdimuon_isoTk_1mm, "isoTk_1mm/F");
435 m_lowdimuon->Branch("plusmatches", &m_lowdimuon_plusmatches, "plusmatches/I");
436 m_lowdimuon->Branch("plushits", &m_lowdimuon_plushits, "plushits/I");
437 m_lowdimuon->Branch("plusnormchi2", &m_lowdimuon_plusnormchi2, "plusnormchi2/F");
438 m_lowdimuon->Branch("minusmatches", &m_lowdimuon_minusmatches, "minusmatches/I");
439 m_lowdimuon->Branch("minushits", &m_lowdimuon_minushits, "minushits/I");
440 m_lowdimuon->Branch("minusnormchi2", &m_lowdimuon_minusnormchi2, "minusnormchi2/F");
441 m_lowdimuon->Branch("lxy", &m_lowdimuon_lxy, "lxy/F");
442 m_lowdimuon->Branch("bbbarlike", &m_lowdimuon_bbbarlike, "bbbarlike/I");
443
444 m_lowdimuon->Branch("trigpt1", &m_lowdimuon_trigpt1, "trigpt1/F");
445 m_lowdimuon->Branch("trigeta1", &m_lowdimuon_trigeta1, "trigeta1/F");
446 m_lowdimuon->Branch("trigpt2", &m_lowdimuon_trigpt2, "trigpt2/F");
447 m_lowdimuon->Branch("trigeta2", &m_lowdimuon_trigeta2, "trigeta2/F");
448
449
450 /*
451 m_lowdimuon->Branch("dphi100", &m_lowdimuon_dphi100, "dphi100/F");
452 m_lowdimuon->Branch("dphi200", &m_lowdimuon_dphi200, "dphi200/F");
453 m_lowdimuon->Branch("dphi300", &m_lowdimuon_dphi300, "dphi300/F");
454 m_lowdimuon->Branch("dphi425", &m_lowdimuon_dphi425, "dphi425/F");
455 m_lowdimuon->Branch("dphi470", &m_lowdimuon_dphi470, "dphi470/F");
456 m_lowdimuon->Branch("dphi510", &m_lowdimuon_dphi510, "dphi510/F");
457 m_lowdimuon->Branch("dphi565", &m_lowdimuon_dphi565, "dphi565/F");
458 m_lowdimuon->Branch("dphi620", &m_lowdimuon_dphi620, "dphi620/F");
459 m_lowdimuon->Branch("dphi670", &m_lowdimuon_dphi670, "dphi670/F");
460 m_lowdimuon->Branch("dphi720", &m_lowdimuon_dphi720, "dphi720/F");
461 m_lowdimuon->Branch("dphi800", &m_lowdimuon_dphi800, "dphi800/F");
462 m_lowdimuon->Branch("dphi900", &m_lowdimuon_dphi900, "dphi900/F");
463 */
464
465 m_lowdimuon->Branch("plus_thits", &m_lowdimuon_plus_thits, "plus_thits/I");
466 m_lowdimuon->Branch("plus_mhits", &m_lowdimuon_plus_mhits, "plus_mhits/I");
467 m_lowdimuon->Branch("plus_phits", &m_lowdimuon_plus_phits, "plus_phits/I");
468
469 m_lowdimuon->Branch("minus_thits", &m_lowdimuon_minus_thits, "minus_thits/I");
470 m_lowdimuon->Branch("minus_mhits", &m_lowdimuon_minus_mhits, "minus_mhits/I");
471 m_lowdimuon->Branch("minus_phits", &m_lowdimuon_minus_phits, "minus_phits/I");
472
473 m_lowdimuon->Branch("plus_dxy", &m_lowdimuon_plus_dxy, "plus_dxy/F");
474 m_lowdimuon->Branch("minus_dxy", &m_lowdimuon_minus_dxy, "minus_dxy/F");
475
476 m_lowdimuon->Branch("plus_pt", &m_lowdimuon_plus_pt, "plus_pt/F");
477 m_lowdimuon->Branch("minus_pt", &m_lowdimuon_minus_pt, "minus_pt/F");
478 m_lowdimuon->Branch("plus_eta", &m_lowdimuon_plus_eta, "plus_eta/F");
479 m_lowdimuon->Branch("minus_eta", &m_lowdimuon_minus_eta, "minus_eta/F");
480 m_lowdimuon->Branch("plus_phi", &m_lowdimuon_plus_phi, "plus_phi/F");
481 m_lowdimuon->Branch("minus_phi", &m_lowdimuon_minus_phi, "minus_phi/F");
482
483 m_lowdimuon->Branch("Deltaphi", &m_lowdimuon_Deltaphi, "Deltaphi/F");
484
485 m_lowdimuon->Branch("dz", &m_lowdimuon_dz, "dz/F");
486
487 m_lowdimuon->Branch("pvsize", &m_primaryverticessize, "pvsize/I");
488 m_lowdimuon->Branch("trackssize", &m_generaltrackssize, "trackssize/I");
489 m_lowdimuon->Branch("maxtrackspervtx", &m_maxtrackspervtx, "maxtrackspervtx/I");
490
491 m_dimuorphan = tFile->make<TTree>("dimuorphan", "dimuorphan");
492 m_dimuorphan->Branch("run", &m_run, "run/I");
493 m_dimuorphan->Branch("lumi", &m_lumi, "lumi/I");
494 m_dimuorphan->Branch("event", &m_event, "event/I");
495 m_dimuorphan->Branch("trigger", &m_trigger, "trigger/I");
496 m_dimuorphan->Branch("trigger15", &m_trigger15, "trigger15/I");
497 m_dimuorphan->Branch("trigger20", &m_trigger20, "trigger20/I");
498 m_dimuorphan->Branch("trigger24", &m_trigger24, "trigger24/I");
499 m_dimuorphan->Branch("trigger30", &m_trigger30, "trigger30/I");
500 m_dimuorphan->Branch("trigger40", &m_trigger40, "trigger40/I");
501 m_dimuorphan->Branch("trigger40eta", &m_trigger40eta, "trigger40eta/I");
502 m_dimuorphan->Branch("muon1pt", &m_muon1pt, "muon1pt/F");
503 m_dimuorphan->Branch("muon1eta", &m_muon1eta, "muon1eta/F");
504 m_dimuorphan->Branch("muon2pt", &m_muon2pt, "muon2pt/F");
505 m_dimuorphan->Branch("muon2eta", &m_muon2eta, "muon2eta/F");
506 m_dimuorphan->Branch("muon3pt", &m_muon3pt, "muon3pt/F");
507 m_dimuorphan->Branch("muon3eta", &m_muon3eta, "muon3eta/F");
508 m_dimuorphan->Branch("muon4pt", &m_muon4pt, "muon4pt/F");
509 m_dimuorphan->Branch("muon4eta", &m_muon4eta, "muon4eta/F");
510 m_dimuorphan->Branch("muontrigpt", &m_muontrigpt, "muontrigpt/F");
511 m_dimuorphan->Branch("muontrigeta", &m_muontrigeta, "muontrigeta/F");
512 m_dimuorphan->Branch("deltaphi", &m_dimuorphan_deltaphi, "deltaphi/F");
513 m_dimuorphan->Branch("orphanpt", &m_dimuorphan_orphanpt, "orphanpt/F");
514 m_dimuorphan->Branch("orphaneta", &m_dimuorphan_orphaneta, "orphaneta/F");
515 m_dimuorphan->Branch("orphanphi", &m_dimuorphan_orphanphi, "orphanphi/F");
516 m_dimuorphan->Branch("orphanisglobal", &m_dimuorphan_orphanisglobal, "orphanisglobal/I");
517 m_dimuorphan->Branch("orphanmatches", &m_dimuorphan_orphanmatches, "orphan/I");
518 m_dimuorphan->Branch("orphanstationmask", &m_dimuorphan_stationmask, "stationmask/I");
519 m_dimuorphan->Branch("orphanhits", &m_dimuorphan_orphanhits, "orphanhits/I");
520 m_dimuorphan->Branch("orphanchi2", &m_dimuorphan_orphanchi2, "orphanchi2/F");
521 m_dimuorphan->Branch("containstrig", &m_dimuorphan_containstrig, "containstrig/I");
522 m_dimuorphan->Branch("containstrig2", &m_dimuorphan_containstrig2, "containstrig2/I");
523 m_dimuorphan->Branch("mass", &m_dimuorphan_mass, "mass/F");
524 m_dimuorphan->Branch("corr_mass", &m_dimuorphan_corr_mass, "corr_mass/F");
525 m_dimuorphan->Branch("recomass", &m_dimuorphan_recomass, "recomass/F");
526 m_dimuorphan->Branch("pt", &m_dimuorphan_pt, "pt/F");
527 m_dimuorphan->Branch("eta", &m_dimuorphan_eta, "eta/F");
528 m_dimuorphan->Branch("phi", &m_dimuorphan_phi, "phi/F");
529 m_dimuorphan->Branch("dphi", &m_dimuorphan_dphi, "dphi/F");
530 m_dimuorphan->Branch("dr", &m_dimuorphan_dr, "dr/F");
531 m_dimuorphan->Branch("vprob", &m_dimuorphan_vprob, "vprob/F");
532 m_dimuorphan->Branch("vnchi2", &m_dimuorphan_vnchi2, "vnchi2/F");
533 m_dimuorphan->Branch("lxy", &m_dimuorphan_lxy, "lxy/F");
534 m_dimuorphan->Branch("lxyz", &m_dimuorphan_lxyz, "lxyz/F");
535 m_dimuorphan->Branch("caloiso", &m_dimuorphan_caloiso, "caloiso/F");
536 m_dimuorphan->Branch("iso", &m_dimuorphan_isoTk, "isoTk/F");
537 m_dimuorphan->Branch("isoTk_3mm", &m_dimuorphan_isoTk_3mm, "isoTk_3mm/F");
538 m_dimuorphan->Branch("isoTk_2mm", &m_dimuorphan_isoTk_2mm, "isoTk_2mm/F");
539 m_dimuorphan->Branch("isoTk_1mm", &m_dimuorphan_isoTk_1mm, "isoTk_1mm/F");
540 m_dimuorphan->Branch("isoPF_1mm", &m_dimuorphan_isoPF_1mm, "isoPF_1mm/F");
541 m_dimuorphan->Branch("orphan_isoTk_3mm", &m_dimuorphan_orphan_isoTk_3mm, "orphan_isoTk_3mm/F");
542 m_dimuorphan->Branch("orphan_isoTk_2mm", &m_dimuorphan_orphan_isoTk_2mm, "orphan_isoTk_2mm/F");
543 m_dimuorphan->Branch("orphan_isoTk_1mm", &m_dimuorphan_orphan_isoTk_1mm, "orphan_isoTk_1mm/F");
544 m_dimuorphan->Branch("orphan_isoPF_1mm", &m_dimuorphan_orphan_isoPF_1mm, "orphan_isoPF_1mm/F");
545 m_dimuorphan->Branch("dimuoncontainstrig", &m_dimuorphan_dimuoncontainstrig, "dimuoncontainstrig/I");
546 m_dimuorphan->Branch("dimuonbbbarlike", &m_dimuorphan_dimuonbbbarlike, "eta/I");
547
548 m_dimuorphan->Branch("trigpt1", &m_dimuorphan_trigpt1, "trigpt1/F");
549 m_dimuorphan->Branch("trigeta1", &m_dimuorphan_trigeta1, "trigeta1/F");
550 m_dimuorphan->Branch("trigpt2", &m_dimuorphan_trigpt2, "trigpt2/F");
551 m_dimuorphan->Branch("trigeta2", &m_dimuorphan_trigeta2, "trigeta2/F");
552
553 m_dimuorphan->Branch("plusmatches", &m_dimuorphan_plusmatches, "plusmatches/I");
554 m_dimuorphan->Branch("plushits", &m_dimuorphan_plushits, "plushits/I");
555 m_dimuorphan->Branch("plusnormchi2", &m_dimuorphan_plusnormchi2, "plusnormchi2/F");
556 m_dimuorphan->Branch("minusmatches", &m_dimuorphan_minusmatches, "minusmatches/I");
557 m_dimuorphan->Branch("minushits", &m_dimuorphan_minushits, "minushits/I");
558 m_dimuorphan->Branch("minusnormchi2", &m_dimuorphan_minusnormchi2, "minusnormchi2/F");
559
560 m_dimuorphan->Branch("plusstationmask", &m_dimuorphan_plusstationmask, "plusstationmask/I");
561 m_dimuorphan->Branch("minusstationmask", &m_dimuorphan_minusstationmask, "minusstationmask/I");
562
563 m_dimuorphan->Branch("plus_qoverpError", &m_dimuorphan_plus_qoverpError, "plus_qoverpError/F");
564 m_dimuorphan->Branch("plus_ptError", &m_dimuorphan_plus_ptError, "plus_ptError/F");
565 m_dimuorphan->Branch("plus_phiError", &m_dimuorphan_plus_phiError, "plus_phiError/F");
566 m_dimuorphan->Branch("plus_etaError", &m_dimuorphan_plus_etaError, "plus_etaError/F");
567 m_dimuorphan->Branch("plus_isGlobal", &m_dimuorphan_plus_isGlobal, "plus_isGlobal/I");
568 m_dimuorphan->Branch("plus_isStandAlone", &m_dimuorphan_plus_isStandAlone, "plus_isStandAlone/I");
569 m_dimuorphan->Branch("plus_GlobalHits", &m_dimuorphan_plus_GlobalHits, "plus_GlobalHits/I");
570 m_dimuorphan->Branch("plus_GlobalChi2", &m_dimuorphan_plus_GlobalChi2, "plus_GlobalChi2/F");
571 m_dimuorphan->Branch("plus_StandAloneHits", &m_dimuorphan_plus_StandAloneHits, "plus_StandAloneHits/I");
572 m_dimuorphan->Branch("plus_StandAloneChi2", &m_dimuorphan_plus_StandAloneChi2, "plus_StandAloneChi2/F");
573
574 m_dimuorphan->Branch("minus_qoverpError", &m_dimuorphan_minus_qoverpError, "minus_qoverpError/F");
575 m_dimuorphan->Branch("minus_ptError", &m_dimuorphan_minus_ptError, "minus_ptError/F");
576 m_dimuorphan->Branch("minus_phiError", &m_dimuorphan_minus_phiError, "minus_phiError/F");
577 m_dimuorphan->Branch("minus_etaError", &m_dimuorphan_minus_etaError, "minus_etaError/F");
578 m_dimuorphan->Branch("minus_isGlobal", &m_dimuorphan_minus_isGlobal, "minus_isGlobal/I");
579 m_dimuorphan->Branch("minus_isStandAlone", &m_dimuorphan_minus_isStandAlone, "minus_isStandAlone/I");
580 m_dimuorphan->Branch("minus_GlobalHits", &m_dimuorphan_minus_GlobalHits, "minus_GlobalHits/I");
581 m_dimuorphan->Branch("minus_GlobalChi2", &m_dimuorphan_minus_GlobalChi2, "minus_GlobalChi2/F");
582 m_dimuorphan->Branch("minus_StandAloneHits", &m_dimuorphan_minus_StandAloneHits, "minus_StandAloneHits/I");
583 m_dimuorphan->Branch("minus_StandAloneChi2", &m_dimuorphan_minus_StandAloneChi2, "minus_StandAloneChi2/F");
584
585 /*
586 m_dimuorphan->Branch("dphi100", &m_dimuorphan_dphi100, "dphi100/F");
587 m_dimuorphan->Branch("dphi200", &m_dimuorphan_dphi200, "dphi200/F");
588 m_dimuorphan->Branch("dphi300", &m_dimuorphan_dphi300, "dphi300/F");
589 m_dimuorphan->Branch("dphi425", &m_dimuorphan_dphi425, "dphi425/F");
590 m_dimuorphan->Branch("dphi470", &m_dimuorphan_dphi470, "dphi470/F");
591 m_dimuorphan->Branch("dphi510", &m_dimuorphan_dphi510, "dphi510/F");
592 m_dimuorphan->Branch("dphi565", &m_dimuorphan_dphi565, "dphi565/F");
593 m_dimuorphan->Branch("dphi620", &m_dimuorphan_dphi620, "dphi620/F");
594 m_dimuorphan->Branch("dphi670", &m_dimuorphan_dphi670, "dphi670/F");
595 m_dimuorphan->Branch("dphi720", &m_dimuorphan_dphi720, "dphi720/F");
596 m_dimuorphan->Branch("dphi800", &m_dimuorphan_dphi800, "dphi800/F");
597 m_dimuorphan->Branch("dphi900", &m_dimuorphan_dphi900, "dphi900/F");
598 */
599
600 m_dimuorphan->Branch("plus_thits", &m_dimuorphan_plus_thits, "plus_thits/I");
601 m_dimuorphan->Branch("plus_mhits", &m_dimuorphan_plus_mhits, "plus_mhits/I");
602 m_dimuorphan->Branch("plus_phits", &m_dimuorphan_plus_phits, "plus_phits/I");
603
604 m_dimuorphan->Branch("minus_thits", &m_dimuorphan_minus_thits, "minus_thits/I");
605 m_dimuorphan->Branch("minus_mhits", &m_dimuorphan_minus_mhits, "minus_mhits/I");
606 m_dimuorphan->Branch("minus_phits", &m_dimuorphan_minus_phits, "minus_phits/I");
607
608 m_dimuorphan->Branch("plus_dxy", &m_dimuorphan_plus_dxy, "plus_dxy/F");
609 m_dimuorphan->Branch("minus_dxy", &m_dimuorphan_minus_dxy, "minus_dxy/F");
610
611 m_dimuorphan->Branch("plus_pt", &m_dimuorphan_plus_pt, "plus_pt/F");
612 m_dimuorphan->Branch("minus_pt", &m_dimuorphan_minus_pt, "minus_pt/F");
613 m_dimuorphan->Branch("plus_eta", &m_dimuorphan_plus_eta, "plus_eta/F");
614 m_dimuorphan->Branch("minus_eta", &m_dimuorphan_minus_eta, "minus_eta/F");
615 m_dimuorphan->Branch("plus_phi", &m_dimuorphan_plus_phi, "plus_phi/F");
616 m_dimuorphan->Branch("minus_phi", &m_dimuorphan_minus_phi, "minus_phi/F");
617
618 m_dimuorphan->Branch("Deltaphi", &m_dimuorphan_Deltaphi, "Deltaphi/F");
619 m_dimuorphan->Branch("Deltaphi_orphan", &m_dimuorphan_Deltaphi_orphan, "Deltaphi_orphan/F");
620
621 m_dimuorphan->Branch("dz1", &m_dimuorphan_dz1, "dz1/F");
622 m_dimuorphan->Branch("dz2", &m_dimuorphan_dz2, "dz2/F");
623 m_dimuorphan->Branch("deltaz", &m_dimuorphan_deltaz, "deltaz/F");
624
625 m_dimuorphan->Branch("muonssize", &m_dimuorphan_muonssize, "muonssize/I");
626
627 m_dimuorphan->Branch("pvsize", &m_primaryverticessize, "pvsize/I");
628 m_dimuorphan->Branch("trackssize", &m_generaltrackssize, "trackssize/I");
629 m_dimuorphan->Branch("maxtrackspervtx", &m_maxtrackspervtx, "maxtrackspervtx/I");
630
631 m_dimuorphan->Branch("dr1", &m_dimuorphan_dr1, "dr1/F");
632 m_dimuorphan->Branch("dr2", &m_dimuorphan_dr2, "dr2/F");
633
634 m_dimuorphan->Branch("trackdensity", &m_dimuorphan_trackdensity, "trackdensity/I");
635
636 m_dimudimu = tFile->make<TTree>("dimudimu", "dimudimu");
637 m_dimudimu->Branch("run", &m_run, "run/I");
638 m_dimudimu->Branch("lumi", &m_lumi, "lumi/I");
639 m_dimudimu->Branch("event", &m_event, "event/I");
640 m_dimudimu->Branch("trigger", &m_trigger, "trigger/I");
641 m_dimudimu->Branch("muon1pt", &m_muon1pt, "muon1pt/F");
642 m_dimudimu->Branch("muon1eta", &m_muon1eta, "muon1eta/F");
643 m_dimudimu->Branch("muon2pt", &m_muon2pt, "muon2pt/F");
644 m_dimudimu->Branch("muon2eta", &m_muon2eta, "muon2eta/F");
645 m_dimudimu->Branch("muon3pt", &m_muon3pt, "muon3pt/F");
646 m_dimudimu->Branch("muon3eta", &m_muon3eta, "muon3eta/F");
647 m_dimudimu->Branch("muon4pt", &m_muon4pt, "muon4pt/F");
648 m_dimudimu->Branch("muon4eta", &m_muon4eta, "muon4eta/F");
649 m_dimudimu->Branch("muontrigpt", &m_muontrigpt, "muontrigpt/F");
650 m_dimudimu->Branch("muontrigeta", &m_muontrigeta, "muontrigeta/F");
651 m_dimudimu->Branch("orphans", &m_dimudimu_orphans, "orphans/I");
652 m_dimudimu->Branch("wholemass", &m_dimudimu_wholemass, "wholemass/F");
653 m_dimudimu->Branch("wholept", &m_dimudimu_wholept, "wholept/F");
654 m_dimudimu->Branch("deltaphi", &m_dimudimu_deltaphi, "deltaphi/F");
655 m_dimudimu->Branch("containstrig", &m_dimudimu_containstrig, "containstrig/I");
656 m_dimudimu->Branch("containstrig2", &m_dimudimu_containstrig2, "containstrig2/I");
657 m_dimudimu->Branch("massC", &m_dimudimu_massC, "massC/F");
658 m_dimudimu->Branch("ptC", &m_dimudimu_ptC, "ptC/F");
659 m_dimudimu->Branch("etaC", &m_dimudimu_etaC, "etaC/F");
660 m_dimudimu->Branch("phiC", &m_dimudimu_phiC, "phiC/F");
661 m_dimudimu->Branch("drC", &m_dimudimu_drC, "drC/F");
662 m_dimudimu->Branch("vprobC", &m_dimudimu_vprobC, "vprobC/F");
663 m_dimudimu->Branch("lxyC", &m_dimudimu_lxyC, "lxyC/F");
664 m_dimudimu->Branch("lxyzC", &m_dimudimu_lxyzC, "lxyzC/F");
665 m_dimudimu->Branch("C_isoTk", &m_dimudimu_C_isoTk, "C_isoTk/F");
666 m_dimudimu->Branch("C_isoTk_3mm", &m_dimudimu_C_isoTk_3mm, "C_isoTk_3mm/F");
667 m_dimudimu->Branch("C_isoTk_2mm", &m_dimudimu_C_isoTk_2mm, "C_isoTk_2mm/F");
668 m_dimudimu->Branch("C_isoTk_1mm", &m_dimudimu_C_isoTk_1mm, "C_isoTk_1mm/F");
669 m_dimudimu->Branch("C_isoPF_1mm", &m_dimudimu_C_isoPF_1mm, "C_isoPF_1mm/F");
670 m_dimudimu->Branch("massF", &m_dimudimu_massF, "massF/F");
671 m_dimudimu->Branch("ptF", &m_dimudimu_ptF, "ptF/F");
672 m_dimudimu->Branch("etaF", &m_dimudimu_etaF, "etaF/F");
673 m_dimudimu->Branch("phiF", &m_dimudimu_phiF, "phiF/F");
674 m_dimudimu->Branch("drF", &m_dimudimu_drF, "drF/F");
675 m_dimudimu->Branch("vprobF", &m_dimudimu_vprobF, "vprobF/F");
676 m_dimudimu->Branch("lxyF", &m_dimudimu_lxyF, "lxyF/F");
677 m_dimudimu->Branch("lxyzF", &m_dimudimu_lxyzF, "lxyzF/F");
678 m_dimudimu->Branch("F_isoTk", &m_dimudimu_F_isoTk, "F_isoTk/F");
679 m_dimudimu->Branch("F_isoTk_3mm", &m_dimudimu_F_isoTk_3mm, "F_isoTk_3mm/F");
680 m_dimudimu->Branch("F_isoTk_2mm", &m_dimudimu_F_isoTk_2mm, "F_isoTk_2mm/F");
681 m_dimudimu->Branch("F_isoTk_1mm", &m_dimudimu_F_isoTk_1mm, "F_isoTk_1mm/F");
682 m_dimudimu->Branch("F_isoPF_1mm", &m_dimudimu_F_isoPF_1mm, "F_isoPF_1mm/F");
683
684 m_dimudimu->Branch("dz1", &m_dimudimu_dz1, "dz1/F");
685 m_dimudimu->Branch("dz2", &m_dimudimu_dz2, "dz2/F");
686 m_dimudimu->Branch("deltaz", &m_dimudimu_deltaz, "deltaz/F");
687
688 m_dimudimu->Branch("pvsize", &m_primaryverticessize, "pvsize/I");
689 m_dimudimu->Branch("trackssize", &m_generaltrackssize, "trackssize/I");
690 m_dimudimu->Branch("maxtrackspervtx", &m_maxtrackspervtx, "maxtrackspervtx/I");
691 }
692
693
694 FitNtuple::~FitNtuple()
695 {
696
697 // do anything here that needs to be done at desctruction time
698 // (e.g. close files, deallocate resources etc.)
699 }
700
701
702 //
703 // member functions
704 //
705
706
707 double scalePt(double pt, double eta, double phi, int charge) {
708
709 double b = -5.03313e-6;
710 double c = -4.41463e-5;
711 double d0 = -0.000148871;
712 double e0 = 1.59501;
713 double d1 = 7.95495e-05;
714 double e1 = -0.364823;
715 double d2 = 0.000152032;
716 double e2 = 0.410195;
717
718 double d = 0;
719 double e = 0;
720
721 double signeta = 0;
722
723 if (fabs(eta) <= 0.9) {
724 d = d0; e = e0;
725 }
726 if (eta > 0.9) {
727 d = d1; e = e1;
728 }
729 if (eta < -0.9) {
730 d = d2; e = e2;
731 }
732
733 if (eta >= 0) signeta = 1;
734 if (eta < 0) signeta = -1;
735
736 return 1+b*pt+c*charge*pt*signeta*eta*eta+charge*d*pt*sin(phi+e);
737
738 }
739
740 // ------------ method called to for each event ------------
741 void FitNtuple::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
742 {
743
744 edm::ESHandle<Propagator> propagator;
745 iSetup.get<TrackingComponentsRecord>().get("SteppingHelixPropagatorAny", propagator);
746
747 edm::ESHandle<MagneticField> magneticField;
748 iSetup.get<IdealMagneticFieldRecord>().get(magneticField);
749
750 edm::ESHandle<GlobalTrackingGeometry> globalGeometry;
751 iSetup.get<GlobalTrackingGeometryRecord>().get(globalGeometry);
752
753 edm::ESHandle<TrackerGeometry> tkgeom;
754 iSetup.get<TrackerDigiGeometryRecord>().get(tkgeom);
755
756 edm::Handle<reco::BeamSpot> theBeamSpot;
757 iEvent.getByLabel("offlineBeamSpot",theBeamSpot);
758
759 // get the run and event number
760 m_run = iEvent.id().run();
761 m_lumi = iEvent.id().luminosityBlock();
762 m_event = iEvent.id().event();
763
764 // mu-jets (muons grouped by mass and vertex compatibility)
765 edm::Handle<pat::MultiMuonCollection> muJets;
766 iEvent.getByLabel(m_muJets, muJets);
767
768 // orphans (muons not found in any group)
769 edm::Handle<pat::MuonCollection> orphans;
770 iEvent.getByLabel(m_muJetOrphans, orphans);
771
772 edm::Handle<reco::TrackCollection> tracks;
773 iEvent.getByLabel("generalTracks", tracks);
774
775 edm::Handle<reco::PFCandidateCollection> pfCandidates;
776 iEvent.getByLabel("particleFlow", pfCandidates);
777
778 edm::Handle<pat::MuonCollection> allmuons;
779 iEvent.getByLabel(m_muons, allmuons);
780
781 // find the top four muons in the event (-1000. if not found)
782 edm::Handle<pat::MuonCollection> muons;
783 iEvent.getByLabel(m_muonsTriggerMatch, muons);
784 m_muon1pt = -1000.; m_muon1eta = -1000.;
785 m_muon2pt = -1000.; m_muon2eta = -1000.;
786 m_muon3pt = -1000.; m_muon3eta = -1000.;
787 m_muon4pt = -1000.; m_muon4eta = -1000.;
788 m_muontrigpt = -1000.; m_muontrigeta = -1000.;
789 std::vector<pat::MuonCollection::const_iterator> hightrigmuons;
790 std::vector<pat::MuonCollection::const_iterator> lowtrigmuons;
791
792 double iso_track_pt_treshold = 0.5;
793
794 /*
795 Cylinder::CylinderPointer m_cylinder100 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 100.);
796 Cylinder::CylinderPointer m_cylinder200 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 200.);
797 Cylinder::CylinderPointer m_cylinder300 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 300.);
798 Cylinder::CylinderPointer m_cylinder425 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 425.);
799 Cylinder::CylinderPointer m_cylinder470 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 470.);
800 Cylinder::CylinderPointer m_cylinder510 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 510.);
801 Cylinder::CylinderPointer m_cylinder565 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 565.);
802 Cylinder::CylinderPointer m_cylinder620 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 620.);
803 Cylinder::CylinderPointer m_cylinder670 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 670.);
804 Cylinder::CylinderPointer m_cylinder720 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 720.);
805 Cylinder::CylinderPointer m_cylinder800 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 800.);
806 Cylinder::CylinderPointer m_cylinder900 = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 900.);
807 */
808
809 for (pat::MuonCollection::const_iterator muon = muons->begin(); muon != muons->end(); ++muon) {
810 // std::cout << "MUON" << std::endl;
811 if ( fabs(muon->eta()) < 2.4
812 && muon->isTrackerMuon()
813 && muon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration) >= 2
814 && muon->innerTrack()->numberOfValidHits() >= 8
815 && muon->innerTrack()->normalizedChi2() < 4.) {
816 // std::cout << "SELECTED MUON" << std::endl;
817 if (muon->pt() > m_muon1pt) {
818 m_muon4pt = m_muon3pt; m_muon4eta = m_muon3eta;
819 m_muon3pt = m_muon2pt; m_muon3eta = m_muon2eta;
820 m_muon2pt = m_muon1pt; m_muon2eta = m_muon1eta;
821 m_muon1pt = muon->pt(); m_muon1eta = muon->eta();
822 } else if (muon->pt() > m_muon2pt) {
823 m_muon4pt = m_muon3pt; m_muon4eta = m_muon3eta;
824 m_muon3pt = m_muon2pt; m_muon3eta = m_muon2eta;
825 m_muon2pt = muon->pt(); m_muon2eta = muon->eta();
826 } else if (muon->pt() > m_muon3pt) {
827 m_muon4pt = m_muon3pt; m_muon4eta = m_muon3eta;
828 m_muon3pt = muon->pt(); m_muon3eta = muon->eta();
829 } else if (muon->pt() > m_muon4pt) {
830 m_muon4pt = muon->pt(); m_muon4eta = muon->eta();
831 }
832
833 // special muon within a more limited |eta| range, to guarantee the trigger
834 if (fabs(muon->eta()) < 2.4) {
835 if (muon->pt() > m_muontrigpt) {
836 m_muontrigpt = muon->pt(); m_muontrigeta = muon->eta();
837 }
838 }
839
840 if (m_dataset == "SingleMu") {
841 if (muon->pt() > m_trigpt && fabs(muon->eta()) < 2.4) {
842 const pat::TriggerObjectStandAlone *mu01 = muon->triggerObjectMatchByPath("HLT_Mu15_v2");
843 const pat::TriggerObjectStandAlone *mu02 = muon->triggerObjectMatchByPath("HLT_Mu20_v1");
844 const pat::TriggerObjectStandAlone *mu03 = muon->triggerObjectMatchByPath("HLT_Mu24_v1");
845 const pat::TriggerObjectStandAlone *mu04 = muon->triggerObjectMatchByPath("HLT_Mu24_v2");
846 const pat::TriggerObjectStandAlone *mu05 = muon->triggerObjectMatchByPath("HLT_Mu30_v1");
847 const pat::TriggerObjectStandAlone *mu06 = muon->triggerObjectMatchByPath("HLT_Mu30_v2");
848 const pat::TriggerObjectStandAlone *mu07 = muon->triggerObjectMatchByPath("HLT_Mu30_v3");
849 const pat::TriggerObjectStandAlone *mu08 = muon->triggerObjectMatchByPath("HLT_Mu40_v1");
850 const pat::TriggerObjectStandAlone *mu09 = muon->triggerObjectMatchByPath("HLT_Mu40_v2");
851 const pat::TriggerObjectStandAlone *mu10 = muon->triggerObjectMatchByPath("HLT_Mu40_v3");
852 const pat::TriggerObjectStandAlone *mu11 = muon->triggerObjectMatchByPath("HLT_Mu40_v5");
853 const pat::TriggerObjectStandAlone *mu12 = muon->triggerObjectMatchByPath("HLT_Mu40_eta2p1_v1");
854 const pat::TriggerObjectStandAlone *mu13 = muon->triggerObjectMatchByPath("HLT_Mu40_eta2p1_v4");
855 const pat::TriggerObjectStandAlone *mu14 = muon->triggerObjectMatchByPath("HLT_Mu40_eta2p1_v5");
856
857 if ((mu01 != NULL && mu01->collection() == std::string("hltL3MuonCandidates::HLT") && mu01->pt() > m_trigpt) ||
858 (mu02 != NULL && mu02->collection() == std::string("hltL3MuonCandidates::HLT") && mu02->pt() > m_trigpt) ||
859 (mu03 != NULL && mu03->collection() == std::string("hltL3MuonCandidates::HLT") && mu03->pt() > m_trigpt) ||
860 (mu04 != NULL && mu04->collection() == std::string("hltL3MuonCandidates::HLT") && mu04->pt() > m_trigpt) ||
861 (mu05 != NULL && mu05->collection() == std::string("hltL3MuonCandidates::HLT") && mu05->pt() > m_trigpt) ||
862 (mu06 != NULL && mu06->collection() == std::string("hltL3MuonCandidates::HLT") && mu06->pt() > m_trigpt) ||
863 (mu07 != NULL && mu07->collection() == std::string("hltL3MuonCandidates::HLT") && mu07->pt() > m_trigpt) ||
864 (mu08 != NULL && mu08->collection() == std::string("hltL3MuonCandidates::HLT") && mu08->pt() > m_trigpt) ||
865 (mu09 != NULL && mu09->collection() == std::string("hltL3MuonCandidates::HLT") && mu09->pt() > m_trigpt) ||
866 (mu10 != NULL && mu10->collection() == std::string("hltL3MuonCandidates::HLT") && mu10->pt() > m_trigpt) ||
867 (mu11 != NULL && mu11->collection() == std::string("hltL3MuonCandidates::HLT") && mu11->pt() > m_trigpt) ||
868 (mu12 != NULL && mu12->collection() == std::string("hltL3MuonCandidates::HLT") && mu12->pt() > m_trigpt) ||
869 (mu13 != NULL && mu13->collection() == std::string("hltL3MuonCandidates::HLT") && mu13->pt() > m_trigpt) ||
870 (mu14 != NULL && mu14->collection() == std::string("hltL3MuonCandidates::HLT") && mu14->pt() > m_trigpt))
871 hightrigmuons.push_back(muon);
872 }
873 }
874 if (m_dataset == "DoubleMu") {
875 if (muon->pt() > m_trigpt && fabs(muon->eta()) < 0.9) {
876 //if (muon->pt() > m_trigpt) {
877 const pat::TriggerObjectStandAlone *mu01 = muon->triggerObjectMatchByPath("HLT_DoubleMu6_v1");
878 const pat::TriggerObjectStandAlone *mu02 = muon->triggerObjectMatchByPath("HLT_DoubleMu7_v1");
879 const pat::TriggerObjectStandAlone *mu03 = muon->triggerObjectMatchByPath("HLT_DoubleMu7_v2");
880 const pat::TriggerObjectStandAlone *mu04 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v22");
881 const pat::TriggerObjectStandAlone *mu05 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v21");
882 const pat::TriggerObjectStandAlone *mu06 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v19");
883 const pat::TriggerObjectStandAlone *mu07 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v18");
884 const pat::TriggerObjectStandAlone *mu08 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v17");
885 const pat::TriggerObjectStandAlone *mu09 = muon->triggerObjectMatchByPath("HLT_Mu13_Mu8_v16");
886 const pat::TriggerObjectStandAlone *mu10 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v22");
887 const pat::TriggerObjectStandAlone *mu11 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v21");
888 const pat::TriggerObjectStandAlone *mu12 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v19");
889 const pat::TriggerObjectStandAlone *mu13 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v18");
890 const pat::TriggerObjectStandAlone *mu14 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v17");
891 const pat::TriggerObjectStandAlone *mu15 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v16");
892
893 if ((mu01 != NULL && mu01->collection() == std::string("hltL3MuonCandidates::HLT") && mu01->pt() > m_trigpt) ||
894 (mu02 != NULL && mu02->collection() == std::string("hltL3MuonCandidates::HLT") && mu02->pt() > m_trigpt) ||
895 (mu03 != NULL && mu03->collection() == std::string("hltL3MuonCandidates::HLT") && mu03->pt() > m_trigpt) ||
896 (mu04 != NULL && mu04->collection() == std::string("hltL3MuonCandidates::HLT") && mu04->pt() > m_trigpt) ||
897 (mu05 != NULL && mu05->collection() == std::string("hltL3MuonCandidates::HLT") && mu05->pt() > m_trigpt) ||
898 (mu06 != NULL && mu06->collection() == std::string("hltL3MuonCandidates::HLT") && mu06->pt() > m_trigpt) ||
899 (mu07 != NULL && mu07->collection() == std::string("hltL3MuonCandidates::HLT") && mu07->pt() > m_trigpt) ||
900 (mu08 != NULL && mu08->collection() == std::string("hltL3MuonCandidates::HLT") && mu08->pt() > m_trigpt) ||
901 (mu09 != NULL && mu09->collection() == std::string("hltL3MuonCandidates::HLT") && mu09->pt() > m_trigpt) ||
902 (mu10 != NULL && mu10->collection() == std::string("hltL3MuonCandidates::HLT") && mu10->pt() > m_trigpt) ||
903 (mu11 != NULL && mu11->collection() == std::string("hltL3MuonCandidates::HLT") && mu11->pt() > m_trigpt) ||
904 (mu12 != NULL && mu12->collection() == std::string("hltL3MuonCandidates::HLT") && mu12->pt() > m_trigpt) ||
905 (mu13 != NULL && mu13->collection() == std::string("hltL3MuonCandidates::HLT") && mu13->pt() > m_trigpt) ||
906 (mu14 != NULL && mu14->collection() == std::string("hltL3MuonCandidates::HLT") && mu14->pt() > m_trigpt) ||
907 (mu15 != NULL && mu15->collection() == std::string("hltL3MuonCandidates::HLT") && mu15->pt() > m_trigpt))
908 hightrigmuons.push_back(muon);
909 }
910 //if (muon->pt() > 8. && fabs(muon->eta()) < 0.9) {
911 if (muon->pt() > 8.) {
912 const pat::TriggerObjectStandAlone *mu01 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v22");
913 const pat::TriggerObjectStandAlone *mu02 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v21");
914 const pat::TriggerObjectStandAlone *mu03 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v19");
915 const pat::TriggerObjectStandAlone *mu04 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v18");
916 const pat::TriggerObjectStandAlone *mu05 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v17");
917 const pat::TriggerObjectStandAlone *mu06 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v16");
918 // const pat::TriggerObjectStandAlone *mu07 = muon->triggerObjectMatchByPath("HLT_Mu17_Mu8_v11");
919
920 if ((mu01 != NULL && mu01->collection() == std::string("hltL3MuonCandidates::HLT") && mu01->pt() > 8.) ||
921 (mu02 != NULL && mu02->collection() == std::string("hltL3MuonCandidates::HLT") && mu02->pt() > 8.) ||
922 (mu03 != NULL && mu03->collection() == std::string("hltL3MuonCandidates::HLT") && mu03->pt() > 8.) ||
923 (mu04 != NULL && mu04->collection() == std::string("hltL3MuonCandidates::HLT") && mu04->pt() > 8.) ||
924 (mu05 != NULL && mu05->collection() == std::string("hltL3MuonCandidates::HLT") && mu05->pt() > 8.) ||
925 (mu06 != NULL && mu06->collection() == std::string("hltL3MuonCandidates::HLT") && mu06->pt() > 8.)
926 // || (mu07 != NULL && mu07->collection() == std::string("hltL3MuonCandidates::HLT") && mu07->pt() > 8.)
927 )
928 lowtrigmuons.push_back(muon);
929 }
930 }
931 }
932 }
933
934 // // all tracker-tracks
935 // edm::Handle<reco::TrackCollection> tracks;
936 // iEvent.getByLabel("generalTracks", tracks);
937
938 // // generator-level 4-vectors
939 // edm::Handle<reco::GenParticleCollection> genParticles;
940 // iEvent.getByLabel("genParticles", genParticles);
941
942 // find the closest primary vertex (in Z) to the first muJet with a valid vertex
943 edm::Handle<reco::VertexCollection> primaryVertices;
944 iEvent.getByLabel("offlinePrimaryVertices", primaryVertices);
945 reco::VertexCollection::const_iterator closestPrimaryVertex = primaryVertices->end();
946 if (muJets->size() > 0) {
947 pat::MultiMuonCollection::const_iterator muJet0 = muJets->end();
948 for (pat::MultiMuonCollection::const_iterator muJet = muJets->begin(); muJet != muJets->end(); ++muJet) {
949 if (muJet->vertexValid()) {
950 muJet0 = muJet;
951 break;
952 }
953 }
954
955 if (muJet0 != muJets->end()) {
956 for (reco::VertexCollection::const_iterator vertex = primaryVertices->begin(); vertex != primaryVertices->end(); ++vertex) {
957 if (vertex->isValid() && !vertex->isFake() && vertex->tracksSize() > 3 && fabs(vertex->z()) < 24.) {
958 if (closestPrimaryVertex == primaryVertices->end() || fabs(vertex->z() - muJet0->vertexPoint().z()) < fabs(closestPrimaryVertex->z() - muJet0->vertexPoint().z())) {
959 closestPrimaryVertex = vertex;
960 }
961 } // end vertex quality cuts
962 } // end loop over primary vertices
963 } // end if muJet0 exists
964 } // end if muJets->size > 0
965
966 //m_primaryverticessize = primaryVertices->size();
967 m_primaryverticessize = 0;
968
969 for (reco::VertexCollection::const_iterator vertex = primaryVertices->begin(); vertex != primaryVertices->end(); ++vertex) {
970 if (vertex->isValid() && !vertex->isFake() && vertex->tracksSize() > 3 && fabs(vertex->z()) < 24.) {
971 m_primaryverticessize++;
972 }
973 }
974
975 m_generaltrackssize = tracks->size();
976 m_maxtrackspervtx = 0;
977
978 for (reco::VertexCollection::const_iterator vertex = primaryVertices->begin(); vertex != primaryVertices->end(); ++vertex) {
979 if (vertex->isValid() && fabs(vertex->z()) < 24.) {
980 if (int(vertex->tracksSize()) > m_maxtrackspervtx) m_maxtrackspervtx = int(vertex->tracksSize());
981 }
982 }
983
984 // find out which trigger bits were fired
985 edm::Handle<pat::TriggerEvent> triggerEvent;
986 iEvent.getByLabel("patTriggerEvent", triggerEvent);
987 m_trigger = 0;
988
989 if (m_dataset == "DoubleMu") {
990 if (triggerEvent->path("HLT_DoubleMu6_v1") && triggerEvent->path("HLT_DoubleMu6_v1")->wasAccept()) m_trigger += 1;
991 if (triggerEvent->path("HLT_DoubleMu7_v1") && triggerEvent->path("HLT_DoubleMu7_v1")->wasAccept()) m_trigger += 2;
992 if (triggerEvent->path("HLT_DoubleMu7_v2") && triggerEvent->path("HLT_DoubleMu7_v2")->wasAccept()) m_trigger += 2;
993 if (triggerEvent->path("HLT_Mu13_Mu8_v2") && triggerEvent->path("HLT_Mu13_Mu8_v2")->wasAccept()) m_trigger += 4;
994 if (triggerEvent->path("HLT_Mu13_Mu8_v3") && triggerEvent->path("HLT_Mu13_Mu8_v3")->wasAccept()) m_trigger += 4;
995 if (triggerEvent->path("HLT_Mu13_Mu8_v4") && triggerEvent->path("HLT_Mu13_Mu8_v4")->wasAccept()) m_trigger += 4;
996 if (triggerEvent->path("HLT_Mu13_Mu8_v6") && triggerEvent->path("HLT_Mu13_Mu8_v6")->wasAccept()) m_trigger += 4;
997 if (triggerEvent->path("HLT_Mu13_Mu8_v7") && triggerEvent->path("HLT_Mu13_Mu8_v7")->wasAccept()) m_trigger += 4;
998 if (triggerEvent->path("HLT_Mu17_Mu8_v2") && triggerEvent->path("HLT_Mu17_Mu8_v2")->wasAccept()) m_trigger += 8;
999 if (triggerEvent->path("HLT_Mu17_Mu8_v22") && triggerEvent->path("HLT_Mu17_Mu8_v22")->wasAccept()) m_trigger += 8;
1000 if (triggerEvent->path("HLT_Mu17_Mu8_v21") && triggerEvent->path("HLT_Mu17_Mu8_v21")->wasAccept()) m_trigger += 8;
1001 if (triggerEvent->path("HLT_Mu17_Mu8_v19") && triggerEvent->path("HLT_Mu17_Mu8_v19")->wasAccept()) m_trigger += 8;
1002 if (triggerEvent->path("HLT_Mu17_Mu8_v18") && triggerEvent->path("HLT_Mu17_Mu8_v18")->wasAccept()) m_trigger += 8;
1003 if (triggerEvent->path("HLT_Mu17_Mu8_v17") && triggerEvent->path("HLT_Mu17_Mu8_v17")->wasAccept()) m_trigger += 8;
1004 if (triggerEvent->path("HLT_Mu17_Mu8_v16") && triggerEvent->path("HLT_Mu17_Mu8_v16")->wasAccept()) m_trigger += 8;
1005 }
1006 if (m_dataset == "SingleMu") {
1007 if (triggerEvent->path("HLT_Mu15_v2") && triggerEvent->path("HLT_Mu15_v2")->wasAccept()) m_trigger += 1;
1008 if (triggerEvent->path("HLT_Mu20_v1") && triggerEvent->path("HLT_Mu20_v1")->wasAccept()) m_trigger += 2;
1009 if (triggerEvent->path("HLT_Mu24_v1") && triggerEvent->path("HLT_Mu24_v1")->wasAccept()) m_trigger += 4;
1010 if (triggerEvent->path("HLT_Mu24_v2") && triggerEvent->path("HLT_Mu24_v2")->wasAccept()) m_trigger += 4;
1011 if (triggerEvent->path("HLT_Mu30_v1") && triggerEvent->path("HLT_Mu30_v1")->wasAccept()) m_trigger += 8;
1012 if (triggerEvent->path("HLT_Mu30_v2") && triggerEvent->path("HLT_Mu30_v2")->wasAccept()) m_trigger += 8;
1013 if (triggerEvent->path("HLT_Mu30_v3") && triggerEvent->path("HLT_Mu30_v3")->wasAccept()) m_trigger += 8;
1014 if (triggerEvent->path("HLT_Mu40_v1") && triggerEvent->path("HLT_Mu40_v1")->wasAccept()) m_trigger += 16;
1015 if (triggerEvent->path("HLT_Mu40_v2") && triggerEvent->path("HLT_Mu40_v2")->wasAccept()) m_trigger += 16;
1016 if (triggerEvent->path("HLT_Mu40_v3") && triggerEvent->path("HLT_Mu40_v3")->wasAccept()) m_trigger += 16;
1017 if (triggerEvent->path("HLT_Mu40_v5") && triggerEvent->path("HLT_Mu40_v5")->wasAccept()) m_trigger += 16;
1018 if (triggerEvent->path("HLT_Mu40_eta2p1_v1") && triggerEvent->path("HLT_Mu40_eta2p1_v1")->wasAccept()) m_trigger += 16;
1019 if (triggerEvent->path("HLT_Mu40_eta2p1_v4") && triggerEvent->path("HLT_Mu40_eta2p1_v4")->wasAccept()) m_trigger += 16;
1020 if (triggerEvent->path("HLT_Mu40_eta2p1_v5") && triggerEvent->path("HLT_Mu40_eta2p1_v5")->wasAccept()) m_trigger += 16;
1021 }
1022
1023
1024
1025 m_trigger15 = 0;
1026 m_trigger20 = 0;
1027 m_trigger24 = 0;
1028 m_trigger30 = 0;
1029 m_trigger40 = 0;
1030 m_trigger40eta = 0;
1031
1032 if (triggerEvent->path("HLT_Mu15_v2") && triggerEvent->path("HLT_Mu15_v2")->wasAccept()) m_trigger15 = 1;
1033 if (triggerEvent->path("HLT_Mu20_v1") && triggerEvent->path("HLT_Mu20_v1")->wasAccept()) m_trigger20 = 1;
1034 if (triggerEvent->path("HLT_Mu24_v2") && triggerEvent->path("HLT_Mu24_v2")->wasAccept()) m_trigger24 = 1;
1035 if (triggerEvent->path("HLT_Mu30_v3") && triggerEvent->path("HLT_Mu30_v3")->wasAccept()) m_trigger30 = 1;
1036 if (triggerEvent->path("HLT_Mu40_v1") && triggerEvent->path("HLT_Mu40_v1")->wasAccept()) m_trigger40 = 1;
1037 if (triggerEvent->path("HLT_Mu40_v2") && triggerEvent->path("HLT_Mu40_v2")->wasAccept()) m_trigger40 = 1;
1038 if (triggerEvent->path("HLT_Mu40_v3") && triggerEvent->path("HLT_Mu40_v3")->wasAccept()) m_trigger40 = 1;
1039 if (triggerEvent->path("HLT_Mu40_v5") && triggerEvent->path("HLT_Mu40_v5")->wasAccept()) m_trigger40 = 1;
1040 if (triggerEvent->path("HLT_Mu40_eta2p1_v1") && triggerEvent->path("HLT_Mu40_eta2p1_v1")->wasAccept()) m_trigger40eta = 1;
1041 if (triggerEvent->path("HLT_Mu40_eta2p1_v4") && triggerEvent->path("HLT_Mu40_eta2p1_v4")->wasAccept()) m_trigger40eta = 1;
1042 if (triggerEvent->path("HLT_Mu40_eta2p1_v5") && triggerEvent->path("HLT_Mu40_eta2p1_v5")->wasAccept()) m_trigger40eta = 1;
1043
1044 ////////////////////////////////////////////////////////// lowdimuon and highdimuon:
1045
1046 //if (muJets->size() == 1 && (*muJets)[0].numberOfDaughters() == 2 && orphans->size() == 0 && m_trigger > 0) {
1047 if (muJets->size() == 1 && (*muJets)[0].numberOfDaughters() == 2 && orphans->size() == 0) {
1048 pat::MultiMuonCollection::const_iterator muJet = muJets->begin();
1049
1050 m_lowdimuon_trigpt1 = -20;
1051 m_lowdimuon_trigeta1 = -20;
1052 m_lowdimuon_trigpt2 = -20;
1053 m_lowdimuon_trigeta2 = -20;
1054
1055 m_lowdimuon_containstrig = 0;
1056 m_lowdimuon_containstrig2 = 0;
1057
1058 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = lowtrigmuons.begin(); iter != lowtrigmuons.end(); ++iter) {
1059 if (muJet->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1060 muJet->sameTrack(&*(muJet->muon(0)->innerTrack()), &*((*iter)->innerTrack()))) {
1061 m_lowdimuon_containstrig2++;
1062 }
1063 if (muJet->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1064 muJet->sameTrack(&*(muJet->muon(1)->innerTrack()), &*((*iter)->innerTrack()))) {
1065 m_lowdimuon_containstrig2++;
1066 }
1067 }
1068
1069 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = hightrigmuons.begin(); iter != hightrigmuons.end(); ++iter) {
1070 if (muJet->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1071 muJet->sameTrack(&*(muJet->muon(0)->innerTrack()), &*((*iter)->innerTrack()))) {
1072 m_lowdimuon_containstrig++;
1073 if (m_lowdimuon_containstrig2 == 2) {
1074 m_lowdimuon_trigpt1 = muJet->muon(0)->pt();
1075 m_lowdimuon_trigeta1 = muJet->muon(0)->eta();
1076 m_lowdimuon_trigpt2 = muJet->muon(1)->pt();
1077 m_lowdimuon_trigeta2 = muJet->muon(1)->eta();
1078 }
1079 }
1080 if (muJet->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1081 muJet->sameTrack(&*(muJet->muon(1)->innerTrack()), &*((*iter)->innerTrack()))) {
1082 m_lowdimuon_containstrig++;
1083 if (m_lowdimuon_containstrig2 == 2) {
1084 m_lowdimuon_trigpt1 = muJet->muon(1)->pt();
1085 m_lowdimuon_trigeta1 = muJet->muon(1)->eta();
1086 m_lowdimuon_trigpt2 = muJet->muon(0)->pt();
1087 m_lowdimuon_trigeta2 = muJet->muon(0)->eta();
1088 }
1089 }
1090 }
1091
1092 // generator-level mass using matched genParticles (for resolution of fit peak)
1093 m_lowdimuon_genmass = -1000.;
1094 if (muJet->muon(0)->genParticlesSize() == 1 && muJet->muon(1)->genParticlesSize() == 1) {
1095 const reco::GenParticle *mu0 = muJet->muon(0)->genParticle();
1096 const reco::GenParticle *mu1 = muJet->muon(1)->genParticle();
1097
1098 double total_energy = mu0->energy() + mu1->energy();
1099 double total_px = mu0->px() + mu1->px();
1100 double total_py = mu0->py() + mu1->py();
1101 double total_pz = mu0->pz() + mu1->pz();
1102 m_lowdimuon_genmass = sqrt(total_energy*total_energy - total_px*total_px - total_py*total_py - total_pz*total_pz);
1103 }
1104
1105 m_lowdimuon_mass = muJet->mass();
1106 m_lowdimuon_corr_mass = -1;
1107 m_lowdimuon_recomass = muJet->mass();
1108 m_lowdimuon_pt = muJet->pt();
1109 m_lowdimuon_eta = muJet->eta();
1110 m_lowdimuon_phi = muJet->phi();
1111 m_lowdimuon_dr = muJet->dRmax();
1112 if (muJet->daughter(0)->charge() > 0) {
1113 m_lowdimuon_pluspx = muJet->daughter(0)->px();
1114 m_lowdimuon_pluspy = muJet->daughter(0)->py();
1115 m_lowdimuon_pluspz = muJet->daughter(0)->pz();
1116 m_lowdimuon_minuspx = muJet->daughter(1)->px();
1117 m_lowdimuon_minuspy = muJet->daughter(1)->py();
1118 m_lowdimuon_minuspz = muJet->daughter(1)->pz();
1119 }
1120 else {
1121 m_lowdimuon_pluspx = muJet->daughter(1)->px();
1122 m_lowdimuon_pluspy = muJet->daughter(1)->py();
1123 m_lowdimuon_pluspz = muJet->daughter(1)->pz();
1124 m_lowdimuon_minuspx = muJet->daughter(0)->px();
1125 m_lowdimuon_minuspy = muJet->daughter(0)->py();
1126 m_lowdimuon_minuspz = muJet->daughter(0)->pz();
1127 }
1128 m_lowdimuon_vprob = -1000.;
1129 m_lowdimuon_vnchi2 = -1000.;
1130 m_lowdimuon_vx = -1000.;
1131 m_lowdimuon_vy = -1000.;
1132 m_lowdimuon_vz = -1000.;
1133
1134 m_lowdimuon_lxy = -1000.;
1135
1136 // replace all values with vertex-updated values if vertex fitting succeeded
1137 if (muJet->vertexValid()) {
1138 m_lowdimuon_mass = muJet->vertexMass();
1139 m_lowdimuon_corr_mass = muJet->vertexMass();
1140 m_lowdimuon_pt = muJet->vertexMomentum().perp();
1141 m_lowdimuon_eta = muJet->vertexMomentum().eta();
1142 m_lowdimuon_phi = muJet->vertexMomentum().phi();
1143 m_lowdimuon_dr = muJet->dRmax(true);
1144 if (muJet->daughter(0)->charge() > 0) {
1145 m_lowdimuon_pluspx = muJet->vertexMomentum(0).x();
1146 m_lowdimuon_pluspy = muJet->vertexMomentum(0).y();
1147 m_lowdimuon_pluspz = muJet->vertexMomentum(0).z();
1148 m_lowdimuon_minuspx = muJet->vertexMomentum(1).x();
1149 m_lowdimuon_minuspy = muJet->vertexMomentum(1).y();
1150 m_lowdimuon_minuspz = muJet->vertexMomentum(1).z();
1151 }
1152 else {
1153 m_lowdimuon_pluspx = muJet->vertexMomentum(1).x();
1154 m_lowdimuon_pluspy = muJet->vertexMomentum(1).y();
1155 m_lowdimuon_pluspz = muJet->vertexMomentum(1).z();
1156 m_lowdimuon_minuspx = muJet->vertexMomentum(0).x();
1157 m_lowdimuon_minuspy = muJet->vertexMomentum(0).y();
1158 m_lowdimuon_minuspz = muJet->vertexMomentum(0).z();
1159 }
1160 m_lowdimuon_vprob = muJet->vertexProb();
1161 m_lowdimuon_vnchi2 = muJet->vertexNormalizedChi2();
1162
1163 if (closestPrimaryVertex != primaryVertices->end()) {
1164 m_lowdimuon_vx = muJet->vertexPoint().x() - closestPrimaryVertex->x();
1165 m_lowdimuon_vy = muJet->vertexPoint().y() - closestPrimaryVertex->y();
1166 m_lowdimuon_vz = muJet->vertexPoint().z() - closestPrimaryVertex->z();
1167 m_lowdimuon_lxy = muJet->lxy(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
1168 }
1169
1170 double scale0 = scalePt(muJet->vertexMomentum(0).perp(),muJet->vertexMomentum(0).eta(),muJet->vertexMomentum(0).phi(),muJet->muon(0)->charge());
1171 double scale1 = scalePt(muJet->vertexMomentum(1).perp(),muJet->vertexMomentum(1).eta(),muJet->vertexMomentum(1).phi(),muJet->muon(1)->charge());
1172
1173 double e1 = pow(pow(muJet->vertexMomentum(0).perp()*scale0,2)+pow(muJet->vertexMomentum(0).z(),2)+pow(muJet->muon(0)->mass(),2),0.5);
1174 double e2 = pow(pow(muJet->vertexMomentum(1).perp()*scale1,2)+pow(muJet->vertexMomentum(1).z(),2)+pow(muJet->muon(1)->mass(),2),0.5);
1175 double e = e1 + e2;
1176
1177 double px = muJet->vertexMomentum(0).x()*scale0 + muJet->vertexMomentum(1).x()*scale1;
1178 double py = muJet->vertexMomentum(0).y()*scale0 + muJet->vertexMomentum(1).y()*scale1;
1179 double pz = muJet->vertexMomentum(0).z() + muJet->vertexMomentum(1).z();
1180
1181 m_lowdimuon_corr_mass = pow(e*e-px*px-py*py-pz*pz,0.5);
1182
1183 } // end of replacements with quantities measured at the vertex
1184
1185 m_lowdimuon_isoTk = muJet->centralTrackIsolation();
1186
1187 m_lowdimuon_isoTk_3mm = 0;
1188 m_lowdimuon_isoTk_2mm = 0;
1189 m_lowdimuon_isoTk_1mm = 0;
1190
1191 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
1192 bool track_is_muon = false;
1193 for (pat::MuonCollection::const_iterator muon = allmuons->begin(); muon != allmuons->end(); ++muon) {
1194 if (muJet->sameTrack(&*track,&*(muon->innerTrack()))) { track_is_muon = true; break; }
1195 }
1196 if (!track_is_muon) {
1197 double dphi = muJet->phi() - track->phi();
1198 if (dphi > M_PI) dphi -= 2.*M_PI;
1199 if (dphi < -M_PI) dphi += 2.*M_PI;
1200 double deta = muJet->eta() - track->eta();
1201 double dR = sqrt(pow(dphi, 2) + pow(deta, 2));
1202 if (dR < 0.4 && track->pt() > iso_track_pt_treshold) {
1203 double dz = fabs(track->dz(theBeamSpot->position())-muJet->dz(theBeamSpot->position()));
1204 if (dz < 0.3) m_lowdimuon_isoTk_3mm += track->pt();
1205 if (dz < 0.2) m_lowdimuon_isoTk_2mm += track->pt();
1206 if (dz < 0.1) m_lowdimuon_isoTk_1mm += track->pt();
1207 }
1208 }
1209 }
1210
1211 std::vector<reco::MuonChamberMatch> plusmatches, minusmatches;
1212 if (muJet->daughter(0)->charge() > 0) {
1213 m_lowdimuon_plusmatches = muJet->muon(0)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1214 m_lowdimuon_minusmatches = muJet->muon(1)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1215 plusmatches = muJet->muon(0)->matches();
1216 minusmatches = muJet->muon(1)->matches();
1217 m_lowdimuon_plushits = muJet->muon(0)->innerTrack()->numberOfValidHits();
1218 m_lowdimuon_minushits = muJet->muon(1)->innerTrack()->numberOfValidHits();
1219 m_lowdimuon_plusnormchi2 = muJet->muon(0)->innerTrack()->normalizedChi2();
1220 m_lowdimuon_minusnormchi2 = muJet->muon(1)->innerTrack()->normalizedChi2();
1221
1222 m_lowdimuon_plus_thits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1223 m_lowdimuon_plus_mhits = 0;
1224 if (muJet->muon(0)->isGlobalMuon()) m_lowdimuon_plus_mhits = muJet->muon(0)->globalTrack()->hitPattern().numberOfValidMuonHits();
1225 m_lowdimuon_plus_phits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidPixelHits();
1226
1227 m_lowdimuon_minus_thits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1228 m_lowdimuon_minus_mhits = 0;
1229 if (muJet->muon(1)->isGlobalMuon()) m_lowdimuon_minus_mhits = muJet->muon(1)->globalTrack()->hitPattern().numberOfValidMuonHits();
1230 m_lowdimuon_minus_phits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidPixelHits();
1231
1232
1233 m_lowdimuon_plus_dxy = muJet->muon(0)->innerTrack()->dxy(theBeamSpot->position());
1234 m_lowdimuon_minus_dxy = muJet->muon(1)->innerTrack()->dxy(theBeamSpot->position());
1235
1236 m_lowdimuon_plus_pt = muJet->muon(0)->pt();
1237 m_lowdimuon_minus_pt = muJet->muon(1)->pt();
1238 m_lowdimuon_plus_eta = muJet->muon(0)->eta();
1239 m_lowdimuon_minus_eta = muJet->muon(1)->eta();
1240 m_lowdimuon_plus_phi = muJet->muon(0)->phi();
1241 m_lowdimuon_minus_phi = muJet->muon(1)->phi();
1242
1243 m_lowdimuon_Deltaphi = m_lowdimuon_plus_phi-m_lowdimuon_minus_phi;
1244 if (m_lowdimuon_Deltaphi > M_PI) m_lowdimuon_Deltaphi -= 2.*M_PI;
1245 if (m_lowdimuon_Deltaphi < -M_PI) m_lowdimuon_Deltaphi += 2.*M_PI;
1246
1247
1248 /*
1249 for (trackingRecHit_iterator hit = muJet->muon(0)->innerTrack()->recHitsBegin(); hit != muJet->muon(0)->innerTrack()->recHitsEnd(); ++hit) {
1250 if ((*hit)->isValid()) {
1251 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1252 m_lowdimuon_plus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1253 break;
1254 }
1255 }
1256 for (trackingRecHit_iterator hit = muJet->muon(1)->innerTrack()->recHitsBegin(); hit != muJet->muon(1)->innerTrack()->recHitsEnd(); ++hit) {
1257 if ((*hit)->isValid()) {
1258 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1259 m_lowdimuon_minus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1260 break;
1261 }
1262 }
1263 */
1264
1265 }
1266 else {
1267 m_lowdimuon_plusmatches = muJet->muon(1)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1268 m_lowdimuon_minusmatches = muJet->muon(0)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1269 plusmatches = muJet->muon(1)->matches();
1270 minusmatches = muJet->muon(0)->matches();
1271 m_lowdimuon_plushits = muJet->muon(1)->innerTrack()->numberOfValidHits();
1272 m_lowdimuon_minushits = muJet->muon(0)->innerTrack()->numberOfValidHits();
1273 m_lowdimuon_plusnormchi2 = muJet->muon(1)->innerTrack()->normalizedChi2();
1274 m_lowdimuon_minusnormchi2 = muJet->muon(0)->innerTrack()->normalizedChi2();
1275
1276 m_lowdimuon_plus_thits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1277 m_lowdimuon_plus_mhits = 0;
1278 if (muJet->muon(1)->isGlobalMuon()) m_lowdimuon_plus_mhits = muJet->muon(1)->globalTrack()->hitPattern().numberOfValidMuonHits();
1279 m_lowdimuon_plus_phits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidPixelHits();
1280
1281 m_lowdimuon_minus_thits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1282 m_lowdimuon_minus_mhits = 0;
1283 if (muJet->muon(0)->isGlobalMuon()) m_lowdimuon_minus_mhits = muJet->muon(0)->globalTrack()->hitPattern().numberOfValidMuonHits();
1284 m_lowdimuon_minus_phits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidPixelHits();
1285
1286 m_lowdimuon_plus_dxy = muJet->muon(1)->innerTrack()->dxy(theBeamSpot->position());
1287 m_lowdimuon_minus_dxy = muJet->muon(0)->innerTrack()->dxy(theBeamSpot->position());
1288
1289 m_lowdimuon_plus_pt = muJet->muon(1)->pt();
1290 m_lowdimuon_minus_pt = muJet->muon(0)->pt();
1291 m_lowdimuon_plus_eta = muJet->muon(1)->eta();
1292 m_lowdimuon_minus_eta = muJet->muon(0)->eta();
1293 m_lowdimuon_plus_phi = muJet->muon(1)->phi();
1294 m_lowdimuon_minus_phi = muJet->muon(0)->phi();
1295
1296 m_lowdimuon_Deltaphi = m_lowdimuon_plus_phi-m_lowdimuon_minus_phi;
1297 if (m_lowdimuon_Deltaphi > M_PI) m_lowdimuon_Deltaphi -= 2.*M_PI;
1298 if (m_lowdimuon_Deltaphi < -M_PI) m_lowdimuon_Deltaphi += 2.*M_PI;
1299 }
1300
1301 /*
1302 const pat::Muon *muplus = NULL;
1303 const pat::Muon *muminus = NULL;
1304
1305 if (muJet->muon(0)->charge() > 0) {
1306 muplus = &*muJet->muon(0);
1307 muminus = &*muJet->muon(1);
1308 }
1309 if (muJet->muon(0)->charge() < 0) {
1310 muplus = &*muJet->muon(1);
1311 muminus = &*muJet->muon(0);
1312 }
1313
1314 FreeTrajectoryState plus_initial(GlobalPoint(muplus->vx(), muplus->vy(), muplus->vz()), GlobalVector(muplus->px(), muplus->py(), muplus->pz()), 1, &*magneticField);
1315 FreeTrajectoryState minus_initial(GlobalPoint(muminus->vx(), muminus->vy(), muminus->vz()), GlobalVector(muminus->px(), muminus->py(), muminus->pz()), -1, &*magneticField);
1316
1317 if (fabs(m_lowdimuon_eta) < 0.9) {
1318 TrajectoryStateOnSurface plus_final = propagator->propagate(plus_initial, *m_cylinder100);
1319 TrajectoryStateOnSurface minus_final = propagator->propagate(minus_initial, *m_cylinder100);
1320 if (plus_final.isValid() && minus_final.isValid()) {
1321 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1322 while (dphi > M_PI) dphi -= 2.*M_PI;
1323 while (dphi < -M_PI) dphi += 2.*M_PI;
1324 m_lowdimuon_dphi100 = dphi;
1325 }
1326 plus_final = propagator->propagate(plus_initial, *m_cylinder200);
1327 minus_final = propagator->propagate(minus_initial, *m_cylinder200);
1328 if (plus_final.isValid() && minus_final.isValid()) {
1329 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1330 while (dphi > M_PI) dphi -= 2.*M_PI;
1331 while (dphi < -M_PI) dphi += 2.*M_PI;
1332 m_lowdimuon_dphi200 = dphi;
1333 }
1334 plus_final = propagator->propagate(plus_initial, *m_cylinder300);
1335 minus_final = propagator->propagate(minus_initial, *m_cylinder300);
1336 if (plus_final.isValid() && minus_final.isValid()) {
1337 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1338 while (dphi > M_PI) dphi -= 2.*M_PI;
1339 while (dphi < -M_PI) dphi += 2.*M_PI;
1340 m_lowdimuon_dphi300 = dphi;
1341 }
1342 plus_final = propagator->propagate(plus_initial, *m_cylinder425);
1343 minus_final = propagator->propagate(minus_initial, *m_cylinder425);
1344 if (plus_final.isValid() && minus_final.isValid()) {
1345 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1346 while (dphi > M_PI) dphi -= 2.*M_PI;
1347 while (dphi < -M_PI) dphi += 2.*M_PI;
1348 m_lowdimuon_dphi425 = dphi;
1349 }
1350 plus_final = propagator->propagate(plus_initial, *m_cylinder470);
1351 minus_final = propagator->propagate(minus_initial, *m_cylinder470);
1352 if (plus_final.isValid() && minus_final.isValid()) {
1353 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1354 while (dphi > M_PI) dphi -= 2.*M_PI;
1355 while (dphi < -M_PI) dphi += 2.*M_PI;
1356 m_lowdimuon_dphi470 = dphi;
1357 }
1358 plus_final = propagator->propagate(plus_initial, *m_cylinder510);
1359 minus_final = propagator->propagate(minus_initial, *m_cylinder510);
1360 if (plus_final.isValid() && minus_final.isValid()) {
1361 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1362 while (dphi > M_PI) dphi -= 2.*M_PI;
1363 while (dphi < -M_PI) dphi += 2.*M_PI;
1364 m_lowdimuon_dphi510 = dphi;
1365 }
1366 plus_final = propagator->propagate(plus_initial, *m_cylinder565);
1367 minus_final = propagator->propagate(minus_initial, *m_cylinder565);
1368 if (plus_final.isValid() && minus_final.isValid()) {
1369 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1370 while (dphi > M_PI) dphi -= 2.*M_PI;
1371 while (dphi < -M_PI) dphi += 2.*M_PI;
1372 m_lowdimuon_dphi565 = dphi;
1373 }
1374 plus_final = propagator->propagate(plus_initial, *m_cylinder620);
1375 minus_final = propagator->propagate(minus_initial, *m_cylinder620);
1376 if (plus_final.isValid() && minus_final.isValid()) {
1377 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1378 while (dphi > M_PI) dphi -= 2.*M_PI;
1379 while (dphi < -M_PI) dphi += 2.*M_PI;
1380 m_lowdimuon_dphi620 = dphi;
1381 }
1382 plus_final = propagator->propagate(plus_initial, *m_cylinder670);
1383 minus_final = propagator->propagate(minus_initial, *m_cylinder670);
1384 if (plus_final.isValid() && minus_final.isValid()) {
1385 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1386 while (dphi > M_PI) dphi -= 2.*M_PI;
1387 while (dphi < -M_PI) dphi += 2.*M_PI;
1388 m_lowdimuon_dphi670 = dphi;
1389 }
1390 plus_final = propagator->propagate(plus_initial, *m_cylinder720);
1391 minus_final = propagator->propagate(minus_initial, *m_cylinder720);
1392 if (plus_final.isValid() && minus_final.isValid()) {
1393 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1394 while (dphi > M_PI) dphi -= 2.*M_PI;
1395 while (dphi < -M_PI) dphi += 2.*M_PI;
1396 m_lowdimuon_dphi720 = dphi;
1397 }
1398 plus_final = propagator->propagate(plus_initial, *m_cylinder800);
1399 minus_final = propagator->propagate(minus_initial, *m_cylinder800);
1400 if (plus_final.isValid() && minus_final.isValid()) {
1401 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1402 while (dphi > M_PI) dphi -= 2.*M_PI;
1403 while (dphi < -M_PI) dphi += 2.*M_PI;
1404 m_lowdimuon_dphi800 = dphi;
1405 }
1406 plus_final = propagator->propagate(plus_initial, *m_cylinder900);
1407 minus_final = propagator->propagate(minus_initial, *m_cylinder900);
1408 if (plus_final.isValid() && minus_final.isValid()) {
1409 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1410 while (dphi > M_PI) dphi -= 2.*M_PI;
1411 while (dphi < -M_PI) dphi += 2.*M_PI;
1412 m_lowdimuon_dphi900 = dphi;
1413 }
1414 }
1415 */
1416
1417 m_lowdimuon_bbbarlike = -1;
1418 if (m_lowdimuon_isoTk > 4.5 || m_lowdimuon_lxy > 0.2) m_lowdimuon_bbbarlike = 1;
1419 if (m_lowdimuon_isoTk < 4.5 && m_lowdimuon_lxy < 0.2) m_lowdimuon_bbbarlike = 0;
1420
1421 m_lowdimuon_dz = muJet->dz(theBeamSpot->position());
1422
1423
1424
1425 if (m_lowdimuon_containstrig > 0 && m_lowdimuon_containstrig2 > 1) {
1426 m_lowdimuon->Fill();
1427 }
1428 }
1429
1430 ////////////////////////////////////////////////////////// dimuorphan
1431
1432 if (muJets->size() == 1 && (*muJets)[0].numberOfDaughters() == 2 && orphans->size() == 1 && m_trigger > 0) {
1433 //if (muJets->size() == 1 && (*muJets)[0].numberOfDaughters() == 2 && orphans->size() == 1) {
1434 pat::MultiMuonCollection::const_iterator muJet = muJets->begin();
1435 pat::MuonCollection::const_iterator orphan = orphans->begin();
1436
1437 m_dimuorphan_deltaphi = muJet->phi() - orphan->phi();
1438 while (m_dimuorphan_deltaphi > M_PI) m_dimuorphan_deltaphi -= 2.*M_PI;
1439 while (m_dimuorphan_deltaphi < -M_PI) m_dimuorphan_deltaphi += 2.*M_PI;
1440
1441 double deta = orphan->innerTrack()->eta() - muJet->muon(0)->innerTrack()->eta();
1442 double dphi = orphan->innerTrack()->phi() - muJet->muon(0)->innerTrack()->phi();
1443 if (dphi > M_PI) dphi -= 2.*M_PI;
1444 if (dphi < -M_PI) dphi += 2.*M_PI;
1445 m_dimuorphan_dr1 = pow(dphi*dphi+deta*deta,0.5);
1446
1447 deta = orphan->innerTrack()->eta() - muJet->muon(1)->innerTrack()->eta();
1448 dphi = orphan->innerTrack()->phi() - muJet->muon(1)->innerTrack()->phi();
1449 if (dphi > M_PI) dphi -= 2.*M_PI;
1450 if (dphi < -M_PI) dphi += 2.*M_PI;
1451 m_dimuorphan_dr2 = pow(dphi*dphi+deta*deta,0.5);
1452
1453 m_dimuorphan_orphanpt = orphan->pt();
1454 m_dimuorphan_orphaneta = orphan->eta();
1455 m_dimuorphan_orphanphi = orphan->phi();
1456 m_dimuorphan_orphanisglobal = (orphan->isGlobalMuon() ? 1 : 0);
1457 m_dimuorphan_orphanmatches = orphan->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1458 m_dimuorphan_stationmask = orphan->stationMask();
1459 m_dimuorphan_orphanhits = (orphan->innerTrack().isAvailable() ? orphan->innerTrack()->numberOfValidHits(): -1);
1460 m_dimuorphan_orphanchi2 = (orphan->innerTrack().isAvailable() ? orphan->innerTrack()->normalizedChi2(): -1.);
1461
1462 m_dimuorphan_trigpt1 = -20;
1463 m_dimuorphan_trigeta1 = -20;
1464 m_dimuorphan_trigpt2 = -20;
1465 m_dimuorphan_trigeta2 = -20;
1466
1467 m_dimuorphan_containstrig = 0;
1468 m_dimuorphan_containstrig2 = 0;
1469
1470 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = hightrigmuons.begin(); iter != hightrigmuons.end(); ++iter) {
1471 if (orphan->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1472 muJet->sameTrack(&*(orphan->innerTrack()), &*((*iter)->innerTrack()))) {
1473 m_dimuorphan_containstrig++;
1474 }
1475 }
1476
1477 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = hightrigmuons.begin(); iter != hightrigmuons.end(); ++iter) {
1478 if (muJet->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1479 muJet->sameTrack(&*(muJet->muon(0)->innerTrack()), &*((*iter)->innerTrack()))) {
1480 m_dimuorphan_containstrig2++;
1481 if (m_dimuorphan_containstrig2 == 2) {
1482 m_dimuorphan_trigpt1 = muJet->muon(0)->pt();
1483 m_dimuorphan_trigeta1 = muJet->muon(0)->eta();
1484 m_dimuorphan_trigpt2 = orphan->pt();
1485 m_dimuorphan_trigeta2 = orphan->eta();
1486 }
1487 }
1488 if (muJet->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
1489 muJet->sameTrack(&*(muJet->muon(1)->innerTrack()), &*((*iter)->innerTrack()))) {
1490 m_dimuorphan_containstrig2++;
1491 if (m_dimuorphan_containstrig2 == 2) {
1492 m_dimuorphan_trigpt1 = muJet->muon(1)->pt();
1493 m_dimuorphan_trigeta1 = muJet->muon(1)->eta();
1494 m_dimuorphan_trigpt2 = orphan->pt();
1495 m_dimuorphan_trigeta2 = orphan->eta();
1496 }
1497 }
1498 }
1499
1500
1501
1502 m_dimuorphan_mass = muJet->mass();
1503 m_dimuorphan_corr_mass = -1;
1504 m_dimuorphan_recomass = muJet->mass();
1505 m_dimuorphan_pt = muJet->pt();
1506 m_dimuorphan_eta = muJet->eta();
1507 m_dimuorphan_phi = muJet->phi();
1508 m_dimuorphan_dr = muJet->dRmax();
1509 m_dimuorphan_vprob = -1.;
1510 m_dimuorphan_vnchi2 = -1.;
1511 m_dimuorphan_lxy = -1000.;
1512 m_dimuorphan_lxyz = -1000.;
1513 m_dimuorphan_caloiso = muJet->centralCaloIsolation();
1514 m_dimuorphan_isoTk = muJet->centralTrackIsolation();
1515
1516 m_dimuorphan_isoTk_3mm = 0;
1517 m_dimuorphan_isoTk_2mm = 0;
1518 m_dimuorphan_isoTk_1mm = 0;
1519
1520 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
1521 bool track_is_muon = false;
1522 if (muJet->sameTrack(&*track,&*(muJet->muon(0)->innerTrack())) || muJet->sameTrack(&*track,&*(muJet->muon(1)->innerTrack()))) track_is_muon = true;
1523 if (!track_is_muon) {
1524 double dphi = muJet->phi() - track->phi();
1525 if (dphi > M_PI) dphi -= 2.*M_PI;
1526 if (dphi < -M_PI) dphi += 2.*M_PI;
1527 double deta = muJet->eta() - track->eta();
1528 double dR = sqrt(pow(dphi, 2) + pow(deta, 2));
1529 if (dR < 0.4 && track->pt() > iso_track_pt_treshold) {
1530 double dz = fabs(track->dz(theBeamSpot->position())-muJet->dz(theBeamSpot->position()));
1531 if (dz < 0.3) m_dimuorphan_isoTk_3mm += track->pt();
1532 if (dz < 0.2) m_dimuorphan_isoTk_2mm += track->pt();
1533 if (dz < 0.1) m_dimuorphan_isoTk_1mm += track->pt();
1534 }
1535 }
1536 }
1537
1538 m_dimuorphan_isoPF_1mm = 0.0;
1539 for (reco::PFCandidateCollection::const_iterator pfCand = pfCandidates->begin(); pfCand != pfCandidates->end(); ++pfCand) {
1540 if ( pfCand->particleId() == reco::PFCandidate::h ) {
1541 double dPhi = muJet->phi() - pfCand->phi();
1542 if (dPhi > M_PI) dPhi -= 2.*M_PI;
1543 if (dPhi < -M_PI) dPhi += 2.*M_PI;
1544 double dEta = muJet->eta() - pfCand->eta();
1545 double dR = sqrt( dPhi*dPhi + dEta*dEta );
1546 // if ( dR < 0.4 && pfCand->pt() > 0.5 ) {
1547 if ( dR < 0.4 ) {
1548 double dz = fabs( pfCand->vertex().z() - theBeamSpot->position().z() - muJet->dz(theBeamSpot->position()) );
1549 if ( dz < 0.1 ) m_dimuorphan_isoPF_1mm += pfCand->pt();
1550 }
1551 }
1552 }
1553
1554 m_dimuorphan_orphan_isoTk_3mm = 0;
1555 m_dimuorphan_orphan_isoTk_2mm = 0;
1556 m_dimuorphan_orphan_isoTk_1mm = 0;
1557
1558 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
1559 if (!muJet->sameTrack(&*track,&*(orphan->innerTrack()))) {
1560 double dphi = orphan->innerTrack()->phi() - track->phi();
1561 if (dphi > M_PI) dphi -= 2.*M_PI;
1562 if (dphi < -M_PI) dphi += 2.*M_PI;
1563 double deta = orphan->innerTrack()->eta() - track->eta();
1564 double dR = sqrt(pow(dphi, 2) + pow(deta, 2));
1565 if (dR < 0.4 && track->pt() > iso_track_pt_treshold) {
1566 double dz = fabs(track->dz(theBeamSpot->position())-orphan->innerTrack()->dz(theBeamSpot->position()));
1567 if (dz < 0.3) m_dimuorphan_orphan_isoTk_3mm += track->pt();
1568 if (dz < 0.2) m_dimuorphan_orphan_isoTk_2mm += track->pt();
1569 if (dz < 0.1) m_dimuorphan_orphan_isoTk_1mm += track->pt();
1570 }
1571 }
1572 }
1573
1574 m_dimuorphan_orphan_isoPF_1mm = 0.0;
1575 for (reco::PFCandidateCollection::const_iterator pfCand = pfCandidates->begin(); pfCand != pfCandidates->end(); ++pfCand) {
1576 if ( pfCand->particleId() == reco::PFCandidate::h ) {
1577 double dPhi = orphan->innerTrack()->phi() - pfCand->phi();
1578 if (dPhi > M_PI) dPhi -= 2.*M_PI;
1579 if (dPhi < -M_PI) dPhi += 2.*M_PI;
1580 double dEta = orphan->innerTrack()->eta() - pfCand->eta();
1581 double dR = sqrt( dPhi*dPhi + dEta*dEta );
1582 // if ( dR < 0.4 && pfCand->pt() > 0.5 ) {
1583 if ( dR < 0.4 ) {
1584 double dz = fabs( pfCand->vertex().z() - theBeamSpot->position().z() - orphan->innerTrack()->dz(theBeamSpot->position()) );
1585 if ( dz < 0.1 ) m_dimuorphan_orphan_isoPF_1mm += pfCand->pt();
1586 }
1587 }
1588 }
1589
1590 m_dimuorphan_trackdensity = 0;
1591
1592 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
1593 if (!muJet->sameTrack(&*track,&*(muJet->muon(0)->innerTrack())) && !muJet->sameTrack(&*track,&*(muJet->muon(1)->innerTrack()))) {
1594 double dphi0 = muJet->muon(0)->innerTrack()->phi() - track->phi();
1595 if (dphi0 > M_PI) dphi0 -= 2.*M_PI;
1596 if (dphi0 < -M_PI) dphi0 += 2.*M_PI;
1597 double deta0 = muJet->muon(0)->innerTrack()->eta() - track->eta();
1598 double dR0 = sqrt(pow(dphi0, 2) + pow(deta0, 2));
1599 double dphi1 = muJet->muon(1)->innerTrack()->phi() - track->phi();
1600 if (dphi1 > M_PI) dphi1 -= 2.*M_PI;
1601 if (dphi1 < -M_PI) dphi1 += 2.*M_PI;
1602 double deta1 = muJet->muon(1)->innerTrack()->eta() - track->eta();
1603 double dR1 = sqrt(pow(dphi1, 2) + pow(deta1, 2));
1604 if ((dR0 < 0.2 || dR1 < 0.2) && track->pt() > 1.) {
1605 double dz0 = fabs(track->dz(theBeamSpot->position())-muJet->muon(0)->innerTrack()->dz(theBeamSpot->position()));
1606 double dz1 = fabs(track->dz(theBeamSpot->position())-muJet->muon(1)->innerTrack()->dz(theBeamSpot->position()));
1607 if (dz0 < 0.2 || dz1 < 0.2) m_dimuorphan_trackdensity++;
1608 }
1609 }
1610 }
1611
1612 if (muJet->vertexValid()) {
1613 m_dimuorphan_deltaphi = muJet->vertexMomentum().phi() - orphan->phi();
1614 while (m_dimuorphan_deltaphi > M_PI) m_dimuorphan_deltaphi -= 2.*M_PI;
1615 while (m_dimuorphan_deltaphi < -M_PI) m_dimuorphan_deltaphi += 2.*M_PI;
1616
1617 m_dimuorphan_mass = muJet->vertexMass();
1618 m_dimuorphan_corr_mass = muJet->vertexMass();
1619 m_dimuorphan_pt = muJet->vertexMomentum().perp();
1620 m_dimuorphan_eta = muJet->vertexMomentum().eta();
1621 m_dimuorphan_phi = muJet->vertexMomentum().phi();
1622 m_dimuorphan_dr = muJet->dRmax(true);
1623 m_dimuorphan_vprob = muJet->vertexProb();
1624 m_dimuorphan_vnchi2 = muJet->vertexNormalizedChi2();
1625 if (closestPrimaryVertex != primaryVertices->end()) {
1626 m_dimuorphan_lxy = muJet->lxy(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
1627 m_dimuorphan_lxyz = muJet->lxyz(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
1628 }
1629
1630 double scale0 = scalePt(muJet->vertexMomentum(0).perp(),muJet->vertexMomentum(0).eta(),muJet->vertexMomentum(0).phi(),muJet->muon(0)->charge());
1631 double scale1 = scalePt(muJet->vertexMomentum(1).perp(),muJet->vertexMomentum(1).eta(),muJet->vertexMomentum(1).phi(),muJet->muon(1)->charge());
1632
1633 double e1 = pow(pow(muJet->vertexMomentum(0).perp()*scale0,2)+pow(muJet->vertexMomentum(0).z(),2)+pow(muJet->muon(0)->mass(),2),0.5);
1634 double e2 = pow(pow(muJet->vertexMomentum(1).perp()*scale1,2)+pow(muJet->vertexMomentum(1).z(),2)+pow(muJet->muon(1)->mass(),2),0.5);
1635 double e = e1 + e2;
1636
1637 double px = muJet->vertexMomentum(0).x()*scale0 + muJet->vertexMomentum(1).x()*scale1;
1638 double py = muJet->vertexMomentum(0).y()*scale0 + muJet->vertexMomentum(1).y()*scale1;
1639 double pz = muJet->vertexMomentum(0).z() + muJet->vertexMomentum(1).z();
1640
1641 m_dimuorphan_corr_mass = pow(e*e-px*px-py*py-pz*pz,0.5);
1642
1643 }
1644
1645 if (muJet->daughter(0)->charge() > 0) {
1646 m_dimuorphan_plusmatches = muJet->muon(0)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1647 m_dimuorphan_minusmatches = muJet->muon(1)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1648 m_dimuorphan_plushits = muJet->muon(0)->innerTrack()->numberOfValidHits();
1649 m_dimuorphan_minushits = muJet->muon(1)->innerTrack()->numberOfValidHits();
1650 m_dimuorphan_plusnormchi2 = muJet->muon(0)->innerTrack()->normalizedChi2();
1651 m_dimuorphan_minusnormchi2 = muJet->muon(1)->innerTrack()->normalizedChi2();
1652
1653 m_dimuorphan_plusstationmask = muJet->muon(0)->stationMask();
1654 m_dimuorphan_minusstationmask = muJet->muon(1)->stationMask();
1655
1656 m_dimuorphan_plus_thits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1657 m_dimuorphan_plus_mhits = 0;
1658 if (muJet->muon(0)->isGlobalMuon()) m_dimuorphan_plus_mhits = muJet->muon(0)->globalTrack()->hitPattern().numberOfValidMuonHits();
1659 m_dimuorphan_plus_phits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidPixelHits();
1660
1661 m_dimuorphan_minus_thits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1662 m_dimuorphan_minus_mhits = 0;
1663 if (muJet->muon(1)->isGlobalMuon()) m_dimuorphan_minus_mhits = muJet->muon(1)->globalTrack()->hitPattern().numberOfValidMuonHits();
1664 m_dimuorphan_minus_phits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidPixelHits();
1665
1666 m_dimuorphan_plus_dxy = muJet->muon(0)->innerTrack()->dxy(theBeamSpot->position());
1667 m_dimuorphan_minus_dxy = muJet->muon(1)->innerTrack()->dxy(theBeamSpot->position());
1668
1669 m_dimuorphan_plus_pt = muJet->muon(0)->pt();
1670 m_dimuorphan_minus_pt = muJet->muon(1)->pt();
1671 m_dimuorphan_plus_eta = muJet->muon(0)->eta();
1672 m_dimuorphan_minus_eta = muJet->muon(1)->eta();
1673 m_dimuorphan_plus_phi = muJet->muon(0)->phi();
1674 m_dimuorphan_minus_phi = muJet->muon(1)->phi();
1675
1676 m_dimuorphan_Deltaphi = m_dimuorphan_plus_phi-m_dimuorphan_minus_phi;
1677 if (m_dimuorphan_Deltaphi > M_PI) m_dimuorphan_Deltaphi -= 2.*M_PI;
1678 if (m_dimuorphan_Deltaphi < -M_PI) m_dimuorphan_Deltaphi += 2.*M_PI;
1679 /*
1680 for (trackingRecHit_iterator hit = muJet->muon(0)->innerTrack()->recHitsBegin(); hit != muJet->muon(0)->innerTrack()->recHitsEnd(); ++hit) {
1681 if ((*hit)->isValid()) {
1682 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1683 m_dimuorphan_plus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1684 break;
1685 }
1686 }
1687 for (trackingRecHit_iterator hit = muJet->muon(1)->innerTrack()->recHitsBegin(); hit != muJet->muon(1)->innerTrack()->recHitsEnd(); ++hit) {
1688 if ((*hit)->isValid()) {
1689 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1690 m_dimuorphan_minus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1691 break;
1692 }
1693 }
1694 */
1695
1696 }
1697 else {
1698 m_dimuorphan_plusmatches = muJet->muon(1)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1699 m_dimuorphan_minusmatches = muJet->muon(0)->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
1700 m_dimuorphan_plushits = muJet->muon(1)->innerTrack()->numberOfValidHits();
1701 m_dimuorphan_minushits = muJet->muon(0)->innerTrack()->numberOfValidHits();
1702 m_dimuorphan_plusnormchi2 = muJet->muon(1)->innerTrack()->normalizedChi2();
1703 m_dimuorphan_minusnormchi2 = muJet->muon(0)->innerTrack()->normalizedChi2();
1704
1705 m_dimuorphan_plusstationmask = muJet->muon(1)->stationMask();
1706 m_dimuorphan_minusstationmask = muJet->muon(0)->stationMask();
1707
1708 m_dimuorphan_plus_thits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1709 m_dimuorphan_plus_mhits = 0;
1710 if (muJet->muon(1)->isGlobalMuon()) m_dimuorphan_plus_mhits = muJet->muon(1)->globalTrack()->hitPattern().numberOfValidMuonHits();
1711 m_dimuorphan_plus_phits = muJet->muon(1)->innerTrack()->hitPattern().numberOfValidPixelHits();
1712
1713 m_dimuorphan_minus_thits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidTrackerHits();
1714 m_dimuorphan_minus_mhits = 0;
1715 if (muJet->muon(0)->isGlobalMuon()) m_dimuorphan_minus_mhits = muJet->muon(0)->globalTrack()->hitPattern().numberOfValidMuonHits();
1716 m_dimuorphan_minus_phits = muJet->muon(0)->innerTrack()->hitPattern().numberOfValidPixelHits();
1717
1718 m_dimuorphan_plus_dxy = muJet->muon(1)->innerTrack()->dxy(theBeamSpot->position());
1719 m_dimuorphan_minus_dxy = muJet->muon(0)->innerTrack()->dxy(theBeamSpot->position());
1720
1721 m_dimuorphan_plus_pt = muJet->muon(1)->pt();
1722 m_dimuorphan_minus_pt = muJet->muon(0)->pt();
1723 m_dimuorphan_plus_eta = muJet->muon(1)->eta();
1724 m_dimuorphan_minus_eta = muJet->muon(0)->eta();
1725 m_dimuorphan_plus_phi = muJet->muon(1)->phi();
1726 m_dimuorphan_minus_phi = muJet->muon(0)->phi();
1727
1728 m_dimuorphan_Deltaphi = m_dimuorphan_plus_phi-m_dimuorphan_minus_phi;
1729 if (m_dimuorphan_Deltaphi > M_PI) m_dimuorphan_Deltaphi -= 2.*M_PI;
1730 if (m_dimuorphan_Deltaphi < -M_PI) m_dimuorphan_Deltaphi += 2.*M_PI;
1731 /*
1732 for (trackingRecHit_iterator hit = muJet->muon(0)->innerTrack()->recHitsBegin(); hit != muJet->muon(0)->innerTrack()->recHitsEnd(); ++hit) {
1733 if ((*hit)->isValid()) {
1734 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1735 m_dimuorphan_minus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1736 break;
1737 }
1738 }
1739 for (trackingRecHit_iterator hit = muJet->muon(1)->innerTrack()->recHitsBegin(); hit != muJet->muon(1)->innerTrack()->recHitsEnd(); ++hit) {
1740 if ((*hit)->isValid()) {
1741 GlobalPoint p = tkgeom->idToDet((*hit)->geographicalId())->toGlobal((*hit)->localPosition());
1742 m_dimuorphan_plus_innerR = pow(p.x()*p.x()+p.y()*p.y(),0.5);
1743 break;
1744 }
1745 }
1746 */
1747
1748 }
1749 m_dimuorphan_Deltaphi_orphan = m_dimuorphan_orphanphi - m_dimuorphan_phi;
1750 if (m_dimuorphan_Deltaphi_orphan > M_PI) m_dimuorphan_Deltaphi_orphan -= 2.*M_PI;
1751 if (m_dimuorphan_Deltaphi_orphan < -M_PI) m_dimuorphan_Deltaphi_orphan += 2.*M_PI;
1752
1753 m_dimuorphan_plus_isGlobal = 0;
1754 m_dimuorphan_plus_isStandAlone = 0;
1755 m_dimuorphan_plus_GlobalHits = -1;
1756 m_dimuorphan_plus_GlobalChi2 = -1;
1757 m_dimuorphan_plus_StandAloneHits = -1;
1758 m_dimuorphan_plus_StandAloneChi2 = -1;
1759
1760 m_dimuorphan_minus_isGlobal = 0;
1761 m_dimuorphan_minus_isStandAlone = 0;
1762 m_dimuorphan_minus_GlobalHits = -1;
1763 m_dimuorphan_minus_GlobalChi2 = -1;
1764 m_dimuorphan_minus_StandAloneHits = -1;
1765 m_dimuorphan_minus_StandAloneChi2 = -1;
1766
1767 if (muJet->daughter(0)->charge() > 0) {
1768 m_dimuorphan_plus_qoverpError = muJet->muon(0)->innerTrack()->qoverpError();
1769 m_dimuorphan_plus_ptError = muJet->muon(0)->innerTrack()->ptError();
1770 m_dimuorphan_plus_phiError = muJet->muon(0)->innerTrack()->phiError();
1771 m_dimuorphan_plus_etaError = muJet->muon(0)->innerTrack()->etaError();
1772 if (muJet->muon(0)->isGlobalMuon()) {
1773 m_dimuorphan_plus_isGlobal = 1;
1774 m_dimuorphan_plus_GlobalHits = muJet->muon(0)->globalTrack()->numberOfValidHits();
1775 m_dimuorphan_plus_GlobalChi2 = muJet->muon(0)->globalTrack()->normalizedChi2();
1776 }
1777 if (muJet->muon(0)->isStandAloneMuon()) {
1778 m_dimuorphan_plus_isStandAlone = 1;
1779 m_dimuorphan_plus_StandAloneHits = muJet->muon(0)->outerTrack()->numberOfValidHits();
1780 m_dimuorphan_plus_StandAloneChi2 = muJet->muon(0)->outerTrack()->normalizedChi2();
1781 }
1782 m_dimuorphan_minus_qoverpError = muJet->muon(1)->innerTrack()->qoverpError();
1783 m_dimuorphan_minus_ptError = muJet->muon(1)->innerTrack()->ptError();
1784 m_dimuorphan_minus_phiError = muJet->muon(1)->innerTrack()->phiError();
1785 m_dimuorphan_minus_etaError = muJet->muon(1)->innerTrack()->etaError();
1786 if (muJet->muon(1)->isGlobalMuon()) {
1787 m_dimuorphan_minus_isGlobal = 1;
1788 m_dimuorphan_minus_GlobalHits = muJet->muon(1)->globalTrack()->numberOfValidHits();
1789 m_dimuorphan_minus_GlobalChi2 = muJet->muon(1)->globalTrack()->normalizedChi2();
1790 }
1791 if (muJet->muon(1)->isStandAloneMuon()) {
1792 m_dimuorphan_minus_isStandAlone = 1;
1793 m_dimuorphan_minus_StandAloneHits = muJet->muon(1)->outerTrack()->numberOfValidHits();
1794 m_dimuorphan_minus_StandAloneChi2 = muJet->muon(1)->outerTrack()->normalizedChi2();
1795 }
1796 }
1797 else {
1798 m_dimuorphan_plus_qoverpError = muJet->muon(1)->innerTrack()->qoverpError();
1799 m_dimuorphan_plus_ptError = muJet->muon(1)->innerTrack()->ptError();
1800 m_dimuorphan_plus_phiError = muJet->muon(1)->innerTrack()->phiError();
1801 m_dimuorphan_plus_etaError = muJet->muon(1)->innerTrack()->etaError();
1802 if (muJet->muon(1)->isGlobalMuon()) {
1803 m_dimuorphan_plus_isGlobal = 1;
1804 m_dimuorphan_plus_GlobalHits = muJet->muon(1)->globalTrack()->numberOfValidHits();
1805 m_dimuorphan_plus_GlobalChi2 = muJet->muon(1)->globalTrack()->normalizedChi2();
1806 }
1807 if (muJet->muon(1)->isStandAloneMuon()) {
1808 m_dimuorphan_plus_isStandAlone = 1;
1809 m_dimuorphan_plus_StandAloneHits = muJet->muon(1)->outerTrack()->numberOfValidHits();
1810 m_dimuorphan_plus_StandAloneChi2 = muJet->muon(1)->outerTrack()->normalizedChi2();
1811 }
1812 m_dimuorphan_minus_qoverpError = muJet->muon(0)->innerTrack()->qoverpError();
1813 m_dimuorphan_minus_ptError = muJet->muon(0)->innerTrack()->ptError();
1814 m_dimuorphan_minus_phiError = muJet->muon(0)->innerTrack()->phiError();
1815 m_dimuorphan_minus_etaError = muJet->muon(0)->innerTrack()->etaError();
1816 if (muJet->muon(0)->isGlobalMuon()) {
1817 m_dimuorphan_minus_isGlobal = 1;
1818 m_dimuorphan_minus_GlobalHits = muJet->muon(0)->globalTrack()->numberOfValidHits();
1819 m_dimuorphan_minus_GlobalChi2 = muJet->muon(0)->globalTrack()->normalizedChi2();
1820 }
1821 if (muJet->muon(0)->isStandAloneMuon()) {
1822 m_dimuorphan_minus_isStandAlone = 1;
1823 m_dimuorphan_minus_StandAloneHits = muJet->muon(0)->outerTrack()->numberOfValidHits();
1824 m_dimuorphan_minus_StandAloneChi2 = muJet->muon(0)->outerTrack()->normalizedChi2();
1825 }
1826 }
1827
1828 m_dimuorphan_dimuonbbbarlike = -1;
1829 if (m_dimuorphan_isoTk > 4.5 || m_dimuorphan_lxy > 0.2) m_dimuorphan_dimuonbbbarlike = 1;
1830 if (m_dimuorphan_isoTk < 4.5 && m_dimuorphan_lxy < 0.2) m_dimuorphan_dimuonbbbarlike = 0;
1831
1832 /*
1833 const pat::Muon *muplus = NULL;
1834 const pat::Muon *muminus = NULL;
1835
1836 if (muJet->muon(0)->charge() > 0) {
1837 muplus = &*muJet->muon(0);
1838 muminus = &*muJet->muon(1);
1839 }
1840 if (muJet->muon(0)->charge() < 0) {
1841 muplus = &*muJet->muon(1);
1842 muminus = &*muJet->muon(0);
1843 }
1844
1845 FreeTrajectoryState plus_initial(GlobalPoint(muplus->vx(), muplus->vy(), muplus->vz()), GlobalVector(muplus->px(), muplus->py(), muplus->pz()), 1, &*magneticField);
1846 FreeTrajectoryState minus_initial(GlobalPoint(muminus->vx(), muminus->vy(), muminus->vz()), GlobalVector(muminus->px(), muminus->py(), muminus->pz()), -1, &*magneticField);
1847
1848 Cylinder::CylinderPointer m_cylinder = Cylinder::build(Cylinder::PositionType(0., 0., 0.), Cylinder::RotationType(), 600.);
1849 Plane::PlanePointer m_plane = Plane::build(Plane::PositionType(0., 0., 800.), Plane::RotationType());
1850
1851
1852 if (fabs(m_dimuorphan_eta) < 1.1) {
1853 TrajectoryStateOnSurface plus_final = propagator->propagate(plus_initial, *m_cylinder);
1854 TrajectoryStateOnSurface minus_final = propagator->propagate(minus_initial, *m_cylinder);
1855 if (plus_final.isValid() && minus_final.isValid()) {
1856 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1857 while (dphi > M_PI) dphi -= 2.*M_PI;
1858 while (dphi < -M_PI) dphi += 2.*M_PI;
1859 m_dimuorphan_dphi = dphi;
1860 }
1861 }
1862 else {
1863 TrajectoryStateOnSurface plus_final = propagator->propagate(plus_initial, *m_plane);
1864 TrajectoryStateOnSurface minus_final = propagator->propagate(minus_initial, *m_plane);
1865 if (plus_final.isValid() && minus_final.isValid()) {
1866 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1867 while (dphi > M_PI) dphi -= 2.*M_PI;
1868 while (dphi < -M_PI) dphi += 2.*M_PI;
1869 m_dimuorphan_dphi = dphi;
1870 }
1871 }
1872
1873 if (fabs(m_dimuorphan_eta) < 0.9) {
1874 TrajectoryStateOnSurface plus_final = propagator->propagate(plus_initial, *m_cylinder100);
1875 TrajectoryStateOnSurface minus_final = propagator->propagate(minus_initial, *m_cylinder100);
1876 if (plus_final.isValid() && minus_final.isValid()) {
1877 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1878 while (dphi > M_PI) dphi -= 2.*M_PI;
1879 while (dphi < -M_PI) dphi += 2.*M_PI;
1880 m_dimuorphan_dphi100 = dphi;
1881 }
1882 plus_final = propagator->propagate(plus_initial, *m_cylinder200);
1883 minus_final = propagator->propagate(minus_initial, *m_cylinder200);
1884 if (plus_final.isValid() && minus_final.isValid()) {
1885 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1886 while (dphi > M_PI) dphi -= 2.*M_PI;
1887 while (dphi < -M_PI) dphi += 2.*M_PI;
1888 m_dimuorphan_dphi200 = dphi;
1889 }
1890 plus_final = propagator->propagate(plus_initial, *m_cylinder300);
1891 minus_final = propagator->propagate(minus_initial, *m_cylinder300);
1892 if (plus_final.isValid() && minus_final.isValid()) {
1893 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1894 while (dphi > M_PI) dphi -= 2.*M_PI;
1895 while (dphi < -M_PI) dphi += 2.*M_PI;
1896 m_dimuorphan_dphi300 = dphi;
1897 }
1898 plus_final = propagator->propagate(plus_initial, *m_cylinder425);
1899 minus_final = propagator->propagate(minus_initial, *m_cylinder425);
1900 if (plus_final.isValid() && minus_final.isValid()) {
1901 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1902 while (dphi > M_PI) dphi -= 2.*M_PI;
1903 while (dphi < -M_PI) dphi += 2.*M_PI;
1904 m_dimuorphan_dphi425 = dphi;
1905 }
1906 plus_final = propagator->propagate(plus_initial, *m_cylinder470);
1907 minus_final = propagator->propagate(minus_initial, *m_cylinder470);
1908 if (plus_final.isValid() && minus_final.isValid()) {
1909 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1910 while (dphi > M_PI) dphi -= 2.*M_PI;
1911 while (dphi < -M_PI) dphi += 2.*M_PI;
1912 m_dimuorphan_dphi470 = dphi;
1913 }
1914 plus_final = propagator->propagate(plus_initial, *m_cylinder510);
1915 minus_final = propagator->propagate(minus_initial, *m_cylinder510);
1916 if (plus_final.isValid() && minus_final.isValid()) {
1917 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1918 while (dphi > M_PI) dphi -= 2.*M_PI;
1919 while (dphi < -M_PI) dphi += 2.*M_PI;
1920 m_dimuorphan_dphi510 = dphi;
1921 }
1922 plus_final = propagator->propagate(plus_initial, *m_cylinder565);
1923 minus_final = propagator->propagate(minus_initial, *m_cylinder565);
1924 if (plus_final.isValid() && minus_final.isValid()) {
1925 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1926 while (dphi > M_PI) dphi -= 2.*M_PI;
1927 while (dphi < -M_PI) dphi += 2.*M_PI;
1928 m_dimuorphan_dphi565 = dphi;
1929 }
1930 plus_final = propagator->propagate(plus_initial, *m_cylinder620);
1931 minus_final = propagator->propagate(minus_initial, *m_cylinder620);
1932 if (plus_final.isValid() && minus_final.isValid()) {
1933 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1934 while (dphi > M_PI) dphi -= 2.*M_PI;
1935 while (dphi < -M_PI) dphi += 2.*M_PI;
1936 m_dimuorphan_dphi620 = dphi;
1937 }
1938 plus_final = propagator->propagate(plus_initial, *m_cylinder670);
1939 minus_final = propagator->propagate(minus_initial, *m_cylinder670);
1940 if (plus_final.isValid() && minus_final.isValid()) {
1941 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1942 while (dphi > M_PI) dphi -= 2.*M_PI;
1943 while (dphi < -M_PI) dphi += 2.*M_PI;
1944 m_dimuorphan_dphi670 = dphi;
1945 }
1946 plus_final = propagator->propagate(plus_initial, *m_cylinder720);
1947 minus_final = propagator->propagate(minus_initial, *m_cylinder720);
1948 if (plus_final.isValid() && minus_final.isValid()) {
1949 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1950 while (dphi > M_PI) dphi -= 2.*M_PI;
1951 while (dphi < -M_PI) dphi += 2.*M_PI;
1952 m_dimuorphan_dphi720 = dphi;
1953 }
1954 plus_final = propagator->propagate(plus_initial, *m_cylinder800);
1955 minus_final = propagator->propagate(minus_initial, *m_cylinder800);
1956 if (plus_final.isValid() && minus_final.isValid()) {
1957 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1958 while (dphi > M_PI) dphi -= 2.*M_PI;
1959 while (dphi < -M_PI) dphi += 2.*M_PI;
1960 m_dimuorphan_dphi800 = dphi;
1961 }
1962 plus_final = propagator->propagate(plus_initial, *m_cylinder900);
1963 minus_final = propagator->propagate(minus_initial, *m_cylinder900);
1964 if (plus_final.isValid() && minus_final.isValid()) {
1965 double dphi = plus_final.globalPosition().phi() - minus_final.globalPosition().phi();
1966 while (dphi > M_PI) dphi -= 2.*M_PI;
1967 while (dphi < -M_PI) dphi += 2.*M_PI;
1968 m_dimuorphan_dphi900 = dphi;
1969 }
1970 }
1971 */
1972
1973 m_dimuorphan_dz1 = muJet->dz(theBeamSpot->position());
1974 m_dimuorphan_dz2 = orphan->innerTrack()->dz(theBeamSpot->position());
1975 m_dimuorphan_deltaz = m_dimuorphan_dz1 - m_dimuorphan_dz2;
1976
1977 m_dimuorphan_muonssize = allmuons->size();
1978
1979 if ((m_dimuorphan_containstrig > 0 || m_dimuorphan_containstrig2 > 0) && fabs(m_dimuorphan_deltaz) < 0.1) m_dimuorphan->Fill();
1980
1981 }
1982
1983 /////////////////////////////////////////////////////////////
1984 // DIMUON + DIMUON //
1985 /////////////////////////////////////////////////////////////
1986
1987 if ( muJets->size() == 2
1988 && (*muJets)[0].numberOfDaughters() == 2
1989 && (*muJets)[1].numberOfDaughters() == 2
1990 && m_trigger > 0) {
1991 const pat::MultiMuon *muJet0 = &((*muJets)[0]);
1992 const pat::MultiMuon *muJet1 = &((*muJets)[1]);
1993
1994 bool muJet0_canBeC = false;
1995 bool muJet1_canBeC = false;
1996 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = hightrigmuons.begin(); iter != hightrigmuons.end(); ++iter) {
1997 if ( muJet0->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable()
1998 && muJet0->sameTrack(&*(muJet0->muon(0)->innerTrack()), &*((*iter)->innerTrack())) ) {
1999 muJet0_canBeC = true;
2000 }
2001 if ( muJet0->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable()
2002 && muJet0->sameTrack(&*(muJet0->muon(1)->innerTrack()), &*((*iter)->innerTrack())) ) {
2003 muJet0_canBeC = true;
2004 }
2005 if ( muJet1->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable()
2006 && muJet1->sameTrack(&*(muJet1->muon(0)->innerTrack()), &*((*iter)->innerTrack())) ) {
2007 muJet1_canBeC = true;
2008 }
2009 if ( muJet1->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable()
2010 && muJet1->sameTrack(&*(muJet1->muon(1)->innerTrack()), &*((*iter)->innerTrack())) ) {
2011 muJet1_canBeC = true;
2012 }
2013 }
2014
2015 for (std::vector<pat::MuonCollection::const_iterator>::const_iterator iter = lowtrigmuons.begin(); iter != lowtrigmuons.end(); ++iter) {
2016 if (muJet0->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
2017 muJet0->sameTrack(&*(muJet0->muon(0)->innerTrack()), &*((*iter)->innerTrack()))) {
2018 m_dimudimu_containstrig2++;
2019 }
2020 if (muJet0->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
2021 muJet0->sameTrack(&*(muJet0->muon(1)->innerTrack()), &*((*iter)->innerTrack()))) {
2022 m_dimudimu_containstrig2++;
2023 }
2024
2025 if (muJet1->muon(0)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
2026 muJet1->sameTrack(&*(muJet1->muon(0)->innerTrack()), &*((*iter)->innerTrack()))) {
2027 m_dimudimu_containstrig2++;
2028 }
2029 if (muJet1->muon(1)->innerTrack().isAvailable() && (*iter)->innerTrack().isAvailable() &&
2030 muJet1->sameTrack(&*(muJet1->muon(1)->innerTrack()), &*((*iter)->innerTrack()))) {
2031 m_dimudimu_containstrig2++;
2032 }
2033 }
2034
2035 const pat::MultiMuon *muJetC = NULL;
2036 const pat::MultiMuon *muJetF = NULL;
2037 if (muJet0_canBeC && muJet1_canBeC) {
2038 m_dimudimu_containstrig = 1;
2039
2040 if (m_trandom3.Integer(2) == 0) {
2041 muJetC = muJet0;
2042 muJetF = muJet1;
2043 }
2044 else {
2045 muJetC = muJet1;
2046 muJetF = muJet0;
2047 }
2048 }
2049 else if (muJet0_canBeC) {
2050 m_dimudimu_containstrig = 1;
2051
2052 muJetC = muJet0;
2053 muJetF = muJet1;
2054 }
2055 else if (muJet1_canBeC) {
2056 m_dimudimu_containstrig = 1;
2057
2058 muJetC = muJet1;
2059 muJetF = muJet0;
2060 }
2061 else {
2062 m_dimudimu_containstrig = 0;
2063 }
2064
2065 if (muJetC != NULL && muJetF != NULL) {
2066 m_dimudimu_orphans = orphans->size();
2067
2068 double total_energy = muJetC->energy() + muJetF->energy();
2069 double total_px = muJetC->px() + muJetF->px();
2070 double total_py = muJetC->py() + muJetF->py();
2071 double total_pz = muJetC->pz() + muJetF->pz();
2072
2073 m_dimudimu_wholemass = sqrt(total_energy*total_energy - total_px*total_px - total_py*total_py - total_pz*total_pz);
2074 m_dimudimu_wholept = sqrt(total_px*total_px + total_py*total_py);
2075
2076 m_dimudimu_deltaphi = muJetC->phi() - muJetF->phi();
2077 while (m_dimudimu_deltaphi > M_PI) m_dimudimu_deltaphi -= 2.*M_PI;
2078 while (m_dimudimu_deltaphi < -M_PI) m_dimudimu_deltaphi += 2.*M_PI;
2079
2080 m_dimudimu_massC = muJetC->mass();
2081 m_dimudimu_ptC = muJetC->pt();
2082 m_dimudimu_etaC = muJetC->eta();
2083 m_dimudimu_phiC = muJetC->phi();
2084 m_dimudimu_drC = muJetC->dRmax();
2085 m_dimudimu_vprobC = -1.;
2086 m_dimudimu_lxyC = -1000.;
2087 m_dimudimu_lxyzC = -1000.;
2088 m_dimudimu_C_isoTk = muJetC->centralTrackIsolation();
2089
2090 m_dimudimu_C_isoTk_3mm = 0;
2091 m_dimudimu_C_isoTk_2mm = 0;
2092 m_dimudimu_C_isoTk_1mm = 0;
2093
2094 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
2095 bool track_is_muon = false;
2096 if (muJetC->sameTrack(&*track,&*(muJetC->muon(0)->innerTrack())) || muJetC->sameTrack(&*track,&*(muJetC->muon(1)->innerTrack()))) track_is_muon = true;
2097 if (!track_is_muon) {
2098 double dphi = muJetC->phi() - track->phi();
2099 if (dphi > M_PI) dphi -= 2.*M_PI;
2100 if (dphi < -M_PI) dphi += 2.*M_PI;
2101 double deta = muJetC->eta() - track->eta();
2102 double dR = sqrt(pow(dphi, 2) + pow(deta, 2));
2103 if (dR < 0.4 && track->pt() > iso_track_pt_treshold) {
2104 double dz = fabs(track->dz(theBeamSpot->position())-muJetC->dz(theBeamSpot->position()));
2105 if (dz < 0.3) m_dimudimu_C_isoTk_3mm += track->pt();
2106 if (dz < 0.2) m_dimudimu_C_isoTk_2mm += track->pt();
2107 if (dz < 0.1) m_dimudimu_C_isoTk_1mm += track->pt();
2108 }
2109 }
2110 }
2111
2112 m_dimudimu_C_isoPF_1mm = 0.0;
2113 for (reco::PFCandidateCollection::const_iterator pfCand = pfCandidates->begin(); pfCand != pfCandidates->end(); ++pfCand) {
2114 if ( pfCand->particleId() == reco::PFCandidate::h ) {
2115 double dPhi = muJetC->phi() - pfCand->phi();
2116 if (dPhi > M_PI) dPhi -= 2.*M_PI;
2117 if (dPhi < -M_PI) dPhi += 2.*M_PI;
2118 double dEta = muJetC->eta() - pfCand->eta();
2119 double dR = sqrt( dPhi*dPhi + dEta*dEta );
2120 // if ( dR < 0.4 && pfCand->pt() > 0.5 ) {
2121 if ( dR < 0.4 ) {
2122 double dz = fabs( pfCand->vertex().z() - theBeamSpot->position().z() - muJetC->dz(theBeamSpot->position()) );
2123 if ( dz < 0.1 ) m_dimudimu_C_isoPF_1mm += pfCand->pt();
2124 }
2125 }
2126 }
2127
2128 m_dimudimu_massF = muJetF->mass();
2129 std::cout << "m_dimudimu_massF = muJetF->mass(): " << m_dimudimu_massF << std::endl;
2130 m_dimudimu_ptF = muJetF->pt();
2131 m_dimudimu_etaF = muJetF->eta();
2132 m_dimudimu_phiF = muJetF->phi();
2133 m_dimudimu_drF = muJetF->dRmax();
2134 m_dimudimu_vprobF = -1.;
2135 m_dimudimu_lxyF = -1000.;
2136 m_dimudimu_lxyzF = -1000.;
2137 m_dimudimu_F_isoTk = muJetF->centralTrackIsolation();
2138
2139 m_dimudimu_F_isoTk_3mm = 0;
2140 m_dimudimu_F_isoTk_2mm = 0;
2141 m_dimudimu_F_isoTk_1mm = 0;
2142
2143 for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
2144 bool track_is_muon = false;
2145 if (muJetF->sameTrack(&*track,&*(muJetF->muon(0)->innerTrack())) || muJetF->sameTrack(&*track,&*(muJetF->muon(1)->innerTrack()))) track_is_muon = true;
2146 if (!track_is_muon) {
2147 double dphi = muJetF->phi() - track->phi();
2148 if (dphi > M_PI) dphi -= 2.*M_PI;
2149 if (dphi < -M_PI) dphi += 2.*M_PI;
2150 double deta = muJetF->eta() - track->eta();
2151 double dR = sqrt(pow(dphi, 2) + pow(deta, 2));
2152 if (dR < 0.4 && track->pt() > iso_track_pt_treshold) {
2153 double dz = fabs(track->dz(theBeamSpot->position())-muJetF->dz(theBeamSpot->position()));
2154 if (dz < 0.3) m_dimudimu_F_isoTk_3mm += track->pt();
2155 if (dz < 0.2) m_dimudimu_F_isoTk_2mm += track->pt();
2156 if (dz < 0.1) m_dimudimu_F_isoTk_1mm += track->pt();
2157 }
2158 }
2159 }
2160
2161 m_dimudimu_F_isoPF_1mm = 0.0;
2162 for (reco::PFCandidateCollection::const_iterator pfCand = pfCandidates->begin(); pfCand != pfCandidates->end(); ++pfCand) {
2163 if ( pfCand->particleId() == reco::PFCandidate::h ) {
2164 double dPhi = muJetF->phi() - pfCand->phi();
2165 if (dPhi > M_PI) dPhi -= 2.*M_PI;
2166 if (dPhi < -M_PI) dPhi += 2.*M_PI;
2167 double dEta = muJetF->eta() - pfCand->eta();
2168 double dR = sqrt( dPhi*dPhi + dEta*dEta );
2169 // if ( dR < 0.4 && pfCand->pt() > 0.5 ) {
2170 if ( dR < 0.4 ) {
2171 double dz = fabs( pfCand->vertex().z() - theBeamSpot->position().z() - muJetF->dz(theBeamSpot->position()) );
2172 if ( dz < 0.1 ) m_dimudimu_F_isoPF_1mm += pfCand->pt();
2173 }
2174 }
2175 }
2176
2177 if (muJetC->vertexValid() && muJetF->vertexValid()) {
2178 total_energy = muJetC->vertexP4().energy() + muJetF->vertexP4().energy();
2179 total_px = muJetC->vertexMomentum().x() + muJetF->vertexMomentum().x();
2180 total_py = muJetC->vertexMomentum().y() + muJetF->vertexMomentum().y();
2181 total_pz = muJetC->vertexMomentum().z() + muJetF->vertexMomentum().z();
2182
2183 m_dimudimu_wholemass = sqrt(total_energy*total_energy - total_px*total_px - total_py*total_py - total_pz*total_pz);
2184 m_dimudimu_wholept = sqrt(total_px*total_px + total_py*total_py);
2185
2186 m_dimudimu_deltaphi = muJetC->vertexMomentum().phi() - muJetF->vertexMomentum().phi();
2187 while (m_dimudimu_deltaphi > M_PI) m_dimudimu_deltaphi -= 2.*M_PI;
2188 while (m_dimudimu_deltaphi < -M_PI) m_dimudimu_deltaphi += 2.*M_PI;
2189
2190 m_dimudimu_massC = muJetC->vertexMass();
2191 m_dimudimu_ptC = muJetC->vertexMomentum().perp();
2192 m_dimudimu_etaC = muJetC->vertexMomentum().eta();
2193 m_dimudimu_phiC = muJetC->vertexMomentum().phi();
2194 m_dimudimu_drC = muJetC->dRmax(true);
2195 m_dimudimu_vprobC = muJetC->vertexProb();
2196 if (closestPrimaryVertex != primaryVertices->end()) {
2197 m_dimudimu_lxyC = muJetC->lxy(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
2198 m_dimudimu_lxyzC = muJetC->lxyz(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
2199 }
2200
2201 m_dimudimu_massF = muJetF->vertexMass();
2202 std::cout << "m_dimudimu_massF = muJetF->vertexMass(): " << m_dimudimu_massF << std::endl;
2203
2204 m_dimudimu_ptF = muJetF->vertexMomentum().perp();
2205 m_dimudimu_etaF = muJetF->vertexMomentum().eta();
2206 m_dimudimu_phiF = muJetF->vertexMomentum().phi();
2207 m_dimudimu_drF = muJetF->dRmax(true);
2208 m_dimudimu_vprobF = muJetF->vertexProb();
2209 if (closestPrimaryVertex != primaryVertices->end()) {
2210 m_dimudimu_lxyF = muJetF->lxy(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
2211 m_dimudimu_lxyzF = muJetF->lxyz(GlobalPoint(closestPrimaryVertex->x(), closestPrimaryVertex->y(), closestPrimaryVertex->z()));
2212 }
2213 }
2214
2215 m_dimudimu_dz1 = muJet0->dz(theBeamSpot->position());
2216 m_dimudimu_dz2 = muJet1->dz(theBeamSpot->position());
2217 m_dimudimu_deltaz = m_dimudimu_dz1 - m_dimudimu_dz2;
2218
2219 if (m_dimudimu_containstrig > 0 && fabs(m_dimudimu_deltaz) < 0.1) m_dimudimu->Fill();
2220 }
2221 }
2222 }
2223
2224
2225 // ------------ method called once each job just before starting event loop ------------
2226 void FitNtuple::beginJob()
2227 {
2228 }
2229
2230 // ------------ method called once each job just after ending the event loop ------------
2231 void FitNtuple::endJob() {
2232 }
2233
2234 //define this as a plug-in
2235 DEFINE_FWK_MODULE(FitNtuple);