00001 #include <algorithm>
00002 #include <iterator>
00003 #include <vector>
00004 #include <map>
00005 #include <set>
00006
00007 #include <Math/SVector.h>
00008 #include <Math/SMatrix.h>
00009 #include <Math/MatrixFunctions.h>
00010
00011 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
00012 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
00013 #include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h"
00014
00015 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00016
00017 #include "TrackingTools/TrajectoryParametrization/interface/GlobalTrajectoryParameters.h"
00018 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00019 #include "TrackingTools/TrajectoryState/interface/FreeTrajectoryState.h"
00020 #include "TrackingTools/GeomPropagators/interface/AnalyticalImpactPointExtrapolator.h"
00021 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00022 #include "TrackingTools/TransientTrack/interface/TransientTrackFromFTSFactory.h"
00023
00024 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
00025 #include "RecoVertex/VertexPrimitives/interface/VertexFitter.h"
00026 #include "RecoVertex/VertexPrimitives/interface/VertexState.h"
00027 #include "RecoVertex/VertexPrimitives/interface/CachingVertex.h"
00028 #include "RecoVertex/VertexPrimitives/interface/LinearizedTrackState.h"
00029 #include "RecoVertex/VertexPrimitives/interface/ConvertError.h"
00030 #include "RecoVertex/VertexPrimitives/interface/ConvertToFromReco.h"
00031 #include "RecoVertex/VertexTools/interface/LinearizedTrackStateFactory.h"
00032 #include "RecoVertex/VertexTools/interface/VertexTrackFactory.h"
00033 #include "RecoVertex/VertexTools/interface/VertexDistance3D.h"
00034 #include "RecoVertex/VertexTools/interface/GeometricAnnealing.h"
00035 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
00036 #include "RecoVertex/AdaptiveVertexFit/interface/AdaptiveVertexFitter.h"
00037
00038 #include "RecoVertex/GhostTrackFitter/interface/GhostTrack.h"
00039 #include "RecoVertex/GhostTrackFitter/interface/GhostTrackState.h"
00040 #include "RecoVertex/GhostTrackFitter/interface/GhostTrackFitter.h"
00041 #include "RecoVertex/GhostTrackFitter/interface/GhostTrackPrediction.h"
00042 #include "RecoVertex/GhostTrackFitter/interface/SequentialGhostTrackFitter.h"
00043 #include "RecoVertex/GhostTrackFitter/interface/KalmanGhostTrackUpdater.h"
00044
00045 #include "RecoVertex/GhostTrackFitter/interface/GhostTrackVertexFinder.h"
00046
00047
00048
00049 #ifdef DEBUG
00050 # include "RecoBTag/SecondaryVertex/interface/TrackKinematics.h"
00051 #endif
00052
00053 using namespace reco;
00054
00055 namespace {
00056 using namespace ROOT::Math;
00057
00058 typedef SVector<double, 3> Vector3;
00059
00060 typedef SMatrix<double, 3, 3, MatRepSym<double, 3> > Matrix3S;
00061 typedef SMatrix<double, 2, 2, MatRepSym<double, 2> > Matrix2S;
00062 typedef SMatrix<double, 2, 3> Matrix23;
00063
00064 typedef ReferenceCountingPointer<VertexTrack<5> >
00065 RefCountedVertexTrack;
00066 typedef ReferenceCountingPointer<LinearizedTrackState<5> >
00067 RefCountedLinearizedTrackState;
00068
00069 struct VtxTrackIs : public std::unary_function<
00070 RefCountedVertexTrack, bool> {
00071 VtxTrackIs(const TransientTrack &track) : track(track) {}
00072 bool operator()(const RefCountedVertexTrack &vtxTrack) const
00073 { return vtxTrack->linearizedTrack()->track() == track; }
00074
00075 const TransientTrack &track;
00076 };
00077
00078 static inline Vector3 conv(const GlobalVector &vec)
00079 {
00080 Vector3 result;
00081 result[0] = vec.x();
00082 result[1] = vec.y();
00083 result[2] = vec.z();
00084 return result;
00085 }
00086
00087 static inline double sqr(double arg) { return arg * arg; }
00088 }
00089
00090 struct GhostTrackVertexFinder::FinderInfo {
00091 FinderInfo(const CachingVertex<5> &primaryVertex,
00092 const GhostTrack &ghostTrack, const BeamSpot &beamSpot,
00093 bool hasBeamSpot, bool hasPrimaries) :
00094 primaryVertex(primaryVertex), pred(ghostTrack.prediction()),
00095 prior(ghostTrack.prior()), states(ghostTrack.states()),
00096 beamSpot(beamSpot), hasBeamSpot(hasBeamSpot),
00097 hasPrimaries(hasPrimaries), field(0) {}
00098
00099 const CachingVertex<5> &primaryVertex;
00100 GhostTrackPrediction pred;
00101 GhostTrackPrediction prior;
00102 std::vector<GhostTrackState> states;
00103 VertexState beamSpot;
00104 bool hasBeamSpot;
00105 bool hasPrimaries;
00106 const MagneticField *field;
00107 TransientTrack ghostTrack;
00108 };
00109
00110 static TransientTrack transientGhostTrack(const GhostTrackPrediction &pred,
00111 const MagneticField *field)
00112 { return TransientTrackFromFTSFactory().build(pred.fts(field)); }
00113
00114 static double vtxErrorLong(const GlobalError &error, const GlobalVector &dir)
00115 {
00116 return ROOT::Math::Similarity(conv(dir), error.matrix_new());
00117 }
00118
00119 static GlobalPoint vtxMean(const GlobalPoint &p1, const GlobalError &e1,
00120 const GlobalPoint &p2, const GlobalError &e2)
00121 {
00122 GlobalVector diff = p2 - p1;
00123
00124 double err1 = vtxErrorLong(e1, diff);
00125 double err2 = vtxErrorLong(e2, diff);
00126
00127 double weight = err1 / (err1 + err2);
00128
00129 return p1 + weight * diff;
00130 }
00131
00132 static VertexState stateMean(const VertexState &v1, const VertexState &v2)
00133 {
00134 return VertexState(vtxMean(v1.position(), v1.error(),
00135 v2.position(), v2.error()),
00136 v1.error() + v2.error());
00137 }
00138
00139 static bool covarianceUpdate(Matrix3S &cov, const Vector3 &residual,
00140 const Matrix3S &error, double &chi2,
00141 double theta, double phi)
00142 {
00143 using namespace ROOT::Math;
00144
00145 Matrix23 jacobian;
00146 jacobian(0, 0) = std::cos(phi) * std::cos(theta);
00147 jacobian(0, 1) = std::sin(phi) * std::cos(theta);
00148 jacobian(0, 2) = -std::sin(theta);
00149 jacobian(1, 0) = -std::sin(phi);
00150 jacobian(1, 1) = std::cos(phi);
00151
00152 Matrix2S measErr = Similarity(jacobian, error);
00153 Matrix2S combErr = Similarity(jacobian, cov) + measErr;
00154 if (!measErr.Invert() || !combErr.Invert())
00155 return false;
00156
00157 cov -= SimilarityT(jacobian * cov, combErr);
00158 chi2 += Similarity(jacobian * residual, measErr);
00159
00160 return true;
00161 }
00162
00163 static CachingVertex<5> vertexAtState(const TransientTrack &ghostTrack,
00164 const GhostTrackPrediction &pred,
00165 const GhostTrackState &state)
00166 {
00167 LinearizedTrackStateFactory linTrackFactory;
00168 VertexTrackFactory<5> vertexTrackFactory;
00169
00170 if (!state.isValid())
00171 return CachingVertex<5>();
00172
00173 GlobalPoint pca1 = pred.position(state.lambda());
00174 GlobalError err1 = pred.positionError(state.lambda());
00175
00176 GlobalPoint pca2 = state.globalPosition();
00177 GlobalError err2 = state.cartesianError();
00178
00179 GlobalPoint point = vtxMean(pca1, err1, pca2, err2);
00180
00181 TransientTrack recTrack = state.track();
00182
00183 RefCountedLinearizedTrackState linState[2] = {
00184 linTrackFactory.linearizedTrackState(point, ghostTrack),
00185 linTrackFactory.linearizedTrackState(point, recTrack)
00186 };
00187 if ( !linState[0]->isValid() || !linState[1]->isValid() )
00188 return CachingVertex<5>();
00189
00190 Matrix3S cov = SMatrixIdentity();
00191 cov *= 10000;
00192 double chi2 = 0.;
00193 if (!covarianceUpdate(cov, conv(pca1 - point), err1.matrix_new(), chi2,
00194 linState[0]->predictedStateParameters()[1],
00195 linState[0]->predictedStateParameters()[2]) ||
00196 !covarianceUpdate(cov, conv(pca2 - point), err2.matrix_new(), chi2,
00197 linState[1]->predictedStateParameters()[1],
00198 linState[1]->predictedStateParameters()[2]))
00199 return CachingVertex<5>();
00200
00201 GlobalError error(cov);
00202 VertexState vtxState(point, error);
00203
00204 std::vector<RefCountedVertexTrack> linTracks(2);
00205 linTracks[0] = vertexTrackFactory.vertexTrack(linState[0], vtxState);
00206 linTracks[1] = vertexTrackFactory.vertexTrack(linState[1], vtxState);
00207
00208 return CachingVertex<5>(point, error, linTracks, chi2);
00209 }
00210
00211 static RefCountedVertexTrack relinearizeTrack(
00212 const RefCountedVertexTrack &track,
00213 const VertexState &state,
00214 const VertexTrackFactory<5> factory)
00215 {
00216 RefCountedLinearizedTrackState linTrack = track->linearizedTrack();
00217 linTrack = linTrack->stateWithNewLinearizationPoint(state.position());
00218 return factory.vertexTrack(linTrack, state);
00219 }
00220
00221 static std::vector<RefCountedVertexTrack> relinearizeTracks(
00222 const std::vector<RefCountedVertexTrack> &tracks,
00223 const VertexState &state)
00224 {
00225 VertexTrackFactory<5> vertexTrackFactory;
00226
00227 std::vector<RefCountedVertexTrack> finalTracks;
00228 finalTracks.reserve(tracks.size());
00229
00230 for(std::vector<RefCountedVertexTrack>::const_iterator iter =
00231 tracks.begin(); iter != tracks.end(); ++iter)
00232 finalTracks.push_back(
00233 relinearizeTrack(*iter, state, vertexTrackFactory));
00234
00235 return finalTracks;
00236 }
00237
00238 double GhostTrackVertexFinder::vertexCompat(const CachingVertex<5> &vtx1,
00239 const CachingVertex<5> &vtx2,
00240 const FinderInfo &info,
00241 double scale1, double scale2)
00242 {
00243 Vector3 diff = conv(vtx2.position() - vtx1.position());
00244 Matrix3S cov = scale1 * vtx1.error().matrix_new() +
00245 scale2 * vtx2.error().matrix_new();
00246
00247 return sqr(ROOT::Math::Dot(diff, diff)) / ROOT::Math::Similarity(cov, diff);
00248 }
00249
00250 static double trackVertexCompat(const CachingVertex<5> &vtx,
00251 const RefCountedVertexTrack &vertexTrack)
00252 {
00253 using namespace ROOT::Math;
00254
00255 TransientTrack track = vertexTrack->linearizedTrack()->track();
00256 GlobalPoint point = vtx.position();
00257 AnalyticalImpactPointExtrapolator extrap(track.field());
00258 TrajectoryStateOnSurface tsos =
00259 extrap.extrapolate(track.impactPointState(), point);
00260
00261 if (!tsos.isValid())
00262 return 1.0e6;
00263
00264 GlobalPoint point1 = vtx.position();
00265 GlobalPoint point2 = tsos.globalPosition();
00266 Vector3 dir = conv(point2 - point1);
00267 Matrix3S error = vtx.error().matrix_new() +
00268 tsos.cartesianError().matrix().Sub<Matrix3S>(0, 0);
00269 if (!error.Invert())
00270 return 1.0e6;
00271
00272 return ROOT::Math::Similarity(error, dir);
00273 }
00274
00275 GhostTrackVertexFinder::GhostTrackVertexFinder() :
00276 maxFitChi2_(10.0),
00277 mergeThreshold_(3.0),
00278 primcut_(2.0),
00279 seccut_(5.0),
00280 fitType_(kAlwaysWithGhostTrack)
00281 {
00282 }
00283
00284 GhostTrackVertexFinder::GhostTrackVertexFinder(
00285 double maxFitChi2, double mergeThreshold, double primcut,
00286 double seccut, FitType fitType) :
00287 maxFitChi2_(maxFitChi2),
00288 mergeThreshold_(mergeThreshold),
00289 primcut_(primcut),
00290 seccut_(seccut),
00291 fitType_(fitType)
00292 {
00293 }
00294
00295 GhostTrackVertexFinder::~GhostTrackVertexFinder()
00296 {
00297 }
00298
00299 GhostTrackFitter &GhostTrackVertexFinder::ghostTrackFitter() const
00300 {
00301 if (!ghostTrackFitter_.get())
00302 ghostTrackFitter_.reset(new GhostTrackFitter);
00303
00304 return *ghostTrackFitter_;
00305 }
00306
00307 VertexFitter<5> &GhostTrackVertexFinder::vertexFitter(bool primary) const
00308 {
00309 std::auto_ptr<VertexFitter<5> > *ptr =
00310 primary ? &primVertexFitter_ : &secVertexFitter_;
00311 if (!ptr->get())
00312 ptr->reset(new AdaptiveVertexFitter(
00313 GeometricAnnealing(primary
00314 ? primcut_ : seccut_)));
00315
00316 return **ptr;
00317 }
00318
00319 #ifdef DEBUG
00320 static void debugTrack(const TransientTrack &track, const TransientVertex *vtx = 0)
00321 {
00322 const FreeTrajectoryState &fts = track.initialFreeState();
00323 GlobalVector mom = fts.momentum();
00324 std::cout << "\tpt = " << mom.perp() << ", eta = " << mom.eta() << ", phi = " << mom.phi() << ", charge = " << fts.charge();
00325 if (vtx && vtx->hasTrackWeight())
00326 std::cout << ", w = " << vtx->trackWeight(track) << std::endl;
00327 else
00328 std::cout << std::endl;
00329 }
00330
00331 static void debugTracks(const std::vector<TransientTrack> &tracks, const TransientVertex *vtx = 0)
00332 {
00333 TrackKinematics kin;
00334 for(std::vector<TransientTrack>::const_iterator iter = tracks.begin();
00335 iter != tracks.end(); ++iter) {
00336 debugTrack(*iter, vtx);
00337 try {
00338 TrackBaseRef track = iter->trackBaseRef();
00339 kin.add(*track);
00340 } catch(...) {
00341
00342 }
00343 }
00344
00345 std::cout << "mass = " << kin.vectorSum().M() << std::endl;
00346 }
00347
00348 static void debugTracks(const std::vector<RefCountedVertexTrack> &tracks)
00349 {
00350 std::vector<TransientTrack> tracks_;
00351 for(std::vector<RefCountedVertexTrack>::const_iterator iter = tracks.begin();
00352 iter != tracks.end(); ++iter) {
00353 std::cout << "+ " << (*iter)->linearizedTrack()->linearizationPoint() << std::endl;
00354 tracks_.push_back((*iter)->linearizedTrack()->track());
00355 }
00356 debugTracks(tracks_);
00357 }
00358
00359 static void debugVertex(const TransientVertex &vtx, const GhostTrackPrediction &pred)
00360 {
00361 double origin = 0.;
00362 std::cout << vtx.position() << "-> " << vtx.totalChiSquared() << " with " << vtx.degreesOfFreedom() << std::endl;
00363
00364 double err = std::sqrt(vtxErrorLong(
00365 vtx.positionError(),
00366 pred.direction())
00367 / pred.rho2());
00368 std::cout << "* " << (pred.lambda(vtx.position()) * pred.rho() - origin)
00369 << " +-" << err
00370 << " => " << ((pred.lambda(vtx.position()) * pred.rho() - origin) / err)
00371 << std::endl;
00372
00373 std::vector<TransientTrack> tracks = vtx.originalTracks();
00374 debugTracks(tracks, &vtx);
00375 }
00376 #endif
00377
00378 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00379 const Vertex &primaryVertex,
00380 const GlobalVector &direction,
00381 double coneRadius,
00382 const std::vector<TransientTrack> &tracks) const
00383 {
00384 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00385 RecoVertex::convertError(primaryVertex.error()),
00386 direction, coneRadius, tracks);
00387 }
00388
00389 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00390 const Vertex &primaryVertex,
00391 const GlobalVector &direction,
00392 double coneRadius,
00393 const reco::BeamSpot &beamSpot,
00394 const std::vector<TransientTrack> &tracks) const
00395 {
00396 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00397 RecoVertex::convertError(primaryVertex.error()),
00398 direction, coneRadius, beamSpot, tracks);
00399 }
00400
00401 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00402 const Vertex &primaryVertex,
00403 const GlobalVector &direction,
00404 double coneRadius,
00405 const reco::BeamSpot &beamSpot,
00406 const std::vector<TransientTrack> &primaries,
00407 const std::vector<TransientTrack> &tracks) const
00408 {
00409 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00410 RecoVertex::convertError(primaryVertex.error()),
00411 direction, coneRadius, beamSpot, primaries, tracks);
00412 }
00413
00414 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00415 const GlobalPoint &primaryPosition,
00416 const GlobalError &primaryError,
00417 const GlobalVector &direction,
00418 double coneRadius,
00419 const std::vector<TransientTrack> &tracks) const
00420 {
00421 GhostTrack ghostTrack =
00422 ghostTrackFitter().fit(primaryPosition, primaryError,
00423 direction, coneRadius, tracks);
00424
00425 CachingVertex<5> primary(primaryPosition, primaryError,
00426 std::vector<RefCountedVertexTrack>(), 0.);
00427
00428 std::vector<TransientVertex> result = vertices(ghostTrack, primary);
00429
00430 #ifdef DEBUG
00431 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00432 iter != result.end(); ++iter)
00433 debugVertex(*iter, ghostTrack.prediction());
00434 #endif
00435
00436 return result;
00437 }
00438
00439 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00440 const GlobalPoint &primaryPosition,
00441 const GlobalError &primaryError,
00442 const GlobalVector &direction,
00443 double coneRadius,
00444 const BeamSpot &beamSpot,
00445 const std::vector<TransientTrack> &tracks) const
00446 {
00447 GhostTrack ghostTrack =
00448 ghostTrackFitter().fit(primaryPosition, primaryError,
00449 direction, coneRadius, tracks);
00450
00451 CachingVertex<5> primary(primaryPosition, primaryError,
00452 std::vector<RefCountedVertexTrack>(), 0.);
00453
00454 std::vector<TransientVertex> result =
00455 vertices(ghostTrack, primary, beamSpot, true);
00456
00457 #ifdef DEBUG
00458 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00459 iter != result.end(); ++iter)
00460 debugVertex(*iter, ghostTrack.prediction());
00461 #endif
00462
00463 return result;
00464 }
00465
00466 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00467 const GlobalPoint &primaryPosition,
00468 const GlobalError &primaryError,
00469 const GlobalVector &direction,
00470 double coneRadius,
00471 const BeamSpot &beamSpot,
00472 const std::vector<TransientTrack> &primaries,
00473 const std::vector<TransientTrack> &tracks) const
00474 {
00475 GhostTrack ghostTrack =
00476 ghostTrackFitter().fit(primaryPosition, primaryError,
00477 direction, coneRadius, tracks);
00478
00479 std::vector<RefCountedVertexTrack> primaryVertexTracks;
00480 if (!primaries.empty()) {
00481 LinearizedTrackStateFactory linTrackFactory;
00482 VertexTrackFactory<5> vertexTrackFactory;
00483
00484 VertexState state(primaryPosition, primaryError);
00485
00486 for(std::vector<TransientTrack>::const_iterator iter =
00487 primaries.begin(); iter != primaries.end(); ++iter) {
00488
00489 RefCountedLinearizedTrackState linState =
00490 linTrackFactory.linearizedTrackState(
00491 primaryPosition, *iter);
00492
00493 primaryVertexTracks.push_back(
00494 vertexTrackFactory.vertexTrack(
00495 linState, state));
00496 }
00497 }
00498
00499 CachingVertex<5> primary(primaryPosition, primaryError,
00500 primaryVertexTracks, 0.);
00501
00502 std::vector<TransientVertex> result =
00503 vertices(ghostTrack, primary, beamSpot, true, true);
00504
00505 #ifdef DEBUG
00506 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00507 iter != result.end(); ++iter)
00508 debugVertex(*iter, ghostTrack.prediction());
00509 #endif
00510
00511 return result;
00512 }
00513
00514 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00515 const GlobalPoint &primaryPosition,
00516 const GlobalError &primaryError,
00517 const GhostTrack &ghostTrack) const
00518 {
00519 CachingVertex<5> primary(primaryPosition, primaryError,
00520 std::vector<RefCountedVertexTrack>(), 0.);
00521
00522 std::vector<TransientVertex> result = vertices(ghostTrack, primary);
00523
00524 #ifdef DEBUG
00525 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00526 iter != result.end(); ++iter)
00527 debugVertex(*iter, ghostTrack.prediction());
00528 #endif
00529
00530 return result;
00531 }
00532
00533
00534 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00535 const GlobalPoint &primaryPosition,
00536 const GlobalError &primaryError,
00537 const BeamSpot &beamSpot,
00538 const GhostTrack &ghostTrack) const
00539 {
00540 CachingVertex<5> primary(primaryPosition, primaryError,
00541 std::vector<RefCountedVertexTrack>(), 0.);
00542
00543 std::vector<TransientVertex> result =
00544 vertices(ghostTrack, primary, beamSpot, true);
00545
00546 #ifdef DEBUG
00547 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00548 iter != result.end(); ++iter)
00549 debugVertex(*iter, ghostTrack.prediction());
00550 #endif
00551
00552 return result;
00553 }
00554
00555 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00556 const GlobalPoint &primaryPosition,
00557 const GlobalError &primaryError,
00558 const BeamSpot &beamSpot,
00559 const std::vector<TransientTrack> &primaries,
00560 const GhostTrack &ghostTrack) const
00561 {
00562 std::vector<RefCountedVertexTrack> primaryVertexTracks;
00563 if (!primaries.empty()) {
00564 LinearizedTrackStateFactory linTrackFactory;
00565 VertexTrackFactory<5> vertexTrackFactory;
00566
00567 VertexState state(primaryPosition, primaryError);
00568
00569 for(std::vector<TransientTrack>::const_iterator iter =
00570 primaries.begin(); iter != primaries.end(); ++iter) {
00571
00572 RefCountedLinearizedTrackState linState =
00573 linTrackFactory.linearizedTrackState(
00574 primaryPosition, *iter);
00575
00576 primaryVertexTracks.push_back(
00577 vertexTrackFactory.vertexTrack(
00578 linState, state));
00579 }
00580 }
00581
00582 CachingVertex<5> primary(primaryPosition, primaryError,
00583 primaryVertexTracks, 0.);
00584
00585 std::vector<TransientVertex> result =
00586 vertices(ghostTrack, primary, beamSpot, true, true);
00587
00588 #ifdef DEBUG
00589 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00590 iter != result.end(); ++iter)
00591 debugVertex(*iter, ghostTrack.prediction());
00592 #endif
00593
00594 return result;
00595 }
00596
00597 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00598 const Vertex &primaryVertex,
00599 const GhostTrack &ghostTrack) const
00600 {
00601 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00602 RecoVertex::convertError(primaryVertex.error()),
00603 ghostTrack);
00604 }
00605
00606 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00607 const Vertex &primaryVertex,
00608 const BeamSpot &beamSpot,
00609 const GhostTrack &ghostTrack) const
00610 {
00611 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00612 RecoVertex::convertError(primaryVertex.error()),
00613 beamSpot, ghostTrack);
00614 }
00615
00616 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00617 const Vertex &primaryVertex,
00618 const BeamSpot &beamSpot,
00619 const std::vector<TransientTrack> &primaries,
00620 const GhostTrack &ghostTrack) const
00621 {
00622 return vertices(RecoVertex::convertPos(primaryVertex.position()),
00623 RecoVertex::convertError(primaryVertex.error()),
00624 beamSpot, primaries, ghostTrack);
00625 }
00626
00627 static GhostTrackPrediction dummyPrediction(
00628 const Vertex &primaryVertex, const Track &ghostTrack)
00629 {
00630 return GhostTrackPrediction(
00631 RecoVertex::convertPos(primaryVertex.position()),
00632 RecoVertex::convertError(primaryVertex.error()),
00633 GlobalVector(ghostTrack.px(), ghostTrack.py(),
00634 ghostTrack.pz()), 0.05);
00635 }
00636
00637 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00638 const Vertex &primaryVertex,
00639 const Track &ghostTrack,
00640 const std::vector<TransientTrack> &tracks,
00641 const std::vector<float> &weights) const
00642 {
00643 GhostTrack ghostTrack_(ghostTrack, tracks, weights,
00644 dummyPrediction(primaryVertex, ghostTrack),
00645 RecoVertex::convertPos(
00646 primaryVertex.position()),
00647 true);
00648
00649 CachingVertex<5> primary(
00650 RecoVertex::convertPos(primaryVertex.position()),
00651 RecoVertex::convertError(primaryVertex.error()),
00652 std::vector<RefCountedVertexTrack>(), 0.);
00653
00654 std::vector<TransientVertex> result = vertices(ghostTrack_, primary);
00655
00656 #ifdef DEBUG
00657 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00658 iter != result.end(); ++iter)
00659 debugVertex(*iter, ghostTrack_.prediction());
00660 #endif
00661
00662 return result;
00663 }
00664
00665 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00666 const Vertex &primaryVertex,
00667 const Track &ghostTrack,
00668 const BeamSpot &beamSpot,
00669 const std::vector<TransientTrack> &tracks,
00670 const std::vector<float> &weights) const
00671 {
00672 GhostTrack ghostTrack_(ghostTrack, tracks, weights,
00673 dummyPrediction(primaryVertex, ghostTrack),
00674 RecoVertex::convertPos(primaryVertex.position()),
00675 true);
00676
00677 CachingVertex<5> primary(
00678 RecoVertex::convertPos(primaryVertex.position()),
00679 RecoVertex::convertError(primaryVertex.error()),
00680 std::vector<RefCountedVertexTrack>(), 0.);
00681
00682 std::vector<TransientVertex> result =
00683 vertices(ghostTrack_, primary, beamSpot, true);
00684
00685 #ifdef DEBUG
00686 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00687 iter != result.end(); ++iter)
00688 debugVertex(*iter, ghostTrack_.prediction());
00689 #endif
00690
00691 return result;
00692 }
00693
00694 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
00695 const Vertex &primaryVertex,
00696 const Track &ghostTrack,
00697 const BeamSpot &beamSpot,
00698 const std::vector<TransientTrack> &primaries,
00699 const std::vector<TransientTrack> &tracks,
00700 const std::vector<float> &weights) const
00701 {
00702 GhostTrack ghostTrack_(ghostTrack, tracks, weights,
00703 dummyPrediction(primaryVertex, ghostTrack),
00704 RecoVertex::convertPos(primaryVertex.position()),
00705 true);
00706
00707 std::vector<RefCountedVertexTrack> primaryVertexTracks;
00708 if (!primaries.empty()) {
00709 LinearizedTrackStateFactory linTrackFactory;
00710 VertexTrackFactory<5> vertexTrackFactory;
00711
00712 VertexState state(
00713 RecoVertex::convertPos(primaryVertex.position()),
00714 RecoVertex::convertError(primaryVertex.error()));
00715
00716 for(std::vector<TransientTrack>::const_iterator iter =
00717 primaries.begin(); iter != primaries.end(); ++iter) {
00718
00719 RefCountedLinearizedTrackState linState =
00720 linTrackFactory.linearizedTrackState(
00721 state.position(), *iter);
00722
00723 primaryVertexTracks.push_back(
00724 vertexTrackFactory.vertexTrack(
00725 linState, state));
00726 }
00727 }
00728
00729 CachingVertex<5> primary(
00730 RecoVertex::convertPos(primaryVertex.position()),
00731 RecoVertex::convertError(primaryVertex.error()),
00732 primaryVertexTracks, 0.);
00733
00734 std::vector<TransientVertex> result =
00735 vertices(ghostTrack_, primary, beamSpot, true, true);
00736
00737 #ifdef DEBUG
00738 for(std::vector<TransientVertex>::const_iterator iter = result.begin();
00739 iter != result.end(); ++iter)
00740 debugVertex(*iter, ghostTrack_.prediction());
00741 #endif
00742
00743 return result;
00744 }
00745
00746 static double fitChi2(const CachingVertex<5> &vtx)
00747 {
00748 return vtx.totalChiSquared() / vtx.degreesOfFreedom();
00749 }
00750
00751 std::vector<CachingVertex<5> > GhostTrackVertexFinder::initialVertices(
00752 const FinderInfo &info) const
00753 {
00754 std::vector<CachingVertex<5> > vertices;
00755 for(std::vector<GhostTrackState>::const_iterator iter =
00756 info.states.begin(); iter != info. states.end(); ++iter) {
00757
00758 if (!iter->isValid())
00759 continue;
00760
00761 GhostTrackState state(*iter);
00762 state.linearize(info.pred);
00763
00764 if (!state.isValid())
00765 continue;
00766
00767 CachingVertex<5> vtx = vertexAtState(info.ghostTrack,
00768 info.pred, state);
00769
00770 if (vtx.isValid())
00771 vertices.push_back(vtx);
00772 }
00773
00774 return vertices;
00775 }
00776
00777 static void mergeTrackHelper(const std::vector<RefCountedVertexTrack> &tracks,
00778 std::vector<RefCountedVertexTrack> &newTracks,
00779 const VertexState &state,
00780 const VtxTrackIs &ghostTrackFinder,
00781 RefCountedVertexTrack &ghostTrack,
00782 const VertexTrackFactory<5> &factory)
00783 {
00784 for(std::vector<RefCountedVertexTrack>::const_iterator iter =
00785 tracks.begin(); iter != tracks.end(); ++iter) {
00786 bool gt = ghostTrackFinder(*iter);
00787 if (gt && ghostTrack)
00788 continue;
00789
00790 RefCountedVertexTrack track =
00791 relinearizeTrack(*iter, state, factory);
00792
00793 if (gt)
00794 ghostTrack = *iter;
00795 else
00796 newTracks.push_back(*iter);
00797 }
00798 }
00799
00800 CachingVertex<5> GhostTrackVertexFinder::mergeVertices(
00801 const CachingVertex<5> &vertex1,
00802 const CachingVertex<5> &vertex2,
00803 const FinderInfo &info,
00804 bool isPrimary) const
00805 {
00806 VertexTrackFactory<5> vertexTrackFactory;
00807
00808 VertexState state = stateMean(vertex1.vertexState(),
00809 vertex2.vertexState());
00810 std::vector<RefCountedVertexTrack> linTracks;
00811 VtxTrackIs isGhostTrack(info.ghostTrack);
00812 RefCountedVertexTrack vtxGhostTrack;
00813
00814 mergeTrackHelper(vertex1.tracks(), linTracks, state,
00815 isGhostTrack, vtxGhostTrack, vertexTrackFactory);
00816 mergeTrackHelper(vertex2.tracks(), linTracks, state,
00817 isGhostTrack, vtxGhostTrack, vertexTrackFactory);
00818
00819 if (vtxGhostTrack &&
00820 (fitType_ == kAlwaysWithGhostTrack || linTracks.size() < 2))
00821 linTracks.push_back(
00822 vertexTrackFactory.vertexTrack(
00823 vtxGhostTrack->linearizedTrack(),
00824 vtxGhostTrack->vertexState()));
00825
00826 try {
00827 CachingVertex<5> vtx;
00828 if (info.hasBeamSpot && isPrimary)
00829 vtx = vertexFitter(true).vertex(linTracks,
00830 info.beamSpot.position(),
00831 info.beamSpot.error());
00832 else
00833 vtx = vertexFitter(isPrimary).vertex(linTracks);
00834 if (vtx.isValid())
00835 return vtx;
00836 } catch(const VertexException &e) {
00837
00838 }
00839
00840 return CachingVertex<5>();
00841 }
00842
00843 bool GhostTrackVertexFinder::recursiveMerge(
00844 std::vector<CachingVertex<5> > &vertices,
00845 const FinderInfo &info) const
00846 {
00847 typedef std::pair<unsigned int, unsigned int> Indices;
00848
00849 std::multimap<float, Indices> compatMap;
00850 unsigned int n = vertices.size();
00851 for(unsigned int i = 0; i < n; i++) {
00852 const CachingVertex<5> &v1 = vertices[i];
00853 for(unsigned int j = i + 1; j < n; j++) {
00854 const CachingVertex<5> &v2 = vertices[j];
00855
00856 float compat = vertexCompat(v1, v2, info,
00857 i == 0 ? primcut_ : seccut_, seccut_);
00858
00859 if (compat > mergeThreshold_)
00860 continue;
00861
00862 compatMap.insert(
00863 std::make_pair(compat, Indices(i, j)));
00864 }
00865 }
00866
00867 bool changesMade = false;
00868 bool repeat = true;
00869 while(repeat) {
00870 repeat = false;
00871 for(std::multimap<float, Indices>::const_iterator iter =
00872 compatMap.begin(); iter != compatMap.end(); ++iter) {
00873 unsigned int v1 = iter->second.first;
00874 unsigned int v2 = iter->second.second;
00875
00876 CachingVertex<5> newVtx =
00877 mergeVertices(vertices[v1], vertices[v2],
00878 info, v1 == 0);
00879 if (!newVtx.isValid() ||
00880 (v1 != 0 && fitChi2(newVtx) > maxFitChi2_))
00881 continue;
00882
00883 std::swap(vertices[v1], newVtx);
00884 vertices.erase(vertices.begin() + v2);
00885 n--;
00886
00887 std::multimap<float, Indices> newCompatMap;
00888 for(++iter; iter != compatMap.end(); ++iter) {
00889 if (iter->second.first == v1 ||
00890 iter->second.first == v2 ||
00891 iter->second.second == v1 ||
00892 iter->second.second == v2)
00893 continue;
00894
00895 Indices indices = iter->second;
00896 indices.first -= indices.first > v2;
00897 indices.second -= indices.second > v2;
00898
00899 newCompatMap.insert(std::make_pair(
00900 iter->first, indices));
00901 }
00902 std::swap(compatMap, newCompatMap);
00903
00904 for(unsigned int i = 0; i < n; i++) {
00905 if (i == v1)
00906 continue;
00907
00908 const CachingVertex<5> &other = vertices[i];
00909 float compat = vertexCompat(
00910 vertices[v1], other, info,
00911 v1 == 0 ? primcut_ : seccut_,
00912 i == 0 ? primcut_ : seccut_);
00913
00914 if (compat > mergeThreshold_)
00915 continue;
00916
00917 compatMap.insert(
00918 std::make_pair(
00919 compat,
00920 Indices(std::min(i, v1),
00921 std::max(i, v1))));
00922 }
00923
00924 changesMade = true;
00925 repeat = true;
00926 break;
00927 }
00928 }
00929
00930 return changesMade;
00931 }
00932
00933 bool GhostTrackVertexFinder::reassignTracks(
00934 std::vector<CachingVertex<5> > &vertices_,
00935 const FinderInfo &info) const
00936 {
00937 std::vector<std::pair<RefCountedVertexTrack,
00938 std::vector<RefCountedVertexTrack> > >
00939 trackBundles(vertices_.size());
00940
00941 VtxTrackIs isGhostTrack(info.ghostTrack);
00942
00943 bool assignmentChanged = false;
00944 for(std::vector<CachingVertex<5> >::const_iterator iter =
00945 vertices_.begin(); iter != vertices_.end(); ++iter) {
00946 std::vector<RefCountedVertexTrack> vtxTracks = iter->tracks();
00947
00948 if (vtxTracks.empty()) {
00949 LinearizedTrackStateFactory linTrackFactory;
00950 VertexTrackFactory<5> vertexTrackFactory;
00951
00952 RefCountedLinearizedTrackState linState =
00953 linTrackFactory.linearizedTrackState(
00954 iter->position(), info.ghostTrack);
00955
00956 trackBundles[iter - vertices_.begin()].first =
00957 vertexTrackFactory.vertexTrack(
00958 linState, iter->vertexState());
00959 }
00960
00961 for(std::vector<RefCountedVertexTrack>::const_iterator track =
00962 vtxTracks.begin(); track != vtxTracks.end(); ++track) {
00963
00964 if (isGhostTrack(*track)) {
00965 trackBundles[iter - vertices_.begin()]
00966 .first = *track;
00967 continue;
00968 }
00969
00970 if ((*track)->weight() < 1e-3) {
00971 trackBundles[iter - vertices_.begin()]
00972 .second.push_back(*track);
00973 continue;
00974 }
00975
00976 unsigned int idx = iter - vertices_.begin();
00977 double best = 1.0e9;
00978 for(std::vector<CachingVertex<5> >::const_iterator
00979 vtx = vertices_.begin();
00980 vtx != vertices_.end(); ++vtx) {
00981 if (info.pred.lambda(vtx->position()) <
00982 info.pred.lambda(vertices_[0].position()))
00983 continue;
00984
00985 double compat =
00986 sqr(trackVertexCompat(*vtx, *track));
00987
00988 compat /= (vtx == vertices_.begin()) ?
00989 primcut_ : seccut_;
00990
00991 if (compat < best) {
00992 best = compat;
00993 idx = vtx - vertices_.begin();
00994 }
00995 }
00996
00997 if ((int)idx != iter - vertices_.begin())
00998 assignmentChanged = true;
00999
01000 trackBundles[idx].second.push_back(*track);
01001 }
01002 }
01003
01004 if (!assignmentChanged)
01005 return false;
01006
01007 VertexTrackFactory<5> vertexTrackFactory;
01008 std::vector<CachingVertex<5> > vertices;
01009 vertices.reserve(vertices_.size());
01010
01011 for(std::vector<CachingVertex<5> >::const_iterator iter =
01012 vertices_.begin(); iter != vertices_.end(); ++iter) {
01013 const std::vector<RefCountedVertexTrack> &tracks =
01014 trackBundles[iter - vertices_.begin()].second;
01015 if (tracks.empty())
01016 continue;
01017
01018 CachingVertex<5> vtx;
01019
01020 if (tracks.size() == 1) {
01021 const TransientTrack &track =
01022 tracks[0]->linearizedTrack()->track();
01023
01024 int idx = -1;
01025 for(std::vector<GhostTrackState>::const_iterator iter =
01026 info.states.begin();
01027 iter != info.states.end(); ++iter) {
01028 if (iter->isTrack() && iter->track() == track) {
01029 idx = iter - info.states.begin();
01030 break;
01031 }
01032 }
01033 if (idx < 0)
01034 continue;
01035
01036 vtx = vertexAtState(info.ghostTrack, info.pred,
01037 info.states[idx]);
01038 if (!vtx.isValid())
01039 continue;
01040 } else {
01041 std::vector<RefCountedVertexTrack> linTracks =
01042 relinearizeTracks(tracks, iter->vertexState());
01043
01044 if (fitType_ == kAlwaysWithGhostTrack)
01045 linTracks.push_back(relinearizeTrack(
01046 trackBundles[iter - vertices_.begin()].first,
01047 iter->vertexState(), vertexTrackFactory));
01048
01049 bool primary = iter == vertices_.begin();
01050 try {
01051 if (primary && info.hasBeamSpot)
01052 vtx = vertexFitter(true).vertex(
01053 linTracks,
01054 info.beamSpot.position(),
01055 info.beamSpot.error());
01056 else
01057 vtx = vertexFitter(primary).vertex(
01058 linTracks);
01059 } catch(const VertexException &e) {
01060
01061 }
01062 if (!vtx.isValid())
01063 return false;
01064 }
01065
01066 vertices.push_back(vtx);
01067 };
01068
01069 std::swap(vertices_, vertices);
01070 return true;
01071 }
01072
01073 void GhostTrackVertexFinder::refitGhostTrack(
01074 std::vector<CachingVertex<5> > &vertices,
01075 FinderInfo &info) const
01076 {
01077 VtxTrackIs isGhostTrack(info.ghostTrack);
01078
01079 std::vector<GhostTrackState> states;
01080 std::vector<unsigned int> oldStates;
01081 oldStates.reserve(info.states.size());
01082
01083 for(std::vector<CachingVertex<5> >::const_iterator iter =
01084 vertices.begin(); iter != vertices.end(); ++iter) {
01085 std::vector<RefCountedVertexTrack> vtxTracks = iter->tracks();
01086
01087 oldStates.clear();
01088 for(std::vector<RefCountedVertexTrack>::const_iterator track =
01089 vtxTracks.begin(); track != vtxTracks.end(); ++track) {
01090
01091 if (isGhostTrack(*track) || (*track)->weight() < 1e-3)
01092 continue;
01093
01094 const TransientTrack &tt =
01095 (*track)->linearizedTrack()->track();
01096
01097 int idx = -1;
01098 for(std::vector<GhostTrackState>::const_iterator iter =
01099 info.states.begin();
01100 iter != info.states.end(); ++iter) {
01101 if (iter->isTrack() && iter->track() == tt) {
01102 idx = iter - info.states.begin();
01103 break;
01104 }
01105 }
01106
01107 if (idx >= 0)
01108 oldStates.push_back(idx);
01109 }
01110
01111 if (oldStates.size() == 1)
01112 states.push_back(info.states[oldStates[0]]);
01113 else
01114 states.push_back(GhostTrackState(iter->vertexState()));
01115 }
01116
01117 KalmanGhostTrackUpdater updater;
01118 SequentialGhostTrackFitter fitter;
01119 double ndof, chi2;
01120 info.pred = fitter.fit(updater, info.prior, states, ndof, chi2);
01121 TransientTrack ghostTrack = transientGhostTrack(info.pred, info.field);
01122
01123 std::swap(info.states, states);
01124 states.clear();
01125
01126 std::vector<CachingVertex<5> > newVertices;
01127 for(std::vector<CachingVertex<5> >::const_iterator iter =
01128 vertices.begin(); iter != vertices.end(); ++iter) {
01129 std::vector<RefCountedVertexTrack> vtxTracks = iter->tracks();
01130
01131 int idx = -1;
01132 bool redo = false;
01133 for(std::vector<RefCountedVertexTrack>::iterator track =
01134 vtxTracks.begin(); track != vtxTracks.end(); ++track) {
01135
01136 if (isGhostTrack(*track)) {
01137 LinearizedTrackStateFactory linTrackFactory;
01138 VertexTrackFactory<5> vertexTrackFactory;
01139
01140 *track = vertexTrackFactory.vertexTrack(
01141 linTrackFactory.linearizedTrackState(
01142 iter->position(), ghostTrack),
01143 iter->vertexState());
01144 redo = true;
01145 continue;
01146 }
01147
01148 const TransientTrack &tt =
01149 (*track)->linearizedTrack()->track();
01150
01151 if (idx >= 0) {
01152 idx = -1;
01153 break;
01154 }
01155
01156 for(std::vector<GhostTrackState>::const_iterator it =
01157 info.states.begin();
01158 it != info.states.end(); ++it) {
01159 if (!it->isTrack())
01160 continue;
01161
01162 if (it->track() == tt) {
01163 idx = it - info.states.begin();
01164 break;
01165 }
01166 }
01167 }
01168
01169 if (idx >= 0) {
01170 CachingVertex<5> vtx =
01171 vertexAtState(ghostTrack, info.pred,
01172 info.states[idx]);
01173 if (vtx.isValid())
01174 newVertices.push_back(vtx);
01175 } else if (redo) {
01176 bool primary = iter == vertices.begin();
01177 CachingVertex<5> vtx;
01178 if (primary && info.hasBeamSpot)
01179 vtx = vertexFitter(true).vertex(
01180 vtxTracks,
01181 info.beamSpot.position(),
01182 info.beamSpot.error());
01183 else
01184 vtx = vertexFitter(primary).vertex(vtxTracks);
01185 if (vtx.isValid())
01186 newVertices.push_back(vtx);
01187 } else
01188 newVertices.push_back(*iter);
01189 }
01190
01191 std::swap(newVertices, vertices);
01192 info.ghostTrack = ghostTrack;
01193 }
01194
01195
01196
01197 std::vector<TransientVertex> GhostTrackVertexFinder::vertices(
01198 const GhostTrack &ghostTrack,
01199 const CachingVertex<5> &primary,
01200 const BeamSpot &beamSpot,
01201 bool hasBeamSpot, bool hasPrimaries) const
01202 {
01203 FinderInfo info(primary, ghostTrack, beamSpot,
01204 hasBeamSpot, hasPrimaries);
01205
01206 std::vector<TransientVertex> result;
01207 if (info.states.empty())
01208 return result;
01209
01210 info.field = info.states[0].track().field();
01211 info.ghostTrack = transientGhostTrack(info.pred, info.field);
01212
01213 std::vector<CachingVertex<5> > vertices = initialVertices(info);
01214 if (primary.isValid()) {
01215 vertices.push_back(primary);
01216 if (vertices.size() > 1)
01217 std::swap(vertices.front(), vertices.back());
01218 }
01219
01220 unsigned int reassigned = 0;
01221 while(reassigned < 3) {
01222 if (vertices.size() < 2)
01223 break;
01224
01225 #ifdef DEBUG
01226 for(std::vector<CachingVertex<5> >::const_iterator iter =
01227 vertices.begin(); iter != vertices.end(); ++iter)
01228
01229 debugVertex(*iter, ghostTrack.prediction());
01230
01231 std::cout << "----- recursive merging: ---------" << std::endl;
01232 #endif
01233
01234 bool changed = recursiveMerge(vertices, info);
01235 if ((!changed && !reassigned) || vertices.size() < 2)
01236 break;
01237 if (changed)
01238 reassigned = 0;
01239
01240 if (fitType_ == kRefitGhostTrackWithVertices) {
01241 refitGhostTrack(vertices, info);
01242 changed = true;
01243 }
01244
01245 try {
01246 #ifdef DEBUG
01247 for(std::vector<CachingVertex<5> >::const_iterator
01248 iter = vertices.begin();
01249 iter != vertices.end(); ++iter)
01250 debugVertex(*iter, ghostTrack.prediction());
01251 std::cout << "----- reassignment: ---------" << std::endl;
01252 #endif
01253 if (reassignTracks(vertices, info)) {
01254 reassigned++;
01255 changed = true;
01256 } else
01257 reassigned = 0;
01258 } catch(const VertexException &e) {
01259
01260 }
01261
01262 if (!changed)
01263 break;
01264 }
01265
01266 for(std::vector<CachingVertex<5> >::const_iterator iter =
01267 vertices.begin(); iter != vertices.end(); ++iter) {
01268 std::vector<RefCountedVertexTrack> tracks = iter->tracks();
01269 std::vector<RefCountedVertexTrack> newTracks;
01270 newTracks.reserve(tracks.size());
01271
01272 std::remove_copy_if(tracks.begin(), tracks.end(),
01273 std::back_inserter(newTracks),
01274 VtxTrackIs(info.ghostTrack));
01275
01276 if (newTracks.empty())
01277 continue;
01278
01279 CachingVertex<5> vtx(iter->vertexState(), newTracks,
01280 iter->totalChiSquared());
01281 result.push_back(vtx);
01282 }
01283
01284 return result;
01285 }
01286
1.7.3