|
|
|
#include <DAClusterizerInZ_vect.h>
Classes | |
| struct | track_t |
| struct | vertex_t |
Public Member Functions | |
| double | beta0 (const double betamax, track_t &tks, vertex_t &y) const |
| std::vector< std::vector < reco::TransientTrack > > | clusterize (const std::vector< reco::TransientTrack > &tracks) const |
| DAClusterizerInZ_vect (const edm::ParameterSet &conf) | |
| void | dump (const double beta, const vertex_t &y, const track_t &tks, const int verbosity=0) const |
| double | Eik (double const &t_z, double const &k_z, double const &t_dz2) const |
| track_t | fill (const std::vector< reco::TransientTrack > &tracks) const |
| double | local_exp (double const &inp) const |
| void | local_exp_list (double *arg_inp, double *arg_out, const int arg_arr_size) const |
| bool | merge (vertex_t &) const |
| bool | merge (vertex_t &y, double &beta) const |
| bool | purge (vertex_t &, track_t &, double &, const double) const |
| bool | split (const double beta, track_t &t, vertex_t &y) const |
| void | splitAll (vertex_t &y) const |
| double | update (double beta, track_t >racks, vertex_t &gvertices, bool useRho0, double &rho0) const |
| std::vector< TransientVertex > | vertices (const std::vector< reco::TransientTrack > &tracks, const int verbosity=0) const |
Private Attributes | |
| float | betamax_ |
| float | betastop_ |
| double | coolingFactor_ |
| double | d0CutOff_ |
| double | dzCutOff_ |
| int | maxIterations_ |
| bool | use_vdt_ |
| bool | useTc_ |
| bool | verbose_ |
| float | vertexSize_ |
Description: separates event tracks into clusters along the beam line
Version which auto-vectorizes with gcc 4.6 or newer
Definition at line 21 of file DAClusterizerInZ_vect.h.
| DAClusterizerInZ_vect::DAClusterizerInZ_vect | ( | const edm::ParameterSet & | conf | ) |
Definition at line 15 of file DAClusterizerInZ_vect.cc.
References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), LogDebug, and builder_last_value_cfg::Tmin.
{
// some defaults to avoid uninitialized variables
verbose_ = conf.getUntrackedParameter<bool> ("verbose", false);
use_vdt_ = conf.getUntrackedParameter<bool> ("use_vdt", false);
betamax_ = 0.1;
betastop_ = 1.0;
coolingFactor_ = 0.6;
maxIterations_ = 100;
vertexSize_ = 0.01; // 0.1 mm
dzCutOff_ = 4.0;
// configure
double Tmin = conf.getParameter<double> ("Tmin");
vertexSize_ = conf.getParameter<double> ("vertexSize");
coolingFactor_ = conf.getParameter<double> ("coolingFactor");
useTc_=true;
if(coolingFactor_<0){
coolingFactor_=-coolingFactor_;
useTc_=false;
}
d0CutOff_ = conf.getParameter<double> ("d0CutOff");
dzCutOff_ = conf.getParameter<double> ("dzCutOff");
maxIterations_ = 100;
if (Tmin == 0) {
LogDebug("DAClusterizerinZ_vectorized") << "DAClusterizerInZ: invalid Tmin" << Tmin
<< " reset do default " << 1. / betamax_ << endl;
} else {
betamax_ = 1. / Tmin;
}
}
Definition at line 301 of file DAClusterizerInZ.cc.
References a, b, i, gen::k, funct::log(), nt, funct::pow(), and w().
{
double T0=0; // max Tc for beta=0
// estimate critical temperature from beta=0 (T=inf)
unsigned int nt=tks.size();
for(vector<vertex_t>::iterator k=y.begin(); k!=y.end(); k++){
// vertex fit at T=inf
double sumwz=0;
double sumw=0;
for(unsigned int i=0; i<nt; i++){
double w=tks[i].pi/tks[i].dz2;
sumwz+=w*tks[i].z;
sumw +=w;
}
k->z=sumwz/sumw;
// estimate Tcrit, eventually do this in the same loop
double a=0, b=0;
for(unsigned int i=0; i<nt; i++){
double dx=tks[i].z-(k->z);
double w=tks[i].pi/tks[i].dz2;
a+=w*pow(dx,2)/tks[i].dz2;
b+=w;
}
double Tc= 2.*a/b; // the critical temperature of this vertex
if(Tc>T0) T0=Tc;
}// vertex loop (normally there should be only one vertex at beta=0)
if (T0>1./betamax){
return betamax/pow(coolingFactor_, int(log(T0*betamax)/log(coolingFactor_))-1 );
}else{
// ensure at least one annealing step
return betamax/coolingFactor_;
}
}
| std::vector<std::vector<reco::TransientTrack> > DAClusterizerInZ_vect::clusterize | ( | const std::vector< reco::TransientTrack > & | tracks | ) | const [virtual] |
Implements TrackClusterizerInZ.
| void DAClusterizerInZ_vect::dump | ( | const double | beta, |
| const vertex_t & | y, | ||
| const track_t & | tks, | ||
| const int | verbosity = 0 |
||
| ) | const |
Definition at line 546 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::track_t::_dz2, DAClusterizerInZ_vect::track_t::_pi, DAClusterizerInZ_vect::vertex_t::_pk, DAClusterizerInZ_vect::track_t::_z, DAClusterizerInZ_vect::vertex_t::_z, DAClusterizerInZ_vect::track_t::_Z_sum, beta, Measurement1D::error(), funct::exp(), F(), DAClusterizerInZ_vect::track_t::GetSize(), DAClusterizerInZ_vect::vertex_t::GetSize(), reco::TrackBase::highPurity, i, funct::log(), LogDebug, nt, AlCaHLTBitMon_ParallelJobs::p, mathSSE::sqrt(), DAClusterizerInZ_vect::track_t::tt, and Measurement1D::value().
{
const unsigned int nv = y.GetSize();
const unsigned int nt = tks.GetSize();
LogDebug("DAClusterizerinZ_vectorized") << "-----DAClusterizerInZ::dump ----" << endl;
LogDebug("DAClusterizerinZ_vectorized") << "beta=" << beta << " betamax= " << betamax_ << endl;
LogDebug("DAClusterizerinZ_vectorized") << " z= ";
LogDebug("DAClusterizerinZ_vectorized") << setprecision(4);
for (unsigned int ivertex = 0; ivertex < nv; ++ ivertex) {
LogDebug("DAClusterizerinZ_vectorized") << setw(8) << fixed << y._z[ivertex];
}
LogDebug("DAClusterizerinZ_vectorized") << endl << "T=" << setw(15) << 1. / beta
<< " Tc= ";
LogDebug("DAClusterizerinZ_vectorized") << endl
<< " pk=";
double sumpk = 0;
for (unsigned int ivertex = 0; ivertex < nv; ++ ivertex) {
LogDebug("DAClusterizerinZ_vectorized") << setw(8) << setprecision(3) << fixed << y._pk[ivertex];
sumpk += y._pk[ivertex];
}
LogDebug("DAClusterizerinZ_vectorized") << endl;
if (verbosity > 0) {
double E = 0, F = 0;
LogDebug("DAClusterizerinZ_vectorized") << endl;
LogDebug("DAClusterizerinZ_vectorized")
<< "---- z +/- dz ip +/-dip pt phi eta weights ----"
<< endl;
LogDebug("DAClusterizerinZ_vectorized") << setprecision(4);
for (unsigned int i = 0; i < nt; i++) {
if (tks._Z_sum[i] > 0) {
F -= log(tks._Z_sum[i]) / beta;
}
double tz = tks._z[i];
LogDebug("DAClusterizerinZ_vectorized") << setw(3) << i << ")" << setw(8) << fixed << setprecision(4)
<< tz << " +/-" << setw(6) << sqrt(tks._dz2[i]);
if (tks.tt[i]->track().quality(reco::TrackBase::highPurity)) {
LogDebug("DAClusterizerinZ_vectorized") << " *";
} else {
LogDebug("DAClusterizerinZ_vectorized") << " ";
}
if (tks.tt[i]->track().hitPattern().hasValidHitInFirstPixelBarrel()) {
LogDebug("DAClusterizerinZ_vectorized") << "+";
} else {
LogDebug("DAClusterizerinZ_vectorized") << "-";
}
LogDebug("DAClusterizerinZ_vectorized") << setw(1)
<< tks.tt[i]->track().hitPattern().pixelBarrelLayersWithMeasurement(); // see DataFormats/TrackReco/interface/HitPattern.h
LogDebug("DAClusterizerinZ_vectorized") << setw(1)
<< tks.tt[i]->track().hitPattern().pixelEndcapLayersWithMeasurement();
LogDebug("DAClusterizerinZ_vectorized") << setw(1) << hex
<< tks.tt[i]->track().hitPattern().trackerLayersWithMeasurement()
- tks.tt[i]->track().hitPattern().pixelLayersWithMeasurement()
<< dec;
LogDebug("DAClusterizerinZ_vectorized") << "=" << setw(1) << hex
<< tks.tt[i]->track().trackerExpectedHitsOuter().numberOfHits()
<< dec;
Measurement1D IP =
tks.tt[i]->stateAtBeamLine().transverseImpactParameter();
LogDebug("DAClusterizerinZ_vectorized") << setw(8) << IP.value() << "+/-" << setw(6) << IP.error();
LogDebug("DAClusterizerinZ_vectorized") << " " << setw(6) << setprecision(2)
<< tks.tt[i]->track().pt() * tks.tt[i]->track().charge();
LogDebug("DAClusterizerinZ_vectorized") << " " << setw(5) << setprecision(2)
<< tks.tt[i]->track().phi() << " " << setw(5)
<< setprecision(2) << tks.tt[i]->track().eta();
double sump = 0.;
for (unsigned int ivertex = 0; ivertex < nv; ++ ivertex) {
if ((tks._pi[i] > 0) && (tks._Z_sum[i] > 0)) {
//double p=pik(beta,tks[i],*k);
double p = y._pk[ivertex] * exp(-beta * Eik(tks._z[i], y._z[ivertex], tks._dz2[i])) / tks._Z_sum[i];
if (p > 0.0001) {
LogDebug("DAClusterizerinZ_vectorized") << setw(8) << setprecision(3) << p;
} else {
LogDebug("DAClusterizerinZ_vectorized") << " . ";
}
E += p * Eik(tks._z[i], y._z[ivertex], tks._dz2[i]);
sump += p;
} else {
LogDebug("DAClusterizerinZ_vectorized") << " ";
}
}
LogDebug("DAClusterizerinZ_vectorized") << endl;
}
LogDebug("DAClusterizerinZ_vectorized") << endl << "T=" << 1 / beta << " E=" << E << " n=" << y.GetSize()
<< " F= " << F << endl << "----------" << endl;
}
}
| double DAClusterizerInZ_vect::Eik | ( | double const & | t_z, |
| double const & | k_z, | ||
| double const & | t_dz2 | ||
| ) | const |
Definition at line 110 of file DAClusterizerInZ_vect.cc.
References funct::pow().
{
return pow(t_z - k_z, 2) / t_dz2;
}
| track_t DAClusterizerInZ_vect::fill | ( | const std::vector< reco::TransientTrack > & | tracks | ) | const |
| double DAClusterizerInZ_vect::local_exp | ( | double const & | inp | ) | const [inline] |
Definition at line 50 of file DAClusterizerInZ_vect.cc.
References funct::exp().
| void DAClusterizerInZ_vect::local_exp_list | ( | double * | arg_inp, |
| double * | arg_out, | ||
| const int | arg_arr_size | ||
| ) | const [inline] |
Definition at line 58 of file DAClusterizerInZ_vect.cc.
References funct::exp(), and i.
| bool DAClusterizerInZ_vect::merge | ( | vertex_t & | y, |
| double & | beta | ||
| ) | const |
Definition at line 261 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::vertex_t::_pk, DAClusterizerInZ_vect::vertex_t::_sw, DAClusterizerInZ_vect::vertex_t::_swE, DAClusterizerInZ_vect::vertex_t::_z, DAClusterizerInZ_vect::vertex_t::GetSize(), gen::k, funct::pow(), DAClusterizerInZ_vect::vertex_t::RemoveItem(), and rho.
{
// merge clusters that collapsed or never separated,
// only merge if the estimated critical temperature of the merged vertex is below the current temperature
// return true if vertices were merged, false otherwise
const unsigned int nv = y.GetSize();
if (nv < 2)
return false;
for (unsigned int k = 0; (k + 1) < nv; k++) {
if (fabs(y._z[k + 1] - y._z[k]) < 2.e-2) {
double rho=y._pk[k] + y._pk[k+1];
double swE=y._swE[k]+y._swE[k+1] - y._pk[k]*y._pk[k+1] /rho *pow(y._z[k+1]-y._z[k],2);
double Tc=2*swE/(y._sw[k]+y._sw[k+1]);
if(Tc*beta<1){
if(rho>0){
y._z[k] = (y._pk[k]*y._z[k] + y._pk[k+1]*y._z[k + 1])/rho;
}else{
y._z[k] = 0.5 * (y._z[k] + y._z[k + 1]);
}
y._pk[k] = rho;
y._sw[k]+=y._sw[k+1];
y._swE[k]=swE;
y.RemoveItem(k+1);
return true;
}
}
}
return false;
}
| bool DAClusterizerInZ_vect::merge | ( | vertex_t & | y | ) | const |
Definition at line 299 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::vertex_t::_pk, DAClusterizerInZ_vect::vertex_t::_z, DAClusterizerInZ_vect::vertex_t::GetSize(), gen::k, LogDebug, and DAClusterizerInZ_vect::vertex_t::RemoveItem().
{
// merge clusters that collapsed or never separated, return true if vertices were merged, false otherwise
const unsigned int nv = y.GetSize();
if (nv < 2)
return false;
for (unsigned int k = 0; (k + 1) < nv; k++) {
//if ((k+1)->z - k->z<1.e-2){ // note, no fabs here, maintains z-ordering (with split()+merge() at every temperature)
if (fabs(y._z[k + 1] - y._z[k]) < 1.e-2) { // with fabs if only called after freeze-out (splitAll() at highter T)
y._pk[k] += y._pk[k + 1];
y._z[k] = 0.5 * (y._z[k] + y._z[k + 1]);
y.RemoveItem(k + 1);
if ( verbose_)
{
LogDebug("DAClusterizerinZ_vectorized") << "Merging vertices k = " << k << std::endl;
}
return true;
}
}
return false;
}
| bool DAClusterizerInZ_vect::purge | ( | vertex_t & | y, |
| track_t & | tks, | ||
| double & | rho0, | ||
| const double | beta | ||
| ) | const |
Definition at line 327 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::track_t::_dz2, DAClusterizerInZ_vect::track_t::_pi, DAClusterizerInZ_vect::vertex_t::_pk, DAClusterizerInZ_vect::track_t::_z, DAClusterizerInZ_vect::vertex_t::_z, DAClusterizerInZ_vect::track_t::_Z_sum, DAClusterizerInZ_vect::track_t::GetSize(), DAClusterizerInZ_vect::vertex_t::GetSize(), i, gen::k, SiStripMonitorClusterAlca_cfi::k0, LogDebug, nt, AlCaHLTBitMon_ParallelJobs::p, and DAClusterizerInZ_vect::vertex_t::RemoveItem().
{
// eliminate clusters with only one significant/unique track
const unsigned int nv = y.GetSize();
const unsigned int nt = tks.GetSize();
if (nv < 2)
return false;
double sumpmin = nt;
unsigned int k0 = nv;
for (unsigned int k = 0; k < nv; k++) {
int nUnique = 0;
double sump = 0;
double pmax = y._pk[k] / (y._pk[k] + rho0 * local_exp(-beta * dzCutOff_
* dzCutOff_));
for (unsigned int i = 0; i < nt; i++) {
if (tks._Z_sum[i] > 0) {
//double p=pik(beta,tks[i],*k);
double p = y._pk[k] * local_exp(-beta * Eik(tks._z[i], y._z[k],
tks._dz2[i])) / tks._Z_sum[i];
sump += p;
if ((p > 0.9 * pmax) && (tks._pi[i] > 0)) {
nUnique++;
}
}
}
if ((nUnique < 2) && (sump < sumpmin)) {
sumpmin = sump;
k0 = k;
}
}
if (k0 != nv) {
if (verbose_) {
LogDebug("DAClusterizerinZ_vectorized") << "eliminating prototype at " << y._z[k0] << " with sump="
<< sumpmin << endl;
}
//rho0+=k0->pk;
y.RemoveItem(k0);
return true;
} else {
return false;
}
}
Definition at line 420 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::track_t::_dz2, DAClusterizerInZ_vect::vertex_t::_pk, DAClusterizerInZ_vect::vertex_t::_sw, DAClusterizerInZ_vect::vertex_t::_swE, DAClusterizerInZ_vect::track_t::_z, DAClusterizerInZ_vect::vertex_t::_z, DAClusterizerInZ_vect::track_t::_Z_sum, hlt_jetmet_dqm_QT_fromfile_cfg::critical, alignCSCRings::e, epsilon, DAClusterizerInZ_vect::track_t::GetSize(), DAClusterizerInZ_vect::vertex_t::GetSize(), i, DAClusterizerInZ_vect::vertex_t::InsertItem(), gen::k, nt, AlCaHLTBitMon_ParallelJobs::p, p1, p2, edm::second(), split, w(), and w2.
{
// split only critical vertices (Tc >~ T=1/beta <==> beta*Tc>~1)
// an update must have been made just before doing this (same beta, no merging)
// returns true if at least one cluster was split
double epsilon=1e-3; // split all single vertices by 10 um
unsigned int nv = y.GetSize();
// avoid left-right biases by splitting highest Tc first
std::vector<std::pair<double, unsigned int> > critical;
for(unsigned int k=0; k<nv; k++){
double Tc= 2*y._swE[k]/y._sw[k];
if (beta*Tc > 1.){
critical.push_back( make_pair(Tc, k));
}
}
if (critical.size()==0) return false;
stable_sort(critical.begin(), critical.end(), std::greater<std::pair<double, unsigned int> >() );
bool split=false;
const unsigned int nt = tks.GetSize();
for(unsigned int ic=0; ic<critical.size(); ic++){
unsigned int k=critical[ic].second;
// estimate subcluster positions and weight
double p1=0, z1=0, w1=0;
double p2=0, z2=0, w2=0;
for(unsigned int i=0; i<nt; i++){
if (tks._Z_sum[i] > 0) {
double p = y._pk[k] * local_exp(-beta * Eik(tks._z[i], y._z[k],
tks._dz2[i])) / tks._Z_sum[i];
double w=p/tks._dz2[i];
if(tks._z[i]<y._z[k]){
p1+=p; z1+=w*tks._z[i]; w1+=w;
}else{
p2+=p; z2+=w*tks._z[i]; w2+=w;
}
}
}
if(w1>0){ z1=z1/w1;} else{z1=y._z[k]-epsilon;}
if(w2>0){ z2=z2/w2;} else{z2=y._z[k]+epsilon;}
// reduce split size if there is not enough room
if( ( k > 0 ) && ( y._z[k-1]>=z1 ) ){ z1=0.5*(y._z[k]+y._z[k-1]); }
if( ( k+1 < nv) && ( y._z[k+1]<=z2 ) ){ z2=0.5*(y._z[k]+y._z[k+1]); }
// split if the new subclusters are significantly separated
if( (z2-z1)>epsilon){
split=true;
double pk1=p1*y._pk[k]/(p1+p2);
double pk2=p2*y._pk[k]/(p1+p2);
y._z[k] = z2;
y._pk[k] = pk2;
y.InsertItem(k, z1, pk1);
nv++;
// adjust remaining pointers
for(unsigned int jc=ic; jc<critical.size(); jc++){
if (critical[jc].second>k) {critical[jc].second++;}
}
}
}
return split;
}
| void DAClusterizerInZ_vect::splitAll | ( | vertex_t & | y | ) | const |
Definition at line 422 of file DAClusterizerInZ.cc.
References alignCSCRings::e, epsilon, gen::k, DAClusterizerInZ_vect::vertex_t::pk, and DAClusterizerInZ_vect::vertex_t::z.
{
double epsilon=1e-3; // split all single vertices by 10 um
double zsep=2*epsilon; // split vertices that are isolated by at least zsep (vertices that haven't collapsed)
vector<vertex_t> y1;
for(vector<vertex_t>::iterator k=y.begin(); k!=y.end(); k++){
if ( ( (k==y.begin())|| (k-1)->z < k->z - zsep) && (((k+1)==y.end() )|| (k+1)->z > k->z + zsep)) {
// isolated prototype, split
vertex_t vnew;
vnew.z = k->z -epsilon;
(*k).z = k->z+epsilon;
vnew.pk= 0.5* (*k).pk;
(*k).pk= 0.5* (*k).pk;
y1.push_back(vnew);
y1.push_back(*k);
}else if( y1.empty() || (y1.back().z < k->z -zsep)){
y1.push_back(*k);
}else{
y1.back().z -=epsilon;
k->z+=epsilon;
y1.push_back(*k);
}
}// vertex loop
y=y1;
}
| double DAClusterizerInZ_vect::update | ( | double | beta, |
| track_t & | gtracks, | ||
| vertex_t & | gvertices, | ||
| bool | useRho0, | ||
| double & | rho0 | ||
| ) | const |
Definition at line 116 of file DAClusterizerInZ_vect.cc.
References DAClusterizerInZ_vect::vertex_t::_ei, DAClusterizerInZ_vect::vertex_t::_ei_cache, DAClusterizerInZ_vect::track_t::_pi, DAClusterizerInZ_vect::vertex_t::_se, DAClusterizerInZ_vect::vertex_t::_sw, DAClusterizerInZ_vect::vertex_t::_swE, DAClusterizerInZ_vect::vertex_t::_swz, DAClusterizerInZ_vect::track_t::_Z_sum, beta, delta, DAClusterizerInZ_vect::track_t::GetSize(), DAClusterizerInZ_vect::vertex_t::GetSize(), gen::k, LogDebug, nt, funct::pow(), testEve_cfg::tracks, and w().
{
//update weights and vertex positions
// mass constrained annealing without noise
// returns the squared sum of changes of vertex positions
const unsigned int nt = gtracks.GetSize();
const unsigned int nv = gvertices.GetSize();
//initialize sums
double sumpi = 0;
// to return how much the prototype moved
double delta = 0;
unsigned int itrack;
unsigned int ivertex;
// intial value of a sum
double Z_init = 0;
// define kernels
auto kernel_calc_exp_arg = [ &beta, nv ] (
const unsigned int itrack,
track_t const& tracks,
vertex_t const& vertices )
{
const double track_z = tracks._z[itrack];
const double track_dz2 = tracks._dz2[itrack];
// auto-vectorized
for ( unsigned int ivertex = 0; ivertex < nv; ++ivertex)
{
double mult_res = ( track_z - vertices._z[ivertex] );
vertices._ei_cache[ivertex] = -beta *
( mult_res * mult_res ) / track_dz2;
}
};
// auto kernel_add_Z = [ nv, &Z_init ] (vertex_t const& vertices) raises a warning
// we declare the return type with " -> double"
auto kernel_add_Z = [ nv, &Z_init ] (vertex_t const& vertices) -> double
{
double ZTemp = Z_init;
for (unsigned int ivertex = 0; ivertex < nv; ++ivertex) {
ZTemp += vertices._pk[ivertex] * vertices._ei[ivertex];
}
return ZTemp;
};
auto kernel_calc_normalization = [ &beta, nv ] (const unsigned int track_num,
track_t & tks_vec,
vertex_t & y_vec )
{
double w;
const double tmp_trk_pi = tks_vec._pi[track_num];
const double tmp_trk_Z_sum = tks_vec._Z_sum[track_num];
const double tmp_trk_dz2 = tks_vec._dz2[track_num];
const double tmp_trk_z = tks_vec._z[track_num];
// auto-vectorized
for (unsigned int k = 0; k < nv; ++k) {
y_vec._se[k] += tmp_trk_pi* y_vec._ei[k] / tmp_trk_Z_sum;
w = y_vec._pk[k] * tmp_trk_pi * y_vec._ei[k] / tmp_trk_Z_sum / tmp_trk_dz2;
y_vec._sw[k] += w;
y_vec._swz[k] += w * tmp_trk_z;
y_vec._swE[k] += w * y_vec._ei_cache[k]/(-beta);
}
};
// not vectorized
for (ivertex = 0; ivertex < nv; ++ivertex) {
gvertices._se[ivertex] = 0.0;
gvertices._sw[ivertex] = 0.0;
gvertices._swz[ivertex] = 0.0;
gvertices._swE[ivertex] = 0.0;
}
// independpent of loop
if ( useRho0 )
{
Z_init = rho0 * local_exp(-beta * dzCutOff_ * dzCutOff_); // cut-off
}
// loop over tracks
for (itrack = 0; itrack < nt; ++itrack) {
kernel_calc_exp_arg(itrack, gtracks, gvertices);
local_exp_list(gvertices._ei_cache, gvertices._ei, nv);
//vdt::cephes_single_exp_vect( y_vec._ei, nv );
gtracks._Z_sum[itrack] = kernel_add_Z(gvertices);
// used in the next major loop to follow
if (!useRho0)
sumpi += gtracks._pi[itrack];
if (gtracks._Z_sum[itrack] > 0) {
kernel_calc_normalization(itrack, gtracks, gvertices);
}
}
// now update z and pk
auto kernel_calc_z = [ &delta, &sumpi, nv, this, useRho0 ] (vertex_t & vertices )
{
// does not vectorizes
for (unsigned int ivertex = 0; ivertex < nv; ++ ivertex )
{
if (vertices._sw[ivertex] > 0)
{
double znew = vertices._swz[ ivertex ] / vertices._sw[ ivertex ];
// prevents from vectorizing
delta += pow( vertices._z[ ivertex ] - znew, 2 );
vertices._z[ ivertex ] = znew;
}
else {
edm::LogInfo("sumw") << "invalid sum of weights in fit: " << vertices._sw[ivertex]
<< endl;
if (this->verbose_) {
LogDebug("DAClusterizerinZ_vectorized") << " a cluster melted away ? pk=" << vertices._pk[ ivertex ] << " sumw="
<< vertices._sw[ivertex] << endl;
}
}
// dont do, if rho cut
if ( ! useRho0 )
{
vertices._pk[ ivertex ] = vertices._pk[ ivertex ] * vertices._se[ ivertex ] / sumpi;
}
}
};
kernel_calc_z(gvertices);
// return how much the prototypes moved
return delta;
}
| std::vector<TransientVertex> DAClusterizerInZ_vect::vertices | ( | const std::vector< reco::TransientTrack > & | tracks, |
| const int | verbosity = 0 |
||
| ) | const |
float DAClusterizerInZ_vect::betamax_ [private] |
Definition at line 206 of file DAClusterizerInZ_vect.h.
float DAClusterizerInZ_vect::betastop_ [private] |
Definition at line 207 of file DAClusterizerInZ_vect.h.
double DAClusterizerInZ_vect::coolingFactor_ [private] |
Definition at line 205 of file DAClusterizerInZ_vect.h.
double DAClusterizerInZ_vect::d0CutOff_ [private] |
Definition at line 209 of file DAClusterizerInZ_vect.h.
double DAClusterizerInZ_vect::dzCutOff_ [private] |
Definition at line 208 of file DAClusterizerInZ_vect.h.
int DAClusterizerInZ_vect::maxIterations_ [private] |
Definition at line 204 of file DAClusterizerInZ_vect.h.
bool DAClusterizerInZ_vect::use_vdt_ [private] |
Definition at line 210 of file DAClusterizerInZ_vect.h.
bool DAClusterizerInZ_vect::useTc_ [private] |
Definition at line 211 of file DAClusterizerInZ_vect.h.
bool DAClusterizerInZ_vect::verbose_ [private] |
Definition at line 202 of file DAClusterizerInZ_vect.h.
float DAClusterizerInZ_vect::vertexSize_ [private] |
Definition at line 203 of file DAClusterizerInZ_vect.h.
1.7.3