ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/codeHeavyBaryons/HeavyBaryonsAnalysis/Lambda/LambdaAnalyzer/src/LambdaAnalyzer.cc
Revision: 1.1.1.1 (vendor branch)
Committed: Mon Dec 13 15:58:16 2010 UTC (14 years, 5 months ago) by mirena
Content type: text/plain
Branch: LambdaBAnalyzer
CVS Tags: START
Changes since 1.1: +0 -0 lines
Log Message:

File Contents

# Content
1 // -*- C++ -*-
2 //
3 // Package: LambdaAnalyzer
4 // Class: LambdaAnalyzer
5 //
6 /**\class LambdaAnalyzer LambdaAnalyzer.cc Lambda/LambdaAnalyzer/src/LambdaAnalyzer.cc
7
8 Description: [one line class summary]
9
10 Implementation:
11 [Notes on implementation]
12 */
13 //
14 // Original Author: Mirena Ivova Rikova,15 R-025,+41227678361,
15 // Created: Wed Aug 4 12:34:06 CEST 2010
16 // $Id$
17 //
18 //
19
20
21 // system include files
22 #include <memory>
23 #include<fstream>
24
25 // user include files
26 #include "FWCore/Framework/interface/Frameworkfwd.h"
27 #include "FWCore/Framework/interface/EDAnalyzer.h"
28
29 #include "FWCore/Framework/interface/Event.h"
30 #include "FWCore/Framework/interface/MakerMacros.h"
31
32 #include "FWCore/ParameterSet/interface/ParameterSet.h"
33
34 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
35
36 #include "DataFormats/TrackReco/interface/Track.h"
37 #include "DataFormats/Candidate/interface/VertexCompositeCandidateFwd.h"
38 #include "DataFormats/Candidate/interface/VertexCompositeCandidate.h"
39 #include "FWCore/ServiceRegistry/interface/Service.h"
40 #include "CommonTools/UtilAlgos/interface/TFileService.h"
41 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
42 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
43 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
44
45 #include "DataFormats/VertexReco/interface/Vertex.h"
46 #include "DataFormats/VertexReco/interface/VertexFwd.h"
47
48 #include "DataFormats/PatCandidates/interface/Muon.h"
49 #include "DataFormats/MuonReco/interface/MuonFwd.h"
50 #include "DataFormats/MuonReco/interface/Muon.h"
51
52 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
53 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
54 //#include "TrackingTools/TransientTrack/interface/TransientTrack.h"
55 //#include "MagneticField/Engine/interface/MagneticField.h"
56
57 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
58 #include "PhysicsTools/CandUtils/interface/AddFourMomenta.h"
59
60 //try trajectory backwards extrapolation
61 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
62 #include "TrackingTools/GeomPropagators/interface/Propagator.h"
63 //triggers
64 #include "DataFormats/Common/interface/TriggerResults.h"
65 #include "FWCore/Common/interface/TriggerNames.h"
66
67 //triggers
68 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
69 #include "CondFormats/L1TObjects/interface/L1GtTriggerMenu.h"
70 #include "CondFormats/DataRecord/interface/L1GtTriggerMenuRcd.h"
71 #include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h"
72
73 //Vertex
74 #include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
75 #include "RecoVertex/VertexTools/interface/VertexDistance3D.h"
76 #include "TVector3.h"
77
78 //Kinematic vertex fit
79
80 #include "RecoVertex/KinematicFit/interface/KinematicParticleVertexFitter.h"
81 #include "RecoVertex/KinematicFitPrimitives/interface/KinematicParticleFactoryFromTransientTrack.h"
82 #include "RecoVertex/KinematicFit/interface/MassKinematicConstraint.h"
83 #include "RecoVertex/KinematicFit/interface/KinematicParticleFitter.h"
84 #include "RecoVertex/KinematicFitPrimitives/interface/MultiTrackKinematicConstraint.h"
85 #include "RecoVertex/KinematicFit/interface/KinematicConstrainedVertexFitter.h"
86 #include "RecoVertex/KinematicFit/interface/TwoTrackMassKinematicConstraint.h"
87
88 //Trigger matching
89
90 #include "MuonAnalysis/MuonAssociators/interface/L1MuonMatcherAlgo.h"
91
92 #include "TH1F.h"
93 #include "TH2F.h"
94 #include "TTree.h"
95
96
97
98 //
99 // class declaration
100 //
101
102 class LambdaAnalyzer : public edm::EDAnalyzer {
103 public:
104 explicit LambdaAnalyzer(const edm::ParameterSet&);
105 ~LambdaAnalyzer();
106
107
108 private:
109 virtual void beginJob() ;
110 virtual void analyze(const edm::Event&, const edm::EventSetup&);
111 virtual void endJob() ;
112
113 virtual bool pass_HLT_Mu0_TkMu0_Jpsi(const edm::Event&);
114 virtual bool pass_HLT_Mu0_TkMu0_OST_Jpsi(const edm::Event&);
115 virtual bool pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v1(const edm::Event&);
116 virtual bool pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v2(const edm::Event&);
117 virtual bool pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v3(const edm::Event&);
118 virtual bool pass_HLT_DoubleMu0_Quarkonium_v1(const edm::Event&);
119
120 virtual void push_anything_back(); //fill tree with default values for _after_fit_ parameters in case of not succeessful kin fit
121
122
123 L1MuonMatcherAlgo matcher_;
124
125
126 unsigned int tmpNvtx;
127
128 TTree* tree_;
129 unsigned int nLambdaB;
130
131 //reco particles info
132 float run;
133 float lumiBlock;
134 float event;
135 // float priVtxX, priVtxY, priVtxZ, priVtxXE, priVtxYE, priVtxZE, priVtxCL;
136 // float priVtxBSx, priVtxBSy, priVtxBSz, priVtxBSxE, priVtxBSyE, priVtxBSzE, priVtxBSCL;
137 float beamSpotX, beamSpotY, beamSpotZ;
138 unsigned int hlt_mu3, hlt_mu5, hlt_mu7, hlt_mu9, hlt_2mu0, hlt_2mu3, hlt_2mu0L2, hlt_2mu3JPsi, hlt_BJPsiMuMu;
139 unsigned int hlt_mu0trk0, hlt_mu3trk0, hlt_mu0trkmu0, hlt_mu3trkmu0, hlt_L1muOpen, hlt_L12muOpen, hlt_L12muOpenTight;
140 unsigned int hlt_mu0_L1muOpen, hlt_mu3_L1muOpen, hlt_mu5_L1muOpen, hlt_mu5trk0, hlt_mu0trkmu0_nocharge;
141 unsigned int hlt_mu0trkmu0_OST, hlt_2mu0_Qv1, hlt_2mu0_QLSv1, hlt_mu0trkmu0_OST_tightV1, hlt_mu0trkmu0_OST_tightV2, hlt_mu0trkmu0_OST_tightV3;
142
143 vector<float> *priVtxX;
144 vector<float> *priVtxY;
145 vector<float> *priVtxZ;
146 vector<float> *priVtxXE;
147 vector<float> *priVtxYE;
148 vector<float> *priVtxZE;
149 vector<float> *priVtxCL;
150 vector<float> *priVtxBSx;
151 vector<float> *priVtxBSy;
152 vector<float> *priVtxBSz;
153 vector<float> *priVtxBSxE;
154 vector<float> *priVtxBSyE;
155 vector<float> *priVtxBSzE;
156 vector<float> *priVtxBSCL;
157
158
159
160
161
162
163
164 //before kinematic fit
165 vector<float> *JpsiVtxX;
166 vector<float> *JpsiVtxY;
167 vector<float> *JpsiVtxZ;
168 vector<float> *LambdaVtxX;
169 vector<float> *LambdaVtxY;
170 vector<float> *LambdaVtxZ;
171 vector<float> *JpsiVtxProb;
172 vector<float> *JpsiVtxCL;
173 vector<float> *LambdaVtxProb;
174 vector<float> *LambdaVtxCL;
175 vector<float> *LambdaBMass;
176 vector<float> *LambdaBPx;
177 vector<float> *LambdaBPy;
178 vector<float> *LambdaBPz;
179 vector<float> *JpsiPx;
180 vector<float> *JpsiPy;
181 vector<float> *JpsiPz;
182 vector<float> *LambdaPx;
183 vector<float> *LambdaPy;
184 vector<float> *LambdaPz;
185 vector<float> *JpsiMass;
186 vector<float> *LambdaMass;
187 vector<float> *mupPx;
188 vector<float> *mupPy;
189 vector<float> *mupPz;
190 vector<float> *mumPx;
191 vector<float> *mumPy;
192 vector<float> *mumPz;
193 vector<float> *trackPPx;
194 vector<float> *trackPPy;
195 vector<float> *trackPPz;
196 vector<float> *trackMPx;
197 vector<float> *trackMPy;
198 vector<float> *trackMPz;
199 vector<float> *JpsiLambdaVtxDistXY;
200 vector<float> *JpsiLambdaVtxDist3D;
201 vector<float> *openingAngleLambdaB;
202 vector<float> *openingAngleLambda;
203 vector<float> *LxyPV;
204 vector<float> *LxyBS;
205 vector<float> *L3DPV;
206 vector<float> *L3DBS;
207 vector<float> *ctauLambdaBxyPV;
208 vector<float> *ctauLambdaBxyBS;
209 vector<float> *ctauLambdaB3DPV;
210 vector<float> *ctauLambdaB3DBS;
211 vector<float> *ctauLambdaXY;
212 vector<float> *dRmup;
213 vector<float> *dRmum;
214 vector<float> *dRtrackp;
215 vector<float> *dRtrackm;
216 vector<float> *dRLambdaVtx;
217 vector<float> *mupChi2;
218 vector<float> *mumChi2;
219 vector<float> *mupNDOF;
220 vector<float> *mumNDOF;
221 vector<float> *mupD0;
222 vector<float> *mumD0;
223 vector<float> *mupNvalidHits;
224 vector<float> *mumNvalidHits;
225 vector<float> *mupIsGlobal;
226 vector<float> *mumIsGlobal;
227 vector<float> *mupIsTracker;
228 vector<float> *mumIsTracker;
229
230 //after kinematic fit
231 vector<float> *JpsiVtxX_afterKinFit;
232 vector<float> *JpsiVtxY_afterKinFit;
233 vector<float> *JpsiVtxZ_afterKinFit;
234 vector<float> *LambdaVtxX_afterKinFit;
235 vector<float> *LambdaVtxY_afterKinFit;
236 vector<float> *LambdaVtxZ_afterKinFit;
237 vector<float> *JpsiVtxProb_afterKinFit;
238 vector<float> *JpsiVtxCL_afterKinFit;
239 vector<float> *LambdaVtxProb_afterKinFit;
240 vector<float> *LambdaVtxCL_afterKinFit;
241 vector<float> *LambdaBMass_afterKinFit;
242 vector<float> *LambdaBPx_afterKinFit;
243 vector<float> *LambdaBPy_afterKinFit;
244 vector<float> *LambdaBPz_afterKinFit;
245 vector<float> *JpsiPx_afterKinFit;
246 vector<float> *JpsiPy_afterKinFit;
247 vector<float> *JpsiPz_afterKinFit;
248 vector<float> *LambdaPx_afterKinFit;
249 vector<float> *LambdaPy_afterKinFit;
250 vector<float> *LambdaPz_afterKinFit;
251 vector<float> *JpsiMass_afterKinFit;
252 vector<float> *LambdaMass_afterKinFit;
253 vector<float> *mupPx_afterKinFit;
254 vector<float> *mupPy_afterKinFit;
255 vector<float> *mupPz_afterKinFit;
256 vector<float> *mumPx_afterKinFit;
257 vector<float> *mumPy_afterKinFit;
258 vector<float> *mumPz_afterKinFit;
259 vector<float> *trackPPx_afterKinFit;
260 vector<float> *trackPPy_afterKinFit;
261 vector<float> *trackPPz_afterKinFit;
262 vector<float> *trackMPx_afterKinFit;
263 vector<float> *trackMPy_afterKinFit;
264 vector<float> *trackMPz_afterKinFit;
265 vector<float> *JpsiLambdaVtxDistXY_afterKinFit;
266 vector<float> *JpsiLambdaVtxDist3D_afterKinFit;
267 vector<float> *openingAngleLambdaB_afterKinFit;
268 vector<float> *openingAngleLambda_afterKinFit;
269 vector<float> *LxyPV_afterKinFit;
270 vector<float> *LxyBS_afterKinFit;
271 vector<float> *L3DPV_afterKinFit;
272 vector<float> *L3DBS_afterKinFit;
273 vector<float> *ctauLambdaBxyPV_afterKinFit;
274 vector<float> *ctauLambdaBxyBS_afterKinFit;
275 vector<float> *ctauLambdaB3DPV_afterKinFit;
276 vector<float> *ctauLambdaB3DBS_afterKinFit;
277 vector<float> *ctauLambdaXY_afterKinFit;
278 vector<float> *LambdaBX_afterKinFit;
279 vector<float> *LambdaBY_afterKinFit;
280 vector<float> *LambdaBZ_afterKinFit;
281 vector<float> *LambdaBVtxX_afterKinFit;
282 vector<float> *LambdaBVtxY_afterKinFit;
283 vector<float> *LambdaBVtxZ_afterKinFit;
284 vector<float> *LambdaBNDOF_afterKinFit;
285 vector<float> *LambdaBchiSquared_afterKinFit;
286 vector<float> *LambdaBVtxProb_afterKinFit;
287 vector<float> *LambdaBPt_afterKinFit;
288 vector<float> *LambdaBEta_afterKinFit;
289 vector<float> *LambdaBPhi_afterKinFit;
290 vector<float> *dRmup_afterKinFit;
291 vector<float> *dRmum_afterKinFit;
292 vector<float> *dRtrackp_afterKinFit;
293 vector<float> *dRtrackm_afterKinFit;
294 vector<float> *dRLambdaVtx_afterKinFit;
295 vector<float> *mupChi2_afterKinFit;
296 vector<float> *mumChi2_afterKinFit;
297 vector<float> *mupNDOF_afterKinFit;
298 vector<float> *mumNDOF_afterKinFit;
299 vector<float> *mupD0_afterKinFit;
300 vector<float> *mumD0_afterKinFit;
301 vector<float> *mupNvalidHits_afterKinFit;
302 vector<float> *mumNvalidHits_afterKinFit;
303 vector<float> *mupIsGlobal_afterKinFit;
304 vector<float> *mumIsGlobal_afterKinFit;
305 vector<float> *mupIsTracker_afterKinFit;
306 vector<float> *mumIsTracker_afterKinFit;
307
308 vector<bool> *mupIsL1Matched;
309 vector<bool> *mumIsL1Matched;
310 vector<bool> *isHLTMatched;
311
312 int genLambdaJpsi;
313
314
315 //MC truth info
316 float truePriVtxX, truePriVtxY, truePriVtxZ;
317 float trueLambdaBPx, trueLambdaBPy,trueLambdaBPz;
318 float trueLambdaBDecayVtxX, trueLambdaBDecayVtxY, trueLambdaBDecayVtxZ;
319 float trueJpsiPx, trueJpsiPy, trueJpsiPz;
320 float trueJpsiDecayVtxX, trueJpsiDecayVtxY, trueJpsiDecayVtxZ;
321 float trueMumPx, trueMumPy, trueMumPz;
322 float trueMupPx, trueMupPy, trueMupPz;
323 float trueLambdaPx, trueLambdaPy, trueLambdaPz;
324 float trueLambdaDecayVtxX,trueLambdaDecayVtxY, trueLambdaDecayVtxZ;
325 float trueTrackPPx, trueTrackPPy, trueTrackPPz;
326 float trueTrackMPx, trueTrackMPy, trueTrackMPz;
327 float trueLxyPV;
328 float trueLxyBS;
329 float trueL3DPV;
330 float trueL3DBS;
331 float truectauLambdaBxyPV;
332 float truectauLambdaB3DPV;
333 float trueLambdaBMass;
334
335 vector<int> truthMatch, truthLambda, truthJpsi;
336
337
338 //parameters from config file
339 bool MCtruth_;
340 bool excludeKshortTracks_;
341 edm::InputTag triggerTag_;
342 edm::InputTag triggerEventTag_;
343
344 std::vector< std::vector<TString> > hltPaths;
345 std::vector< std::vector<std::string> > hltFilters;
346 std::vector< std::vector<unsigned int> > hltRuns;
347
348 // ----------member data ---------------------------
349 };
350
351 //
352 // constants, enums and typedefs
353 //
354
355 //
356 // static data member definitions
357 //
358
359 //
360 // constructors and destructor
361 //
362 LambdaAnalyzer::LambdaAnalyzer(const edm::ParameterSet& iConfig)
363 :
364 matcher_(iConfig),
365 tree_(0),
366 nLambdaB(0),
367 run(0),
368 lumiBlock(0),
369 event(0),
370 beamSpotX(0),
371 beamSpotY(0),
372 beamSpotZ(0),
373 hlt_mu3(0),
374 hlt_mu5(0),
375 hlt_mu7(0),
376 hlt_mu9(0),
377 hlt_2mu0(0),
378 hlt_2mu3(0),
379 hlt_2mu0L2(0),
380 hlt_2mu3JPsi(0),
381 hlt_BJPsiMuMu(0),
382 hlt_mu0trk0(0),
383 hlt_mu3trk0(0),
384 hlt_mu0trkmu0(0),
385 hlt_mu3trkmu0(0),
386 hlt_L1muOpen(0),
387 hlt_L12muOpen(0),
388 hlt_L12muOpenTight(0),
389 hlt_mu0_L1muOpen(0),
390 hlt_mu3_L1muOpen(0),
391 hlt_mu5_L1muOpen(0),
392 hlt_mu5trk0(0),
393 hlt_mu0trkmu0_nocharge(0),
394 hlt_mu0trkmu0_OST(0),
395 hlt_2mu0_Qv1(0),
396 hlt_2mu0_QLSv1(0),
397 hlt_mu0trkmu0_OST_tightV1(0),
398 hlt_mu0trkmu0_OST_tightV2(0),
399 hlt_mu0trkmu0_OST_tightV3(0),
400 priVtxX(0),
401 priVtxY(0),
402 priVtxZ(0),
403 priVtxXE(0),
404 priVtxYE(0),
405 priVtxZE(0),
406 priVtxCL(0),
407 priVtxBSx(0),
408 priVtxBSy(0),
409 priVtxBSz(0),
410 priVtxBSxE(0),
411 priVtxBSyE(0),
412 priVtxBSzE(0),
413 priVtxBSCL(0),
414 JpsiVtxX(0),
415 JpsiVtxY(0),
416 JpsiVtxZ(0),
417 LambdaVtxX(0),
418 LambdaVtxY(0),
419 LambdaVtxZ(0),
420 JpsiVtxProb(0),
421 JpsiVtxCL(0),
422 LambdaVtxProb(0),
423 LambdaVtxCL(0),
424 LambdaBMass(0),
425 LambdaBPx(0),
426 LambdaBPy(0),
427 LambdaBPz(0),
428 JpsiPx(0),
429 JpsiPy(0),
430 JpsiPz(0),
431 LambdaPx(0),
432 LambdaPy(0),
433 LambdaPz(0),
434 JpsiMass(0),
435 LambdaMass(0),
436 mupPx(0),
437 mupPy(0),
438 mupPz(0),
439 mumPx(0),
440 mumPy(0),
441 mumPz(0),
442 trackPPx(0),
443 trackPPy(0),
444 trackPPz(0),
445 trackMPx(0),
446 trackMPy(0),
447 trackMPz(0),
448 JpsiLambdaVtxDistXY(0),
449 JpsiLambdaVtxDist3D(0),
450 openingAngleLambdaB(0),
451 openingAngleLambda(0),
452 LxyPV(0),
453 LxyBS(0),
454 L3DPV(0),
455 L3DBS(0),
456 ctauLambdaBxyPV(0),
457 ctauLambdaBxyBS(0),
458 ctauLambdaB3DPV(0),
459 ctauLambdaB3DBS(0),
460 ctauLambdaXY(0),
461 dRmup(0),
462 dRmum(0),
463 dRtrackp(0),
464 dRtrackm(0),
465 dRLambdaVtx(0),
466 mupChi2(0),
467 mumChi2(0),
468 mupNDOF(0),
469 mumNDOF(0),
470 mupD0(0),
471 mumD0(0),
472 mupNvalidHits(0),
473 mumNvalidHits(0),
474 mupIsGlobal(0),
475 mumIsGlobal(0),
476 mupIsTracker(0),
477 mumIsTracker(0),
478 JpsiVtxX_afterKinFit(0),
479 JpsiVtxY_afterKinFit(0),
480 JpsiVtxZ_afterKinFit(0),
481 LambdaVtxX_afterKinFit(0),
482 LambdaVtxY_afterKinFit(0),
483 LambdaVtxZ_afterKinFit(0),
484 JpsiVtxProb_afterKinFit(0),
485 JpsiVtxCL_afterKinFit(0),
486 LambdaVtxProb_afterKinFit(0),
487 LambdaVtxCL_afterKinFit(0),
488 LambdaBMass_afterKinFit(0),
489 LambdaBPx_afterKinFit(0),
490 LambdaBPy_afterKinFit(0),
491 LambdaBPz_afterKinFit(0),
492 JpsiPx_afterKinFit(0),
493 JpsiPy_afterKinFit(0),
494 JpsiPz_afterKinFit(0),
495 LambdaPx_afterKinFit(0),
496 LambdaPy_afterKinFit(0),
497 LambdaPz_afterKinFit(0),
498 JpsiMass_afterKinFit(0),
499 LambdaMass_afterKinFit(0),
500 mupPx_afterKinFit(0),
501 mupPy_afterKinFit(0),
502 mupPz_afterKinFit(0),
503 mumPx_afterKinFit(0),
504 mumPy_afterKinFit(0),
505 mumPz_afterKinFit(0),
506 trackPPx_afterKinFit(0),
507 trackPPy_afterKinFit(0),
508 trackPPz_afterKinFit(0),
509 trackMPx_afterKinFit(0),
510 trackMPy_afterKinFit(0),
511 trackMPz_afterKinFit(0),
512 JpsiLambdaVtxDistXY_afterKinFit(0),
513 JpsiLambdaVtxDist3D_afterKinFit(0),
514 openingAngleLambdaB_afterKinFit(0),
515 openingAngleLambda_afterKinFit(0),
516 LxyPV_afterKinFit(0),
517 LxyBS_afterKinFit(0),
518 L3DPV_afterKinFit(0),
519 L3DBS_afterKinFit(0),
520 ctauLambdaBxyPV_afterKinFit(0),
521 ctauLambdaBxyBS_afterKinFit(0),
522 ctauLambdaB3DPV_afterKinFit(0),
523 ctauLambdaB3DBS_afterKinFit(0),
524 ctauLambdaXY_afterKinFit(0),
525 LambdaBX_afterKinFit(0),
526 LambdaBY_afterKinFit(0),
527 LambdaBZ_afterKinFit(0),
528 LambdaBVtxX_afterKinFit(0),
529 LambdaBVtxY_afterKinFit(0),
530 LambdaBVtxZ_afterKinFit(0),
531 LambdaBNDOF_afterKinFit(0),
532 LambdaBchiSquared_afterKinFit(0),
533 LambdaBVtxProb_afterKinFit(0),
534 LambdaBPt_afterKinFit(0),
535 LambdaBEta_afterKinFit(0),
536 LambdaBPhi_afterKinFit(0),
537 dRmup_afterKinFit(0),
538 dRmum_afterKinFit(0),
539 dRtrackp_afterKinFit(0),
540 dRtrackm_afterKinFit(0),
541 dRLambdaVtx_afterKinFit(0),
542 mupChi2_afterKinFit(0),
543 mumChi2_afterKinFit(0),
544 mupNDOF_afterKinFit(0),
545 mumNDOF_afterKinFit(0),
546 mupD0_afterKinFit(0),
547 mumD0_afterKinFit(0),
548 mupNvalidHits_afterKinFit(0),
549 mumNvalidHits_afterKinFit(0),
550 mupIsGlobal_afterKinFit(0),
551 mumIsGlobal_afterKinFit(0),
552 mupIsTracker_afterKinFit(0),
553 mumIsTracker_afterKinFit(0),
554 mupIsL1Matched(0),
555 mumIsL1Matched(0),
556 isHLTMatched(0),
557 genLambdaJpsi(0),
558 truePriVtxX(0),
559 truePriVtxY(0),
560 truePriVtxZ(0),
561 trueLambdaBPx(0),
562 trueLambdaBPy(0),
563 trueLambdaBPz(0),
564 trueLambdaBDecayVtxX(0),
565 trueLambdaBDecayVtxY(0),
566 trueLambdaBDecayVtxZ(0),
567 trueJpsiPx(0),
568 trueJpsiPy(0),
569 trueJpsiPz(0),
570 trueJpsiDecayVtxX(0),
571 trueJpsiDecayVtxY(0),
572 trueJpsiDecayVtxZ(0),
573 trueMumPx(0),
574 trueMumPy(0),
575 trueMumPz(0),
576 trueMupPx(0),
577 trueMupPy(0),
578 trueMupPz(0),
579 trueLambdaPx(0),
580 trueLambdaPy(0),
581 trueLambdaPz(0),
582 trueLambdaDecayVtxX(0),
583 trueLambdaDecayVtxY(0),
584 trueLambdaDecayVtxZ(0),
585 trueTrackPPx(0),
586 trueTrackPPy(0),
587 trueTrackPPz(0),
588 trueTrackMPx(0),
589 trueTrackMPy(0),
590 trueTrackMPz(0),
591 trueLxyPV(0),
592 trueLxyBS(0),
593 trueL3DPV(0),
594 trueL3DBS(0),
595 truectauLambdaBxyPV(0),
596 truectauLambdaB3DPV(0),
597 trueLambdaBMass(0),
598 truthMatch(0), truthLambda(0), truthJpsi(0)
599
600 {
601 //now do what ever initialization is needed
602 edm::Service<TFileService> fs;
603 //then declare any histograms needed
604 MCtruth_ = iConfig.getParameter<bool>("MCtruth");
605 excludeKshortTracks_ = iConfig.getParameter<bool>("excludeKshortTracks");
606 triggerTag_ = iConfig.getParameter<edm::InputTag>("TriggerTag");
607 triggerEventTag_ = iConfig.getParameter<edm::InputTag>("TriggerEventTag");
608
609 ifstream hltFile("trigger_20101203.txt");
610 int nVersions = 0;
611 hltFile >> nVersions;
612 for (int iVer=0 ; iVer<nVersions ; iVer++){
613 std::vector<TString> hltPaths_tmp;
614 std::vector<std::string> hltFilters_tmp;
615 std::vector<unsigned int> hltRuns_tmp;
616 int nRuns = 0;
617 hltFile >> nRuns;
618 for (int iRun=0 ; iRun<nRuns ; iRun++){
619 int run_tmp = 0;
620 hltFile >> run_tmp;
621 hltRuns_tmp.push_back(run_tmp);
622 }
623 hltRuns.push_back(hltRuns_tmp);
624 int nPaths = 0;
625 hltFile >> nPaths;
626 for (int iPath=0 ; iPath<nPaths ; iPath++){
627 TString path_tmp = "";
628 std::string filter_tmp = "";
629 hltFile >> path_tmp >> filter_tmp;
630 hltPaths_tmp.push_back(path_tmp);
631 hltFilters_tmp.push_back(filter_tmp);
632 }
633 hltPaths.push_back(hltPaths_tmp);
634 hltFilters.push_back(hltFilters_tmp);
635 }
636
637
638 }
639
640 LambdaAnalyzer::~LambdaAnalyzer()
641 {
642
643 // do anything here that needs to be done at desctruction time
644 // (e.g. close files, deallocate resources etc.)
645
646 }
647
648
649 //
650 // member functions
651 //
652
653 // ------------ method called to for each event ------------
654 void
655 LambdaAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
656 {
657
658 run = iEvent.id().run();
659 lumiBlock = iEvent.luminosityBlock();
660 event = iEvent.id().event();
661
662
663 float pi = 3.14159265;
664
665 using namespace edm;
666 using namespace std;
667 using namespace reco;
668
669 Handle<reco::VertexCompositeCandidateCollection> theKshorts;
670 iEvent.getByLabel( InputTag("generalV0Candidates:Kshort"), theKshorts );
671
672 Handle<reco::VertexCompositeCandidateCollection> theLambdas;
673 iEvent.getByLabel( InputTag("generalV0Candidates:Lambda"), theLambdas );
674
675 Handle<reco::MuonCollection> theMuons;
676 iEvent.getByLabel( InputTag("muons"), theMuons);
677
678 Handle<reco::BeamSpot> theBeamSpot;
679 iEvent.getByLabel(string("offlineBeamSpot"), theBeamSpot );
680
681 //get primary vertex
682 Handle<reco::VertexCollection> thePrimaryVertex;
683 iEvent.getByLabel("offlinePrimaryVertices", thePrimaryVertex);
684
685
686
687 reco::VertexCollection::const_iterator bestVtx = thePrimaryVertex->begin();
688
689
690 //get primary vertex with beam spot constraint
691 Handle<reco::VertexCollection>thePrimaryVertexBS;
692 iEvent.getByLabel("offlinePrimaryVerticesWithBS",thePrimaryVertexBS);
693
694
695 reco::VertexCollection::const_iterator bestVtxBS = thePrimaryVertexBS->begin();
696
697
698
699 //cout<<"event has "<<theLambdas->size()<<" Lambdas, "<<theKshorts->size()<<" Kshorts and "<<theMuons->size()<<" muons"<<endl;
700
701 //HLT results
702 edm::Handle<edm::TriggerResults> hltresults;
703 iEvent.getByLabel(triggerTag_, hltresults);
704 const edm::TriggerNames & triggerNames_ = iEvent.triggerNames(*hltresults);
705 int ntrigs=hltresults->size();
706 if (!MCtruth_) {
707 for (int itrig=0; itrig< ntrigs; itrig++) {
708 TString trigName = triggerNames_.triggerName(itrig);
709 int hltflag = (*hltresults)[itrig].accept();
710 if (trigName=="HLT_DoubleMu3_BJPsi") hlt_BJPsiMuMu = hltflag;
711 if (trigName=="HLT_Mu3") hlt_mu3 = hltflag;
712 if (trigName=="HLT_Mu5") hlt_mu5 = hltflag;
713 if (trigName=="HLT_Mu7") hlt_mu7 = hltflag;
714 if (trigName=="HLT_Mu9") hlt_mu9 = hltflag;
715 if (trigName=="HLT_DoubleMu0") hlt_2mu0 = hltflag;
716 if (trigName=="HLT_DoubleMu3") hlt_2mu3 = hltflag;
717 if (trigName=="HLT_DoubleMu3_Jpsi") hlt_2mu3JPsi = hltflag;
718 if (trigName=="HLT_Mu0_Track0_Jpsi") hlt_mu0trk0 = hltflag;
719 if (trigName=="HLT_Mu3_Track0_Jpsi") hlt_mu3trk0 = hltflag;
720 if (trigName=="HLT_Mu0_TkMu0_Jpsi") hlt_mu0trkmu0 = hltflag;
721 if (trigName=="HLT_Mu3_TkMu0_Jpsi") hlt_mu3trkmu0 = hltflag;
722 if (trigName=="HLT_L1MuOpen") hlt_L1muOpen = hltflag;
723 if (trigName=="HLT_L1DoubleMuOpen") hlt_L12muOpen = hltflag;
724 if (trigName=="HLT_L1DoubleMuOpen_Tight") hlt_L12muOpenTight = hltflag;
725 if (trigName=="HLT_L2DoubleMu0") hlt_2mu0L2 = hltflag;
726 if (trigName=="HLT_Mu0_L1MuOpen") hlt_mu0_L1muOpen = hltflag;
727 if (trigName=="HLT_Mu3_L1MuOpen") hlt_mu3_L1muOpen = hltflag;
728 if (trigName=="HLT_Mu5_L1MuOpen") hlt_mu5_L1muOpen = hltflag;
729 if (trigName=="HLT_Mu5_Track0_Jpsi") hlt_mu5trk0 = hltflag;
730 if (trigName=="HLT_Mu0_TkMu0_Jpsi_NoCharge") hlt_mu0trkmu0_nocharge = hltflag;
731
732 if (trigName=="HLT_Mu0_TkMu0_OST_Jpsi") hlt_mu0trkmu0_OST = hltflag;
733 if (trigName=="HLT_DoubleMu0_Quarkonium_v1") hlt_2mu0_Qv1 = hltflag;
734 if (trigName=="HLT_DoubleMu0_Quarkonium_LS_v1") hlt_2mu0_QLSv1 = hltflag;
735 if (trigName=="HLT_Mu0_TkMu0_OST_Jpsi_Tight_v1") hlt_mu0trkmu0_OST_tightV1 = hltflag;
736 if (trigName=="HLT_Mu0_TkMu0_OST_Jpsi_Tight_v2") hlt_mu0trkmu0_OST_tightV2 = hltflag;
737 if (trigName=="HLT_Mu0_TkMu0_OST_Jpsi_Tight_v3") hlt_mu0trkmu0_OST_tightV3 = hltflag;
738
739 }
740
741 } else {
742 for (int itrig=0; itrig< ntrigs; itrig++) {
743 TString trigName = triggerNames_.triggerName(itrig);
744 int hltflag = (*hltresults)[itrig].accept();
745 if (trigName=="HLT_L1DoubleMuOpen") hlt_L12muOpen = hltflag;
746 if (trigName=="HLT_Mu0_Track0_Jpsi") hlt_mu0trk0 = hltflag;
747 if (trigName=="HLT_DoubleMu0") hlt_2mu0 = hltflag;
748 }
749
750 edm::ESHandle<L1GtTriggerMenu> l1GtMenu;
751 iSetup.get<L1GtTriggerMenuRcd>().get(l1GtMenu);
752 const L1GtTriggerMenu* m_l1GtMenu = l1GtMenu.product();
753 (const_cast<L1GtTriggerMenu*> (m_l1GtMenu))->buildGtConditionMap(); //...ugly
754
755 const AlgorithmMap& algorithmAliasMap = l1GtMenu->gtAlgorithmAliasMap();
756
757 edm::InputTag m_l1GtReadoutRecordTag("gtDigis");
758 edm::Handle<L1GlobalTriggerReadoutRecord> gtReadoutRecord;
759 iEvent.getByLabel(m_l1GtReadoutRecordTag, gtReadoutRecord);
760 const L1GlobalTriggerReadoutRecord* gtReadoutRecordPtr = gtReadoutRecord.product();
761 const std::vector<bool>& gtDecisionWord = gtReadoutRecordPtr->decisionWord();
762
763 bool pass_L1SingleMuOpen = gtDecisionWord[(algorithmAliasMap.find("L1_SingleMuOpen")->second).algoBitNumber()];
764 bool pass_L1SingleMu0 = gtDecisionWord[(algorithmAliasMap.find("L1_SingleMu0")->second).algoBitNumber()];
765 bool pass_L1SingleMu3 = gtDecisionWord[(algorithmAliasMap.find("L1_SingleMu3")->second).algoBitNumber()];
766
767 hlt_mu0trkmu0 = (pass_L1SingleMuOpen||pass_L1SingleMu0||pass_L1SingleMu3) && pass_HLT_Mu0_TkMu0_Jpsi(iEvent);
768 hlt_mu0trkmu0_OST = hlt_mu0trkmu0;
769 hlt_mu0trkmu0_OST_tightV1 = (pass_L1SingleMuOpen||pass_L1SingleMu0||pass_L1SingleMu3) && pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v1(iEvent);
770 hlt_mu0trkmu0_OST_tightV2 = pass_L1SingleMu0 && hlt_mu0trkmu0_OST_tightV1;
771 hlt_mu0trkmu0_OST_tightV3 = pass_L1SingleMu0 && hlt_mu0trkmu0_OST_tightV1;
772 hlt_2mu0_Qv1 = pass_L1SingleMuOpen && pass_HLT_DoubleMu0_Quarkonium_v1(iEvent);
773 }
774
775
776
777 int runVer = 9999;
778
779 if(MCtruth_)
780 {runVer = 6;}
781 else{
782 for (unsigned int iVer=0 ; iVer<hltRuns.size() && runVer == 9999; iVer++){
783 for (unsigned int iRun=0 ; iRun<hltRuns[iVer].size() ; iRun++)
784 if (hltRuns[iVer][iRun] == iEvent.id().run()) {
785 runVer = iVer;
786 continue;
787 }
788 }
789 }
790 matcher_.init(iSetup);
791
792 edm::Handle<std::vector<l1extra::L1MuonParticle> > L1Muons;
793 iEvent.getByLabel(InputTag("l1extraParticles"),L1Muons);
794
795 edm::Handle<trigger::TriggerEvent> hltevent;
796 iEvent.getByLabel(triggerEventTag_, hltevent);
797 const trigger::TriggerObjectCollection& hltobjects(hltevent->getObjects());
798
799 //loop on Lambdas and muons
800 //Lambdas
801 if(theLambdas->size()>0 && theMuons->size()>=2){
802 for(reco::VertexCompositeCandidateCollection::const_iterator iLambda=theLambdas->begin(); iLambda!=theLambdas->end();++iLambda){
803
804 const Candidate * dau0Cand = iLambda->daughter(0);
805 const Candidate * dau1Cand = iLambda->daughter(1);
806 //check for opposite charge of daughters
807 if(dau0Cand->charge()==dau1Cand->charge()){
808 //cout<<"Lambda0 daughters have same sign. Continue."<<endl;
809 continue;
810 }
811 CompositeCandidate LambdaCand;
812 LambdaCand.addDaughter(*dau0Cand);
813 LambdaCand.addDaughter(*dau1Cand);
814
815 AddFourMomenta add;
816 add.set(LambdaCand);
817
818 const reco::RecoChargedCandidate* daughter0 = static_cast<const reco::RecoChargedCandidate*>(iLambda->daughter(0));
819 const reco::RecoChargedCandidate* daughter1 = static_cast<const reco::RecoChargedCandidate*>(iLambda->daughter(1));
820 vector<RecoChargedCandidate> LambdaDaughters;
821 vector<TrackRef> LambdaDauTracks;
822
823 //check momentum and sort: pion first, proton second
824 if( daughter0->momentum().mag2() < daughter1->momentum().mag2() ){
825 LambdaDaughters.push_back(*(dynamic_cast<const reco::RecoChargedCandidate *> (daughter0)) );
826 LambdaDaughters.push_back(*(dynamic_cast<const reco::RecoChargedCandidate *> (daughter1)) );
827 }else{
828 LambdaDaughters.push_back(*(dynamic_cast<const reco::RecoChargedCandidate *> (daughter1)) );
829 LambdaDaughters.push_back(*(dynamic_cast<const reco::RecoChargedCandidate *> (daughter0)) );
830
831 }
832 for(unsigned int j=0; j<LambdaDaughters.size(); ++j){
833 LambdaDauTracks.push_back(LambdaDaughters[j].track() );
834 }
835
836 //loop on all Kshors in the event to check for common tracks with Lambda
837 if(excludeKshortTracks_){
838 for(reco::VertexCompositeCandidateCollection::const_iterator iKshort=theKshorts->begin(); iKshort!=theKshorts->end(); ++iKshort){
839 for(unsigned int idau=0; idau<LambdaDauTracks.size();++idau){
840 if(iKshort->daughter(0)->charge()==LambdaDauTracks[idau]->charge() &&
841 iKshort->daughter(0)->momentum()==LambdaDauTracks[idau]->momentum() ){
842 //cout<<"Lambda0 dau track matched with Ks dau track. Continue."<<endl;
843 continue;
844 }
845 if(iKshort->daughter(1)->charge()==LambdaDauTracks[idau]->charge() &&
846 iKshort->daughter(1)->momentum()==LambdaDauTracks[idau]->momentum() ){
847 //cout<<"Lambda0 dau track matched with Ks dau track. Continue."<<endl;
848 continue;
849 }
850 }
851 }//end of loop on Kshorts
852 }// end of if(excludeKshortTracks)
853
854
855 //now loop on muons
856
857 float deltaRs = 999.;
858 float deltaPhis = 999.;
859
860 //first the positive muon
861 for(reco::MuonCollection::const_iterator iMuonP=theMuons->begin(); iMuonP!=theMuons->end(); ++iMuonP){
862 if(iMuonP->charge() == 1){
863 const Candidate &muPCand = *iMuonP;
864 TrackRef trackP = iMuonP->track();
865 if(trackP.isNull()){continue;}
866 bool match = false;
867 for(unsigned int j=0;j<LambdaDauTracks.size();++j){
868 if(trackP->charge()==LambdaDauTracks[j]->charge() &&
869 trackP->momentum() == LambdaDauTracks[j]->momentum() ){
870 //cout<<"Match found between muP and Lambda track with P = "<<trackP->momentum()<<" and "<<LambdaDauTracks[j]->momentum()<<endl;
871 match = true;
872 }
873 if(match) break;
874 }
875 if(match){// track is already used in making the Lambda
876 //cout<<"Match found between muP and Lambda track"<<endl;
877 continue;
878 }
879 //now the negative muon
880 for(reco::MuonCollection::const_iterator iMuonM=theMuons->begin(); iMuonM!=theMuons->end(); ++iMuonM){
881 if(iMuonM->charge() == -1){
882 const Candidate &muMCand = *iMuonM;
883 TrackRef trackM = iMuonM->track();
884 if(trackM.isNull()){continue;}
885 bool match = false;
886 for(unsigned int j=0;j<LambdaDauTracks.size();++j){
887 if(trackM->charge()==LambdaDauTracks[j]->charge() &&
888 trackM->momentum() == LambdaDauTracks[j]->momentum() ){
889 //cout<<"Match found between muM and Lambda track with P = "<<trackM->momentum()<<" and "<<LambdaDauTracks[j]->momentum()<<endl;
890 match = true;
891 }
892 if(match) break;
893 }
894 if(match){// track is already used in making the Lambda
895 //cout<<"Match found between muM and Lambda track"<<endl;
896 continue;
897 }
898 //four good tracks found
899
900 //KalmanVertexFitter to fit the four particles
901 edm::ESHandle<TransientTrackBuilder> theB;
902 iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder",theB);
903
904 vector<TransientTrack> MuMuTT;
905 MuMuTT.push_back((*theB).build(&trackP));
906 MuMuTT.push_back((*theB).build(&trackM));
907
908 vector<TransientTrack> PiProtonTT;
909 PiProtonTT.push_back((*theB).build(&LambdaDauTracks[0]));
910 PiProtonTT.push_back((*theB).build(&LambdaDauTracks[1]));
911
912
913 KalmanVertexFitter kvf;
914 TransientVertex tvJpsi = kvf.vertex(MuMuTT);
915 TransientVertex tvLambda = kvf.vertex(PiProtonTT);
916 if(tvJpsi.isValid() && tvLambda.isValid()){
917 nLambdaB++;
918 mupChi2->push_back(trackP->chi2());
919 mumChi2->push_back(trackM->chi2());
920 mupNDOF->push_back(trackP->ndof());
921 mumNDOF->push_back(trackM->ndof());
922 mupD0->push_back(trackP->d0());
923 mumD0->push_back(trackM->d0());
924 mupNvalidHits->push_back(trackP->numberOfValidHits());
925 mumNvalidHits->push_back(trackM->numberOfValidHits());
926 mupIsGlobal->push_back(iMuonP->isGlobalMuon());
927 mumIsGlobal->push_back(iMuonM->isGlobalMuon());
928 mupIsTracker->push_back(iMuonP->isTrackerMuon());
929 mumIsTracker->push_back(iMuonM->isTrackerMuon());
930
931
932 CompositeCandidate JpsiCand;
933 JpsiCand.addDaughter(muPCand);
934 JpsiCand.addDaughter(muMCand);
935
936 AddFourMomenta add;
937 add.set(JpsiCand);
938 //create a composite candidate
939 CompositeCandidate LambdaB;
940 LambdaB.addDaughter(*dau0Cand);
941 LambdaB.addDaughter(*dau1Cand);
942 LambdaB.addDaughter(muPCand);
943 LambdaB.addDaughter(muMCand);
944
945 //AddFourMomenta add;
946 add.set(LambdaB);
947
948 //fill candidates information in vectors
949 LambdaBMass->push_back(LambdaB.mass());
950 LambdaBPx->push_back(LambdaB.px());
951 LambdaBPy->push_back(LambdaB.py());
952 LambdaBPz->push_back(LambdaB.pz());
953 JpsiPx->push_back(muPCand.px() + muMCand.px());
954 JpsiPy->push_back(muPCand.py() + muMCand.py());
955 JpsiPz->push_back(muPCand.pz() + muMCand.pz());
956 LambdaPx->push_back(dau0Cand->px() + dau1Cand->px());
957 LambdaPy->push_back(dau0Cand->py() + dau1Cand->py());
958 LambdaPz->push_back(dau0Cand->pz() + dau1Cand->pz());
959 JpsiMass->push_back(JpsiCand.mass());
960 LambdaMass->push_back(LambdaCand.mass());
961 mupPx->push_back(muPCand.px());
962 mupPy->push_back(muPCand.py());
963 mupPz->push_back(muPCand.pz());
964 mumPx->push_back(muMCand.px());
965 mumPy->push_back(muMCand.py());
966 mumPz->push_back(muMCand.pz());
967 for(uint k=0; k<iLambda->numberOfDaughters(); ++k){
968 if(iLambda->daughter(k)->charge()>0){
969
970 trackPPx->push_back(iLambda->daughter(k)->px());
971 trackPPy->push_back(iLambda->daughter(k)->py());
972 trackPPz->push_back(iLambda->daughter(k)->pz());
973 }
974 if(iLambda->daughter(k)->charge()<0){
975 trackMPx->push_back(iLambda->daughter(k)->px());
976 trackMPy->push_back(iLambda->daughter(k)->py());
977 trackMPz->push_back(iLambda->daughter(k)->pz());
978 }
979 }
980
981 ///////////////////////////////////////////////////////////
982 //primary vertex selection
983
984 //select best PV from the min opening angle for LambdaB
985 float minOpeningAngle = -999;
986 reco::VertexCollection::const_iterator bestVtx = thePrimaryVertex->begin();
987 for (reco::VertexCollection::const_iterator vtx = thePrimaryVertex->begin();
988 vtx != thePrimaryVertex->end(); ++vtx){
989 //opening for LambdaB, before kin fit
990 float openingAngleLambdaB_tmp = ( ( (tvJpsi.position().x()-vtx->x())*LambdaB.px()
991 + (tvJpsi.position().y()-vtx->y())*LambdaB.py()
992 + (tvJpsi.position().z()-vtx->z())*LambdaB.pz() ) /
993 sqrt( ( (tvJpsi.position().x()-vtx->x())*(tvJpsi.position().x()-vtx->x()) +
994 (tvJpsi.position().y()-vtx->y())*(tvJpsi.position().y()-vtx->y()) +
995 (tvJpsi.position().z()-vtx->z())*(tvJpsi.position().z()-vtx->z()) ) *
996 ( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() ) )
997 );
998 if(openingAngleLambdaB_tmp > minOpeningAngle){
999 minOpeningAngle = openingAngleLambdaB_tmp;
1000 bestVtx = vtx;
1001 }
1002
1003 }
1004
1005 //fill primary vertex info
1006 //primary vertex
1007 priVtxX->push_back(bestVtx->x());
1008 priVtxY->push_back(bestVtx->y());
1009 priVtxZ->push_back(bestVtx->z());
1010 priVtxXE->push_back(bestVtx->xError());
1011 priVtxYE->push_back(bestVtx->yError());
1012 priVtxZE->push_back(bestVtx->zError());
1013 priVtxCL->push_back(ChiSquaredProbability((double)(bestVtx->chi2()),(double)(bestVtx->ndof())) );
1014
1015 //select best primary vertex with BS constraint from min opening angle for LambdaB
1016 minOpeningAngle = -999;
1017 reco::VertexCollection::const_iterator bestVtxBS = thePrimaryVertexBS->begin();
1018 for (reco::VertexCollection::const_iterator vtxBS = thePrimaryVertexBS->begin();
1019 vtxBS != thePrimaryVertexBS->end(); ++vtxBS){
1020 //opening for LambdaB, before kin fit
1021 float openingAngleLambdaB_tmp = ( ( (tvJpsi.position().x()-vtxBS->x())*LambdaB.px()
1022 + (tvJpsi.position().y()-vtxBS->y())*LambdaB.py()
1023 + (tvJpsi.position().z()-vtxBS->z())*LambdaB.pz() ) /
1024 sqrt( ( (tvJpsi.position().x()-vtxBS->x())*(tvJpsi.position().x()-vtxBS->x()) +
1025 (tvJpsi.position().y()-vtxBS->y())*(tvJpsi.position().y()-vtxBS->y()) +
1026 (tvJpsi.position().z()-vtxBS->z())*(tvJpsi.position().z()-vtxBS->z()) ) *
1027 ( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() ) )
1028 );
1029 if(openingAngleLambdaB_tmp > minOpeningAngle){
1030 minOpeningAngle = openingAngleLambdaB_tmp;
1031 bestVtxBS = vtxBS;
1032 }
1033 }
1034
1035 //fill primary vertex with BS constraint info
1036 priVtxBSx->push_back(bestVtxBS->x());
1037 priVtxBSy->push_back(bestVtxBS->y());
1038 priVtxBSz->push_back(bestVtxBS->z());
1039 priVtxBSxE->push_back(bestVtxBS->xError());
1040 priVtxBSyE->push_back(bestVtxBS->yError());
1041 priVtxBSzE->push_back(bestVtxBS->zError());
1042 priVtxBSCL->push_back(ChiSquaredProbability((double)(bestVtxBS->chi2()),(double)(bestVtxBS->ndof())) );
1043
1044
1045
1046 //end of primary vertex selection
1047 //////////////////////////////////////////////
1048
1049
1050 TVector3 LambdaVtx;
1051 TVector3 JpsiVtx;
1052 VertexDistanceXY vdistXY;
1053 VertexDistance3D vdist3D;
1054
1055 LambdaVtx.SetXYZ(iLambda->vx(), iLambda->vy(), iLambda->vz());
1056 JpsiVtx.SetXYZ(tvJpsi.position().x(), tvJpsi.position().y(), tvJpsi.position().z());
1057
1058 Measurement1D JpsiLambdaVtxDistXYMeasure = vdistXY.distance(tvJpsi, tvLambda);
1059 Measurement1D JpsiLambdaVtxDist3DMeasure = vdist3D.distance(tvJpsi, tvLambda);
1060
1061 JpsiLambdaVtxDistXY->push_back(JpsiLambdaVtxDistXYMeasure.value());
1062 JpsiLambdaVtxDist3D->push_back(JpsiLambdaVtxDist3DMeasure.value());
1063
1064 Vertex vertexJpsi = tvJpsi;
1065 JpsiVtxProb->push_back(TMath::Prob(vertexJpsi.chi2(),(int)vertexJpsi.ndof()));
1066 JpsiVtxCL->push_back(ChiSquaredProbability( (double)(vertexJpsi.chi2()),(double)(vertexJpsi.ndof()) ));
1067
1068
1069
1070 JpsiVtxX->push_back(tvJpsi.position().x());
1071 JpsiVtxY->push_back(tvJpsi.position().y());
1072 JpsiVtxZ->push_back(tvJpsi.position().z());
1073 LambdaVtxX->push_back(iLambda->vx());
1074 LambdaVtxY->push_back(iLambda->vy());
1075 LambdaVtxZ->push_back(iLambda->vz());
1076
1077 Vertex vertexLambda = tvLambda;
1078 LambdaVtxProb->push_back(TMath::Prob(vertexLambda.chi2(),(int)vertexLambda.ndof()));
1079 LambdaVtxCL->push_back(ChiSquaredProbability( (double)(vertexLambda.chi2()),(double)(vertexLambda.ndof()) ) );
1080
1081 //cos of angle between momentum and decay length
1082 //opening for LambdaB, before kin fit
1083 openingAngleLambdaB->push_back( ( (tvJpsi.position().x()-bestVtx->x())*LambdaB.px()
1084 + (tvJpsi.position().y()-bestVtx->y())*LambdaB.py()
1085 + (tvJpsi.position().z()-bestVtx->z())*LambdaB.pz() ) /
1086 sqrt( ( (tvJpsi.position().x()-bestVtx->x())*(tvJpsi.position().x()-bestVtx->x()) +
1087 (tvJpsi.position().y()-bestVtx->y())*(tvJpsi.position().y()-bestVtx->y()) +
1088 (tvJpsi.position().z()-bestVtx->z())*(tvJpsi.position().z()-bestVtx->z()) ) *
1089 ( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() ) )
1090 );
1091
1092 openingAngleLambda->push_back( ( (iLambda->vx()-tvJpsi.position().x())*(dau0Cand->px() + dau1Cand->px()) +
1093 (iLambda->vy()-tvJpsi.position().y())*(dau0Cand->py() + dau1Cand->py()) +
1094 (iLambda->vz()-tvJpsi.position().z())*(dau0Cand->pz() + dau1Cand->pz()) ) /
1095 sqrt( (
1096 (iLambda->vx()-tvJpsi.position().x())*(iLambda->vx()-tvJpsi.position().x()) +
1097 (iLambda->vy()-tvJpsi.position().y())*(iLambda->vy()-tvJpsi.position().y()) +
1098 (iLambda->vz()-tvJpsi.position().z())*(iLambda->vz()-tvJpsi.position().z()) ) *
1099 (
1100 (dau0Cand->px() + dau1Cand->px())*(dau0Cand->px() + dau1Cand->px()) +
1101 (dau0Cand->py() + dau1Cand->py())*(dau0Cand->py() + dau1Cand->py()) +
1102 (dau0Cand->pz() + dau1Cand->pz())*(dau0Cand->pz() + dau1Cand->pz())
1103 )
1104 )
1105 );
1106
1107 //decay lendth L and proper decay length ctau=L*m/p with beam spot and PV
1108 LxyPV->push_back( ( (tvJpsi.position().x()-bestVtx->x())*LambdaB.px() + (tvJpsi.position().y()-bestVtx->y())*LambdaB.py() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py()) );
1109 LxyBS->push_back( ( (tvJpsi.position().x()-theBeamSpot->position().x())*LambdaB.px() + (tvJpsi.position().y()-theBeamSpot->position().y())*LambdaB.py() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() ) );
1110 L3DPV->push_back( ( (tvJpsi.position().x()-bestVtx->x())*LambdaB.px() + (tvJpsi.position().y()-bestVtx->y())*LambdaB.py() + (tvJpsi.position().z()-bestVtx->z())*LambdaB.pz()) /sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() ) );
1111 L3DBS->push_back( ( (tvJpsi.position().x()-theBeamSpot->position().x())*LambdaB.px() + (tvJpsi.position().y()-theBeamSpot->position().y())*LambdaB.py() +(tvJpsi.position().z()-theBeamSpot->position().z())*LambdaB.pz() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() ) );
1112
1113 float LxyPVtmp = ( (tvJpsi.position().x()-bestVtx->x())*LambdaB.px() + (tvJpsi.position().y()-bestVtx->y())*LambdaB.py() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py());
1114 float LxyBStmp = ( (tvJpsi.position().x()-theBeamSpot->position().x())*LambdaB.px() + (tvJpsi.position().y()-theBeamSpot->position().y())*LambdaB.py() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() );
1115 float L3DPVtmp = ( (tvJpsi.position().x()-bestVtx->x())*LambdaB.px() + (tvJpsi.position().y()-bestVtx->y())*LambdaB.py() + (tvJpsi.position().z()-bestVtx->z())*LambdaB.pz()) /sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() );
1116 float L3DBStmp = ( (tvJpsi.position().x()-theBeamSpot->position().x())*LambdaB.px() + (tvJpsi.position().y()-theBeamSpot->position().y())*LambdaB.py() +(tvJpsi.position().z()-theBeamSpot->position().z())*LambdaB.pz() )/sqrt( LambdaB.px()*LambdaB.px() + LambdaB.py()*LambdaB.py() + LambdaB.pz()*LambdaB.pz() );
1117
1118 ctauLambdaBxyPV->push_back(LxyPVtmp*LambdaB.mass()/LambdaB.pt());
1119 ctauLambdaBxyBS->push_back(LxyBStmp*LambdaB.mass()/LambdaB.pt());
1120 ctauLambdaB3DPV->push_back(L3DPVtmp*LambdaB.mass()/LambdaB.p());
1121 ctauLambdaB3DBS->push_back(L3DBStmp*LambdaB.mass()/LambdaB.p());
1122
1123 float LambdaLxytmp = ( (iLambda->vx()-tvJpsi.position().x())*(dau0Cand->px()+dau1Cand->px()) +
1124 (iLambda->vy()-tvJpsi.position().y())*(dau0Cand->py()+dau1Cand->py()) ) /
1125 sqrt( (dau0Cand->px()+dau1Cand->px())*(dau0Cand->px()+dau1Cand->px()) + (dau0Cand->py()+dau1Cand->py())*(dau0Cand->py()+dau1Cand->py()) );
1126 ctauLambdaXY->push_back(LambdaLxytmp*LambdaCand.mass()/LambdaCand.pt() );
1127
1128 ///////////////////////////////////////////////
1129 // Trigger Matching
1130 TrajectoryStateOnSurface propagated;
1131 int match1 = matcher_.match(*iMuonP,*L1Muons,deltaRs,deltaPhis, propagated);
1132 int match2 = matcher_.match(*iMuonM,*L1Muons,deltaRs,deltaPhis, propagated);
1133
1134 mupIsL1Matched->push_back(match1!=-1);
1135 mumIsL1Matched->push_back(match2!=-1);
1136
1137 bool mupHLTmatched = false;
1138 bool mumHLTmatched = false;
1139 bool hltMatched = false;
1140
1141 for (unsigned int iPath=0 ; iPath<hltPaths[runVer].size() && !hltMatched ; iPath++){
1142 for (int itrig=0; itrig< ntrigs; itrig++) {
1143 TString trigName = triggerNames_.triggerName(itrig);
1144 if (trigName == hltPaths[runVer][iPath] && (*hltresults)[itrig].accept()) {
1145
1146 trigger::size_type index(0);
1147 std::string hltObjTag = hltFilters[runVer][iPath]+"::HLT";
1148 edm::InputTag it_track(hltObjTag);
1149
1150 index=hltevent->filterIndex(it_track);
1151 if (index<hltevent->sizeFilters()) {
1152 const trigger::Keys& KEYS (hltevent->filterKeys(index));
1153 const trigger::size_type nTO(KEYS.size());
1154 for (trigger::size_type i=0; i!=nTO; ++i) {
1155 const trigger::TriggerObject& TO( hltobjects.at(KEYS[i]) );
1156 if (deltaR(iMuonP->eta(),iMuonP->phi(),TO.eta(),TO.phi())<0.02) mupHLTmatched=true;
1157 if (deltaR(iMuonM->eta(),iMuonM->phi(),TO.eta(),TO.phi())<0.02) mumHLTmatched=true;
1158 }
1159 }
1160 }
1161 }
1162
1163 TString filterName(hltFilters[runVer][iPath]);
1164 /*
1165 if (filterName.Contains("Single")) {
1166 if ((mupHLTmatched || mumHLTmatched)) hltMatched = true;
1167 } else {
1168 if ((mupHLTmatched && mumHLTmatched)) hltMatched = true;
1169 }
1170 */
1171 if (filterName.Contains("Double") || filterName.Contains("DiMuon") || filterName.Contains("Jpsi")) {
1172 if ((mupHLTmatched && mumHLTmatched)) hltMatched = true;
1173 } else {
1174 if ((mupHLTmatched || mumHLTmatched)) hltMatched = true;
1175 }
1176
1177 }
1178
1179 isHLTMatched->push_back(hltMatched);
1180
1181 ///////////////////////////////////////////////
1182 // Kinematic fit
1183
1184 KinematicParticleFactoryFromTransientTrack pFactory;
1185 KinematicParticleVertexFitter Fitter;
1186 KinematicParticleFitter constFitter;
1187 KinematicConstrainedVertexFitter constVertexFitter;
1188
1189 const ParticleMass muonMass(0.1056583);
1190 float muonSigma = 1E-10;
1191
1192 const ParticleMass protonMass(0.938272013);
1193 float protonSigma = 1E-10;
1194
1195 const ParticleMass pionMass(0.139570);
1196 float pionSigma = 0.000016;
1197
1198 vector<RefCountedKinematicParticle> allParticlesMu;
1199 allParticlesMu.push_back(pFactory.particle (MuMuTT[0], muonMass, float(0), float(0), muonSigma));
1200 allParticlesMu.push_back(pFactory.particle (MuMuTT[1], muonMass, float(0), float(0), muonSigma));
1201 RefCountedKinematicTree JpsiTree = Fitter.fit(allParticlesMu);
1202 if(JpsiTree->isEmpty()) {push_anything_back(); continue;}
1203
1204 double nominalJpsiMass = 3.096916;
1205 float jpsiMsigma = 0.00004;
1206 KinematicConstraint *jpsi_const = new MassKinematicConstraint( nominalJpsiMass, jpsiMsigma);
1207
1208 JpsiTree = constFitter.fit(jpsi_const,JpsiTree);
1209
1210 if(JpsiTree->isEmpty()) {push_anything_back();continue;}
1211
1212 JpsiTree->movePointerToTheTop();
1213 RefCountedKinematicParticle Jpsi_branch = JpsiTree->currentParticle();
1214
1215 if (!Jpsi_branch->currentState().isValid()) {push_anything_back();continue;}
1216
1217 vector<RefCountedKinematicParticle> allLambdaDaughters;
1218 allLambdaDaughters.push_back(pFactory.particle (PiProtonTT[0], pionMass, float(0), float(0), pionSigma));
1219 allLambdaDaughters.push_back(pFactory.particle (PiProtonTT[1], protonMass, float(0), float(0), protonSigma));
1220 RefCountedKinematicTree LambdaTree = Fitter.fit(allLambdaDaughters);
1221 if(LambdaTree->isEmpty()) {push_anything_back();continue;}
1222
1223 double nominalLambdaMass = 1.115683;
1224 float lambdaMsigma = 0.000006;
1225 KinematicConstraint *lambda_const = new MassKinematicConstraint( nominalLambdaMass, lambdaMsigma);
1226
1227 LambdaTree = constFitter.fit(lambda_const,LambdaTree);
1228
1229 if(LambdaTree->isEmpty()) {push_anything_back();continue;}
1230
1231 LambdaTree->movePointerToTheTop();
1232 RefCountedKinematicParticle Lambda_branch = LambdaTree->currentParticle();
1233 if (!Lambda_branch->currentState().isValid()) {push_anything_back();continue;}
1234
1235 LambdaTree->movePointerToTheFirstChild();
1236 RefCountedKinematicParticle Pion_branch = LambdaTree->currentParticle();
1237 if (!Pion_branch->currentState().isValid()) {push_anything_back();continue;}
1238
1239 LambdaTree->movePointerToTheNextChild();
1240 RefCountedKinematicParticle Proton_branch = LambdaTree->currentParticle();
1241 if (!Proton_branch->currentState().isValid()) {push_anything_back();continue;}
1242
1243
1244 vector<RefCountedKinematicParticle> allLambdaBDaughters;
1245 allLambdaBDaughters.push_back(pFactory.particle (MuMuTT[0], muonMass, float(0), float(0), muonSigma));
1246 allLambdaBDaughters.push_back(pFactory.particle (MuMuTT[1], muonMass, float(0), float(0), muonSigma));
1247 allLambdaBDaughters.push_back(Lambda_branch);
1248
1249 MultiTrackKinematicConstraint *jpsi_mtc = new TwoTrackMassKinematicConstraint(nominalJpsiMass);
1250 RefCountedKinematicTree LambdaBTree = constVertexFitter.fit(allLambdaBDaughters,jpsi_mtc);
1251 if(LambdaBTree->isEmpty()) {push_anything_back();continue;}
1252
1253
1254 LambdaBTree->movePointerToTheTop();
1255 RefCountedKinematicParticle LambdaB_branch = LambdaBTree->currentParticle();
1256 if (!LambdaB_branch->currentState().isValid()) {push_anything_back();continue;}
1257 RefCountedKinematicVertex LambdaB_DecayVtx = LambdaBTree->currentDecayVertex();
1258
1259 AlgebraicVector7 LambdaB_par = LambdaB_branch->currentState().kinematicParameters().vector();
1260 AlgebraicMatrix77 LambdaB_err = LambdaB_branch->currentState().kinematicParametersError().matrix();
1261 GlobalVector LambdaBvector(LambdaB_par[3], LambdaB_par[4], LambdaB_par[5]); // the fitted momentum vector of the LambdaB
1262
1263
1264 LambdaBTree->movePointerToTheFirstChild();
1265 RefCountedKinematicParticle MuP_branch = LambdaBTree->currentParticle();
1266 if (!MuP_branch->currentState().isValid()) {push_anything_back();continue;}
1267
1268 LambdaBTree->movePointerToTheNextChild();
1269 RefCountedKinematicParticle MuM_branch = LambdaBTree->currentParticle();
1270 if (!MuM_branch->currentState().isValid()) {push_anything_back();continue;}
1271
1272 LambdaBTree->movePointerToTheNextChild();
1273 Lambda_branch = LambdaBTree->currentParticle();
1274 if (!Lambda_branch->currentState().isValid()) {push_anything_back();continue;}
1275
1276 //fill candidates information in vectors
1277 LambdaBMass_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().mass());
1278 LambdaBPx_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().momentum().x());
1279 LambdaBPy_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().momentum().y());
1280 LambdaBPz_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().momentum().z());
1281 LambdaBX_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().position().x());
1282 LambdaBY_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().position().y());
1283 LambdaBZ_afterKinFit->push_back(LambdaB_branch->currentState().kinematicParameters().position().z());
1284 LambdaBEta_afterKinFit->push_back(LambdaBvector.eta());
1285 LambdaBPhi_afterKinFit->push_back(LambdaBvector.phi());
1286 LambdaBPt_afterKinFit->push_back(LambdaBvector.perp());
1287 LambdaBNDOF_afterKinFit->push_back((int)LambdaB_branch->degreesOfFreedom());
1288 LambdaBchiSquared_afterKinFit->push_back(LambdaB_branch->chiSquared());
1289 LambdaBVtxProb_afterKinFit->push_back(TMath::Prob(LambdaB_branch->chiSquared(), (int)LambdaB_branch->degreesOfFreedom()) );
1290 LambdaBVtxX_afterKinFit->push_back(LambdaB_DecayVtx->position().x());
1291 LambdaBVtxY_afterKinFit->push_back(LambdaB_DecayVtx->position().y());
1292 LambdaBVtxZ_afterKinFit->push_back(LambdaB_DecayVtx->position().z());
1293
1294 mumPx_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[3]);
1295 mumPy_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[4]);
1296 mumPz_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[5]);
1297
1298 mupPx_afterKinFit->push_back((MuP_branch->currentState().kinematicParameters().vector())[3]);
1299 mupPy_afterKinFit->push_back((MuP_branch->currentState().kinematicParameters().vector())[4]);
1300 mupPz_afterKinFit->push_back((MuP_branch->currentState().kinematicParameters().vector())[5]);
1301
1302 JpsiPx_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[3]+(MuP_branch->currentState().kinematicParameters().vector())[3]);
1303 JpsiPy_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[4]+(MuP_branch->currentState().kinematicParameters().vector())[4]);
1304 JpsiPz_afterKinFit->push_back((MuM_branch->currentState().kinematicParameters().vector())[5]+(MuP_branch->currentState().kinematicParameters().vector())[5]);
1305
1306 LambdaPx_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[3]);
1307 LambdaPy_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[4]);
1308 LambdaPz_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[5]);
1309
1310 JpsiVtxX_afterKinFit->push_back((LambdaB_branch->currentState().kinematicParameters().vector())[0]);
1311 JpsiVtxY_afterKinFit->push_back((LambdaB_branch->currentState().kinematicParameters().vector())[1]);
1312 JpsiVtxZ_afterKinFit->push_back((LambdaB_branch->currentState().kinematicParameters().vector())[2]);
1313
1314 LambdaVtxX_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[0]);
1315 LambdaVtxY_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[1]);
1316 LambdaVtxZ_afterKinFit->push_back((Lambda_branch->currentState().kinematicParameters().vector())[2]);
1317
1318 LambdaVtxProb_afterKinFit->push_back(TMath::Prob(Lambda_branch->chiSquared(),(int)Lambda_branch->degreesOfFreedom()));
1319 JpsiVtxProb_afterKinFit->push_back(TMath::Prob(Jpsi_branch->chiSquared(),(int)Jpsi_branch->degreesOfFreedom()));
1320
1321 float LambdaBVtxX_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[0];
1322 float LambdaBVtxY_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[1];
1323 float LambdaBVtxZ_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[2];
1324
1325 float LambdaBPx_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[3];
1326 float LambdaBPy_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[4];
1327 float LambdaBPz_fit_tmp = (LambdaB_branch->currentState().kinematicParameters().vector())[5];
1328
1329 float PVx_tmp = bestVtx->x();
1330 float PVy_tmp = bestVtx->y();
1331 float PVz_tmp = bestVtx->z();
1332
1333 float LambdaBMass_fit_tmp = LambdaB_branch->currentState().kinematicParameters().mass();
1334 //float LambdaMass_fit_tmp = Lambda_branch->currentState().kinematicParameters().mass();
1335
1336
1337 float oAngleLambdaB_fit_tmp =
1338 ( (LambdaBVtxX_fit_tmp - PVx_tmp)*LambdaBPx_fit_tmp
1339 + (LambdaBVtxY_fit_tmp - PVy_tmp)*LambdaBPy_fit_tmp
1340 + (LambdaBVtxZ_fit_tmp - PVz_tmp)*LambdaBPz_fit_tmp ) /
1341 sqrt( ((LambdaBVtxX_fit_tmp - PVx_tmp)*(LambdaBVtxX_fit_tmp - PVx_tmp) +
1342 (LambdaBVtxY_fit_tmp - PVy_tmp)*(LambdaBVtxY_fit_tmp - PVy_tmp) +
1343 (LambdaBVtxZ_fit_tmp - PVz_tmp)*(LambdaBVtxZ_fit_tmp - PVz_tmp)) * (
1344 (LambdaBPx_fit_tmp * LambdaBPx_fit_tmp) +
1345 (LambdaBPy_fit_tmp * LambdaBPy_fit_tmp) +
1346 (LambdaBPz_fit_tmp * LambdaBPz_fit_tmp) ) );
1347
1348
1349 // float LambdaVtxX_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[0];
1350 // float LambdaVtxY_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[1];
1351 // float LambdaVtxZ_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[2];
1352
1353 // float LambdaPx_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[3];
1354 // float LambdaPy_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[4];
1355 // float LambdaPz_fit_tmp = (Lambda_branch->currentState().kinematicParameters().vector())[5];
1356
1357
1358
1359 // float oAngleLambda_fit_tmp =
1360 // ( (LambdaVtxX_fit_tmp - LambdaBVtxX_fit_tmp)*LambdaPx_fit_tmp
1361 // + (LambdaVtxY_fit_tmp - LambdaBVtxY_fit_tmp)*LambdaPy_fit_tmp
1362 // + (LambdaVtxZ_fit_tmp - LambdaBVtxZ_fit_tmp)*LambdaPz_fit_tmp )/
1363 // sqrt( ((LambdaVtxX_fit_tmp - LambdaBVtxX_fit_tmp)*(LambdaVtxX_fit_tmp - LambdaBVtxX_fit_tmp) +
1364 // (LambdaVtxY_fit_tmp - LambdaBVtxY_fit_tmp)*(LambdaVtxY_fit_tmp - LambdaBVtxY_fit_tmp) +
1365 // (LambdaVtxZ_fit_tmp - LambdaBVtxZ_fit_tmp)*(LambdaVtxZ_fit_tmp - LambdaBVtxZ_fit_tmp)) * (
1366 // (LambdaPx_fit_tmp*LambdaPx_fit_tmp) +
1367 // (LambdaPy_fit_tmp*LambdaPy_fit_tmp) +
1368 // (LambdaPz_fit_tmp*LambdaPz_fit_tmp) ) );
1369
1370
1371
1372
1373 openingAngleLambdaB_afterKinFit->push_back(oAngleLambdaB_fit_tmp);
1374 //openingAngleLambda_afterKinFit->push_back(oAngleLambda_fit_tmp);
1375
1376
1377
1378 //cout<<"two muons have valid vtx"<<endl;
1379 tmpNvtx++;
1380
1381
1382 float L3DPV_fit_tmp = ( (LambdaBVtxX_fit_tmp-PVx_tmp)*LambdaBPx_fit_tmp + (LambdaBVtxY_fit_tmp-PVy_tmp)*LambdaBPy_fit_tmp + (LambdaBVtxZ_fit_tmp-PVz_tmp)*LambdaBPz_fit_tmp )/
1383 sqrt(LambdaBPx_fit_tmp*LambdaBPx_fit_tmp + LambdaBPy_fit_tmp*LambdaBPy_fit_tmp + LambdaBPz_fit_tmp*LambdaBPz_fit_tmp);
1384 ctauLambdaB3DPV_afterKinFit->push_back(L3DPV_fit_tmp*LambdaBMass_fit_tmp/sqrt(LambdaBPx_fit_tmp*LambdaBPx_fit_tmp+LambdaBPy_fit_tmp*LambdaBPy_fit_tmp+LambdaBPz_fit_tmp*LambdaBPz_fit_tmp) );
1385
1386 float LxyPV_fit_tmp = ( (LambdaBVtxX_fit_tmp-PVx_tmp)*LambdaBPx_fit_tmp + (LambdaBVtxY_fit_tmp-PVy_tmp)*LambdaBPy_fit_tmp )/
1387 sqrt(LambdaBPx_fit_tmp*LambdaBPx_fit_tmp + LambdaBPy_fit_tmp*LambdaBPy_fit_tmp);
1388 ctauLambdaBxyPV_afterKinFit->push_back(LxyPV_fit_tmp*LambdaBMass_fit_tmp/sqrt(LambdaBPx_fit_tmp*LambdaBPx_fit_tmp+LambdaBPy_fit_tmp*LambdaBPy_fit_tmp) );
1389
1390 /*
1391 float LambdaLxy_fit_tmp = ( (LambdaVtxX_fit_tmp-LambdaBVtxX_fit_tmp)*LambdaPx_fit_tmp + (LambdaVtxY_fit_tmp-LambdaBVtxY_fit_tmp)*LambdaPy_fit_tmp )/
1392 sqrt(LambdaPx_fit_tmp*LambdaPx_fit_tmp + LambdaPy_fit_tmp*LambdaPy_fit_tmp);
1393
1394 ctauLambdaXY_afterKinFit->push_back(LambdaLxy_fit_tmp*LambdaMass_fit_tmp/sqrt(LambdaPx_fit_tmp*LambdaPx_fit_tmp+LambdaPy_fit_tmp*LambdaPy_fit_tmp) );
1395
1396 */
1397
1398
1399 }//end of if MuMu vertex is valid
1400 }//end of if(iMuonP->charge() == -1)
1401 }//end of loop of iMuonP
1402 }//end of if(iMuonP->charge() == 1)
1403 }//end of loop of iMuonP
1404 }//end of loop on Lambdas
1405 }//end of loop on Lambdas and muons combined (if ...)
1406
1407 // Beam Spot
1408 beamSpotX = theBeamSpot->position().x();
1409 beamSpotY = theBeamSpot->position().y();
1410 beamSpotZ = theBeamSpot->position().z();
1411
1412
1413
1414 ///////////////////////////////
1415 /// MC part
1416 if(nLambdaB > 0 && MCtruth_){
1417
1418 // Beam Spot
1419 beamSpotX = theBeamSpot->position().x();
1420 beamSpotY = theBeamSpot->position().y();
1421 beamSpotZ = theBeamSpot->position().z();
1422
1423 Handle<GenParticleCollection>theGenParticles;
1424 iEvent.getByLabel("genParticlePlusGEANT", theGenParticles);
1425
1426 genLambdaJpsi = -1;
1427 for(size_t iGenParticle=0; iGenParticle<theGenParticles->size();++iGenParticle){
1428 const Candidate & genLambdab = (*theGenParticles)[iGenParticle];
1429 if(abs( genLambdab.pdgId() )==5122){//if Lambda^0_b is generated
1430 //cout<<"found gen Lambda_b";
1431 int iJpsi(-1), iLambda(-1);
1432 bool wrong = false;
1433 for(uint i=0; i<genLambdab.numberOfDaughters();++i){
1434 //check Lambda_b for Jpsi and Lambda daughters
1435 const Candidate * genDau = genLambdab.daughter(i);
1436 //cout<<" ="<<genDau->pdgId();
1437 int imu1(-1), imu2(-1), ipi(-1), ip(-1);
1438 for(uint j=0; j<genDau->numberOfDaughters();++j){
1439 //cout<<" =="<<genDau->daughter(j)->pdgId();
1440 if(genDau->daughter(j)->pdgId()==13 && genDau->pdgId()==443)imu1=j;
1441 if(genDau->daughter(j)->pdgId()==-13 && genDau->pdgId()==443)imu2=j;
1442 if(abs(genDau->daughter(j)->pdgId())==211 && abs(genDau->pdgId())==3122 )ipi=j;
1443 if(abs(genDau->daughter(j)->pdgId())==2212 && abs(genDau->pdgId())==3122)ip=j;
1444
1445 if(genDau->pdgId()==443 && abs(genDau->daughter(j)->pdgId())!=13 && genDau->daughter(j)->pdgId()!=22)wrong=true;
1446 if(abs(genDau->pdgId())==3122 && abs(genDau->daughter(j)->pdgId())!=211 && abs(genDau->daughter(j)->pdgId())!=2212 && genDau->daughter(j)->pdgId()!=22)wrong = true;
1447 }//end of loop on daughters of daughters
1448 if(genDau->pdgId()!=443 && abs(genDau->pdgId())!=3122 && genDau->pdgId()!=22)wrong = true;
1449 if(imu1!=-1 && imu2!=-1 && !wrong)iJpsi = i;
1450 if(ipi!=-1 && ip!=-1 && !wrong)iLambda = i;
1451 }//end of loop on Lambda^0_b MC daughters
1452 //cout<<endl;
1453 if(iJpsi!=-1 && iLambda!=-1 && !wrong){
1454 //cout<<"setting genLambdaJpsi = 1"<<endl;
1455 genLambdaJpsi = 1;
1456
1457 //fill MC truth info
1458 trueLambdaBMass = genLambdab.p4().M();
1459
1460 //write out info from daughters
1461 const Candidate * genJpsi = genLambdab.daughter(iJpsi);
1462 const Candidate * genLambda = genLambdab.daughter(iLambda);
1463
1464 truePriVtxX = genLambdab.vx();
1465 truePriVtxY = genLambdab.vy();
1466 truePriVtxZ = genLambdab.vz();
1467
1468 trueLambdaBPx = genLambdab.px();
1469 trueLambdaBPy = genLambdab.py();
1470 trueLambdaBPz = genLambdab.pz();
1471
1472 trueLambdaBDecayVtxX = genJpsi->vx();
1473 trueLambdaBDecayVtxY = genJpsi->vy();
1474 trueLambdaBDecayVtxZ = genJpsi->vz();
1475
1476 trueJpsiPx = genJpsi->px();
1477 trueJpsiPy = genJpsi->py();
1478 trueJpsiPz = genJpsi->pz();
1479
1480 for(uint j=0; j<genJpsi->numberOfDaughters(); ++j){
1481 if(genJpsi->daughter(j)->pdgId()==13){//13 is a mu-
1482 trueMumPx = genJpsi->daughter(j)->px();
1483 trueMumPy = genJpsi->daughter(j)->py();
1484 trueMumPz = genJpsi->daughter(j)->pz();
1485
1486 trueJpsiDecayVtxX = genJpsi->daughter(j)->vx();
1487 trueJpsiDecayVtxY = genJpsi->daughter(j)->vy();
1488 trueJpsiDecayVtxZ = genJpsi->daughter(j)->vz();
1489 }
1490 if(genJpsi->daughter(j)->pdgId()==-13){//-13 is a mu+
1491 trueMupPx = genJpsi->daughter(j)->px();
1492 trueMupPy = genJpsi->daughter(j)->py();
1493 trueMupPz = genJpsi->daughter(j)->pz();
1494 }
1495 }
1496 trueLambdaPx = genLambda->px();
1497 trueLambdaPy = genLambda->py();
1498 trueLambdaPz = genLambda->pz();
1499
1500 //cout<<"starting daughters loop with genLambda id = "<<genLambda->pdgId()<<" and with "<<genLambda->numberOfDaughters()<<" daughters"<<endl;
1501 for(uint j=0; j<genLambda->numberOfDaughters(); ++j){
1502 //cout<<"genLambda daughter has id = "<<genLambda->daughter(j)->pdgId()<<" with ndau "<<genLambda->daughter(j)->numberOfDaughters()<<endl;
1503 //cout<<"no granddaughters, filling genLambda info"<<endl;
1504 if(genLambda->daughter(j)->charge()>0){
1505 trueTrackPPx = genLambda->daughter(j)->px();
1506 trueTrackPPy = genLambda->daughter(j)->py();
1507 trueTrackPPz = genLambda->daughter(j)->pz();
1508 trueLambdaDecayVtxX = genLambda->daughter(j)->vx();
1509 trueLambdaDecayVtxY = genLambda->daughter(j)->vy();
1510 trueLambdaDecayVtxZ = genLambda->daughter(j)->vz();
1511 }
1512 if(genLambda->daughter(j)->charge()<0){
1513 trueTrackMPx = genLambda->daughter(j)->px();
1514 trueTrackMPy = genLambda->daughter(j)->py();
1515 trueTrackMPz = genLambda->daughter(j)->pz();
1516 }
1517 }
1518
1519 //true decay lengths
1520 trueLxyPV = ( (trueLambdaBDecayVtxX-truePriVtxX)*trueLambdaBPx +
1521 (trueLambdaBDecayVtxY-truePriVtxY)*trueLambdaBPy ) /
1522 sqrt(trueLambdaBPx*trueLambdaBPx + trueLambdaBPy*trueLambdaBPy);
1523
1524 trueL3DPV = ( (trueLambdaBDecayVtxX-truePriVtxX)*trueLambdaBPx +
1525 (trueLambdaBDecayVtxY-truePriVtxY)*trueLambdaBPy +
1526 (trueLambdaBDecayVtxZ-truePriVtxZ)*trueLambdaBPz) /
1527 sqrt(trueLambdaBPx*trueLambdaBPx + trueLambdaBPy*trueLambdaBPy + trueLambdaBPz*trueLambdaBPz);
1528
1529 truectauLambdaBxyPV = trueLxyPV*trueLambdaBMass/sqrt(trueLambdaBPx*trueLambdaBPx + trueLambdaBPy*trueLambdaBPy );
1530 truectauLambdaB3DPV = trueL3DPV*trueLambdaBMass/sqrt(trueLambdaBPx*trueLambdaBPx + trueLambdaBPy*trueLambdaBPy + trueLambdaBPz*trueLambdaBPz);
1531
1532
1533 //determine MC truth
1534
1535 //calculate true eta and phi for all tracks
1536 float trueMupPhi = atan(trueMupPy/trueMupPx);
1537 if(trueMupPx<0 && trueMupPy<0) trueMupPhi -=pi;
1538 if(trueMupPx<0 && trueMupPy>0) trueMupPhi +=pi;
1539 float trueMupP = sqrt(trueMupPx*trueMupPx + trueMupPy*trueMupPy + trueMupPz*trueMupPz);
1540 float trueMupEta = 0.5*(log( (trueMupP + trueMupPz)/(trueMupP - trueMupPz)) );
1541
1542 float trueMumPhi = atan(trueMumPy/trueMumPx);
1543 if(trueMumPx<0 && trueMumPy<0) trueMumPhi -=pi;
1544 if(trueMumPx<0 && trueMumPy>0) trueMumPhi +=pi;
1545 float trueMumP = sqrt( trueMumPx*trueMumPx + trueMumPy*trueMumPy + trueMumPz*trueMumPz);
1546 float trueMumEta = 0.5*(log( (trueMumP + trueMumPz)/(trueMumP - trueMumPz)) );
1547
1548 float trueTrackPPhi = atan(trueTrackPPy/trueTrackPPx);
1549 if ( trueTrackPPx < 0 && trueTrackPPy < 0 ) trueTrackPPhi -= pi;
1550 if ( trueTrackPPx < 0 && trueTrackPPy > 0 ) trueTrackPPhi += pi;
1551 float trueTrackPP = sqrt( trueTrackPPx*trueTrackPPx + trueTrackPPy*trueTrackPPy + trueTrackPPz*trueTrackPPz );
1552 float trueTrackPEta = 0.5*log( (trueTrackPP + trueTrackPPz)/(trueTrackPP - trueTrackPPz) );
1553
1554 float trueTrackMPhi = atan(trueTrackMPy/trueTrackMPx);
1555 if ( trueTrackMPx < 0 && trueTrackMPy < 0 ) trueTrackMPhi -= pi;
1556 if ( trueTrackMPx < 0 && trueTrackMPy > 0 ) trueTrackMPhi += pi;
1557 float trueTrackMP = sqrt( trueTrackMPx*trueTrackMPx + trueTrackMPy*trueTrackMPy + trueTrackMPz*trueTrackMPz );
1558 float trueTrackMEta = 0.5*log( (trueTrackMP + trueTrackMPz)/(trueTrackMP - trueTrackMPz) );
1559
1560 float RcutMu = 100;
1561 float RcutTrack = 100;
1562 float RcutVtx = 100;
1563
1564 truthMatch.clear(); truthLambda.clear(); truthJpsi.clear();
1565
1566 //loop to check all LambdaB candidates found
1567 for (uint i=0; i<mupPx->size(); ++i){
1568
1569 bool istrueMup = false;
1570 bool istrueMum = false;
1571 bool istrueTrackP = false;
1572 bool istrueTrackM = false;
1573 bool istrueLambda = false;
1574 bool istrueJpsi = false;
1575 bool istrueLambdaB = false;
1576
1577 //calculate eta and phi for all tracks in B candidate
1578 float mupPhi = atan(mupPy->at(i)/mupPx->at(i));
1579 if( mupPx->at(i)<0 && mupPy->at(i)<0 ) mupPhi -=pi;
1580 if( mupPx->at(i)<0 && mupPy->at(i)>0 ) mupPhi +=pi;
1581 float mupP = sqrt( mupPx->at(i)*mupPx->at(i) + mupPy->at(i)*mupPy->at(i) + mupPz->at(i)*mupPz->at(i) );
1582 float mupEta = 0.5*log( (mupP + mupPz->at(i))/(mupP - mupPz->at(i)) );
1583
1584 float mumPhi = atan(mumPy->at(i)/mumPx->at(i));
1585 if ( mumPx->at(i) < 0 && mumPy->at(i) < 0 ) mumPhi -= pi;
1586 if ( mumPx->at(i) < 0 && mumPy->at(i) > 0 ) mumPhi += pi;
1587 float mumP = sqrt( mumPx->at(i)*mumPx->at(i) + mumPy->at(i)*mumPy->at(i) + mumPz->at(i)*mumPz->at(i) );
1588 float mumEta = 0.5*log( (mumP + mumPz->at(i))/(mumP - mumPz->at(i)) );
1589
1590 float trackPPhi = atan(trackPPy->at(i)/trackPPx->at(i));
1591 if ( trackPPx->at(i) < 0 && trackPPy->at(i) < 0 ) trackPPhi -= pi;
1592 if ( trackPPx->at(i) < 0 && trackPPy->at(i) > 0 ) trackPPhi += pi;
1593 float trackPP = sqrt( trackPPx->at(i)*trackPPx->at(i) + trackPPy->at(i)*trackPPy->at(i) + trackPPz->at(i)*trackPPz->at(i) );
1594 float trackPEta = 0.5*log( (trackPP + trackPPz->at(i))/(trackPP - trackPPz->at(i)) );
1595
1596 float trackMPhi = atan(trackMPy->at(i)/trackMPx->at(i));
1597 if ( trackMPx->at(i) < 0 && trackMPy->at(i) < 0 ) trackMPhi -= pi;
1598 if ( trackMPx->at(i) < 0 && trackMPy->at(i) > 0 ) trackMPhi += pi;
1599 float trackMP = sqrt( trackMPx->at(i)*trackMPx->at(i) + trackMPy->at(i)*trackMPy->at(i) + trackMPz->at(i)*trackMPz->at(i) );
1600 float trackMEta = 0.5*log( (trackMP + trackMPz->at(i))/(trackMP - trackMPz->at(i)) );
1601
1602
1603
1604 float deltaRmup = sqrt( (mupEta-trueMupEta)*(mupEta-trueMupEta) +
1605 (mupPhi-trueMupPhi)*(mupPhi-trueMupPhi) );
1606 if(deltaRmup<RcutMu)istrueMup = true;
1607 float deltaRmum = sqrt( (mumEta-trueMumEta)*(mumEta-trueMumEta) +
1608 (mumPhi-trueMumPhi)*(mumPhi-trueMumPhi) );
1609 if(deltaRmum<RcutMu)istrueMum = true;
1610
1611 float deltaRtrackp =sqrt( (trackPEta-trueTrackPEta)*(trackPEta-trueTrackPEta) +
1612 (trackPPhi-trueTrackPPhi)*(trackPPhi-trueTrackPPhi) );
1613 if(deltaRtrackp<RcutTrack)istrueTrackP = true;
1614
1615 float deltaRtrackm =sqrt( (trackMEta-trueTrackMEta)*(trackMEta-trueTrackMEta) +
1616 (trackMPhi-trueTrackMPhi)*(trackMPhi-trueTrackMPhi) );
1617 if(deltaRtrackm<RcutTrack)istrueTrackM = true;
1618
1619 float deltaRLambdaVtx = sqrt( (trueLambdaDecayVtxX-LambdaVtxX->at(i))*(trueLambdaDecayVtxX-LambdaVtxX->at(i)) +
1620 (trueLambdaDecayVtxY-LambdaVtxY->at(i))*(trueLambdaDecayVtxY-LambdaVtxY->at(i)) +
1621 (trueLambdaDecayVtxZ-LambdaVtxZ->at(i))*(trueLambdaDecayVtxZ-LambdaVtxZ->at(i)) );
1622
1623 dRmup->push_back(deltaRmup);
1624 dRmum->push_back(deltaRmum);
1625 dRtrackp->push_back(deltaRtrackp);
1626 dRtrackm->push_back(deltaRtrackm);
1627 dRLambdaVtx->push_back(deltaRLambdaVtx);
1628
1629 if(istrueMup && istrueMum)istrueJpsi = true;
1630 if(istrueTrackP && istrueTrackM && (deltaRLambdaVtx<RcutVtx) )istrueLambda = true;
1631 if(istrueJpsi && istrueLambda)istrueLambdaB = true;
1632
1633 if(istrueJpsi){
1634 truthJpsi.push_back(1);
1635 }else{
1636 truthJpsi.push_back(-1);
1637 }
1638 if(istrueLambda){
1639 truthLambda.push_back(1);
1640 }else{
1641 truthLambda.push_back(-1);
1642 }
1643 if(istrueLambdaB){
1644 truthMatch.push_back(1);
1645 }else{
1646 truthMatch.push_back(-1);
1647 }
1648
1649
1650 //cout<<"istrueMup = "<<istrueMup<<", istrueMum = "<<istrueMum<<", istrueTrackP = "<<istrueTrackP<<", istrueTrackM = "<<istrueTrackM<<endl;
1651 //cout<<"istrueJpsi = "<<istrueJpsi<<", istrueLambda = "<<istrueLambda<<", istrueLambdaB = "<<istrueLambdaB<<endl;
1652 //cout<<"trueLambdaBMass = "<<trueLambdaBMass<<endl;
1653
1654 }//end of loop to check the reco Lambda B with MC
1655 }//end of if the gen Lambdab is correct
1656 }//end of if Lambda^0_b
1657 }//end of loop on genParticles
1658 ///
1659 ///////////////////////////////
1660
1661 //fill truthMatch with zeros if we didn't fidn a truth match
1662 if ( truthMatch.size()==0 ) { //if no signal event has been found, fill with zeros
1663 for (uint i = 0; i<mupPx->size(); i++) {
1664 truthJpsi.push_back(0);
1665 truthLambda.push_back(0);
1666 truthMatch.push_back(0);
1667 }
1668 }
1669
1670
1671 }//end of if nLambdaB candidates > 0
1672
1673 //fill tree and clear vectors
1674 if(nLambdaB > 0){
1675 tree_->Fill();
1676 }
1677 nLambdaB = 0;
1678 run = 0;
1679 lumiBlock = 0;
1680 event = 0;
1681
1682 priVtxX->clear();
1683 priVtxY->clear();
1684 priVtxZ->clear();
1685 priVtxXE->clear();
1686 priVtxYE->clear();
1687 priVtxZE->clear();
1688 priVtxCL->clear();
1689 priVtxBSx->clear();
1690 priVtxBSy->clear();
1691 priVtxBSz->clear();
1692 priVtxBSxE->clear();
1693 priVtxBSyE->clear();
1694 priVtxBSzE->clear();
1695 priVtxBSCL->clear();
1696 beamSpotX = 0;
1697 beamSpotY = 0;
1698 beamSpotZ = 0;
1699 hlt_mu3 = 0;
1700 hlt_mu5 = 0;
1701 hlt_mu7 = 0;
1702 hlt_mu9 = 0;
1703 hlt_2mu0 = 0;
1704 hlt_2mu3 = 0;
1705 hlt_2mu0L2 = 0;
1706 hlt_2mu3JPsi = 0;
1707 hlt_BJPsiMuMu = 0;
1708 hlt_mu0trk0 = 0;
1709 hlt_mu3trk0 = 0;
1710 hlt_mu0trkmu0 = 0;
1711 hlt_mu3trkmu0 = 0;
1712 hlt_L1muOpen = 0;
1713 hlt_L12muOpen = 0;
1714 hlt_L12muOpenTight = 0;
1715 hlt_mu0_L1muOpen= 0;
1716 hlt_mu3_L1muOpen= 0;
1717 hlt_mu5_L1muOpen= 0;
1718 hlt_mu5trk0 = 0;
1719 hlt_mu0trkmu0_nocharge = 0;
1720 hlt_mu0trkmu0_OST = 0;
1721 hlt_2mu0_Qv1 = 0;
1722 hlt_2mu0_QLSv1 = 0;
1723 hlt_mu0trkmu0_OST_tightV1 = 0;
1724 hlt_mu0trkmu0_OST_tightV2 = 0;
1725 hlt_mu0trkmu0_OST_tightV3 = 0;
1726 JpsiVtxX->clear();
1727 JpsiVtxY->clear();
1728 JpsiVtxZ->clear();
1729 LambdaVtxX->clear();
1730 LambdaVtxY->clear();
1731 LambdaVtxZ->clear();
1732 JpsiVtxProb->clear();
1733 JpsiVtxCL->clear();
1734 LambdaVtxProb->clear();
1735 LambdaVtxCL->clear();
1736 LambdaBMass->clear();
1737 LambdaBPx->clear();
1738 LambdaBPy->clear();
1739 LambdaBPz->clear();
1740 JpsiPx->clear();
1741 JpsiPy->clear();
1742 JpsiPz->clear();
1743 LambdaPx->clear();
1744 LambdaPy->clear();
1745 LambdaPz->clear();
1746 JpsiMass->clear();
1747 LambdaMass->clear();
1748 mupPx->clear();
1749 mupPy->clear();
1750 mupPz->clear();
1751 mumPx->clear();
1752 mumPy->clear();
1753 mumPz->clear();
1754 trackPPx->clear();
1755 trackPPy->clear();
1756 trackPPz->clear();
1757 trackMPx->clear();
1758 trackMPy->clear();
1759 trackMPz->clear();
1760 JpsiLambdaVtxDistXY->clear();
1761 JpsiLambdaVtxDist3D->clear();
1762 openingAngleLambdaB->clear();
1763 openingAngleLambda->clear();
1764 LxyPV->clear();
1765 LxyBS->clear();
1766 L3DPV->clear();
1767 L3DBS->clear();
1768 ctauLambdaBxyPV->clear();
1769 ctauLambdaBxyBS->clear();
1770 ctauLambdaB3DPV->clear();
1771 ctauLambdaB3DBS->clear();
1772 dRmup->clear();
1773 dRmum->clear();
1774 dRtrackp->clear();
1775 dRtrackm->clear();
1776 dRLambdaVtx->clear();
1777 mupChi2->clear();
1778 mumChi2->clear();
1779 mupNDOF->clear();
1780 mumNDOF->clear();
1781 mupD0->clear();
1782 mumD0->clear();
1783 mupNvalidHits->clear();
1784 mumNvalidHits->clear();
1785 mupIsGlobal->clear();
1786 mumIsGlobal->clear();
1787 mupIsTracker->clear();
1788 mumIsTracker->clear();
1789 ctauLambdaXY->clear();
1790
1791 JpsiVtxX_afterKinFit->clear();
1792 JpsiVtxY_afterKinFit->clear();
1793 JpsiVtxZ_afterKinFit->clear();
1794 LambdaVtxX_afterKinFit->clear();
1795 LambdaVtxY_afterKinFit->clear();
1796 LambdaVtxZ_afterKinFit->clear();
1797 JpsiVtxProb_afterKinFit->clear();
1798 JpsiVtxCL_afterKinFit->clear();
1799 LambdaVtxProb_afterKinFit->clear();
1800 LambdaVtxCL_afterKinFit->clear();
1801 LambdaBMass_afterKinFit->clear();
1802 LambdaBPx_afterKinFit->clear();
1803 LambdaBPy_afterKinFit->clear();
1804 LambdaBPz_afterKinFit->clear();
1805 JpsiPx_afterKinFit->clear();
1806 JpsiPy_afterKinFit->clear();
1807 JpsiPz_afterKinFit->clear();
1808 LambdaPx_afterKinFit->clear();
1809 LambdaPy_afterKinFit->clear();
1810 LambdaPz_afterKinFit->clear();
1811 JpsiMass_afterKinFit->clear();
1812 LambdaMass_afterKinFit->clear();
1813 mupPx_afterKinFit->clear();
1814 mupPy_afterKinFit->clear();
1815 mupPz_afterKinFit->clear();
1816 mumPx_afterKinFit->clear();
1817 mumPy_afterKinFit->clear();
1818 mumPz_afterKinFit->clear();
1819 trackPPx_afterKinFit->clear();
1820 trackPPy_afterKinFit->clear();
1821 trackPPz_afterKinFit->clear();
1822 trackMPx_afterKinFit->clear();
1823 trackMPy_afterKinFit->clear();
1824 trackMPz_afterKinFit->clear();
1825 JpsiLambdaVtxDistXY_afterKinFit->clear();
1826 JpsiLambdaVtxDist3D_afterKinFit->clear();
1827 openingAngleLambdaB_afterKinFit->clear();
1828 openingAngleLambda_afterKinFit->clear();
1829 LxyPV_afterKinFit->clear();
1830 LxyBS_afterKinFit->clear();
1831 L3DPV_afterKinFit->clear();
1832 L3DBS_afterKinFit->clear();
1833 ctauLambdaBxyPV_afterKinFit->clear();
1834 ctauLambdaBxyBS_afterKinFit->clear();
1835 ctauLambdaB3DPV_afterKinFit->clear();
1836 ctauLambdaB3DBS_afterKinFit->clear();
1837 dRmup_afterKinFit->clear();
1838 dRmum_afterKinFit->clear();
1839 dRtrackp_afterKinFit->clear();
1840 dRtrackm_afterKinFit->clear();
1841 dRLambdaVtx_afterKinFit->clear();
1842 mupChi2_afterKinFit->clear();
1843 mumChi2_afterKinFit->clear();
1844 mupNDOF_afterKinFit->clear();
1845 mumNDOF_afterKinFit->clear();
1846 mupD0_afterKinFit->clear();
1847 mumD0_afterKinFit->clear();
1848 mupNvalidHits_afterKinFit->clear();
1849 mumNvalidHits_afterKinFit->clear();
1850 mupIsGlobal_afterKinFit->clear();
1851 mumIsGlobal_afterKinFit->clear();
1852 mupIsTracker_afterKinFit->clear();
1853 mumIsTracker_afterKinFit->clear();
1854 ctauLambdaXY_afterKinFit->clear();
1855 LambdaBX_afterKinFit->clear();
1856 LambdaBY_afterKinFit->clear();
1857 LambdaBZ_afterKinFit->clear();
1858 LambdaBVtxX_afterKinFit->clear();
1859 LambdaBVtxY_afterKinFit->clear();
1860 LambdaBVtxZ_afterKinFit->clear();
1861 LambdaBNDOF_afterKinFit->clear();
1862 LambdaBchiSquared_afterKinFit->clear();
1863 LambdaBVtxProb_afterKinFit->clear();
1864 LambdaBPt_afterKinFit->clear();
1865 LambdaBEta_afterKinFit->clear();
1866 LambdaBPhi_afterKinFit->clear();
1867
1868 mupIsL1Matched->clear();
1869 mumIsL1Matched->clear();
1870 isHLTMatched->clear();
1871
1872 genLambdaJpsi = 0;
1873
1874
1875
1876 //MC truth
1877 truePriVtxX = 0;
1878 truePriVtxY = 0;
1879 truePriVtxZ = 0;
1880 trueLambdaBPx = 0;
1881 trueLambdaBPy = 0;
1882 trueLambdaBPz = 0;
1883 trueLambdaBDecayVtxX = 0;
1884 trueLambdaBDecayVtxY = 0;
1885 trueLambdaBDecayVtxZ = 0;
1886 trueJpsiPx = 0;
1887 trueJpsiPy = 0;
1888 trueJpsiPz = 0;
1889 trueJpsiDecayVtxX = 0;
1890 trueJpsiDecayVtxY = 0;
1891 trueJpsiDecayVtxZ = 0;
1892 trueMumPx = 0;
1893 trueMumPy = 0;
1894 trueMumPz = 0;
1895 trueMupPx = 0;
1896 trueMupPy = 0;
1897 trueMupPz = 0;
1898 trueLambdaPx = 0;
1899 trueLambdaPy = 0;
1900 trueLambdaPz = 0;
1901 trueLambdaDecayVtxX = 0;
1902 trueLambdaDecayVtxY = 0;
1903 trueLambdaDecayVtxZ = 0;
1904 trueTrackPPx = 0;
1905 trueTrackPPy = 0;
1906 trueTrackPPz = 0;
1907 trueTrackMPx = 0;
1908 trueTrackMPy = 0;
1909 trueTrackMPz = 0;
1910 trueLxyPV = 0;
1911 trueLxyBS = 0;
1912 trueL3DPV = 0;
1913 trueL3DBS = 0;
1914 truectauLambdaBxyPV = 0;
1915 truectauLambdaB3DPV = 0;
1916 trueLambdaBMass = 0;
1917
1918 truthMatch.clear();
1919 truthLambda.clear();
1920 truthJpsi.clear();
1921
1922
1923
1924
1925
1926
1927
1928 #ifdef THIS_IS_AN_EVENT_EXAMPLE
1929 Handle<ExampleData> pIn;
1930 iEvent.getByLabel("example",pIn);
1931 #endif
1932
1933 #ifdef THIS_IS_AN_EVENTSETUP_EXAMPLE
1934 ESHandle<SetupData> pSetup;
1935 iSetup.get<SetupRecord>().get(pSetup);
1936 #endif
1937 }
1938
1939
1940 // ------------ method called once each job just before starting event loop ------------
1941 void
1942 LambdaAnalyzer::beginJob()
1943 {
1944 tmpNvtx = 0;
1945 cout<<"begin job"<<endl;
1946 edm::Service<TFileService> fs;
1947 tree_ = fs->make<TTree>("ntpl","LambdaB ntuple");
1948
1949 //branches for reco particles
1950 tree_->Branch("nLambdaB",&nLambdaB,"nBLambda/i");
1951 tree_->Branch("run",&run,"run/f");
1952 tree_->Branch("lumiBlock",&lumiBlock,"lumiBlock/f");
1953 tree_->Branch("event",&event,"event/f");
1954
1955 tree_->Branch("priVtxX",&priVtxX);
1956 tree_->Branch("priVtxY",&priVtxY);
1957 tree_->Branch("priVtxZ",&priVtxZ);
1958 tree_->Branch("priVtxXE",&priVtxXE);
1959 tree_->Branch("priVtxYE",&priVtxYE);
1960 tree_->Branch("priVtxZE",&priVtxZE);
1961 tree_->Branch("priVtxCL",&priVtxCL);
1962 tree_->Branch("priVtxBSx",&priVtxBSx);
1963 tree_->Branch("priVtxBSy",&priVtxBSy);
1964 tree_->Branch("priVtxBSz",&priVtxBSz);
1965 tree_->Branch("priVtxBSxE",&priVtxBSxE);
1966 tree_->Branch("priVtxBSyE",&priVtxBSyE);
1967 tree_->Branch("priVtxBSzE",&priVtxBSzE);
1968 tree_->Branch("priVtxBSCL",&priVtxBSCL);
1969 tree_->Branch("beamSpotX",&beamSpotX, "beamSpotX/f");
1970 tree_->Branch("beamSpotY",&beamSpotY, "beamSpotY/f");
1971 tree_->Branch("beamSpotZ",&beamSpotZ, "beamSpotZ/f");
1972 tree_->Branch("hlt_mu3",&hlt_mu3,"hlt_mu3/i");
1973 tree_->Branch("hlt_mu5",&hlt_mu5,"hlt_mu5/i");
1974 tree_->Branch("hlt_mu7",&hlt_mu7,"hlt_mu7/i");
1975 tree_->Branch("hlt_mu9",&hlt_mu9,"hlt_mu9/i");
1976 tree_->Branch("hlt_2mu0",&hlt_2mu0,"hlt_2mu0/i");
1977 tree_->Branch("hlt_2mu3",&hlt_2mu3,"hlt_2mu3/i");
1978 tree_->Branch("hlt_2mu0L2",&hlt_2mu0L2,"hlt_2mu0L2/i");
1979 tree_->Branch("hlt_2mu3JPsi",&hlt_2mu3JPsi,"hlt_2mu3JPsi/i");
1980 tree_->Branch("hlt_BJPsiMuMu",&hlt_BJPsiMuMu,"hlt_BJPsiMuMu/i");
1981 tree_->Branch("hlt_mu0trk0",&hlt_mu0trk0,"hlt_mu0trk0/i");
1982 tree_->Branch("hlt_mu3trk0",&hlt_mu3trk0,"hlt_mu3trk0/i");
1983 tree_->Branch("hlt_mu0trkmu0",&hlt_mu0trkmu0,"hlt_mu0trkmu0/i");
1984 tree_->Branch("hlt_mu3trkmu0",&hlt_mu3trkmu0,"hlt_mu3trkmu0/i");
1985 tree_->Branch("hlt_L1muOpen",&hlt_L1muOpen,"hlt_L1muOpen/i");
1986 tree_->Branch("hlt_L12muOpen",&hlt_L12muOpen,"hlt_L12muOpen/i");
1987 tree_->Branch("hlt_L12muOpenTight",&hlt_L12muOpenTight,"hlt_L12muOpenTight/i");
1988 tree_->Branch("hlt_mu0_L1muOpen",&hlt_mu0_L1muOpen,"hlt_mu0_L1muOpen/i");
1989 tree_->Branch("hlt_mu3_L1muOpen",&hlt_mu3_L1muOpen,"hlt_mu3_L1muOpen/i");
1990 tree_->Branch("hlt_mu5_L1muOpen",&hlt_mu5_L1muOpen,"hlt_mu5_L1muOpen/i");
1991 tree_->Branch("hlt_mu5trk0",&hlt_mu5trk0,"hlt_mu5trk0/i");
1992 tree_->Branch("hlt_mu0trkmu0_nocharge",&hlt_mu0trkmu0_nocharge,"hlt_mu0trkmu0_nocharge/i");
1993 tree_->Branch("hlt_mu0trkmu0_OST",&hlt_mu0trkmu0_OST,"hlt_mu0trkmu0_OST/i");
1994 tree_->Branch("hlt_2mu0_Qv1",&hlt_2mu0_Qv1,"hlt_2mu0_Qv1/i");
1995 tree_->Branch("hlt_2mu0_QLSv1",&hlt_2mu0_QLSv1,"hlt_2mu0_QLSv1/i");
1996 tree_->Branch("hlt_mu0trkmu0_OST_tightV1",&hlt_mu0trkmu0_OST_tightV1,"hlt_mu0trkmu0_OST_tightV1/i");
1997 tree_->Branch("hlt_mu0trkmu0_OST_tightV2",&hlt_mu0trkmu0_OST_tightV2,"hlt_mu0trkmu0_OST_tightV2/i");
1998 tree_->Branch("hlt_mu0trkmu0_OST_tightV3",&hlt_mu0trkmu0_OST_tightV3,"hlt_mu0trkmu0_OST_tightV3/i");
1999 tree_->Branch("JpsiVtxX",&JpsiVtxX);
2000 tree_->Branch("JpsiVtxY",&JpsiVtxY);
2001 tree_->Branch("JpsiVtxZ",&JpsiVtxZ);
2002 tree_->Branch("LambdaVtxX",&LambdaVtxX);
2003 tree_->Branch("LambdaVtxY",&LambdaVtxY);
2004 tree_->Branch("LambdaVtxZ",&LambdaVtxZ);
2005 tree_->Branch("JpsiVtxProb",&JpsiVtxProb);
2006 tree_->Branch("JpsiVtxCL",&JpsiVtxCL);
2007 tree_->Branch("LambdaVtxProb",&LambdaVtxProb);
2008 tree_->Branch("LambdaVtxCL",&LambdaVtxCL);
2009 tree_->Branch("LambdaBMass",&LambdaBMass);
2010 tree_->Branch("LambdaBPx",&LambdaBPx);
2011 tree_->Branch("LambdaBPy",&LambdaBPy);
2012 tree_->Branch("LambdaBPz",&LambdaBPz);
2013 tree_->Branch("JpsiPx",&JpsiPx);
2014 tree_->Branch("JpsiPy",&JpsiPy);
2015 tree_->Branch("JpsiPz",&JpsiPz);
2016 tree_->Branch("LambdaPx",&LambdaPx);
2017 tree_->Branch("LambdaPy",&LambdaPy);
2018 tree_->Branch("LambdaPz",&LambdaPz);
2019 tree_->Branch("JpsiMass",&JpsiMass);
2020 tree_->Branch("LambdaMass",&LambdaMass);
2021 tree_->Branch("mupPx",&mupPx);
2022 tree_->Branch("mupPy",&mupPy);
2023 tree_->Branch("mupPz",&mupPz);
2024 tree_->Branch("mumPx",&mumPx);
2025 tree_->Branch("mumPy",&mumPy);
2026 tree_->Branch("mumPz",&mumPz);
2027 tree_->Branch("trackPPx",&trackPPx);
2028 tree_->Branch("trackPPy",&trackPPy);
2029 tree_->Branch("trackPPz",&trackPPz);
2030 tree_->Branch("trackMPx",&trackMPx);
2031 tree_->Branch("trackMPy",&trackMPy);
2032 tree_->Branch("trackMPz",&trackMPz);
2033 tree_->Branch("JpsiLambdaVtxDistXY",&JpsiLambdaVtxDistXY);
2034 tree_->Branch("JpsiLambdaVtxDist3D",&JpsiLambdaVtxDist3D);
2035 tree_->Branch("openingAngleLambdaB",&openingAngleLambdaB);
2036 tree_->Branch("openingAngleLambda",&openingAngleLambda);
2037 tree_->Branch("LxyPV",&LxyPV);
2038 tree_->Branch("LxyBS",&LxyBS);
2039 tree_->Branch("L3DPV",&L3DPV);
2040 tree_->Branch("L3DBS",&L3DBS);
2041 tree_->Branch("ctauLambdaBxyPV",&ctauLambdaBxyPV);
2042 tree_->Branch("ctauLambdaBxyBS",&ctauLambdaBxyBS);
2043 tree_->Branch("ctauLambdaB3DPV",&ctauLambdaB3DPV);
2044 tree_->Branch("ctauLambdaB3DBS",&ctauLambdaB3DBS);
2045 tree_->Branch("dRmup",&dRmup);
2046 tree_->Branch("dRmum",&dRmum);
2047 tree_->Branch("dRtrackp",&dRtrackp);
2048 tree_->Branch("dRtrackm",&dRtrackm);
2049 tree_->Branch("dRLambdaVtx",&dRLambdaVtx);
2050 tree_->Branch("mupChi2",&mupChi2);
2051 tree_->Branch("mumChi2",&mumChi2);
2052 tree_->Branch("mupNDOF",&mupNDOF);
2053 tree_->Branch("mumNDOF",&mumNDOF);
2054 tree_->Branch("mupD0",&mupD0);
2055 tree_->Branch("mumD0",&mumD0);
2056 tree_->Branch("mupNvalidHits",&mupNvalidHits);
2057 tree_->Branch("mumNvalidHits",&mumNvalidHits);
2058 tree_->Branch("mupIsGlobal",&mupIsGlobal);
2059 tree_->Branch("mumIsGlobal",&mumIsGlobal);
2060 tree_->Branch("mupIsTracker",&mupIsTracker);
2061 tree_->Branch("mumIsTracker",&mumIsTracker);
2062 tree_->Branch("ctauLambdaXY",&ctauLambdaXY);
2063
2064 //after kinematic fit
2065 tree_->Branch("JpsiVtxX_afterKinFit",&JpsiVtxX_afterKinFit);
2066 tree_->Branch("JpsiVtxY_afterKinFit",&JpsiVtxY_afterKinFit);
2067 tree_->Branch("JpsiVtxZ_afterKinFit",&JpsiVtxZ_afterKinFit);
2068 tree_->Branch("LambdaVtxX_afterKinFit",&LambdaVtxX_afterKinFit);
2069 tree_->Branch("LambdaVtxY_afterKinFit",&LambdaVtxY_afterKinFit);
2070 tree_->Branch("LambdaVtxZ_afterKinFit",&LambdaVtxZ_afterKinFit);
2071 tree_->Branch("JpsiVtxProb_afterKinFit",&JpsiVtxProb_afterKinFit);
2072 tree_->Branch("JpsiVtxCL_afterKinFit",&JpsiVtxCL_afterKinFit);
2073 tree_->Branch("LambdaVtxProb_afterKinFit",&LambdaVtxProb_afterKinFit);
2074 tree_->Branch("LambdaVtxCL_afterKinFit",&LambdaVtxCL_afterKinFit);
2075 tree_->Branch("LambdaBMass_afterKinFit",&LambdaBMass_afterKinFit);
2076 tree_->Branch("LambdaBPx_afterKinFit",&LambdaBPx_afterKinFit);
2077 tree_->Branch("LambdaBPy_afterKinFit",&LambdaBPy_afterKinFit);
2078 tree_->Branch("LambdaBPz_afterKinFit",&LambdaBPz_afterKinFit);
2079 tree_->Branch("JpsiPx_afterKinFit",&JpsiPx_afterKinFit);
2080 tree_->Branch("JpsiPy_afterKinFit",&JpsiPy_afterKinFit);
2081 tree_->Branch("JpsiPz_afterKinFit",&JpsiPz_afterKinFit);
2082 tree_->Branch("LambdaPx_afterKinFit",&LambdaPx_afterKinFit);
2083 tree_->Branch("LambdaPy_afterKinFit",&LambdaPy_afterKinFit);
2084 tree_->Branch("LambdaPz_afterKinFit",&LambdaPz_afterKinFit);
2085 tree_->Branch("JpsiMass_afterKinFit",&JpsiMass_afterKinFit);
2086 tree_->Branch("LambdaMass_afterKinFit",&LambdaMass_afterKinFit);
2087 tree_->Branch("mupPx_afterKinFit",&mupPx_afterKinFit);
2088 tree_->Branch("mupPy_afterKinFit",&mupPy_afterKinFit);
2089 tree_->Branch("mupPz_afterKinFit",&mupPz_afterKinFit);
2090 tree_->Branch("mumPx_afterKinFit",&mumPx_afterKinFit);
2091 tree_->Branch("mumPy_afterKinFit",&mumPy_afterKinFit);
2092 tree_->Branch("mumPz_afterKinFit",&mumPz_afterKinFit);
2093 tree_->Branch("trackPPx_afterKinFit",&trackPPx_afterKinFit);
2094 tree_->Branch("trackPPy_afterKinFit",&trackPPy_afterKinFit);
2095 tree_->Branch("trackPPz_afterKinFit",&trackPPz_afterKinFit);
2096 tree_->Branch("trackMPx_afterKinFit",&trackMPx_afterKinFit);
2097 tree_->Branch("trackMPy_afterKinFit",&trackMPy_afterKinFit);
2098 tree_->Branch("trackMPz_afterKinFit",&trackMPz_afterKinFit);
2099 tree_->Branch("JpsiLambdaVtxDistXY_afterKinFit",&JpsiLambdaVtxDistXY_afterKinFit);
2100 tree_->Branch("JpsiLambdaVtxDist3D_afterKinFit",&JpsiLambdaVtxDist3D_afterKinFit);
2101 tree_->Branch("openingAngleLambdaB_afterKinFit",&openingAngleLambdaB_afterKinFit);
2102 tree_->Branch("openingAngleLambda_afterKinFit",&openingAngleLambda_afterKinFit);
2103 tree_->Branch("LxyPV_afterKinFit",&LxyPV_afterKinFit);
2104 tree_->Branch("LxyBS_afterKinFit",&LxyBS_afterKinFit);
2105 tree_->Branch("L3DPV_afterKinFit",&L3DPV_afterKinFit);
2106 tree_->Branch("L3DBS_afterKinFit",&L3DBS_afterKinFit);
2107 tree_->Branch("ctauLambdaBxyPV_afterKinFit",&ctauLambdaBxyPV_afterKinFit);
2108 tree_->Branch("ctauLambdaBxyBS_afterKinFit",&ctauLambdaBxyBS_afterKinFit);
2109 tree_->Branch("ctauLambdaB3DPV_afterKinFit",&ctauLambdaB3DPV_afterKinFit);
2110 tree_->Branch("ctauLambdaB3DBS_afterKinFit",&ctauLambdaB3DBS_afterKinFit);
2111 tree_->Branch("dRmup_afterKinFit",&dRmup_afterKinFit);
2112 tree_->Branch("dRmum_afterKinFit",&dRmum_afterKinFit);
2113 tree_->Branch("dRtrackp_afterKinFit",&dRtrackp_afterKinFit);
2114 tree_->Branch("dRtrackm_afterKinFit",&dRtrackm_afterKinFit);
2115 tree_->Branch("dRLambdaVtx_afterKinFit",&dRLambdaVtx_afterKinFit);
2116 tree_->Branch("mupChi2_afterKinFit",&mupChi2_afterKinFit);
2117 tree_->Branch("mumChi2_afterKinFit",&mumChi2_afterKinFit);
2118 tree_->Branch("mupNDOF_afterKinFit",&mupNDOF_afterKinFit);
2119 tree_->Branch("mumNDOF_afterKinFit",&mumNDOF_afterKinFit);
2120 tree_->Branch("mupD0_afterKinFit",&mupD0_afterKinFit);
2121 tree_->Branch("mumD0_afterKinFit",&mumD0_afterKinFit);
2122 tree_->Branch("mupNvalidHits_afterKinFit",&mupNvalidHits_afterKinFit);
2123 tree_->Branch("mumNvalidHits_afterKinFit",&mumNvalidHits_afterKinFit);
2124 tree_->Branch("mupIsGlobal_afterKinFit",&mupIsGlobal_afterKinFit);
2125 tree_->Branch("mumIsGlobal_afterKinFit",&mumIsGlobal_afterKinFit);
2126 tree_->Branch("mupIsTracker_afterKinFit",&mupIsTracker_afterKinFit);
2127 tree_->Branch("mumIsTracker_afterKinFit",&mumIsTracker_afterKinFit);
2128 tree_->Branch("ctauLambdaXY_afterKinFit",&ctauLambdaXY_afterKinFit);
2129 tree_->Branch("LambdaBX_afterKinFit",&LambdaBX_afterKinFit);
2130 tree_->Branch("LambdaBY_afterKinFit",&LambdaBY_afterKinFit);
2131 tree_->Branch("LambdaBZ_afterKinFit",&LambdaBZ_afterKinFit);
2132 tree_->Branch("LambdaBVtxX_afterKinFit",&LambdaBVtxX_afterKinFit);
2133 tree_->Branch("LambdaBVtxY_afterKinFit",&LambdaBVtxY_afterKinFit);
2134 tree_->Branch("LambdaBVtxZ_afterKinFit",&LambdaBVtxZ_afterKinFit);
2135 tree_->Branch("LambdaBNDOF_afterKinFit",&LambdaBNDOF_afterKinFit);
2136 tree_->Branch("LambdaBchiSquared_afterKinFit",&LambdaBchiSquared_afterKinFit);
2137 tree_->Branch("LambdaBVtxProb_afterKinFit",&LambdaBVtxProb_afterKinFit);
2138 tree_->Branch("LambdaBPt_afterKinFit",&LambdaBPt_afterKinFit);
2139 tree_->Branch("LambdaBEta_afterKinFit",&LambdaBEta_afterKinFit);
2140 tree_->Branch("LambdaBPhi_afterKinFit",&LambdaBPhi_afterKinFit);
2141 tree_->Branch("mupIsL1Matched",&mupIsL1Matched);
2142 tree_->Branch("mumIsL1Matched",&mumIsL1Matched);
2143 tree_->Branch("isHLTMatched",&isHLTMatched);
2144
2145 tree_->Branch("genLambdaJpsi", &genLambdaJpsi, "genLambdaJpsi/I");
2146
2147 // branches for MC truth
2148 tree_->Branch("truePriVtxX", &truePriVtxX, "truePriVtxX/f");
2149 tree_->Branch("truePriVtxY", &truePriVtxY, "truePriVtxY/f");
2150 tree_->Branch("truePriVtxZ", &truePriVtxZ, "truePriVtxZ/f");
2151 tree_->Branch("trueLambdaBPx",&trueLambdaBPx, "trueLambdaBPx/f");
2152 tree_->Branch("trueLambdaBPy",&trueLambdaBPy, "trueLambdaBPy/f");
2153 tree_->Branch("trueLambdaBPz",&trueLambdaBPz, "trueLambdaBPz/f");
2154 tree_->Branch("trueLambdaBDecayVtxX",&trueLambdaBDecayVtxX, "trueLambdaBDecayVtxX/f");
2155 tree_->Branch("trueLambdaBDecayVtxY",&trueLambdaBDecayVtxY, "trueLambdaBDecayVtxY/f");
2156 tree_->Branch("trueLambdaBDecayVtxZ",&trueLambdaBDecayVtxZ, "trueLambdaBDecayVtxZ/f");
2157 tree_->Branch("trueJpsiPx",&trueJpsiPx, "trueJpsiPx/f");
2158 tree_->Branch("trueJpsiPy",&trueJpsiPy, "trueJpsiPy/f");
2159 tree_->Branch("trueJpsiPz",&trueJpsiPz, "trueJpsiPz/f");
2160 tree_->Branch("trueJpsiDecayVtxX",&trueJpsiDecayVtxX, "trueJpsiDecayVtxX/f");
2161 tree_->Branch("trueJpsiDecayVtxY",&trueJpsiDecayVtxY, "trueJpsiDecayVtxY/f");
2162 tree_->Branch("trueJpsiDecayVtxZ",&trueJpsiDecayVtxZ, "trueJpsiDecayVtxZ/f");
2163 tree_->Branch("trueMumPx",&trueMumPx, "trueMumPx/f");
2164 tree_->Branch("trueMumPy",&trueMumPy, "trueMumPy/f");
2165 tree_->Branch("trueMumPz",&trueMumPz, "trueMumPz/f");
2166 tree_->Branch("trueMupPx",&trueMupPx, "trueMupPx/f");
2167 tree_->Branch("trueMupPy",&trueMupPy, "trueMupPy/f");
2168 tree_->Branch("trueMupPz",&trueMupPz, "trueMupPz/f");
2169 tree_->Branch("trueLambdaPx",&trueLambdaPx, "trueLambdaPx/f");
2170 tree_->Branch("trueLambdaPy",&trueLambdaPy, "trueLambdaPy/f");
2171 tree_->Branch("trueLambdaPz",&trueLambdaPz, "trueLambdaPz/f");
2172 tree_->Branch("trueLambdaDecayVtxX",&trueLambdaDecayVtxX, "trueLambdaDecayVtxX/f");
2173 tree_->Branch("trueLambdaDecayVtxY",&trueLambdaDecayVtxY, "trueLambdaDecayVtxY/f");
2174 tree_->Branch("trueLambdaDecayVtxZ",&trueLambdaDecayVtxZ, "trueLambdaDecayVtxZ/f");
2175 tree_->Branch("trueTrackPPx",&trueTrackPPx, "trueTrackPPx/f");
2176 tree_->Branch("trueTrackPPy",&trueTrackPPy, "trueTrackPPy/f");
2177 tree_->Branch("trueTrackPPz",&trueTrackPPz, "trueTrackPPz/f");
2178 tree_->Branch("trueTrackMPx",&trueTrackMPx, "trueTrackMPx/f");
2179 tree_->Branch("trueTrackMPy",&trueTrackMPy, "trueTrackMPy/f");
2180 tree_->Branch("trueTrackMPz",&trueTrackMPz, "trueTrackMPz/f");
2181 tree_->Branch("trueLxyPV",&trueLxyPV, "trueLxyPV/f");
2182 tree_->Branch("trueLxyBS",&trueLxyBS, "trueLxyBS/f");
2183 tree_->Branch("trueL3DPV",&trueL3DPV, "trueL3DPV/f");
2184 tree_->Branch("trueL3DBS",&trueL3DBS, "trueL3DBS/f");
2185 tree_->Branch("truectauLambdaBxyPV",&truectauLambdaBxyPV, "truectauLambdaBxyPV/f");
2186 tree_->Branch("truectauLambdaB3DPV",&truectauLambdaB3DPV, "truectauLambdaB3DPV/f");
2187 tree_->Branch("trueLambdaBMass",&trueLambdaBMass, "trueLambdaBMass/f");
2188
2189 tree_->Branch("truthMatch",&truthMatch);
2190 tree_->Branch("truthLambda",&truthLambda);
2191 tree_->Branch("truthJpsi",&truthJpsi);
2192
2193
2194
2195 }
2196
2197 // ------------ method called once each job just after ending the event loop ------------
2198 void
2199 LambdaAnalyzer::endJob() {
2200 cout<<"number of valid vertices: "<<tmpNvtx<<endl;
2201 // tree_->GetDirectory()->cd();
2202 // tree_->Write();
2203 }
2204
2205
2206
2207 bool LambdaAnalyzer::pass_HLT_Mu0_TkMu0_Jpsi(const edm::Event& iEvent){
2208
2209 edm::Handle<trigger::TriggerEvent> hltevent;
2210 iEvent.getByLabel(triggerEventTag_, hltevent);
2211 const trigger::TriggerObjectCollection&
2212 hltobjects(hltevent->getObjects());
2213
2214 trigger::size_type index(0);
2215 edm::InputTag it_track("hltMu0TrackJpsiTrackMassFiltered::REDIGI36X");
2216
2217 index=hltevent->filterIndex(it_track);
2218 if (index<hltevent->sizeFilters()) {
2219 const trigger::Keys& KEYS (hltevent->filterKeys(index));
2220 const trigger::size_type nTO(KEYS.size());
2221 for (trigger::size_type i=0; i!=nTO; ++i) {
2222 const trigger::TriggerObject& TO1( hltobjects.at(KEYS[i]) );
2223 double pt1 = TO1.pt(); double eta1 = TO1.eta(); double phi1 = TO1.phi();
2224 i++;
2225 const trigger::TriggerObject& TO2( hltobjects.at(KEYS[i]) );
2226 double pt2 = TO2.pt(); double eta2 = TO2.eta(); double phi2 = TO2.phi();
2227 double massJpsi = sqrt(2*pt1*pt2*(TMath::CosH(eta1)*TMath::CosH(eta2)-TMath::SinH(eta1)*TMath::SinH(eta2)-TMath::Sin(phi1)*TMath::Sin(phi2)-TMath::Cos(phi1)*TMath::Cos(phi2)));
2228 if (massJpsi>2.5 && massJpsi<4.1) return true;
2229 }
2230 }
2231
2232 return false;
2233 }
2234
2235 bool LambdaAnalyzer::pass_HLT_Mu0_TkMu0_OST_Jpsi(const edm::Event& iEvent) {
2236 return pass_HLT_Mu0_TkMu0_Jpsi(iEvent);
2237 }
2238
2239 bool LambdaAnalyzer::pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v1(const edm::Event& iEvent)
2240 {
2241
2242 edm::Handle<reco::MuonCollection> theMuons;
2243 iEvent.getByLabel( edm::InputTag("muons"), theMuons);
2244
2245 edm::Handle<trigger::TriggerEvent> hltevent;
2246 iEvent.getByLabel(triggerEventTag_, hltevent);
2247 const trigger::TriggerObjectCollection&
2248 hltobjects(hltevent->getObjects());
2249
2250 trigger::size_type index(0);
2251 edm::InputTag it_track("hltMu0TrackJpsiTrackMassFiltered::REDIGI36X");
2252
2253 index=hltevent->filterIndex(it_track);
2254 if (index<hltevent->sizeFilters()) {
2255 const trigger::Keys& KEYS (hltevent->filterKeys(index));
2256 const trigger::size_type nTO(KEYS.size());
2257 for (trigger::size_type i=0; i!=nTO; ++i) {
2258 const trigger::TriggerObject& TO1( hltobjects.at(KEYS[i]) );
2259 double pt1 = TO1.pt(); double eta1 = TO1.eta(); double phi1 = TO1.phi();
2260 i++;
2261 const trigger::TriggerObject& TO2( hltobjects.at(KEYS[i]) );
2262 double pt2 = TO2.pt(); double eta2 = TO2.eta(); double phi2 = TO2.phi();
2263 double minDR1 = 999.;
2264 int q1 = 0;
2265 for(reco::MuonCollection::const_iterator iMuon=theMuons->begin(); iMuon!=theMuons->end(); ++iMuon){
2266 double dR = deltaR(eta1,phi1,iMuon->eta(),iMuon->phi());
2267 if(dR<minDR1){
2268 minDR1 = dR;
2269 q1 = iMuon->charge();
2270 }
2271 }
2272 double minDR2 = 999.;
2273 int q2 = 0;
2274 for(reco::MuonCollection::const_iterator iMuon=theMuons->begin(); iMuon!=theMuons->end(); ++iMuon){
2275 double dR = deltaR(eta2,phi2,iMuon->eta(),iMuon->phi());
2276 if(dR<minDR2){
2277 minDR2 = dR;
2278 q2 = iMuon->charge();
2279 }
2280 }
2281 double massJpsi = sqrt(2*pt1*pt2*(TMath::CosH(eta1)*TMath::CosH(eta2)-TMath::SinH(eta1)*TMath::SinH(eta2)-TMath::Sin(phi1)*TMath::Sin(phi2)-TMath::Cos(phi1)*TMath::Cos(phi2)));
2282 bool cowboy=false;
2283 if(minDR1<minDR2) cowboy = q1*deltaPhi(phi1,phi2) > 0;
2284 else cowboy = q2*deltaPhi(phi2,phi1) > 0;
2285 if (massJpsi>2.5 && massJpsi<4.1 && !cowboy) return true;
2286 }
2287 }
2288
2289
2290 return false;
2291 }
2292
2293 bool LambdaAnalyzer::pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v2(const
2294 edm::Event& iEvent)
2295 {
2296 return pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v1(iEvent);
2297 }
2298
2299 bool LambdaAnalyzer::pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v3(const
2300 edm::Event& iEvent)
2301 {
2302 return pass_HLT_Mu0_TkMu0_OST_Jpsi_Tight_v2(iEvent);
2303 }
2304
2305
2306
2307 bool LambdaAnalyzer::pass_HLT_DoubleMu0_Quarkonium_v1(const edm::Event& iEvent)
2308 {
2309
2310 edm::Handle<trigger::TriggerEvent> hltevent;
2311 iEvent.getByLabel(triggerEventTag_, hltevent);
2312 const trigger::TriggerObjectCollection& hltobjects(hltevent->getObjects());
2313
2314 trigger::size_type index(0);
2315 edm::InputTag it_track("hltDiMuonL3PreFiltered0::REDIGI36X");
2316
2317 index=hltevent->filterIndex(it_track);
2318 if (index<hltevent->sizeFilters()) {
2319 const trigger::Keys& KEYS (hltevent->filterKeys(index));
2320 const trigger::size_type nTO(KEYS.size());
2321 for (trigger::size_type i=0; i!=nTO; ++i) {
2322 const trigger::TriggerObject& TO1( hltobjects.at(KEYS[i]) );
2323 double pt1 = TO1.pt(); double eta1 = TO1.eta(); double phi1 = TO1.phi();
2324 for (trigger::size_type j=i+1; j!=nTO; ++j) {
2325 const trigger::TriggerObject& TO2( hltobjects.at(KEYS[j]) );
2326 double pt2 = TO2.pt(); double eta2 = TO2.eta(); double phi2 = TO2.phi();
2327 double massJpsi = sqrt(2*pt1*pt2*(TMath::CosH(eta1)*TMath::CosH(eta2)-TMath::SinH(eta1)*TMath::SinH(eta2)-TMath::Sin(phi1)*TMath::Sin(phi2)-TMath::Cos(phi1)*TMath::Cos(phi2)));
2328 if (massJpsi>1.5 && massJpsi<14.5) return true;
2329 }
2330 }
2331 }
2332
2333 return false;
2334 }
2335
2336
2337 void LambdaAnalyzer::push_anything_back()
2338 {
2339
2340 LambdaBMass_afterKinFit->push_back(-99999);
2341 LambdaBPx_afterKinFit->push_back(-99999);
2342 LambdaBPy_afterKinFit->push_back(-99999);
2343 LambdaBPz_afterKinFit->push_back(-99999);
2344 LambdaBX_afterKinFit->push_back(-99999);
2345 LambdaBY_afterKinFit->push_back(-99999);
2346 LambdaBZ_afterKinFit->push_back(-99999);
2347 LambdaBEta_afterKinFit->push_back(-99999);
2348 LambdaBPhi_afterKinFit->push_back(-99999);
2349 LambdaBPt_afterKinFit->push_back(-99999);
2350 LambdaBNDOF_afterKinFit->push_back(-99999);
2351 LambdaBchiSquared_afterKinFit->push_back(-99999);
2352 LambdaBVtxProb_afterKinFit->push_back(-99999);
2353 LambdaBVtxX_afterKinFit->push_back(-99999);
2354 LambdaBVtxY_afterKinFit->push_back(-99999);
2355 LambdaBVtxZ_afterKinFit->push_back(-99999);
2356
2357 mumPx_afterKinFit->push_back(-99999);
2358 mumPy_afterKinFit->push_back(-99999);
2359 mumPz_afterKinFit->push_back(-99999);
2360
2361 mupPx_afterKinFit->push_back(-99999);
2362 mupPy_afterKinFit->push_back(-99999);
2363 mupPz_afterKinFit->push_back(-99999);
2364
2365 JpsiPx_afterKinFit->push_back(-99999);
2366 JpsiPy_afterKinFit->push_back(-99999);
2367 JpsiPz_afterKinFit->push_back(-99999);
2368
2369 LambdaPx_afterKinFit->push_back(-99999);
2370 LambdaPy_afterKinFit->push_back(-99999);
2371 LambdaPz_afterKinFit->push_back(-99999);
2372
2373 JpsiVtxX_afterKinFit->push_back(-99999);
2374 JpsiVtxY_afterKinFit->push_back(-99999);
2375 JpsiVtxZ_afterKinFit->push_back(-99999);
2376
2377 LambdaVtxX_afterKinFit->push_back(-99999);
2378 LambdaVtxY_afterKinFit->push_back(-99999);
2379 LambdaVtxZ_afterKinFit->push_back(-99999);
2380
2381 LambdaVtxProb_afterKinFit->push_back(-99999);
2382 JpsiVtxProb_afterKinFit->push_back(-99999);
2383
2384 openingAngleLambdaB_afterKinFit->push_back(-99999);
2385 ctauLambdaB3DPV_afterKinFit->push_back(-99999);
2386 ctauLambdaBxyPV_afterKinFit->push_back(-99999);
2387 ctauLambdaXY_afterKinFit->push_back(-99999);
2388
2389 return;
2390
2391 }
2392
2393
2394 //define this as a plug-in
2395 DEFINE_FWK_MODULE(LambdaAnalyzer);