epsproc.sphPlot module

ePSproc spherical polar plotting functions

Collection of functions for plotting 2D spherical polar data \(I(\theta, \phi)\).

Main plotters are:

28/11/19
  • Changed plotTypeSelector to dictionary method (and added more methods).

13/09/19

  • Added sphFromBLMPlot().
  • Some additional plot type fuctionality added.
  • Some tidying up and reorganising… hopefully nothing broken…

14/08/19 v1

epsproc.sphPlot.plotTypeSelector(dataPlot, pType='a', axisUW='Eke', returnDict=False)[source]

Set plotting data type.

Parameters:
  • dataPlot (np.array, Xarray) – Data for plotting, output converted to type defined by pType.
  • pType (char, optional, default 'a') –

    Set type of plot.

    • ’a’ (abs) = np.abs(dataPlot)
    • ’a2’ (abs^2) = np.abs(dataPlot**2)
    • ’r’ (real) = np.real(dataPlot)
    • ’i’ (imag) = np.imag(dataPlot)
    • ’p’ (product) = dataPlot * np.conj(dataPlot)
    • ’phase’ = np.angle(dataPlot)
    • ’phaseUW’ = np.unwrap(np.angle(dataPlot), axis = axisUW)
  • axisUW (str, optional, default 'Eke') – Axis to use for phase unwrapping (for pType = ‘phaseUW’). Axis name must be in passed Xarray.
  • returnDict (optional, default = False) – If true, return dictionary of types & methods instead of data array.
Returns:

  • dataPlot (Xarray or np.array) – Input data structure converted to pType.
  • pTypeDict (dictionary) – Structure of valid types, formula and functions.

epsproc.sphPlot.sphFromBLMPlot(BLMXin, res=50, pType='r', plotFlag=False, facetDim=None, backend='mpl', fnType=None)[source]

Calculate spherical harmonic expansions from BLM parameters and plot.

Surfaces calculated as:

\[I(\theta,\phi)=\sum_{L,M}\beta_{L,M}Y_{L,M}(\theta,\phi)\]
Parameters:
  • dataIn (Xarray) – Input set of BLM parameters, or other (L,M) exapansion dataType.
  • res (int, optional, default 50) – Resolution for output (theta,phi) grids.
  • pType (char, optional, default 'r' (real part)) – Set (data) type to plot. See plotTypeSelector().
  • plotFlag (bool, optional, default False) – Set plotting True/False. Note that this will plot for all facetDim.
  • facetDim (str, optional, default None) – Dimension to use for subplots. Currently set for a single dimension only. For matplotlib backend: one figure per surface. For plotly backend: subplots per surface.
  • backend (str, optional, default 'mpl' (matplotlib)) – Set backend used for plotting. See sphSumPlotX() for details.
  • fnType (str, optional, default = 'sph') – Set function for expansion parameters, default is YLM from scipy.special.sph_harm. See ep.sphCalc() for details.
Returns:

  • Xarray – Xarray containing the calcualted surfaces I(theta,phi,…)
  • fig – List of figure handles.

epsproc.sphPlot.sphPlotHV(dataIn)[source]
epsproc.sphPlot.sphPlotMPL(dataPlot, theta, phi, convention='phys')[source]

Plot spherical polar function (R,theta,phi) to a Cartesian grid, using Matplotlib.

Parameters:
  • dataPlot (np.array or Xarray) – Values to plot, single surface only, with dims (theta,phi).
  • phi (theta,) – Angles defining spherical polar grid, 2D arrays.
Returns:

Handle to matplotlib figure.

Return type:

fig

epsproc.sphPlot.sphPlotPL(dataPlot, theta, phi, facetDim='Eke', rc=None)[source]

Plot spherical polar function (R,theta,phi) to a Cartesian grid, using Plotly.

Parameters:
  • dataPlot (np.array or Xarray) – Values to plot, single surface only, with dims (theta,phi).
  • phi (theta,) – Angles defining spherical polar grid, 2D arrays.
  • facetDim (str, default 'Eke') – Dimension to use for faceting (subplots), currently set for single dim only.
Returns:

Handle to figure.

Return type:

fig

epsproc.sphPlot.sphSumPlotX(dataIn, pType='a', facetDim='Eke', backend='mpl', convention='phys')[source]

Plot sum of spherical harmonics from an Xarray.

Parameters:
  • dataIn (Xarray) –

    Input structure can be

    • Set of precalculated Ylms, dims (theta,phi) or (theta,phi,LM).
    • Set of precalculated mfpads, dims (theta,phi), (theta,phi,LM) or (theta,phi,LM,facetDim).
    • If (LM) dimension is present, it is summed over before plotting.
    • If facetDim is present this is used for subplots, currently only one facetDim is supported here.
  • pType (char, optional, default 'a' (abs value)) – Set (data) type of plot. See plotTypeSelector().
  • facetDim (str, optional, default Eke) –

    Dimension to use for subplots.

    • Currently set for a single dimension only.
    • For matplotlib backend: one figure per surface.
    • For plotly backend: subplots per surface.
  • backend (str, optional, default 'mpl') –

    Set backend used for plotting.

    • mpl matplotlib: basic 3D plotting, one figure per surface.
    • pl plotly: fancier 3D plotting, interactive in Jupyter but may fail at console.
      Subplots for surfaces.
    • hv holoviews: fancier plotting with additional back-end options.
      Can facet on specific data types.
Returns:

List of figure handles.

Return type:

fig

Examples

>>> YlmX = sphCalc(2, res = 50)
>>> sphSumPlotX(YlmX)

Note

Pretty basic functionality here, should add more colour mapping options and multiple plots, alternative back-ends, support for more dimensions etc.

epsproc.sphPlot.sphToCart(R, theta, phi, convention='phys')[source]

Convert spherical polar coords \((R,\theta,\phi)\) to Cartesian \((X,Y,Z)\).

Parameters:
  • theta, phi (R,) – Spherical polar coords \((R,\theta,\phi)\).
  • convention (str, optional, default = 'phys') – Specify choice of Spherical Polar coordinate system, ‘phys’ or ‘maths’ (see note below).
Returns:

X, Y, Z – Cartesian coords \((X,Y,Z)\).

Return type:

np.arrays

Conversion defined with the usual physics convention , where:

  • \(R\) is the radial distance from the origin
  • \(\theta\) is the polar angle (defined relative to the z-axis), \(0\leq\theta\leq\pi\)
  • \(\phi\) is the azimuthal angle (defined relative to the x-axis), \(0\leq\theta\leq2\pi\)
  • \(X = R * np.sin(phi) * np.cos(theta)\)
  • \(Y = R * np.sin(phi) * np.sin(theta)\)
  • \(Z = R * np.cos(phi)\)

Specify convention = ‘maths’ to use alternative definition with theta and phi swapped.