Skip to content
Snippets Groups Projects
Commit f2d83250 authored by Sylvester Joosten's avatar Sylvester Joosten
Browse files

Resolve "Move closer to EDM4hep (part 2)"

parent f184edc8
Branches
Tags
1 merge request!75Resolve "Move closer to EDM4hep (part 2)"
Pipeline #28234 failed
Stage: config
Stage: build
Stage: deploy
Show whitespace changes
Inline Side-by-side
eic_data.yaml
+ 58
− 138
View file @ f2d83250
......@@ -45,39 +45,6 @@ components:
explicit operator bool() const {return valid();}
"
## simple weight that defaults to 1 if not set
eic::Weight:
Members:
- float value
ExtraCode:
declaration: "
Weight() : value{1.} {}\n
Weight(double w) : value {static_cast<float>(w)} {}\n
Weight& operator=(double w) {value = static_cast<float>(w); return *this;}\n
operator float() const {return value;}
"
eic::VectorXY:
Members:
- float x // [mm] or [GeV]
- float y //
ExtraCode:
includes: "#include <cmath>\n"
declaration: "
VectorXY() : x{0}, y{0} {}\n
VectorXY(double xx, double yy) : x{static_cast<float>(xx)}, y{static_cast<float>(yy)} {}\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);}\n
float r() const {return mag();}\n
float phi() const {return atan2(y, x);}\n
operator std::pair<float, float>() const {return {x, y};}\n
float dot(const VectorXY& rhs) const {return x*rhs.x + y*rhs.y;}\n
VectorXY add(const VectorXY& rhs) const {return {x+rhs.x, y+rhs.y};}\n
VectorXY subtract(const VectorXY& rhs) const {return {x-rhs.x, y-rhs.y};}\n
VectorXY scale(double f) const {return {f*x, f*y};}\n
"
# Copy from EDM4hep with conversion operators added
# TODO: uncomment Vector2f conversions when we migrate to new EDM4hep version
eicd::Vector2f :
......@@ -117,94 +84,47 @@ components:
includes: "
#include <edm4hep/Vector3f.h>\n
"
eic::VectorXYZ:
# Copy from EDM4hep with conversion operators added
eicd::Vector3d :
Members:
- float x // [mm] or [GeV]
- float y //
- float z //
- double x
- double y
- double z
ExtraCode:
includes: "#include <cmath>\n#include<tuple>"
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 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
float r() const {return mag();}\n
float theta() const {return acos(z/mag());}\n
float phi() const {return atan2(y,x);}\n
float eta() const {return -log(tan(0.5*theta()));}\n
operator std::tuple<float, float, float>() {return {x, y, z};}\n
float dot(const VectorXYZ& rhs) const {return x*rhs.x + y*rhs.y + z*rhs.z;}\n
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
Vector3d() : x(0),y(0),z(0) {}\n
Vector3d(double xx, double yy, double zz) : x(xx),y(yy),z(zz) {}\n
Vector3d(const double* v) : x(v[0]),y(v[1]),z(v[2]) {}\n
bool operator==(const Vector3d& v) const { return (x==v.x&&y==v.y&&z==v.z) ; }\n
double operator[](unsigned i) const { return *( &x + i ) ; }\n
Vector3d(const edm4hep::Vector3d& v) : x{v.x}, y{v.y}, z{v.z} {}
operator edm4hep::Vector3d() const {return {x, y, z};}\n
"
eic::VectorXYZT:
Members:
- double x // [mm] or [GeV]
- double y //
- double z //
- double t // [ns] or [GeV]
ExtraCode:
includes: "#include <cmath>\n#include <tuple>"
declaration: "
VectorXYZT() : x{0}, y{0}, z{0}, t{0} {}\n
VectorXYZT(double xx, double yy, double zz, double tt) : x{xx}, y{yy}, z{zz}, t{tt} {}\n
double& operator[](unsigned i) {return *(&x + i);}\n
const double& operator[](unsigned i) const {return *(&x + i);}\n
double mag() const {return std::hypot(x, y, z);}\n
double r() const {return mag();}\n
double theta() const {return acos(z/mag());}\n
double phi() const {return atan2(y,x);}\n
double eta() const {return -log(tan(0.5*theta()));}
double energy() const {return t;}\n
double mass() const {return sqrt(t*t - x*x - y*y - z*z);}\n
operator std::tuple<double, double, double, double>() {return {x, y, z, t};}\n
double dot(const VectorXYZT& rhs) const {return t*rhs.t - x*rhs.x - y*rhs.y - z*rhs.z;}\n
VectorXYZT add(const VectorXYZT& rhs) const {return {x+rhs.x, y+rhs.y, z+rhs.z, t+rhs.t};}\n
VectorXYZT subtract(const VectorXYZT& rhs) const {return {x-rhs.x, y-rhs.y, z-rhs.z, t-rhs.t};}\n
VectorXYZT scale(double f) const {return {f*x, f*y, f*z, f*t};}\n
includes: "
#include <edm4hep/Vector3d.h>\n
"
eic::CovDiagXYZ:
eicd::CovDiag3f:
Members:
- float xx
- float yy
- float zz
ExtraCode:
declaration: "
CovDiagXYZ() : xx{0}, yy{0}, zz{0} {}\n
CovDiagXYZ(double x, double y, double z) : xx{static_cast<float>(x)}, yy{static_cast<float>(y)}, zz{static_cast<float>(z)} {}\n
CovDiag3f() : xx{0}, yy{0}, zz{0} {}\n
CovDiag3f(double x, double y, double z) : xx{static_cast<float>(x)}, yy{static_cast<float>(y)}, zz{static_cast<float>(z)} {}\n
float operator()(unsigned i, unsigned j) const {return (i == j) ? *(&xx + i) : 0.;}\n
"
eic::CovDiagXYZT:
Members:
- double xx
- double yy
- double zz
- double tt
ExtraCode:
declaration: "
CovDiagXYZT() : xx{0}, yy{0}, zz{0}, tt{0} {}\n
CovDiagXYZT(double x, double y, double z, double t) : xx{x}, yy{y}, zz{z}, tt{t} {}\n
double operator()(unsigned i, unsigned j) const {return (i == j) ? *(&xx + i) : 0.;}\n
"
eic::CovXY:
eicd::Cov2f:
Members:
- float xx
- float yy
- float xy
ExtraCode:
declaration: "
CovXY() : xx{0}, yy{0}, xy{0} {}\n
CovXY(double vx, double vy, double vxy = 0)\n
Cov2f() : xx{0}, yy{0}, xy{0} {}\n
Cov2f(double vx, double vy, double vxy = 0)\n
: xx{static_cast<float>(vx)}, yy{static_cast<float>(vy)}, xy{static_cast<float>(vxy)} {}\n
float operator()(unsigned i, unsigned j) const {\n
// diagonal\n
......@@ -222,9 +142,9 @@ components:
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
- float weight // weight of the hit
eic::CovXYZ:
eicd::Cov3f:
Members:
- float xx
- float yy
......@@ -234,12 +154,12 @@ components:
- float yz
ExtraCode:
declaration: "
CovXYZ() : xx{0}, yy{0}, zz{0}, xy{0}, xz{0}, yz{0} {}\n
CovXYZ(double vx, double vy, double vz, double vxy, double vxz, double vyz)\n
Cov3f() : xx{0}, yy{0}, zz{0}, xy{0}, xz{0}, yz{0} {}\n
Cov3f(double vx, double vy, double vz, double vxy, double vxz, double vyz)\n
: xx{static_cast<float>(vx)}, yy{static_cast<float>(vy)}, zz{static_cast<float>(vz)},\n
xy{static_cast<float>(vxy)}, xz{static_cast<float>(vxz)}, yz{static_cast<float>(vyz)} {}\n
CovXYZ(double vx, double vy, double vz)\n
: CovXYZ(vx, vy, vz, 0, 0, 0) {} \n
Cov3f(double vx, double vy, double vz)\n
: Cov3f(vx, vy, vz, 0, 0, 0) {} \n
float operator()(unsigned i, unsigned j) const {\n
// diagonal\n
if (i == j) {\n
......@@ -255,15 +175,15 @@ components:
## A point along a track
eic::TrackPoint:
Members:
- eic::VectorXYZ position // Position of the trajectory point [mm]
- eic::CovXYZ positionError // Error on the position
- eic::VectorXYZ momentum // 3-momentum at the point [GeV]
- eic::CovXYZ momentumError // Error on the 3-momentum
- eicd::Vector3f position // Position of the trajectory point [mm]
- eicd::Cov3f positionError // Error on the position
- eicd::Vector3f momentum // 3-momentum at the point [GeV]
- eicd::Cov3f momentumError // Error on the 3-momentum
- float time // Time at this point [ns]
- float timeError // Error on the time at this point
- float theta // polar direction of the track at the surface [rad]
- float phi // azimuthal direction of the track at the surface [rad]
- eic::CovXY directionError // Error on the polar and azimuthal angles
- eicd::Cov2f directionError // Error on the polar and azimuthal angles
- float pathlength // Pathlength from the origin to this point
- float pathlengthError // Error on the pathlenght
......@@ -282,7 +202,7 @@ datatypes:
- int32_t type // Event type identifier (TBD).
- int32_t proc // Process identifier (TBD).
- int32_t source // Source/identifier (TBD)
- eic::Weight weight // Optional event weight (useful for MC)
- float weight // Optional event weight (useful for MC)
## ==========================================================================
......@@ -295,25 +215,25 @@ datatypes:
Members:
- eic::Index ID // Unique particle index
- eicd::Vector3f p // Momentum [GeV]
- eicd::Vector3f v // Vertex [mm]
- eicd::Vector3d v // Vertex [mm]
- float time // Time in [ns]
- int32_t pid // Particle PDG code
- int16_t status // Status code
- int16_t charge // Particle charge (or sign)
- eic::Weight weight // Particle weight, e.g. from PID algorithm [0-1]
- float weight // Particle weight, e.g. from PID algorithm [0-1]
eic::ReconstructedParticle:
Description: "EIC Reconstructed Particle"
Author: "W. Armstrong, S. Joosten"
Members:
- eic::Index ID // Unique particle index
- eic::VectorXYZ p // Momentum vector [GeV]
- eic::VectorXYZ v // Vertex [mm]
- eicd::Vector3f p // Momentum vector [GeV]
- eicd::Vector3d v // Vertex [mm]
- float time // Time in [ns]
- int32_t pid // PID of reconstructed particle.
- int16_t status // Status code
- int16_t charge // Particle charge (or sign)
- eic::Weight weight // Particle weight, e.g. from PID algorithm [0-1]
- float weight // Particle weight, e.g. from PID algorithm [0-1]
- float theta // Polar angle of this particle [rad]
- float phi // Azimuthal angle of this particle [rad]
- float momentum // particle 3-momentum magnitude [GeV]
......@@ -347,11 +267,11 @@ datatypes:
- 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 dimension // The dimension information of the cell [mm].
- eicd::Vector3f position // The global position of the hit in world coordinates [mm].
- eicd::Vector3f dimension // The dimension information of the cell [mm].
- 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].
- eicd::Vector3f local // The local coordinates of the hit in the detector segment [mm].
## ==========================================================================
## Clustering
......@@ -376,11 +296,11 @@ datatypes:
- 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).
- eicd::Vector3f position // Global position of the cluster [mm].
- eicd::Cov3f positionError // Covariance matrix of the position (6 Parameters).
- 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
- eicd::Cov2f 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()
......@@ -412,9 +332,9 @@ datatypes:
- float npe // Estimated number of photo-electrons [#]
- float time // Time [ns]
- float timeError // Error on the time [ns]
- eic::VectorXYZ position // PMT hit position [mm]
- eic::VectorXYZ local // The local position of the hit in detector coordinates [mm]
- eic::VectorXYZ dimension // The dimension information of the pixel [mm].
- eicd::Vector3f position // PMT hit position [mm]
- eicd::Vector3f local // The local position of the hit in detector coordinates [mm]
- eicd::Vector3f dimension // The dimension information of the pixel [mm].
eic::RingImage:
Description: "EIC Ring Image Cluster"
......@@ -422,8 +342,8 @@ datatypes:
Members:
- eic::Index ID // Unique cluster ID
- float npe // Number of photo-electrons [#]
- eic::VectorXYZ position // Global position of the cluster [mm]
- eic::VectorXYZ positionError // Error on the position
- eicd::Vector3f position // Global position of the cluster [mm]
- eicd::Vector3f positionError // Error on the position
- float theta // Opening angle of the ring [rad, 0->pi]
- float thetaError // Error on the opening angle
- float radius // Radius of the best fit ring [mm]
......@@ -448,8 +368,8 @@ datatypes:
Members:
- eic::Index ID // Unique hit ID, same as one of the involved raw hits.
- uint64_t cellID // The detector specific (geometrical) cell id.
- eic::VectorXYZ position // Hit (cell) position and time [mm, ns]
- eic::CovDiagXYZ positionError // Covariance Matrix
- eicd::Vector3f position // Hit (cell) position and time [mm, ns]
- eicd::CovDiag3f positionError // Covariance Matrix
- float time // Hit time
- float timeError // Error on the time
- float edep // Energy deposit in this hit [GeV]
......@@ -464,7 +384,7 @@ datatypes:
## Members:
## eic::Index ID // Unique identifier
## eic::Index seedID // Index of corresponding initial track parameters
## eic::Weight weight // prototrack weight, in case we share pixels [0-1]
## float weight // prototrack weight, in case we share pixels [0-1]
## VectorMembers:
## int32_t hits // tracker hit indicies
......@@ -492,12 +412,12 @@ datatypes:
Author: "W. Armstrong, S. Joosten"
Members:
- int32_t type // Type of track parameters (-1/seed, 0/head, ...)
- eic::VectorXY loc // 2D location on surface
- eic::CovXY locError // Covariance on loc
- eicd::Vector2f loc // 2D location on surface
- eicd::Cov2f locError // Covariance on loc
- float theta // Track polar angle [rad]
- float phi // Track azimuthal angle [rad]
- float qOverP // [e/GeV]
- eic::CovXYZ momentumError // Covariance on theta, phi and qOverP
- eicd::Cov3f momentumError // Covariance on theta, phi and qOverP
- float time // Track time [ns]
- float timeError // Error on the time
- float charge // Particle charge
......@@ -511,8 +431,8 @@ datatypes:
- int32_t type // Flag that defines the type of track
- float chi2 // Total chi2 (sum) of the track fit
- int32_t ndf // Numbers of degrees of freedom of the track fit
- eic::VectorXYZ momentum // Track 3-momentum at the vertex [GeV]
- eic::CovXYZ momentumError // Covariance matrix on the momentum
- eicd::Vector3f momentum // Track 3-momentum at the vertex [GeV]
- eicd::Cov3f momentumError // Covariance matrix on the momentum
- float time // Track time at the vertex [ns]
- float timeError // Error on the track vertex time
- float charge // Particle charge
......@@ -541,7 +461,7 @@ datatypes:
Author: "W. Armstrong, S. Joosten"
Members:
- eic::Index ID // Unique vertex ID
- eic::VectorXYZ position // Postion of vertex [mm]
- eicd::Vector3f position // Postion of vertex [mm]
- float time // Time of vertex [ns]
- float chi2 // Chi squared of the vertex fit.
- float probability // Probability of the vertex fit
......
utils/include/eicd/vector_utils.h
+ 40
− 18
View file @ f2d83250
......@@ -21,14 +21,21 @@ 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 VectorHasA = requires(V v) { v.a; };
template <class V> concept VectorHasB = requires(V v) { v.b; };
template <class V> concept ClassVector = requires(V v) { v.x(); };
template <class V>
concept Vector2D =
concept Vector2D_XY =
VectorHasX<V>&& VectorHasY<V> && !VectorHasZ<V> && !ClassVector<V>;
template <class V>
concept Vector2D_AB =
VectorHasA<V>&& VectorHasB<V> && !VectorHasZ<V> && !ClassVector<V>;
template <class V> concept Vector2D = Vector2D_XY<V> || Vector2D_AB<V>;
template <class V>
concept Vector3D =
VectorHasX<V>&& VectorHasY<V>&& VectorHasZ<V> && !ClassVector<V>;
template <class V> concept VectorND = Vector2D<V> || Vector3D<V>;
template <class V> concept VectorND_XYZ = Vector2D_XY<V> || Vector3D<V>;
inline double etaToAngle(const double eta) {
return std::atan(std::exp(-eta)) * 2.;
......@@ -37,13 +44,20 @@ inline double angleToEta(const double theta) {
return -std::log(std::tan(0.5 * theta));
}
// Utility getters to accomodate different vector types
template <VectorND_XYZ V> auto vector_x(const V& v) { return v.x; }
template <VectorND_XYZ V> auto vector_y(const V& v) { return v.y; }
template <Vector3D V> auto vector_z(const V& v) { return v.z; }
template <Vector2D_AB V> auto vector_x(const V& v) { return v.a; }
template <Vector2D_AB V> auto vector_y(const V& v) { return v.b; }
// 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);
using FloatType = decltype(vector_x(V()));
const double sth = sin(theta);
const double cth = cos(theta);
const double sph = sin(phi);
......@@ -55,24 +69,26 @@ V sphericalToVector(const double r, const double theta, const double phi) {
}
template <Vector3D V> double anglePolar(const V& v) {
return std::atan2(std::hypot(v.x, v.y), v.z);
return std::atan2(std::hypot(vector_x(v), vector_y(v)), vector_z(v));
}
template <VectorND V> double angleAzimuthal(const V& v) {
return std::atan2(v.y, v.x);
return std::atan2(vector_y(v), vector_x(v));
}
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);
return std::hypot(vector_x(v), vector_y(v));
}
template <Vector3D V> double magnitude(const V& v) {
return std::hypot(v.x, v.y, v.z);
return std::hypot(vector_x(v), vector_y(v), vector_z(v));
}
template <Vector3D V> double magnitudeTransverse(const V& v) {
return std::hypot(v.x, v.y);
return std::hypot(vector_x(v), vector_y(v));
}
template <Vector3D V> double magnitudeLongitudinal(const V& v) {
return vector_z(v);
}
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) {
......@@ -80,8 +96,12 @@ template <VectorND V> V normalizeVector(const V& v, double norm = 1.) {
}
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}; }
template <Vector3D V> V vectorTransverse(const V& v) {
return {vector_x(v), vector_y(v), 0};
}
template <Vector3D V> V vectorLongitudinal(const V& v) {
return {0, 0, vector_z(v)};
}
// Two vector functions
template <VectorND V> double angleBetween(const V& v1, const V& v2) {
const double dot = v1 * v2;
......@@ -101,28 +121,30 @@ template <Vector3D V> double projection(const V& v, const V& v1) {
} // namespace eicd
template <eicd::Vector2D V> V operator+(const V& v1, const V& v2) {
return {v1.x + v2.x, v1.y + v2.y};
return {vector_x(v1) + vector_x(v2), vector_y(v1) + vector_y(v2)};
}
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};
return {vector_x(v1) + vector_x(v2), vector_y(v1) + vector_y(v2),
vector_z(v1) + vector_z(v2)};
}
template <eicd::Vector2D V> double operator*(const V& v1, const V& v2) {
return v1.x * v2.x + v1.y * v2.y;
return vector_x(v1) * vector_x(v2) + vector_y(v1) * vector_y(v2);
}
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;
return vector_x(v1) * vector_x(v2) + vector_y(v1) * vector_y(v2) +
vector_z(v1) * vector_z(v2);
}
template <eicd::Vector2D V> V operator*(const double d, const V& v) {
return {d * v.x, d * v.y};
return {d * vector_x(v), d * vector_y(v)};
}
template <eicd::Vector3D V> V operator*(const double d, const V& v) {
return {d * v.x, d * v.y, d * v.z};
return {d * vector_x(v), d * vector_y(v), d * vector_z(v)};
}
template <eicd::Vector2D V> V operator*(const V& v, const double d) {
return {d * v.x, d * v.y};
return {d * vector_x(v), d * vector_y(v)};
}
template <eicd::Vector3D V> V operator*(const V& v, const double d) {
return {d * v.x, d * v.y, d * v.z};
return {d * vector_x(v), d * vector_y(v), d * vector_z(v)};
}
template <eicd::VectorND V> V operator-(const V& v1, const V& v2) {
return v1 + (-1. * v2);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment