pyfar.samplings¶
Collection of sampling schemes and related functionality. For information on
the used coordinate systems refer to the
concepts
.
Classes:
|
Voronoi diagrams on the surface of a sphere. |
Functions:
|
Calculate sampling weights for numeric integration. |
|
Create a cuboid sampling with equidistant spacings in x, y, and z. |
|
Generate a sampling based on the center points of the twelve dodecahedron faces. |
|
Generate sampling of the sphere with equally spaced angles. |
|
Sampling based on partitioning into faces with equal area. |
|
Generate an equiangular sampling of the sphere. |
|
Return a Hyperinterpolation sampling grid. |
|
Return Fliege-Maier spherical sampling grid. |
|
Generate sampling of the sphere based on the Gaussian quadrature. |
|
Spherical sampling grid according to the great circle distance criterion. |
|
Generate a sampling from the center points of the twenty icosahedron faces. |
|
Return Lebedev spherical sampling grid. |
|
Return spherical t-design sampling grid. |
- class pyfar.samplings.SphericalVoronoi(sampling, round_decimals=12, center=0.0)[source]¶
Bases:
scipy.spatial._spherical_voronoi.SphericalVoronoi
Voronoi diagrams on the surface of a sphere. Note that
calculate_sph_voronoi_weights
can be used directly, if only the sampling weights are needed.Methods:
__init__
(sampling[, round_decimals, center])Calculate a Voronoi diagram on the sphere for the given samplings points.
copy
()Return a copy of the Voronoi object.
- __init__(sampling, round_decimals=12, center=0.0)[source]¶
Calculate a Voronoi diagram on the sphere for the given samplings points.
- Parameters
sampling (Coordinates) – Spherical sampling.
round_decimals (int) – Number of decimals to be rounded for checking for equal radius. The default is
12
.center (double) – Center point of the voronoi diagram. The default is
0
.
- Returns
voronoi – Spherical voronoi diagram as implemented in
scipy.spatial
.- Return type
See also
- pyfar.samplings.calculate_sph_voronoi_weights(sampling, normalize=True, center=[0, 0, 0], round_decimals=12)[source]¶
Calculate sampling weights for numeric integration.
Uses the class method
calculate_areas
fromSphericalVoronoi
to calculate the weights. It requires a spherical sampling grid with a single radius and usesscipy.spatial.SphericalVoronoi
in the background.- Parameters
sampling (Coordinates) – Sampling points on a sphere, i.e., all points must have the same radius.
normalize (boolean, optional) – Normalize the samplings weights to
sum(weights)=1
. Otherwise the weights sum to. The default is
True
.center (list) – Center of the spherical sampling grid. The default is
[0, 0, 0]
.round_decimals (int, optional) – Round to round_decimals digits to check for equal radius. The default is
12
.
- Returns
weigths – Sampling weights of size samplings.csize.
- Return type
ndarray, double
- pyfar.samplings.cart_equidistant_cube(n_points)[source]¶
Create a cuboid sampling with equidistant spacings in x, y, and z.
The cube will have dimensions 1 x 1 x 1.
- Parameters
n_points (int, tuple) – Number of points in the sampling. If a single value is given, the number of sampling positions will be the same in every axis. If a tuple is given, the number of points will be set as
(n_x, n_y, n_z)
.- Returns
sampling – Sampling positions. Does not contain sampling weights.
- Return type
- pyfar.samplings.sph_dodecahedron(radius=1.0)[source]¶
Generate a sampling based on the center points of the twelve dodecahedron faces.
- Parameters
radius (number, optional) – Radius of the sampling grid. The default is
1
.- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
- pyfar.samplings.sph_equal_angle(delta_angles, radius=1.0)[source]¶
Generate sampling of the sphere with equally spaced angles.
This sampling contain points at the North and South Pole. See
sph_equiangular
,sph_gaussian
, andsph_great_circle
for samplings that do not contain points at the poles.- Parameters
delta_angles (tuple, number) – Tuple that gives the angular spacing in azimuth and colatitude in degrees. If a number is provided, the same spacing is applied in both dimensions.
radius (number, optional) – Radius of the sampling grid. The default is
1
.
- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
- pyfar.samplings.sph_equal_area(n_points, radius=1.0)[source]¶
Sampling based on partitioning into faces with equal area.
For detailed information, see 1.
- Parameters
n_points (int) – Number of points corresponding to the number of partitions of the sphere.
radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
References
- 1
P. Leopardi, “A partition of the unit sphere into regions of equal area and small diameter,” Electronic Transactions on Numerical Analysis, vol. 25, no. 12, pp. 309–327, 2006.
- pyfar.samplings.sph_equiangular(n_points=None, sh_order=None, radius=1.0)[source]¶
Generate an equiangular sampling of the sphere.
For detailed information, see 2, Chapter 3.2. This sampling does not contain points at the North and South Pole and is typically used for spherical harmonics processing. See
sph_equal_angle
andsph_great_circle
for samplings containing points at the poles.- Parameters
n_points (int, tuple of two ints) – Number of sampling points in azimuth and elevation. Either n_points or sh_order must be provided. The default is
None
.sh_order (int) – Maximum applicable spherical harmonic order. If this is provided, ‘n_points’ is set to
2 * sh_order + 1
. Either n_points or sh_order must be provided. The default isNone
.radius (number, optional) – Radius of the sampling grid. The default is
1
.
- Returns
sampling – Sampling positions including sampling weights.
- Return type
References
- 2
B. Rafaely, Fundamentals of spherical array processing, 1st ed. Berlin, Heidelberg, Germany: Springer, 2015.
- pyfar.samplings.sph_extremal(n_points=None, sh_order=None, radius=1.0)[source]¶
Return a Hyperinterpolation sampling grid.
After Sloan and Womersley 3. The samplings are available for 1 <= sh_order <= 200 (
n_points = (sh_order + 1)^2
).- Parameters
n_points (int) – Number of sampling points in the grid. Related to the spherical harmonic order by
n_points = (sh_order + 1)**2
. Either n_points or sh_order must be provided. The default isNone
.sh_order (int) – Maximum applicable spherical harmonic order. Related to the number of points by
sh_order = np.sqrt(n_points) - 1
. Either n_points or sh_order must be provided. The default isNone
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions including sampling weights.
- Return type
Notes
This implementation uses precalculated sets of points from 4. The data up to
sh_order = 99
are loaded the first time this function is called. The remaining data is loaded upon request.References
- pyfar.samplings.sph_fliege(n_points=None, sh_order=None, radius=1.0)[source]¶
Return Fliege-Maier spherical sampling grid.
For detailed information, see 5. Call
sph_fliege
for a list of possible values for n_points and sh_order.- Parameters
n_points (int, optional) – Number of sampling points in the grid. Related to the spherical harmonic order by
n_points = (sh_order + 1)**2
. Either n_points or sh_order must be provided. The default isNone
.sh_order (int, optional) – Maximum applicable spherical harmonic order. Related to the number of points by
sh_order = np.sqrt(n_points) - 1
. Either n_points or sh_order must be provided. The default isNone
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions including sampling weights.
- Return type
Notes
This implementation uses pre-calculated points from the SOFiA toolbox 6. Possible combinations of n_points and sh_order are:
n_points
sh_order
4
1
9
2
16
3
25
4
36
5
49
6
64
7
81
8
100
9
121
10
144
11
169
12
196
13
225
14
256
15
289
16
324
17
361
18
400
19
441
20
484
21
529
22
576
23
625
24
676
25
729
26
784
27
841
28
900
29
References
- pyfar.samplings.sph_gaussian(n_points=None, sh_order=None, radius=1.0)[source]¶
Generate sampling of the sphere based on the Gaussian quadrature.
For detailed information, see 7 (Section 3.3). This sampling does not contain points at the North and South Pole and is typically used for spherical harmonics processing. See
sph_equal_angle
andsph_great_circle
for samplings containing points at the poles.- Parameters
n_points (int, tuple of two ints) – Number of sampling points in azimuth and elevation. Either n_points or sh_order must be provided. The default is
None
.sh_order (int) – Maximum applicable spherical harmonic order. If this is provided, n_points is set to
(2 * (sh_order + 1), sh_order + 1)
. Either n_points or sh_order must be provided. The default isNone
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions including sampling weights.
- Return type
References
- 7
B. Rafaely, Fundamentals of spherical array processing, 1st ed. Berlin, Heidelberg, Germany: Springer, 2015.
- pyfar.samplings.sph_great_circle(elevation=array([- 90.0, - 80.0, - 70.0, - 60.0, - 50.0, - 40.0, - 30.0, - 20.0, - 10.0, 0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0]), gcd=10, radius=1, azimuth_res=1, match=360)[source]¶
Spherical sampling grid according to the great circle distance criterion.
Sampling grid where neighboring points of the same elevation have approx. the same great circle distance across elevations 8.
- Parameters
elevation (array like, optional) – Contains the elevation from wich the sampling grid is generated, with
(
: North Pole,
: South Pole). The default is
np.linspace(-90, 90, 19)
.gcd (number, optional) – Desired great circle distance (GCD). Note that the actual GCD of the sampling grid is equal or smaller then the desired GCD and that the GCD may vary across elevations. The default is
10
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.azimuth_res (number, optional) – Minimum resolution of the azimuth angle in degree. The default is
1
.match (number, optional) – Forces azimuth entries to appear with a period of match degrees. E.g., if
match=90
, the grid includes the azimuth angles 0, 90, 180, and 270 degrees. The default is360
.
- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
References
- 8
B. P. Bovbjerg, F. Christensen, P. Minnaar, and X. Chen, “Measuring the head-related transfer functions of an artificial head with a high directional resolution,” 109th AES Convention, Los Angeles, USA, Sep. 2000.
- pyfar.samplings.sph_icosahedron(radius=1.0)[source]¶
Generate a sampling from the center points of the twenty icosahedron faces.
- Parameters
radius (number, optional) – Radius of the sampling grid. The default is
1
.- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
- pyfar.samplings.sph_lebedev(n_points=None, sh_order=None, radius=1.0)[source]¶
Return Lebedev spherical sampling grid.
For detailed information, see 9. For a list of available values for n_points and sh_order call
sph_lebedev
.- Parameters
n_points (int, optional) – Number of sampling points in the grid. Related to the spherical harmonic order by
n_points = (sh_order + 1)**2
. Either n_points or sh_order must be provided. The default isNone
.sh_order (int, optional) – Maximum applicable spherical harmonic order. Related to the number of points by
sh_order = np.sqrt(n_points) - 1
. Either n_points or sh_order must be provided. The default isNone
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions including sampling weights.
- Return type
Notes
This is a Python port of the Matlab Code written by Rob Parrish 10.
References
- 9
V.I. Lebedev, and D.N. Laikov “A quadrature formula for the sphere of the 131st algebraic order of accuracy” Doklady Mathematics, Vol. 59, No. 3, 1999, pp. 477-481.
- 10
https://de.mathworks.com/matlabcentral/fileexchange/27097- getlebedevsphere
- pyfar.samplings.sph_t_design(degree=None, sh_order=None, criterion='const_energy', radius=1.0)[source]¶
Return spherical t-design sampling grid.
For detailed information, see 11. For a spherical harmonic order
, a t-Design of degree
for constant energy or
additionally ensuring a constant angular spread of energy is required 12. For a given degree t
points will be generated, except for t = 3, 5, 7, 9, 11, 13, and 15. T-designs allow for an inverse spherical harmonic transform matrix calculated as
with
being the hermitian transpose of the spherical harmonics matrix.
- Parameters
degree (int) – T-design degree between
1
and180
. Either degree or sh_order must be provided. The default isNone
.sh_order (int) – Maximum applicable spherical harmonic order. Related to the degree by
degree = 2 * sh_order
(const_energy
) anddegree = 2 * sh_order + 1
(const_angular_spread
). Either degree or sh_order must be provided. The default isNone
.criterion (
const_energy
,const_angular_spread
) – Design criterion ensuring only a constant energy or additionally constant angular spread of energy. The default isconst_energy
.radius (number, optional) – Radius of the sampling grid in meters. The default is
1
.
- Returns
sampling – Sampling positions. Sampling weights can be obtained from
calculate_sph_voronoi_weights
.- Return type
Notes
This function downloads a pre-calculated set of points from 13 . The data up to
degree = 99
are loaded the first time this function is called. The remaining data is loaded upon request.References
- 11
C. An, X. Chen, I. H. Sloan, and R. S. Womersley, “Well Conditioned Spherical Designs for Integration and Interpolation on the Two-Sphere,” SIAM Journal on Numerical Analysis, vol. 48, no. 6, pp. 2135–2157, Jan. 2010.
- 12
F. Zotter, M. Frank, and A. Sontacchi, “The Virtual T-Design Ambisonics-Rig Using VBAP,” in Proceedings on the Congress on Sound and Vibration, 2010.
- 13