Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-04-01 23:54:10

0001 #ifndef __L1Trigger_VertexFinder_VertexFinder_h__
0002 #define __L1Trigger_VertexFinder_VertexFinder_h__
0003 
0004 #include "DataFormats/Common/interface/Ptr.h"
0005 #include "DataFormats/L1Trigger/interface/VertexWord.h"
0006 #include "DataFormats/Math/interface/deltaPhi.h"
0007 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "L1Trigger/VertexFinder/interface/AlgoSettings.h"
0010 #include "L1Trigger/VertexFinder/interface/RecoVertex.h"
0011 
0012 #include <algorithm>
0013 #include <iterator>
0014 #include <vector>
0015 
0016 namespace l1tVertexFinder {
0017 
0018   // Returns the number of bits needed to represent and integer decimal
0019   static constexpr unsigned BitsToRepresent(unsigned x) { return x < 2 ? 1 : 1 + BitsToRepresent(x >> 1); }
0020 
0021   typedef std::vector<L1Track> FitTrackCollection;
0022   typedef std::vector<RecoVertex<>> RecoVertexCollection;
0023 
0024   class VertexFinder {
0025   public:
0026     /// Constructor and destructor
0027     VertexFinder(FitTrackCollection& fitTracks, const AlgoSettings& settings) {
0028       fitTracks_ = fitTracks;
0029       settings_ = &settings;
0030     }
0031     ~VertexFinder() {}
0032 
0033     /// Helper structs/classes
0034     struct SortTracksByZ0 {
0035       inline bool operator()(const L1Track track0, const L1Track track1) { return (track0.z0() < track1.z0()); }
0036     };
0037 
0038     struct SortTracksByPt {
0039       inline bool operator()(const L1Track track0, const L1Track track1) {
0040         return (std::abs(track0.pt()) > std::abs(track1.pt()));
0041       }
0042     };
0043 
0044     /// Accessors
0045 
0046     /// Storage for tracks out of the L1 Track finder
0047     const FitTrackCollection& fitTracks() const { return fitTracks_; }
0048     /// Number of iterations
0049     unsigned int iterationsPerTrack() const { return double(iterations_) / double(fitTracks_.size()); }
0050     /// Storage for tracks out of the L1 Track finder
0051     unsigned int numInputTracks() const { return fitTracks_.size(); }
0052     /// Number of iterations
0053     unsigned int numIterations() const { return iterations_; }
0054     /// Number of reconstructed vertices
0055     unsigned int numVertices() const { return vertices_.size(); }
0056     /// Number of emulation vertices
0057     unsigned int numVerticesEmulation() const { return verticesEmulation_.size(); }
0058     /// Reconstructed primary vertex
0059     template <typename T>
0060     T PrimaryVertex() const {
0061       if ((settings_->vx_precision() == Precision::Simulation) && (pv_index_ < vertices_.size()))
0062         return vertices_[pv_index_];
0063       else if ((settings_->vx_precision() == Precision::Emulation) && (pv_index_ < vertices_.size()))
0064         return verticesEmulation_[pv_index_];
0065       else {
0066         edm::LogWarning("VertexFinder") << "PrimaryVertex::No primary vertex has been found.";
0067         return RecoVertex<>();
0068       }
0069     }
0070     /// Reconstructed Primary Vertex Id
0071     unsigned int primaryVertexId() const { return pv_index_; }
0072     /// Returns the z positions of the reconstructed primary vertices
0073     const std::vector<RecoVertex<>>& vertices() const { return vertices_; }
0074     /// Returns the emulation primary vertices
0075     const l1t::VertexWordCollection& verticesEmulation() const { return verticesEmulation_; }
0076 
0077     /// Find the primary vertex
0078     void findPrimaryVertex();
0079     /// Associate the primary vertex with the real one
0080     void associatePrimaryVertex(double trueZ0);
0081     /// Gap Clustering Algorithm
0082     void GapClustering();
0083     /// Find maximum distance in two clusters of tracks
0084     float maxDistance(RecoVertex<> cluster0, RecoVertex<> cluster1);
0085     /// Find minimum distance in two clusters of tracks
0086     float minDistance(RecoVertex<> cluster0, RecoVertex<> cluster1);
0087     /// Find average distance in two clusters of tracks
0088     float meanDistance(RecoVertex<> cluster0, RecoVertex<> cluster1);
0089     /// Find distance between centres of two clusters
0090     float centralDistance(RecoVertex<> cluster0, RecoVertex<> cluster1);
0091     /// Simple Merge Algorithm
0092     void agglomerativeHierarchicalClustering();
0093     /// Adaptive Vertex Reconstruction algorithm
0094     void adaptiveVertexReconstruction();
0095     /// High pT Vertex Algorithm
0096     void fastHistoLooseAssociation();
0097     /// Histogramming algorithm
0098     void fastHisto(const TrackerTopology* tTopo);
0099     /// Histogramming algorithm (emulation)
0100     void fastHistoEmulation();
0101 
0102     /// Sort vertices in pT
0103     void sortVerticesInPt();
0104     /// Sort vertices in z
0105     void sortVerticesInZ0();
0106 
0107     /// Print an ASCII histogram
0108     template <class data_type, typename stream_type = std::ostream>
0109     void printHistogram(stream_type& stream,
0110                         std::vector<data_type> data,
0111                         int width = 80,
0112                         int minimum = 0,
0113                         int maximum = -1,
0114                         std::string title = "",
0115                         std::string color = "");
0116 
0117     template <typename ForwardIterator, typename T>
0118     void strided_iota(ForwardIterator first, ForwardIterator last, T value, T stride) {
0119       while (first != last) {
0120         *first++ = value;
0121         value += stride;
0122       }
0123     }
0124 
0125     /// Vertexing algorithms
0126 
0127     /// Compute the vertex parameters
0128     void computeAndSetVertexParameters(RecoVertex<>& vertex,
0129                                        const std::vector<float>& bin_centers,
0130                                        const std::vector<unsigned int>& counts);
0131     /// DBSCAN algorithm
0132     void DBSCAN();
0133     /// High pT Vertex Algorithm
0134     void HPV();
0135     /// Kmeans Algorithm
0136     void Kmeans();
0137     /// Find maximum distance in two clusters of tracks
0138     void PVR();
0139 
0140   private:
0141     const AlgoSettings* settings_;
0142     RecoVertexCollection vertices_;
0143     l1t::VertexWordCollection verticesEmulation_;
0144     unsigned int numMatchedVertices_;
0145     FitTrackCollection fitTracks_;
0146     unsigned int pv_index_;
0147     unsigned int iterations_;
0148   };
0149 
0150 }  // end namespace l1tVertexFinder
0151 
0152 #endif