spaudiopy.sph
Spherical Harmonics.
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['axes.grid'] = True
import spaudiopy as spa
spa.plot.sh_coeffs_subplot([np.sqrt(4*np.pi) * np.array([1, 0, 0, 0]),
np.sqrt(4/3*np.pi) * np.array([0, 1, 0, 0]),
np.sqrt(4/3*np.pi) * np.array([0, 0, 1, 0]),
np.sqrt(4/3*np.pi) * np.array([0, 0, 0, 1])],
titles=["$Y_{0, 0}$", "$Y_{1, -1}$",
"$Y_{1, 0}$", "$Y_{1, 1}$"])
Functions
|
Convert B-format input to first order SH signals. |
|
Order N spatially bandlimited Dirac pulse at central angle d. |
|
Spectral equalization gain G(kr)|N for diffuse field of order N_sph. |
|
Modal weights for spatial butterworth filter / beamformer. |
|
Modal weights for beamformer resulting in a cardioid. |
|
Check if condition number for a least-squares SHT(N_sph) is high. |
|
Design spherical filter bank analysis and reconstruction matrix. |
|
Modal weights for beamformer resulting in a hyper-cardioid. |
|
Perform the inverse spherical harmonics transform. |
|
Return max-rE modal weight coefficients for spherical harmonics order N. |
|
Modal weights for beamformer resulting with max-rE weighting. |
|
Mode strength b_n(kr) for an incident plane wave on sphere. |
|
Convert N3D (orthonormal) to SN3D (Schmidt semi-normalized) signals. |
|
Return coordinates of shape='tetrahedron' only, yet. |
|
Calculate the diffuse field pressure frequency response of a spherical scatterer, up to SH order N. |
|
Little helper that projects x, y, z onto unit sphere. |
|
Calculate r_E vector and magnitude from loudspeaker gains. |
Repeat each coefficient in 'c' m times per spherical order n. |
|
|
Rotate spherical harmonics coefficients. |
|
Evaluates the spherical harmonics up to order N_sph for given angles. |
|
Multiply SH vector a_nm and b_nm in (discrete) space. |
|
Computes a Wigner-D matrix for the rotation of spherical harmonics. |
|
Convert first order SH input signals to B-format. |
|
Spherical harmonics transform of f for appropriate point sets. |
|
Spherical harmonics transform of f as least-squares solution. |
|
Convert SN3D (Schmidt semi-normalized) to N3D (orthonormal) signals. |
|
Convert soundfield tetraeder mic input signals to B-format by SHT. |
|
Reconstruction factor for restoring amplitude/energy preservation. |
|
Spherical Hankel function of the second kind. |
|
Get B format signal channels for source in direction azi/zen. |
|
Source signal(s) plane wave encoded in spherical harmonics. |
|
Make modal weighting / tapering unit amplitude in steering direction. |
- spaudiopy.sph.sh_matrix(N_sph, azi, zen, sh_type='real')[source]
Evaluates the spherical harmonics up to order N_sph for given angles.
Matrix returns spherical harmonics up to order \(N\) evaluated at the given angles/grid.
\[\begin{split}\mathbf{Y} = \left[ \begin{array}{ccccc} Y_0^0(\theta[0], \phi[0]) & Y_1^{-1}(\theta[0], \phi[0]) & Y_1^0(\theta[0], \phi[0]) & \dots & Y_N^N(\theta[0], \phi[0]) \\ Y_0^0(\theta[1], \phi[1]) & Y_1^{-1}(\theta[1], \phi[1]) & Y_1^0(\theta[1], \phi[1]) & \dots & Y_N^N(\theta[1], \phi[1]) \\ \vdots & \vdots & \vdots & \vdots & \vdots \\ Y_0^0(\theta[Q-1], \phi[Q-1]) & Y_1^{-1}(\theta[Q-1], \phi[Q-1]) & Y_1^0(\theta[Q-1], \phi[Q-1]) & \dots & Y_N^N(\theta[Q-1], \phi[Q-1]) \end{array} \right]\end{split}\]where
\[Y_n^m(\theta, \phi) = \sqrt{\frac{2n + 1}{4 \pi} \frac{(n-m)!}{(n+m)!}} P_n^m(\cos \theta) e^{i m \phi}\]When using sh_type=’real’, the real spherical harmonics \(Y_{n,m}(\theta, \phi)\) are implemented as a relation to \(Y_n^m(\theta, \phi)\).
- Parameters:
N_sph (int) – Maximum SH order.
azi ((Q,) array_like) – Azimuth.
zen ((Q,) array_like) – Colatitude / Zenith.
sh_type (‘complex’ or ‘real’ spherical harmonics.)
- Returns:
Ymn ((Q, (N+1)**2) numpy.ndarray) – Matrix of spherical harmonics.
Notes
The convention used here is also known as N3D-ACN (for sh_type=’real’).
- spaudiopy.sph.sht(f, N_sph, azi, zen, sh_type, weights=None, Y_nm=None)[source]
Spherical harmonics transform of f for appropriate point sets.
If f is a QxS matrix then the transform is applied to each column of f, and returns the coefficients at each column of F_nm respectively.
- Parameters:
f ((Q, S)) – The spherical function(S) evaluated at Q directions ‘azi/zen’.
N_sph (int) – Maximum SH order.
azi ((Q,) array_like) – Azimuth.
zen ((Q,) array_like) – Colatitude / Zenith.
sh_type (‘complex’ or ‘real’ spherical harmonics.)
weights ((Q,) array_like, optional) – Quadrature weights.
Y_nm ((Q, (N+1)**2) numpy.ndarray, optional) – Matrix of spherical harmonics.
- Returns:
F_nm (((N+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
- spaudiopy.sph.sht_lstsq(f, N_sph, azi, zen, sh_type, Y_nm=None)[source]
Spherical harmonics transform of f as least-squares solution.
If f is a QxS matrix then the transform is applied to each column of f, and returns the coefficients at each column of F_nm respectively.
- Parameters:
f ((Q, S)) – The spherical function(S) evaluated at Q directions ‘azi/zen’.
N_sph (int) – Maximum SH order.
azi ((Q,) array_like) – Azimuth.
zen ((Q,) array_like) – Colatitude / Zenith.
sh_type (‘complex’ or ‘real’ spherical harmonics.)
Y_nm ((Q, (N+1)**2) numpy.ndarray, optional) – Matrix of spherical harmonics.
- Returns:
F_nm (((N+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
- spaudiopy.sph.inverse_sht(F_nm, azi, zen, sh_type, N_sph=None, Y_nm=None)[source]
Perform the inverse spherical harmonics transform.
- Parameters:
F_nm (((N+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
azi ((Q,) array_like) – Azimuth.
zen ((Q,) array_like) – Colatitude / Zenith.
sh_type (‘complex’ or ‘real’ spherical harmonics.)
N_sph (int, optional) – Maximum SH order.
Y_nm ((Q, (N+1)**2) numpy.ndarray, optional) – Matrix of spherical harmonics.
- Returns:
f ((Q, S)) – The spherical function(S) evaluated at Q directions ‘azi/zen’.
- spaudiopy.sph.sh_rotation_matrix(N_sph, yaw, pitch, roll, sh_type='real', return_as_blocks=False)[source]
Computes a Wigner-D matrix for the rotation of spherical harmonics.
- Parameters:
N_sph (int) – Maximum SH order.
yaw (float) – Rotation around Z axis.
pitch (float) – Rotation around Y axis.
roll (float) – Rotation around X axis.
sh_type (‘complex’ or ‘real’ spherical harmonics.) – Currently only ‘real’ is supported.
return_as_blocks (optional, default is False.) – Return full block diagonal matrix, or a list of blocks if True.
- Returns:
R ((…, (N_sph+1)**2, (N_sph+1)**2) numpy.ndarray) – A block diagonal matrix R with shape (…, (N_sph+1)**2, (N_sph+1)**2), or a list r_blocks of numpy arrays [r(n) for n in range(N_sph)], where the shape of r is (…, 2*n-1, 2*n-1).
See also
spaudiopy.sph.rotate_sh()
Apply rotation to SH signals.
References
Implemented according to: Ivanic, Joseph, and Klaus Ruedenberg. “Rotation matrices for real spherical harmonics. Direct determination by recursion.” The Journal of Physical Chemistry 100.15 (1996): 6342-6347.
Ported from https://git.iem.at/audioplugins/IEMPluginSuite .
- spaudiopy.sph.rotate_sh(F_nm, yaw, pitch, roll, sh_type='real')[source]
Rotate spherical harmonics coefficients.
- Parameters:
F_nm ((…, (N_sph+1)**2) numpy.ndarray) – Spherical harmonics coefficients
yaw (float) – Rotation around Z axis.
pitch (float) – Rotation around Y axis.
roll (float) – Rotation around X axis.
sh_type (‘complex’ or ‘real’ spherical harmonics.) – Currently only ‘real’ is supported.
- Returns:
F_nm_rot ((…, (N_sph+1)**2) numpy.ndarray) – Rotated spherical harmonics coefficients.
Examples
N_sph = 5 y1 = spa.sph.sh_matrix(N_sph, 0, np.pi/2, 'real') y2 = 0.5 * spa.sph.sh_matrix(N_sph, np.pi/3, np.pi/2, 'real') y_org = (4 * np.pi) / (N_sph + 1)**2 * (y1 + y2) y_rot = spa.sph.rotate_sh(y_org, -np.pi/2, np.pi/8, np.pi/4) spa.plot.sh_coeffs_subplot([y_org, y_rot], titles=['before', 'after'], cbar=False)
- spaudiopy.sph.check_cond_sht(N_sph, azi, zen, sh_type, lim=None)[source]
Check if condition number for a least-squares SHT(N_sph) is high.
- spaudiopy.sph.n3d_to_sn3d(F_nm, sh_axis=0)[source]
Convert N3D (orthonormal) to SN3D (Schmidt semi-normalized) signals.
- Parameters:
F_nm (((N_sph+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
sh_axis (int, optional) – SH axis. The default is 0.
- Returns:
F_nm (((N_sph+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
- spaudiopy.sph.sn3d_to_n3d(F_nm, sh_axis=0)[source]
Convert SN3D (Schmidt semi-normalized) to N3D (orthonormal) signals.
- Parameters:
F_nm (((N_sph+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
sh_axis (int, optional) – SH axis. The default is 0.
- Returns:
F_nm (((N_sph+1)**2, S) numpy.ndarray) – Matrix of spherical harmonics coefficients of spherical function(S).
- spaudiopy.sph.sh_to_b(F_nm, W_weight=None)[source]
Convert first order SH input signals to B-format.
- Parameters:
F_nm ((4, S) numpy.ndarray) – First order spherical harmonics function coefficients.
W-weight (float) – Weight on W channel.
- Returns:
B_nm ((4, S) numpy.ndarray) – B-format signal(S).
- spaudiopy.sph.b_to_sh(B, W_weight=None)[source]
Convert B-format input to first order SH signals.
- Parameters:
B_nm ((4, S) numpy.ndarray) – B-format signal(S).
W-weight (float) – Weight on W channel.
- Returns:
F_nm ((4, S) numpy.ndarray) – First order spherical harmonics function coefficients.
- spaudiopy.sph.soundfield_to_b(sig, W_weight=None)[source]
Convert soundfield tetraeder mic input signals to B-format by SHT.
- Parameters:
sig ((4, S)) – Tetraeder mic signals(S).
W-weight (float) – Weight on W channel.
- Returns:
B_nm ((4, S) numpy.ndarray) – B-format signal(S).
- spaudiopy.sph.src_to_b(sig, src_azi, src_zen)[source]
Get B format signal channels for source in direction azi/zen.
- spaudiopy.sph.src_to_sh(sig, src_azi, src_zen, N_sph, sh_type='real')[source]
Source signal(s) plane wave encoded in spherical harmonics.
- Parameters:
sig ((num_src, S) numpy.ndarray) – Source signal(s).
src_azi (array_like)
src_zen (array_like)
N_sph (int)
sh_type (‘real’ (default) or ‘complex’, optional)
- Returns:
((N_sph+1)**2, S) numpy.ndarray – Source signal(s) in SHD.
Examples
src_sig = np.array([1, 0.5], ndmin=2).T * np.random.randn(2, 1000) N_sph = 3 src_azi = [np.pi/2, -np.pi/3] src_zen = [np.pi/4, np.pi/2] sig_nm = spa.sph.src_to_sh(src_sig, src_azi, src_zen, N_sph) spa.plot.sh_rms_map(sig_nm)
- spaudiopy.sph.bandlimited_dirac(N_sph, d, w_n=None)[source]
Order N spatially bandlimited Dirac pulse at central angle d.
- Parameters:
N_sph (int) – SH order.
d ((Q,) array_like) – Central angle in rad.
w_n ((N_sph,) array_like, optional. Default is None.) – Tapering window w_n.
- Returns:
dirac ((Q,) array_like) – Amplitude at central angle d.
Notes
Normalize with
\[\sum^N \frac{2N + 1}{4 \pi} = \frac{(N+1)^2}{4 \pi} .\]References
Rafaely, B. (2015). Fundamentals of Spherical Array Processing. Springer., eq. (1.60).
Examples
dirac_azi = np.deg2rad(0) dirac_zen = np.deg2rad(90) N_sph = 5 # cross section azi = np.linspace(0, 2 * np.pi, 720, endpoint=True) # Bandlimited Dirac pulse dirac_bandlim = 4 * np.pi / (N_sph + 1) ** 2 * \ spa.sph.bandlimited_dirac(N_sph, azi - dirac_azi) spa.plot.polar(azi, dirac_bandlim)
- spaudiopy.sph.max_rE_weights(N_sph)[source]
Return max-rE modal weight coefficients for spherical harmonics order N.
See also
spaudiopy.sph.unity_gain()
Unit amplitude compensation.
References
Zotter, F., & Frank, M. (2012). All-Round Ambisonic Panning and Decoding. Journal of Audio Engineering Society, eq. (10).
Examples
dirac_azi = np.deg2rad(45) dirac_zen = np.deg2rad(45) N = 5 # cross section azi = np.linspace(0, 2 * np.pi, 720, endpoint=True) # Bandlimited Dirac pulse, with max r_E tapering window w_n = spa.sph.max_rE_weights(N) w_n = spa.sph.unity_gain(w_n) dirac_tapered = spa.sph.bandlimited_dirac(N, azi - dirac_azi, w_n=w_n) spa.plot.polar(azi, dirac_tapered)
- spaudiopy.sph.r_E(p, g)[source]
Calculate r_E vector and magnitude from loudspeaker gains.
- Parameters:
p ((Q, 3) numpy.ndarray) – Q loudspeaker position vectors.
g ((S, Q) numpy.ndarray) – Q gain vectors per source S.
- Returns:
rE ((S, 3) numpy.ndarray) – rE vector.
rE_mag ((S,) array_like) – rE magnitude (radius).
References
Zotter, F., & Frank, M. (2012). All-Round Ambisonic Panning and Decoding. Journal of Audio Engineering Society, eq. (16).
- spaudiopy.sph.project_on_sphere(x, y, z)[source]
Little helper that projects x, y, z onto unit sphere.
- spaudiopy.sph.repeat_per_order(c)[source]
Repeat each coefficient in ‘c’ m times per spherical order n.
- Parameters:
c ((N,) array_like) – Coefficients up to SH order N.
- Returns:
c_reshaped (((N+1)**2,) array like) – Reshaped input coefficients.
- spaudiopy.sph.spherical_hn2(n, z, derivative=False)[source]
Spherical Hankel function of the second kind.
- Parameters:
n (int, array_like) – Order of the spherical Hankel function (n >= 0).
z (complex or float, array_like) – Argument of the spherical Hankel function.
derivative (bool, optional) – If True, the value of the derivative (rather than the function itself) is returned.
- Returns:
hn2 (array_like)
References
http://mathworld.wolfram.com/SphericalHankelFunctionoftheSecondKind.html
- spaudiopy.sph.mode_strength(n, kr, sphere_type='rigid')[source]
Mode strength b_n(kr) for an incident plane wave on sphere.
- Parameters:
n (int) – Degree.
kr (array_like) – kr vector, product of wavenumber k and radius r_0.
sphere_type (‘rigid’ or ‘open’)
- Returns:
b_n (array_like) – Mode strength b_n(kr).
References
Rafaely, B. (2015). Fundamentals of Spherical Array Processing. Springer. eq. (4.4) and (4.5).
- spaudiopy.sph.pressure_on_sphere(N_sph, kr, weights=None)[source]
Calculate the diffuse field pressure frequency response of a spherical scatterer, up to SH order N.
- Parameters:
N_sph (int) – SH order.
kr (array_like) – kr vector, product of wavenumber k and radius r_0.
weights ((N_sph+1,) array_like) – SH order weights.
- Returns:
p (array_like) – Pressure p(kr)|N.
References
Ben-Hur, Z., Brinkmann, F., Sheaffer, J., et.al. (2017). Spectral equalization in binaural signals represented by order-truncated spherical harmonics. The Journal of the Acoustical Society of America, eq. (11).
- spaudiopy.sph.binaural_coloration_compensation(N_sph, f, r_0=0.0875, w_taper=None)[source]
Spectral equalization gain G(kr)|N for diffuse field of order N_sph. This filter compensates the high frequency roll of that occurs for order truncated SH signals. It models the human head as a rigid sphere of radius r_0 (e.g. 0.0875m) and compensates the binaural signals.
- Parameters:
N_sph (int) – SH order.
f (array_like) – Time-frequency in Hz.
r_0 (radius) – Rigid sphere radius (approx. human head).
w_taper ((N+1,) array_like) – SH order weights for tapering. See e.g. ‘process.half_sided_Hann’.
- Returns:
gain (array_like) – Filter gain(kr).
See also
spaudiopy.process.gain_clipping()
Limit maximum gain.
References
Hold, C., Gamper, H., Pulkki, V., Raghuvanshi, N., & Tashev, I. J. (2019). Improving Binaural Ambisonics Decoding by Spherical Harmonics Domain Tapering and Coloration Compensation. In IEEE International Conference on Acoustics, Speech and Signal Processing.
Examples
fs = 48000 f = np.linspace(0, fs / 2, 1000) # target spherical harmonics order N (>= 3) N_sph = 5 # tapering window w_rE = spa.sph.max_rE_weights(N) compensation_tapered = spa.sph.binaural_coloration_compensation( N_sph, f, w_taper=w_rE) compensation_tapered_lim = spa.process.gain_clipping( compensation_tapered, spa.utils.from_db(12)) spa.plot.freq_resp(f, [compensation_tapered, compensation_tapered_lim], ylim=(-5, 25), labels=[r'$N=5, max_{rE}$', 'with soft lim'])
- spaudiopy.sph.unity_gain(w_n)[source]
Make modal weighting / tapering unit amplitude in steering direction.
- Parameters:
w_n ((N+1,) array_like) – Modal weighting factors.
- Returns:
w_n ((N+1,) array_like) – Modal weighting factors, adjusted for unit amplitude.
Examples
- spaudiopy.sph.hypercardioid_modal_weights(N_sph)[source]
Modal weights for beamformer resulting in a hyper-cardioid.
- Parameters:
N_sph (int) – SH order.
- Returns:
w_n ((N+1,) array_like) – Modal weighting factors.
Notes
Also called max-DI or normalized PWD.
Examples
N_sph = 5 w_n = spa.sph.hypercardioid_modal_weights(N_sph) w_nm = spa.sph.repeat_per_order(w_n) * spa.sph.sh_matrix(N_sph, np.pi/4, np.pi/4, 'real') spa.plot.sh_coeffs(w_nm)
- spaudiopy.sph.cardioid_modal_weights(N_sph)[source]
Modal weights for beamformer resulting in a cardioid.
- Parameters:
N_sph (int) – SH order.
- Returns:
w_n ((N+1,) array_like) – Modal weighting factors.
Notes
Also called in-phase weights, where \(w_{n, cardioid} = 4\pi*w_{n, inphase} / (N+1)\) .
Examples
N_sph = 5 w_n = spa.sph.cardioid_modal_weights(N_sph) w_nm = spa.sph.repeat_per_order(w_n) * spa.sph.sh_matrix(N_sph, np.pi/4, np.pi/4, 'real') spa.plot.sh_coeffs(w_nm)
- spaudiopy.sph.maxre_modal_weights(N_sph, UNITAMP=True)[source]
Modal weights for beamformer resulting with max-rE weighting.
- Parameters:
N_sph (int) – SH order.
UNITAMP (bool, optional (default:True))
- Returns:
w_n ((N+1,) array_like) – Modal weighting factors.
Notes
Can be compensated for unit amplitude.
Examples
N_sph = 5 w_n = spa.sph.maxre_modal_weights(N_sph) w_nm = spa.sph.repeat_per_order(w_n) * spa.sph.sh_matrix(N_sph, np.pi/4, np.pi/4, 'real') spa.plot.sh_coeffs(w_nm)
- spaudiopy.sph.butterworth_modal_weights(N_sph, k, n_c, UNITAMP=True)[source]
Modal weights for spatial butterworth filter / beamformer.
- Parameters:
N_sph (int) – SH order.
k (int (float)) – Filter order
n_c (int (float)) – Cut-on SH order.
UNITAMP (bool, optional (default:True))
- Returns:
w_n ((N+1,) array_like) – Modal weighting factors.
Notes
Can be compensated for unit amplitude.
References
Devaraju, B. (2015). Understanding filtering on the sphere.
Examples
N_sph = 5 w_n = spa.sph.butterworth_modal_weights(N_sph, 5, 3) w_nm = spa.sph.repeat_per_order(w_n) * spa.sph.sh_matrix(N_sph, np.pi/4, np.pi/4, 'real') spa.plot.sh_coeffs(w_nm)
- spaudiopy.sph.sph_filterbank_reconstruction_factor(w_nm, num_secs, mode=None)[source]
Reconstruction factor for restoring amplitude/energy preservation.
- Parameters:
w_nm (((N+1)**2,), array_like) – SH beam coefficients.
num_secs (int) – Number of SH beamformers.
mode (‘amplitude’ or ‘energy’)
- Raises:
ValueError – If mode is not specified.
- Returns:
beta (float) – Reconstruction factor.
References
Hold, C., Politis, A., Mc Cormack, L., & Pulkki, V. (2021). Spatial Filter Bank Design in the Spherical Harmonic Domain. EUSIPCO 2021.
- spaudiopy.sph.design_sph_filterbank(N_sph, sec_azi, sec_zen, c_n, sh_type, mode)[source]
Design spherical filter bank analysis and reconstruction matrix.
- Parameters:
N_sph (int) – SH order.
sec_azi ((J,) array_like) – Sector azimuth steering directions.
sec_zen ((J,) array_like) – Sector zenith/colatitude steering directions.
c_n ((N,) array_like) – SH Modal weights, describing (axisymmetric) pattern.
sh_type (‘real’ or ‘complex’)
mode (‘perfect’ or ‘energy’) – Design achieves perfect reconstruction or energy reconstruction.
- Raises:
ValueError – If mode not specified.
- Returns:
A ((J, (N+1)**2) numpy.ndarray) – Analysis matrix.
B (((N+1)**2, J) numpy.ndarray) – Resynthesis matrix.
References
Hold, C., Schlecht, S. J., Politis, A., & Pulkki, V. (2021). Spatial Filter Bank in the Spherical Harmonic Domain : Reconstruction and Application. WASPAA 2021.
Examples
N_sph = 3 sec_dirs = spa.utils.cart2sph(*spa.grids.load_t_design(2*N_sph).T) c_n = spa.sph.maxre_modal_weights(N_sph) [A, B] = spa.sph.design_sph_filterbank(N_sph, sec_dirs[0], sec_dirs[1], c_n, 'real', 'perfect') # diffuse input SH signal in_nm = np.random.randn((N_sph+1)**2, 1000) # Sector signals (Analysis) s_sec = A @ in_nm # Reconstruction to SH domain out_nm = B @ s_sec # Test perfect reconstruction print(spa.utils.test_diff(in_nm, out_nm))
- spaudiopy.sph.sh_mult(a_nm, b_nm, sh_type)[source]
Multiply SH vector a_nm and b_nm in (discrete) space.
- Parameters:
a_nm (((N1+1)**2,), array_like) – SH coefficients.
b_nm (((N2+1)**2,), array_like) – SH coefficients.
sh_type (‘real’ or ‘complex’)
- Returns:
c_nm (((N1+N2+1)**2,), array_like) – SH coefficients.
Examples
spa.plot.sh_coeffs(4*spa.sph.sh_mult([1, 0, 1, 0], [0, 1, 0, 0], 'real'))