Restrict kinematics functions to position arrays in cartesian coordinates

[sfmig]

Is your feature request related to a problem? Please describe.
Right now we only check if the position array has a time dimension. So in theory we can pass position arrays in polar coordinates to compute_displacement, compute_velocity, compute_acceleration.

This could be confusing for users because the result of passing a position array in polar coordinates to compute_displacement is not the same as the displacement vector in polar coordinates. The same applies to compute_velocity, compute_acceleration and their corresponding vectors.

Additionally, the quantities computed by passing a position array in polar coordinates to the current implementation of the kinematics functions would rarely be useful as a vector. E.g. for compute_displacement, the computed displacement array would have the time derivative of rho along the first spatial dimension and the time derivative of phi along the second spatial dimension. I think it would be very rare if a user wants these two quantities paired up as a vector in their analysis.

If a user explicitly wants these quantities (e.g., the time derivative of phi may be useful by itself), it would be more transparent to directly take the time derivative of that quantity (using xarray's differentiate) rather than to use the displacement wrapper.

Describe the solution you'd like
In the kinematics functions, validate that the input position data is in cartesian coordinates.

Describe alternatives you've considered
\

Additional context
\

[lochhh]

I have a rather different opinion regarding this. I don't see a need to restrict what these functions can do - we expect users to know what's in their data array and the current docstrings already clarify "The input data array containing position vectors in cartesian coordinates, with time as a dimension." ← perhaps we should even remove the mention of "cartesian" here 🤔

Having these kinds of restrictions will mean users cannot use these functions when they have 1D or 3D data, even if their data is in cartesian space.

In terms of position in polar coordinates, keeping the current functions as is allows users to compute radial/angular velocity/acceleration; if they wanted, they could also input "rho" or "phi" separately (e.g. pos.sel(space_pol="phi")) into these functions to say, compute tangential velocity.

[niksirbi]

I see the merits of both arguments for and against restricting the kinematic functions to Cartesian space input. I have a kind of compromise to propose:

• Restrict compute_velocity , compute_displacement and compute_acceleration to cartesian space as per @sfmig 's suggestion. The main reason for me is although the derivative of the phi angle is useful, taking the derivative of rho is a bit non-sensical and confusing.
• Expose the current _compute_approximate_time_derivative function as a public function named compute_time_derivative. This function would only check for the existence of the "time" coordinate. I removed "approximate" from the name to make it less wordy and because as long as we are working with discrete computers, the derivative will always be approximate (and the docstring anyway specifies the nature of that approcimation). The idea would be that users could use this for any custon derivation, e.g. compute_time_derivative(pos.sel(space_pol="phi")). In general this function would just broadcast over all dimensions other than time, so it would allow maximum freedom.

Would you both be ok with such a solution?

[sfmig]

@niksirbi I was gonna suggest that same thing so yeah, I like this option 😁

[sfmig]

I was writing a longer explanation to this, I paste it below mostly for future reference.

It includes the expressions for the radial and tangential velocity, and the radial and tangential acceleration, as a function of the rho and theta coordinates of the position vector.

I agree users should be responsible for the operations they apply to their data, but I think in this case the names of the functions could be genuinely misleading for an input vector in polar coordinates. Below I explain why.

Just to clarify, when I say using a cartesian coordinate system in movement, I mean expressing any vectors (position, displacement, velocity, acceleration,...) as x-y coordinates in 2D, or as x-y-z coordinates in 3D.

Case for displacement

• for an input array of position vectors in time in polar coordinates (rho, theta), applying the compute_displacement function would return an array (\delta rho, \delta theta), with \delta representing the difference between consecutive frames.

• this output is NOT:

  • the displacement vector in polar coordinates. Instead, the output (for a given frame) is a vector, with its two components (in 2D) being the difference in rho and the difference in theta between consecutive frames.

  • Example 1: it is easy to see how this output is Not the displacement vector in polar coordinates in the case of linear motion (see diagram below). If we assume r_1 = 2*r_2, the displacement vector between t=0 and t=1 should be the vector difference between those two positions (in green, offset for visibility). Its rho coordinate will be r_2 and its theta coordinate theta_0 = theta_1 = theta. But if we pass these two vectors in polar coordinates to our compute_displacement function we would get(r_2, 0), which is just the difference in their polar coordinates.

    • 

  • Example 2: another example would be a particle moving in circular motion around the origin. The position vectors between consecutive frames will have constant rho in polar coordinates, so the computed "displacement vector" will have rho=0. But the actual displacement vector (the vector from the previous to the current position) does have some length, and therefore a rho!=0

• I think many users could assume (based on the name of the function) that under the hood we are doing the appropriate transformation to return the displacement vector in polar coordinates.

Case for velocity

• for an input array of position vectors in time in polar coordinates (rho, theta), applying the compute_velocity function would return an array (rho_dot, theta_dot), with dot representing the time-derivative of the variable (either rho or theta).

• This output is NOT:

  • the velocity vector expressed in polar coordinates (meaning, the rho and theta coordinates of the velocity vector), nor

  • the velocity vector expressed as its radial and tangential components. Only the first component of the output would be the radial velocity; see the equation for this below from proofwiki (our rho_dot is their dr/dt, and our theta_dot is their d(theta)/dt, vectors are represented in bold):

  • 

• Like with the displacement case above, I think users may assume that we are doing the appropriate transformation to return the velocity vector in polar coordinates. Confusingly enough, the decomposition of the velocity into radial and tangential components shown above is sometimes called "velocity in polar coordinates".

• Additionally, I think having the pair of values (rho_dot, theta_dot) together in an array would rarely be useful. It may be for some people, but in that case it would be more transparent to compute them using a compute_time_derivative function.

Case for acceleration

• A very similar situation to the velocity would happen with the acceleration, except the output of compute_acceleration would return us rho_dot_dot and theta_dot_dot (the second time-derivative of rho and theta respectively).

• This output is NOT:

  • the acceleration vector expressed in polar coordinates, nor

  • the acceleration vector expressed as its radial and tangential components. See the equation below for this from proofwiki (again our dot is their d()/dt):

  • 

• Like with the velocity, users may expect the result in the equation above when passing a position vector in polar coordinates to compute_acceleration (again, the decomposition above is sometimes called "acceleration in polar coordinates").

Because of these reasons I think we should:

• keep a reference to cartesian in the docstring of these functions, and restrict them to cartesian (2D or 3D) only
• expose the compute_time_derivative, so that users can compute it if they want to (of any quantity)

I realise I should have given more background to this issue, specifically should have also cross-linked to this other closely related issue.

[niksirbi]

Thanks for this detailed explanation Sofia, it will be very useful for when/if we implement the tangential/radial velocity.
I also agree that some users may wrongly assume we are doing something clever, when taking derivatives in polar coordinates (which we are so far not doing).

[lochhh]

Thanks both for the detailed explanations! I see now how this may be confusing.