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)[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)
Returns: ADMX – ADMs in Xarray format, dims as per
epsproc.utils.ADMdimList()
Return type: Xarray
Examples
>>> # Default case >>> ADMX = setADMs() >>> ADMX
>>> # 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.
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')[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, optional, default None) – (Theta, Phi) grid resolution, outputs will be of dim [res,res].
- angs (list of 2D np.arrays, [thetea, phi], optional, default None) – If passed, use these grids for calculation
- 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.
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)[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.
- 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)