126 constexpr double c_cm_ns = CLHEP::c_light * CLHEP::ns / CLHEP::cm;
131 const auto&
tracks = *tracksH;
135 const auto& t0In = *t0H;
139 const auto& tmtdIn = *tmtdH;
143 const auto& sigmat0In = *sigmat0H;
147 const auto& sigmatmtdIn = *sigmatmtdH;
151 const auto& pathLengthIn = *pathLengthH;
155 const auto& pIn = *pH;
159 const auto& vtxs = *vtxsH;
162 std::vector<float> t0OutRaw;
163 std::vector<float> sigmat0OutRaw;
164 std::vector<float> t0safeOutRaw;
165 std::vector<float> sigmat0safeOutRaw;
166 std::vector<float> probPiOutRaw;
167 std::vector<float> probKOutRaw;
168 std::vector<float> probPOutRaw;
171 for (
unsigned int itrack = 0; itrack <
tracks.size(); ++itrack) {
174 float t0 = t0In[trackref];
176 float sigmat0safe = sigmat0In[trackref];
177 float sigmatmtd = sigmatmtdIn[trackref];
178 float sigmat0 = sigmatmtd;
186 double rsigmat = 1. / sigmatmtd;
190 int vtxidxmindz = -1;
191 int vtxidxminchisq = -1;
195 for (
unsigned int ivtx = 0; ivtx < vtxs.size(); ++ivtx) {
209 double dtsig = dt * rsigmat;
210 double chisq = dz * dz * rsigmazsq + dtsig * dtsig;
213 vtxidxminchisq = ivtx;
222 if (vtxidxmindz >= 0 && !(vtxs[vtxidxmindz].tError() > 0. && vtxs[vtxidxmindz].tError() <
vtxMaxSigmaT_)) {
223 vtxidx = vtxidxmindz;
224 }
else if (vtxidxminchisq >= 0) {
225 vtxidx = vtxidxminchisq;
226 }
else if (vtxidxmindz >= 0) {
227 vtxidx = vtxidxmindz;
232 if (vtxidx >= 0 && vtxs[vtxidx].tError() > 0. && vtxs[vtxidx].tError() <
vtxMaxSigmaT_) {
236 double dtnom =
std::abs(t0 - vtxnom.
t());
237 double dtsignom = dtnom * rsigmat;
238 double chisqnom = dznom * dznom * rsigmazsq + dtsignom * dtsignom;
245 sigmat0safe = sigmatmtd;
248 double tmtd = tmtdIn[trackref];
249 double pathlength = pathLengthIn[trackref];
250 double magp = pIn[trackref];
252 double gammasq_k = 1. + magp * magp / m_k / m_k;
253 double beta_k =
std::sqrt(1. - 1. / gammasq_k);
254 double t0_k = tmtd - pathlength / beta_k * c_inv;
256 double gammasq_p = 1. + magp * magp / m_p / m_p;
257 double beta_p =
std::sqrt(1. - 1. / gammasq_p);
258 double t0_p = tmtd - pathlength / beta_p * c_inv;
260 double chisqmin = chisqnom;
262 double chisqmin_pi = chisqnom;
276 double chisqdz = dz * dz * rsigmazsq;
279 double dtsig_k = dt_k * rsigmat;
280 double chisq_k = chisqdz + dtsig_k * dtsig_k;
283 chisqmin_k = chisq_k;
287 double dtsig_p = dt_p * rsigmat;
288 double chisq_p = chisqdz + dtsig_p * dtsig_p;
291 chisqmin_p = chisq_p;
298 sigmat0safe = sigmatmtd;
304 sigmat0safe = sigmatmtd;
310 double rawprob_pi =
exp(-0.5 * chisqmin_pi);
311 double rawprob_k =
exp(-0.5 * chisqmin_k);
312 double rawprob_p =
exp(-0.5 * chisqmin_p);
314 double normprob = 1. / (rawprob_pi + rawprob_k + rawprob_p);
316 prob_pi = rawprob_pi * normprob;
317 prob_k = rawprob_k * normprob;
318 prob_p = rawprob_p * normprob;
320 double prob_heavy = 1. - prob_pi;
328 t0OutRaw.push_back(t0);
329 sigmat0OutRaw.push_back(sigmat0);
330 t0safeOutRaw.push_back(t0safe);
331 sigmat0safeOutRaw.push_back(sigmat0safe);
332 probPiOutRaw.push_back(prob_pi);
333 probKOutRaw.push_back(prob_k);
334 probPOutRaw.push_back(prob_p);
edm::EDGetTokenT< reco::VertexCollection > vtxsToken_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
edm::EDGetTokenT< edm::ValueMap< float > > sigmat0Token_
edm::EDGetTokenT< edm::ValueMap< float > > pathLengthToken_
const Point & position() const
position
double maxDtSignificance_
Abs< T >::type abs(const T &t)
float trackWeight(const TREF &r) const
returns the weight with which a Track has contributed to the vertex-fit.
edm::EDGetTokenT< edm::ValueMap< float > > tmtdToken_
static char sigmat0safeName[]
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
double dzError() const
error on dz
edm::EDGetTokenT< edm::ValueMap< float > > t0Token_
edm::EDGetTokenT< edm::ValueMap< float > > sigmatmtdToken_
edm::EDGetTokenT< edm::ValueMap< float > > pToken_
static char sigmat0Name[]
edm::EDGetTokenT< reco::TrackCollection > tracksToken_
void fillValueMap(edm::Event &iEvent, const edm::Handle< H > &handle, const std::vector< T > &vec, const std::string &name) const
double tError() const
error on t
double t() const
t coordinate