Coordinate systems

Introduction

This notebooks provides examples of curvilinear coordinates. Particularly, it presents the coordinate surfaces for three common coordinate systems.

Use of PyVista in this notebook

We use PyVista to represent the coordinate surfaces and ìtkwidgets to render them interactively in the notebook.

If installing things manually, the following would be required

conda install -c conda-forge pyvista
conda install -c conda-forge itkwidgets


More information: https://docs.pyvista.org/examples/index.html


In [1]:
import numpy as np
import pyvista as pv
from itkwidgets import view

In [2]:
red = (0.9, 0.1, 0.1)
blue = (0.2, 0.5, 0.7)
green = (0.3, 0.7, 0.3)

Cylindrical coordinates

The ISO standard 80000-2 recommends the use of $(ρ, φ, z)$, where $ρ$ is the radial coordinate, $\varphi$ the azimuth, and $z$ the height.

For the conversion between cylindrical and Cartesian coordinates, it is convenient to assume that the reference plane of the former is the Cartesian $xy$-plane (with equation $z=0$, and the cylindrical axis is the Cartesian $z$-axis. Then the $z$-coordinate is the same in both systems, and the correspondence between cylindrical $(\rho, \varphi)$ and Cartesian $(x, y)$ are the same as for polar coordinates, namely

\begin{align} x &= \rho \cos \varphi \\ y &= \rho \sin \varphi \end{align}

in one direction, and

\begin{align} \rho &= \sqrt{x^2+y^2} \\ \varphi &= \begin{cases} 0 & \mbox{if } x = 0 \mbox{ and } y = 0\\ \arcsin\left(\frac{y}{\rho}\right) & \mbox{if } x \geq 0 \\ \arctan\left(\frac{y}{x}\right) & \mbox{if } x > 0 \\ -\arcsin\left(\frac{y}{\rho}\right) + \pi & \mbox{if } x < 0 \end{cases} \end{align}

in the other. It is also common to use $\varphi = \mathrm{atan2}(y, x)$.


In [3]:
# Cylinder
phi, z = np.mgrid[0:2*np.pi:31j, -2:2*np.pi:31j]
x = 2*np.cos(phi)
y = 2*np.sin(phi)
cylinder = pv.StructuredGrid(x, y, z)

In [4]:
# Vertical plane
rho, z = np.mgrid[0:3:31j, -2:2*np.pi:31j]
x = rho*np.cos(np.pi/4)
y = rho*np.sin(np.pi/4)
vert_plane = pv.StructuredGrid(x, y, z)

In [5]:
# Horizontal plane
rho, phi = np.mgrid[0:3:31j, 0:2*np.pi:31j]
x = rho*np.cos(phi)
y = rho*np.sin(phi)
z = np.ones_like(x)
hor_plane = pv.StructuredGrid(x, y, z)

In [6]:
view(geometries=[cylinder, vert_plane, hor_plane],
     geometry_colors=[blue, red, green])


Spherical coordinates

The use of $(r, θ, φ)$ to denote radial distance, inclination (or elevation), and azimuth, respectively, is common practice in physics, and is specified by ISO standard 80000-2.

The spherical coordinates of a point can be obtained from its Cartesian coordinate system $(x, y, z)$

\begin{align} r&=\sqrt{x^2 + y^2 + z^2} \\ \theta &= \arccos\frac{z}{\sqrt{x^2 + y^2 + z^2}} = \arccos\frac{z}{r} \\ \varphi &= \arctan \frac{y}{x} \end{align}

The inverse tangent denoted in $\varphi = \arctan\left(\frac{y}{x}\right)$ must be suitably defined , taking into account the correct quadrant of $(x, y)$ (using the function atan2).

Conversely, the Cartesian coordinates may be retrieved from the spherical coordinates (radius $r$, inclination $\theta$, azimuth $\varphi$), where $r \in [0, \infty)$, $\theta \in [0, \pi]$ and $\varphi \in [0, 2\pi)$, by:

\begin{align} x&=r \, \sin\theta \, \cos\varphi \\ y&=r \, \sin\theta \, \sin\varphi \\ z&=r \, \cos\theta \end{align}

In [7]:
theta, phi = np.mgrid[0:np.pi:21j, 0:np.pi:21j]

In [8]:
# Sphere
x = np.sin(phi) * np.cos(theta)
y = np.sin(phi) * np.sin(theta)
z = np.cos(phi)
sphere = pv.StructuredGrid(x, y, z)
sphere2 =  pv.StructuredGrid(-x, -y, z)

In [9]:
# Cone
x = theta/3 * np.cos(phi)
y = theta/3 * np.sin(phi)
z = theta/3
cone = pv.StructuredGrid(x, y, z)
cone2 = pv.StructuredGrid(-x, -y, z)

In [10]:
# Plane
x = theta/np.pi
y = theta/np.pi
z = phi - np.pi/2
plane = pv.StructuredGrid(x, y, z)

In [11]:
view(geometries=[sphere, sphere2, cone, cone2, plane],
     geometry_colors=[blue, blue, red, red, green],
     geometry_opacities=[1.0, 0.5, 1.0, 0.5, 1.0])


Ellipsoidal coordinates


In [34]:
v, theta = np.mgrid[0:np.pi/2:21j, 0:np.pi:21j]
a = 3
b = 2
c = 1

In [30]:
# Ellipsoid
lam = 3
x =  np.sqrt(a**2  + lam) * np.cos(v) * np.cos(theta)
y = np.sqrt(b**2  + lam)* np.cos(v) * np.sin(theta)
z = np.sqrt(c**2 + lam) * np.sin(v)
ellipsoid = pv.StructuredGrid(x, y, z)
ellipsoid2 = pv.StructuredGrid(x, y, -z)
ellipsoid3 = pv.StructuredGrid(x, -y, z)
ellipsoid4 = pv.StructuredGrid(x, -y, -z)

In [31]:
# Hyperboloid of one sheet
mu = 2
x = np.sqrt(a**2 + mu) * np.cosh(v) * np.cos(theta)
y = np.sqrt(b**2 + mu) * np.cosh(v) * np.sin(theta)
z = np.sqrt(c**2 + mu) * np.sinh(v)
z = np.sqrt(c**2 + mu) * np.sinh(v)
hyper = pv.StructuredGrid(x, y, z)
hyper2 = pv.StructuredGrid(x, y, -z)
hyper3 = pv.StructuredGrid(x, -y, z)
hyper4 = pv.StructuredGrid(x, -y, -z)

In [32]:
# Hyperboloid of two sheets
nu = 1
x = np.sqrt(a**2 + nu) * np.cosh(v)
y = np.sqrt(c**2 + nu) * np.sinh(v) * np.sin(theta)
z = np.sqrt(b**2 + nu) * np.sinh(v) * np.cos(theta)
hyper_up = pv.StructuredGrid(x, y, z)
hyper_down = pv.StructuredGrid(-x, y, z)
hyper_up2 = pv.StructuredGrid(x, -y, z)
hyper_down2 = pv.StructuredGrid(-x, -y, z)

In [33]:
view(geometries=[ellipsoid, ellipsoid2, ellipsoid3, ellipsoid4,
                 hyper, hyper2, hyper3, hyper4,
                 hyper_up, hyper_down, hyper_up2, hyper_down2],
     geometry_colors=[blue]*4 + [red]*4 + [green]*4)


References

  1. Wikipedia contributors. "Cylindrical coordinate system." Wikipedia, The Free Encyclopedia. Wikipedia, The Free Encyclopedia, 12 Dec. 2016. Web. 9 Feb. 2017.

  2. Wikipedia contributors. "Spherical coordinate system." Wikipedia, The Free Encyclopedia. Wikipedia, The Free Encyclopedia, 29 Jan. 2017. Web. 9 Feb. 2017.

  3. Sullivan et al., (2019). PyVista: 3D plotting and mesh analysis through a streamlined interface for the Visualization Toolkit (VTK). Journal of Open Source Software, 4(37), 1450, https://doi.org/10.21105/joss.01450


In [17]:
from IPython.core.display import HTML
def css_styling():
    styles = open('./styles/custom_barba.css', 'r').read()
    return HTML(styles)
css_styling()


Out[17]:

In [ ]: