epsproc package

Subpackages

Submodules

Module contents

epsproc.ADMdimList(sType='stacked')[source]

Return standard list of dimensions for frame definitions, from epsproc.sphCalc.setADMs().

Parameters

sType (string, optional, default = 'stacked') – Selected ‘stacked’ or ‘unstacked’ dimensions. Set ‘sDict’ to return a dictionary of unstacked <> stacked dims mappings for use with xr.stack({dim mapping}).

Returns

list

Return type

set of dimension labels.

epsproc.AFBLMCalcLoop(matE, AKQS=numpy.array, eAngs=[0, 0, 0], thres=1e-06, p=0, R=0, verbose=1)[source]

Calculate inner loop for MFBLMs, based on passed set of matrix elements (Xarray).

Loop code based on Matlab code ePSproc_MFBLM.m

This works, but not very clean or efficient - should be possible to parallelize & vectorise… but this is OK for testing, benchmarking and verification purposes.

Parameters
  • matE (Xarray) – Contains one set of matrix elements to use for calculation. Currently assumes these are a 1D list, with associated (l,m,mu) parameters, as set by mfblm().

  • AKQS (Xarray or np.array, optional, default = np.array([0,0,0,1], ndmin = 2)) – Array containing alignment parameters (axis distribution moments), $A_{Q,S}^{K}(t)$. Format is [K,Q,S,value…] if np.array. For Xarray, (K,Q,S) coords are set as multindex dimension ADM, as per epsproc.setADMs(). Default value corresponds to an isotropic axis distribution.

  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.

  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.

  • p (int, optional, default p = 0) – LF polarization term. Currently only valid for p = 0

  • R (int, optional, default R = 0) – LF polarization term (from tensor contraction). Currently only valid for p = 0

  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.

    • 1: Print basic info.

    • 2: Print intermediate C parameter array to terminal when running.

Returns

  • BLMX (Xarray) – Set of B(L,M; eAngs, Eke) terms for supplied matrix elements, in an Xarray. For cases where no values are calculated (below threshold), return an array with B00 = 0 only.

  • Limitations & To Do

  • ——————–

  • * Currently set with (p,R) values passed, but only valid for (0,0) (not full sum over R terms as shown in formalism above.)

  • * Set to accept a single set of matrix elements (single E), assuming looping over E (and other parameters) elsewhere.

  • * Not explicitly parallelized here, should be done by calling function.

  • (Either via Xarray methods, or numba/dask…? http (//xarray.pydata.org/en/stable/computation.html#wrapping-custom-computation))

  • * Coded for ease, not efficiency - there will be lots of repeated ang. mom. calcs. when run over many sets of matrix elements.

  • * Scale factor currently not propagated.

class epsproc.Arrow3D(*args: Any, **kwargs: Any)[source]

Bases: FancyArrowPatch

Define Arrow3D plotting class

Code verbatim from StackOverflow post https://stackoverflow.com/a/22867877 Thanks to CT Zhu, https://stackoverflow.com/users/2487184/ct-zhu

draw(renderer)[source]
epsproc.BLMdimList(sType='stacked')[source]

Return standard list of dimensions for calculated BLM.

Parameters

sType (string, optional, default = 'stacked') – Selected ‘stacked’ or ‘unstacked’ dimensions. Set ‘sDict’ to return a dictionary of unstacked <> stacked dims mappings for use with xr.stack({dim mapping}).

Returns

list

Return type

set of dimension labels.

epsproc.BLMplot(BLM, thres=0.01, thresType='abs', selDims=None, XS=False, xDim='Eke', pType='r', col='Labels', row=None, pStyle='line', backend='xr', returnPlot=False, renderPlot=True, **kwargs)[source]

Plotting routines for BLM values from Xarray. Plot line or surface plot, with various backends available.

Parameters
  • BLM (Xarray) – Input data for plotting, dims assumed to be as per ePSproc.util.BLMdimList()

  • thres (float, optional, default 1e-2) – Value used for thresholding results, only values > thres will be included in the plot. Either abs or relative (%) value, according to thresType setting. Set thres = None to skip.

  • thresType (str, optional, default = 'abs') – Set to ‘abs’ or ‘pc’ for absolute threshold, or relative value (%age of max value in dataset)

  • selDims (dict, optional, default = None) – Dimensions to select from input Xarray.

  • XS (bool, optional, default = False) – Use absolute XS value for l=0 if true (from BLM.XS); otherwise show current values (usually normalised to =1).

  • xDim (str, optional, default = 'Eke') – Dimension to use for x-axis, also used for thresholding. Default plots (Eke, BLM) surfaces with subplots for (Euler). Change to ‘Euler’ to plot (Euler, BLM) with (Eke) subplots.

  • col (strings, optional, default = 'Labels', None) – Set for dim handling in plot. For xr backend these define plot layout.

  • row (strings, optional, default = 'Labels', None) – Set for dim handling in plot. For xr backend these define plot layout.

  • pType (char, optional, default 'a' (abs values)) – Set (data) type to plot. See plotTypeSelector().

  • pStyle (string, optional, default = 'line') – NOT YET PROPERLY IMPLEMENTED For XR backend can also set ‘surf’.

  • backend (str, optional, default = 'xr') – Plotter to use. Currently supports ‘xr’ = Xarray, ‘hv’ = Holoviews Default is ‘xr’ for Xarray internal plotting. May be switched according to plot type in future… **kwargs are currently passed to the plotter backend.

  • returnXR (bool, optional, default = False) – Return plot data to caller if true.

  • returnPlot (bool, optional, default = False) – For hv backend, return object to caller?

  • renderPlot (bool, optional, default = True) – For hv backend, render plot (if running in a notebook)?

Notes

Updates in issue #27: https://github.com/phockett/ePSproc/issues/27 Demo: https://epsproc.readthedocs.io/en/dev/tests/hvPlotters_fn_tests_150720_v250122.html

  • Proper dim handling to be implemented. See code elsewhere(?).

    XR plotter currently uses ‘cDims’ in some places too.

  • Add XR and HV data return, currently only return objects.

27/01/22: Updated with better selection & basic XR dim handling, but still needs work.

Added basic HV plotter routines.

epsproc.EDCSFileParse(fileName, verbose=1)[source]

Parse an ePS file for EDCS segments.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool or int, optional) – If true, print segment info.

Returns

  • list – [lineStart, lineStop], ints for line #s found from start and end phrases.

  • list – dumpSegs, list of lines read from file.

  • Lists contain entries for each dumpIdy segment found in the file.

epsproc.EDCSSegParse(dumpSeg, verbose=False)[source]

Extract values from EDCS file segments.

Parameters
  • dumpSeg (list) – One EDCS segment, from dumpSegs[], as returned by epsproc.IO.EDCSFileParse()

  • verbose (bool, int, default = False) – Print additional info during run.

Returns

  • np.array – EDCS, array of scattering XS, [theta, Cross Section (Angstrom^2)]

  • list – attribs, list [Label, value, units]

Notes

Only reads (theta,I) information from an “EDCS - differential cross section program” segment. Currently this is a bit messy, and relies on fixed EDCS format. No error checking as yet. Not yet reading all attribs.

Example

>>> EDCS, attribs = EDCSSegParse(dumpSegs[0])
epsproc.EDCSSegsParseX(dumpSegs)[source]

Extract data from ePS EDCS segments into usable form.

Parameters

dumpSegs (list) – Set of dumpIdy segments, i.e. dumpSegs, as returned by epsproc.IO.EDCSFileParse()

Returns

  • xr.array – Xarray data array, containing cross sections. Dimensions (Eke, theta)

  • int – Number of blank segments found. (CURRENTLY not implemented.)

Example

>>> data = EDCSSegsParseX(dumpSegs)

Notes

A rather cut-down version of epsproc.IO.dumpIdySegsParseX(), no error checking currently implemented.

epsproc.EPR(QNs=None, p=None, ep=None, nonzeroFlag=True, form='2d', dlist=['l', 'lp', 'P', 'p', 'R-p', 'R'], phaseConvention='S')[source]

Define polarization tensor (LF) for 1-photon case.

Define field terms (from QM book, corrected vs. original S&U version - see beta-general-forms\_rewrite\_290917.lyx):

\[\begin{split}\begin{equation} E_{PR}(\hat{e})=[e\otimes e^{*}]_{R}^{P}=[P]^{\frac{1}{2}}\sum_{p}(-1)^{R}\left(\begin{array}{ccc} 1 & 1 & P\\ p & R-p & -R \end{array}\right)e_{p}e_{R-p}^{*} \end{equation}\end{split}\]
Parameters
  • QNs (np.array, optional, default = None) – List of QNs [l, lp, P, m, mp, R] to compute 3j terms for. If not supplied all allowed terms for the one-photon case, l=lp=1, will be set.

  • p (list or array, optional, default = None) – Specify polarization terms p. If set to None, all allowed values for the one-photon case will be set, p=[-1,0,1]

  • ep (list or array, optional, default = None) – Relative strengths for the fields ep. If set to None, all terms will be set to unity, ep = 1

  • nonzeroFlag (bool, optional, default = True) – Drop null terms before returning values if true.

  • form (string, optional, default = '2d') – For options see ep.w3jTable()

  • phaseConvention (optional, str, default = 'S') –

    Set phase conventions:

    • ’S’ : Standard derivation.

    • ’R’ : Reduced form geometric tensor derivation.

    • ’E’ : ePolyScat, may have additional changes in numerics, e.g. conjugate Wigner D.

    See setPhaseConventions() for more details.

Examples

# Generate full EPR list with defaults >>> EPRtable = EPR()

# Return as Xarray >>> EPRtable = EPR(form = ‘xarray’)

Note

Currently not handling ep correctly! Should implement as passed Xarray for correct (p,p’) assignment.

class epsproc.Efield(Edef=None, units=None, verbose=True)[source]

Bases: object

Basic class for handling E-field generation, data and related functions.

Currently set for field creation based on Gaussian defined in time-domain. Defaults to a single (E,t) point if no parameters are supplied - this can be used with epsproc.geomFunc.EPR()

Starting to look a bit top-heavy… passing directly to Xarray might save some structural overhead here?

Define Gaussian field in the time-domain, with amplitude $E_0$, width $sigma$ and carrier frequency $omega_0$:

\[ \begin{align}\begin{aligned}E(\omega_0,t)=E_{0}e^{-(t^2/2\sigma^{2})}e^{i\omega_0 t}\\Carrier frequency in rad/s:\end{aligned}\end{align} \]
\[ \begin{align}\begin{aligned}\omega_0 = 2\pi \Omega_0\\FWHM and Gaussian width (see https://en.wikipedia.org/wiki/Gaussian_function for full definitions):\end{aligned}\end{align} \]
\[ \begin{align}\begin{aligned}FWHM = 2\sigma\sqrt{2ln2}\\Parameters ---------- If no parameters are passed, return single point with E0 = 1 For a time-dependent field, pass at least sigma.\\E0 : float, optional, default = None Set (relative) field strength. If None, set E0 = 1\\sigma : float, optional, default = None Set Gaussian width :math:\sigma, where :math:FWHM=2\sigma\sqrt{2ln2}\\t : array, optional, default = None Set time axis for E-field. If None, defaults to np.arange(-(5*sigma), 5*sigma, self.dt), and self.dt = 0.01*sigma\\dt : float, optional, default = None Set time-step. Default will set from t-axis, or self.dt = 0.01*sigma if sigma only passed. Pass here to override.\\l0 : float, optional, default = None Central wavelength in wavelength working units. If supplied this will be used to set f0 and w0.\\f0, w0 : float, optional, default = None Carrier frequency in freq. units, and :math:\omega_0 in (rad/s). If only f0 supplied, w0 will be calculated. If f0 and w0 supplied, w0 will be used directly. If None, this will be ignored.\\units : dictionary, optional, default = None Set units. If not passed, default dict will be created with units (THz, nm, fs) set.\\ TODO: - best way to handle unit conversions? Did some of this for FROG code already? Yes, see python est_codes\E_fields_redux_231118_class.py - Consider other structures here, dataclass or namedtuple? https://stackoverflow.com/questions/354883/how-do-i-return-multiple-values-from-a-function - Xarray conversion? Use Xarray directly in class? - Bug with FFT/iFFT code? If shift=True and ishift=True then get extra freq components on iFFT - must be bug in axis assignment/conversion? OK if shift = False.\end{aligned}\end{align} \]
ECalc()[source]
EFFT()[source]
EiFFT(Emod=None, Ninput=None, comment='')[source]
FWHMGau()[source]
Gau()[source]
calcSpectrogram(signal=None, gate=None, fType='blind', N=None, verbose=False)[source]
chirp(A, resetPhase=False, comment=None)[source]

Add quadratic phase (chirp) in spectral domain. Requires an E-field object, and chirp parameter A.

TODO: check Diels for nomenclature here.

fConv()[source]

Convert freq. domain axis to lenght & energy units. and set fields.

finst()[source]

Calculate f(t) as instantaneous freq.

phaseMask(domain=None, thres=None)[source]

Mask pulse phase away from main features. Should convert to use masked arrays for reversibility, or just return a mask (currently overwrites Ef).

plot(plotType='phaseUW', Nplot=None, domainList=None, thres=0.01)[source]

Basic plot with Matplotlib

plotSpectrogram()[source]

VERY basic spectrogram plotter from old code - for quick testing only.

TODO: fix axes for with/withou fft shift.

BETTER: use plotSpectrogramHV() instead, this uses full axes as set.

plotSpectrogramHV(N=None)[source]
printDef()[source]
removePhase(domain=None)[source]

Remove imaginary field terms (reset phase).

setEf()[source]

Calculate pulse properties and E-fields based on Edef.

setPhase(phaseVec=None, domain=None)[source]

Set imaginary field terms (set phase) for specified domain & recalculated FFT.

toXarrayDS(N=None)[source]
epsproc.MFBLMCalcLoop(matE, eAngs=[0, 0, 0], thres=1e-06, p=0, R=0, verbose=1)[source]

Calculate inner loop for MFBLMs, based on passed set of matrix elements (Xarray).

Loop code based on Matlab code ePSproc_MFBLM.m

This works, but not very clean or efficient - should be possible to parallelize & vectorise… but this is OK for testing, benchmarking and verification purposes.

Parameters
  • matE (Xarray) – Contains one set of matrix elements to use for calculation. Currently assumes these are a 1D list, with associated (l,m,mu) parameters, as set by epsproc.MFBLM.mfblm().

  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.

  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.

  • p (int, optional, default p = 0) – LF polarization term. Currently only valid for p = 0

  • R (int, optional, default R = 0) – LF polarization term (from tensor contraction). Currently only valid for p = 0

  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.

    • 1: Print basic info.

    • 2: Print intermediate C parameter array to terminal when running.

Returns

  • BLMX (Xarray) – Set of B(L,M; eAngs, Eke) terms for supplied matrix elements, in an Xarray. For cases where no values are calculated (below threshold), return an array with B00 = 0 only.

  • Limitations & To Do

  • ——————–

  • * Currently set with (p,R) values passed, but only valid for (0,0) (not full sum over R terms as shown in formalism above.)

  • * Set to accept a single set of matrix elements (single E), assuming looping over E (and other parameters) elsewhere.

  • * Not explicitly parallelized here, should be done by calling function.

  • (Either via Xarray methods, or numba/dask…? http (//xarray.pydata.org/en/stable/computation.html#wrapping-custom-computation))

  • * Coded for ease, not efficiency - there will be lots of repeated ang. mom. calcs. when run over many sets of matrix elements.

  • .. rst-class:: strike

  • * Scale factor currently not propagated.

  • * (SF now propagated via Xarrays and implemented in main calling function )

epsproc.MFproj(QNs=None, RX=None, nonzeroFlag=True, form='2d', dlist=['l', 'lp', 'P', 'mu', 'mup', 'Rp', 'R'], eNames=['Ph', 'Th', 'Ch'], phaseConvention='S')[source]

Define MF projection term, \(\Lambda_{R',R}(R_{\hat{n}})\):

\[\begin{split}\begin{equation} \Lambda_{R',R}(R_{\hat{n}})=(-1)^{(R')}\left(\begin{array}{ccc} 1 & 1 & P\\ \mu & -\mu' & R' \end{array}\right)D_{-R',-R}^{P}(R_{\hat{n}}) \end{equation}\end{split}\]

Then…

\[\begin{eqnarray} \beta_{L,-M}^{\mu_{i},\mu_{f}} & = & \sum_{P,R',R}{\color{red}E_{P-R}(\hat{e};\mu_{0})}\sum_{l,m,\mu}\sum_{l',m',\mu'}(-1)^{(\mu'-\mu_{0})}{\color{red}\Lambda_{R',R}(R_{\hat{n}};\mu,P,R,R')B_{L,-M}(l,l',m,m')}I_{l,m,\mu}^{p_{i}\mu_{i},p_{f}\mu_{f}}(E)I_{l',m',\mu'}^{p_{i}\mu_{i},p_{f}\mu_{f}*}(E) \end{eqnarray}\]
Parameters
  • phaseConvention (optional, str, default = 'S') – Set phase conventions: - ‘S’ : Standard derivation. - ‘R’ : Reduced form geometric tensor derivation. - ‘E’ : ePolyScat, conjugate Wigner D. See setPhaseConventions() for more details.

  • eNames (optional, list, default = ['Ph','Th','Ch']) –

    Set names for Euler angles in output. Note:

    • eNames = [‘P’,’T’,’C’] matches defaults in epsproc.sphCalcs.setPolGeoms() and epsproc.sphCalcs.wDcalc(), but conflicts with QNs ‘P’.

    • eNames = [‘Phi’,’Theta’,’Chi’] may give issues elsewhere, e.g. when multiplying by Ylms with same dim names.

Notes

This is very similar to \(E_{PR}\) term.

12/08/22 Added option for Euler angle dim names to avoid conflicts elsewhere.

Examples

>>> lTerm, lambdaTable, lambdaD, QNs = MFproj(form = 'xarray')
epsproc.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.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

Formalism

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.Wigner3jCached(j_1, j_2, j_3, m_1, m_2, m_3)[source]

Wrapper for 3j caching with functools.lru_cache

epsproc.Wigner3jQNs(QNs)[source]
epsproc.Wigner_D_element_Cached(pRot, P, Rp, R)[source]

Wrapper for WignerD caching with functools.lru_cache

epsproc.afblm(daIn, selDims={'Type': 'L'}, AKQS=numpy.array, eAngs=[0, 0, 0], thres=0.0001, sumDims=('l', 'm', 'mu', 'Cont', 'Targ', 'Total', 'it'), SFflag=True, verbose=1)[source]

Calculate MFBLMs for a range of (E, sym) cases. Default is to calculated for all symmetries at each energy.

Parameters
  • da (Xarray) – Contains matrix elements to use for calculation. Matrix elements will be sorted by energy and BLMs calculated for each set.

  • selDims (dict, optional, default = {'Type':'L'}) – Additional sub-selection to be applied to matrix elements before BLM calculation. Default selects just length gauge results.

  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.

  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.

  • sumDims (tuple, optional, default = ('l','m','mu','Cont','Targ','Total','it')) – Defines which terms are summed over (coherent) in the MFBLM calculation. (These are used to flatten the Xarray before calculation.) Default includes sum over (l,m), symmetries and degeneracies (but not energies).

  • SFflag (bool, default = True) – Normalise by scale factor to give X-sections (B00) in Mb

  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.

    • 1: Print basic info.

    • 2: Print intermediate C parameter array to terminal when running.

Returns

  • Xarray – Calculation results BLM, dims (Euler, Eke, l,m). Some global attributes are also appended.

  • Limitations

  • ———–

  • Currently set to loop calcualtions over energy only, and all symmetries.

  • Pass single {‘Cont’ (‘sym’} to calculated for only one symmetry group.)

  • TODO (In future this will be more elegant.)

  • TODO (Setting selDims in output structure needs more thought for netCDF save compatibility.)

epsproc.afblmXprod(matEin, QNs=None, AKQS=None, EPRX=None, p=[0], BLMtable=None, BLMtableResort=None, lambdaTerm=None, polProd=None, AFterm=None, thres=0.01, thresDims='Eke', selDims={'Type': 'L'}, sqThres=True, dropThres=True, sumDims=['mu', 'mup', 'l', 'lp', 'm', 'mp', 'S-Rp'], sumDimsPol=['P', 'R', 'Rp', 'p'], symSum=True, degenDrop=True, SFflag=False, SFflagRenorm=False, BLMRenorm=1, squeeze=False, phaseConvention='E', basisReturn='BLM', verbose=0, **kwargs)[source]

Implement \(\beta_{LM}^{AF}\) calculation as product of tensors.

\[\begin{eqnarray} \beta_{L,-M}^{\mu_{i},\mu_{f}} & =(-1)^{M} & \sum_{P,R',R}{[P]^{\frac{1}{2}}}{E_{P-R}(\hat{e};\mu_{0})}\sum_{l,m,\mu}\sum_{l',m',\mu'}(-1)^{(\mu'-\mu_{0})}{\Lambda_{R'}(\mu,P,R')B_{L,-M}(l,l',m,m')}I_{l,m,\mu}^{p_{i}\mu_{i},p_{f}\mu_{f}}(E)I_{l',m',\mu'}^{p_{i}\mu_{i},p_{f}\mu_{f}*}(E)\sum_{K,Q,S}\Delta_{L,M}(K,Q,S)A_{Q,S}^{K}(t) \end{eqnarray}\]

Where each component is defined by fns. in :py:module:`epsproc.geomFunc.geomCalc` module.

24/11/22 Added sqThres = True, dropThres = True for use with INITIAL matEleSelector() call only - in some cases may get issues with dim drops here otherwise.

04/05/22 Added **kwargs, unused but allows for arb basis dict unpack and passing from other functions. May want to pipe back to Full basis return however.

06/08/21 Added basic handling for degenerate states, including degenDrop option.

Updated docs, but still rather messy!

27/07/21 Removed eulerAngs & RX input options, since these are redundant (and lead to confusion here!).

For cases where E-field and alignment distribution are rotated, set AKQS in rotated frame, see https://epsproc.readthedocs.io/en/latest/tests/ePSproc_frame_rotation_tests_Dec2019.html Also removed selDims={‘it’:1}, which can result in issues! In current code, adding ‘it’ to sumDims doesn’t work (dim issues somewhere), so may be best to treat independently…?

03/05/21 Tidying up a bit & improving/wrapping for fitting use (inc. basis function reuse).

10/09/20 Verified (but messy) version, with updated defaults.

01/09/20 Verified (but messy) version, including correct renormalisation routines.

15/06/20 In progress! Using mfblmXprod() as template, with just modified lambda term, and new alignment term, to change.

For basic use see the docs: https://epsproc.readthedocs.io/en/dev/demos/ePSproc_class_demo_161020.html#Compute-LF/AF-beta_{LM}-and-PADs

Dev code:

Parameters
  • matE (Xarray) – Xarray containing matrix elements, with QNs (l,m), as created by readMatEle()

  • settings (*** Optional calculation) –

  • selDims (dict, default = {'Type':'L'}) – Selection parameters for calculations, may be be checked and appened herein.

  • sumDims (list, default = ['mu', 'mup', 'l','lp','m','mp','S-Rp']) – Main summation dims, will be checked herein.

  • sumDimsPol (list, default = ['P','R','Rp','p']) – Additional polarization summation dims.

  • symSum (bool, default = True) – Sum over symmetries sets (=={Cont, Targ, Total}) if true.

  • degenDrop (bool) – Flag to set dropping of degenerate components.

  • thres (float, default = 1e-2) – Set threshold value, used to set input matrix elements and again for outputs.

  • thresDims (str, default = 'Eke') – Set threshold dimension (set to be contiguous).

  • verbose (bool or int) – Print output?

*** Optional renormalisation settings (mainly for testing only)

SFflagbool, default = False

Multiply input matrix elements by complex scale-factor if true.

SFflagRenormbool, default = False

Renorm output BLMs by complex scale-factor if true.

BLMRenormint, default = 1

Set different BLM renorm conventions. If 1 renorm by B00. See code for further details.

sqThresbool, default = True

Used for initial matrix element threholding call only. Default = True as previous code, but can cause issues in custom cases (usually accidentally dropping singleton mu coord). Added 24/11/22

dropThresbool, default = True

Used for initial matrix element threholding call only. Default = True as previous code, but can cause issues in custom cases (usually accidentally dropping singleton mu coord). Added 24/11/22

squeezebool, default = False

Squeeze output array after thresholding? Note: this may cause dim issues if True.

*** Optional input data/basis functions (mainly for fitting routine use)

QNsnp.array, optional, default = None

List of QNs as generated by genllpMatE(). Will be generated if not passed.

AKQSXarray, optional, default = None

Alignment parameters, as set by setADMs(). Defaults to isotropic case if not passed.

EPRXXarray, optional, default = None

E-field parameters, as generated by EPR(). Defaults to normalised/unity field, pol = p (below).

plist or array, optional, default = [0]

Specify polarization terms p. Possibly currently only valid for p=0, TBC See https://epsproc.readthedocs.io/en/latest/methods/geometric_method_dev_260220_090420_tidy.html#E_{P,R}-tensor

BLMtable, BLMtableResortXarrays, optional, default = None

Beta calculation parameters, as defined by geomCalc.betaTerm(). BLMtableResort includes phase settings & param renaming as herein.

lambdaTermXarray, optional, default = None

Lambda term parameters, as defined by geomCalc.MFproj()

AFtermXarray, optional, default = None

Alignment term as defined by geomCalc.deltaLMKQS()

polProdXarray, optional, default = None

Polarization tensor as defined by EPRXresort * lambdaTermResort * AFterm

phaseConventionoptional, str, default = ‘E’

Set phase conventions with epsproc.geomCalc.setPhaseConventions(). To use preset phase conventions, pass existing dictionary.

basisReturnoptional, str, default = “BLM”
  • ‘BLM’ return Xarray of results only.

  • ‘Full’ return Xarray of results + basis set dictionary as set during the run.

  • ‘Product’, as full, but minimal basis set with products only.

  • ‘Results’ or ‘Legacy’ direct return of various calc. results Xarrays.

**kwargs, unused but allows for arb basis dict unpack and passing from other functions.

Returns

  • Xarray – Set of AFBLM calculation results

  • dict – Dictionary of basis functions, only if basisReturn != ‘BLM’ (see basisReturn paramter notes).

Notes

Cross-section outputs now set as:

  • XSraw = direct AF calculation output.

  • XSrescaled = XSraw * sqrt(4pi)

  • XSiso = direct sum over matrix elements

Where XSrescaled == XSiso == ePS GetCro output for isotropic distribution. Optionally set SFflag = True to multiply by (complex) scale-factor.

epsproc.arraySort2D(a, col)[source]

Sort np.array a by specified column col. From https://thispointer.com/sorting-2d-numpy-array-by-column-or-row-in-python/

epsproc.betaTerm(QNs=None, Lmin=0, Lmax=10, nonzeroFlag=True, form='2d', dlist=['l', 'lp', 'L', 'm', 'mp', 'M'], phaseConvention='S')[source]

Define BLM coupling tensor

Define field terms (from QM book, corrected vs. original S&U version - see beta-general-forms\_rewrite\_290917.lyx):

\[\begin{split}\begin{equation} B_{L,M}=(-1)^{m}\left(\frac{(2l+1)(2l'+1)(2L+1)}{4\pi}\right)^{1/2}\left(\begin{array}{ccc} l & l' & L\\ 0 & 0 & 0 \end{array}\right)\left(\begin{array}{ccc} l & l' & L\\ -m & m' & M \end{array}\right) \end{equation}\end{split}\]
Parameters
  • QNs (np.array, optional, default = None) – List of QNs [l, lp, L, m, mp, M] to compute 3j terms for. If not supplied all allowed terms for {Lmin, Lmax} will be set. (If supplied, values for Lmin, Lmax and mFlag are not used.)

  • Lmin (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • Lmax (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • mFlag (bool, optional, default = True) – m, mp take all values -l…+l if mFlag=True, or =0 only if mFlag=False

  • nonzeroFlag (bool, optional, default = True) – Drop null terms before returning values if true.

  • form (string, optional, default = '2d') – For options see ep.w3jTable()

  • dlist (list of labels, optional, default ['l','lp','L','m','mp','M']) – Used to label array for Xarray output case.

  • phaseConvention (optional, str, default = 'S') – Set phase conventions: - ‘S’ : Standard derivation. - ‘R’ : Reduced form geometric tensor derivation. - ‘E’ : ePolyScat, may have additional changes in numerics, e.g. conjugate Wigner D. See setPhaseConventions() for more details.

Examples

>>> Lmax = 2
>>> BLMtable = betaTerm(Lmax = Lmax, form = 'xds')
>>> BLMtable = betaTerm(Lmax = Lmax, form = 'xdaLM')
epsproc.blmXarray(BLM, Eke)[source]

Create Xarray from BLM list, format BLM = [L, M, Beta], at a single energy Eke.

Array sorting only valid for 2D BLM array, for B00=0 case pass BLM = None.

epsproc.conv_ev_atm(data, to='ev')[source]

Convert eV <> Hartree (atomic units)

Parameters
  • data (int, float, np.array) – Values to convert.

  • to (str, default = 'ev') –

    • ‘ev’ to convert H > eV

    • ’H’ to convert eV > H

Return type

data converted in converted units.

epsproc.conv_ev_nm(data)[source]

Convert E(eV) <> lambda(nm).

epsproc.dataGroupSel(data, dInd)[source]
epsproc.dataTypesList()[source]

Return a dict of allowed dataTypes, corresponding to epsproc processed data.

Each dataType lists ‘source’, ‘desc’ and ‘recordType’ fields.

  • ‘source’ fields correspond to ePS functions which get or generate the data.

  • ‘desc’ brief description of the dataType.

  • ‘recordType’ gives the required segment in ePS files (and associated parser). If the segment is not present in the source file, then the dataType will not be available.

  • ‘def’ provides definition function handle (if applicable).

  • ‘dims’ lists results from def(sType = ‘sDict’)

TODO: best choice of data structure here? Currently nested dictionary.

epsproc.deltaLMKQS(EPRX, AKQS, phaseConvention='S')[source]

Calculate aligned-frame “alignment” term:

\[\begin{equation} \sum_{K,Q,S}\Delta_{L,M}(K,Q,S)A_{Q,S}^{K}(t) \end{equation}\]
\[\begin{split}\begin{equation} \Delta_{L,M}(K,Q,S)=(2K+1)^{1/2}(-1)^{K+Q}\left(\begin{array}{ccc} P & K & L\\ R & -Q & -M \end{array}\right)\left(\begin{array}{ccc} P & K & L\\ R' & -S & S-R' \end{array}\right) \end{equation}\end{split}\]

15/06/20 IN PROGRESS

Parameters
  • EPRX (Xarray) – Polarization terms in an Xarray, as set by epsproc.geomCalc.EPR()

  • AKQS (Xarray) – Alignement terms in an Xarray, as set by epsproc.setADMs()

Returns

  • AFterm (Xarray) – Full term, including multiplication and sum over (K,Q,S) (note S-Rp term is retained).

  • DeltaKQS (Xarray) – Alignment term \(\Delta_{L,M}(K,Q,S)\).

  • To do

  • —–

  • - Add optional inputs.

  • - Add error checks.

  • See other similar functions for schemes.

epsproc.dumpIdyFileParse(fileName, verbose=1)[source]

Parse an ePS file for dumpIdy segments.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool or int, optional) – If true, print segment info.

Returns

  • list – [lineStart, lineStop], ints for line #s found from start and end phrases.

  • list – dumpSegs, list of lines read from file.

  • Lists contain entries for each dumpIdy segment found in the file.

epsproc.dumpIdySegParse(dumpSeg)[source]

Extract values from dumpIdy file segments.

Parameters

dumpSeg (list) – One dumpIdy segment, from dumpSegs[], as returned by epsproc.IO.dumpIdyFileParse()

Returns

  • np.array – rawIdy, array of matrix elements, [m,l,mu,ip,it,Re,Im]

  • list – attribs, list [Label, value, units]

Notes

Currently this is a bit messy, and relies on fixed DumpIDY format. No error checking as yet. Not yet reading all attribs.

Example

>>> matEle, attribs = dumpIdySegParse(dumpSegs[0])
epsproc.dumpIdySegsParseX(dumpSegs, ekeListUn, symSegs, verbose=1)[source]

Extract data from ePS dumpIdy segments into usable form.

Parameters
  • dumpSegs (list) – Set of dumpIdy segments, i.e. dumpSegs, as returned by epsproc.IO.dumpIdyFileParse()

  • ekeListUn (list) – List of energies, used for error-checking and Xarray rearraging, as returned by epsproc.IO.scatEngFileParse()

  • verbose (bool, default True) – Print job info from file header if true.

Returns

  • xr.array – Xarray data array, containing matrix elements etc. Dimensions (LM, Eke, Sym, mu, it, Type)

  • int – Number of blank segments found.

Example

>>> data = dumpIdySegsParseX(dumpSegs)
epsproc.eulerDimList(sType='stacked')[source]

Return standard list of dimensions for frame definitions, from epsproc.sphCalc.setPolGeoms().

Parameters

sType (string, optional, default = 'stacked') – Selected ‘stacked’ or ‘unstacked’ dimensions. Set ‘sDict’ to return a dictionary of unstacked <> stacked dims mappings for use with xr.stack({dim mapping}).

Returns

list

Return type

set of dimension labels.

epsproc.fileParse(fileName, startPhrase=None, endPhrase=None, comment=None, verbose=0)[source]

Parse a file, return segment(s) from startPhrase:endPhase, excluding comments.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • startPhrase (str, optional) – Phrase denoting start of section to read. Default = None

  • endPhase (str or list, optional) – Phrase denoting end of section to read. Default = None

  • comment (str, optional) – Phrase denoting comment lines, which are skipped. Default = None

  • verbose (int, optional, default = 1) – Level of verbosity in output. - 0 no printed output - 1 print summary info only - 2 print detailed info

Returns

  • list – [lineStart, lineStop], ints for line #s found from start and end phrases.

  • list – segments, list of lines read from file.

  • All lists can contain multiple entries, if more than one segment matches the search criteria.

epsproc.genKQSterms(Kmin=0, Kmax=2, mFlag=True)[source]
epsproc.genKQStermsFromTensors(EPR, AKQS, uniqueFlag=True, phaseConvention='S')[source]

Generate all QNs for \(\Delta_{L,M}(K,Q,S)\) from existing tensors (Xarrays) \(E_{P,R}\) and \(A^K_{Q,S}\).

Cf. epsproc.geomFunc.genllpMatE(), code adapted from there.

Parameters
  • matE (Xarray) – Xarray containing matrix elements, with QNs (l,m), as created by readMatEle()

  • uniqueFlag (bool, default = True) – Check for duplicates and remove (can occur with some forms of matrix elements).

  • mFlag (bool, optional, default = True) – m, mp take all passed values if mFlag=True, or =0 only if mFlag=False

  • phaseConvention (optional, str, default = 'S') – Set phase conventions with epsproc.geomCalc.setPhaseConventions(). To use preset phase conventions, pass existing dictionary. If matE.attrs[‘phaseCons’] is already set, this will be used instead of passed args.

Returns

  • QNs1, QNs2 (two 2D np.arrays) – Values take all allowed combinations [‘P’,’K’,’L’,’R’,’Q’,’M’] and [‘P’,’K’,’L’,’Rp’,’S’,’S-Rp’] from supplied matE. Note phase conventions not applied to QN lists as yet.

  • To do

  • —–

  • - Implement output options (see dev. function w3jTable).

epsproc.genLM(Lmax, allM=True)[source]

Return array of (L,M) up to supplied Lmax

If allM=False only M=0 terms will be set.

TODO: add return type options, include conversion to SHtools types.

epsproc.genllL(Lmin=0, Lmax=10, mFlag=True)[source]

Generate quantum numbers for angular momentum contractions (l, lp, L)

Parameters
  • Lmin (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • Lmax (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • mFlag (bool, optional, default = True) – m, mp take all values -l…+l if mFlag=True, or =0 only if mFlag=False

Returns

QNs – Values take all allowed combinations [‘l’,’lp’,’L’,’m’,’mp’,’M’] up to l=lp=Lmax, one set per row.

Return type

2D np.array

Examples

>>> # Calculate up to Lmax = 2
>>> QNs = genllL(Lmax=2)
>>> # Use with w3jTable function to calculate Wigner 3j terms
>>> w3j = w3jTable(QNs = QNs)

To do

  • Implement output options (see dev. function w3jTable).

epsproc.genllLList(Llist, uniqueFlag=True, mFlag=True)[source]

Generate quantum numbers for angular momentum contractions (l, lp, L) from a passed list, (m, mp, M)=0 or all allowed terms.

Parameters
  • Llist (list) – Values [l, lp, L] to use for calculations. Note this needs to be 2D array in current form of function, i.e. defined as np.array([[L1,L2,L3],…])

  • uniqueFlag (bool, optional, default = True) – Drop duplicate [l,lp,L] sets from list.

  • mFlag (bool, optional, default = True) – m, mp take all values -l…+l if mFlag=True, or =0 only if mFlag=False

Returns

QNs – Values take all allowed combinations [‘l’,’lp’,’L’,’m’,’mp’,’M’] up to l=lp=Lmax, one set per row.

Return type

2D np.array

Examples

>>> # Set from an array
>>> QNs = genllLList(np.array([[1,1,2],[1,3,2],[1,1,2]]), mFlag = True)
>>> # Use with w3jTable function to calculate Wigner 3j terms
>>> w3j = w3jTable(QNs = QNs)

To do

  • Implement output options (see dev. function w3jTable).

epsproc.genllpMatE(matE, uniqueFlag=True, mFlag=True, phaseConvention='S')[source]

Generate quantum numbers for angular momentum contractions (l, lp, L, m, mp, M) from sets of matrix elements.

Parameters
  • matE (Xarray) – Xarray containing matrix elements, with QNs (l,m), as created by readMatEle()

  • uniqueFlag (bool, default = True) – Check for duplicates and remove (can occur with some forms of matrix elements).

  • mFlag (bool, optional, default = True) – m, mp take all passed values if mFlag=True, or =0 only if mFlag=False

  • phaseConvention (optional, str, default = 'S') – Set phase conventions with epsproc.geomCalc.setPhaseConventions(). To use preset phase conventions, pass existing dictionary. If matE.attrs[‘phaseCons’] is already set, this will be used instead of passed args.

Returns

QNs – Values take all allowed combinations [‘l’,’lp’,’L’,’m’,’mp’,’M’] from supplied matE

Return type

2D np.array

Examples

>>> # From demo matrix elements
>>> dataFile = os.path.join(dataPath, 'n2_3sg_0.1-50.1eV_A2.inp.out')  # Set for sample N2 data for testing
>>> # Scan data file
>>> dataSet = ep.readMatEle(fileIn = dataFile)
>>> QNs = genllpMatE(dataSet[0])
>>> # Use with w3jTable function to calculate Wigner 3j terms
>>> w3j = w3jTable(QNs = QNs)

To do

  • Implement output options (see dev. function w3jTable).

epsproc.getCroFileParse(fileName, verbose=1)[source]

Parse an ePS file for GetCro/CrossSection segments.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool or int, optional) – If true, print segment info.

Returns

  • list – [lineStart, lineStop], ints for line #s found from start and end phrases.

  • list – dumpSegs, list of lines read from file.

  • Lists contain entries for each dumpIdy segment found in the file.

epsproc.getCroSegParse(dumpSeg)[source]

Extract values from GetCro/CrossSection file segments.

Parameters

dumpSeg (list) – One CrossSection segment, from dumpSegs[], as returned by epsproc.IO.getCroFileParse()

Returns

  • np.array – CrossSections, table of results vs. energy.

  • list – attribs, list [Label, value, units]

Notes

Currently this is a bit messy, and relies on fixed CrossSection output format. No error checking as yet. Not yet reading all attribs.

Example

>>> XS, attribs = getCroSegParse(dumpSegs[0])
epsproc.getCroSegsParseX(dumpSegs, symSegs, ekeList)[source]

Extract data from ePS getCro/CrossSecion segments into usable form.

Parameters

dumpSegs (list) – Set of dumpIdy segments, i.e. dumpSegs, as returned by epsproc.IO.getCroFileParse()

Returns

  • xr.array – Xarray data array, containing cross sections. Dimensions (Eke, theta)

  • int – Number of blank segments found. (CURRENTLY not implemented.)

Example

>>> data = getCroSegsParseX(dumpSegs)

Notes

A rather cut-down version of epsproc.IO.dumpIdySegsParseX(), no error checking currently implemented.

epsproc.getFiles(fileIn=None, fileBase=None, fType='.out', verbose=True)[source]

Scan dir for ePS (or other) file(s) and return results as a list. File endings specified by fType, default .out. Files are checked for existence with

Parameters
  • fileIn (str, list of strs, optional, default = None) – File(s) to read (file in working dir, or full path). Defaults to current working dir if only a file name is supplied. If None, fileBase dir will be scanned for files. If a list, items will be tested for validity. For consistent results, pass raw strings, e.g. fileIn = r"C:\share\code\ePSproc\python_dev\no2_demo_ePS.out"

  • fileBase (str, optional.) – Dir to scan for files. Currently only accepts a single dir. Defaults to current working dir if no other parameters are passed.

  • fType (str, optional) – File ending for ePS output files, default ‘.out’

  • verbose (bool, optional) – Print output details, default True.

Returns

List of Xarray data arrays, containing matrix elements etc. from each file scanned.

Return type

list

Note: scans only a single dir, no subdirs, using os.listdir. See also classes.multiJob.ePSmultiJob.scanDirs() for alternative with Glob & subdir checks.

epsproc.headerFileParse(fileName, verbose=True)[source]

Parse an ePS file for header & input job info.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool, default True) – Print job info from file header if true.

Returns

  • jobInfo (dict) – Dictionary generated from job details.

  • TO DO

  • —–

  • - Tidy up methods - maybe with parseDigits?

  • - Tidy up dict output.

epsproc.jobSummary(jobInfo=None, molInfo=None, tolConv=0.01)[source]

Print some jobInfo stuff & plot molecular structure. (Currently very basic.)

Parameters
  • jobInfo (dict, default = None) – Dictionary of job data, as generated by :py:function:`epsproc.IO.headerFileParse()` from source ePS output file.

  • molInfo (dict, default = None) – Dictionary of molecule data, as generated by epsproc.IO.molInfoParse() from source ePS output file.

  • tolConv (float, default = 1e-2) – Used to check for convergence in ExpOrb outputs, which defines single-center expansion of orbitals.

Returns

  • JobInfo (list)

  • orbInfo (dict) – Properties of ionizing orbital, as determined from (jobInfo, molInfo).

History

20/09/20 v2 Added orbInfo dict, and use this to hold all orbital related outputs for return. May break old codes (pre v1.2.6-dev).

Moved orbInfo to a separate function.

epsproc.lmPlot(data, pType='a', thres=0.01, thresType='abs', SFflag=True, logFlag=False, eulerGroup=True, selDims=None, sumDims=None, plotDims=None, squeeze=False, fillna=False, xDim='Eke', backend='sns', cmap=None, figsize=None, verbose=False, mMax=10, titleString=None, titleDetails=True, labelRound=3, labelCols=[2, 2], dimMaps={}, lDimLabel=None, mDimLabel=None, catLegend=True)[source]

Plotting routine for ePS matrix elements & BLMs.

First pass - based on new codes + util functions from sphPlot.py, and matE sorting codes.

Parameters
  • data (Xarray, data to plot.) – Should work for any Xarray, but optimised for dataTypes: - matE, matrix elements - BLM paramters - ADMs

  • pType (char, optional, default 'a' (abs values)) – Set (data) type to plot. See plotTypeSelector().

  • thres (float, optional, default 1e-2) – Value used for thresholding results, only values > thres will be included in the plot. Either abs or relative (%) value, according to thresType setting.

  • thresType (str, optional, default = 'abs') – Set to ‘abs’ or ‘pc’ for absolute threshold, or relative value (%age of max value in dataset)

  • SFflag (bool, optional, default = True) – For dataType = matE: Multiply by E-dependent scale factor. For dataType = BLM: Multiply by cross-section (XS) (I.e. if False, normalised BLMs are plotted, with B00 = 1)

  • logFlag (bool, optional, default = False) – Plot values on log10 scale.

  • eulerGroup (bool, optional, default = True) – Group Euler angles by set and use labels (currenly a bit flakey…)

  • selDims (dict, optional, default = {'Type':'L'}) – Dimensions to select from input Xarray.

  • sumDims (tuple, optional, default = None) – Dimensions to sum over from the input Xarray.

  • plotDims (tuple, optional, default = None) – Dimensions to stack for plotting, also controls order of stacking (hence sorting and plotting). Default case plots all dims not included in (xDim, sumDims, selDims). E.g. for matrix elements this will be (‘l’,’m’,’mu’,’Cont’,’Targ’,’Total’,’it’,’Type’) TO DO: auto generation for different dataType, also based on selDims and sumDims selections. 11/03/20: partially fixed, now add any missing dims to plot automatically. NOTE: this currently fails for ‘Labels’ case, not sure why.

  • squeeze (bool, optional, default = True) – Drop singleton dimensions from plot.

  • fillna (bool, optional, default = False) – Fill NaN values with 0 if True.

  • xDim (str, optional, default = 'Eke') – Dimension to use for x-axis, also used for thresholding. Default plots (Eke, LM) surfaces. NOTE: if xDim is a MultiIndex, pass as a dictionary mapping, otherwise it may be unstacked during data prep. E.g. for plotting stacked (L,M), set xDim = {‘LM’:[‘L’,’M’]} NOTE: this is currently also passed to matEleSelector(), so stacked dims must exist in inital Xarray. THIS SHOULD BE FIXED, since it’s a pain.

  • backend (str, optional, default = 'sns') – Plotter to use. Default is ‘sns’ for Seaborn clustermap plot. Set to ‘xr’ for Xarray internal plotting (not all passed args will be used in this case). May be switched according to plot type in future…

  • cmap (str, optional, default = None) – Cmap option to pass to sns clustermap plot.

  • figsize (tuple, optional, default None) – Tuple for Seaborn figure size (ratio), e.g. figsize = (15,5). Useful for setting a long axis explicitly in cases with large dimensional disparity. Default results in a square (ish) aspect.

  • verbose (bool, optional, default False) – Print debug info.

  • mMax (int, optional, default = 10) – Used as default for m colour mapping, in cases where coords are not parsed. TODO: fix this, it’s ugly.

  • titleString (str, optional, default = None) – Set main title for plot. If not set, filename or data.attr[‘jobLabel’] will be used if present. Some plot settings labels will also be appended to this.

  • titleDetails (bool, optional, default = True) – If True add plot details to title string.

  • labelRound (int, optional, default = 3) – Round to N d.p. for plot labels (legend text). This is also passed to epsproc.multiDimXrToPD() for column labels (for E, t or Euler dims only)

  • labelCols (list, optional, default = [2,2]) – Number of columns for 1st and 2nd key maps in SNS Clustermap.

  • dimMaps (dict, optional, default = {}) –

    Pass to override default dim mappings lDims and mDims. Default cases:

    lDims = list(set([‘l’, ‘K’, ‘lp’, ‘L’])-set(xDim)) mDims = [‘m’, ‘mp’, ‘mu’, ‘mup’, ‘Q’, ‘S’, ‘mu0’, ‘Lambda’, ‘M’]

    TODO: Should update to allow for arb arg passing to function + general unpacking?

  • lDimLabel (str or list, optional, default = None) – Dim(s) to use for l-dimension labels. - Set dim name to use specific dim. - Set ‘all’ to use all dims (old behaviour) - First found dim from lDims will be used if None.

  • mDimLabel (str or list, optional, default = None) – Dim(s) to use for m-dimension labels. - Set dim name to use specific dim. - Set ‘all’ to use all dims (old behaviour) - First found dim from mDims will be used if None.

  • catLegend (bool, optional, default = True) – Include 2nd “Categories” legend entries? 08/08/22 Added option, but should also do this automatically for empty cases - needs dim handling consolidation first!

Returns

  • daPlot (Xarray) – Data subset as plotted.

  • legendList (list) – Labels & colour maps

  • g (figure object)

Notes

To do

  • Improved dim handling, maybe use epsproc.util.matEdimList() (and related functions) to avoid hard-coding multiple cases here.
    • Partially implemented in dimMap.

    • Currently throws an error for None in symmetry plotting cases for unrecognised dataType. TO FIX!
      • 14/01/21 fixed, but possibly introduced other issues!

    • Consider consolidating mapping styles further, and keeping a global list of parent and child dims? Or otherwise pulling from stacked Xarray?

  • Improved handling of sets of polarization geometries (angles).

  • CONSOLIDATE stacked/unstacked dim handling. At the moment some functions use stacked, some unstacked, which is leading to annoying bugs.

  • Should update to allow for arb arg passing to function + general unpacking? 08/08/22: added optional dimMaps, fully manual at the moment, only for lDims and mDims currently.

History

  • 08/08/22 - Added some additional optional args dimMaps and labelCols for plot formatting.

  • 14/01/21 - Fixed bug in handling nested dims for unknown datatypes (hopefully), specifically for Sym dims.

  • 27/02/20 - handling for MultiIndex cases for Wigner 3j type data (function of (l,lp,L,m,mp,M))
    • Fixed issue with conversion to Pandas table - should now handle dims better.

    • Improved multidim dropna() using Xarray version (previously used Pandas version).

    • Added xDim as dict mapping capability.

  • 05/12/19 - debug

  • 29/11/19 - initial version.

Examples

(See https://github.com/phockett/ePSproc/blob/master/epsproc/tests/ePSproc_demo_matE_plotting_Nov2019.ipynb)

epsproc.lmSymSummary(data)[source]

Display summary info data tables.

Works nicely in a notebook cell, with Pandas formatted table… but not from function?

For a more sophisticated Pandas conversion, see epsproc.util.conversion.multiDimXrToPD()

epsproc.matEdimList(sType='stacked')[source]

Return standard list of dimensions for matrix elements.

Parameters

sType (string, optional, default = 'stacked') – Selected ‘stacked’ or ‘unstacked’ dimensions. Set ‘sDict’ to return a dictionary of unstacked <> stacked dims mappings for use with xr.stack({dim mapping}).

Returns

list

Return type

set of dimension labels.

epsproc.matEleGroupDim(data, dimGroups=[3, 4, 2])[source]

Group ePS matrix elements by redundant labels.

Default is to group by [‘ip’, ‘it’, ‘mu’] terms, all have only a few values.

TODO: better ways to do this? Shoud be possible at Xarray level.

Parameters

data (list) – Sections from dumpIdy segment, as created in dumpIdySegsParseX() Ordering is [labels, matElements, attribs].

epsproc.matEleGroupDimX(daIn)[source]

Group ePS matrix elements by redundant labels (Xarray version).

Group by [‘ip’, ‘it’, ‘mu’] terms, all have only a few values. Rename ‘ip’:1,2 as ‘Type’:’L’,’V’

TODO: better ways to do this? Via Stack/Unstack? http://xarray.pydata.org/en/stable/api.html#id16 See also tests in funcTests_210819.py for more versions/tests.

Parameters

data (Xarray) – Data array with matrix elements to be split and recombined by dims.

Returns

  • data (Xarray) – Data array with reordered matrix elements (dimensions).

  • NOTE Oct 2022 (this is currently failing at ‘it’ restack in XR >2022.3, likely due to change in selectors?)

  • See https (//github.com/phockett/ePSproc/issues/64)

  • Should rewrite in any case!

epsproc.matEleGroupDimXnested(da)[source]

Group ePS matrix elements by redundant labels (Xarray version).

Group by [‘ip’, ‘it’, ‘mu’] terms, all have only a few values.

TODO: better ways to do this? See also tests in funcTests_210819.py for more versions/tests.

Parameters

data (Xarray) – Data array with matrix elements to be split and recombined by dims.

epsproc.matEleSelector(da, thres=None, inds=None, dims=None, sq=False, drop=True)[source]

Select & threshold raw matrix elements in an Xarray. Wraps Xarray.sel(), plus some additional options.

See Xarray docs for more: http://xarray.pydata.org/en/stable/user-guide/indexing.html

Parameters
  • da (Xarray) – Set of matrix elements to sub-select

  • thres (float, optional, default None) – Threshold value for abs(matElement), keep only elements > thres. This is element-wise.

  • inds (dict, optional, default None) – Dicitonary of additional selection criteria, in name:value format. These correspond to parameter dimensions in the Xarray structure. E.g. inds = {‘Type’:’L’,’Cont’:’A2’} Slices are also acceptable, e.g. inds = {‘Eke’:slice(1,5,4)}

  • dims (str or list of strs, dimensions to look for max & threshold, default None) – Set for dimension-wise thresholding. If set, this is used instead of element-wise thresholding. List of dimensions, which will be checked vs. threshold for max value, according to abs(dim.max) > threshold This allows for consistent selection of continuous parameters over a dimension, by a threshold.

  • sq (bool, optional, default False) – Squeeze output singleton dimensions.

  • drop (bool, optional, default True) – Passed to da.where() for thresholding, drop coord labels for values below threshold.

Returns

Xarray structure of selected matrix elements. Note that Nans are dropped if possible.

Return type

daOut

Example

>>> daOut = matEleSelector(da, inds = {'Type':'L','Cont':'A2'})

Notes

xr.sel(inds) is used here. For single values xr.sel({name:[value]}) or xr.sel({name:value}) is different! Automatically squeeze out dim in latter case. (Tested on xr v0.15)

E.g., for selecting a single Eke value: da.sel({‘Eke’:[1.1]}) # Keeps Eke dim da.sel({‘Eke’:1.1}) # Drops Eke to non-dimension coord. da.sel({‘Eke’:1.1}, drop=True) # Drops Eke completely da.sel({‘Eke’:[1.1]}, drop=True) # Keeps Eke da.sel({‘Eke’:[1.1]}, drop=True).squeeze() # Drops Eke to non-dim coord

epsproc.mfblm(daIn, selDims={'Type': 'L'}, eAngs=[0, 0, 0], thres=0.0001, sumDims=('l', 'm', 'mu', 'Cont', 'Targ', 'Total', 'it'), SFflag=True, verbose=1)[source]

Calculate MFBLMs for a range of (E, sym) cases. Default is to calculated for all symmetries at each energy.

Parameters
  • da (Xarray) – Contains matrix elements to use for calculation. Matrix elements will be sorted by energy and BLMs calculated for each set.

  • selDims (dict, optional, default = {'Type':'L'}) – Additional sub-selection to be applied to matrix elements before BLM calculation. Default selects just length gauge results.

  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Single set of Euler angles defining polarization geometry.

  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.

  • sumDims (tuple, optional, default = ('l','m','mu','Cont','Targ','Total','it')) – Defines which terms are summed over (coherent) in the MFBLM calculation. (These are used to flatten the Xarray before calculation.) Default includes sum over (l,m), symmetries and degeneracies (but not energies).

  • SFflag (bool, default = True) – Normalise by scale factor to give X-sections (B00) in Mb

  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.

    • 1: Print basic info.

    • 2: Print intermediate C parameter array to terminal when running.

Returns

  • Xarray – Calculation results BLM, dims (Euler, Eke, l,m). Some global attributes are also appended.

  • Limitations

  • ———–

  • Currently set to loop calcualtions over energy only, and all symmetries.

  • Pass single {‘Cont’ (‘sym’} to calculated for only one symmetry group.)

  • TODO (In future this will be more elegant.)

  • TODO (Setting selDims in output structure needs more thought for netCDF save compatibility.)

epsproc.mfblmEuler(da, selDims={'Type': 'L'}, eAngs=[0, 0, 0], thres=0.0001, sumDims=('l', 'm', 'mu', 'Cont', 'Targ', 'Total', 'it'), SFflag=True, verbose=1)[source]

Wrapper for epsproc.MFBLM.mfblm() for a set of Euler angles. All other parameters are simply passed to mfblm(). Calculate MFBLMs for a range of (E, sym) cases. Default is to calculated for all symmetries at each energy.

Parameters
  • da (Xarray) – Contains matrix elements to use for calculation. Matrix elements will be sorted by energy and BLMs calculated for each set.

  • selDims (dict, optional, default = {'Type':'L'}) – Additional sub-selection to be applied to matrix elements before BLM calculation. Default selects just length gauge results.

  • eAngs ([phi,theta,chi], optional, default = [0,0,0]) – Set of Euler angles defining polarization geometry. List or np.array, dims(N, 3). Basic Xarray support also in place, as generated by epsproc.sphCalc.setPolGeoms()

  • thres (float, optional, default = 1e-4) – Threshold value for testing significance of terms. Terms < thres will be dropped.

  • sumDims (tuple, optional, default = ('l','m','mu','Cont','Targ','Total','it')) – Defines which terms are summed over (coherent) in the MFBLM calculation. (These are used to flatten the Xarray before calculation.) Default includes sum over (l,m), symmetries and degeneracies (but not energies).

  • SFflag (bool, default = True) – Normalise by scale factor to give X-sections (B00) in Mb

  • verbose (int, optional, default 1) –

    Verbosity level:

    • 0: Silent run.

    • 1: Print basic info.

    • 2: Print intermediate C parameter array to terminal when running.

Returns

  • Xarray – Calculation results BLM, dims (Euler(P,T,C), Eke, BLM(l,m)) - as per epsproc.util.BLMdimList(). Some global attributes are also appended.

  • Limitations

  • ———–

  • Currently set to loop calcualtions over energy only, and all symmetries.

  • Pass single {‘Cont’ (‘sym’} to calculated for only one symmetry group.)

  • TODO (In future this will be more elegant.)

  • TODO (Setting selDims in output structure needs more thought for netCDF save compatibility.)

epsproc.mfblmXprod(matEin, QNs=None, EPRX=None, p=[0], BLMtable=None, BLMtableResort=None, lambdaTerm=None, RX=None, eulerAngs=None, polProd=None, thres=0.01, thresDims='Eke', selDims={'Type': 'L', 'it': 1}, sumDims=['mu', 'mup', 'l', 'lp', 'm', 'mp'], sumDimsPol=['P', 'R', 'Rp', 'p'], symSum=True, SFflag=False, squeeze=False, phaseConvention='E', basisReturn='BLM', verbose=0, **kwargs)[source]

Implement \(\beta_{LM}^{MF}\) calculation as product of tensors.

\[\begin{split}\begin{eqnarray} \beta_{L,-M}^{\mu_{i},\mu_{f}} & = & \sum_{l,m,\mu}\sum_{l',m',\mu'}(-1)^{(\mu'-\mu_{0})}{B_{L,-M}}\nonumber \\ & \times & \sum_{P,R',R}{E_{P-R}(\hat{e})\Lambda_{R',R}(R_{\hat{n}})}I_{l,m,\mu}^{p_{i}\mu_{i},p_{f}\mu_{f}}(E)I_{l',m',\mu'}^{p_{i}\mu_{i},p_{f}\mu_{f}*}(E) \end{eqnarray}\end{split}\]

Where each component is defined by fns. in :py:module:`epsproc.geomFunc.geomCalc` module.

12/08/22 Updating from afblmXprod routine for use with fitting functions.

Added ProductBasis return as per afblmGeom case, for use in fitting. Added **kwargs, unused but allows for arb basis dict unpack and passing from other functions. May want to pipe back to Full basis return however. Updated docs as per afblmXprod. TODO: still needs a tidy up and update, including BLM renorm options, again see afblmGeom code.

16/03/20 In progress!

Dev code:

http://localhost:8888/lab/tree/dev/ePSproc/geometric_method_dev_Betas_090320.ipynb D:codeePSprocpython_devePSproc_MFBLM_Numba_dev_tests_120220.PY

Parameters
  • matE (Xarray) – Xarray containing matrix elements, with QNs (l,m), as created by readMatEle()

  • settings (*** Optional calculation) –

  • selDims (dict, default = {'it':1, 'Type':'L'}) – Selection parameters for calculations, may be be checked and appened herein.

  • sumDims (list, default = ['mu', 'mup', 'l','lp','m','mp']) – Main summation dims, will be checked herein.

  • sumDimsPol (list, default = ['P','R','Rp','p']) – Additional polarization summation dims.

  • symSum (bool, default = True) – Sum over symmetries sets (=={Cont, Targ, Total}) if true.

  • degenDrop (bool) – Flag to set dropping of degenerate components. NOT IMPLEMENTED FOR MF CASE - see afblmXprod.

  • thres (float, default = 1e-2) – Set threshold value, used to set input matrix elements and again for outputs.

  • thresDims (str, default = 'Eke') – Set threshold dimension (set to be contiguous).

  • verbose (bool or int) – Print output?

*** Optional renormalisation settings (mainly for testing only)

SFflagbool, default = False

Multiply input matrix elements by complex scale-factor if true.

SFflagRenormbool, default = False

Renorm output BLMs by complex scale-factor if true. NOT IMPLEMENTED FOR MF CASE - see afblmXprod.

BLMRenormint, default = 1

Set different BLM renorm conventions. If 1 renorm by B00. See code for further details. NOT IMPLEMENTED FOR MF CASE - see afblmXprod.

squeezebool, default = False

Squeeze output array after thresholding? Note: this may cause dim issues if True.

*** Optional input data/basis functions (mainly for fitting routine use)

QNsnp.array, optional, default = None

List of QNs as generated by genllpMatE(). Will be generated if not passed.

EPRXXarray, optional, default = None

E-field parameters, as generated by EPR(). Defaults to normalised/unity field, pol = p (below).

plist or array, optional, default = [0]

Specify polarization terms p. Possibly currently only valid for p=0, TBC See https://epsproc.readthedocs.io/en/latest/methods/geometric_method_dev_260220_090420_tidy.html#E_{P,R}-tensor

BLMtable, BLMtableResortXarrays, optional, default = None

Beta calculation parameters, as defined by geomCalc.betaTerm(). BLMtableResort includes phase settings & param renaming as herein.

lambdaTermXarray, optional, default = None

Lambda term parameters, as defined by geomCalc.MFproj()

RXXarray, optional, default = None

Polarization geometries as defined by epsproc.sphCalc.setPolGeoms(). If not set, defaults are used by epsproc.geomFunc.geomCalc.MFproj(). If not set, but Euler angles are set, then these will be used.

eulerAngslist or np.array of Euler angles (p(hi), t(heta), c(hi)), optional.

Alternative definition for polarization geometries, as used by epsproc.sphCalc.setPolGeoms(). List or array [p,t,c…], shape (Nx3). List or array including set labels, [label,p,t,c…], shape (Nx4)

polProdXarray, optional, default = None

Polarization tensor as defined by EPRXresort * lambdaTermResort

phaseConventionoptional, str, default = ‘E’

Set phase conventions with epsproc.geomCalc.setPhaseConventions(). To use preset phase conventions, pass existing dictionary.

basisReturnoptional, str, default = “BLM”
  • ‘BLM’ return Xarray of results only.

  • ‘Full’ return Xarray of results + basis set dictionary as set during the run.

  • ‘Product’, as full, but minimal basis set with products only.

  • ‘Results’ or ‘Legacy’ direct return of various calc. results Xarrays.

**kwargs, unused but allows for arb basis dict unpack and passing from other functions.

Returns

  • Xarray – Set of AFBLM calculation results

  • dict – Dictionary of basis functions, only if basisReturn != ‘BLM’ (see basisReturn paramter notes).

Notes

Cross-section outputs currently defined as XS = direct MF calculation output.

Optionally set SFflag = True to multiply by (complex) scale-factor.

OTHER RENORM options not implemented as yet, see afblmXprod for details.

epsproc.mfpad(dataIn, thres=0.01, inds={'Type': 'L', 'it': 1}, res=50, R=None, p=0)[source]
Parameters
  • dataIn (Xarray) – Contains set(s) of matrix elements to use, as output by epsproc.readMatEle().

  • thres (float, optional, default 1e-2) – Threshold value for matrix elements to use in calculation.

  • ind (dictionary, optional.) – Used for sub-selection of matrix elements from Xarrays. Default set for length gauage, single it component only, inds = {‘Type’:’L’,’it’:’1’}.

  • res (int, optional, default 50) – Resolution for output (theta,phi) grids.

  • R (list or Xarray of Euler angles or quaternions, optional.) – Define LF > MF polarization geometry/rotations. For default case (R = None), calls epsproc.setPolGeoms(), where 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, [0 pi/2 0] for x-pol, [pi/2 pi/2 0] for y-pol.

  • p (int, optional.) – Defines LF polarization state, p = -1…1, default p = 0 (linearly pol light along z-axis). TODO: add summation over p for multiple pol states in LF.

Returns

  • Ta – Xarray (theta, phi, E, Sym) of MFPADs, summed over (l,m)

  • Tlm – Xarray (theta, phi, E, Sym, lm) of MFPAD components, expanded over all (l,m)

epsproc.molInfoParse(fileName, verbose=True)[source]

Parse an ePS file for input molecule info.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool, default True) – Print job info from file header if true.

Returns

molInfo – Dictionary with atom & orbital details.

Return type

dict

Notes

Only tested for Molden input (MoldenCnv2006).

epsproc.molPlot(molInfo)[source]

Basic 3D scatter plot from molInfo data.

epsproc.multiDimXrToPD(da, colDims=None, rowDims=None, thres=None, squeeze=True, dropna=True, fillna=False, colRound=2, verbose=False)[source]

Convert multidim Xarray to stacked Pandas 2D array, (rowDims, colDims)

Parameters
  • da (Xarray) – Array for conversion.

  • colDims (list of dims for columns, default = None) –

  • rowDims (specifiy both colDims and) –

  • NOTE (if xDim is a MultiIndex, pass as a dictionary mapping, otherwise it may be unstacked during data prep.) –

  • ordering (For full control over dim stack) –

  • rowDims

  • NOTE

  • (L (E.g. for plotting stacked) –

  • M) (['L','M']}) –

  • {'LM' (set xDim =) –

  • thres (float, optional, default = None) – Threshold values in output (pd table only) TODO: generalise this and use matEleSelector() for input?

  • squeeze (bool, optional, default = True) – Drop singleton dimensions.

  • dropna (bool, optional, default = True) – Drop all NaN dimensions from output pd data frame (columnwise and rowise).

  • fillna (bool, optional, default = False) – Fill any NaN values with 0.0. Useful for plotting/making data contiguous.

  • colRound (int, optional, default = True) – Round column values to colRound dp. Only applied for Eke, Ehv, Euler or t dimensions.

Returns

  • daRestackpd (pandas data frame (2D) with sorted data.)

  • daRestack (Xarray with restacked data.)

Method

Restack Xarray by specified dims, including basic dims checking, then use da.to_pandas().

12/03/20 Function adapted from lmPlot() code.

Note

This might casue epsproc.lmPlot() to fail for singleton x-dimensions if squeeze = True. TO do: add work-around, see lines 114-122.

epsproc.parseLineDigits(testLine)[source]

Use regular expressions to extract digits from a string. https://stackoverflow.com/questions/4289331/how-to-extract-numbers-from-a-string-in-python

epsproc.plotTypeSelector(dataPlot=None, 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. Note: this parameter is required, unless returnDict = True.

  • 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.readMatEle(fileIn=None, fileBase=None, fType='.out', recordType='DumpIdy', verbose=1, stackE=True, stackDim='Eke')[source]

Read ePS file(s) and return results as Xarray data structures. File endings specified by fType, default *.out.

Parameters
  • fileIn (str, list of strs, optional.) – File(s) to read (file in working dir, or full path). Defaults to current working dir if only a file name is supplied. For consistent results, pass raw strings, e.g. fileIn = r"C:\share\code\ePSproc\python_dev\no2_demo_ePS.out"

  • fileBase (str, optional.) – Dir to scan for files. Currently only accepts a single dir. Defaults to current working dir if no other parameters are passed.

  • fType (str, optional) – File ending for ePS output files, default ‘.out’

  • recordType (str, optional, default 'DumpIdy') – Type of record to scan for, currently set for ‘DumpIdy’, ‘EDCS’ or ‘CrossSection’. For a full list of descriptions, types and sources, run: >>> epsproc.util.dataTypesList()

  • verbose (int, optional, default = 1) – Level of verbosity in output. - 0 no printed output - 1 print summary info only - 2 print detailed info

  • stackE (bool, optional, default = True) – Identify and stack multi-part jobs to single array (by E) if True.

  • stackDim (bool, optional, default = 'Eke') – Dim to stack. Note if stackE=True, any dim can be set here (not just E).

Returns

  • list – List of Xarray data arrays, containing matrix elements etc. from each file scanned.

  • To do

  • —–

  • - Change to pathlib paths.

  • - Implement outputType options…?

  • 20/07/22 Added dict return for dumpIdy methods for troubleshooting if getCroSegsParseX fails. – Note this only returns if stackE = False is also set.

  • 13/10/20 Adapted to use grouped lists for multi-file jobs, should be back-compatible if stackE = False set.

Examples

>>> dataSet = readMatEle()  # Scan current dir
>>> fileIn = r'C:\share\code\ePSproc\python_dev\no2_demo_ePS.out'
>>> dataSet = readMatEle(fileIn)  # Scan single file
>>> dataSet = readMatEle(fileBase = r'C:\share\code\ePSproc\python_dev') # Scan dir

Note

  • Files are scanned for matrix element output flagged by “DumpIdy” headers.

  • Each segment found is parsed for attributes and data (set of matrix elements).

  • Matrix elements and attributes are combined and output as an Xarray array.

epsproc.readOrb3D(fileIn=None, fileBase=None, fType='_Orb.dat', verbose=True)[source]

Read ePS 3D data file(s) and return results. File endings specified by fType, default *_Orb.dat.

fileInstr, list of strs, optional.

File(s) to read (file in working dir, or full path). Defaults to current working dir if only a file name is supplied. For consistent results, pass raw strings, e.g. fileIn = r”C:sharecodeePSprocpython_dev

o2_demo_ePS.out”

fileBasestr, optional.

Dir to scan for files. Currently only accepts a single dir. Defaults to current working dir if no other parameters are passed.

fTypestr, optional

File ending for ePS output files, default ‘_Orb.dat’

verbosebool, optional

Print output details, default True.

list

List of data arrays, containing matrix elements etc. from each file scanned.

# TODO: Change output to Xarray?

>>> dataSet = readOrb3D()  # Scan current dir
>>> fileIn = r'C:\share\code\ePSproc\python_dev\DABCOSA2PPCA2PP_10.5eV_Orb.dat'
>>> dataSet = readOrb3D(fileIn)  # Scan single file
>>> dataSet = readOrb3D(fileBase = r'C:\share\code\ePSproc\python_dev') # Scan dir
epsproc.readOrbCoords(f, headerLines)[source]
epsproc.readOrbData(f, headerLines)[source]
epsproc.readOrbElements(f, n)[source]
epsproc.readOrbHeader(f)[source]
epsproc.readXarray(fileName, filePath=None, engine='h5netcdf', forceComplex=False, forceArray=True, **kwargs)[source]

Wrapper for backend Xarray file readers.

epsproc.remapllpL(dataIn, QNs, form='dict', method='sel', dlist=['l', 'lp', 'L', 'm', 'mp', 'M'], verbose=0)[source]

Remap Wigner 3j table, with QNs (l,lp,L,m,mp,M) > tensor forms.

Tensors are stored by (l,lp,L) triples, with corresponding arrays [m,mp,M], as a dictionary, or Xarray dataarray or dataset.

Parameters
  • dataIn (np.array) – Array of data values corresponding to rows (coords) in QNs.

  • QNs (np.array) – List of QNs [l, lp, L, m, mp, M] to compute 3j terms for. If not supplied all allowed terms for {Lmin, Lmax} will be set. (If supplied, values for Lmin, Lmax and mFlag are not used.)

  • form (str, optional, default = 'dict') – Type of output structure. - ‘dict’ : dictionary with keys (l,lp,L), coordinate tables - ‘3d’ : dictionary with keys (l,lp,L), 3D arrays indexed by [l+m, lp+mp, L+M]; this case also sets (0,0,0) term as ‘w3j0’. - ‘xdaLM’ : Xarray dataarray, with stacked dims [‘lSet’,’mSet’] - ‘xds’ : Xarray dataset, with one array per (l,lp,L) - ‘xdalist’ : List of Xarray dataarrays, one per (l,lp,L)

  • method (str, optional, default = 'sel') – Method for selection from input array. - ‘sel’ : use full selection routine, epsproc.selQNsRow(), should be most robust. - ‘n’ : sort based on np.unique reindexing. Should be faster, but possibly not robust…?

  • dlist (list, optional, default = ['l','lp','L','m','mp','M']) – Full list of dimension names to use/set.

Returns

w3j – Wigner3j(l,lp,L,m,mp,M) values corresponding to rows (coords) in QNs, type according to input form. Data structure sorted by (l,lp,L) triples.

Return type

Xarray, dictionary

epsproc.scatEngFileParse(fileName, verbose=1)[source]

Parse an ePS file for ScatEng list.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool or int, optional) – If true, print segment info.

Returns

  • list – ekeList, np array of energies set in the ePS file.

  • Lists contain entries for each dumpIdy segment found in the file.

epsproc.selQNsRow(QNs, QNmask, fields=None, verbose=1)[source]

Basic routine for selecting/filtering valid rows from an input array of QNs.

This is similar to methods used in Matlab version, but likely rather slow.

Parameters
  • QNs (np.array) – Values to be filtered.

  • QNmask (list or np.array) – Values to be matched. Must be same dims as QNs, unless fields is also set. If set to None field will be skipped (i.e. match all values)

  • fields (list or np.array, optional, default = None) – Which fields to match in QNs.

  • verbose (int, optional, default = 1) – Print info to terminal.

Returns

  • QNs[boolMask] (np.array of selected values.)

  • boolMask (np.array corresponding to selected values.)

Examples

>>> QNsMask, mask = selQNsRow(QNs,[None, None, None, 0, 0, 0], verbose = False)
>>> # With fields
>>> QNsMask, mask = selQNsRow(QNs,[0, 0, 0], fields = [3,4,5], verbose = False)
epsproc.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.setPhaseConventions(phaseConvention='S', typeList=False)[source]

Set phase convention/choices for geometric functions.

20/03/20 - first attempt. Aim to centralise all phase choices here to keep things clean and easy to debug/change.

Set as dictionary for each term, to be appended to results Xarray.

Parameters
  • phaseConvention (optional, str, default = 'S') – Set phase conventions: - ‘S’ : Standard derivation. - ‘R’ : Reduced form geometric tensor derivation. - ‘E’ : ePolyScat, may have additional changes in numerics, e.g. conjugate Wigner D. If a dict of phaseConventions is passed they will simply be returned - this is for transparency/consistency over multiple fns which call setPhaseConventions()… although may be an issue in some cases.

  • typeList (optional, bool, default = False) – If true, return list of supported options instead of list of phase choices.

Note

If a dict of phaseConventions is passed they will simply be returned - this is for transparency/consistency over multiple fns which call setPhaseConventions()… although may be an issue in some cases.

epsproc.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(Lmax=2, Lmin=0, res=None, angs=None, XFlag=True, fnType='sph', convention='phys', conj=False, realCSphase=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. 23/09/22 - hacked in ‘complex’ and ‘real’ here too, but TODO: split by kind and backend. 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.

  • realCSphase (bool, optional, default = False) – Apply (-1)^m term in real harmonics calc. Generally will already be included in complex spherical harmonics defn., so can be left as False.

Note that either res OR angs needs to be passed.

Outputs

  • 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
.. math::

begin{equation} Y_{l,m}(theta,phi) = (-1)^msqrt{frac{2l+1}{4pi} frac{(l-m)!}{(l+m)!}}e^{i m phi} P^m_l(cos(theta)) end{equation}

And for real harmonics:
.. math::

begin{equation} begin{aligned} Y_{ell m}&={begin{cases}{sqrt {2}},Im [{Y_{ell }^{|m|}}]&{text{if}}mlt0 \Y_{ell }^{0}&{text{if}}m=0 \{sqrt {2}},Re [{Y_{ell }^{m}}]&{text{if}}mgt0.end{cases}} end{aligned} end{equation}

See https://en.wikipedia.org/wiki/Spherical_harmonics#Real_form.
Demos:
- Complex: https://epsproc.readthedocs.io/en/dev/special_topics/ePSproc_docs_working_with_spherical_harmonics_200922.html
- Real: https://epsproc.readthedocs.io/en/dev/special_topics/ePSproc_docs_working_with_real_harmonics_220922.html
- For scipy backend, realCSphase = False case matches SHtools defn, see https://epsproc.readthedocs.io/en/dev/special_topics/ePSproc_docs_working_with_real_harmonics_220922.html#Converting-coefficients-real-to-complex-form

Example

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

Notes

TODO:

  • additional backends.

  • remove hard-coded dim names for flexibility (or set from ref YLMdimList())

  • remove ‘fnType’ and/or use ‘kind’ instead (SHtool notation)

epsproc.sphFromBLMPlot(BLMXin, res=50, pType='a', plotFlag=False, facetDim=None, backend='mpl', fnType=None, conj=False)[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.

  • conj (bool, optional, default = False) – Use conjugate harmonics.

Returns

  • Xarray – Xarray containing the calcualted surfaces I(theta,phi,…)

  • fig – List of figure handles.

epsproc.sphPlotHV(dataIn)[source]
epsproc.sphPlotMPL(dataPlot, theta, phi, convention='phys', tString=None)[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).

  • theta (np.arrays) – Angles defining spherical polar grid, 2D arrays.

  • phi (np.arrays) – Angles defining spherical polar grid, 2D arrays.

  • convention (str, optional, default = 'phys') – Set spherical polar coord convention, see epsproc.sphToCart().

  • tString (str, optional, default = None) – Text to be used for plot title. This will be appended with other data info, if set. If facetDim is passed here, this will be used to set the label.

Returns

Handle to matplotlib figure.

Return type

fig

epsproc.sphPlotPL(dataPlot, theta, phi, facetDim='Eke', surfMap=None, rc=None, norm='global', convention='phys', plotFlag=True, verbose=False)[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).

  • theta (np.arrays) – Angles defining spherical polar grid, 2D arrays.

  • phi (np.arrays) – Angles defining spherical polar grid, 2D arrays.

  • facetDim (str, default 'Eke') – Dimension to use for faceting (subplots), currently set for single dim only.

  • surfMap (str, optional, default = None) – Additional specification to use for Plotly surfacecolor specification. If None not specified == z surf map. If R == radial value surf map.

  • rc (array, optional, default = None) – If set, use to define layout grid [rows, columns]. If not set, this will be set for nCols = 5

  • norm (str, optional, default = 'global') – Set for plot normalisation. - ‘global’ : use same (x,y,z) limits for all plots. - ‘local’ : auto set (x,y,z) limits for each plot.

  • convention (str, optional, default = 'phys') – Spherical polar coord convention, see epsproc.sphToCart()

  • plotFlag (bool, optional, default = True) – Set plotFlag=False bypass for plotter object return only.

Returns

Handle to figure.

Return type

fig

Notes

For additional dev notes, see: - http://localhost:8889/notebooks/dev/ePSproc/plottingDev/plotly_subplot_tests_051020.ipynb - http://localhost:8889/notebooks/dev/ePSproc/classDev/ePSproc_multijob_class_tests_N2O_011020_Stimpy_PLrenderTESTING_local.ipynb

For Jupyter use: currently (Oct. 2020) working only for Notebook Export to HTML (not nbsphinx), and max of 12 subplots. Possible issues with rendering in Firefox (80.0.1, Oct. 2020). Should be fixable with some effort/testing, see https://plotly.com/python/renderers/

TODO:

For JupyterLab, need additional extensions - see https://plotly.com/python/getting-started/#jupyterlab-support: - conda install -c conda-forge -c plotly jupyter-dash - jupyter labextension install jupyterlab-plotly

In some cases may get partially working installation with, e.g., blank surface plots, or plots via HV only. This usually means JupyterLab needs a restart (and maybe a rebuild). For more see https://plotly.com/python/troubleshooting/

epsproc.sphSumPlotX(dataIn, pType='a', facetDim='Eke', surfMap='R', backend='mpl', convention='phys', titleString=None, plotFlag=True, verbose=True, axisUW=None)[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.

  • surfMap (str, optional, default = None) – Additional specification to use for surface colour specification. NOTE: currently only implemented for Plotly backend (via surfacecolor). If None not specified == z surf map. If R == radial value surf map.

  • 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.

  • convention (str, optional, default = 'phys') – Spherical polar coord convention, see epsproc.sphToCart()

  • titleString (str, optional, default = None) – Additional info to use for plot title.

  • plotFlag (bool, optional, default = True) – Set plotFlag=False bypass for plotter object return only. (Required for Plotly backend only.)

  • axisUW (str, optional, default = None) – Unwrap axis for phase maps, only used if pType = ‘phaseUW’

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.

More sophisticated dim handling is implemented in padPlot(), but for class only.

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

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

Parameters
  • R (np.arrays) – Spherical polar coords \((R,\theta,\phi)\).

  • theta (np.arrays) – Spherical polar coords \((R,\theta,\phi)\).

  • phi (np.arrays) – 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

Definitions

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.

epsproc.stringRepMap(string, replacements, ignore_case=False)[source]

Given a string and a replacement map, it returns the replaced string. :param str string: string to execute replacements on :param dict replacements: replacement dictionary {value to find: value to replace} :param bool ignore_case: whether the match should be case insensitive :rtype: str

CODE from: https://gist.github.com/bgusach/a967e0587d6e01e889fd1d776c5f3729 https://stackoverflow.com/questions/6116978/how-to-replace-multiple-substrings-of-a-string … more or less verbatim.

Thanks to bgusach for the Gist.

epsproc.symFileParse(fileName, verbose=1)[source]

Parse an ePS file for scattering symmetries.

Parameters
  • fileName (str) – File to read (file in working dir, or full path)

  • verbose (bool or int, optional) – If true, print segment info.

Returns

  • list – symSegs, raw lines from the ePS file.

  • Lists contain entries for each ScatSym setting found in file header (job input).

epsproc.symListGen(data)[source]
epsproc.w3jTable(Lmin=0, Lmax=10, QNs=None, mFlag=True, nonzeroFlag=False, form='2d', dlist=['l', 'lp', 'L', 'm', 'mp', 'M'], backend='par', verbose=0)[source]

Calculate/tabulate all wigner 3j terms for a given problem/set of QNs.

\[\begin{split}\begin{equation} \begin{array}{ccc} l & l' & L\\ m & m' & M \end{array} \end{equation}\end{split}\]

Where l, l’ take values Lmin…Lmax (default 0…10). \(\l-lp\<=L<=l+lp\) m, mp take values -l…+l if mFlag=True, or =0 only if mFlag=False

Parameters
  • Lmin (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • Lmax (int, optional, default 0, 10) – Integer values for Lmin and Lmax respectively.

  • QNs (np.array, optional, default = None) – List of QNs [l, lp, L, m, mp, M] to compute 3j terms for. If not supplied all allowed terms for {Lmin, Lmax} will be set. (If supplied, values for Lmin, Lmax and mFlag are not used.) NOTE: some return types will convert QNs to int when constructing Xarray output, unless half-int values present. Functions using epsproc.geomFunc.geomCalc.remapllpL() support half-int values in Xarray.

  • mFlag (bool, optional, default = True) – m, mp take all values -l…+l if mFlag=True, or =0 only if mFlag=False

  • nonzeroFlag (bool, optional, default = False) – Drop null terms before returning values if true.

  • form (string, optional, default = '2d') –

    Defines return format. Options are:
    • 2d, return 2D np.array, rows [l, lp, L, m, mp, M, 3j]

    • xarray, return xarray

      This is nice for easy selection/indexing, but may be problematic for large Lmax if unstacked (essentailly similar to nd case).

    • nd, return ND np.array, dims indexed as [l, lp, L, l+m, lp+mp, L+M], with values 3j.

      This is suitable for direct indexing, but will have a lot of zero entries and may be large.

    • ndsparse, return ND sparse array, dims indexed as [l, lp, L, l+m, lp+mp, L+M], with values 3j.

    Additional options are set via remapllpL(). This additionally sorts values by (l,lp,L) triples, which is useful in some cases.
    • ’dict’ : dictionary with keys (l,lp,L), coordinate tables

    • ’3d’ : dictionary with keys (l,lp,L), 3D arrays indexed by [l+m, lp+mp, L+M]; this case also sets (0,0,0) term as ‘w3j0’.

    • ’xdaLM’ : Xarray dataarray, with stacked dims [‘lSet’,’mSet’]

    • ’xds’ : Xarray dataset, with one array per (l,lp,L)

    • ’xdalist’ : List of Xarray dataarrays, one per (l,lp,L)

  • dlist (list of labels, optional, default ['l','lp','L','m','mp','M']) – Used to label array for Xarray output case.

  • backend (str, optional, default = 'par') – See Implementation note below.

Returns

w3j – Wigner3j(l,lp,L,m,mp,M) values corresponding to rows (coords) in QNs, type according to input form.

Return type

np.array, Xarray, dictionary

Implementation

Currently set to run:
  • ‘vec’: w3jguVecCPU(), which uses sf.Wigner3j on the back-end, with additional vectorisation over supplied QNs via Numba’s @guvectorize.

  • ‘par’: w3jprange(), which uses sf.Wigner3j on the back-end, with parallelization over QNs via Numba’s @njit with a prange loop.

epsproc.w3jguVecCPU(QNs, w3j_QNs)

Wrapper for 3j with vectorization via @numba.guvectorize([“void(int32[:,:], float64[:])”], ‘(n,m)->(n)’, target = ‘parallel’).

Parameters
  • QNs (np.array) – Array of QNs to calculated Wigner 3j terms for, columns [l,lp,L,m,mp,M].

  • w3j_QNs (np.array) – Empty array to hold results (no return from @guvectorize). Create as w3j_QNs = np.zeros(QNs.shape[0])

epsproc.w3jprange(QNs)

Wrapper for 3j with @numba.njit(parallel=True), using prange parallel loop.

In testing (Feb 2020) on an AMD Threadripper 1950X (16 core) this was (usually) fastest case, and beat vectorised version.

Parameters

QNs (np.array) – Array of QNs to calculated Wigner 3j terms for, columns [l,lp,L,m,mp,M].

Returns

w3j_QNs – Array of Wigner 3j results, one per row of input QNs.

Return type

np.array

epsproc.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.

  • only) (Options for setting angles (use one) –

  • 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).

XFlagbool, optional, default True

Flag for output. If true, output is Xarray. If false, np.arrays

dlistlist, optional, default [‘lp’,’mu’,’mu0’]

Labels for Xarray QN dims.

eNameslist, optional, default [‘P’,’T’,’C’]

Labels for Xarray Euler dims.

conjFlagbool, optional, default = False

If true, return complex conjuage values.

sfErrorbool, 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

Outputs

  • 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)
epsproc.writeOrb3Dvtk(dataSet)[source]

Write ePS 3D data file(s) to vtk format. This can be opened in, e.g., Paraview.

Parameters
  • dataSet (list) – List of data arrays, containing matrix elements etc. from each file scanned. Assumes format as output by readOrb3D(), [fileName, headerLines, coords, data]

  • TODO (#) –

Returns

List of output files.

Return type

list

epsproc.writeXarray(dataIn, fileName=None, filePath=None, engine='h5netcdf', forceComplex=False, **kwargs)[source]

Wrapper for backend Xarray file writers.