Cylindrical transformations with well-defined phi angle

[jngrad]

Follow-up to #4014

Description of changes:

• add overloads to convert from Cartesian to cylindrical coordinates with custom longitudinal axis and phi angle

[RudolfWeeber]

I have a hard time understanding the approach with the phi_0.

I'd hae preferred the approach with a basis transform, analogous to Python code below, since, at least for me, it is clearer how and that that works.

If you'd like to stay with the current implementation, please add an explanation in the comments explaining the approach, particularly the role of phi_0.

Could you please also consider the tests from test_to_cylindrical_basis_change from the below Python code. They are constructed such that the expected results can be inferred easily by a human.

[RudolfWeeber]

import numpy as np

def basis_change(b1,b2,b3, v):
    M = np.linalg.inv(np.array((b1/np.linalg.norm(b1),b2/np.linalg.norm(b2),b3/np.linalg.norm(b3))).T)
    return np.dot(M, v)

def to_cylindrical(x):
    return (np.sqrt(x[0]**2 +x[1]**2), np.arctan2(x[1], x[0]), x[2])


def to_cylindrical_basis_change(x, axis=None, orientation=None):
    x_t=basis_change(orientation, np.cross(axis,orientation), axis, x)
    return to_cylindrical(x_t)

def test_to_cylindrical():
    t = to_cylindrical(np.array((np.sqrt(4.5),np.sqrt(4.5),1.2)))
    np.testing.assert_allclose(t, [3,np.pi/4, 1.2])

b_x = np.array((1,0,0))
b_y = np.array((0,1,0))
b_z = np.array((0,0,1))

def test_basis_transform():
    v = [1,2,3]
   
    # identity transform
    np.testing.assert_allclose(basis_change(b_x,b_y,b_z,v), v)

    # swap coordinates
    np.testing.assert_allclose(basis_change(b_z,b_y,b_x,[v[2],v[1],v[0]]), v)
    
    v1 = np.array(([2,2,2])) /np.sqrt(12)
    v2 = np.array(([3,3,-6])) /np.sqrt(9+9+36)
    c = np.cross(v1,v2)
    v3 = c /np.linalg.norm(c)

    np.testing.assert_allclose(
        basis_change(v1,v2,v3, 0.1*v1+0.2*v2 -0.3*v3),
        [0.1,0.2,-0.3],atol=1E-10)


def test_to_cylindrical_basis_change():
    v1 = np.array((1,3,4))
    v2 = np.array((-3.0,7,2))
    axis = np.cross(v1,v2)
    axis /= np.linalg.norm(axis)


    np.testing.assert_allclose(
        to_cylindrical_basis_change(v1,axis=axis,orientation=v1),
        [np.linalg.norm(v1),0,0],atol=1E-10)
    
    angle_v1_v2 = np.arccos(np.dot(v1,v2)/np.linalg.norm(v1)/np.linalg.norm(v2))
    np.testing.assert_allclose(
        to_cylindrical_basis_change(v2+2*axis,axis=axis,orientation=v1),
        [np.linalg.norm(v2),angle_v1_v2,2],atol=1E-10)
    np.testing.assert_allclose(
        to_cylindrical_basis_change(v1+2*axis,axis=axis,orientation=v2),
        [np.linalg.norm(v1),-angle_v1_v2,2],atol=1E-10)

[jngrad]

I've added the right-handedness test case and updated the documentation of the phi angle.

In the python code sample, I could not reproduce the following snippet in the C++ unit test:

np.testing.assert_allclose(
        to_cylindrical_basis_change(v2+2*axis,axis=axis,orientation=v1),
        [np.linalg.norm(v2),angle_v1_v2,2],atol=1E-10)

When adding 2*axis to the input vector v2, I get random values for z instead of 2.

[KaiSzuttor]

can this be closed because of #4114?