Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Compare vector results with ROOT's. #36

Draft
wants to merge 10 commits into
base: main
Choose a base branch
from
Draft
10 changes: 5 additions & 5 deletions .github/workflows/ci.yml
Original file line number Diff line number Diff line change
Expand Up @@ -75,11 +75,11 @@ jobs:
# steps:
# - uses: actions/checkout@v2

# - name: Get Conda
# uses: conda-incubator/setup-miniconda@v2
# with:
# environment-file: environment.yml
# activate-environment: vector
- name: Get Conda
uses: conda-incubator/setup-miniconda@v2
with:
environment-file: environment.yml
activate-environment: vector

# - name: Run tests
# shell: "bash -l {0}"
Expand Down
13 changes: 8 additions & 5 deletions src/vector/_compute/spatial/dot.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@

# specialized
def xy_z_xy_z(lib, x1, y1, z1, x2, y2, z2):
return x1 * x2 + y1 * y2 + z1 * z2
return lib.nan_to_num(x1 * x2 + y1 * y2 + z1 * z2, nan=0.0)


def xy_z_xy_theta(lib, x1, y1, z1, x2, y2, theta2):
Expand Down Expand Up @@ -277,7 +277,7 @@ def rhophi_z_xy_eta(lib, rho1, phi1, z1, x2, y2, eta2):

# specialized
def rhophi_z_rhophi_z(lib, rho1, phi1, z1, rho2, phi2, z2):
return rho1 * rho2 * lib.cos(phi1 - phi2) + z1 * z2
return lib.nan_to_num(rho1 * rho2 * lib.cos(phi1 - phi2) + z1 * z2, nan=0.0)


def rhophi_z_rhophi_theta(lib, rho1, phi1, z1, rho2, phi2, theta2):
Expand Down Expand Up @@ -336,8 +336,9 @@ def rhophi_theta_rhophi_z(lib, rho1, phi1, theta1, rho2, phi2, z2):

# specialized
def rhophi_theta_rhophi_theta(lib, rho1, phi1, theta1, rho2, phi2, theta2):
return (
rho1 * rho2 * (lib.cos(phi1 - phi2) + 1 / (lib.tan(theta1) * lib.tan(theta2)))
return lib.nan_to_num(
rho1 * rho2 * (lib.cos(phi1 - phi2) + 1 / (lib.tan(theta1) * lib.tan(theta2))),
nan=0,
)


Expand Down Expand Up @@ -407,7 +408,9 @@ def rhophi_eta_rhophi_eta(lib, rho1, phi1, eta1, rho2, phi2, eta2):
expmeta2 = lib.exp(-eta2)
invtantheta1 = 0.5 * (1 - expmeta1 ** 2) / expmeta1
invtantheta2 = 0.5 * (1 - expmeta2 ** 2) / expmeta2
return rho1 * rho2 * (lib.cos(phi1 - phi2) + invtantheta1 * invtantheta2)
return lib.nan_to_num(
rho1 * rho2 * (lib.cos(phi1 - phi2) + invtantheta1 * invtantheta2), nan=0.0
)


dispatch_map = {
Expand Down
2 changes: 1 addition & 1 deletion src/vector/_compute/spatial/eta.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,7 @@ def rhophi_z(lib, rho, phi, z):


def rhophi_theta(lib, rho, phi, theta):
return -lib.log(lib.tan(0.5 * theta))
return lib.nan_to_num(-lib.log(lib.tan(0.5 * theta)), nan=0.0)


def rhophi_eta(lib, rho, phi, eta):
Expand Down
81 changes: 81 additions & 0 deletions tests/root/test_EulerAngles.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,81 @@
# Copyright (c) 2019-2021, Jonas Eschle, Jim Pivarski, Eduardo Rodrigues, and Henry Schreiner.
#
# Distributed under the 3-clause BSD license, see accompanying file LICENSE
# or https://github.com/scikit-hep/vector for details.

import pytest

# If ROOT is not available, skip these tests.
ROOT = pytest.importorskip("ROOT")

# 4D constructor arguments to get all the weird cases.
constructor = [
(0, 0, 0, 0),
(0, 0, 1, 0), # theta == 0.0
(0, 0, -1, 0),
(0, 0, 1, 0),
(0, 0, 0, 4294967296),
(0, 4294967296, 0, 0),
(0, 0, 0, 10),
(0, 0, 0, -10),
(1, 2, 3, 0),
(1, 2, 3, 10),
(1, 2, 3, -10),
(1.0, 2.0, 3.0, 2.5),
(1, 2, 3, 2.5),
(1, 2, 3, -2.5),
]

# Coordinate conversion methods to apply to the VectorObject4D.
coordinate_list = [
"to_xyzt",
"to_xythetat", # may fail for constructor2
"to_xyetat",
"to_rhophizt",
"to_rhophithetat",
"to_rhophietat",
"to_xyztau",
"to_xythetatau",
"to_xyetatau",
"to_rhophiztau",
"to_rhophithetatau",
"to_rhophietatau",
]


@pytest.fixture(scope="module", params=coordinate_list)
def coordinates(request):
return request.param


angle_list = [
0,
0.0,
0.7853981633974483,
-0.7853981633974483,
1.5707963267948966,
-1.5707963267948966,
3.141592653589793,
-3.141592653589793,
6.283185307179586,
-6.283185307179586,
]


@pytest.fixture(scope="module", params=angle_list)
def angle(request):
return request.param


scalar_list = [
0,
-1,
1.0,
100000.0000,
-100000.0000,
]


@pytest.fixture(scope="module", params=scalar_list)
def scalar(request):
return request.param
246 changes: 246 additions & 0 deletions tests/root/test_Polar2DVector.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,246 @@
# Copyright (c) 2019-2021, Jonas Eschle, Jim Pivarski, Eduardo Rodrigues, and Henry Schreiner.
#
# Distributed under the 3-clause BSD license, see accompanying file LICENSE
# or https://github.com/scikit-hep/vector for details.

import numpy as np
import pytest

import vector

# If ROOT is not available, skip these tests.
ROOT = pytest.importorskip("ROOT")

# ROOT.Math.Polar2DVector constructor arguments to get all the weird cases.
# r > 0 and phi from 0 to 360 deg?
constructor = [
(0, 0),
(0, 10),
(1, 0),
(10, 0),
(1, 10),
(10.0, 10),
(1.0, 2.5),
(1, 2.5),
(1, 6.283185307179586),
]

# Coordinate conversion methods to apply to the VectorObject2D.
coordinate_list = [
"to_xy",
"to_rhophi",
]


@pytest.fixture(scope="module", params=coordinate_list)
def coordinates(request):
return request.param


angle_list = [
0,
0.0,
0.25 * np.pi,
-0.25 * np.pi,
0.5 * np.pi,
-0.5 * np.pi,
np.pi,
-np.pi,
2 * np.pi,
-2 * np.pi,
]


@pytest.fixture(scope="module", params=angle_list)
def angle(request):
return request.param


scalar_list = [
0,
-1,
1.0,
100000.0000,
-100000.0000,
]


@pytest.fixture(scope="module", params=scalar_list)
def scalar(request):
return request.param


# Run a test that compares ROOT's 'Dot()' with vector's 'dot' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Dot(constructor, coordinates):
assert ROOT.Math.Polar2DVector(*constructor).Dot(
ROOT.Math.Polar2DVector(*constructor)
) == pytest.approx(
getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().dot(
getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
),
1.0e-6,
1.0e-6,
)


# Run a test that compares ROOT's 'Mag2()' with vector's 'rho2' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Mag2(constructor, coordinates):
assert ROOT.Math.Polar2DVector(*constructor).Mag2() == pytest.approx(
getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().rho2
)


# Run a test that compares ROOT's 'R()' with vector's 'rho' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Mag(constructor, coordinates):
assert ROOT.Math.sqrt(
ROOT.Math.Polar2DVector(*constructor).Mag2()
) == pytest.approx(
getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)().rho
)


# Run a test that compares ROOT's 'Phi()' with vector's 'rho' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Phi(constructor, coordinates):
assert ROOT.Math.Polar2DVector(*constructor).Phi() == pytest.approx(
getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)().phi
)


# Run a test that compares ROOT's 'Rotate()' with vector's 'rotateZ' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Rotate(constructor, angle, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor)
ref_vec.Rotate(angle)
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().rotateZ(angle)
res_vec = vec.rotateZ(angle)
assert ref_vec.R() == pytest.approx(
res_vec.rho,
1.0e-6,
1.0e-6,
)
assert ref_vec.Phi() == pytest.approx(
res_vec.phi,
1.0e-6,
1.0e-6,
)


# Run a test that compares ROOT's 'Unit()' with vector's 'unit' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_Unit(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).Unit()
vec = getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
res_vec = vec.unit
assert ref_vec.R() == pytest.approx(res_vec().rho)
assert ref_vec.Phi() == pytest.approx(res_vec().phi)


# Run a test that compares ROOT's 'X()' and 'Y()' with vector's 'x' and 'y' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_X_and_Y(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor)
vec = getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
assert ref_vec.X() == pytest.approx(vec.x) and ref_vec.Y() == pytest.approx(vec.y)


# Run a test that compares ROOT's '__add__' with vector's 'add' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_add(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).__add__(
ROOT.Math.Polar2DVector(*constructor)
)
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().add(
getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
)
assert ref_vec.R() == pytest.approx(
vec.rho,
1.0e-6,
1.0e-6,
)
assert ref_vec.Phi() == pytest.approx(
vec.phi,
1.0e-6,
1.0e-6,
)


# Run a test that compares ROOT's '__sub__' with vector's 'subtract' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_sub(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).__sub__(
ROOT.Math.Polar2DVector(*constructor)
)
vec1 = getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
vec2 = getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
res_vec = vec1.subtract(vec2)
assert ref_vec.R() == pytest.approx(
res_vec.rho,
1.0e-6,
1.0e-6,
)
assert ref_vec.Phi() == pytest.approx(
res_vec.phi,
1.0e-6,
1.0e-6,
)


# Run a test that compares ROOT's '__neg__' with vector's '__neg__' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_neg(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).__neg__()
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().__neg__
assert ref_vec.R() == pytest.approx(vec().rho)
assert ref_vec.Phi() == pytest.approx(vec().phi)


# Run a test that compares ROOT's '__mul__' with vector's 'mul' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_mul(constructor, scalar, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).__mul__(scalar)
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().__mul__(scalar)
assert ref_vec.R() == pytest.approx(vec.rho)
assert ref_vec.Phi() == pytest.approx(vec.phi)


# Run a test that compares ROOT's '__truediv__' with vector's '__truediv__' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_truediv(constructor, scalar, coordinates):
# FIXME:
if scalar != 0:
ref_vec = ROOT.Math.Polar2DVector(*constructor).__truediv__(scalar)
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().__truediv__(scalar)
assert ref_vec.R() == pytest.approx(vec.rho)
assert ref_vec.Phi() == pytest.approx(vec.phi)


# Run a test that compares ROOT's '__eq__' with vector's 'isclose' for all cases.
@pytest.mark.parametrize("constructor", constructor)
def test_eq(constructor, coordinates):
ref_vec = ROOT.Math.Polar2DVector(*constructor).__eq__(
ROOT.Math.Polar2DVector(*constructor)
)
vec = getattr(
vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates
)().isclose(
getattr(vector.obj(**dict(zip(["rho", "phi"], constructor))), coordinates)()
)
assert ref_vec == vec
Loading