File indexing completed on 2024-04-06 12:29:15
0001 #ifndef RecoVertex_PrimaryVertexProducer_DAClusterizerInZ_vect_h
0002 #define RecoVertex_PrimaryVertexProducer_DAClusterizerInZ_vect_h
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "RecoVertex/PrimaryVertexProducer/interface/TrackClusterizerInZ.h"
0013 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0014 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include <vector>
0017 #include "DataFormats/Math/interface/Error.h"
0018 #include "RecoVertex/VertexTools/interface/VertexDistanceXY.h"
0019 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0020
0021 class DAClusterizerInZ_vect final : public TrackClusterizerInZ {
0022 public:
0023 static void fillPSetDescription(edm::ParameterSetDescription &desc);
0024
0025
0026 struct track_t {
0027 std::vector<double> zpca_vec;
0028 std::vector<double> dz2_vec;
0029 std::vector<double> sum_Z_vec;
0030 std::vector<double> tkwt_vec;
0031 std::vector<unsigned int> kmin;
0032 std::vector<unsigned int> kmax;
0033 std::vector<const reco::TransientTrack *> tt;
0034
0035 double osumtkwt;
0036
0037 void addItemSorted(double new_zpca, double new_dz2, const reco::TransientTrack *new_tt, double new_tkwt) {
0038
0039 unsigned int i = 0;
0040 for (i = 0; i < zpca_vec.size(); i++) {
0041 if (new_dz2 > dz2_vec[i])
0042 break;
0043 }
0044 insertItem(i, new_zpca, new_dz2, new_tt, new_tkwt);
0045 }
0046
0047 void insertItem(
0048 unsigned int i, double new_zpca, double new_dz2, const reco::TransientTrack *new_tt, double new_tkwt) {
0049 zpca_vec.insert(zpca_vec.begin() + i, new_zpca);
0050 dz2_vec.insert(dz2_vec.begin() + i, new_dz2);
0051 tt.insert(tt.begin() + i, new_tt);
0052 tkwt_vec.insert(tkwt_vec.begin() + i, new_tkwt);
0053 sum_Z_vec.insert(sum_Z_vec.begin() + i, 1.0);
0054 kmin.insert(kmin.begin() + i, 0);
0055 kmax.insert(kmax.begin() + i, 0);
0056 }
0057
0058 unsigned int getSize() const { return zpca_vec.size(); }
0059
0060
0061 void extractRaw() {
0062 zpca = &zpca_vec.front();
0063 dz2 = &dz2_vec.front();
0064 tkwt = &tkwt_vec.front();
0065 sum_Z = &sum_Z_vec.front();
0066 }
0067
0068
0069 double *__restrict__ zpca;
0070 double *__restrict__ dz2;
0071 double *__restrict__ tkwt;
0072 double *__restrict__ sum_Z;
0073 };
0074
0075
0076 struct vertex_t {
0077 std::vector<double> zvtx_vec;
0078 std::vector<double> rho_vec;
0079
0080 std::vector<double> exp_arg_vec;
0081 std::vector<double> exp_vec;
0082 std::vector<double> sw_vec;
0083 std::vector<double> swz_vec;
0084 std::vector<double> se_vec;
0085 std::vector<double> swE_vec;
0086
0087 unsigned int getSize() const { return zvtx_vec.size(); }
0088
0089 void addItem(double new_zvtx, double new_rho) {
0090 zvtx_vec.push_back(new_zvtx);
0091 rho_vec.push_back(new_rho);
0092 exp_arg_vec.push_back(0.0);
0093 exp_vec.push_back(0.0);
0094 sw_vec.push_back(0.0);
0095 swz_vec.push_back(0.0);
0096 se_vec.push_back(0.0);
0097 swE_vec.push_back(0.0);
0098
0099 extractRaw();
0100 }
0101
0102 void insertItem(unsigned int k, double new_zvtx, double new_rho, track_t &tks) {
0103 zvtx_vec.insert(zvtx_vec.begin() + k, new_zvtx);
0104 rho_vec.insert(rho_vec.begin() + k, new_rho);
0105
0106 exp_arg_vec.insert(exp_arg_vec.begin() + k, 0.0);
0107 exp_vec.insert(exp_vec.begin() + k, 0.0);
0108 sw_vec.insert(sw_vec.begin() + k, 0.0);
0109 swz_vec.insert(swz_vec.begin() + k, 0.0);
0110 se_vec.insert(se_vec.begin() + k, 0.0);
0111 swE_vec.insert(swE_vec.begin() + k, 0.0);
0112
0113
0114 for (unsigned int i = 0; i < tks.getSize(); i++) {
0115 if (tks.kmin[i] > k) {
0116 tks.kmin[i]++;
0117 }
0118 if ((tks.kmax[i] >= k) || (tks.kmax[i] == tks.kmin[i])) {
0119 tks.kmax[i]++;
0120 }
0121 }
0122
0123 extractRaw();
0124 }
0125
0126 void removeItem(unsigned int k, track_t &tks) {
0127 zvtx_vec.erase(zvtx_vec.begin() + k);
0128 rho_vec.erase(rho_vec.begin() + k);
0129
0130 exp_arg_vec.erase(exp_arg_vec.begin() + k);
0131 exp_vec.erase(exp_vec.begin() + k);
0132 sw_vec.erase(sw_vec.begin() + k);
0133 swz_vec.erase(swz_vec.begin() + k);
0134 se_vec.erase(se_vec.begin() + k);
0135 swE_vec.erase(swE_vec.begin() + k);
0136
0137
0138 for (unsigned int i = 0; i < tks.getSize(); i++) {
0139 if (tks.kmax[i] > k) {
0140 tks.kmax[i]--;
0141 }
0142 if ((tks.kmin[i] > k) || (((tks.kmax[i] < (tks.kmin[i] + 1)) && (tks.kmin[i] > 0)))) {
0143 tks.kmin[i]--;
0144 }
0145 }
0146
0147 extractRaw();
0148 }
0149
0150
0151 double *__restrict__ zvtx;
0152 double *__restrict__ rho;
0153 double *__restrict__ exp_arg;
0154 double *__restrict__ exp;
0155 double *__restrict__ sw;
0156 double *__restrict__ swz;
0157 double *__restrict__ se;
0158 double *__restrict__ swE;
0159
0160
0161 void extractRaw() {
0162 zvtx = &zvtx_vec.front();
0163 rho = &rho_vec.front();
0164 exp = &exp_vec.front();
0165 sw = &sw_vec.front();
0166 swz = &swz_vec.front();
0167 se = &se_vec.front();
0168 swE = &swE_vec.front();
0169 exp_arg = &exp_arg_vec.front();
0170 }
0171 };
0172
0173 DAClusterizerInZ_vect(const edm::ParameterSet &conf);
0174
0175 std::vector<std::vector<reco::TransientTrack> > clusterize(
0176 const std::vector<reco::TransientTrack> &tracks) const override;
0177
0178 std::vector<TransientVertex> vertices(const std::vector<reco::TransientTrack> &tracks) const override;
0179 std::vector<TransientVertex> vertices_no_blocks(const std::vector<reco::TransientTrack> &tracks) const;
0180 std::vector<TransientVertex> vertices_in_blocks(const std::vector<reco::TransientTrack> &tracks) const;
0181 std::vector<TransientVertex> fill_vertices(double beta, double rho0, track_t &tracks, vertex_t &vertices) const;
0182
0183 track_t fill(const std::vector<reco::TransientTrack> &tracks) const;
0184
0185 void set_vtx_range(double beta, track_t >racks, vertex_t &gvertices) const;
0186
0187 void clear_vtx_range(track_t >racks, vertex_t &gvertices) const;
0188
0189 unsigned int thermalize(
0190 double beta, track_t >racks, vertex_t &gvertices, const double delta_max, const double rho0 = 0.) const;
0191
0192 double update(
0193 double beta, track_t >racks, vertex_t &gvertices, const double rho0 = 0, const bool updateTc = false) const;
0194
0195 void dump(
0196 const double beta, const vertex_t &y, const track_t &tks, const int verbosity = 0, const double rho0 = 0.) const;
0197 bool merge(vertex_t &y, track_t &tks, double &beta) const;
0198 bool purge(vertex_t &, track_t &, double &, const double) const;
0199 bool split(const double beta, track_t &t, vertex_t &y, double threshold = 1.) const;
0200
0201 double beta0(const double betamax, track_t const &tks, vertex_t const &y) const;
0202 void verify(const vertex_t &v, const track_t &tks, unsigned int nv = 999999, unsigned int nt = 999999) const;
0203
0204 private:
0205 double zdumpcenter_;
0206 double zdumpwidth_;
0207
0208 double vertexSize_;
0209 unsigned int maxIterations_;
0210 double coolingFactor_;
0211 double betamax_;
0212 double betastop_;
0213 double dzCutOff_;
0214 double d0CutOff_;
0215
0216 double mintrkweight_;
0217 double uniquetrkweight_;
0218 double uniquetrkminp_;
0219 double zmerge_;
0220 double betapurge_;
0221
0222 unsigned int convergence_mode_;
0223 double delta_highT_;
0224 double delta_lowT_;
0225
0226 double sel_zrange_;
0227 const double zrange_min_ = 0.1;
0228
0229 bool runInBlocks_;
0230 unsigned int block_size_;
0231 double overlap_frac_;
0232 };
0233
0234
0235 #endif