diff --git a/CMakeLists.txt b/CMakeLists.txt
index d61767cd27059897de71b342b0a5ae64cb6b7c71..52ac9db021b8950221e1c903844cc0a56e322fa0 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -1,8 +1,11 @@
-cmake_minimum_required(VERSION 3.8)
+cmake_minimum_required(VERSION 3.12)
project(EICD
VERSION 1.1.0
LANGUAGES CXX)
+set(CMAKE_CXX_STANDARD 20)
+set(CMAKE_CXX_STANDARD_REQUIRED ON)
+
#--- Declare options -----------------------------------------------------------
option(BUILD_DATA_MODEL "Run podio class generator yaml file" ON)
@@ -47,7 +50,7 @@ PODIO_GENERATE_DICTIONARY(eicd ${headers}
set_target_properties(eicd-dictgen PROPERTIES EXCLUDE_FROM_ALL TRUE)
target_sources(eicd PRIVATE eicd.cxx)
-add_subdirectory(source)
+add_subdirectory(utils)
install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/eicd
diff --git a/eic_data.yaml b/eic_data.yaml
index 53cdde0c6125e4e773acb768fe9d38bfc7bdfc08..7932b84e9c481ee339af4b108f32b1fc646626e0 100644
--- a/eic_data.yaml
+++ b/eic_data.yaml
@@ -41,7 +41,7 @@ components:
bool operator==(const eic::Index& rhs) const {return equals(rhs);}\n
bool operator!=(const eic::Index& rhs) const {return !equals(rhs);}\n
bool operator<(const eic::Index& rhs) const {return long_form() < rhs.long_form();}\n
- int64_t long_form() const {int64_t l = static_cast<int64_t>(source) << 32 | value; return l;}\n
+ uint64_t long_form() const {int64_t l = static_cast<int64_t>(source) << 32 | value; return l;}\n
explicit operator bool() const {return valid();}
"
@@ -57,24 +57,6 @@ components:
operator float() const {return value;}
"
-
- ## first-second pair of float s
- eic::FloatPair:
- Members:
- - float first
- - float second
- ExtraCode:
- includes: "#include <tuple>"
- declaration: "
- FloatPair() : first{0}, second{0} {}\n
- FloatPair(double a, double b) : first{static_cast<float>(a)}, second{static_cast<float>(b)} {}\n
- FloatPair(const std::pair<float, float>& p) : first{p.first}, second{p.second} {}\n
- FloatPair& operator=(const std::pair<float, float>& p) {first = p.first; second = p.second; return *this;}\n
- float& operator[](unsigned i) {return *(&first + i);}\n
- const float& operator[](unsigned i) const {return *(&first + i);}\n
- operator std::pair<float, float>() const {return {first, second};}\n
- "
-
eic::VectorXY:
Members:
- float x // [mm] or [GeV]
@@ -96,23 +78,44 @@ components:
VectorXY scale(double f) const {return {f*x, f*y};}\n
"
- eic::Direction:
+ # Copy from EDM4hep with conversion operators added
+ # TODO: uncomment Vector2f conversions when we migrate to new EDM4hep version
+ eicd::Vector2f :
Members:
- - float theta // [rad, 0->pi]
- - float phi // [rad, -pi->pi]
+ - float a
+ - float b
+ ExtraCode:
+ declaration: "
+ Vector2f() : a(0),b(0) {}\n
+ Vector2f(float aa,float bb) : a(aa),b(bb) {}\n
+ Vector2f(const float* v) : a(v[0]), b(v[1]) {}\n
+ bool operator==(const Vector2f& v) const { return (a==v.a&&b==v.b) ; }\n
+ float operator[](unsigned i) const { return *( &a + i ) ; }\n
+ //Vector2f(const edm4hep::Vector2f& v) : a{v.a}, b{v.b} {}\n
+ //operator edm4hep::Vector2f() const {return {a, b};}\n
+ "
+ includes: "
+ // #include <edm4hep/Vector2f.h>\n
+ "
+
+ # Copy from EDM4hep with conversion operators added
+ eicd::Vector3f :
+ Members:
+ - float x
+ - float y
+ - float z
ExtraCode:
- includes: "#include <cmath>\n#include <tuple>"
declaration: "
- Direction() : theta{0}, phi{0} {}\n
- Direction(double th, double ph) : theta{static_cast<float>(th)}, phi{static_cast<float>(ph)} {}\n
- Direction(double x, double y, double z)\n
- : theta{static_cast<float>(acos(z/std::hypot(x,y,z)))}\n
- , phi{static_cast<float>(atan2(y,x))} {}\n
- template <class VectorType> Direction(const VectorType& v) : Direction(v.theta(), v.phi()) {}\n
- operator std::pair<float, float>() const {return {theta, phi};}\n
- float eta() const {return -log(tan(0.5*theta));}\n
- Direction add(const Direction& rhs) const {return {theta+rhs.theta, phi+rhs.phi};}\n
- Direction subtract(const Direction& rhs) const {return {theta-rhs.theta, phi-rhs.phi};}\n
+ Vector3f() : x(0),y(0),z(0) {}\n
+ Vector3f(float xx, float yy, float zz) : x(xx),y(yy),z(zz) {}\n
+ Vector3f(const float* v) : x(v[0]),y(v[1]),z(v[2]) {}\n
+ bool operator==(const Vector3f& v) const { return (x==v.x&&y==v.y&&z==v.z) ; }\n
+ float operator[](unsigned i) const { return *( &x + i ) ; }\n
+ Vector3f(const edm4hep::Vector3f& v) : x{v.x}, y{v.y}, z{v.z} {}
+ operator edm4hep::Vector3f() const {return {x, y, z};}\n
+ "
+ includes: "
+ #include <edm4hep/Vector3f.h>\n
"
eic::VectorXYZ:
@@ -125,7 +128,8 @@ components:
declaration: "
VectorXYZ() : x{0}, y{0}, z{0} {}\n
VectorXYZ(double xx, double yy, double zz) : x{static_cast<float>(xx)}, y{static_cast<float>(yy)}, z{static_cast<float>(zz)} {}\n
- template<class VectorPolarType> VectorXYZ(const VectorPolarType& v) : VectorXYZ(v.x(), v.y(), v.z()) {}\n
+ template <class V>
+ VectorXYZ(const V& v) : x{v.x}, y{v.y}, z{v.z} {}\n
float& operator[](unsigned i) {return *(&x + i);}\n
const float& operator[](unsigned i) const {return *(&x + i);}\n
float mag() const {return std::hypot(x, y, z);}\n
@@ -138,25 +142,7 @@ components:
VectorXYZ add(const VectorXYZ& rhs) const {return {x+rhs.x, y+rhs.y, z+rhs.z};}\n
VectorXYZ subtract(const VectorXYZ& rhs) const {return {x-rhs.x, y-rhs.y, z-rhs.z};}\n
VectorXYZ scale(double f) const {return {f*x, f*y, f*z};}\n
- "
- eic::VectorPolar:
- Members:
- - float r // [mm] or [GeV]
- - float theta // [rad, 0->pi]
- - float phi // [rad, -pi->pi]
- ExtraCode:
- includes: "#include <cmath>\n#include<tuple>"
- declaration: "
- VectorPolar() : r{0}, theta{0}, phi{0} {}\n
- VectorPolar(double rr, double th, double ph) : r{static_cast<float>(rr)}, theta{static_cast<float>(th)}, phi{static_cast<float>(ph)} {}\n
- template<class VectorXYZType> VectorPolar(const VectorXYZType& v) : VectorPolar(v.r(), v.theta(), v.phi()) {}\n
- float mag() const {return r;}\n
- float x() const {return r * cos(phi) * sin(theta);}\n
- float y() const {return r * sin(phi) * sin(theta);}\n
- float z() const {return r * cos(theta);}\n
- float eta() const {return -log(tan(0.5*theta));}\n
- operator std::tuple<float, float, float>() {return {r, theta, phi};}\n
- "
+ "
eic::VectorXYZT:
Members:
@@ -265,12 +251,6 @@ components:
return *(&xy + i + j - 1);\n
}\n
"
- ## ProtoCluster hit relation
- eic::ProtoClusterHit:
- Members:
- - eic::Index ID // ID of the hit
- - uint32_t index // Raw index of the hit in the relevant array
- - eic::Weight weight // weight of the hit
## A point along a track
eic::TrackPoint:
@@ -314,8 +294,8 @@ datatypes:
Author: "W. Armstrong, S. Joosten"
Members:
- eic::Index ID // Unique particle index
- - eic::VectorXYZ p // Momentum [GeV]
- - eic::VectorXYZ v // Vertex [mm]
+ - eicd::Vector3f p // Momentum [GeV]
+ - eicd::Vector3f v // Vertex [mm]
- float time // Time in [ns]
- int32_t pid // Particle PDG code
- int16_t status // Status code
@@ -334,7 +314,8 @@ datatypes:
- int16_t status // Status code
- int16_t charge // Particle charge (or sign)
- eic::Weight weight // Particle weight, e.g. from PID algorithm [0-1]
- - eic::Direction direction // Direction (theta/phi of this particle [rad])
+ - float theta // Polar angle of this particle [rad]
+ - float phi // Azimuthal angle of this particle [rad]
- float momentum // particle 3-momentum magnitude [GeV]
- float energy // Particle energy, consistent with PID assigment [GeV]
- float mass // The mass of the particle in [GeV]
@@ -353,80 +334,61 @@ datatypes:
Description: "Raw (digitized) calorimeter hit"
Author: "W. Armstrong, S. Joosten"
Members:
- - eic::Index ID // Unique hit ID. For MC, the value equals the Geant4 hit index.
- - int64_t cellID // The detector specific (geometrical) cell id.
- - int64_t amplitude // The amplitude of the hit in ADC counts.
- - int64_t time // Timing in TDC
+ - uint64_t cellID // The detector specific (geometrical) cell id.
+ - uint64_t amplitude // The amplitude of the hit in ADC counts.
+ - uint64_t timeStamp // Timing in TDC
eic::CalorimeterHit:
Description: "Calorimeter hit"
Author: "W. Armstrong, S. Joosten"
Members:
- - eic::Index ID // Unique hit ID, same as one of the involved raw hits.
- - int64_t cellID // The detector specific (geometrical) cell id.
- - int32_t layer // Layer for this hit
- - int32_t sector // Sector for this hit
+ - uint64_t cellID // The detector specific (geometrical) cell id.
- float energy // The energy for this hit in [GeV].
- float energyError // Error on energy [GeV].
- float time // The time of the hit in [ns].
+ - float timeError // Error on the time
- eic::VectorXYZ position // The global position of the hit in world coordinates [mm].
- - eic::VectorXYZ local // The local position of the hit in detector coordinates [mm].
- eic::VectorXYZ dimension // The dimension information of the cell [mm].
- ## TODO: deterimine if position/dimension is really what we want
- ## compared to position/positionError for consistency with TrackerHit
-
+ - int32_t sector // Sector that this hit occured in
+ - int32_t layer // Layer that the hit occured in
+ - eic::VectorXYZ local // The local coordinates of the hit in the detector segment [mm].
## ==========================================================================
## Clustering
## ==========================================================================
eic::ProtoCluster:
- Description: "Relational info linking hits to their associated cluster"
+ Description: "Collection of hits identified by the clustering algorithm to belong together"
Author: "S. Joosten"
- Members:
- - eic::Index ID // ID of the cluster
+ OneToManyRelations:
+ - eic::CalorimeterHit hits // Hits associated with this cluster
VectorMembers:
- - eic::ProtoClusterHit hits // List of hits associated with the cluster
+ - float weights // Weight for each of the hits, mirrors hits array
eic::Cluster:
- Description: "EIC cluster"
+ Description: "EIC hit cluster, reworked to more closely resemble EDM4hep"
Author: "W. Armstrong, S. Joosten, C.Peng"
Members:
- - eic::Index ID // Unique ID for this cluster, value identical to ProtoCluster ID
+ # main variables
+ - int32_t type // Flagword that defines the type of the cluster
- float energy // Reconstructed energy of the cluster [GeV].
- float energyError // Error on the cluster energy [GeV]
- float time // [ns]
+ - float timeError // Error on the cluster time
- uint32_t nhits // Number of hits in the cluster.
- eic::VectorXYZ position // Global position of the cluster [mm].
- eic::CovXYZ positionError // Covariance matrix of the position (6 Parameters).
- - float radius // Shower radius [mm]
- - float skewness // Shower skewness [unitless]
- - eic::VectorPolar polar // Cluster position polar information
- - float eta // Cluster pseudorapidity
- - eic::Direction direction // Intrinsic direction of the cluster propagation [rad, 0->pi, -pi->pi]
- - eic::Index mcID // For MC only - associated MC particle
-
- eic::ClusterLayer:
- Description: "2D Cluster in a single layer for a multi-layer detector"
- Author: "S. Joosten, C. Peng"
- Members:
- - eic::Index ID // Unique layer ID
- - eic::Index clusterID // Associated full 3D cluster, -1 if none
- - int32_t layer // Layer number for this cluster layer
- - uint32_t nhits // Number of hits
- - float energy // Energy in this cluster layer [GeV]
- - float energyError // Error on energy [GeV]
- - float radius // Shower radius [mm]
- - float skewness // Skewness of hits distribution
- - eic::VectorXYZ position // Global center position. [mm]
-
- eic::MergedClusterRelations:
- Description: "Relational info between a merged cluster and its parents"
- Author: "S. Joosten"
- Members:
- - eic::Index clusterID // Associated cluster ID
- - uint32_t size // Number of valid parents
- - std::array<eic::Index, 4> parent // (Up to 4) parents for this cluster
+ - float intrinsicTheta // Intrinsic cluster propagation direction polar angle [rad]
+ - float intrinsicPhi // Intrinsic cluster propagation direction azimuthal angle [rad]
+ - eic::CovXY intrinsicDirectionError // Error on the intrinsic cluster propagation direction
+ VectorMembers:
+ - float shapeParameters // Should be set in metadata, for now radius/skewness
+ - float hitContributions // Energy contributions of the hits. Runs parallel to ::hits()
+ - float subdetectorEnergies // Energies observed in each subdetector used for this cluster.
+ OneToManyRelations:
+ - eic::Cluster clusters // Clusters that have been combined to form this cluster
+ - eic::CalorimeterHit hits // Hits that have been combined to form this cluster
+ #- eic::ParticleID particleIDs // Particle IDs sorted by likelihood, TODO
## ==========================================================================
## RICH/Cherenkov data structures
@@ -437,7 +399,7 @@ datatypes:
Author: "S. Joosten, C. Peng"
Members:
- eic::Index ID // Unique hit ID. For MC, the value equals the Geant4 hit index.
- - int64_t cellID // The detector specific (geometrical) cell id.
+ - uint64_t cellID // The detector specific (geometrical) cell id.
- uint32_t amplitude // PMT signal amplitude [ADC]
- uint32_t time // PMT signal time [TDC]
@@ -446,7 +408,7 @@ datatypes:
Author: "S. Joosten, C. Peng"
Members:
- eic::Index ID // Unique hit ID, same as one of the involved raw hits.
- - int64_t cellID // The detector specific (geometrical) cell id.
+ - uint64_t cellID // The detector specific (geometrical) cell id.
- float npe // Estimated number of photo-electrons [#]
- float time // Time [ns]
- float timeError // Error on the time [ns]
@@ -476,7 +438,7 @@ datatypes:
Author: "W. Armstrong, S. Joosten"
Members:
- eic::Index ID // Unique hit ID. For MC, the value equals the Geant4 hit index.
- - int64_t cellID // The detector specific (geometrical) cell id.
+ - uint64_t cellID // The detector specific (geometrical) cell id.
- int32_t time // TDC value.
- int32_t charge // ADC value
@@ -485,7 +447,7 @@ datatypes:
Author: "W. Armstrong, S. Joosten"
Members:
- eic::Index ID // Unique hit ID, same as one of the involved raw hits.
- - int64_t cellID // The detector specific (geometrical) cell id.
+ - uint64_t cellID // The detector specific (geometrical) cell id.
- eic::VectorXYZ position // Hit (cell) position and time [mm, ns]
- eic::CovDiagXYZ positionError // Covariance Matrix
- float time // Hit time
diff --git a/source/include/eicd/helpers.h b/source/include/eicd/helpers.h
deleted file mode 100644
index 0fb2e9d6d31ed0760c029f6ccc5d4d9469b6ea1d..0000000000000000000000000000000000000000
--- a/source/include/eicd/helpers.h
+++ /dev/null
@@ -1,43 +0,0 @@
-#ifndef EICD_HELPERS_HH
-#define EICD_HELPERS_HH
-
-
-#include <algorithm>
-#include <cmath>
-#include <exception>
-#include <limits>
-#include <string>
-#include <vector>
-
-#include <Math/Vector4D.h>
-
-#include "eicd/TrackParametersCollection.h"
-#include "eicd/ReconstructedParticleCollection.h"
-#include "eicd/ReconstructedParticleData.h"
-
-namespace eicd::helpers {
-
- /** Four momentum from track and mass.
- * Get a vector of 4-momenta from raw tracking info, using an externally
- * provided particle mass assumption.
- */
- inline auto momenta_from_tracking(const std::vector<eic::TrackParametersData>& tracks,
- const double mass)
- {
- std::vector<ROOT::Math::PxPyPzMVector> momenta{tracks.size()};
- // transform our raw tracker info into proper 4-momenta
- std::transform(tracks.begin(), tracks.end(), momenta.begin(), [mass](const auto& track) {
- // make sure we don't divide by zero
- if (fabs(track.qOverP) < 1e-9) {
- return ROOT::Math::PxPyPzMVector{};
- }
- const double p = fabs(1. / track.qOverP);
- const double px = p * cos(track.phi) * sin(track.theta);
- const double py = p * sin(track.phi) * sin(track.theta);
- const double pz = p * cos(track.theta);
- return ROOT::Math::PxPyPzMVector{px, py, pz, mass};
- });
- return momenta;
- }
-} // namespace eicd::helpers
-#endif
diff --git a/source/src/helpers.cxx b/source/src/helpers.cxx
deleted file mode 100644
index d52ef6cd0ba1ad3d121246913a2e65dfbad11e3b..0000000000000000000000000000000000000000
--- a/source/src/helpers.cxx
+++ /dev/null
@@ -1,3 +0,0 @@
-#include "eicd/helpers.h"
-
-
diff --git a/source/CMakeLists.txt b/utils/CMakeLists.txt
similarity index 80%
rename from source/CMakeLists.txt
rename to utils/CMakeLists.txt
index f820ec3874635437379a348f469ff18afaf36464..dba4b55269ca897529d0144b4fd30d8eb64ffee7 100644
--- a/source/CMakeLists.txt
+++ b/utils/CMakeLists.txt
@@ -10,35 +10,36 @@ find_package(ROOT REQUIRED COMPONENTS GenVector MathCore)
# )
#add_custom_target(G__NPDetGeoCad_ROOTDICT DEPENDS ${CMAKE_CURRENT_BINARY_DIR}/G__NPDetGeoCad.cxx)
-add_library(eicd_core SHARED
- src/helpers.cxx
+add_library(eicd_utils SHARED
+ src/utils.cxx
)
-target_compile_features(eicd_core
+target_compile_features(eicd_utils
PUBLIC cxx_auto_type
PUBLIC cxx_trailing_return_types
PUBLIC cxx_std_17
PRIVATE cxx_variadic_templates
)
-target_compile_options(eicd_core PRIVATE
+target_compile_options(eicd_utils PRIVATE
-Wno-extra
-Wno-ignored-qualifiers
-Wno-overloaded-virtual
-Wno-shadow)
-target_link_libraries(eicd_core
+target_link_libraries(eicd_utils
PUBLIC eicd
PUBLIC ROOT::GenVector ROOT::MathCore)
-target_include_directories(eicd_core
+target_include_directories(eicd_utils
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}>
PUBLIC $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>
PUBLIC $<INSTALL_INTERFACE:include>
)
install(FILES
- include/eicd/helpers.h
+ include/eicd/analysis_utils.h
+ include/eicd/vector_utils.h
DESTINATION include/eicd
)
#install(FILES
@@ -47,7 +48,7 @@ install(FILES
# DESTINATION lib)
-install(TARGETS eicd_core
+install(TARGETS eicd_utils
EXPORT ${PROJECT_NAME}Targets
LIBRARY DESTINATION lib
ARCHIVE DESTINATION lib
diff --git a/utils/include/eicd/analysis_utils.h b/utils/include/eicd/analysis_utils.h
new file mode 100644
index 0000000000000000000000000000000000000000..8a55a27cbe2d58c951be5194ebdbaaf83ebfc7f2
--- /dev/null
+++ b/utils/include/eicd/analysis_utils.h
@@ -0,0 +1,43 @@
+#ifndef EICD_UTILS_ANALYSIS_HH
+#define EICD_UTILS_ANALYSIS_HH
+
+#include <algorithm>
+#include <cmath>
+#include <exception>
+#include <limits>
+#include <string>
+#include <vector>
+
+#include <Math/Vector4D.h>
+
+#include <eicd/ReconstructedParticleCollection.h>
+#include <eicd/ReconstructedParticleData.h>
+#include <eicd/TrackParametersCollection.h>
+
+namespace eicd {
+
+/** Four momentum from track and mass.
+ * Get a vector of 4-momenta from raw tracking info, using an externally
+ * provided particle mass assumption.
+ */
+inline auto
+momenta_from_tracking(const std::vector<eic::TrackParametersData>& tracks,
+ const double mass) {
+ std::vector<ROOT::Math::PxPyPzMVector> momenta{tracks.size()};
+ // transform our raw tracker info into proper 4-momenta
+ std::transform(tracks.begin(), tracks.end(), momenta.begin(),
+ [mass](const auto& track) {
+ // make sure we don't divide by zero
+ if (fabs(track.qOverP) < 1e-9) {
+ return ROOT::Math::PxPyPzMVector{};
+ }
+ const double p = fabs(1. / track.qOverP);
+ const double px = p * cos(track.phi) * sin(track.theta);
+ const double py = p * sin(track.phi) * sin(track.theta);
+ const double pz = p * cos(track.theta);
+ return ROOT::Math::PxPyPzMVector{px, py, pz, mass};
+ });
+ return momenta;
+}
+} // namespace eicd
+#endif
diff --git a/utils/include/eicd/vector_utils.h b/utils/include/eicd/vector_utils.h
new file mode 100644
index 0000000000000000000000000000000000000000..1c816282ece9f3b1c51447c2f3aeb87103747fb5
--- /dev/null
+++ b/utils/include/eicd/vector_utils.h
@@ -0,0 +1,133 @@
+#ifndef EICD_UTILS_VECTOR_HH
+#define EICD_UTILS_VECTOR_HH
+
+#include <algorithm>
+#include <cmath>
+#include <exception>
+#include <limits>
+#include <string>
+#include <vector>
+
+#include <Math/Vector4D.h>
+
+#include <edm4hep/Vector3f.h>
+
+#include <eicd/ReconstructedParticleCollection.h>
+#include <eicd/ReconstructedParticleData.h>
+#include <eicd/TrackParametersCollection.h>
+
+namespace eicd {
+
+template <class V> concept VectorHasX = requires(V v) { v.x; };
+template <class V> concept VectorHasY = requires(V v) { v.y; };
+template <class V> concept VectorHasZ = requires(V v) { v.z; };
+template <class V> concept ClassVector = requires(V v) { v.x(); };
+template <class V>
+concept Vector2D =
+ VectorHasX<V>&& VectorHasY<V> && !VectorHasZ<V> && !ClassVector<V>;
+template <class V>
+concept Vector3D =
+ VectorHasX<V>&& VectorHasY<V>&& VectorHasZ<V> && !ClassVector<V>;
+template <class V> concept VectorND = Vector2D<V> || Vector3D<V>;
+
+inline double etaToAngle(const double eta) {
+ return std::atan(std::exp(-eta)) * 2.;
+}
+inline double angleToEta(const double theta) {
+ return -std::log(std::tan(0.5 * theta));
+}
+
+// inline edm4hep::Vector2f VectorFromPolar(const double r, const double theta)
+// {
+// return {r * sin(theta), r * cos(theta)};
+//}
+template <Vector3D V = edm4hep::Vector3f>
+V sphericalToVector(const double r, const double theta, const double phi) {
+ using FloatType = decltype(V().x);
+ const double sth = sin(theta);
+ const double cth = cos(theta);
+ const double sph = sin(phi);
+ const double cph = cos(phi);
+ const FloatType x = r * sth * cph;
+ const FloatType y = r * sth * sph;
+ const FloatType z = r * cth;
+ return {x, y, z};
+}
+
+template <Vector3D V> double anglePolar(const V& v) {
+ return std::atan2(std::hypot(v.x, v.y), v.z);
+}
+template <VectorND V> double angleAzimuthal(const V& v) {
+ return std::atan2(v.y, v.x);
+}
+template <Vector3D V> double eta(const V& v) {
+ return angleToEta(anglePolar(v));
+}
+template <Vector2D V> double magnitude(const V& v) {
+ return std::hypot(v.x, v.y);
+}
+template <Vector3D V> double magnitude(const V& v) {
+ return std::hypot(v.x, v.y, v.z);
+}
+template <Vector3D V> double magnitudeTransverse(const V& v) {
+ return std::hypot(v.x, v.y);
+}
+template <Vector3D V> double magnitudeLongitudinal(const V& v) { return v.z; }
+template <VectorND V> V normalizeVector(const V& v, double norm = 1.) {
+ const double old = magnitude(v);
+ if (old == 0) {
+ return v;
+ }
+ return (norm > 0) ? v * norm / old : v * 0;
+}
+template <Vector3D V> V vectorTransverse(const V& v) { return {v.x, v.y, 0}; }
+template <Vector3D V> V vectorLongitudinal(const V& v) { return {0, 0, v.z}; }
+// Two vector functions
+template <VectorND V> double angleBetween(const V& v1, const V& v2) {
+ const double dot = v1 * v2;
+ if (dot == 0) {
+ return 0.;
+ }
+ return acos(dot / (magnitude(v1) * magnitude(v2)));
+}
+// Project v onto v1
+template <Vector3D V> double projection(const V& v, const V& v1) {
+ const double norm = magnitude(v1);
+ if (norm == 0) {
+ return magnitude(v);
+ }
+ return v * v1 / norm;
+}
+
+} // namespace eicd
+template <eicd::Vector2D V> V operator+(const V& v1, const V& v2) {
+ return {v1.x + v2.x, v1.y + v2.y};
+}
+template <eicd::Vector3D V> V operator+(const V& v1, const V& v2) {
+ return {v1.x + v2.x, v1.y + v2.y, v1.z + v2.z};
+}
+template <eicd::Vector2D V> double operator*(const V& v1, const V& v2) {
+ return v1.x * v2.x + v1.y * v2.y;
+}
+template <eicd::Vector3D V> double operator*(const V& v1, const V& v2) {
+ return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z;
+}
+template <eicd::Vector2D V> V operator*(const double d, const V& v) {
+ return {d * v.x, d * v.y};
+}
+template <eicd::Vector3D V> V operator*(const double d, const V& v) {
+ return {d * v.x, d * v.y, d * v.z};
+}
+template <eicd::Vector2D V> V operator*(const V& v, const double d) {
+ return {d * v.x, d * v.y};
+}
+template <eicd::Vector3D V> V operator*(const V& v, const double d) {
+ return {d * v.x, d * v.y, d * v.z};
+}
+template <eicd::VectorND V> V operator-(const V& v1, const V& v2) {
+ return v1 + (-1. * v2);
+}
+template <eicd::VectorND V> V operator/(const V& v, const double d) {
+ return (1. / d) * v;
+}
+#endif
diff --git a/utils/src/utils.cxx b/utils/src/utils.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..2e1256b5e6315132d0206281b00bda86924bcdea
--- /dev/null
+++ b/utils/src/utils.cxx
@@ -0,0 +1,2 @@
+#include "eicd/analysis_utils.h"
+#include "eicd/vector_utils.h"