epsproc.sphCalc module

ePSproc spherical function calculations.

Collection of functions for calculating Spherical Tensors: Ylm, wignerD etc.

For spherical harmonics, currently using scipy.special.sph_harm

For other functions, using Moble’s spherical_functions package https://github.com/moble/spherical_functions

See tests/Spherical function testing Aug 2019.ipynb

04/12/19 Added setPolGeoms() to define frames as Xarray.
Added setADMs() to define ADMs as Xarray

02/12/19 Added basic TKQ multipole frame rotation routine. 27/08/19 Added wDcalc for Wigner D functions. 14/08/19 v1 Implmented sphCalc

epsproc.sphCalc.TKQarrayRot(TKQ, eAngs)[source]

Frame rotation for multipoles $T_{K,Q}$.

Basic frame rotation code, see https://github.com/phockett/Quantum-Metrology-with-Photoelectrons/blob/master/Alignment/Alignment-1.ipynb for examples.

Parameters:
  • TKQ (np.array) – Values defining the initial distribution, [K,Q,TKQ]
  • eAngs (list or np.array) – List of Euler angles (theta,phi,chi) defining rotated frame.
Returns:

  • TKQRot (np.array) – Multipoles $T’_{K,Q}$ in rotated frame, as an np.array [K,Q,TKQ].
  • TODO (redo with Moble’s functions, and Xarray input & output.)
  • Formalism
  • ———-
  • For the state multipoles, frame rotations are fairly straightforward
  • (Eqn. 4.41 in Blum)
  • .. math:: – begin{equation} leftlangle T(J’,J)_{KQ}^{dagger}rightrangle =sum_{q}leftlangle T(J’,J)_{Kq}^{dagger}rightrangle D(Omega)_{qQ}^{K*} end{equation}
  • Where $D(Omega)_{qQ}^{K}$ is a Wigner rotation operator, for a*
  • rotation defined by a set of Euler angles $Omega={theta,phi,chi}$.
  • Hence the multipoles transform, as expected, as irreducible tensors,
  • i.e. components $q$ are mixed by rotation, but terms of different
  • rank $K$ are not.

epsproc.sphCalc.TKQarrayRotX(TKQin, RX, form=2)[source]

Frame rotation for multipoles $T_{K,Q}$.

Basic frame rotation code, see https://github.com/phockett/Quantum-Metrology-with-Photoelectrons/blob/master/Alignment/Alignment-1.ipynb for examples.

Parameters:
  • TKQin (Xarray) – Values defining the initial distribution, [K,Q,TKQ]. Other dimensions will be propagated.
  • RX (Xarray defining frame rotations, from epsproc.setPolGeoms()) – List of Euler angles (theta,phi,chi) and corresponding quaternions defining rotated frame.
Returns:

TKQRot – Multipoles $T’_{K,Q}$ in rotated frame, as an np.array [K,Q,TKQ].

Return type:

Xarray

For the state multipoles, frame rotations are fairly straightforward (Eqn. 4.41 in Blum):

\[\begin{equation} \left\langle T(J',J)_{KQ}^{\dagger}\right\rangle =\sum_{q}\left\langle T(J',J)_{Kq}^{\dagger}\right\rangle D(\Omega)_{qQ}^{K*} \end{equation}\]

Where $D(Omega)_{qQ}^{K*}$ is a Wigner rotation operator, for a rotation defined by a set of Euler angles $Omega={theta,phi,chi}$. Hence the multipoles transform, as expected, as irreducible tensors, i.e. components $q$ are mixed by rotation, but terms of different rank $K$ are not.

Examples

>>> vFlag = 2
>>> RX = ep.setPolGeoms(vFlag = vFlag)  # Package version
>>> RX
>>> testADMX = ep.setADMs(ADMs=[[0,0,0,1],[2,0,0,0.5]])
>>> testADMX
>>> testADMrot, wDX, wDXre = TKQarrayRotX(testADMX, RX)
>>> testADMrot
>>> testADMrot.attrs['dataType'] = 'ADM'
>>> sph, _ = sphFromBLMPlot(testADMrot, facetDim = 'Euler', plotFlag = True)
epsproc.sphCalc.setADMs(ADMs=[0, 0, 0, 1], KQSLabels=None, t=None, addS=False, name=None, tUnits='ps')[source]

Create Xarray from ADMs, or create default case ADM(K,Q,S) = [0,0,0,1].

Parameters:
  • ADMs (list or np.array, default = [0,0,0,1]) – Set of ADMs = [K, Q, S, ADM]. If multiple ADMs are provided per (K,Q,S) index, they are set to the t axis (if provided), or indexed numerically.
  • KQSLabels (list or np.array, optional, default = None) – If passed, assume ADMs are unabelled, and use (K,Q,S) indicies provided here.
  • t (list or np.array, optional, default = None) – If passed, use for dimension defining ADM sets (usually time). Defaults to numerical label if not passed, t = np.arange(0,ADMs.shape[1])
  • addS (bool, default = False) – If set, append S = 0 to ADMs. This allows for passing of [K,Q,ADM] type values (e.g. for symmetric top case)
  • name (str, optional, default = None) – Set a name for the array. If None, will be set to ‘ADM’ (same as dataType attrib)
  • tUnits (str, optional, default = 'ps') – Units for temporal axis, if set.
Returns:

ADMX – ADMs in Xarray format, dims as per epsproc.utils.ADMdimList()

Return type:

Xarray

Examples

>>> # Default case
>>> ADMX = setADMs()
>>> ADMX
>>> # Set with ranges (as an array [K,Q,S, t(0), t(1)...]), with values per t
>>> tPoints = 10
>>> ADMX = setADMs(ADMs = [[0,0,0, *np.ones(10)], [2,0,0, *np.linspace(0,1,tPoints)], [4,0,0, *np.linspace(0,0.5,tPoints)]])
>>> # With full N2 rotational wavepacket ADM set from demo data (ePSproc\datalignment), where modPath defines root...
>>> # Load ADMs for N2
>>> from scipy.io import loadmat
>>> ADMdataFile = os.path.join(modPath, 'data', 'alignment', 'N2_ADM_VM_290816.mat')
>>> ADMs = loadmat(ADMdataFile)
>>> ADMX = setADMs(ADMs = ADMs['ADM'], KQSLabels = ADMs['ADMlist'], addS = True)
>>> ADMX
epsproc.sphCalc.setBLMs(BLMs=[0, 0, 1], LMLabels=None, t=None, name=None, tUnits='ps', conformDims=False)[source]

Create Xarray from BLMs, or create default case BLM = [0,0,1].

Parameters:
  • BLMs (list or np.array, default = [0,0,1]) – Set of BLMs = [L, M, BLM]. If multiple BLMs are provided per (L,M) index, they are set to the E or t axis (if provided), or indexed numerically.
  • LMLabels (list or np.array, optional, default = None) – If passed, assume BLMs are unabelled, and use (L,M) indicies provided here.
  • t (list or np.array, optional, default = None) – If passed, use for dimension defining BLM sets (usually time or energy). Defaults to numerical label if not passed, t = np.arange(0,ADMs.shape[1])
  • name (str, optional, default = None) – Set a name for the array. If None, will be set to ‘ADM’ (same as dataType attrib)
  • tUnits (str, optional, default = 'ps') – Units for t axis, if set.
  • conformDims (bool, optional, default = False) – Add any missing dims to match ep.listFuncs.BLMdimList. NOT YET IMPLEMENTED - see PEMtk.toePSproc() for method.
Returns:

  • BLMX (Xarray) – BLMs in Xarray format, dims as per epsproc.utils.BLMdimList()
  • TODO
  • - Implemnt conformDims, see PEMtk.toePSproc() for method.
  • - General handling for additional dims, currently only ‘t’ supported. Should change to dict mapping.

Examples

>>> # Default case
>>> BLMX = setBLMs()
>>> BLMX
>>> # Set with ranges (as an array [L,M, t(0), t(1)...]), with values per t
>>> tPoints = 10
>>> BLMX = setBLMs(ADMs = [[0,0, *np.ones(10)], [2,0, *np.linspace(0,1,tPoints)], [4,0, *np.linspace(0,0.5,tPoints)]])
epsproc.sphCalc.setPolGeoms(eulerAngs=None, quat=None, labels=None, vFlag=2)[source]

Generate Xarray containing polarization geometries as Euler angles and corresponding quaternions.

Define LF > MF polarization geometry/rotations. Provide either eulerAngs or quaternions, but not both (supplied quaternions only will be used in this case).

For default case (eulerAngs = None, quat = None), 3 geometries are calculated, corresponding to z-pol, x-pol and y-pol cases. Defined by Euler angles: (p,t,c) = [0 0 0] for z-pol, (p,t,c) = [0 pi/2 0] for x-pol, (p,t,c) = [pi/2 pi/2 0] for y-pol.

Parameters:
  • eulerAngs (list or np.array of Euler angles (p(hi), t(heta), c(hi)), optional.) – List or array [p,t,c…], shape (Nx3). List or array including set labels, [label,p,t,c…], shape (Nx4)
  • quat (list or np.array of quaternions, optional.) –
  • labels (list of labels, one per set of angles. Optional.) – If not set, states will be labelled numerically.
  • vFlag (version of routine to use, optional, default = 2) – Options: - 1, use labels as sub-dimensional coord. - 2, set labels as non-dimensional coord.
Returns:

  • RX (Xarray of quaternions, with Euler angles as dimensional params.)
  • To do
  • —–
  • - Better label handling, as dictionary? With mixed-type array may get issues later. – (sf.quaternion doesn’t seem to have an issue however.)
  • - Xarray MultiIndex with mixed types? – Tested with pd - not supported: >>> eulerInd = pd.MultiIndex.from_arrays([eulerAngs[:,0].T, eulerAngs[:,1:].T.astype(‘float’)], names = [‘Label’,’P’,’T’,’C’]) # Gives error: # NotImplementedError: > 1 ndim Categorical are not supported at this time

Examples

>>> # Defaults
>>> RXdefault = setPolGeoms()
>>> print(RXdefault)
>>> # Pass Eulers, no labels
>>> pRot = [1.666, 0, np.pi/2]
>>> tRot = [0, np.pi/2, np.pi/2]
>>> cRot = [-1.5, 0, 0]
>>> eulerAngs = np.array([pRot, tRot, cRot]).T
>>> RXePass = setPolGeoms(eulerAngs = eulerAngs)
>>> print(RXePass)
>>> # Pass labels separately
>>> RXePass = setPolGeoms(eulerAngs = eulerAngs, labels = ['1','23','ff'])
>>> print(RXePass)
>>> # Pass Eulers with existing labels
>>> labels = ['A','B','C']
>>> eulerAngs = np.array([labels, pRot, tRot, cRot]).T
>>> RXePass = setPolGeoms(eulerAngs = eulerAngs)
>>> print(RXePass)
>>> # Pass Quaternions and labels
>>> RXqPass = setPolGeoms(quat = RXePass, labels = labels)
>>> print(RXqPass)
>>> # Pass both - only quaternions will be used in this case, and warning displayed.
>>> RXqeTest = setPolGeoms(eulerAngs = eulerAngs, quat = RXePass, labels = labels)
>>> print(RXqeTest)
epsproc.sphCalc.sphCalc(Lmax, Lmin=0, res=None, angs=None, XFlag=True, fnType='sph', convention='phys', conj=False)[source]

Calculate set of spherical harmonics Ylm(theta,phi) on a grid.

Parameters:
  • Lmax (int) – Maximum L for the set. Ylm calculated for Lmin:Lmax, all m.
  • Lmin (int, optional, default 0) – Min L for the set. Ylm calculated for Lmin:Lmax, all m.
  • res (int or list, optional, default None) – (Theta, Phi) grid resolution, outputs will be of dim [res,res] (if int) or [res[0], res[1]] (if list).
  • angs (list of 2D np.arrays, [thetea, phi], optional, default None) – If passed, use these grids for calculation. NOTE: if ‘maths’ convention set, this array will be assumed to be [phi, theta]
  • XFlag (bool, optional, default True) – Flag for output. If true, output is Xarray. If false, np.arrays
  • fnType (str, optional, default = 'sph') – Currently can set to ‘sph’ for SciPy spherical harmonics, or ‘lg’ for SciPy Legendre polynomials. More backends to follow.
  • convention (str, optional, default = 'phys') – Set to ‘phys’ (theta from z-axis) or ‘maths’ (phi from z-axis) spherical polar conventions. (Might need some work!)
  • conj (bool, optional, default = False) – Return complex conjugate.

Note that either res OR angs needs to be passed.

  • if XFlag -
YlmX
3D Xarray, dims (lm,theta,phi)
  • else -
Ylm, lm
3D np.array of values, dims (lm,theta,phi), plus list of lm pairs
Currently set for scipy.special.sph_harm as calculation routine. Note (theta, phi) definition, and normalisation.
https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.sph_harm.html

Example

>>> YlmX = sphCalc(2, res = 50)
epsproc.sphCalc.wDcalc(Lrange=[0, 1], Nangs=None, eAngs=None, R=None, XFlag=True, QNs=None, dlist=['lp', 'mu', 'mu0'], eNames=['P', 'T', 'C'], conjFlag=False, sfError=True, verbose=False)[source]

Calculate set of Wigner D functions D(l,m,mp; R) on a grid.

Parameters:
  • Lrange (list, optional, default [0, 1]) – Range of L to calculate parameters for. If len(Lrange) == 2 assumed to be of form [Lmin, Lmax], otherwise list is used directly. For a given l, all (m, mp) combinations are calculated.
  • QNs (np.array, optional, default = None) – List of QNs [l,m,mp] to compute Wigner D terms for. If supplied, use this instead of Lrange setting.
  • for setting angles (use one only) (Options) –
  • Nangs (int, optional, default None) – If passed, use this to define Euler angles sampled. Ranges will be set as (theta, phi, chi) = (0:pi, 0:pi/2, 0:pi) in Nangs steps.
  • eAngs (np.array, optional, default None) – If passed, use this to define Euler angles sampled. Array of angles, [theta,phi,chi], in radians
  • R (np.array, optional, default None) – If passed, use this to define Euler angles sampled. Array of quaternions, as given by quaternion.from_euler_angles(eAngs).
XFlag : bool, optional, default True
Flag for output. If true, output is Xarray. If false, np.arrays
dlist : list, optional, default [‘lp’,’mu’,’mu0’]
Labels for Xarray QN dims.
eNames : list, optional, default [‘P’,’T’,’C’]
Labels for Xarray Euler dims.
conjFlag : bool, optional, default = False
If true, return complex conjuage values.
sfError : bool, optional, default = None
If not None, set sf.error_on_bad_indices = sfError If True (default case), this will raise a value error on bad QNs. If False, set = 0. See code at https://github.com/moble/spherical_functions/blob/master/spherical_functions/WignerD/__init__.py 05/05/21 - added, but CURRENTLY NOT WORKING
  • if XFlag -
wDX
Xarray, dims (lmmp,Euler)
  • else -
wD, R, lmmp
np.arrays of values, dims (lmmp,Euler), plus list of angles and lmmp sets.
Uses Moble's spherical_functions package for wigner D function.
https://github.com/moble/spherical_functions
Moble's quaternion package for angles and conversions.
https://github.com/moble/quaternion
For testing, see https://epsproc.readthedocs.io/en/latest/tests/Spherical_function_testing_Aug_2019.html

Examples

>>> wDX1 = wDcalc(eAngs = np.array([0,0,0]))
>>> wDX2 = wDcalc(Nangs = 10)