Skip to content

Commit

Permalink
Merge pull request #19935 from lgray/topic_vectorize_zt
Browse files Browse the repository at this point in the history
Vectorized 2D simulated annealing vertexing
  • Loading branch information
cmsbuild authored Sep 12, 2017
2 parents dfad66c + 163dda9 commit 8dd9a2c
Show file tree
Hide file tree
Showing 10 changed files with 1,336 additions and 829 deletions.
6 changes: 3 additions & 3 deletions RecoVertex/Configuration/python/RecoVertex_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -33,14 +33,14 @@
)

#timing
from RecoVertex.PrimaryVertexProducer.TkClusParameters_cff import DA2DParameters
from RecoVertex.PrimaryVertexProducer.TkClusParameters_cff import DA2D_vectParameters
unsortedOfflinePrimaryVertices1D = unsortedOfflinePrimaryVertices.clone()
unsortedOfflinePrimaryVertices1D.TkFilterParameters.minPt = cms.double(0.7)
offlinePrimaryVertices1D=sortedPrimaryVertices.clone(vertices="unsortedOfflinePrimaryVertices1D", particles="trackRefsForJetsBeforeSorting")
offlinePrimaryVertices1DWithBS=sortedPrimaryVertices.clone(vertices="unsortedOfflinePrimaryVertices1D:WithBS", particles="trackRefsForJetsBeforeSorting")
DA2DParameters.TkDAClusParameters.verbose = cms.untracked.bool(False)
DA2D_vectParameters.TkDAClusParameters.verbose = cms.untracked.bool(False)
unsortedOfflinePrimaryVertices4D = unsortedOfflinePrimaryVertices.clone( verbose = cms.untracked.bool(False),
TkClusParameters = DA2DParameters )
TkClusParameters = DA2D_vectParameters )
unsortedOfflinePrimaryVertices4D.TkFilterParameters.minPt = cms.double(0.7)
unsortedOfflinePrimaryVertices4D.TrackTimesLabel = cms.InputTag("trackTimeValueMapProducer:generalTracksConfigurableFlatResolutionModel")
unsortedOfflinePrimaryVertices4D.TrackTimeResosLabel = cms.InputTag("trackTimeValueMapProducer:generalTracksConfigurableFlatResolutionModelResolution")
Expand Down
102 changes: 0 additions & 102 deletions RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZT.h

This file was deleted.

240 changes: 240 additions & 0 deletions RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZT_vect.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,240 @@
#ifndef DAClusterizerInZT_vect_h
#define DAClusterizerInZT_vect_h

/**\class DAClusterizerInZT_vect
Description: separates event tracks into clusters along the beam line
Version which auto-vectorizes with gcc 4.6 or newer
*/

#include "RecoVertex/PrimaryVertexProducer/interface/TrackClusterizerInZ.h"
#include "TrackingTools/TransientTrack/interface/TransientTrack.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include <vector>
#include "DataFormats/Math/interface/Error.h"
#include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
#include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"

#include <memory>

class DAClusterizerInZT_vect final : public TrackClusterizerInZ {

public:

// Internal data structure to
struct track_t {

void addItem( double new_z, double new_t, double new_dz2, double new_dt2, const reco::TransientTrack* new_tt, double new_pi )
{
z.push_back( new_z );
t.push_back( new_t );
dz2.push_back( new_dz2 );
dt2.push_back( new_dt2 );
errsum.push_back( 1./(1./new_dz2 + 1./new_dt2) );
tt.push_back( new_tt );

pi.push_back( new_pi ); // track weight
Z_sum.push_back( 1.0 ); // Z[i] for DA clustering, initial value as done in ::fill
}

unsigned int getSize() const
{
return z.size();
}

// has to be called everytime the items are modified
void extractRaw()
{
z_ = &z.front();
t_ = &t.front();
dz2_ = &dz2.front();
dt2_ = &dt2.front();
errsum_ = &errsum.front();
Z_sum_ = &Z_sum.front();
pi_ = &pi.front();
}

double * z_; // z-coordinate at point of closest approach to the beamline
double * t_; // t-coordinate at point of closest approach to the beamline
double * pi_; // track weight

double * dz2_; // square of the error of z(pca)
double * dt2_; // square of the error of t(pca)
double * errsum_; // sum of squares of the pca errors
double * Z_sum_; // Z[i] for DA clustering

std::vector<double> z; // z-coordinate at point of closest approach to the beamline
std::vector<double> t; // t-coordinate at point of closest approach to the beamline
std::vector<double> dz2; // square of the error of z(pca)
std::vector<double> dt2; // square of the error of t(pca)
std::vector<double> errsum; // sum of squares of the pca errors
std::vector<double> Z_sum; // Z[i] for DA clustering
std::vector<double> pi; // track weight
std::vector< const reco::TransientTrack* > tt; // a pointer to the Transient Track
};

struct vertex_t {

void addItem( double new_z, double new_t, double new_pk )
{
z.push_back( new_z);
t.push_back( new_t);
pk.push_back( new_pk);

ei_cache.push_back( 0.0 );
ei.push_back( 0.0 );
sw.push_back( 0.0 );
swz.push_back( 0.0);
swt.push_back( 0.0);
se.push_back( 0.0);
swE.push_back( 0.0);

extractRaw();
}

unsigned int getSize() const
{
return z.size();
}

// has to be called everytime the items are modified
void extractRaw()
{
z_ = &z.front();
t_ = &t.front();
pk_ = &pk.front();

ei_ = &ei.front();
sw_ = &sw.front();
swz_ = &swz.front();
swt_ = &swt.front();
se_ = &se.front();
swE_ = &swE.front();
ei_cache_ = &ei_cache.front();

}

void insertItem( unsigned int i, double new_z, double new_t, double new_pk )
{
z.insert(z.begin() + i, new_z);
t.insert(t.begin() + i, new_t);
pk.insert(pk.begin() + i, new_pk);

ei_cache.insert(ei_cache.begin() + i, 0.0 );
ei.insert( ei.begin() + i, 0.0 );
sw.insert( sw.begin() + i, 0.0 );
swz.insert(swz.begin() + i, 0.0 );
swt.insert(swt.begin() + i, 0.0 );
se.insert( se.begin() + i, 0.0 );
swE.insert(swE.begin() + i, 0.0 );

extractRaw();
}

void removeItem( unsigned int i )
{
z.erase( z.begin() + i );
t.erase( t.begin() + i );
pk.erase( pk.begin() + i );

ei_cache.erase( ei_cache.begin() + i);
ei.erase( ei.begin() + i);
sw.erase( sw.begin() + i);
swz.erase( swz.begin() + i);
swt.erase( swt.begin() + i);
se.erase(se.begin() + i);
swE.erase(swE.begin() + i);

extractRaw();
}

void debugOut()
{
std::cout << "vertex_t size: " << getSize() << std::endl;

for ( unsigned int i =0; i < getSize(); ++ i)
{
std::cout << " z = " << z_[i] << " t = " << t_[i] << " pk = " << pk_[i] << std::endl;
}
}

std::vector<double> z; // z coordinate
std::vector<double> t; // t coordinate
std::vector<double> pk; // vertex weight for "constrained" clustering

double * z_;
double * t_;
double * pk_;

double * ei_cache_;
double * ei_;
double * sw_;
double * swz_;
double * swt_;
double * se_;
double * swE_;

// --- temporary numbers, used during update
std::vector<double> ei_cache;
std::vector<double> ei;
std::vector<double> sw;
std::vector<double> swz;
std::vector<double> swt;
std::vector<double> se;
std::vector<double> swE;
};

DAClusterizerInZT_vect(const edm::ParameterSet& conf);

std::vector<std::vector<reco::TransientTrack> >
clusterize(const std::vector<reco::TransientTrack> & tracks) const override;

std::vector<TransientVertex>
vertices(const std::vector<reco::TransientTrack> & tracks,
const int verbosity = 0) const ;

track_t fill(const std::vector<reco::TransientTrack> & tracks) const;

double update(double beta, track_t & gtracks,
vertex_t & gvertices, bool useRho0, const double & rho0) const;

void dump(const double beta, const vertex_t & y,
const track_t & tks, const int verbosity = 0) const;
bool merge(vertex_t & y, double & beta)const;
bool purge(vertex_t &, track_t &, double &,
const double) const;
void splitAll( vertex_t & y) const;
bool split(const double beta, track_t &t, vertex_t & y, double threshold = 1. ) const;

double beta0(const double betamax, track_t const & tks, vertex_t const & y) const;


private:
bool verbose_;
double zdumpcenter_;
double zdumpwidth_;

double vertexSize_;
double vertexSizeTime_;
int maxIterations_;
double coolingFactor_;
double betamax_;
double betastop_;
double dzCutOff_;
double d0CutOff_;
double dtCutOff_;
bool useTc_;

double mintrkweight_;
double uniquetrkweight_;
double zmerge_;
double tmerge_;
double betapurge_;

};


//#ifndef DAClusterizerInZT_new_h
#endif
Original file line number Diff line number Diff line change
Expand Up @@ -34,12 +34,13 @@
#include "RecoVertex/PrimaryVertexProducer/interface/TrackFilterForPVFindingBase.h"
#include "RecoVertex/PrimaryVertexProducer/interface/TrackClusterizerInZ.h"
#include "RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZ_vect.h"
#include "RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZT_vect.h"


#include "RecoVertex/PrimaryVertexProducer/interface/TrackFilterForPVFinding.h"
#include "RecoVertex/PrimaryVertexProducer/interface/HITrackFilterForPVFinding.h"
#include "RecoVertex/PrimaryVertexProducer/interface/GapClusterizerInZ.h"
#include "RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZ.h"
#include "RecoVertex/PrimaryVertexProducer/interface/DAClusterizerInZT.h"
#include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
#include "RecoVertex/AdaptiveVertexFit/interface/AdaptiveVertexFitter.h"
//#include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
Expand All @@ -54,10 +55,10 @@

class PrimaryVertexProducer : public edm::stream::EDProducer<> {
public:
explicit PrimaryVertexProducer(const edm::ParameterSet&);
~PrimaryVertexProducer();
PrimaryVertexProducer(const edm::ParameterSet&);
~PrimaryVertexProducer() override;

virtual void produce(edm::Event&, const edm::EventSetup&);
void produce(edm::Event&, const edm::EventSetup&) override;

// access to config
edm::ParameterSet config() const { return theConfig; }
Expand Down
Loading

0 comments on commit 8dd9a2c

Please sign in to comment.