Method development for geometric functions pt 3: \(\beta\) aligned-frame (AF) parameters with geometric functions.
01/09/20 v2 with verified AF code.
09/06/20 v1
Aims:
Develop \(\beta_{L,M}\) formalism for AF, using geometric tensor formalism as already applied to MF case.
Develop corresponding numerical methods - see pt 1 notebook.
Analyse geometric terms - see pt 1 notebook.
\(\beta_{L,M}^{AF}\) rewrite
The various terms defined in pt 1 can be used to redefine the full AF observables, expressed as a set of \(\beta_{L,M}\) coefficients (with the addition of another tensor to define the alignment terms).
The original (full) form for the AF equations, as implemented in ``ePSproc.afblm` <https://epsproc.readthedocs.io/en/dev/modules/epsproc.AFBLM.html>`__ (note, however, that the previous implementation is not fully tested, since it was s…l…o…w… the geometric version should avoid this issue):
Where \(I_{l,m,\mu}^{p_{i}\mu_{i},p_{f}\mu_{f}}(E)\) are the energy-dependent dipole matrix elements, and \(A_{Q,S}^{K}(t)\) define the alignment parameters.
In terms of the geometric parameters, this can be rewritten as:
Where there’s a new alignment tensor:
And the the \(\Lambda_{R',R}\) term is a simplified form of the previously derived MF term:
All phase conventions should be as the MF case, and the numerics for all the ternsors can be used as is… hopefully…
Further notes:
Note \(B_{L,S-R'}(l,l',m,m')}\) instead of \(B_{L,-M}(l,l',m,m')}\) for MF case. This allows for all MF projections to contribute.
Refs for the full AF-PAD formalism above:
Reid, Katharine L., and Jonathan G. Underwood. “Extracting Molecular Axis Alignment from Photoelectron Angular Distributions.” The Journal of Chemical Physics 112, no. 8 (2000): 3643. https://doi.org/10.1063/1.480517.
Underwood, Jonathan G., and Katharine L. Reid. “Time-Resolved Photoelectron Angular Distributions as a Probe of Intramolecular Dynamics: Connecting the Molecular Frame and the Laboratory Frame.” The Journal of Chemical Physics 113, no. 3 (2000): 1067. https://doi.org/10.1063/1.481918.
Stolow, Albert, and Jonathan G. Underwood. “Time-Resolved Photoelectron Spectroscopy of Non-Adiabatic Dynamics in Polyatomic Molecules.” In Advances in Chemical Physics, edited by Stuart A. Rice, 139:497–584. Advances in Chemical Physics. Hoboken, NJ, USA: John Wiley & Sons, Inc., 2008. https://doi.org/10.1002/9780470259498.ch6.
Where [3] has the version as per the full form above (full asymmetric top alignment distribution expansion).
To consider
Normalisation for ADMs? Will matter in cases where abs cross-sections are valid (but not for PADs generally).
Status
10/09/20 Now working - see AF/LF verification notebook for more details. May still have function definition issues in some cases?
TODO: more careful comparison with experimental data (see old processing notebooks…)
Setup
[1]:
# Imports
import numpy as np
import pandas as pd
import xarray as xr
# Special functions
# from scipy.special import sph_harm
import spherical_functions as sf
import quaternion
# Performance & benchmarking libraries
# from joblib import Memory
# import xyzpy as xyz
import numba as nb
# Timings with ttictoc or time
# https://github.com/hector-sab/ttictoc
from ttictoc import TicToc
import time
# Package fns.
# For module testing, include path to module here
import sys
import os
modPath = r'D:\code\github\ePSproc' # Win test machine
# modPath = r'/home/femtolab/github/ePSproc/' # Linux test machine
sys.path.append(modPath)
import epsproc as ep
# TODO: tidy this up!
from epsproc.util import matEleSelector
from epsproc.geomFunc import geomCalc
* pyevtk not found, VTK export not available.
[2]:
from IPython.core.interactiveshell import InteractiveShell
InteractiveShell.ast_node_interactivity = "all"
[3]:
# Plotters
import matplotlib.pyplot as plt
from epsproc.plot import hvPlotters
hvPlotters.setPlotters()
Alignment terms
Axis distribution moments
These are already set up by setADMs().
[4]:
# set default alignment terms - single term A(0,0,0)=1, which corresponds to isotropic distributions
AKQS = ep.setADMs()
AKQS
[4]:
<xarray.DataArray (ADM: 1, t: 1)> array([[1]]) Coordinates: * ADM (ADM) MultiIndex - K (ADM) int64 0 - Q (ADM) int64 0 - S (ADM) int64 0 * t (t) int32 0 Attributes: dataType: ADM
[5]:
AKQSpd,_ = ep.util.multiDimXrToPD(AKQS, colDims='t', squeeze=False)
AKQSpd
[5]:
t | 0 | ||
---|---|---|---|
K | Q | S | |
0 | 0 | 0 | 1 |
[6]:
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(AKQS, xDim = 't', pType = 'r', squeeze = False) # Note squeeze = False required for 1D case (should add this to code!)
daPlotpd
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[6]:
t | 0 | ||
---|---|---|---|
K | Q | S | |
0 | 0 | 0 | 1.0 |
[7]:
# Test multiple t points
AKQS = ep.setADMs(ADMs = [0,0,0,1,1], t=[0,1])
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(AKQS, xDim = 't', pType = 'r', squeeze = False) # Note squeeze = False required for 1D case (should add this to code!)
daPlotpd
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[7]:
t | 0 | 1 | ||
---|---|---|---|---|
K | Q | S | ||
0 | 0 | 0 | 1.0 | 1.0 |
[8]:
# Test additional time-dependent values
AKQS = ep.setADMs(ADMs = [[0,0,0,1,1],[2,0,0,0,0.5]], t=[0,1]) # Nested list or np.array OK.
# AKQS = ep.setADMs(ADMs = np.array([[0,0,0,1,1],[2,0,0,0,0.5]]), t=[0,1])
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(AKQS, xDim = 't', pType = 'r', squeeze = False) # Note squeeze = False required for 1D case (should add this to code!)
daPlotpd
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[8]:
t | 0 | 1 | ||
---|---|---|---|---|
K | Q | S | ||
0 | 0 | 0 | 1.0 | 1.0 |
2 | 0 | 0 | NaN | 0.5 |
Alignment tensor
As previously defined:
Function dev
Based on existing MF functions, see ep.geomFunc.geomCalc and ep.geomFunc.mfblmGeom.
Most similar to existing betaTerm() function, which also has product of 3j terms in a similar manner.
[9]:
# Generate QNs - code adapted from ep.geomFunc.geomUtils.genllL(Lmin = 0, Lmax = 10, mFlag = True):
# Generate QNs for deltaKQS term - 3j product term
def genKQSterms(Kmin = 0, Kmax = 2, mFlag = True):
# Set QNs for calculation, (l,m,mp)
QNs = []
for P in np.arange(0, 2+1): # HARD-CODED R for testing - should get from EPR tensor defn. in full calcs.
for K in np.arange(Kmin, Kmax+1): # Eventually this will come from alignment term
for L in np.arange(np.abs(P-K), P+K+1): # Allowed L given P and K defined
if mFlag: # Include "m" (or equivalent) terms?
mMax = L
RMax = P
QMax = K
else:
mMax = 0
RMax = 0
QMax = 0
for R in np.arange(-RMax, RMax+1):
for Q in np.arange(-QMax, QMax+1):
#for M in np.arange(np.abs(l-lp), l+lp+1):
# for M in np.arange(-mMax, mMax+1):
# Set M - note this implies specific phase choice.
# M = -(m+mp)
# M = (-m+mp)
# if np.abs(M) <= L: # Skip terms with invalid M
# QNs.append([l, lp, L, m, mp, M])
# Run for all possible M
for M in np.arange(-L, L+1):
QNs.append([P, K, L, R, Q, M])
return np.array(QNs)
# Generate QNs from EPR + AKQS tensors
def genKQStermsFromTensors(EPR, AKQS, uniqueFlag = True, phaseConvention = 'S'):
'''
Generate all QNs for :math:`\Delta_{L,M}(K,Q,S)` from existing tensors (Xarrays) :math:`E_{P,R}` and :math:`A^K_{Q,S}`.
Cf. :py:func:`epsproc.geomFunc.genllpMatE`, code adapted from there.
Parameters
----------
matE : Xarray
Xarray containing matrix elements, with QNs (l,m), as created by :py:func:`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 :py:func:`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).
'''
# Local import.
from epsproc.geomFunc.geomCalc import setPhaseConventions
# For transparency/consistency with subfunctions, str/dict now set in setPhaseConventions()
if 'phaseCons' in EPR.attrs.keys():
phaseCons = EPR.attrs['phaseCons']
else:
phaseCons = setPhaseConventions(phaseConvention = phaseConvention)
# Get QNs from inputs
KQScoords = AKQS.unstack().coords # Use unstack here, or np.unique(matE.l), to avoid duplicates
PRcoords = EPR.unstack().coords
# Use passed (m,mp) values, or run for m=mp=0 only.
# if mFlag:
# mList = matE.unstack().m.values
# else:
# mList = 0
# Set QNs for calculation, one set for each 3j term
QNs1 = []
QNs2 = []
for P in PRcoords['P'].values: # Note dictionary syntax for coords objects
for K in KQScoords['K'].values:
for L in np.arange(np.abs(P-K), P+K+1): # Allowed L given P and K defined
# if mFlag: # Include "m" (or equivalent) terms?
# mMax = L
# RMax = P
# QMax = K
# else:
# mMax = 0
# RMax = 0
# QMax = 0
for R in PRcoords['R'].values:
for Q in KQScoords['Q'].values:
# Set M, with +/- phase convention - TBC MAY BE INCORRECT IN THIS CASE/CONTEXT?
# Note that setting phaseCons['afblmCons']['negM'] = phaseCons['genMatEcons']['negm'] is current default case, but doesn't have to be!
if phaseCons['genMatEcons']['negm']:
M = (-R+Q) # Case for M -> -M switch
else:
M = -(R+Q) # Usual phase convention.
if (abs(R)<=P) and (abs(Q)<=K) and (abs(M)<=L): # Check term is valid.
QNs1.append([P, K, L, R, Q, M])
# Set Rp and S - these are essentially independent of R,Q,M, but keep nested for full dim output.
for Rp in PRcoords['R'].values:
for S in KQScoords['S'].values:
if phaseCons['genMatEcons']['negm']:
SRp = (-Rp+S) # Case for M -> -M switch
else:
SRp = -(Rp+S) # Usual phase convention.
# SRp = S-Rp # Set final 3j term, S-Rp
if (abs(Rp)<=P) and (abs(S)<=K) and (abs(SRp)<=L): # Check term is valid.
QNs2.append([P, K, L, Rp, S, SRp])
#for M in np.arange(np.abs(l-lp), l+lp+1):
# for M in np.arange(-mMax, mMax+1):
# Set M - note this implies specific phase choice.
# M = -(m+mp)
# M = (-m+mp)
# if np.abs(M) <= L: # Skip terms with invalid M
# QNs.append([l, lp, L, m, mp, M])
# Run for all possible M
# for M in np.arange(-L, L+1):
# QNs.append([P, K, L, R, Q, M])
if uniqueFlag:
return np.unique(QNs1, axis = 0), np.unique(QNs2, axis = 0)
else:
return np.array(QNs1), np.array(QNs2)
[10]:
test = genKQSterms()
# test.tolist() # Use tolist() to get full array output (not truncated by np)
test.shape
[10]:
(1225, 6)
[11]:
# Can also just use existing fn.
test2 = ep.geomFunc.genllL(Lmax=2)
test2.shape
[11]:
(1225, 6)
[12]:
# Also QN list fn - checks and removes duplicates
QNs = ep.geomFunc.genllL(Lmax=2)
test3 = ep.geomFunc.geomUtils.genllLList(QNs, uniqueFlag = True, mFlag = True)
test3.shape
[12]:
(1225, 6)
All QNs
[13]:
# Then calc 3js.... as per betaTerm
form = 'xdaLM' # xds
dlist1 = ['P', 'K', 'L', 'R', 'Q', 'M']
dlist2 = ['P', 'K', 'L', 'Rp', 'S', 'S-Rp']
QNs = ep.geomFunc.genllL(Lmax=2)
[14]:
# Set phase conventions for this case, extending existing structure
phaseCons = ep.geomFunc.setPhaseConventions('E')
phaseCons['afblmCons'] = {}
# (+/-)M phase selection, set as per existing code, betaCons['negM'] = genMatEcons['negm'] # Use -M term in 3j? Should be anti-correlated with genMatEcons['negm']...? 31/03/20 NOW correlated with mfblmCons['Mphase']
# Note this is correlated with QN generation in genllpMatE() - should set equivalent fn for alignment terms.
# In existing case this arises from M = (-m+mp) or M = -(m+mp) choice.
phaseCons['afblmCons']['negM'] = phaseCons['genMatEcons']['negm']
phaseCons['afblmCons']['negQ'] = True
phaseCons['afblmCons']['negS'] = True
# Apply phase conventions to input QNs
# if phaseCons['mfblmCons']['BLMmPhase']:
# QNsBLMtable[:,3] *= -1
# QNsBLMtable[:,5] *= -1
[15]:
# Calculate two 3j terms, with respective QN sets
thrj1 = ep.geomFunc.w3jTable(QNs = QNs, nonzeroFlag = True, form = form, dlist = dlist1)
thrj2 = ep.geomFunc.w3jTable(QNs = QNs, nonzeroFlag = True, form = form, dlist = dlist2)
[16]:
# Multiplication term...
testMult = thrj1.unstack() * thrj2.unstack()
# This can get large quickly - for Kmax=2 already have 2e6 terms
[17]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, xDim=xDim, pType = 'r')
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
Reduced QNs set
[18]:
# Then calc 3js.... as per betaTerm
form = 'xdaLM' # xds
dlist1 = ['P', 'K', 'L', 'R', 'Q', 'M']
dlist2 = ['P', 'K', 'L', 'Rp', 'S', 'S-Rp']
# Calculate two 3j terms, with respective QN sets
thrj1 = ep.geomFunc.w3jTable(QNs = QNs, nonzeroFlag = True, form = form, dlist = dlist1)
thrj2 = ep.geomFunc.w3jTable(QNs = QNs, nonzeroFlag = True, form = form, dlist = dlist2)
[19]:
# With subselection (currently don't have QN generation fnc of correct form - see geomUtils)
# from epsproc.util import matEleSelector
# thrj1Sel = ep.util.matEleSelector(thrj1, inds = {'Q':0}, sq=False)
[20]:
# Check 3j terms
xDim = {'LM':['L','M']}
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(thrj1, xDim=xDim, pType = 'r')
xDim = {'LM':['L','S-Rp']}
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(thrj2.sel({'S':0}), xDim=xDim, pType = 'r')
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[21]:
# Multiplication term...
testMultSub = thrj1.unstack().sel({'Q':0}) * thrj2.unstack().sel({'S':0}) # SLOW, 91125 elements. ACTUALLY - much faster after a reboot. Sigh.
testMultSub2 = thrj1.sel({'Q':0}).unstack() * thrj2.sel({'S':0}).unstack() # FAST, 28125 elements. Possible issue with which dims are dropped - check previous notes!
# This can get large quickly - for Kmax=2 already have 2e6 terms
testMultSub.notnull().sum() # == 194 for Kmax = 2, Q=S=0
testMultSub2.notnull().sum() # == 194 for Kmax = 2, Q=S=0 OK
[21]:
<xarray.DataArray 'w3jStacked' ()> array(194) Coordinates: Q int64 0 S int64 0
[21]:
<xarray.DataArray 'w3jStacked' ()> array(194)
[22]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMultSub, xDim=xDim, pType = 'r')
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[23]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMultSub2, xDim=xDim, pType = 'r')
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
Plots appear to be the same, aside from ordering of M terms(?)
QNs from existing tensors
[24]:
phaseConvention = 'E'
# Set polarisation term
p=[0]
EPRX = geomCalc.EPR(form = 'xarray', p = p, nonzeroFlag = True, phaseConvention = phaseConvention).unstack().sel({'R-p':0}).drop('R-p') # Set for R-p = 0 for p=0 case (redundant coord) - need to fix in e-field mult term!
EPRXresort = EPRX.squeeze(['l','lp']).drop(['l','lp']) # Safe squeeze & drop of selected singleton dims only.
# Set alignment terms (inc. time-dependence)
AKQS = ep.setADMs(ADMs = [[0,0,0,1,1],[2,0,0,0,0.5]], t=[0,1]) # Nested list or np.array OK.
[25]:
# Calculate alignment term - this cell should form core function, cf. betaTerm() etc.
# Set QNs
QNs1, QNs2 = genKQStermsFromTensors(EPRXresort, AKQS, uniqueFlag = True, phaseConvention = phaseConvention)
# Then calc 3js.... as per betaTerm
form = 'xdaLM' # xds
dlist1 = ['P', 'K', 'L', 'R', 'Q', 'M']
dlist2 = ['P', 'K', 'L', 'Rp', 'S', 'S-Rp']
# Copy QNs and apply any additional phase conventions
QNs1DeltaTable = QNs1.copy()
QNs2DeltaTable = QNs2.copy()
# Set additional phase cons here - these will be set in master function eventually!
# NOTE - only testing for Q=S=0 case initially.
phaseCons['afblmCons']['negM'] = phaseCons['genMatEcons']['negm'] # IF SET TO TRUE THIS KNOCKS OUT M!=0 terms - not sure if this is correct here, depends also on phase cons in genKQStermsFromTensors().
# Yeah, looks like phase error in current case, get terms with R=M, instead of R=-M
# Confusion is due to explicit assignment of +/-M terms in QN generation (only allowed terms), which *already* enforces this phase convention.
phaseCons['afblmCons']['negQ'] = True
phaseCons['afblmCons']['negS'] = True
# Switch signs (m,M) before 3j calcs.
if phaseCons['afblmCons']['negQ']:
QNs1DeltaTable[:,4] *= -1
# Switch sign Q > -Q before 3j calcs.
if phaseCons['afblmCons']['negM']:
QNs1DeltaTable[:,5] *= -1
# Switch sign S > -S before 3j calcs.
if phaseCons['afblmCons']['negS']:
QNs2DeltaTable[:,4] *= -1
# Calculate two 3j terms, with respective QN sets
thrj1 = ep.geomFunc.w3jTable(QNs = QNs1DeltaTable, nonzeroFlag = True, form = form, dlist = dlist1)
thrj2 = ep.geomFunc.w3jTable(QNs = QNs2DeltaTable, nonzeroFlag = True, form = form, dlist = dlist2)
# Multiply
thrjMult = thrj1.unstack() * thrj2.unstack()
# Additional terms & multiplications
Kdegen = np.sqrt(2*thrjMult.K + 1)
KQphase = np.power(-1, np.abs(thrjMult.K + thrjMult.Q))
DeltaKQS = Kdegen * KQphase * thrjMult
# AF term
AFterm = (DeltaKQS * AKQS.unstack()).sum({'K','Q','S'})
[26]:
thrjMult.notnull().sum() # == 69 for test case with all phase switches on, and same for no phase switches (in test case Q=S=0 in any case!)
[26]:
<xarray.DataArray 'w3jStacked' ()> array(69)
[27]:
# Plot
xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(DeltaKQS, xDim=xDim, pType = 'r')
daPlotpd
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[27]:
L | 0 | 1 | 2 | 3 | 4 | ||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
M | 0 | -1 | 0 | 1 | -1 | 0 | 1 | -1 | 0 | 1 | -1 | 0 | 1 | ||||||
K | P | Q | R | Rp | S | S-Rp | |||||||||||||
0 | 0 | 0 | 0 | 0 | 0 | 0 | 1.000000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
1 | 0 | -1 | -1 | 0 | 1 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
0 | 0 | 0 | NaN | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
0 | -1 | 0 | 1 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | 1 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
2 | 0 | -1 | -1 | 0 | 1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | |
0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | ||||
0 | -1 | 0 | 1 | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | 1 | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
2 | 0 | 0 | 0 | 0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
1 | 0 | -1 | -1 | 0 | 1 | NaN | NaN | NaN | 0.074536 | NaN | NaN | 0.223607 | NaN | NaN | 0.127775 | NaN | NaN | NaN | |
0 | 0 | 0 | NaN | NaN | NaN | 0.149071 | NaN | NaN | NaN | NaN | NaN | -0.156492 | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | 0.074536 | NaN | NaN | -0.223607 | NaN | NaN | 0.127775 | NaN | NaN | NaN | ||||
0 | -1 | 0 | 1 | NaN | NaN | 0.149071 | NaN | NaN | NaN | NaN | NaN | -0.156492 | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | NaN | 0.298142 | NaN | NaN | NaN | NaN | NaN | 0.191663 | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | NaN | 0.149071 | NaN | NaN | NaN | NaN | NaN | -0.156492 | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | 1 | NaN | 0.074536 | NaN | NaN | -0.223607 | NaN | NaN | 0.127775 | NaN | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | NaN | 0.149071 | NaN | NaN | NaN | NaN | NaN | -0.156492 | NaN | NaN | NaN | NaN | NaN | ||||
1 | 0 | -1 | NaN | 0.074536 | NaN | NaN | 0.223607 | NaN | NaN | 0.127775 | NaN | NaN | NaN | NaN | NaN | ||||
2 | 0 | -1 | -1 | 0 | 1 | NaN | NaN | NaN | 0.223607 | NaN | NaN | 0.031944 | NaN | NaN | 0.063888 | NaN | NaN | 0.106479 | |
0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.063888 | NaN | NaN | NaN | NaN | NaN | -0.116642 | ||||
1 | 0 | -1 | NaN | NaN | NaN | -0.223607 | NaN | NaN | 0.031944 | NaN | NaN | -0.063888 | NaN | NaN | 0.106479 | ||||
0 | -1 | 0 | 1 | NaN | NaN | NaN | NaN | NaN | -0.063888 | NaN | NaN | NaN | NaN | NaN | -0.116642 | NaN | |||
0 | 0 | 0 | 0.447214 | NaN | NaN | NaN | NaN | 0.127775 | NaN | NaN | NaN | NaN | NaN | 0.127775 | NaN | ||||
1 | 0 | -1 | NaN | NaN | NaN | NaN | NaN | -0.063888 | NaN | NaN | NaN | NaN | NaN | -0.116642 | NaN | ||||
1 | -1 | 0 | 1 | NaN | -0.223607 | NaN | NaN | 0.031944 | NaN | NaN | -0.063888 | NaN | NaN | 0.106479 | NaN | NaN | |||
0 | 0 | 0 | NaN | NaN | NaN | NaN | -0.063888 | NaN | NaN | NaN | NaN | NaN | -0.116642 | NaN | NaN | ||||
1 | 0 | -1 | NaN | 0.223607 | NaN | NaN | 0.031944 | NaN | NaN | 0.063888 | NaN | NaN | 0.106479 | NaN | NaN |
[28]:
# Plot
xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(AFterm, xDim=xDim, pType = 'r')
daPlotpd
Set dataType (No dataType)
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[28]:
L | 0 | 1 | 2 | 3 | 4 | ||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
M | 0 | -1 | 0 | 1 | -1 | 0 | 1 | -1 | 0 | 1 | -1 | 0 | 1 | ||||
P | R | Rp | S-Rp | t | |||||||||||||
0 | 0 | 0 | 0 | 0 | 1.000000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
1 | 1.000000 | NaN | NaN | NaN | NaN | 0.223607 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||
1 | -1 | -1 | 1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
1 | NaN | NaN | NaN | 0.370601 | NaN | NaN | 0.111803 | NaN | NaN | 0.063888 | NaN | NaN | NaN | ||||
0 | 0 | 0 | NaN | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | -0.258798 | NaN | NaN | NaN | NaN | NaN | -0.078246 | NaN | NaN | NaN | ||||
1 | -1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | 0.370601 | NaN | NaN | -0.111803 | NaN | NaN | 0.063888 | NaN | NaN | NaN | ||||
0 | -1 | 1 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
1 | NaN | NaN | -0.258798 | NaN | NaN | NaN | NaN | NaN | -0.078246 | NaN | NaN | NaN | NaN | ||||
0 | 0 | 0 | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | 0.482405 | NaN | NaN | NaN | NaN | NaN | 0.095831 | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | -0.258798 | NaN | NaN | NaN | NaN | NaN | -0.078246 | NaN | NaN | NaN | NaN | ||||
1 | -1 | 1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
1 | NaN | 0.370601 | NaN | NaN | -0.111803 | NaN | NaN | 0.063888 | NaN | NaN | NaN | NaN | NaN | ||||
0 | 0 | 0 | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | -0.258798 | NaN | NaN | NaN | NaN | NaN | -0.078246 | NaN | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | 0.370601 | NaN | NaN | 0.111803 | NaN | NaN | 0.063888 | NaN | NaN | NaN | NaN | NaN | ||||
2 | -1 | -1 | 1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN |
1 | NaN | NaN | NaN | 0.111803 | NaN | NaN | 0.215972 | NaN | NaN | 0.031944 | NaN | NaN | 0.053240 | ||||
0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | NaN | -0.231944 | NaN | NaN | NaN | NaN | NaN | -0.058321 | ||||
1 | -1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | -0.111803 | NaN | NaN | 0.215972 | NaN | NaN | -0.031944 | NaN | NaN | 0.053240 | ||||
0 | -1 | 1 | 0 | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
1 | NaN | NaN | NaN | NaN | NaN | -0.231944 | NaN | NaN | NaN | NaN | NaN | -0.058321 | NaN | ||||
0 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | 0.223607 | NaN | NaN | NaN | NaN | 0.263888 | NaN | NaN | NaN | NaN | NaN | 0.063888 | NaN | ||||
1 | -1 | 0 | NaN | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | -0.231944 | NaN | NaN | NaN | NaN | NaN | -0.058321 | NaN | ||||
1 | -1 | 1 | 0 | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
1 | NaN | -0.111803 | NaN | NaN | 0.215972 | NaN | NaN | -0.031944 | NaN | NaN | 0.053240 | NaN | NaN | ||||
0 | 0 | 0 | NaN | NaN | NaN | NaN | -0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | -0.231944 | NaN | NaN | NaN | NaN | NaN | -0.058321 | NaN | NaN | ||||
1 | -1 | 0 | NaN | NaN | NaN | NaN | 0.200000 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
1 | NaN | 0.111803 | NaN | NaN | 0.215972 | NaN | NaN | 0.031944 | NaN | NaN | 0.053240 | NaN | NaN |
Note here that there are blocks of non-zero terms with M=+/-1, these should drop out later (?) by sums over symmetry and/or other 3j terms… TBC…
Lambda term redux
Use existing function and force/sub-select terms…
[29]:
# Code adapted from mfblmXprod()
eulerAngs = np.array([0,0,0], ndmin=2)
# RX = ep.setPolGeoms(eulerAngs = eulerAngs) # This throws error in geomCalc.MFproj???? Something to do with form of terms passed to wD, line 970 vs. 976 in geomCalc.py
RX = ep.setPolGeoms() # (0,0,0) term in geomCalc.MFproj OK.
lambdaTerm, lambdaTable, lambdaD, QNsLambda = geomCalc.MFproj(RX = RX, form = 'xarray', phaseConvention = phaseConvention)
# lambdaTermResort = lambdaTerm.squeeze().drop('l').drop('lp') # This removes photon (l,lp) dims fully.
lambdaTermResort = lambdaTerm.squeeze(['l','lp']).drop(['l','lp']).sel({'Labels':'z'}) # Safe squeeze & drop of selected singleton dims only.
[30]:
# Plot
xDim = {'PRp':['P','Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(lambdaTermResort, xDim=xDim, pType = 'r')
daPlotpd
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[30]:
P | 0 | 1 | 2 | ||||||||
---|---|---|---|---|---|---|---|---|---|---|---|
Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | ||
R | mu | mup | |||||||||
-2 | 1 | 1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN |
-1 | 0 | 1 | NaN | -0.408248 | NaN | NaN | NaN | 0.316228 | NaN | NaN | NaN |
1 | 0 | NaN | 0.408248 | NaN | NaN | NaN | 0.316228 | NaN | NaN | NaN | |
0 | -1 | 1 | 0.57735 | NaN | -0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN |
0 | 0 | -0.57735 | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | |
1 | -1 | 0.57735 | NaN | 0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | |
1 | -1 | 0 | NaN | NaN | NaN | -0.408248 | NaN | NaN | NaN | 0.316228 | NaN |
0 | -1 | NaN | NaN | NaN | 0.408248 | NaN | NaN | NaN | 0.316228 | NaN | |
2 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
Check component terms are as expected
Should have R = Rp for z case, and wD terms = 1 or 0.
[31]:
lambdaTablepd, _ = ep.util.multiDimXrToPD(lambdaTable, colDims=xDim, dropna=True)
lambdaTablepd
[31]:
P | 0 | 1 | 2 | ||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | ||||
R | l | lp | mu | mup | |||||||||
-2 | 1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | ||||
1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.182574 | NaN | NaN | ||||
0 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.182574 | NaN | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN | ||||
-1 | 1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | ||||
1 | NaN | NaN | -0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | ||||
0 | -1 | NaN | NaN | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | ||||
1 | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | -1 | NaN | NaN | 0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | |||
0 | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN | ||||
0 | 1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | ||||
1 | 0.57735 | NaN | -0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | ||||
0 | -1 | NaN | NaN | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | |||
0 | -0.57735 | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | ||||
1 | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | -1 | 0.57735 | NaN | 0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | |||
0 | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN | ||||
1 | 1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | ||||
1 | NaN | NaN | -0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | ||||
0 | -1 | NaN | NaN | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | ||||
1 | NaN | 0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | -1 | NaN | NaN | 0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | |||
0 | NaN | -0.408248 | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN | ||||
2 | 1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | ||||
1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.182574 | NaN | NaN | ||||
0 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.182574 | NaN | NaN | |||
0 | NaN | NaN | NaN | NaN | NaN | -0.316228 | NaN | NaN | NaN | ||||
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN |
[32]:
lambdaTermResortpd, _ = ep.util.multiDimXrToPD(lambdaTermResort, colDims=xDim, dropna=True)
lambdaTermResortpd
# Think this is as per lambdaTable terms, just different ordering - because lambdaTable doesn't include some phase switches?
[32]:
P | 0 | 1 | 2 | ||||||||
---|---|---|---|---|---|---|---|---|---|---|---|
Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | ||
R | mu | mup | |||||||||
-2 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j |
0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
0 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | 0.447214+0.000000j | NaN | NaN | NaN | NaN | ||
-1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j |
0 | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | ||
1 | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
0 | -1 | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
1 | NaN | -0.408248+0.000000j | NaN | NaN | NaN | 0.316228+0.000000j | NaN | NaN | NaN | ||
1 | -1 | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | |
0 | NaN | 0.408248+0.000000j | NaN | NaN | NaN | 0.316228+0.000000j | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | NaN | ||
0 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j |
0 | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | ||
1 | 0.577350+0.000000j | NaN | -0.408248+0.000000j | NaN | NaN | NaN | 0.182574+0.000000j | NaN | NaN | ||
0 | -1 | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | |
0 | -0.577350+0.000000j | NaN | NaN | NaN | NaN | NaN | 0.365148+0.000000j | NaN | NaN | ||
1 | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | -1 | 0.577350+0.000000j | NaN | 0.408248+0.000000j | NaN | NaN | NaN | 0.182574+0.000000j | NaN | NaN | |
0 | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | NaN | ||
1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j |
0 | NaN | NaN | NaN | -0.408248+0.000000j | NaN | NaN | NaN | 0.316228+0.000000j | NaN | ||
1 | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
0 | -1 | NaN | NaN | NaN | 0.408248+0.000000j | NaN | NaN | NaN | 0.316228+0.000000j | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
1 | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | -1 | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | |
0 | NaN | 0.000000+0.000000j | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | NaN | ||
2 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214+0.000000j |
0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
0 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | ||
1 | NaN | NaN | NaN | NaN | 0.000000+0.000000j | NaN | NaN | NaN | NaN |
[33]:
# wigner D term looks good.
lambdaDpd, _ = ep.util.multiDimXrToPD(lambdaD.sel({'Labels':'z'}), colDims=xDim, dropna=True)
lambdaDpd
[33]:
P | 0 | 1 | 2 | ||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Rp | 2 | 1 | 0 | -1 | -2 | 2 | 1 | 0 | -1 | -2 | 2 | 1 | 0 | -1 | -2 |
R | |||||||||||||||
2 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j |
1 | NaN | NaN | NaN | NaN | NaN | 0.000000-0.000000j | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j |
0 | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j |
-1 | NaN | NaN | NaN | NaN | NaN | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000+0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000+0.000000j | 0.000000-0.000000j |
-2 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 0.000000-0.000000j | 1.000000+0.000000j |
[34]:
# NOTE - accidental correlation R==Rp here unless dim is dropped
# xDim = {'PRp':['P','Rp']}
# lambdaDpd, _ = ep.util.multiDimXrToPD(lambdaD.sel({'Labels':'z'}).unstack().sum('R'), colDims=xDim, dropna=True, squeeze=False) # Throwing output errors... not sure why, has values and dims when tested at cmd???
# lambdaDpd
[35]:
# NOTE - accidental correlation R==Rp here unless dim is dropped
xDim = {'PRp':['P','Rp']}
lambdaTermResortpd, _ = ep.util.multiDimXrToPD(lambdaTermResort.sum('R'), colDims=xDim, dropna=True)
lambdaTermResortpd
[35]:
P | 0 | 1 | 2 | |||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Rp | -2 | -1 | 0 | 1 | 2 | -2 | -1 | 0 | 1 | 2 | -2 | -1 | 0 | 1 | 2 | |
mu | mup | |||||||||||||||
-1 | -1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.447214+0.000000j |
0 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | -0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.316228+0.000000j | 0.000000+0.000000j | |
1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.577350+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | -0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.182574+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | |
0 | -1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.316228+0.000000j | 0.000000+0.000000j |
0 | 0.000000+0.000000j | 0.000000+0.000000j | -0.577350+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.365148+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | |
1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | -0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.316228+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | |
1 | -1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.577350+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.182574+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j |
0 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.408248+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.316228+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | |
1 | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.447214+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j | 0.000000+0.000000j |
Build full calculation from functions
Use mfblmGeom.py as template: basically just need modified lambda term as above, and new alignment term, and rest of calculation should be identical.
[36]:
# Check polProd term - incorporate alignment term here...?
# Existing terms
# polProd = (EPRXresort * lambdaTermResort)
# sumDimsPol = ['P','R','Rp','p']
# polProd = polProd.sum(sumDimsPol)
# Test with alignment term
polProd = (EPRXresort * lambdaTermResort * AFterm)
sumDimsPol = ['P','R','Rp','p', 'S-Rp']
polProd = polProd.sum(sumDimsPol)
# Looks OK - keeps correct dims in test case! 270 terms.
NOW implemented in ep.geomFunc.afblmGeom.py
Testing…
16/06/20
Test code adapted from previous round of AF tests, http://localhost:8888/lab/tree/dev/ePSproc/ePSproc_AFBLM_calcs_bench_100220.ipynb
See also MFBLM test code, http://localhost:8888/lab/tree/dev/ePSproc/geometric_method_dev_2020/geometric_method_dev_pt2_170320_v090620.ipynb
Load data
[37]:
# Load data from modPath\data
dataPath = os.path.join(modPath, 'data', 'photoionization')
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)
dataXS = ep.readMatEle(fileIn = dataFile, recordType = 'CrossSection') # XS info currently not set in NO2 sample file.
*** ePSproc readMatEle(): scanning files for DumpIdy segments.
*** Scanning file(s)
['D:\\code\\github\\ePSproc\\data\\photoionization\\n2_3sg_0.1-50.1eV_A2.inp.out']
*** Reading ePS output file: D:\code\github\ePSproc\data\photoionization\n2_3sg_0.1-50.1eV_A2.inp.out
Expecting 51 energy points.
Expecting 2 symmetries.
Scanning CrossSection segments.
Expecting 102 DumpIdy segments.
Found 102 dumpIdy segments (sets of matrix elements).
Processing segments to Xarrays...
Processed 102 sets of DumpIdy file segments, (0 blank)
*** ePSproc readMatEle(): scanning files for CrossSection segments.
*** Scanning file(s)
['D:\\code\\github\\ePSproc\\data\\photoionization\\n2_3sg_0.1-50.1eV_A2.inp.out']
*** Reading ePS output file: D:\code\github\ePSproc\data\photoionization\n2_3sg_0.1-50.1eV_A2.inp.out
Expecting 51 energy points.
Expecting 2 symmetries.
Scanning CrossSection segments.
Expecting 3 CrossSection segments.
Found 3 CrossSection segments (sets of results).
Processed 3 sets of CrossSection file segments, (0 blank)
Plot ePS results (isotropic case)
[38]:
# Plot cross sections using Xarray functionality
# Set here to plot per file - should add some logic to combine files.
for data in dataXS:
daPlot = data.sel(XC='SIGMA')
daPlot.plot.line(x='Eke', col='Type')
[38]:
<xarray.plot.facetgrid.FacetGrid at 0x1da1a0574e0>
[39]:
# Repeat for betas
for data in dataXS:
daPlot = data.sel(XC='BETA')
daPlot.plot.line(x='Eke', col='Type')
[39]:
<xarray.plot.facetgrid.FacetGrid at 0x1da1838cdd8>
Try new AF calculation - isotropic (default) case
[40]:
# Tabulate & plot matrix elements vs. Eke
selDims = {'it':1, 'Type':'L'}
# selDims = {'Type':'L'}
matE = dataSet[0].sel(selDims) # Set for N2 case, length-gauge results only.
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(matE, xDim = 'Eke', pType = 'r', fillna = True)
daPlotpd
Plotting data n2_3sg_0.1-50.1eV_A2.inp.out, pType=r, thres=0.01, with Seaborn
[40]:
Eke | 0.1 | 1.1 | 2.1 | 3.1 | 4.1 | 5.1 | 6.1 | 7.1 | 8.1 | 9.1 | ... | 41.1 | 42.1 | 43.1 | 44.1 | 45.1 | 46.1 | 47.1 | 48.1 | 49.1 | 50.1 | |||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Cont | Targ | Total | l | m | mu | |||||||||||||||||||||
PU | SG | PU | 1 | -1 | 1 | -6.203556 | 7.496908 | 3.926892 | 1.071093 | -0.335132 | -1.043391 | -1.396885 | -1.554984 | -1.598812 | -1.573327 | ... | -0.396306 | -0.422225 | -0.448003 | -0.473466 | -0.498463 | -0.522867 | -0.546573 | -0.569493 | -0.591558 | -0.612712 |
1 | -1 | -6.203556 | 7.496908 | 3.926892 | 1.071093 | -0.335132 | -1.043391 | -1.396885 | -1.554984 | -1.598812 | -1.573327 | ... | -0.396306 | -0.422225 | -0.448003 | -0.473466 | -0.498463 | -0.522867 | -0.546573 | -0.569493 | -0.591558 | -0.612712 | ||||
3 | -1 | 1 | -2.090641 | -1.723467 | -3.571018 | -1.703191 | 0.232612 | 1.862713 | 3.201022 | 4.298839 | 5.198215 | 5.929303 | ... | 1.200463 | 1.040074 | 0.892501 | 0.757189 | 0.633556 | 0.521004 | 0.418928 | 0.326726 | 0.243800 | 0.169566 | |||
1 | -1 | -2.090641 | -1.723467 | -3.571018 | -1.703191 | 0.232612 | 1.862713 | 3.201022 | 4.298839 | 5.198215 | 5.929303 | ... | 1.200463 | 1.040074 | 0.892501 | 0.757189 | 0.633556 | 0.521004 | 0.418928 | 0.326726 | 0.243800 | 0.169566 | ||||
5 | -1 | 1 | 0.000000 | 0.013246 | 0.000000 | -0.013096 | -0.024367 | -0.033455 | -0.042762 | -0.053765 | -0.067354 | -0.083946 | ... | -0.325544 | -0.318791 | -0.312332 | -0.306211 | -0.300456 | -0.295091 | -0.290130 | -0.285584 | -0.281453 | -0.277737 | |||
1 | -1 | 0.000000 | 0.013246 | 0.000000 | -0.013096 | -0.024367 | -0.033455 | -0.042762 | -0.053765 | -0.067354 | -0.083946 | ... | -0.325544 | -0.318791 | -0.312332 | -0.306211 | -0.300456 | -0.295091 | -0.290130 | -0.285584 | -0.281453 | -0.277737 | ||||
7 | -1 | 1 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.027469 | 0.028083 | 0.028678 | 0.029257 | 0.029823 | 0.030377 | 0.030923 | 0.031461 | 0.031994 | 0.032523 | |||
1 | -1 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.027469 | 0.028083 | 0.028678 | 0.029257 | 0.029823 | 0.030377 | 0.030923 | 0.031461 | 0.031994 | 0.032523 | ||||
SU | SG | SU | 1 | 0 | 0 | 6.246520 | -10.081763 | -8.912214 | -5.342049 | -3.160539 | -1.795783 | -0.867993 | -0.174154 | 0.400643 | 0.926420 | ... | -0.854886 | -0.940017 | -1.021577 | -1.099776 | -1.174787 | -1.246759 | -1.315814 | -1.382058 | -1.445582 | -1.506467 |
3 | 0 | 0 | 2.605768 | 2.775108 | 4.733473 | 1.677812 | -1.450830 | -4.187395 | -6.554491 | -8.601759 | -10.343109 | -11.747922 | ... | -0.293053 | -0.441416 | -0.579578 | -0.708075 | -0.827415 | -0.938075 | -1.040506 | -1.135134 | -1.222361 | -1.302566 | |||
5 | 0 | 0 | 0.000000 | -0.020864 | 0.000000 | 0.024263 | 0.042694 | 0.059142 | 0.077308 | 0.099330 | 0.126227 | 0.157745 | ... | 0.490499 | 0.511546 | 0.531810 | 0.551281 | 0.569951 | 0.587815 | 0.604868 | 0.621107 | 0.636528 | 0.651132 | |||
7 | 0 | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | -0.039566 | -0.041422 | -0.043265 | -0.045092 | -0.046901 | -0.048689 | -0.050453 | -0.052190 | -0.053899 | -0.055577 |
12 rows × 51 columns
Check BetaTerm - QN switch
[41]:
# Set parameters
SFflag = False
phaseConvention = 'E'
thres = 1e-1
selDims = {'it':1, 'Type':'L'}
thresDims = 'Eke'
# sumDims = ['mu', 'mup', 'l','lp','m','mp','S-Rp']
sumDims = ['mu', 'mup', 'l','lp','m','mp']
sumDimsPol = ['P','R','Rp','p']
# Set matrix elements
matE = dataSet[0].copy()
# Setting this here gives usual problem (non-matching XS values to getCro)
# Now testing multiplication at end of calc...
# if SFflag:
# matE.values = matE * matE.SF
matEthres = ep.matEleSelector(matE, thres = thres, inds = selDims, dims = thresDims, sq = True, drop = True)
# Code as per afblmGeom.py for testing - phase issues somewhere here...
QNs = ep.geomFunc.genllpMatE(matEthres, phaseConvention = phaseConvention)
QNsBLMtable = QNs.copy()
phaseCons = ep.geomFunc.setPhaseConventions(phaseConvention = phaseConvention)
# Switch signs (m,M) before 3j calcs.
# if phaseCons['mfblmCons']['BLMmPhase']:
# QNsBLMtable[:,3] *= -1
# QNsBLMtable[:,5] *= -1
#**** TESTING
# QNsBLMtable[:,3] *= -1 # m > -m # IF set, get S-Rp = -2....+2, otherwise only -1....+1
# QNsBLMtable[:,5] *= -1
BLMtable = geomCalc.betaTerm(QNs = QNsBLMtable, form = 'xdaLM', phaseConvention = phaseConvention)
# if BLMtableResort is None:
# Apply additional phase convention
BLMtableResort = BLMtable.copy().unstack()
# if phaseCons['mfblmCons']['negMcoordSwap']:
# BLMtableResort['M'] *= -1
# if phaseCons['mfblmCons']['Mphase']:
# BLMtableResort *= np.power(-1, np.abs(BLMtableResort.M)) # Associated phase term
# if phaseCons['mfblmCons']['negmCoordSwap']:
# BLMtableResort['m'] *= -1 # 25/06/20 Currently sends M>0 terms to 0. But suspicious this is incorrect at matE selection...
# Full phase cons testing....
if phaseCons['mfblmCons']['negMcoordSwap']:
BLMtableResort['M'] *= -1
if phaseCons['mfblmCons']['Mphase']:
BLMtableResort *= np.power(-1, np.abs(BLMtableResort.M)) # Associated phase term
if phaseCons['mfblmCons']['negmCoordSwap']:
BLMtableResort['m'] *= -1
if phaseCons['mfblmCons']['mPhase']:
BLMtableResort *= np.power(-1, np.abs(BLMtableResort.m)) # Associated phase term
# RENAME, M > (S-R') for AF case - this correctly allows for all MF projections!!!
# Some MF phase cons as applied above may also be incorrect?
BLMtableResort = BLMtableResort.rename({'M':'S-Rp'})
# BLMtableResort['m'] *= -1
[42]:
# Plot
xDim = {'LSRp':['L','S-Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(BLMtableResort, xDim=xDim, pType = 'r')
# pd.options.display.max_rows = 50
pd.set_option('display.max_rows', 50)
daPlotpd[0:20]
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[42]:
L | 0 | 2 | 4 | ... | 8 | 10 | ||||||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
S-Rp | 0 | -2 | -1 | 0 | 1 | 2 | -2 | -1 | 0 | 1 | ... | -2 | -1 | 0 | 1 | 2 | -2 | -1 | 0 | 1 | 2 | |||
l | lp | m | mp | |||||||||||||||||||||
1 | 1 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | -0.309019 | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
0 | NaN | NaN | NaN | NaN | 0.218510 | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | 0.282095 | NaN | NaN | -0.126157 | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | -1 | NaN | NaN | NaN | NaN | -0.218510 | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
0 | 0.282095 | NaN | NaN | 0.252313 | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | NaN | NaN | -0.218510 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | -1 | 0.282095 | NaN | NaN | -0.126157 | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
0 | NaN | NaN | 0.218510 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | NaN | -0.309019 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
3 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | 0.082589 | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
0 | NaN | NaN | NaN | NaN | -0.143048 | NaN | NaN | NaN | NaN | 0.194664 | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | NaN | NaN | NaN | 0.202301 | NaN | NaN | NaN | NaN | -0.150786 | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
0 | -1 | NaN | NaN | NaN | NaN | -0.233597 | NaN | NaN | NaN | NaN | -0.238414 | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
0 | NaN | NaN | NaN | 0.247767 | NaN | NaN | NaN | NaN | 0.246233 | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | NaN | NaN | -0.233597 | NaN | NaN | NaN | NaN | -0.238414 | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | -1 | NaN | NaN | NaN | 0.202301 | NaN | NaN | NaN | NaN | -0.150786 | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||
0 | NaN | NaN | -0.143048 | NaN | NaN | NaN | NaN | 0.194664 | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | NaN | 0.082589 | NaN | NaN | NaN | NaN | -0.238414 | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
5 | -1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.155288 | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
20 rows × 26 columns
[43]:
# EPR
EPRX = geomCalc.EPR(form = 'xarray', p = p).unstack().sel({'R-p':0}).drop('R-p')
EPRXresort = EPRX.squeeze(['l','lp']).drop(['l','lp']) # Safe squeeze & drop of selected singleton dims only.
if phaseCons['mfblmCons']['negRcoordSwap']:
EPRXresort['R'] *= -1
# Lambda
RX = ep.setPolGeoms()
# *** Lambda term
lambdaTerm, lambdaTable, lambdaD, _ = geomCalc.MFproj(RX = RX, form = 'xarray', phaseConvention = phaseConvention)
# lambdaTermResort = lambdaTerm.squeeze().drop('l').drop('lp') # This removes photon (l,lp) dims fully.
# lambdaTermResort = lambdaTerm.squeeze(['l','lp']).drop(['l','lp']) # Safe squeeze & drop of selected singleton dims only.
lambdaTermResort = lambdaTerm.squeeze(['l','lp']).drop(['l','lp']).sel({'Labels':'z'}).sum('R') # Safe squeeze & drop of selected singleton dims only, select (0,0,0) term only for pol. geometry.
# AF term
AKQS = ep.setADMs() # If not passed, set to defaults - A(0,0,0)=1 term only, i.e. isotropic distribution.
AFterm, DeltaKQS = geomCalc.deltaLMKQS(EPRXresort, AKQS, phaseConvention = phaseConvention)
[44]:
EPRX.unstack().coords
[44]:
Coordinates:
* l (l) int64 1
* lp (lp) int64 1
* P (P) int64 0 1 2
* p (p) int64 0
* R (R) int64 -1 0 1
[45]:
lambdaTerm.dims
[45]:
('Rp', 'l', 'lp', 'P', 'mu', 'mup', 'R', 'Labels')
[46]:
lambdaTermResort.dims
[46]:
('Rp', 'P', 'mu', 'mup')
[47]:
EPRX.squeeze()
[47]:
<xarray.DataArray (P: 3, R: 3)> array([[ nan, -0.57735027, nan], [ nan, nan, nan], [ nan, 0.81649658, nan]]) Coordinates: l int64 1 lp int64 1 * P (P) int64 0 1 2 p int64 0 * R (R) int64 -1 0 1 Attributes: dataType: EPR phaseCons: {'phaseConvention': 'S', 'genMatEcons': {'negm': False}, 'EPR...
[48]:
plotDimsRed = ['l', 'p', 'lp']
xDim = {'PR':['P','R']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(EPRXresort, plotDims=plotDimsRed, xDim=xDim, pType = 'r', squeeze = False)
daPlotpd
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[48]:
P | 0 | 2 |
---|---|---|
R | 0 | 0 |
p | ||
0 | -0.57735 | 0.816497 |
[49]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LSRp':['P','Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(lambdaTerm.sel({'Labels':'z'}), xDim=xDim, pType = 'r', squeeze = False)
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(lambdaTermResort, xDim=xDim, pType = 'r', squeeze = False)
daPlotpd[0:20]
Set dataType (No dataType)
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[49]:
P | 0 | 1 | 2 | |||||||
---|---|---|---|---|---|---|---|---|---|---|
Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | |
mu | mup | |||||||||
-1 | -1 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.447214 |
0 | NaN | NaN | NaN | -0.408248 | NaN | NaN | NaN | 0.316228 | NaN | |
1 | 0.57735 | NaN | -0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN | |
0 | -1 | NaN | NaN | NaN | 0.408248 | NaN | NaN | NaN | 0.316228 | NaN |
0 | -0.57735 | NaN | NaN | NaN | NaN | NaN | 0.365148 | NaN | NaN | |
1 | NaN | -0.408248 | NaN | NaN | NaN | 0.316228 | NaN | NaN | NaN | |
1 | -1 | 0.57735 | NaN | 0.408248 | NaN | NaN | NaN | 0.182574 | NaN | NaN |
0 | NaN | 0.408248 | NaN | NaN | NaN | 0.316228 | NaN | NaN | NaN | |
1 | NaN | NaN | NaN | NaN | 0.447214 | NaN | NaN | NaN | NaN |
[50]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LSRp':['L','S-Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(DeltaKQS, xDim=xDim, pType = 'r')
daPlotpd
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[50]:
L | 0 | 1 | 2 | ||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
S-Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | ||||||
K | M | P | Q | R | Rp | S | |||||||||
0 | -1 | 1 | 0 | 1 | -1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN |
0 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||||
1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||||
2 | 0 | 1 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | ||
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | |||||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | |||||
1 | 0 | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | |||||
2 | 0 | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | NaN | |||||
0 | 0 | 0 | 0 | 0 | 0 | 1.0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
1 | 0 | 0 | -1 | 0 | NaN | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | ||
0 | 0 | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||||
1 | 0 | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||||
2 | 0 | 0 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | ||
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | |||||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | |||||
1 | 0 | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | |||||
2 | 0 | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | NaN | |||||
1 | 1 | 0 | -1 | -1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | |
0 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||||
1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||||
2 | 0 | -1 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | ||
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | |||||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | |||||
1 | 0 | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | |||||
2 | 0 | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | NaN |
[51]:
# plotDimsRed = ['l', 'm', 'lp', 'mp']
xDim = {'LSRp':['L','S-Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(AFterm, xDim=xDim, pType = 'r')
daPlotpd
Set dataType (No dataType)
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[51]:
L | 0 | 1 | 2 | ||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
S-Rp | 0 | -1 | 0 | 1 | -2 | -1 | 0 | 1 | 2 | ||||
M | P | R | Rp | t | |||||||||
-1 | 1 | 1 | -1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN |
0 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
2 | 1 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | |
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | |||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | |||
1 | 0 | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | |||
2 | 0 | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | NaN | |||
0 | 0 | 0 | 0 | 0 | 1.0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
1 | 0 | -1 | 0 | NaN | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | |
0 | 0 | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | 0 | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
2 | 0 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | |
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | |||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | |||
1 | 0 | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | |||
2 | 0 | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | NaN | |||
1 | 1 | -1 | -1 | 0 | NaN | NaN | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN |
0 | 0 | NaN | NaN | -0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | |||
1 | 0 | NaN | 0.333333 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||
2 | -1 | -2 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | |
-1 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | |||
0 | 0 | NaN | NaN | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | |||
1 | 0 | NaN | NaN | NaN | NaN | NaN | 0.2 | NaN | NaN | NaN | |||
2 | 0 | NaN | NaN | NaN | NaN | -0.2 | NaN | NaN | NaN | NaN |
[52]:
#*** Products
# Matrix element pair-wise multiplication by (l,m,mu) dims
matEconj = matEthres.copy().conj()
# matEconj = matEconj.unstack().rename({'l':'lp','m':'mp','mu':'mup'}) # Full unstack
# matEmult = matEthres.unstack() * matEconj
matEconj = matEconj.unstack('LM').rename({'l':'lp','m':'mp','mu':'mup'}) # Unstack LM only.
matEmult = matEthres.unstack('LM') * matEconj
matEmult.attrs['dataType'] = 'multTest'
# Threshold product and drop dims.
# matEmult = ep.util.matEleSelector(matEmult, thres = thres, dims = thresDims)
matEmult = matEleSelector(matEmult, thres = thres, dims = thresDims)
# Apply additional phase conventions?
if phaseCons['afblmCons']['llpPhase']:
matEmult *= np.power(-1, np.abs(matEmult.l - matEmult.lp))
# matEmult['mp'] *= -1 # Quick test of mp phase
# Product terms with similar dims
BLMprod = matEmult * BLMtableResort # Unstacked case with phase correction - THIS DROPS SYM TERMS? Takes intersection of das - http://xarray.pydata.org/en/stable/computation.html#automatic-alignment
# polProd = (EPRXresort * lambdaTermResort).sum(sumDimsPol) # Sum polarization terms here to keep total dims minimal in product. Here dims = (mu,mup,Euler/Labels)
# polProd = (EPRXresort * lambdaTermResort) # Without polarization terms sum to allow for mupPhase below (reqs. p)
# Test with alignment term
# polProd = (EPRXresort * lambdaTermResort * AFterm)
# With additional phase change(s)
lTermTest = lambdaTermResort.copy()
AFtermTest = AFterm.copy()
# AFtermTest['Rp'] *= -1
# lTermTest['mup'] *=-1
polProd = (EPRXresort * lTermTest * AFtermTest)
# polProd = (lTermTest * AFtermTest) # 25/06/20 this multiplication seems to restrict S-Rp to -1...+1, not sure why. mu/mup phase? Rp phase? Neither (nor both) seems to fix the issue.
# Set additional phase term, (-1)^(mup-p) **** THIS MIGHT BE SPURIOUS FOR GENERAL EPR TENSOR CASE??? Not sure... but definitely won't work if p summed over above!
if phaseCons['mfblmCons']['mupPhase']:
mupPhaseTerm = np.power(-1, np.abs(polProd.mup - polProd.p))
polProd *= mupPhaseTerm
# Additional [P]^1/2 degen term, NOT included in EPR defn.
# Added 09/04/20
polProd *= np.sqrt(2*polProd.P+1)
[53]:
# Check terms...
# xDim = {'LM':['L','M']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(polProd, xDim=xDim, pType = 'r') #, thres=None)
daPlotpd
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[53]:
L | 0 | 2 | |||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
S-Rp | 0 | -2 | -1 | 0 | 1 | 2 | |||||||
M | P | R | Rp | mu | mup | p | t | ||||||
0 | 0 | 0 | 0 | -1 | 1 | 0 | 0 | 0.333333 | NaN | NaN | NaN | NaN | NaN |
0 | 0 | 0 | 0 | 0.333333 | NaN | NaN | NaN | NaN | NaN | ||||
1 | -1 | 0 | 0 | 0.333333 | NaN | NaN | NaN | NaN | NaN | ||||
2 | 0 | -2 | 1 | 1 | 0 | 0 | NaN | NaN | NaN | NaN | NaN | -0.163299 | |
-1 | 0 | 1 | 0 | 0 | NaN | NaN | NaN | NaN | 0.11547 | NaN | |||
1 | 0 | 0 | 0 | NaN | NaN | NaN | NaN | -0.11547 | NaN | ||||
0 | -1 | 1 | 0 | 0 | NaN | NaN | NaN | -0.066667 | NaN | NaN | |||
0 | 0 | 0 | 0 | NaN | NaN | NaN | 0.133333 | NaN | NaN | ||||
1 | -1 | 0 | 0 | NaN | NaN | NaN | -0.066667 | NaN | NaN | ||||
1 | -1 | 0 | 0 | 0 | NaN | NaN | -0.11547 | NaN | NaN | NaN | |||
0 | -1 | 0 | 0 | NaN | NaN | 0.11547 | NaN | NaN | NaN | ||||
2 | -1 | -1 | 0 | 0 | NaN | -0.163299 | NaN | NaN | NaN | NaN |
[54]:
polProd = polProd.sum(sumDimsPol)
# polProd = matEleSelector(polProd, thres = thres) # Select over dims for reduction.
# Test big mult...
# mTerm = polProd.sel({'R':0,'Labels':'z'}) * BLMprod.sum(['Total']) # With selection of z geom. # BLMprod.sum(['Cont', 'Targ', 'Total'])
# mTerm = polProd.sel({'R':0}) * BLMprod # BLMprod.sum(['Cont', 'Targ', 'Total'])
mTerm = polProd * BLMprod
[55]:
# Check terms... BLMprod
xDim = {'LSRp':['L','S-Rp']}
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(testMult, plotDims=plotDimsRed, xDim=xDim, pType = 'r')
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(BLMprod.sel({'Eke':0.1}), xDim=xDim, pType = 'r')
daPlotpd
Set dataType (No dataType)
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\core\nputils.py:215: RuntimeWarning:
All-NaN slice encountered
Plotting data (No filename), pType=r, thres=0.01, with Seaborn
[55]:
L | 0 | 2 | 4 | 6 | ||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
S-Rp | 0 | -2 | 0 | 2 | -2 | 0 | 2 | -2 | 0 | 2 | ||||||||
Cont | Targ | Total | l | lp | m | mp | mu | mup | ||||||||||
PU | SG | PU | 1 | 1 | -1 | -1 | 1 | 1 | NaN | NaN | NaN | -1.098886 | NaN | NaN | NaN | NaN | NaN | NaN |
1 | 1 | -1 | 1.003141 | NaN | -0.448618 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||||
1 | -1 | -1 | 1 | 1.003141 | NaN | -0.448618 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |||||
1 | -1 | -1 | NaN | -1.098886 | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ||||||
3 | -1 | -1 | 1 | 1 | NaN | NaN | NaN | 0.020515 | NaN | NaN | -0.059223 | NaN | NaN | NaN | ||||
1 | 1 | -1 | NaN | NaN | 0.050252 | NaN | NaN | -0.037456 | NaN | NaN | NaN | NaN | ||||||
1 | -1 | -1 | 1 | NaN | NaN | 0.050252 | NaN | NaN | -0.037456 | NaN | NaN | NaN | NaN | |||||
1 | -1 | -1 | NaN | 0.020515 | NaN | NaN | -0.059223 | NaN | NaN | NaN | NaN | NaN | ||||||
3 | 1 | -1 | -1 | 1 | 1 | NaN | NaN | NaN | 0.020515 | NaN | NaN | -0.059223 | NaN | NaN | NaN | |||
1 | 1 | -1 | NaN | NaN | 0.050252 | NaN | NaN | -0.037456 | NaN | NaN | NaN | NaN | ||||||
1 | -1 | -1 | 1 | NaN | NaN | 0.050252 | NaN | NaN | -0.037456 | NaN | NaN | NaN | NaN | |||||
1 | -1 | -1 | NaN | 0.020515 | NaN | NaN | -0.059223 | NaN | NaN | NaN | NaN | NaN | ||||||
3 | -1 | -1 | 1 | 1 | NaN | NaN | NaN | -0.075872 | NaN | NaN | -0.059734 | NaN | NaN | -0.089473 | ||||
1 | 1 | -1 | 0.103893 | NaN | 0.046462 | NaN | NaN | NaN | NaN | NaN | -0.065488 | NaN | ||||||
1 | -1 | -1 | 1 | 0.103893 | NaN | 0.046462 | NaN | NaN | NaN | NaN | NaN | -0.065488 | NaN | |||||
1 | -1 | -1 | NaN | -0.075872 | NaN | NaN | -0.059734 | NaN | NaN | -0.089473 | NaN | NaN | ||||||
SU | SG | SU | 1 | 1 | 0 | 0 | 0 | 0 | 2.114595 | NaN | 1.891351 | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
3 | 0 | 0 | 0 | 0 | NaN | NaN | -0.098066 | NaN | NaN | -0.097458 | NaN | NaN | NaN | NaN | ||||
3 | 1 | 0 | 0 | 0 | 0 | NaN | NaN | -0.098066 | NaN | NaN | -0.097458 | NaN | NaN | NaN | NaN | |||
3 | 0 | 0 | 0 | 0 | 0.186727 | NaN | 0.111343 | NaN | NaN | 0.101851 | NaN | NaN | 0.156936 | NaN |
[56]:
xDim = {'LM':['L','M']}
mTermSum = mTerm.sum(sumDims)
mTermSum.attrs['dataType'] = 'matE' # Set matE here to allow for correct plotting of sym dims.
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermSum, xDim='Eke', sumDims=None, pType = 'r', thres = None, fillna = True, SFflag=False)
daPlotpd[0:50]
Plotting data (No filename), pType=r, thres=None, with Seaborn
No handles with labels found to put in legend.
[56]:
Eke | 0.1 | 1.1 | 2.1 | 3.1 | 4.1 | 5.1 | 6.1 | 7.1 | 8.1 | 9.1 | ... | 41.1 | 42.1 | 43.1 | 44.1 | 45.1 | 46.1 | 47.1 | 48.1 | 49.1 | 50.1 | ||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Cont | L | M | S-Rp | Targ | Total | t | |||||||||||||||||||||
PU | 0 | -1 | -2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 |
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | -2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ||
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.738022 | 0.720190 | 0.713936 | 0.715632 | 0.723239 | 0.735319 | 0.750596 | 0.767824 | 0.785799 | 0.803414 | ... | 0.389743 | 0.371662 | 0.354243 | 0.337490 | 0.321402 | 0.305977 | 0.291206 | 0.277080 | 0.263586 | 0.250712 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | -2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ||
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | -1 | -2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | -2 | SG | PU | 0 | 0.185154 | 0.207022 | 0.214972 | 0.218305 | 0.219175 | 0.218591 | 0.217069 | 0.214871 | 0.212125 | 0.208893 | ... | 0.061546 | 0.058851 | 0.056280 | 0.053827 | 0.051487 | 0.049257 | 0.047132 | 0.045107 | 0.043178 | 0.041342 | ||
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.040253 | 0.091965 | 0.110662 | 0.116569 | 0.114907 | 0.108336 | 0.098461 | 0.086362 | 0.072831 | 0.058489 | ... | -0.051207 | -0.048509 | -0.045862 | -0.043276 | -0.040761 | -0.038323 | -0.035968 | -0.033701 | -0.031523 | -0.029438 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.185154 | 0.207022 | 0.214972 | 0.218305 | 0.219175 | 0.218591 | 0.217069 | 0.214871 | 0.212125 | 0.208893 | ... | 0.061546 | 0.058851 | 0.056280 | 0.053827 | 0.051487 | 0.049257 | 0.047132 | 0.045107 | 0.043178 | 0.041342 | |||
1 | -2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ||
-1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | PU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
SU | 0 | -1 | -2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 |
-1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | -2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ||
-1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | SU | 0 | 0.767107 | 0.794111 | 0.835793 | 0.889138 | 0.952871 | 1.026468 | 1.109505 | 1.201144 | 1.299572 | 1.401236 | ... | 0.022910 | 0.023531 | 0.024251 | 0.025035 | 0.025854 | 0.026686 | 0.027512 | 0.028319 | 0.029095 | 0.029830 | |||
1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | -2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ||
-1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | -1 | -2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |
-1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
0 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
1 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | |||
2 | SG | SU | 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | ... | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 0.000000 |
50 rows × 51 columns
[57]:
# mTermSum.XS.real.squeeze().plot.line(x='Eke', col='Sym');
# ep.util.matEleSelector(mTermSum, thres = None, dims='Eke').real.squeeze().plot.line(x='Eke', col='Sym');
ep.util.matEleSelector(mTermSum.sel({'t':0}), thres = 0.01, dims='Eke').real.squeeze().plot.line(x='Eke', row='S-Rp', col='Sym');
[58]:
# mTermSum.XS.real.squeeze().plot.line(x='Eke', col='Sym');
# ep.util.matEleSelector(mTermSum, thres = None, dims='Eke').real.squeeze().plot.line(x='Eke', col='Sym');
ep.util.matEleSelector(mTermSum.sum('S-Rp'), thres = 0.01, dims='Eke').real.squeeze().plot.line(x='Eke', col='Sym');
Testing functional form
[59]:
phaseConvention = 'E' # Set phase conventions used in the numerics - for ePolyScat matrix elements, set to 'E', to match defns. above.
symSum = False # Sum over symmetry groups, or keep separate?
SFflag = False # Include scaling factor to Mb in calculation?
thres = 1e-4
RX = ep.setPolGeoms() # Set default pol geoms (z,x,y), or will be set by mfblmXprod() defaults - FOR AF case this is only used to set 'z' geom for unity wigner D's - should rationalise this!
start = time.time()
mTermST, mTermS, mTermTest, BetaNormX = ep.geomFunc.afblmXprod(dataSet[0], QNs = None, RX=RX, thres = thres, selDims = {'it':1, 'Type':'L'}, thresDims='Eke', symSum=symSum, phaseConvention=phaseConvention,
SFflag=False, SFflagRenorm = False, BLMRenorm = 1,sumDims = ['mu', 'mup', 'l','lp','m','mp','S-Rp'])
end = time.time()
print('Elapsed time = {0} seconds, for {1} energy points, {2} polarizations, threshold={3}.'.format((end-start), mTermST.Eke.size, RX.size, thres))
# Elapsed time = 3.3885273933410645 seconds, for 51 energy points, 3 polarizations, threshold=0.01.
# Elapsed time = 5.059587478637695 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
Elapsed time = 2.4853365421295166 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
[60]:
# Full results (before summation)
mTermST.attrs['dataType'] = 'matE' # Set matE here to allow for correct plotting of sym dims.
plotDimsRed = ['Labels','L','M'] # Set plotDims to fix dim ordering in plot
if not symSum:
plotDimsRed.extend(['Cont','Targ','Total'])
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, plotDims=plotDimsRed, xDim='Eke', sumDims=None, pType = 'r', thres = None, fillna = True, SFflag=False)
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True) # If plotDims is not passed use default ordering.
Plotting data n2_3sg_0.1-50.1eV_A2.inp.out, pType=r, thres=None, with Seaborn
No handles with labels found to put in legend.
[61]:
mTermST.XSraw.real.squeeze().plot.line(x='Eke', col='Sym');
(mTermST.XSrescaled).real.squeeze().plot.line(x='Eke', col='Sym');
(mTermST.XSiso).real.squeeze().plot.line(x='Eke', col='Sym');
[62]:
# Plot cross sections using Xarray functionality
dataXS[0].sel({'Type':'L', 'XC':'SIGMA'}).plot.line(x='Eke');
# (2*mTermST.XSrescaled).real.squeeze().plot.line(x='Eke', linestyle='dashed');
(mTermST.XSrescaled).pipe(np.abs).squeeze().plot.line(x='Eke', linestyle='dashed');
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\plot\plot.py:86: FutureWarning:
This DataArray contains multi-dimensional coordinates. In the future, these coordinates will be transposed as well unless you specify transpose_coords=False.
[63]:
# Test normalisation as per CG version (see lfblmGeom.py) - NOW SET IN FUNCTION
# XSmatE = (matEthres * matEthres.conj()).sum(['LM','mu']) # Use selected & thresholded matE.
# # NOTE - this may fail in current form if dims are missing.
# normBeta = 3/5 * (1/XSmatE) # Normalise by sum over matrix elements squared.
# normBeta.sel({'Cont':'PU'}).real.plot.line(x='Eke')
# ep.util.matEleSelector(mTermST*normBeta, thres = None, dims='Eke').real.squeeze().plot.line(x='Eke', col='Sym');
[64]:
# Try converting to BL values...
BLvals = ep.util.conversion.conv_BL_BLM(mTermST, to='lg')
ep.util.matEleSelector(BLvals, thres = None, dims='Eke').real.squeeze().plot.line(x='Eke', col='Sym');
plt.figure();
dataXS[0].sel({'Type':'L', 'XC':'BETA'}).plot.line(x='Eke', linestyle = 'dashed'); # Ref. values from GetCro
[65]:
(BLvals.XSrescaled).pipe(np.abs).squeeze().plot.line(x='Eke', linestyle='dashed');
C:\Users\femtolab\.conda\envs\ePSdev\lib\site-packages\xarray\plot\plot.py:86: FutureWarning:
This DataArray contains multi-dimensional coordinates. In the future, these coordinates will be transposed as well unless you specify transpose_coords=False.
Try sym summation…
[66]:
phaseConvention = 'E' # Set phase conventions used in the numerics - for ePolyScat matrix elements, set to 'E', to match defns. above.
symSum = True # Sum over symmetry groups, or keep separate?
SFflag = False # Include scaling factor to Mb in calculation?
thres = 1e-4
RX = ep.setPolGeoms() # Set default pol geoms (z,x,y), or will be set by mfblmXprod() defaults - FOR AF case this is only used to set 'z' geom for unity wigner D's - should rationalise this!
start = time.time()
mTermST, mTermS, mTermTest, BetaNormX = ep.geomFunc.afblmXprod(dataSet[0], QNs = None, RX=RX, thres = thres, selDims = {'it':1, 'Type':'L'}, thresDims='Eke', symSum=symSum, SFflag=SFflag, phaseConvention=phaseConvention)
end = time.time()
print('Elapsed time = {0} seconds, for {1} energy points, {2} polarizations, threshold={3}.'.format((end-start), mTermST.Eke.size, RX.size, thres))
# Elapsed time = 3.3885273933410645 seconds, for 51 energy points, 3 polarizations, threshold=0.01.
# Elapsed time = 5.059587478637695 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
Elapsed time = 2.3291077613830566 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
[67]:
# Full results (before summation)
mTermST.attrs['dataType'] = 'matE' # Set matE here to allow for correct plotting of sym dims.
plotDimsRed = ['Labels','L','M'] # Set plotDims to fix dim ordering in plot
if not symSum:
plotDimsRed.extend(['Cont','Targ','Total'])
daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, plotDims=plotDimsRed, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True, SFflag=False)
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True) # If plotDims is not passed use default ordering.
Plotting data n2_3sg_0.1-50.1eV_A2.inp.out, pType=r, thres=0.01, with Seaborn
No handles with labels found to put in legend.
[68]:
mTermST.XSraw.real.squeeze().plot.line(x='Eke');
mTermST.XSrescaled.real.squeeze().plot.line(x='Eke');
[69]:
# 09/09/20 - updating here for BetaNormX values...
plotThres = 1e-2
ep.util.matEleSelector(BetaNormX['XSraw'], thres = plotThres, dims='Eke', sq=True, drop=True).real.plot.line(x='Eke'); # Values from AF code - note these are small, with complex terms
ep.util.matEleSelector(BetaNormX['XSrescaled'], thres = plotThres, dims='Eke', sq=True, drop=True).real.plot.line(x='Eke'); # Values from AF code - real valued
ep.util.matEleSelector(BetaNormX['XSiso'], thres = plotThres, dims='Eke', sq=True, drop=True).real.plot.line(x='Eke'); # Values from sum over matrix elements
# Ref. values from LF code
# ep.util.matEleSelector(BetaNormX_sph['XS'], thres = plotThres, dims='Eke', sq=True, drop=True).real.plot.line(x='Eke', linestyle='dashed'); # Values from LF code.
[70]:
ep.util.matEleSelector(mTermST, thres = 0.1, dims='Eke').real.squeeze().plot.line(x='Eke');
First attempt… with sym summation
Looks quite different - might be SF issue with combining continua, and/or degen factor…?
Looks same as previous AF code attempt (http://localhost:8888/lab/tree/dev/ePSproc/ePSproc_AFBLM_calcs_bench_100220.ipynb), suggesting something in degen/normalisation in formalism amiss?
EDIT 10/09/20: now fixed, issues with normalisation and Spherical Harmonic and Legendre function defintions - see LF/AF test notebook for details
[71]:
# ep.util.matEleSelector(mTermST, thres = 0.1, dims='Eke').real.squeeze().plot.line(x='Eke');
ep.util.matEleSelector(mTermST, thres = 0.1, dims='Eke').imag.squeeze().plot.line(x='Eke');
Test compared to experimental N2 AF results…
In this case MF PADs look good, so also expect good agreement with AF…
[72]:
# Adapted from ePSproc_AFBLM_testing_010519_300719.m
# Load ADMs for N2
from scipy.io import loadmat
ADMdataFile = os.path.join(modPath, 'data', 'alignment', 'N2_ADM_VM_290816.mat')
ADMs = loadmat(ADMdataFile)
[73]:
# Set tOffset for calcs, 3.76ps!!!
# This is because this is 2-pulse case, and will set t=0 to 2nd pulse (and matches defn. in N2 experimental paper)
tOffset = -3.76
ADMs['time'] = ADMs['time'] + tOffset
[74]:
# Plot
import matplotlib.pyplot as plt
plt.plot(ADMs['time'].T, np.real(ADMs['ADM'].T))
plt.legend(ADMs['ADMlist'])
plt.xlabel('t/ps')
plt.ylabel('ADM')
plt.show()
[74]:
[<matplotlib.lines.Line2D at 0x1da183b19e8>,
<matplotlib.lines.Line2D at 0x1da183b17f0>,
<matplotlib.lines.Line2D at 0x1da183b12e8>,
<matplotlib.lines.Line2D at 0x1da183b1668>]
[74]:
<matplotlib.legend.Legend at 0x1da1838a1d0>
[74]:
Text(0.5, 0, 't/ps')
[74]:
Text(0, 0.5, 'ADM')
[75]:
# Selection & downsampling
trange=[1, 9.5] # Set range in ps for calc
tStep=5 # Set tStep for downsampling
tMask = (ADMs['time']>trange[0]) * (ADMs['time']<trange[1])
ind = np.nonzero(tMask)[1][0::tStep]
At = ADMs['time'][:,ind].squeeze()
ADMin = ADMs['ADM'][:,ind]
print(f"Selecting {ind.size} points")
plt.plot(At, np.real(ADMin.T))
plt.legend(ADMs['ADMlist'])
plt.xlabel('t/ps')
plt.ylabel('ADM')
plt.show()
# Set in Xarray
ADMX = ep.setADMs(ADMs = ADMs['ADM'][:,ind], t=At, KQSLabels = ADMs['ADMlist'], addS = True)
# ADMX
Selecting 87 points
[75]:
[<matplotlib.lines.Line2D at 0x1da19862518>,
<matplotlib.lines.Line2D at 0x1da198622b0>,
<matplotlib.lines.Line2D at 0x1da19862828>,
<matplotlib.lines.Line2D at 0x1da19862748>]
[75]:
<matplotlib.legend.Legend at 0x1da19eb1630>
[75]:
Text(0.5, 0, 't/ps')
[75]:
Text(0, 0.5, 'ADM')
[76]:
# Run with sym summation...
phaseConvention = 'E' # Set phase conventions used in the numerics - for ePolyScat matrix elements, set to 'E', to match defns. above.
symSum = True # Sum over symmetry groups, or keep separate?
SFflag = False # Include scaling factor to Mb in calculation?
SFflagRenorm = False # Renorm terms
BLMRenorm = 1
thres = 1e-4
RX = ep.setPolGeoms() # Set default pol geoms (z,x,y), or will be set by mfblmXprod() defaults - FOR AF case this is only used to set 'z' geom for unity wigner D's - should rationalise this!
start = time.time()
mTermST, mTermS, mTermTest, BetaNormX = ep.geomFunc.afblmXprod(dataSet[0], QNs = None, AKQS=ADMX, RX=RX, thres = thres, selDims = {'it':1, 'Type':'L'}, thresDims='Eke',
symSum=symSum, SFflag=SFflag, BLMRenorm = BLMRenorm,
phaseConvention=phaseConvention)
end = time.time()
print('Elapsed time = {0} seconds, for {1} energy points, {2} polarizations, threshold={3}.'.format((end-start), mTermST.Eke.size, RX.size, thres))
# Elapsed time = 3.3885273933410645 seconds, for 51 energy points, 3 polarizations, threshold=0.01.
# Elapsed time = 5.059587478637695 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
Elapsed time = 6.940988779067993 seconds, for 51 energy points, 3 polarizations, threshold=0.0001.
[77]:
# Full results (before summation)
mTermST.attrs['dataType'] = 'matE' # Set matE here to allow for correct plotting of sym dims.
plotDimsRed = ['Labels','L','M'] # Set plotDims to fix dim ordering in plot
if not symSum:
plotDimsRed.extend(['Cont','Targ','Total'])
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, plotDims=plotDimsRed, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True, SFflag=False)
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True) # If plotDims is not passed use default ordering.
[78]:
# daPlot, daPlotpd, legendList, gFig = ep.lmPlot(mTermST.XS, plotDims=plotDimsRed, xDim='Eke', sumDims=None, pType = 'r', thres = 0.01, fillna = True, SFflag=False)
[79]:
mTermST.XSrescaled.real.plot()
mTermST.where(mTermST.L > 0).real.plot(col='LM')
[79]:
<matplotlib.collections.QuadMesh at 0x1da190f54a8>
[79]:
<xarray.plot.facetgrid.FacetGrid at 0x1da1cfca048>
[80]:
# Plot single energy XS
E = 8.1
mTermST.XSraw.sel({'Eke':E}).real.squeeze().plot.line(x='t');
# mTermST.XSrescaled.sel({'Eke':E}).real.squeeze().plot.line(x='t');
[81]:
# LM with renorm
mTermST.sel({'Eke':E}).where(mTermST['L']>0).real.squeeze().plot.line(x='t');
# Lineshapes look good, but abs. values are low compared to expt. (see below)
[82]:
# # With further renorm by B00(t), BLMRenorm =2 or 3
# # This looks pretty good for shape... but about x2 lower than ref. experimental values????
# # 01/09/20: now looking a bit better with updated renorm routine (BLMRenorm =3 case)
# print(BLMRenorm)
# ep.util.matEleSelector(mTermST, thres = 1e-4, dims='t').sel({'Eke':E}).real.squeeze().plot.line(x='t');
# plt.figure()
# # ep.util.matEleSelector((mTermST/mTermST.XS), thres = 1e-4, dims='t').sel({'Eke':E}).where(mTermST.L > 0).real.squeeze().plot.line(x='t');
# ep.util.matEleSelector((mTermST/mTermST.XSraw), thres = 1e-4, dims='t').sel({'Eke':E}).real.squeeze().plot.line(x='t');
# ep.util.matEleSelector((mTermST/mTermST.XS), thres = 1e-4, dims='t').sel({'Eke':E}).pipe(np.abs).squeeze().plot.line(x='t', linestyle="dashed");
# plt.figure()
# ep.util.matEleSelector(mTermST/mTermST.sel({'L':0,'M':0}), thres = 1e-4, dims='t').sel({'Eke':E}).real.squeeze().plot.line(x='t');
[83]:
# Check Lg converted values... does this fix possible magnitude issue?
ep.util.matEleSelector(ep.conversion.conv_BL_BLM(mTermST, to='lg'), thres = 1e-6, dims='t').sel({'Eke':E}).where(mTermST['L']>0).real.squeeze().plot.line(x='t');
# ep.util.matEleSelector(ep.conversion.conv_BL_BLM(mTermST/mTermST.sel({'L':0,'M':0}), to='lg'), thres = 1e-6, dims='t').sel({'Eke':E}).where(mTermST['L']>0).real.squeeze().plot.line(x='t');
01/09/20: FEELS like progress, but still not quite correct - values a bit small (Sph), or a bit large (Lg) - might be an issue with experimental renorm, since B00 ~ 1.2?
Experimental results
(Fig 2 from Marceau, Claude, Varun Makhija, Dominique Platzer, A. Yu. Naumov, P. B. Corkum, Albert Stolow, D. M. Villeneuve, and Paul Hockett. “Molecular Frame Reconstruction Using Time-Domain Photoionization Interferometry.” Physical Review Letters 119, no. 8 (August 2017): 083401. https://doi.org/10.1103/PhysRevLett.119.083401. Full data at https://doi.org/10.6084/m9.figshare.4480349.v10)
Currently:
Relative signs/phases of betas look good at 4.1eV… but not at 7.1 or 8.1eV! (Experimental MFPAD comparison energy) - suggests still phase and/or renorm issue somewhere?
Values still off, same issue probably.
Might be phase issue with sum over rotation matrix elements, which are conj() in MF case, but not in original derivation here…? This could give 3j sign-flip, hence change calculated alignment response.
UPDATE 07/07/20 (after LF testing and degen partial fix):
Now looking reasonable at 7.1eV by shape, but off at 8.1eV (i.e. transition zone has shifted up in E).
Still have some issue(s) with abs. values, and XS vs. E. This is the case if SF is included, and also if not!
Version info
[84]:
import scooby
scooby.Report(additional=['epsproc', 'xarray', 'holoviews'])
[84]:
Thu Sep 10 18:31:17 2020 Eastern Daylight Time | |||||
OS | Windows | CPU(s) | 32 | Machine | AMD64 |
Architecture | 64bit | RAM | 63.9 GB | Environment | Jupyter |
Python 3.7.3 (default, Apr 24 2019, 15:29:51) [MSC v.1915 64 bit (AMD64)] | |||||
epsproc | 1.2.5-dev | xarray | 0.15.0 | holoviews | 1.13.3 |
numpy | 1.19.1 | scipy | 1.3.0 | IPython | 7.12.0 |
matplotlib | 3.3.1 | scooby | 0.5.6 | ||
Intel(R) Math Kernel Library Version 2020.0.0 Product Build 20191125 for Intel(R) 64 architecture applications |