Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:21:28

0001 #ifndef PixelTrackFitting_RZLine_H
0002 #define PixelTrackFitting_RZLine_H
0003 
0004 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0005 #include "DataFormats/GeometryCommonDetAlgo/interface/GlobalError.h"
0006 #include "CommonTools/Utils/interface/DynArray.h"
0007 
0008 #include "CommonTools/Statistics/interface/LinearFit.h"
0009 
0010 #include <vector>
0011 
0012 class RZLine {
0013 public:
0014   struct ErrZ2_tag {};
0015 
0016   /**
0017    * Constructor for containers of GlobalPoint, GlobalError, and bool
0018    *
0019    * @tparam P  Container of GlobalPoint
0020    * @tparam E  Container of GlobalError
0021    * @tparam B  Container of bool
0022    *
0023    * Container can be e.g. std::vector, std::array, or DynArray.
0024    *
0025    * Although for std::array use this constructor could be specialized
0026    * to use std::array instead of DynArray for temporary storage.
0027    */
0028   template <typename P, typename E, typename B>
0029   RZLine(const P& points, const E& errors, const B& isBarrel) {
0030     const size_t n = points.size();
0031     const size_t nSafe = n > 0 ? n : 1;
0032     float r[nSafe];
0033     float z[nSafe];
0034     float errZ2[nSafe];
0035     for (size_t i = 0; i < n; ++i) {
0036       const GlobalPoint& p = points[i];
0037       r[i] = p.perp();
0038       z[i] = p.z();
0039     }
0040 
0041     float simpleCot2 = n > 1 ? sqr((z[n - 1] - z[0]) / (r[n - 1] - r[0])) : 0.f;
0042     for (size_t i = 0; i < n; ++i) {
0043       errZ2[i] = (isBarrel[i]) ? errors[i].czz() : errors[i].rerr(points[i]) * simpleCot2;
0044     }
0045 
0046     calculate(r, z, errZ2, n);
0047   }
0048 
0049   /**
0050    * Constructor for std::vector of r, z, and z standard deviation
0051    */
0052   RZLine(const std::vector<float>& r, const std::vector<float>& z, const std::vector<float>& errZ) {
0053     const size_t n = errZ.size();
0054     float errZ2[n > 0 ? n : n + 1];
0055     for (size_t i = 0; i < n; ++i)
0056       errZ2[i] = sqr(errZ[i]);
0057     calculate(r.data(), z.data(), errZ2, n);
0058   }
0059 
0060   /**
0061    * Constructor for std::array of r, z, and z standard deviation
0062    */
0063   template <size_t N>
0064   RZLine(const std::array<float, N>& r, const std::array<float, N>& z, const std::array<float, N>& errZ) {
0065     std::array<float, N> errZ2;
0066     for (size_t i = 0; i < N; ++i)
0067       errZ2[i] = sqr(errZ[i]);
0068     calculate(r.data(), z.data(), errZ2.data(), N);
0069   }
0070 
0071   /**
0072    * Constructor for container of r, z, and z variance
0073    *
0074    * @tparam T  Container of float
0075    *
0076    * Container can be e.g. std::vector, std::array, or DynArray.
0077    *
0078    * The ErrZ2_tag parameter is used to distinguish this constructor
0079    * from other 3-parameter constructors.
0080    *
0081    * Passing variance is useful in cases where it is already available
0082    * to avoid making a square of a square root.
0083    */
0084   template <typename T>
0085   RZLine(const T& r, const T& z, const T& errZ2, ErrZ2_tag) {
0086     calculate(r.data(), z.data(), errZ2.data(), r.size());
0087   }
0088 
0089   float cotTheta() const { return cotTheta_; }
0090   float intercept() const { return intercept_; }
0091   float covss() const { return covss_; }
0092   float covii() const { return covii_; }
0093   float covsi() const { return covsi_; }
0094 
0095   float chi2() const { return chi2_; }
0096 
0097 private:
0098   template <typename T>
0099   void calculate(const T* r, const T* z, const T* errZ2, size_t n) {
0100     //Have 0 and 1 cases return same value
0101     if (n < 2) [[unlikely]] {
0102       n = 0;
0103     }
0104     linearFit(r, z, n, errZ2, cotTheta_, intercept_, covss_, covii_, covsi_);
0105     chi2_ = 0.f;
0106     for (size_t i = 0; i < n; ++i) {
0107       chi2_ += sqr(((z[i] - intercept_) - cotTheta_ * r[i])) / errZ2[i];
0108     }
0109   }
0110 
0111   template <typename T>
0112   static constexpr T sqr(T t) {
0113     return t * t;
0114   }
0115 
0116   float cotTheta_;
0117   float intercept_;
0118   float covss_;
0119   float covii_;
0120   float covsi_;
0121   float chi2_;
0122 };
0123 #endif