ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/codeHeavyBaryons/LambdaB/src/LambdaAnalyzer.cc
Revision: 1.2
Committed: Fri Dec 10 10:29:56 2010 UTC (14 years, 5 months ago) by mirena
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +1 -1 lines
State: FILE REMOVED
Log Message:
*** empty log message ***

File Contents

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