File indexing completed on 2024-08-21 04:46:46
0001
0002
0003
0004 #include "RecoHGCal/TICL/interface/SuperclusteringDNNInputs.h"
0005
0006 #include <algorithm>
0007 #include <numeric>
0008 #include <cmath>
0009
0010 #include <Math/Vector2D.h>
0011 #include <Math/Vector3D.h>
0012 #include <Math/Rotation3D.h>
0013 #include <Math/VectorUtil.h>
0014
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016 #include "DataFormats/HGCalReco/interface/Trackster.h"
0017
0018 namespace ticl {
0019
0020 std::vector<float> SuperclusteringDNNInputV1::computeVector(Trackster const& ts_base, Trackster const& ts_toCluster) {
0021
0022
0023
0024
0025 return {{
0026 std::abs(ts_toCluster.barycenter().Eta()) - std::abs(ts_base.barycenter().Eta()),
0027 ts_toCluster.barycenter().Phi() - ts_base.barycenter().phi(),
0028 ts_toCluster.raw_energy(),
0029 ts_toCluster.barycenter().Eta(),
0030 (ts_toCluster.raw_energy() * std::sin(ts_toCluster.barycenter().Theta())),
0031 ts_base.barycenter().Eta(),
0032 ts_base.barycenter().Phi(),
0033 ts_base.raw_energy(),
0034 (ts_base.raw_energy() * std::sin(ts_toCluster.barycenter().Theta())),
0035 }};
0036 }
0037
0038
0039 template <class Vector1, class Vector2>
0040 float CosTheta2D(const Vector1& v1, const Vector2& v2) {
0041 float arg;
0042 float v1_r2 = v1.X() * v1.X() + v1.Y() * v1.Y();
0043 float v2_r2 = v2.X() * v2.X() + v2.Y() * v2.Y();
0044 float ptot2 = v1_r2 * v2_r2;
0045 if (ptot2 <= 0) {
0046 arg = 0.0;
0047 } else {
0048 float pdot = v1.X() * v2.X() + v1.Y() * v2.Y();
0049 using std::sqrt;
0050 arg = pdot / sqrt(ptot2);
0051 if (arg > 1.0)
0052 arg = 1.0;
0053 if (arg < -1.0)
0054 arg = -1.0;
0055 }
0056 return arg;
0057 }
0058 template <class Vector1, class Vector2>
0059 inline float Angle2D(const Vector1& v1, const Vector2& v2) {
0060 return static_cast<float>(std::acos(CosTheta2D(v1, v2)));
0061 }
0062
0063 std::vector<float> SuperclusteringDNNInputV2::computeVector(Trackster const& ts_base, Trackster const& ts_toCluster) {
0064 using ROOT::Math::XYVectorF;
0065 using ROOT::Math::XYZVectorF;
0066 using ROOT::Math::VectorUtil::Angle;
0067 XYZVectorF const& pca_seed_cmsFrame(ts_base.eigenvectors(0));
0068 XYZVectorF const& pca_cand_cmsFrame(ts_toCluster.eigenvectors(0));
0069 XYZVectorF xs(pca_seed_cmsFrame.Cross(XYZVectorF(0, 0, 1)).Unit());
0070 ROOT::Math::Rotation3D rot(xs, xs.Cross(pca_seed_cmsFrame).Unit(), pca_seed_cmsFrame);
0071
0072 XYZVectorF pca_cand_seedFrame = rot(pca_cand_cmsFrame);
0073
0074 float explVar_denominator = std::accumulate(
0075 std::begin(ts_toCluster.eigenvalues()), std::end(ts_toCluster.eigenvalues()), 0.f, std::plus<float>());
0076 float explVarRatio = 0.;
0077 if (explVar_denominator != 0.) {
0078 explVarRatio = ts_toCluster.eigenvalues()[0] / explVar_denominator;
0079 } else {
0080 edm::LogWarning("HGCalTICLSuperclustering")
0081 << "Sum of eigenvalues was zero for trackster. Could not compute explained variance ratio.";
0082 }
0083
0084 return {{
0085 std::abs(ts_toCluster.barycenter().Eta()) - std::abs(ts_base.barycenter().Eta()),
0086 ts_toCluster.barycenter().Phi() - ts_base.barycenter().phi(),
0087 ts_toCluster.raw_energy(),
0088 ts_toCluster.barycenter().Eta(),
0089 (ts_toCluster.raw_energy() * std::sin(ts_toCluster.barycenter().Theta())),
0090 ts_base.barycenter().Eta(),
0091 ts_base.barycenter().Phi(),
0092 ts_base.raw_energy(),
0093 (ts_base.raw_energy() * std::sin(ts_toCluster.barycenter().Theta())),
0094 static_cast<float>(Angle(pca_cand_cmsFrame, pca_seed_cmsFrame)),
0095 Angle2D(XYVectorF(pca_cand_seedFrame.x(), pca_cand_seedFrame.z()), XYVectorF(0, 1)),
0096 Angle2D(XYVectorF(pca_cand_seedFrame.y(), pca_cand_seedFrame.z()), XYVectorF(0, 1)),
0097 Angle2D(XYVectorF(pca_cand_cmsFrame.x(), pca_cand_cmsFrame.y()),
0098 XYVectorF(pca_seed_cmsFrame.x(), pca_seed_cmsFrame.y())),
0099 Angle2D(XYVectorF(pca_cand_cmsFrame.y(), pca_cand_cmsFrame.z()),
0100 XYVectorF(pca_seed_cmsFrame.y(), pca_seed_cmsFrame.z())),
0101 Angle2D(XYVectorF(pca_cand_cmsFrame.x(), pca_cand_cmsFrame.z()),
0102 XYVectorF(pca_seed_cmsFrame.x(), pca_seed_cmsFrame.z())),
0103 ts_toCluster.eigenvalues()[0],
0104 explVarRatio
0105 }};
0106 }
0107
0108 std::unique_ptr<AbstractSuperclusteringDNNInput> makeSuperclusteringDNNInputFromString(std::string dnnInputVersion) {
0109 if (dnnInputVersion == "v1")
0110 return std::make_unique<SuperclusteringDNNInputV1>();
0111 else if (dnnInputVersion == "v2")
0112 return std::make_unique<SuperclusteringDNNInputV2>();
0113 assert(false);
0114 }
0115
0116 }