51 |
|
#include <SimDataFormats/Track/interface/SimTrackContainer.h> |
52 |
|
// BeamSpot |
53 |
|
#include "DataFormats/BeamSpot/interface/BeamSpot.h" |
54 |
– |
// Trigger |
55 |
– |
#include "CondFormats/L1TObjects/interface/L1GtTriggerMenu.h" |
56 |
– |
#include "CondFormats/DataRecord/interface/L1GtTriggerMenuRcd.h" |
57 |
– |
#include "CondFormats/L1TObjects/interface/L1GtPrescaleFactors.h" |
58 |
– |
#include "CondFormats/DataRecord/interface/L1GtPrescaleFactorsAlgoTrigRcd.h" |
59 |
– |
#include "CondFormats/DataRecord/interface/L1GtPrescaleFactorsTechTrigRcd.h" |
60 |
– |
#include "CondFormats/L1TObjects/interface/L1GtTriggerMask.h" |
61 |
– |
#include "CondFormats/DataRecord/interface/L1GtTriggerMaskAlgoTrigRcd.h" |
62 |
– |
#include "CondFormats/DataRecord/interface/L1GtTriggerMaskTechTrigRcd.h" |
63 |
– |
#include "CondFormats/DataRecord/interface/L1GtTriggerMaskVetoAlgoTrigRcd.h" |
64 |
– |
#include "CondFormats/DataRecord/interface/L1GtTriggerMaskVetoTechTrigRcd.h" |
65 |
– |
#include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetupFwd.h" |
66 |
– |
#include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutSetup.h" |
67 |
– |
#include "DataFormats/L1GlobalTrigger/interface/L1GlobalTriggerReadoutRecord.h" |
68 |
– |
|
54 |
|
|
55 |
|
//root |
56 |
|
#include <TROOT.h> |
80 |
|
realData_ = iConfig.getUntrackedParameter<bool>("realData",false); |
81 |
|
analyze_ = iConfig.getUntrackedParameter<bool>("analyzeOnTheFly",false); |
82 |
|
saventuple_ = iConfig.getUntrackedParameter<bool>("saventuple",false); |
83 |
< |
outputfilename_ = iConfig.getUntrackedParameter<string>("OutputFileName"); |
83 |
> |
outputfilename_ = iConfig.getUntrackedParameter<string>("OutputFileName"); |
84 |
> |
histoFileName_ = iConfig.getUntrackedParameter<std::string> ("histoFileName"); |
85 |
|
ntrkdiffcut_ = iConfig.getUntrackedParameter<double>("ntrkdiffcut"); |
86 |
|
nTrkMin_ = iConfig.getUntrackedParameter<int>("nTrkMin"); |
87 |
|
nTrkMax_ = iConfig.getUntrackedParameter<int>("nTrkMax"); |
88 |
|
zsigncut_ = iConfig.getUntrackedParameter<int>("zsigncut"); |
89 |
|
bsSrc = iConfig.getParameter< edm::InputTag >("beamSpot"); |
104 |
– |
applyEvtClean_ = iConfig.getUntrackedParameter<bool>("applyEvtClean",false); |
105 |
– |
histoFileName_ = iConfig.getUntrackedParameter<std::string> ("histoFileName"); |
90 |
|
|
91 |
|
//now do what ever initialization is needed |
92 |
|
pvinfo.clear(); |
93 |
|
|
110 |
– |
|
94 |
|
// Specify the data mode vector |
95 |
|
if(realData_) datamode.push_back(0); |
96 |
|
else { |
110 |
|
|
111 |
|
// pvtxtree_ analyzing the pvtxs ootb |
112 |
|
pvtxtree_ = new TTree("pvtxtree","pvtxtree"); |
130 |
– |
//pvtx using all tracks |
131 |
– |
pvtxtree_->Branch("nrecPV",&nrecPV_,"nrecPV/I"); |
132 |
– |
pvtxtree_->Branch("nTrkPV",&nTrkPV_,"nTrkPV[nrecPV]/I"); |
133 |
– |
pvtxtree_->Branch("sumptsq",&sumptsq_,"sumptsq[nrecPV]/D"); |
134 |
– |
pvtxtree_->Branch("isValid",&isValid_,"isValid/I"); |
135 |
– |
pvtxtree_->Branch("isFake",&isFake_,"isFake/I"); |
136 |
– |
pvtxtree_->Branch("recx",&recx_,"recx[nrecPV]/D"); |
137 |
– |
pvtxtree_->Branch("recy",&recy_,"recy[nrecPV]/D"); |
138 |
– |
pvtxtree_->Branch("recz",&recz_,"recz[nrecPV]/D"); |
139 |
– |
pvtxtree_->Branch("recx_err",&recx_err_,"recx_err[nrecPV]/D"); |
140 |
– |
pvtxtree_->Branch("recy_err",&recy_err_,"recy_err[nrecPV]/D"); |
141 |
– |
pvtxtree_->Branch("recz_err",&recz_err_,"recz_err[nrecPV]/D"); |
113 |
|
|
114 |
+ |
// Event information for the data |
115 |
+ |
pvtxtree_->Branch("glob_runno",&glob_runno_,"glob_runno/I"); |
116 |
+ |
pvtxtree_->Branch("glob_evtno",&glob_evtno_,"glob_evtno/I"); |
117 |
+ |
pvtxtree_->Branch("glob_ls",&glob_ls_,"glob_ls/I"); |
118 |
+ |
pvtxtree_->Branch("glob_bx",&glob_bx_,"glob_bx/I"); |
119 |
+ |
|
120 |
+ |
|
121 |
+ |
pvtxtree_->Branch("nrecPV",&nrecPV_,"nrecPV/I"); |
122 |
+ |
pvtxtree_->Branch("nrecPV1_spl",&nrecPV1_spl_,"nrecPV1_spl/I"); |
123 |
+ |
pvtxtree_->Branch("nrecPV2_spl",&nrecPV2_spl_,"nrecPV2_spl/I"); |
124 |
+ |
// Event level information |
125 |
|
pvtxtree_->Branch("min_zsep",&min_zsep_,"min_zsep/D"); |
126 |
|
pvtxtree_->Branch("min_ntrksep",&min_ntrksep_,"min_ntrksep/D"); |
127 |
|
pvtxtree_->Branch("min_sumpt2sep",&min_sumpt2sep_,"min_sumpt2sep/D"); |
128 |
|
|
129 |
< |
//pvtx using splittracks1 |
130 |
< |
pvtxtree_->Branch("nrecPV1_spl",&nrecPV1_spl_,"nrecPV1_spl/I"); |
131 |
< |
pvtxtree_->Branch("nTrkPV1_spl",&nTrkPV1_spl_,"nTrkPV1_spl[nrecPV1_spl]/I"); |
132 |
< |
pvtxtree_->Branch("sumptsq1_spl",&sumptsq1_spl_,"sumptsq1_spl[nrecPV1_spl]/D"); |
151 |
< |
pvtxtree_->Branch("isValid1_spl",&isValid1_spl_,"isValid1_spl/I"); |
152 |
< |
pvtxtree_->Branch("isFake1_spl",&isFake1_spl_,"isFake1_spl/I"); |
153 |
< |
pvtxtree_->Branch("recx1_spl",&recx1_spl_,"recx1_spl[nrecPV1_spl]/D"); |
154 |
< |
pvtxtree_->Branch("recy1_spl",&recy1_spl_,"recy1_spl[nrecPV1_spl]/D"); |
155 |
< |
pvtxtree_->Branch("recz1_spl",&recz1_spl_,"recz1_spl[nrecPV1_spl]/D"); |
156 |
< |
pvtxtree_->Branch("recx1_err_spl",&recx1_err_spl_,"recx1_err_spl[nrecPV1_spl]/D"); |
157 |
< |
pvtxtree_->Branch("recy1_err_spl",&recy1_err_spl_,"recy1_err_spl[nrecPV1_spl]/D"); |
158 |
< |
pvtxtree_->Branch("recz1_err_spl",&recz1_err_spl_,"recz1_err_spl[nrecPV1_spl]/D"); |
159 |
< |
|
160 |
< |
//pvtx using splittracks2 |
161 |
< |
pvtxtree_->Branch("nrecPV2_spl",&nrecPV2_spl_,"nrecPV2_spl/I"); |
162 |
< |
pvtxtree_->Branch("nTrkPV2_spl",&nTrkPV2_spl_,"nTrkPV2_spl[nrecPV2_spl]/I"); |
163 |
< |
pvtxtree_->Branch("sumptsq2_spl",&sumptsq2_spl_,"sumptsq2_spl[nrecPV2_spl]/D"); |
164 |
< |
pvtxtree_->Branch("isValid2_spl",&isValid2_spl_,"isValid2_spl/I"); |
165 |
< |
pvtxtree_->Branch("isFake2_spl",&isFake2_spl_,"isFake2_spl/I"); |
166 |
< |
pvtxtree_->Branch("recx2_spl",&recx2_spl_,"recx2_spl[nrecPV2_spl]/D"); |
167 |
< |
pvtxtree_->Branch("recy2_spl",&recy2_spl_,"recy2_spl[nrecPV2_spl]/D"); |
168 |
< |
pvtxtree_->Branch("recz2_spl",&recz2_spl_,"recz2_spl[nrecPV2_spl]/D"); |
169 |
< |
pvtxtree_->Branch("recx2_err_spl",&recx2_err_spl_,"recx2_err_spl[nrecPV2_spl]/D"); |
170 |
< |
pvtxtree_->Branch("recy2_err_spl",&recy2_err_spl_,"recy2_err_spl[nrecPV2_spl]/D"); |
171 |
< |
pvtxtree_->Branch("recz2_err_spl",&recz2_err_spl_,"recz2_err_spl[nrecPV2_spl]/D"); |
172 |
< |
|
173 |
< |
//pixeVertices |
174 |
< |
pvtxtree_->Branch("nrecPV_pxlpvtx",&nrecPV_pxlpvtx_,"nrecPV_pxlpvtx/I"); |
175 |
< |
pvtxtree_->Branch("nTrkPV_pxlpvtx",&nTrkPV_pxlpvtx_,"nTrkPV_pxlpvtx[nrecPV_pxlpvtx]/I"); |
176 |
< |
pvtxtree_->Branch("sumptsq_pxlpvtx",&sumptsq_pxlpvtx_,"sumptsq_pxlpvtx[nrecPV_pxlpvtx]/D"); |
177 |
< |
pvtxtree_->Branch("isValid_pxlpvtx",&isValid_pxlpvtx_,"isValid_pxlpvtx/I"); |
178 |
< |
pvtxtree_->Branch("isFake_pxlpvtx",&isFake_pxlpvtx_,"isFake_pxlpvtx/I"); |
179 |
< |
pvtxtree_->Branch("recx_pxlpvtx",&recx_pxlpvtx_,"recx_pxlpvtx[nrecPV_pxlpvtx]/D"); |
180 |
< |
pvtxtree_->Branch("recy_pxlpvtx",&recy_pxlpvtx_,"recy_pxlpvtx[nrecPV_pxlpvtx]/D"); |
181 |
< |
pvtxtree_->Branch("recz_pxlpvtx",&recz_pxlpvtx_,"recz_pxlpvtx[nrecPV_pxlpvtx]/D"); |
182 |
< |
pvtxtree_->Branch("recx_err_pxlpvtx",&recx_err_pxlpvtx_,"recx_err_pxlpvtx[nrecPV_pxlpvtx]/D"); |
183 |
< |
pvtxtree_->Branch("recy_err_pxlpvtx",&recy_err_pxlpvtx_,"recy_err_pxlpvtx[nrecPV_pxlpvtx]/D"); |
184 |
< |
pvtxtree_->Branch("recz_err_pxlpvtx",&recz_err_pxlpvtx_,"recz_err_pxlpvtx[nrecPV_pxlpvtx]/D"); |
129 |
> |
//Fill the variables in the twovtx pair (recvtx1, recvtx2) |
130 |
> |
// Information related to the analyzing the two-vertex method |
131 |
> |
|
132 |
> |
|
133 |
|
|
186 |
– |
//Fill the variables in the twovtx pair (recvtx1, recvtx2) |
134 |
|
pvtxtree_->Branch("nrecPV_twovtx",&nrecPV_twovtx_,"nrecPV_twovtx/I"); |
188 |
– |
pvtxtree_->Branch("nTrkPV1_spl_twovtx",&nTrkPV1_spl_twovtx_,"nTrkPV1_spl_twovtx[nrecPV_twovtx]/I"); |
189 |
– |
pvtxtree_->Branch("nTrkPV2_spl_twovtx",&nTrkPV2_spl_twovtx_,"nTrkPV2_spl_twovtx[nrecPV_twovtx]/I"); |
190 |
– |
pvtxtree_->Branch("sumptsq1_spl_twovtx",&sumptsq1_spl_twovtx_,"sumptsq1_spl_twovtx[nrecPV_twovtx]/D"); |
191 |
– |
pvtxtree_->Branch("sumptsq2_spl_twovtx",&sumptsq2_spl_twovtx_,"sumptsq2_spl_twovtx[nrecPV_twovtx]/D"); |
135 |
|
pvtxtree_->Branch("nTrkPV_twovtx",&nTrkPV_twovtx_,"nTrkPV_twovtx[nrecPV_twovtx]/I"); |
136 |
|
pvtxtree_->Branch("deltax_twovtx",&deltax_twovtx_,"deltax_twovtx[nrecPV_twovtx]/D"); |
137 |
|
pvtxtree_->Branch("deltay_twovtx",&deltay_twovtx_,"deltay_twovtx[nrecPV_twovtx]/D"); |
142 |
|
pvtxtree_->Branch("pullx_twovtx",&pullx_twovtx_,"pullx_twovtx[nrecPV_twovtx]/D"); |
143 |
|
pvtxtree_->Branch("pully_twovtx",&pully_twovtx_,"pully_twovtx[nrecPV_twovtx]/D"); |
144 |
|
pvtxtree_->Branch("pullz_twovtx",&pullz_twovtx_,"pullz_twovtx[nrecPV_twovtx]/D"); |
145 |
+ |
|
146 |
+ |
// Information for the splitVertexColl1 |
147 |
+ |
pvtxtree_->Branch("nTrkPV1_spl_twovtx",&nTrkPV1_spl_twovtx_,"nTrkPV1_spl_twovtx[nrecPV_twovtx]/I"); |
148 |
+ |
pvtxtree_->Branch("ndofPV1_spl_twovtx",&ndofPV1_spl_twovtx_,"ndofPV1_spl_twovtx[nrecPV_twovtx]/D"); |
149 |
+ |
pvtxtree_->Branch("normchi2PV1_spl_twovtx",&normchi2PV1_spl_twovtx_,"normchi2PV1_spl_twovtx[nrecPV_twovtx]/D"); |
150 |
+ |
pvtxtree_->Branch("avgPtPV1_spl_twovtx",&avgPtPV1_spl_twovtx_,"avgPtPV1_spl_twovtx[nrecPV_twovtx]/D"); |
151 |
+ |
pvtxtree_->Branch("errx1_spl_twovtx",&errx1_spl_twovtx_,"errx1_spl_twovtx[nrecPV_twovtx]/D"); |
152 |
+ |
pvtxtree_->Branch("erry1_spl_twovtx",&erry1_spl_twovtx_,"erry1_spl_twovtx[nrecPV_twovtx]/D"); |
153 |
+ |
pvtxtree_->Branch("errz1_spl_twovtx",&errz1_spl_twovtx_,"errz1_spl_twovtx[nrecPV_twovtx]/D"); |
154 |
+ |
|
155 |
+ |
|
156 |
+ |
// Information for the splitVertexColl2 |
157 |
+ |
pvtxtree_->Branch("nTrkPV2_spl_twovtx",&nTrkPV2_spl_twovtx_,"nTrkPV2_spl_twovtx[nrecPV_twovtx]/I"); |
158 |
+ |
pvtxtree_->Branch("ndofPV2_spl_twovtx",&ndofPV2_spl_twovtx_,"ndofPV2_spl_twovtx[nrecPV_twovtx]/D"); |
159 |
+ |
pvtxtree_->Branch("normchi2PV2_spl_twovtx",&normchi2PV2_spl_twovtx_,"normchi2PV2_spl_twovtx[nrecPV_twovtx]/D"); |
160 |
+ |
pvtxtree_->Branch("avgPtPV2_spl_twovtx",&avgPtPV2_spl_twovtx_,"avgPtPV2_spl_twovtx[nrecPV_twovtx]/D"); |
161 |
+ |
pvtxtree_->Branch("errx2_spl_twovtx",&errx2_spl_twovtx_,"errx2_spl_twovtx[nrecPV_twovtx]/D"); |
162 |
+ |
pvtxtree_->Branch("erry2_spl_twovtx",&erry2_spl_twovtx_,"erry2_spl_twovtx[nrecPV_twovtx]/D"); |
163 |
+ |
pvtxtree_->Branch("errz2_spl_twovtx",&errz2_spl_twovtx_,"errz2_spl_twovtx[nrecPV_twovtx]/D"); |
164 |
+ |
|
165 |
|
|
166 |
|
// MC variables |
167 |
|
if(!realData_) { |
179 |
|
pvtxtree_->Branch("deltaz_mct",&deltaz_mct_,"deltaz_mct[nrecPV_mct]/D"); |
180 |
|
pvtxtree_->Branch("pullx_mct",&pullx_mct_,"pullx_mct[nrecPV_mct]/D"); |
181 |
|
pvtxtree_->Branch("pully_mct",&pully_mct_,"pully_mct[nrecPV_mct]/D"); |
182 |
< |
pvtxtree_->Branch("pullz_mct",&pullz_mct_,"pullz_mct[nrecPV_mct]/D"); |
182 |
> |
pvtxtree_->Branch("pullz_mct",&pullz_mct_,"pullz_mct[nrecPV_mct]/D"); |
183 |
> |
pvtxtree_->Branch("errx_mct",&errx_mct_,"errx_mct[nrecPV_mct]/D"); |
184 |
> |
pvtxtree_->Branch("erry_mct",&erry_mct_,"erry_mct[nrecPV_mct]/D"); |
185 |
> |
pvtxtree_->Branch("errz_mct",&errz_mct_,"errz_mct[nrecPV_mct]/D"); |
186 |
> |
pvtxtree_->Branch("nTrkPV_mct",&nTrkPV_mct_,"nTrkPV_mct[nrecPV_mct]/I"); |
187 |
> |
pvtxtree_->Branch("ndofPV_mct",&ndofPV_mct_,"ndofPV_mct[nrecPV_mct]/D"); |
188 |
> |
pvtxtree_->Branch("normchi2PV_mct",&normchi2PV_mct_,"normchi2PV_mct[nrecPV_mct]/D"); |
189 |
> |
pvtxtree_->Branch("avgPtPV_mct",&avgPtPV_mct_,"avgPtPV_mct[nrecPV_mct]/D"); |
190 |
> |
|
191 |
|
// For pvtxs with splittracks1 |
192 |
|
pvtxtree_->Branch("nrecPV_spl1_mct",&nrecPV_spl1_mct_,"nrecPV_spl1_mct/I"); |
193 |
|
pvtxtree_->Branch("deltax_spl1_mct",&deltax_spl1_mct_,"deltax_spl1_mct[nrecPV_spl1_mct]/D"); |
196 |
|
pvtxtree_->Branch("pullx_spl1_mct",&pullx_spl1_mct_,"pullx_spl1_mct[nrecPV_spl1_mct]/D"); |
197 |
|
pvtxtree_->Branch("pully_spl1_mct",&pully_spl1_mct_,"pully_spl1_mct[nrecPV_spl1_mct]/D"); |
198 |
|
pvtxtree_->Branch("pullz_spl1_mct",&pullz_spl1_mct_,"pullz_spl1_mct[nrecPV_spl1_mct]/D"); |
199 |
+ |
pvtxtree_->Branch("errx_spl1_mct",&errx_spl1_mct_,"errx_spl1_mct[nrecPV_spl1_mct]/D"); |
200 |
+ |
pvtxtree_->Branch("erry_spl1_mct",&erry_spl1_mct_,"erry_spl1_mct[nrecPV_spl1_mct]/D"); |
201 |
+ |
pvtxtree_->Branch("errz_spl1_mct",&errz_spl1_mct_,"errz_spl1_mct[nrecPV_spl1_mct]/D"); |
202 |
+ |
pvtxtree_->Branch("nTrkPV_spl1_mct",&nTrkPV_spl1_mct_,"nTrkPV_spl1_mct[nrecPV_spl1_mct]/I"); |
203 |
+ |
pvtxtree_->Branch("ndofPV_spl1_mct",&ndofPV_spl1_mct_,"ndofPV_spl1_mct[nrecPV_spl1_mct]/D"); |
204 |
+ |
pvtxtree_->Branch("normchi2PV_spl1_mct",&normchi2PV_spl1_mct_,"normchi2PV_spl1_mct[nrecPV_spl1_mct]/D"); |
205 |
+ |
pvtxtree_->Branch("avgPtPV_spl1_mct",&avgPtPV_spl1_mct_,"avgPtPV_spl1_mct[nrecPV_spl1_mct]/D"); |
206 |
+ |
|
207 |
|
// For pvtxs with splittracks1 |
208 |
|
pvtxtree_->Branch("nrecPV_spl2_mct",&nrecPV_spl2_mct_,"nrecPV_spl2_mct/I"); |
209 |
|
pvtxtree_->Branch("deltax_spl2_mct",&deltax_spl2_mct_,"deltax_spl2_mct[nrecPV_spl2_mct]/D"); |
211 |
|
pvtxtree_->Branch("deltaz_spl2_mct",&deltaz_spl2_mct_,"deltaz_spl2_mct[nrecPV_spl2_mct]/D"); |
212 |
|
pvtxtree_->Branch("pullx_spl2_mct",&pullx_spl2_mct_,"pullx_spl2_mct[nrecPV_spl2_mct]/D"); |
213 |
|
pvtxtree_->Branch("pully_spl2_mct",&pully_spl2_mct_,"pully_spl2_mct[nrecPV_spl2_mct]/D"); |
214 |
< |
pvtxtree_->Branch("pullz_spl2_mct",&pullz_spl2_mct_,"pullz_spl2_mct[nrecPV_spl2_mct]/D"); |
214 |
> |
pvtxtree_->Branch("pullz_spl2_mct",&pullz_spl2_mct_,"pullz_spl2_mct[nrecPV_spl2_mct]/D"); |
215 |
> |
pvtxtree_->Branch("errx_spl2_mct",&errx_spl2_mct_,"errx_spl2_mct[nrecPV_spl2_mct]/D"); |
216 |
> |
pvtxtree_->Branch("erry_spl2_mct",&erry_spl2_mct_,"erry_spl2_mct[nrecPV_spl2_mct]/D"); |
217 |
> |
pvtxtree_->Branch("errz_spl2_mct",&errz_spl2_mct_,"errz_spl2_mct[nrecPV_spl2_mct]/D"); |
218 |
> |
pvtxtree_->Branch("nTrkPV_spl2_mct",&nTrkPV_spl2_mct_,"nTrkPV_spl2_mct[nrecPV_spl2_mct]/I"); |
219 |
> |
pvtxtree_->Branch("ndofPV_spl2_mct",&ndofPV_spl2_mct_,"ndofPV_spl2_mct[nrecPV_spl2_mct]/D"); |
220 |
> |
pvtxtree_->Branch("normchi2PV_spl2_mct",&normchi2PV_spl2_mct_,"normchi2PV_spl2_mct[nrecPV_spl2_mct]/D"); |
221 |
> |
pvtxtree_->Branch("avgPtPV_spl2_mct",&avgPtPV_spl2_mct_,"avgPtPV_spl2_mct[nrecPV_spl2_mct]/D"); |
222 |
> |
|
223 |
|
} |
224 |
|
} |
225 |
|
|
243 |
|
|
244 |
|
// Book Histograms: |
245 |
|
h_pvtrk = new PVHistograms(); |
259 |
– |
h_pixvtx = new PVHistograms(); |
246 |
|
h_misc = new PVHistograms(); |
247 |
|
h_summary = new PVHistograms(); |
248 |
|
h_others = new PVHistograms(); |
255 |
|
h_pvtrk->Init("pvTrk", ss.str(),"spl"); |
256 |
|
} |
257 |
|
} |
272 |
– |
h_pixvtx->Init("pixVtx"); |
258 |
|
h_misc->Init("misc"); |
259 |
|
|
260 |
|
// Book MC only plots |
303 |
|
theFile->cd(); |
304 |
|
theFile->cd("Summary"); |
305 |
|
h_pvtrk->Save(); |
321 |
– |
h_pixvtx->Save(); |
306 |
|
h_misc->Save(); |
307 |
|
h_summary->Save(); |
308 |
|
if (!realData_) |
340 |
|
gStyle->SetTitleSize(0.055, ""); // size for pad title; default is 0.02 |
341 |
|
gStyle->SetLabelSize(0.03, "XYZ"); // size for axis labels; default is 0.04 |
342 |
|
gStyle->SetStatFontSize(0.08); // size for stat. box |
343 |
< |
gStyle->SetTitleFont(32, "XYZ"); // times-bold-italic font (p. 153) for axes |
344 |
< |
gStyle->SetTitleFont(32, ""); // same for pad title |
345 |
< |
gStyle->SetLabelFont(32, "XYZ"); // same for axis labels |
346 |
< |
gStyle->SetStatFont(32); // same for stat. box |
343 |
> |
gStyle->SetTitleFont(42, "XYZ"); // times-bold-italic font (p. 153) for axes |
344 |
> |
gStyle->SetTitleFont(42, ""); // same for pad title |
345 |
> |
gStyle->SetLabelFont(42, "XYZ"); // same for axis labels |
346 |
> |
gStyle->SetStatFont(42); // same for stat. box |
347 |
|
gStyle->SetLabelOffset(0.006, "Y"); // default is 0.005 |
348 |
|
// |
349 |
|
return; |
530 |
|
using namespace edm; |
531 |
|
using namespace reco; |
532 |
|
|
533 |
+ |
//======================================================================== |
534 |
+ |
// Step 0: Prepare root variables and get information from the Event |
535 |
+ |
//======================================================================== |
536 |
+ |
|
537 |
|
edm::LogInfo("Debug")<<"[PVStudy]"<<endl; |
550 |
– |
//======================================================= |
538 |
|
// Initialize Root-tuple variables if needed |
539 |
< |
//======================================================= |
540 |
< |
//if(saventuple_) |
541 |
< |
SetVarToZero(); |
555 |
< |
|
556 |
< |
//======================================================= |
557 |
< |
// Track accessors |
558 |
< |
//======================================================= |
559 |
< |
//trackCollection |
539 |
> |
SetVarToZero(); |
540 |
> |
|
541 |
> |
// ====== TrackCollection |
542 |
|
static const reco::TrackCollection s_empty_trackColl; |
543 |
|
const reco::TrackCollection *trackColl = &s_empty_trackColl; |
544 |
|
edm::Handle<reco::TrackCollection> trackCollectionHandle; |
548 |
|
} else { |
549 |
|
edm::LogInfo("Debug") << "[PVStudy] trackCollection cannot be found -> using empty collection of same type." <<endl; |
550 |
|
} |
551 |
< |
//splitTrackCollection1 |
551 |
> |
// ====== splitTrackCollection1 |
552 |
|
static const reco::TrackCollection s_empty_splitTrackColl1; |
553 |
|
const reco::TrackCollection *splitTrackColl1 = &s_empty_splitTrackColl1; |
554 |
|
edm::Handle<reco::TrackCollection> splitTrackCollection1Handle; |
558 |
|
} else { |
559 |
|
edm::LogInfo("Debug") << "[PVStudy] splitTrackCollection1 cannot be found -> using empty collection of same type." <<endl; |
560 |
|
} |
561 |
< |
//splitTrackCollection2 |
561 |
> |
// ====== splitTrackCollection2 |
562 |
|
static const reco::TrackCollection s_empty_splitTrackColl2; |
563 |
|
const reco::TrackCollection *splitTrackColl2 = &s_empty_splitTrackColl2; |
564 |
|
edm::Handle<reco::TrackCollection> splitTrackCollection2Handle; |
569 |
|
edm::LogInfo("Debug") << "[PVStudy] splitTrackCollection2 cannot be found -> using empty collection of same type." <<endl; |
570 |
|
} |
571 |
|
|
572 |
< |
//======================================================= |
591 |
< |
// PVTX accessors |
592 |
< |
//======================================================= |
593 |
< |
//vertexCollection |
572 |
> |
// ======= PrimaryVertexCollection |
573 |
|
static const reco::VertexCollection s_empty_vertexColl; |
574 |
|
const reco::VertexCollection *vertexColl = &s_empty_vertexColl; |
575 |
|
edm::Handle<reco::VertexCollection> vertexCollectionHandle; |
579 |
|
} else { |
580 |
|
edm::LogInfo("Debug") << "[PVStudy] vertexCollection cannot be found -> using empty collection of same type." <<endl; |
581 |
|
} |
582 |
< |
//splitVertexCollection1 |
582 |
> |
// ====== splitVertexCollection1 |
583 |
|
static const reco::VertexCollection s_empty_splitVertexColl1; |
584 |
|
const reco::VertexCollection *splitVertexColl1 = &s_empty_splitVertexColl1; |
585 |
|
edm::Handle<reco::VertexCollection> splitVertexCollection1Handle; |
589 |
|
} else { |
590 |
|
edm::LogInfo("Debug") << "[PVStudy] splitVertexCollection1 cannot be found -> using empty collection of same type." <<endl; |
591 |
|
} |
592 |
< |
//splitVertexCollection2 |
592 |
> |
// ====== splitVertexCollection2 |
593 |
|
static const reco::VertexCollection s_empty_splitVertexColl2; |
594 |
|
const reco::VertexCollection *splitVertexColl2 = &s_empty_splitVertexColl2; |
595 |
|
edm::Handle<reco::VertexCollection> splitVertexCollection2Handle; |
601 |
|
} |
602 |
|
|
603 |
|
|
604 |
< |
//======================================================= |
626 |
< |
// GET pixelVertices |
627 |
< |
//======================================================= |
628 |
< |
static const reco::VertexCollection s_empty_pixelVertexColl; |
629 |
< |
const reco::VertexCollection *pixelVertexColl = &s_empty_pixelVertexColl; |
630 |
< |
edm::Handle<reco::VertexCollection> pixelVertexCollectionHandle; |
631 |
< |
iEvent.getByLabel(pixelVertexCollectionTag_, pixelVertexCollectionHandle); |
632 |
< |
if( iEvent.getByLabel(pixelVertexCollectionTag_, pixelVertexCollectionHandle)) { |
633 |
< |
pixelVertexColl = pixelVertexCollectionHandle.product(); |
634 |
< |
} else { |
635 |
< |
edm::LogInfo("Debug") << "[PVStudy] pixelVertexCollection cannot be found. -> using empty collection of same type." <<endl; |
636 |
< |
} |
637 |
< |
|
638 |
< |
//======================================================= |
639 |
< |
// BeamSpot accessors |
640 |
< |
//======================================================= |
641 |
< |
|
604 |
> |
// ======== BeamSpot accessors |
605 |
|
edm::Handle<reco::BeamSpot> recoBeamSpotHandle; |
606 |
|
iEvent.getByLabel(bsSrc,recoBeamSpotHandle); |
607 |
|
reco::BeamSpot bs = *recoBeamSpotHandle; |
608 |
|
const Point beamSpot = recoBeamSpotHandle.isValid() ? Point(recoBeamSpotHandle->x0(), recoBeamSpotHandle->y0(), recoBeamSpotHandle->z0()) : Point(0, 0, 0); |
609 |
|
|
610 |
< |
edm::LogInfo("Debug")<<"[PVStudy] End accessing the track, beamSpot, primary vertex and pixelvertices collections"<<endl; |
610 |
> |
edm::LogInfo("Debug")<<"[PVStudy] End accessing the track, beamSpot, primary vertex collections"<<endl; |
611 |
|
|
612 |
< |
//======================================================= |
650 |
< |
// MC simvtx accessor |
651 |
< |
//======================================================= |
612 |
> |
// ========== MC simvtx accessor |
613 |
|
if (!realData_) { |
614 |
|
edm::Handle<SimVertexContainer> simVtxs; |
615 |
|
iEvent.getByLabel( simG4_, simVtxs); |
618 |
|
iEvent.getByLabel( simG4_, simTrks); |
619 |
|
} |
620 |
|
|
621 |
< |
//======================================================= |
661 |
< |
// GET PDT |
662 |
< |
//======================================================= |
621 |
> |
// ========== GET PDT |
622 |
|
try{ |
623 |
|
iSetup.getData(pdt); |
624 |
|
}catch(const Exception&){ |
625 |
|
edm::LogInfo("Debug") << "[PVStudy] Some problem occurred with the particle data table. This may not work !." <<endl; |
626 |
|
} |
627 |
< |
|
628 |
< |
//======================================================= |
629 |
< |
// Apply event cleaning for firstcoll data and MC |
630 |
< |
//======================================================= |
631 |
< |
if(applyEvtClean_) { |
632 |
< |
// =====Select on the trigger bits |
633 |
< |
edm::ESHandle<L1GtTriggerMenu> menuRcd; |
634 |
< |
iSetup.get<L1GtTriggerMenuRcd>().get(menuRcd) ; |
635 |
< |
const L1GtTriggerMenu* menu = menuRcd.product(); |
636 |
< |
edm::Handle< L1GlobalTriggerReadoutRecord > gtRecord; |
637 |
< |
iEvent.getByLabel( edm::InputTag("gtDigis"), gtRecord); |
638 |
< |
|
680 |
< |
bool isTechBit0 = false; |
681 |
< |
bool isTechBit40 = false; |
682 |
< |
bool isBeamHalo = false; |
683 |
< |
|
684 |
< |
TechnicalTriggerWord tw = gtRecord->technicalTriggerWord(); |
685 |
< |
if ( ! tw.empty() ) { |
686 |
< |
// loop over dec. bit to get total rate (no overlap) |
687 |
< |
for ( int i = 0; i < 64; ++i ) { |
688 |
< |
if ( tw[i] ) { |
689 |
< |
//cout<<"technical number "<<i<<" "<<endl; |
690 |
< |
if (i == 0) isTechBit0 = true; |
691 |
< |
if (i < 40 && i > 35) isBeamHalo = true; // The beamHalo bits are 36-39 |
692 |
< |
if (i == 40 || i == 41) isTechBit40 = true; |
693 |
< |
} |
694 |
< |
} |
695 |
< |
} // =====End select on the trigger bits |
696 |
< |
// =====Loop over the trackColl to tet the fraction of HighPurity tracks |
697 |
< |
int nTracks = 0; |
698 |
< |
int nHighPurityTracks=0; |
699 |
< |
double fHighPurity=0; |
700 |
< |
|
701 |
< |
for (unsigned int i=0; i<trackCollectionHandle->size(); i++, nTracks++) { |
702 |
< |
TrackRef tkref(trackCollectionHandle,i); |
703 |
< |
if( (tkref->qualityMask() & 4 ) == 4) ++nHighPurityTracks; |
704 |
< |
} |
705 |
< |
if(nTracks>0) |
706 |
< |
fHighPurity = double(nHighPurityTracks)/double(nTracks); |
707 |
< |
if(verbose_) |
708 |
< |
cout<<"fraction of HighPurity track is "<<fHighPurity<<endl; |
709 |
< |
|
710 |
< |
// Apply the event selection for MC events |
711 |
< |
if(realData_) { |
712 |
< |
if(!isTechBit40 || isBeamHalo || (fHighPurity<0.2&&nTracks>10) || vertexColl->begin()->isFake()) { |
713 |
< |
glob_runno_ = iEvent.id().run(); |
714 |
< |
glob_evtno_ = iEvent.id().event(); |
715 |
< |
glob_ls_ = iEvent.luminosityBlock(); |
716 |
< |
glob_bx_ = iEvent.bunchCrossing(); |
717 |
< |
return; |
718 |
< |
} |
719 |
< |
else { |
720 |
< |
// TechBit 0, beamHalo are not simulted in MC |
721 |
< |
if ( !isTechBit0 || !isTechBit40 || isBeamHalo || (fHighPurity<0.2&&nTracks>10) || vertexColl->begin()->isFake()) return; |
722 |
< |
} |
627 |
> |
|
628 |
> |
// ======= Get MC information if needed |
629 |
> |
bool MC=false; |
630 |
> |
Handle<HepMCProduct> evtMC; |
631 |
> |
if (!realData_) { |
632 |
> |
iEvent.getByLabel("generator",evtMC); |
633 |
> |
if (!evtMC.isValid()) { |
634 |
> |
MC=false; |
635 |
> |
edm::LogInfo("Debug") << "[PVStudy] no HepMCProduct found"<< endl; |
636 |
> |
} else { |
637 |
> |
edm::LogInfo("Debug") << "[PVStudy] generator HepMCProduct found"<< endl; |
638 |
> |
MC=true; |
639 |
|
} |
640 |
< |
} // End of Apply event cleaning cuts for firstcoll data |
725 |
< |
|
640 |
> |
} |
641 |
|
|
642 |
+ |
//======================================================================== |
643 |
+ |
// Step 1: Apply event cleaning for data and MC |
644 |
+ |
// WARNING: event selection cut are hard coded!! |
645 |
+ |
//======================================================================== |
646 |
+ |
// =====Loop over the trackColl to get the fraction of HighPurity tracks |
647 |
+ |
int nTracks = 0; |
648 |
+ |
int nHighPurityTracks=0; |
649 |
+ |
double fHighPurity=0; |
650 |
+ |
|
651 |
+ |
for (unsigned int i=0; i<trackCollectionHandle->size(); i++, nTracks++) { |
652 |
+ |
TrackRef tkref(trackCollectionHandle,i); |
653 |
+ |
if(tkref->quality(reco::TrackBase::highPurity)) ++nHighPurityTracks; |
654 |
+ |
} |
655 |
+ |
if(nTracks>0) |
656 |
+ |
fHighPurity = double(nHighPurityTracks)/double(nTracks); |
657 |
+ |
if(verbose_) |
658 |
+ |
cout<<"fraction of HighPurity track is "<<fHighPurity<<endl; |
659 |
+ |
|
660 |
+ |
if( (fHighPurity<0.2 && nTracks>10) || vertexColl->begin()->isFake()) { |
661 |
+ |
glob_runno_ = iEvent.id().run(); |
662 |
+ |
glob_evtno_ = iEvent.id().event(); |
663 |
+ |
glob_ls_ = iEvent.luminosityBlock(); |
664 |
+ |
glob_bx_ = iEvent.bunchCrossing(); |
665 |
+ |
return; |
666 |
+ |
} |
667 |
|
|
668 |
< |
//======================================================= |
669 |
< |
// Fill trackparameters of the input tracks to pvtx fitter |
670 |
< |
//======================================================= |
668 |
> |
//======================================================================== |
669 |
> |
// Step 2: Fill histograms for the splitting consistency checks |
670 |
> |
//======================================================================== |
671 |
> |
|
672 |
> |
// === Fill trackparameters of the input tracks to pvtx fitter |
673 |
|
edm::LogInfo("Debug")<<"[PVStudy] Start filling track parameters of the input tracks to pvtx fitter."<<endl; |
674 |
|
//fillTrackHisto(const reco::TrackCollection *trackColl, int datatype, const Point & bs) |
675 |
|
// datatype: unsplittracks (0); splittracks1 (1); splittracks2 (2); |
677 |
|
fillTrackHisto(splitTrackColl1, 1, beamSpot); |
678 |
|
fillTrackHisto(splitTrackColl2, 2, beamSpot); |
679 |
|
edm::LogInfo("Debug")<<"[PVStudy] End filling track parameters of the input tracks to pvtx fitter."<<endl; |
680 |
< |
|
681 |
< |
|
682 |
< |
//======================================================= |
741 |
< |
// Fill pixelVertices related histograms |
742 |
< |
//======================================================= |
743 |
< |
nrecPV_pxlpvtx_ = int (pixelVertexColl->size()); |
744 |
< |
if(pixelVertexColl->size()>0 && pixelVertexColl->begin()->isValid() && !(pixelVertexColl->begin()->isFake())) { |
745 |
< |
//fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple, const Point & bs) { |
746 |
< |
fillTrackHistoInPV(pixelVertexColl, 4, false, true, beamSpot); |
747 |
< |
h_pixvtx->Fill1d("dzErr_pxlpvtx", pixelVertexColl->begin()->zError()); |
748 |
< |
// Get the dZ error of the tracks in the leading pixelVertexColl |
749 |
< |
for(reco::Vertex::trackRef_iterator t = (pixelVertexColl->begin())->tracks_begin(); |
750 |
< |
t!= (pixelVertexColl->begin())->tracks_end(); t++) { |
751 |
< |
|
752 |
< |
if ( (**t).charge() < -1 || (**t).charge() > 1 ) { |
753 |
< |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
754 |
< |
} |
755 |
< |
else |
756 |
< |
h_pixvtx->Fill1d("trkdzErr_pxlpvtx", (**t).dzError()); |
757 |
< |
} |
758 |
< |
// Fill the track->dz() in the PV difference from first pixelpvtx |
759 |
< |
for(reco::Vertex::trackRef_iterator t = (vertexColl->begin())->tracks_begin(); |
760 |
< |
t!= (vertexColl->begin())->tracks_end(); t++) { |
761 |
< |
// illegal charge |
762 |
< |
if ( (**t).charge() < -1 || (**t).charge() > 1 ) { |
763 |
< |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
764 |
< |
} |
765 |
< |
else { |
766 |
< |
if(pixelVertexColl->size()>0) { |
767 |
< |
h_pixvtx->Fill1d("trkdz_pxlpvtxdz", (**t).dz() - pixelVertexColl->begin()->z()); |
768 |
< |
h_pixvtx->Fill1d("trkdz_pxlpvtxdz_pxlpvtxdzerr", fabs((**t).dz() - pixelVertexColl->begin()->z())/pixelVertexColl->begin()->zError()); |
769 |
< |
h_pixvtx->Fill1d("trkdz_pxlpvtxdz_trkdzerr", fabs((**t).dz() - pixelVertexColl->begin()->z())/(**t).dzError()); |
770 |
< |
h_pixvtx->Fill1d("trkdzErr_pvtx", (**t).dzError()); |
771 |
< |
} |
772 |
< |
} |
773 |
< |
} |
774 |
< |
} |
775 |
< |
|
776 |
< |
//======================================================= |
777 |
< |
// Fill number of reconstructed vertices |
778 |
< |
//======================================================= |
779 |
< |
|
680 |
> |
|
681 |
> |
|
682 |
> |
// ==== Fill number of reconstructed vertices |
683 |
|
edm::LogInfo("Debug")<<"[PVStudy] Printing vertexCollection: "<<endl; |
684 |
|
edm::LogInfo("Debug")<<"[PVStudy] Printing splitVertexCollection1: "<<endl; |
685 |
|
edm::LogInfo("Debug")<<"[PVStudy] Printing splitVertexCollection2: "<<endl; |
689 |
|
printRecVtxs(splitVertexCollection1Handle); |
690 |
|
printRecVtxs(splitVertexCollection2Handle); |
691 |
|
} |
692 |
+ |
|
693 |
|
nrecPV_ = int (vertexColl->size()); |
694 |
|
nrecPV1_spl_ = int (splitVertexColl1->size()); |
695 |
|
nrecPV2_spl_ = int (splitVertexColl2->size()); |
700 |
|
h_misc->Fill1d("nrecPVDiff", double(nrecPV1_spl_)-double(nrecPV2_spl_)); |
701 |
|
|
702 |
|
|
703 |
< |
//================================================================= |
704 |
< |
// Fill track parameter ntuple/hist for tracks used in recoVertices |
801 |
< |
//================================================================= |
802 |
< |
|
803 |
< |
//fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple, const Point & bs) { |
703 |
> |
// ======= Fill track parameter ntuple/hist for tracks used in recoVertices |
704 |
> |
//fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, const Point & bs) { |
705 |
|
if(vertexColl->size() > 0 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake())) |
706 |
< |
fillTrackHistoInPV(vertexColl, 0, true, true, beamSpot); |
806 |
< |
|
807 |
< |
if(splitVertexColl1->size() > 0 && splitVertexColl1->begin()->isValid() && !(splitVertexColl1->begin()->isFake())) |
808 |
< |
fillTrackHistoInPV(splitVertexColl1, 1, false, true, beamSpot); |
809 |
< |
|
810 |
< |
if(splitVertexColl2->size() > 0 && splitVertexColl2->begin()->isValid() && !(splitVertexColl2->begin()->isFake())) |
811 |
< |
fillTrackHistoInPV(splitVertexColl2, 2, false, true, beamSpot); |
812 |
< |
|
813 |
< |
//======================================================= |
814 |
< |
// Compare offlinePrimaryVertices with pixelVertices |
815 |
< |
//======================================================= |
816 |
< |
if( (pixelVertexColl->size()>0 && pixelVertexColl->begin()->isValid() && !(pixelVertexColl->begin()->isFake())) |
817 |
< |
&& (vertexColl->size()>0 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake())) ) { |
818 |
< |
h_pixvtx->Fill1d("nrecPV_minus_nrecPxlPV", double (nrecPV_ - nrecPV_pxlpvtx_)); |
819 |
< |
// difference in reconstructed position of the leading pvtx |
820 |
< |
//edm::LogInfo("Debug")<<"recx_[0] = "<< recx_[0] << "recx_pxlpvtx_[0] = "<< recx_pxlpvtx_[0]<<endl; |
821 |
< |
//edm::LogInfo("Debug")<<"recy_[0] = "<< recy_[0] << "recy_pxlpvtx_[0] = "<< recy_pxlpvtx_[0]<<endl; |
822 |
< |
h_pixvtx->Fill1d("recxPV_minus_recxPxlPV", recx_[0] - recx_pxlpvtx_[0]); |
823 |
< |
h_pixvtx->Fill1d("recyPV_minus_recyPxlPV", recy_[0] - recy_pxlpvtx_[0]); |
824 |
< |
h_pixvtx->Fill1d("reczPV_minus_reczPxlPV", recz_[0] - recz_pxlpvtx_[0]); |
825 |
< |
} |
826 |
< |
|
827 |
< |
|
828 |
< |
//========================================================== |
829 |
< |
// Fill secondary/primary min separations for multi-vertices |
830 |
< |
//========================================================== |
706 |
> |
fillTrackHistoInPV(vertexColl, 0, beamSpot); |
707 |
|
|
708 |
+ |
// ======= Fill secondary/primary min separations for multi-vertices |
709 |
|
if(vertexColl->size()>1 && vertexColl->begin()->isValid() && !(vertexColl->begin()->isFake()) ) { |
833 |
– |
|
710 |
|
double min_xsep = 9999.0; |
711 |
|
double min_ysep = 9999.0; |
712 |
|
double min_zsep = 9999.0; |
763 |
|
min_zsep_ = min_zsep; |
764 |
|
min_ntrksep_ = min_ntrksep; |
765 |
|
min_sumpt2sep_ = min_sumpt2sep; |
890 |
– |
|
891 |
– |
|
766 |
|
} // End of if(vertexColl->size()>1) { |
767 |
|
|
768 |
|
edm::LogInfo("Debug")<<"[PVStudy] End filling pvtx parameters."<<endl; |
769 |
|
|
770 |
< |
|
771 |
< |
|
772 |
< |
//======================================================= |
773 |
< |
// PrimaryVertex Matching |
774 |
< |
// 1. In z |z1-z2|< zsigncut * max(sigmaz1, sigmaz2) |
775 |
< |
// 2. |(nTrkPV1 - nTrkPV2)/(nTrkPV1+nTrkPV2)|<ntrkdiffcut_ |
902 |
< |
// Assume the first match is the primary vertex, |
903 |
< |
// since vertexColl are sorted in decreasing order of sum(pT) |
904 |
< |
//======================================================= |
770 |
> |
//======================================================================== |
771 |
> |
// Step 3: PrimaryVertex Matching |
772 |
> |
// 1. In z |z1-z2|< zsigncut * max(sigmaz1, sigmaz2) |
773 |
> |
// 2. |(nTrkPV1 - nTrkPV2)/(nTrkPV1+nTrkPV2)|<ntrkdiffcut_ |
774 |
> |
// 3. The first match is the primary vertex, which has largest sum(pT^2) |
775 |
> |
//======================================================================== |
776 |
|
|
777 |
|
edm::LogInfo("Debug")<<"[PVStudy] matching pvtxs "<<endl; |
778 |
|
reco::VertexCollection s_empty_matchedVertexColl1; |
826 |
|
} |
827 |
|
} |
828 |
|
edm::LogInfo("Debug")<<"[PVStudy] End matching pvtxs"<<endl; |
829 |
< |
|
830 |
< |
//======================================================= |
831 |
< |
// Analyze matchedVertexColls |
832 |
< |
//======================================================= |
829 |
> |
|
830 |
> |
//======================================================================== |
831 |
> |
// Step 4: Analyze matchedVertexColls |
832 |
> |
//======================================================================== |
833 |
|
edm::LogInfo("Debug")<<"[PVStudy] Begin analyzing the matchedVertexColl with size = "<< matchedVertexColl1->size()<< endl; |
834 |
|
|
835 |
+ |
// ==== If it is MC, analyze the res/pull of the unsplit vertexColl first |
836 |
+ |
if(MC){ |
837 |
+ |
// make a list of primary vertices: |
838 |
+ |
std::vector<simPrimaryVertex> simpv; |
839 |
+ |
simpv=getSimPVs(evtMC,""); |
840 |
+ |
// simpv=getSimPVs(evtMC, simVtxs, simTrks); |
841 |
+ |
h_gen->Fill1d("nsimPV", simpv.size()); |
842 |
+ |
nsimPV_ = simpv.size(); |
843 |
+ |
int isimpv = 0; |
844 |
+ |
for(std::vector<simPrimaryVertex>::iterator vsim=simpv.begin(); |
845 |
+ |
vsim!=simpv.end(); vsim++, isimpv++){ |
846 |
+ |
nsimTrkPV_[isimpv] =vsim->nGenTrk; |
847 |
+ |
simx_[isimpv] = vsim->x; |
848 |
+ |
simy_[isimpv] = vsim->y; |
849 |
+ |
simz_[isimpv] = vsim->y; |
850 |
+ |
simptsq_[isimpv] = vsim->ptsq; |
851 |
+ |
//fillMCHisto((std::vector<simPrimaryVertex>::iterator vsim, int isimpv, const reco::VertexCollection *vtxColl, int datatype) { |
852 |
+ |
fillMCHisto(vsim, isimpv, vertexColl, 3); |
853 |
+ |
} |
854 |
+ |
} |
855 |
+ |
|
856 |
|
// Compare the reconstructed vertex position and calculate resolution/pulls |
857 |
|
if(matchedVertexColl1->size() != 0 && matchedVertexColl2->size() != 0) { |
858 |
< |
//fillPvtxHisto(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple) |
859 |
< |
fillTrackHistoInPV(matchedVertexColl1, 1, true, false, beamSpot); |
860 |
< |
fillTrackHistoInPV(matchedVertexColl2, 2, true, false, beamSpot); |
858 |
> |
// ==== Analyze the splitted vtx using MC method |
859 |
> |
if(MC){ |
860 |
> |
// make a list of primary vertices: |
861 |
> |
std::vector<simPrimaryVertex> simpv; |
862 |
> |
simpv=getSimPVs(evtMC,""); |
863 |
> |
// simpv=getSimPVs(evtMC, simVtxs, simTrks); |
864 |
> |
int isimpv = 0; |
865 |
> |
for(std::vector<simPrimaryVertex>::iterator vsim=simpv.begin(); |
866 |
> |
vsim!=simpv.end(); vsim++, isimpv++){ |
867 |
> |
//fillMCHisto((std::vector<simPrimaryVertex>::iterator vsim, int isimpv, const reco::VertexCollection *vtxColl, int datatype) { |
868 |
> |
fillMCHisto(vsim, isimpv, matchedVertexColl1, 1); |
869 |
> |
fillMCHisto(vsim, isimpv, matchedVertexColl2, 2); |
870 |
> |
} |
871 |
> |
} |
872 |
> |
|
873 |
> |
// ==== Analyze res/pull two-vertex method |
874 |
> |
//fillPvtxHisto(const reco::VertexCollection *vertexColl, int datatype) |
875 |
> |
fillTrackHistoInPV(matchedVertexColl1, 1, beamSpot); |
876 |
> |
fillTrackHistoInPV(matchedVertexColl2, 2, beamSpot); |
877 |
|
|
878 |
|
reco::VertexCollection::const_iterator v1; |
879 |
|
reco::VertexCollection::const_iterator v2; |
880 |
+ |
|
881 |
|
nrecPV_twovtx_ = 0; |
882 |
|
for(v1 = matchedVertexColl1->begin(), v2 = matchedVertexColl2->begin(); |
883 |
|
v1!=matchedVertexColl1->end(), v2 != matchedVertexColl2->end(); |
908 |
|
fillHisto(res(fres_[0], fres_[1], fres_[2]), |
909 |
|
error(ferror_[0], ferror_[1], ferror_[2]), |
910 |
|
fntrk_,0); |
911 |
+ |
|
912 |
|
// Fill the ntuple variables |
1003 |
– |
nTrkPV1_spl_twovtx_[nrecPV_twovtx_] = v1->tracksSize(); |
1004 |
– |
nTrkPV2_spl_twovtx_[nrecPV_twovtx_] = v2->tracksSize(); |
1005 |
– |
sumptsq1_spl_twovtx_[nrecPV_twovtx_] = sumPtSquared(*v1); |
1006 |
– |
sumptsq2_spl_twovtx_[nrecPV_twovtx_] = sumPtSquared(*v2); |
913 |
|
nTrkPV_twovtx_[nrecPV_twovtx_] = fntrk_; |
914 |
|
deltax_twovtx_[nrecPV_twovtx_] = fres_[0]; |
915 |
|
deltay_twovtx_[nrecPV_twovtx_] = fres_[1]; |
920 |
|
pullx_twovtx_[nrecPV_twovtx_] = fres_[0]/ferror_[0]; |
921 |
|
pully_twovtx_[nrecPV_twovtx_] = fres_[1]/ferror_[1]; |
922 |
|
pullz_twovtx_[nrecPV_twovtx_] = fres_[2]/ferror_[2]; |
1017 |
– |
nrecPV_twovtx_++; |
1018 |
– |
} // End of analyzing res/pull |
923 |
|
|
924 |
< |
//======================================================= |
925 |
< |
// Fill simulated vertices |
1022 |
< |
//======================================================= |
1023 |
< |
if (!realData_) { |
1024 |
< |
bool MC=false; |
1025 |
< |
Handle<HepMCProduct> evtMC; |
1026 |
< |
iEvent.getByLabel("generator",evtMC); |
1027 |
< |
if (!evtMC.isValid()) { |
1028 |
< |
MC=false; |
1029 |
< |
edm::LogInfo("Debug") << "[PVStudy] no HepMCProduct found"<< endl; |
1030 |
< |
} else { |
1031 |
< |
edm::LogInfo("Debug") << "[PVStudy] generator HepMCProduct found"<< endl; |
1032 |
< |
MC=true; |
1033 |
< |
} |
924 |
> |
|
925 |
> |
//SplittedVertex |
926 |
|
|
927 |
< |
if(MC){ |
928 |
< |
// make a list of primary vertices: |
929 |
< |
std::vector<simPrimaryVertex> simpv; |
930 |
< |
simpv=getSimPVs(evtMC,""); |
931 |
< |
// simpv=getSimPVs(evtMC, simVtxs, simTrks); |
932 |
< |
h_gen->Fill1d("nsimPV", simpv.size()); |
933 |
< |
nsimPV_ = simpv.size(); |
934 |
< |
|
935 |
< |
int isimpv = 0; |
936 |
< |
for(std::vector<simPrimaryVertex>::iterator vsim=simpv.begin(); |
937 |
< |
vsim!=simpv.end(); vsim++, isimpv++){ |
938 |
< |
nsimTrkPV_[isimpv] =vsim->nGenTrk; |
939 |
< |
simx_[isimpv] = vsim->x; |
940 |
< |
simy_[isimpv] = vsim->y; |
941 |
< |
simz_[isimpv] = vsim->y; |
942 |
< |
simptsq_[isimpv] = vsim->ptsq; |
943 |
< |
|
1052 |
< |
//fillMCHisto((std::vector<simPrimaryVertex>::iterator vsim, int isimpv, const reco::VertexCollection *vtxColl, int datatype) { |
1053 |
< |
fillMCHisto(vsim, isimpv, splitVertexColl1, 1); |
1054 |
< |
fillMCHisto(vsim, isimpv, splitVertexColl2, 2); |
1055 |
< |
fillMCHisto(vsim, isimpv, vertexColl, 3); |
1056 |
< |
} |
1057 |
< |
} // End of for(std::vector<simPrimaryVertex>::iterator vsim=simpv.begin() |
1058 |
< |
} // End of if (!realData_) { |
927 |
> |
nTrkPV1_spl_twovtx_[nrecPV_twovtx_] = v1->tracksSize(); |
928 |
> |
ndofPV1_spl_twovtx_[nrecPV_twovtx_] = v1->ndof(); |
929 |
> |
normchi2PV1_spl_twovtx_[nrecPV_twovtx_] = v1->normalizedChi2(); |
930 |
> |
avgPtPV1_spl_twovtx_[nrecPV_twovtx_] = avgPtRecVtx(*v1); |
931 |
> |
errx1_spl_twovtx_[nrecPV_twovtx_] = v1->xError(); |
932 |
> |
erry1_spl_twovtx_[nrecPV_twovtx_] = v1->yError(); |
933 |
> |
errz1_spl_twovtx_[nrecPV_twovtx_] = v1->zError(); |
934 |
> |
|
935 |
> |
nTrkPV2_spl_twovtx_[nrecPV_twovtx_] = v2->tracksSize(); |
936 |
> |
ndofPV2_spl_twovtx_[nrecPV_twovtx_] = v2->ndof(); |
937 |
> |
normchi2PV2_spl_twovtx_[nrecPV_twovtx_] = v2->normalizedChi2(); |
938 |
> |
avgPtPV2_spl_twovtx_[nrecPV_twovtx_] = avgPtRecVtx(*v2); |
939 |
> |
errx2_spl_twovtx_[nrecPV_twovtx_] = v2->xError(); |
940 |
> |
erry2_spl_twovtx_[nrecPV_twovtx_] = v2->yError(); |
941 |
> |
errz2_spl_twovtx_[nrecPV_twovtx_] = v2->zError(); |
942 |
> |
|
943 |
> |
nrecPV_twovtx_++; |
944 |
|
|
945 |
+ |
|
946 |
+ |
|
947 |
+ |
// Print some information of the two tracks events |
948 |
+ |
if(verbose_ && fntrk_ < 4) { |
949 |
+ |
cout<<"Printing matchedVertexColl1 for low ntrack bins"<<endl; |
950 |
+ |
printRecVtx(*v1); |
951 |
+ |
cout<<"Printing matchedVertexColl1 for low ntrack bins"<<endl; |
952 |
+ |
printRecVtx(*v2); |
953 |
+ |
} |
954 |
+ |
} // End of analyzing res/pull |
955 |
|
} // End of if(matchedVertexColl1->size() == 1 && matchedVertexColl2->size() == 1 ) { |
956 |
|
else |
957 |
|
edm::LogInfo("Debug")<<"[PVStudy] WARNING: Cannot find matching pvtxs"<<endl; |
959 |
|
edm::LogInfo("Debug")<<"[PVStudy] End analyzing the matchedVertexColl"<<endl; |
960 |
|
|
961 |
|
|
962 |
< |
|
963 |
< |
// Finally fill the ftree_ |
962 |
> |
//======================================================================== |
963 |
> |
// Step 5: Fill ntuple if saventuple_ is on |
964 |
> |
//======================================================================== |
965 |
|
if(saventuple_) { |
966 |
|
ftree_->Fill(); |
967 |
|
pvtxtree_->Fill(); |
1028 |
|
{ |
1029 |
|
if(vrec.isFake() || !vrec.isValid()) return false; |
1030 |
|
else |
1031 |
< |
return (fabs(vsim.z-vrec.z())<0.0500); // =500um |
1031 |
> |
return true; |
1032 |
> |
//return (fabs(vsim.z-vrec.z())<0.0500); // =500um // remove this hard cut |
1033 |
|
} |
1034 |
|
|
1035 |
|
// Match two reconstructed vertices |
1096 |
|
} |
1097 |
|
} |
1098 |
|
|
1099 |
+ |
void PVStudy::printRecVtx(const reco::Vertex & v){ |
1100 |
+ |
|
1101 |
+ |
cout << "#trk " << std::setw(3) << v.tracksSize() |
1102 |
+ |
<< " chi2 " << std::setw(4) << v.chi2() |
1103 |
+ |
<< " ndof " << std::setw(3) << v.ndof() << endl |
1104 |
+ |
<< " x " << std::setw(8) <<std::fixed << std::setprecision(4) << v.x() |
1105 |
+ |
<< " dx " << std::setw(8) << v.xError()<< endl |
1106 |
+ |
<< " y " << std::setw(8) << v.y() |
1107 |
+ |
<< " dy " << std::setw(8) << v.yError()<< endl |
1108 |
+ |
<< " z " << std::setw(8) << v.z() |
1109 |
+ |
<< " dz " << std::setw(8) << v.zError() |
1110 |
+ |
<< endl; |
1111 |
+ |
} |
1112 |
+ |
|
1113 |
|
void PVStudy::printSimTrks(const Handle<SimTrackContainer> simTrks){ |
1114 |
|
cout << " simTrks type, (momentum), vertIndex, genpartIndex" << endl; |
1115 |
|
int i=1; |
1281 |
|
} |
1282 |
|
} |
1283 |
|
|
1284 |
< |
void PVStudy::fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, bool fillHisto, bool fillNtuple, const Point & bs) { |
1284 |
> |
void PVStudy::fillTrackHistoInPV(const reco::VertexCollection *vertexColl, int datatype, const Point & bs) { |
1285 |
|
string suffix; |
1286 |
|
suffix = ""; // for unsplit tracks |
1287 |
|
if(datatype == 1) suffix = "1_spl"; // for splittracks 1 |
1289 |
|
int ivtx = 0; |
1290 |
|
for(reco::VertexCollection::const_iterator v=vertexColl->begin(); |
1291 |
|
v!=vertexColl->end(); ++v, ++ivtx) { |
1292 |
< |
if(v->isFake()) continue; |
1293 |
< |
if(fillNtuple) { |
1294 |
< |
if(datatype == 0) { |
1295 |
< |
nTrkPV_[ivtx] = v->tracksSize(); |
1296 |
< |
sumptsq_[ivtx] = sumPtSquared(*v); |
1297 |
< |
isValid_[ivtx] = v->isValid(); |
1298 |
< |
isFake_[ivtx] = v->isFake(); |
1299 |
< |
recx_[ivtx] = v->x(); |
1300 |
< |
recy_[ivtx] = v->y(); |
1301 |
< |
recz_[ivtx] = v->z(); |
1302 |
< |
recx_err_[ivtx] = v->xError(); |
1303 |
< |
recy_err_[ivtx] = v->yError(); |
1304 |
< |
recz_err_[ivtx] = v->zError(); |
1305 |
< |
} |
1306 |
< |
if(datatype == 1) { |
1307 |
< |
nTrkPV1_spl_[ivtx] = v->tracksSize(); |
1308 |
< |
sumptsq1_spl_[ivtx] = sumPtSquared(*v); |
1309 |
< |
isValid1_spl_[ivtx] = v->isValid(); |
1310 |
< |
isFake1_spl_[ivtx] = v->isFake(); |
1311 |
< |
recx1_spl_[ivtx] = v->x(); |
1312 |
< |
recy1_spl_[ivtx] = v->y(); |
1402 |
< |
recz1_spl_[ivtx] = v->z(); |
1403 |
< |
recx1_err_spl_[ivtx] = v->xError(); |
1404 |
< |
recy1_err_spl_[ivtx] = v->yError(); |
1405 |
< |
recz1_err_spl_[ivtx] = v->zError(); |
1406 |
< |
|
1407 |
< |
} |
1408 |
< |
if(datatype == 2) { |
1409 |
< |
nTrkPV2_spl_[ivtx] = v->tracksSize(); |
1410 |
< |
sumptsq2_spl_[ivtx] = sumPtSquared(*v); |
1411 |
< |
isValid2_spl_[ivtx] = v->isValid(); |
1412 |
< |
isFake2_spl_[ivtx] = v->isFake(); |
1413 |
< |
recx2_spl_[ivtx] = v->x(); |
1414 |
< |
recy2_spl_[ivtx] = v->y(); |
1415 |
< |
recz2_spl_[ivtx] = v->z(); |
1416 |
< |
recx2_err_spl_[ivtx] = v->xError(); |
1417 |
< |
recy2_err_spl_[ivtx] = v->yError(); |
1418 |
< |
recz2_err_spl_[ivtx] = v->zError(); |
1419 |
< |
} |
1420 |
< |
if(datatype == 4) { // for pixelVertices |
1421 |
< |
nTrkPV_pxlpvtx_[ivtx] = v->tracksSize(); |
1422 |
< |
sumptsq_pxlpvtx_[ivtx] = sumPtSquared(*v); |
1423 |
< |
isValid_pxlpvtx_[ivtx] = v->isValid(); |
1424 |
< |
isFake_pxlpvtx_[ivtx] = v->isFake(); |
1425 |
< |
recx_pxlpvtx_[ivtx] = v->x(); |
1426 |
< |
recy_pxlpvtx_[ivtx] = v->y(); |
1427 |
< |
recz_pxlpvtx_[ivtx] = v->z(); |
1428 |
< |
recx_err_pxlpvtx_[ivtx] = v->xError(); |
1429 |
< |
recy_err_pxlpvtx_[ivtx] = v->yError(); |
1430 |
< |
recz_err_pxlpvtx_[ivtx] = v->zError(); |
1431 |
< |
} |
1432 |
< |
} |
1433 |
< |
} |
1434 |
< |
// For histogram only fill the leading pvtx parameters |
1435 |
< |
if(fillHisto) { |
1436 |
< |
h_pvtrk->Fill1d("nTrkPV"+suffix, vertexColl->begin()->tracksSize()); |
1437 |
< |
h_pvtrk->Fill1d("ndofPV"+suffix, vertexColl->begin()->ndof()); |
1438 |
< |
int nHWTrkPV_ = 0; |
1439 |
< |
|
1440 |
< |
try { |
1441 |
< |
for(reco::Vertex::trackRef_iterator t = vertexColl->begin()->tracks_begin(); |
1442 |
< |
t!=vertexColl->begin()->tracks_end(); t++) { |
1443 |
< |
// illegal charge |
1444 |
< |
if ( (**t).charge() < -1 || (**t).charge() > 1 ) { |
1445 |
< |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
1446 |
< |
} |
1447 |
< |
else { |
1448 |
< |
h_pvtrk->Fill1d("trkPtPV"+suffix, (**t).pt()); |
1449 |
< |
h_pvtrk->Fill1d("trkDxyPV"+suffix, (**t).dxy()); |
1450 |
< |
h_pvtrk->Fill1d("trkDxyCorrPV"+suffix, (**t).dxy(bs)); |
1451 |
< |
h_pvtrk->Fill1d("trkDzPV"+suffix, (**t).dz()); |
1452 |
< |
h_pvtrk->Fill1d("trkEtaPV"+suffix, (**t).eta()); |
1453 |
< |
h_pvtrk->Fill1d("trkPhiPV"+suffix, (**t).phi()); |
1454 |
< |
|
1455 |
< |
if(vertexColl->begin()->trackWeight(*t) > 0.5) |
1456 |
< |
nHWTrkPV_++; |
1292 |
> |
if(!v->isFake()) { |
1293 |
> |
h_pvtrk->Fill1d("nTrkPV"+suffix, vertexColl->begin()->tracksSize()); |
1294 |
> |
h_pvtrk->Fill1d("ndofPV"+suffix, vertexColl->begin()->ndof()); |
1295 |
> |
int nHWTrkPV_ = 0; |
1296 |
> |
try { |
1297 |
> |
for(reco::Vertex::trackRef_iterator t = vertexColl->begin()->tracks_begin(); |
1298 |
> |
t!=vertexColl->begin()->tracks_end(); t++) { |
1299 |
> |
// illegal charge |
1300 |
> |
if ( (**t).charge() < -1 || (**t).charge() > 1 ) { |
1301 |
> |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
1302 |
> |
} |
1303 |
> |
else { |
1304 |
> |
h_pvtrk->Fill1d("trkPtPV"+suffix, (**t).pt()); |
1305 |
> |
h_pvtrk->Fill1d("trkDxyPV"+suffix, (**t).dxy()); |
1306 |
> |
h_pvtrk->Fill1d("trkDxyCorrPV"+suffix, (**t).dxy(bs)); |
1307 |
> |
h_pvtrk->Fill1d("trkDzPV"+suffix, (**t).dz()); |
1308 |
> |
h_pvtrk->Fill1d("trkEtaPV"+suffix, (**t).eta()); |
1309 |
> |
h_pvtrk->Fill1d("trkPhiPV"+suffix, (**t).phi()); |
1310 |
> |
if(vertexColl->begin()->trackWeight(*t) > 0.5) |
1311 |
> |
nHWTrkPV_++; |
1312 |
> |
} |
1313 |
|
} |
1314 |
|
} |
1315 |
+ |
catch (...) { |
1316 |
+ |
// exception thrown when trying to use linked track |
1317 |
+ |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
1318 |
+ |
} |
1319 |
+ |
h_pvtrk->Fill1d("nHWTrkPV"+suffix, nHWTrkPV_); |
1320 |
|
} |
1460 |
– |
catch (...) { |
1461 |
– |
// exception thrown when trying to use linked track |
1462 |
– |
// h_pvtrk->Fill1d("trkPtPV", 0.); |
1463 |
– |
} |
1464 |
– |
h_pvtrk->Fill1d("nHWTrkPV"+suffix, nHWTrkPV_); |
1321 |
|
} |
1466 |
– |
|
1467 |
– |
|
1322 |
|
} |
1323 |
|
|
1324 |
|
void PVStudy::fillMCHisto(std::vector<simPrimaryVertex>::iterator vsim, int isimpv, const reco::VertexCollection *vtxColl, int datatype) |
1344 |
|
|
1345 |
|
for(reco::VertexCollection::const_iterator vrec=vtxColl->begin(); |
1346 |
|
vrec!=vtxColl->end(); ++vrec){ |
1493 |
– |
int nrectrk = vrec->tracksSize(); |
1347 |
|
vsim->recVtx=NULL; |
1495 |
– |
|
1348 |
|
edm::LogInfo("Debug") << "[fillMCHisto] sim primary vertex x = " << vsim->x << "; y = " << vsim->y << "; z = " << vsim->z << endl; |
1349 |
|
edm::LogInfo("Debug") << "[fillMCHisto] Is matched? " << (matchVertex(*vsim,*vrec)?"Yes":"No") << endl; |
1350 |
|
|
1360 |
|
|
1361 |
|
double fres_mct[3]; |
1362 |
|
double ferror_mct[3]; |
1363 |
< |
|
1363 |
> |
int nrectrk = int(vsim->recVtx->tracksSize()); |
1364 |
> |
double ndofPV = vsim->recVtx->ndof(); |
1365 |
> |
double normchi2PV = vsim->recVtx->normalizedChi2(); |
1366 |
> |
double avgPtPV = avgPtRecVtx(*(vsim->recVtx)); |
1367 |
> |
|
1368 |
|
fres_mct[0] = vsim->recVtx->x()-vsim->x; |
1369 |
|
fres_mct[1] = vsim->recVtx->y()-vsim->y; |
1370 |
|
fres_mct[2] = vsim->recVtx->z()-vsim->z; |
1390 |
|
nrectrk, datatype); |
1391 |
|
|
1392 |
|
if(saventuple_) { |
1393 |
< |
//Fill the values for variables in ftree_ |
1393 |
> |
//Fill the values for variables in pvtxtree_ |
1394 |
|
if(datatype == 1) { |
1395 |
|
nTrkPV_spl1_mct_[nrecPV_spl1_mct_] = nrectrk; |
1396 |
+ |
ndofPV_spl1_mct_[nrecPV_spl1_mct_] = ndofPV; |
1397 |
+ |
normchi2PV_spl1_mct_[nrecPV_spl1_mct_] = normchi2PV; |
1398 |
+ |
avgPtPV_spl1_mct_[nrecPV_spl1_mct_] = avgPtPV; |
1399 |
|
deltax_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[0]; |
1400 |
< |
deltay_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[0]; |
1401 |
< |
deltaz_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[0]; |
1400 |
> |
deltay_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[1]; |
1401 |
> |
deltaz_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[2]; |
1402 |
|
pullx_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[0]/ferror_mct[0]; |
1403 |
|
pully_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[1]/ferror_mct[1]; |
1404 |
|
pullz_spl1_mct_[nrecPV_spl1_mct_] = fres_mct[2]/ferror_mct[2]; |
1405 |
+ |
errx_spl1_mct_[nrecPV_spl1_mct_] = ferror_mct[0]; |
1406 |
+ |
erry_spl1_mct_[nrecPV_spl1_mct_] = ferror_mct[1]; |
1407 |
+ |
errz_spl1_mct_[nrecPV_spl1_mct_] = ferror_mct[2]; |
1408 |
|
nrecPV_spl1_mct_++; |
1409 |
|
} |
1548 |
– |
|
1410 |
|
if(datatype == 2) { |
1411 |
< |
nTrkPV_spl2_mct_[nrecPV_spl2_mct_] = nrectrk; |
1411 |
> |
nTrkPV_spl2_mct_[nrecPV_spl2_mct_] = nrectrk; |
1412 |
> |
ndofPV_spl2_mct_[nrecPV_spl2_mct_] = ndofPV; |
1413 |
> |
normchi2PV_spl2_mct_[nrecPV_spl2_mct_] = normchi2PV; |
1414 |
> |
avgPtPV_spl2_mct_[nrecPV_spl2_mct_] = avgPtPV; |
1415 |
|
deltax_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[0]; |
1416 |
< |
deltay_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[0]; |
1417 |
< |
deltaz_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[0]; |
1416 |
> |
deltay_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[1]; |
1417 |
> |
deltaz_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[2]; |
1418 |
|
pullx_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[0]/ferror_mct[0]; |
1419 |
|
pully_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[1]/ferror_mct[1]; |
1420 |
< |
pullz_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[2]/ferror_mct[2]; |
1420 |
> |
pullz_spl2_mct_[nrecPV_spl2_mct_] = fres_mct[2]/ferror_mct[2]; |
1421 |
> |
errx_spl2_mct_[nrecPV_spl2_mct_] = ferror_mct[0]; |
1422 |
> |
erry_spl2_mct_[nrecPV_spl2_mct_] = ferror_mct[1]; |
1423 |
> |
errz_spl2_mct_[nrecPV_spl2_mct_] = ferror_mct[2]; |
1424 |
|
nrecPV_spl2_mct_++; |
1425 |
|
} |
1426 |
|
if(datatype == 3) { |
1427 |
< |
nTrkPV_mct_[nrecPV_mct_] = nrectrk; |
1427 |
> |
nTrkPV_mct_[nrecPV_mct_] = nrectrk; |
1428 |
> |
ndofPV_mct_[nrecPV_mct_] = ndofPV; |
1429 |
> |
normchi2PV_mct_[nrecPV_mct_] = normchi2PV; |
1430 |
> |
avgPtPV_mct_[nrecPV_mct_] = avgPtPV; |
1431 |
|
deltax_mct_[nrecPV_mct_] = fres_mct[0]; |
1432 |
|
deltay_mct_[nrecPV_mct_] = fres_mct[0]; |
1433 |
|
deltaz_mct_[nrecPV_mct_] = fres_mct[0]; |
1434 |
|
pullx_mct_[nrecPV_mct_] = fres_mct[0]/ferror_mct[0]; |
1435 |
|
pully_mct_[nrecPV_mct_] = fres_mct[1]/ferror_mct[1]; |
1436 |
< |
pullz_mct_[nrecPV_mct_] = fres_mct[2]/ferror_mct[2]; |
1436 |
> |
pullz_mct_[nrecPV_mct_] = fres_mct[2]/ferror_mct[2]; |
1437 |
> |
errx_mct_[nrecPV_mct_] = ferror_mct[0]; |
1438 |
> |
erry_mct_[nrecPV_mct_] = ferror_mct[1]; |
1439 |
> |
errz_mct_[nrecPV_mct_] = ferror_mct[2]; |
1440 |
|
nrecPV_mct_++; |
1441 |
|
} |
1442 |
|
} // End of if(saventuple_) { |
1448 |
|
} |
1449 |
|
|
1450 |
|
void PVStudy::SetVarToZero() { |
1578 |
– |
//Greg's variables |
1451 |
|
fntrk_ = 0; |
1452 |
|
//pvtx position (x,y,z) residual and error |
1453 |
|
for(int i = 0; i<3;i++) { |
1455 |
|
ferror_[i] = 0; |
1456 |
|
} |
1457 |
|
|
1458 |
< |
//Yanyan's variables |
1459 |
< |
// Number of reconstructed vertices |
1458 |
> |
// Event level information |
1459 |
> |
glob_runno_ = 0; |
1460 |
> |
glob_evtno_ = 0; |
1461 |
> |
glob_ls_ = 0; |
1462 |
> |
glob_bx_ = 0; |
1463 |
|
nrecPV_ = 0; |
1464 |
|
nrecPV1_spl_ = 0; |
1465 |
|
nrecPV2_spl_ = 0; |
1475 |
|
min_ntrksep_ = 9999.0; |
1476 |
|
min_sumpt2sep_ = -9999.0; |
1477 |
|
|
1478 |
< |
|
1478 |
> |
// Variables filled per Vertex |
1479 |
|
for (int i = 0; i < nMaxPVs_; i++) { |
1605 |
– |
// recoVertices with all tracks |
1606 |
– |
nTrkPV_[i] = 0; // Number of tracks in the pvtx |
1607 |
– |
sumptsq_[i] = 0; |
1608 |
– |
isValid_[i] = -1; |
1609 |
– |
isFake_[i] = -1; |
1610 |
– |
recx_[i] = 0; |
1611 |
– |
recy_[i] = 0; |
1612 |
– |
recz_[i] = 0; |
1613 |
– |
recx_err_[i] = 0; |
1614 |
– |
recy_err_[i] = 0; |
1615 |
– |
recz_err_[i] = 0; |
1616 |
– |
|
1617 |
– |
// recoVertices with splitTrack1 |
1618 |
– |
nTrkPV1_spl_[i] = 0; // Number of tracks in the pvtx |
1619 |
– |
sumptsq1_spl_[i] = 0; |
1620 |
– |
isValid1_spl_[i] = -1; |
1621 |
– |
isFake1_spl_[i] = -1; |
1622 |
– |
recx1_spl_[i] = 0; |
1623 |
– |
recy1_spl_[i] = 0; |
1624 |
– |
recz1_spl_[i] = 0; |
1625 |
– |
recx1_err_spl_[i] = 0; |
1626 |
– |
recy1_err_spl_[i] = 0; |
1627 |
– |
recz1_err_spl_[i] = 0; |
1628 |
– |
|
1629 |
– |
// recoVertices with splitTrack2 |
1630 |
– |
nTrkPV2_spl_[i] = 0; // Number of tracks in the pvtx |
1631 |
– |
sumptsq2_spl_[i] = 0; |
1632 |
– |
isValid2_spl_[i] = -1; |
1633 |
– |
isFake2_spl_[i] = -1; |
1634 |
– |
recx2_spl_[i] = 0; |
1635 |
– |
recy2_spl_[i] = 0; |
1636 |
– |
recz2_spl_[i] = 0; |
1637 |
– |
recx2_err_spl_[i] = 0; |
1638 |
– |
recy2_err_spl_[i] = 0; |
1639 |
– |
recz2_err_spl_[i] = 0; |
1640 |
– |
|
1641 |
– |
//pixelVertices |
1642 |
– |
nTrkPV_pxlpvtx_[i] = 0; // Number of tracks in the pvtx |
1643 |
– |
sumptsq_pxlpvtx_[i] = 0; |
1644 |
– |
isValid_pxlpvtx_[i] = -1; |
1645 |
– |
isFake_pxlpvtx_[i] = -1; |
1646 |
– |
recx_pxlpvtx_[i] = 0; |
1647 |
– |
recy_pxlpvtx_[i] = 0; |
1648 |
– |
recz_pxlpvtx_[i] = 0; |
1649 |
– |
recx_err_pxlpvtx_[i] = 0; |
1650 |
– |
recy_err_pxlpvtx_[i] = 0; |
1651 |
– |
recz_err_pxlpvtx_[i] = 0; |
1652 |
– |
|
1480 |
|
// matched two-vertices |
1481 |
|
nTrkPV1_spl_twovtx_[i] = 0; |
1482 |
+ |
ndofPV1_spl_twovtx_[i] = 0; |
1483 |
+ |
normchi2PV1_spl_twovtx_[i] = 0; |
1484 |
+ |
avgPtPV1_spl_twovtx_[i] = 0; |
1485 |
+ |
errx1_spl_twovtx_[i] = 0; |
1486 |
+ |
erry1_spl_twovtx_[i] = 0; |
1487 |
+ |
errz1_spl_twovtx_[i] = 0; |
1488 |
+ |
|
1489 |
|
nTrkPV2_spl_twovtx_[i] = 0; |
1490 |
< |
nTrkPV_twovtx_[i] = 0; |
1491 |
< |
sumptsq1_spl_twovtx_[i] = 0; // Sum of pT squred in the pvtx |
1492 |
< |
sumptsq2_spl_twovtx_[i] = 0; // Sum of pT squred in the pvtx |
1490 |
> |
ndofPV2_spl_twovtx_[i] = 0; |
1491 |
> |
normchi2PV2_spl_twovtx_[i] = 0; |
1492 |
> |
avgPtPV2_spl_twovtx_[i] = 0; |
1493 |
> |
errx2_spl_twovtx_[i] = 0; |
1494 |
> |
erry2_spl_twovtx_[i] = 0; |
1495 |
> |
errz2_spl_twovtx_[i] = 0; |
1496 |
> |
|
1497 |
|
deltax_twovtx_[i] = 0; |
1498 |
|
deltay_twovtx_[i] = 0; |
1499 |
|
deltaz_twovtx_[i] = 0; |
1518 |
|
pullx_mct_[i] = 0; |
1519 |
|
pully_mct_[i] = 0; |
1520 |
|
pullz_mct_[i] = 0; |
1521 |
< |
|
1521 |
> |
nTrkPV_mct_[i] = 0; |
1522 |
> |
ndofPV_mct_[i] = 0; |
1523 |
> |
normchi2PV_mct_[i] = 0; |
1524 |
> |
avgPtPV_mct_[i] = 0; |
1525 |
> |
errx_mct_[i] = 0; |
1526 |
> |
erry_mct_[i] = 0; |
1527 |
> |
errz_mct_[i] = 0; |
1528 |
> |
|
1529 |
|
deltax_spl1_mct_[i] = 0; |
1530 |
|
deltay_spl1_mct_[i] = 0; |
1531 |
|
deltaz_spl1_mct_[i] = 0; |
1532 |
|
pullx_spl1_mct_[i] = 0; |
1533 |
|
pully_spl1_mct_[i] = 0; |
1534 |
|
pullz_spl1_mct_[i] = 0; |
1535 |
< |
|
1535 |
> |
nTrkPV_spl1_mct_[i] = 0; |
1536 |
> |
ndofPV_spl1_mct_[i] = 0; |
1537 |
> |
normchi2PV_spl1_mct_[i] = 0; |
1538 |
> |
avgPtPV_spl1_mct_[i] = 0; |
1539 |
> |
errx_spl1_mct_[i] = 0; |
1540 |
> |
erry_spl1_mct_[i] = 0; |
1541 |
> |
errz_spl1_mct_[i] = 0; |
1542 |
> |
|
1543 |
|
deltax_spl2_mct_[i] = 0; |
1544 |
|
deltay_spl2_mct_[i] = 0; |
1545 |
|
deltaz_spl2_mct_[i] = 0; |
1546 |
|
pullx_spl2_mct_[i] = 0; |
1547 |
|
pully_spl2_mct_[i] = 0; |
1548 |
|
pullz_spl2_mct_[i] = 0; |
1549 |
+ |
nTrkPV_spl2_mct_[i] = 0; |
1550 |
+ |
ndofPV_spl2_mct_[i] = 0; |
1551 |
+ |
normchi2PV_spl2_mct_[i] = 0; |
1552 |
+ |
avgPtPV_spl2_mct_[i] = 0; |
1553 |
+ |
errx_spl2_mct_[i] = 0; |
1554 |
+ |
erry_spl2_mct_[i] = 0; |
1555 |
+ |
errz_spl2_mct_[i] = 0; |
1556 |
+ |
|
1557 |
|
} |
1558 |
|
} |
1559 |
|
|
1560 |
|
double PVStudy::sumPtSquared(const reco::Vertex & v) { |
1561 |
< |
double sum = 0.; |
1562 |
< |
double pT; |
1563 |
< |
for (reco::Vertex::trackRef_iterator it = v.tracks_begin(); it != v.tracks_end(); it++) { |
1564 |
< |
pT = (**it).pt(); |
1565 |
< |
sum += pT*pT; |
1561 |
> |
if(v.isFake() || (!v.isValid()) || v.tracksSize() == 0) return 0.0; |
1562 |
> |
else { |
1563 |
> |
double sum = 0.; |
1564 |
> |
double pT; |
1565 |
> |
for (reco::Vertex::trackRef_iterator it = v.tracks_begin(); it != v.tracks_end(); it++) { |
1566 |
> |
pT = (**it).pt(); |
1567 |
> |
sum += pT*pT; |
1568 |
> |
} |
1569 |
> |
return sum; |
1570 |
> |
} |
1571 |
> |
} |
1572 |
> |
|
1573 |
> |
|
1574 |
> |
|
1575 |
> |
double PVStudy::avgPtRecVtx(const reco::Vertex & v) { |
1576 |
> |
|
1577 |
> |
if(v.isFake() || !v.isValid() || v.tracksSize()==0 ) return 0; |
1578 |
> |
else { |
1579 |
> |
double sumpT = 0.; |
1580 |
> |
for (reco::Vertex::trackRef_iterator it = v.tracks_begin(); it != v.tracks_end(); it++) { |
1581 |
> |
sumpT += (**it).pt(); |
1582 |
> |
} |
1583 |
> |
return sumpT/double(v.tracksSize()); |
1584 |
|
} |
1707 |
– |
return sum; |
1585 |
|
} |
1586 |
|
|
1587 |
|
|