The pj_rhealpix Module¶
This Python 3.11 module implements the rHEALPix map projection.
Alexander Raichev (AR), 2013-01-26: Refactored code from release 0.3.
NOTE:
All lengths are measured in meters and all angles are measured in radians unless indicated otherwise. By ‘ellipsoid’ below, I mean an oblate ellipsoid of revolution.
- rhealpixdggs.pj_rhealpix.combine_triangles(x: float, y: float, north_square: int = 0, south_square: int = 0, inverse: bool = False) tuple[float, float] ¶
Rearrange point (x, y) in the HEALPix projection by combining the polar triangles into two polar squares. Put the north polar square in position north_square and the south polar square in position south_square. If inverse = True, uncombine the polar triangles.
INPUT:
x, y - Coordinates in the HEALPix projection of the unit sphere.
north_square, south_square - Integers between 0 and 3 indicating the positions of the north_square polar square and south_square polar square respectively. See rhealpix_sphere() docstring for a diagram.
inverse - (Optional; default = False) Boolean. If False, then compute forward function. If True, then compute inverse function.
EXAMPLES:
>>> u, v = -pi/4, pi/3 >>> x, y = combine_triangles(u, v) >>> print(tuple(x.tolist() for x in my_round((x, y), 15))) (-1.832595714594046, 1.570796326794896) >>> print(tuple(x.tolist() for x in my_round(combine_triangles(x, y, inverse=True), 15))) (-0.785398163397448, 1.047197551196598) >>> print(my_round((u, v), 15)) (-0.785398163397448, 1.047197551196598)
- rhealpixdggs.pj_rhealpix.in_rhealpix_image(x: float, y: float, north_square: int = 0, south_square: int = 0) bool ¶
Return True if and only if the point (x, y) lies in the image of the rHEALPix projection of the unit sphere.
EXAMPLES:
>>> eps = 0 # Test boundary points. >>> north_square, south_square = 0, 0 >>> rhp = [ ... (-pi - eps, pi/4 + eps), ... (-pi + north_square*pi/2 - eps, pi/4 + eps), ... (-pi + north_square*pi/2 - eps, 3*pi/4 + eps), ... (-pi + (north_square + 1)*pi/2 + eps, 3*pi/4 + eps), ... (-pi + (north_square + 1)*pi/2 + eps, pi/4 + eps), ... (pi + eps, pi/4 + eps), ... (pi + eps,-pi/4 - eps), ... (-pi + (south_square + 1)*pi/2 + eps,-pi/4 - eps), ... (-pi + (south_square + 1)*pi/2 + eps,-3*pi/4 - eps), ... (-pi + south_square*pi/2 - eps,-3*pi/4 - eps), ... (-pi + south_square*pi/2 -eps,-pi/4 - eps), ... (-pi - eps,-pi/4 - eps) ... ] >>> for p in rhp: ... if not in_rhealpix_image(*p): ... print('Fail') ... >>> print(in_rhealpix_image(0, 0)) True >>> print(in_rhealpix_image(0, pi/4 + 0.1)) False
- rhealpixdggs.pj_rhealpix.rhealpix(a: float = 1, e: float = 0, north_square: int = 0, south_square: int = 0, region: str = 'none') Callable[[float, float, bool, bool], tuple[float, float]] ¶
Return a function object that wraps the rHEALPix projection and its inverse of an ellipsoid with major radius a and eccentricity e.
EXAMPLES:
>>> f = rhealpix(a=2, e=0, north_square=1, south_square=2) >>> print(tuple(x.tolist() for x in my_round(f(0, pi/3, radians=True), 15))) (-0.574951359778215, 2.145747686573111) >>> p = (0, 60) >>> q = f(*p, radians=False) >>> print(tuple(x.tolist() for x in my_round(q, 15))) (-0.574951359778215, 2.145747686573111) >>> print(tuple(x.tolist() for x in my_round(f(*q, radians=False, inverse=True), 15))) (6e-15, 59.999999999999986) >>> print(my_round(p, 15)) (0, 60)
OUTPUT:
A function object of the form f(u, v, radians=False, inverse=False).
- rhealpixdggs.pj_rhealpix.rhealpix_ellipsoid(lam: float, phi: float, e: float = 0, north_square: int = 0, south_square: int = 0, region: str = 'none') tuple[float, float] ¶
Compute the signature functions of the rHEALPix map projection of an oblate ellipsoid with eccentricity e whose authalic sphere is the unit sphere. The north polar square is put in position north_square, and the south polar square is put in position south_square. Works when e = 0 (spherical case) too.
INPUT:
lam, phi - Geographic longitude-latitude coordinates in radian. Assume -pi <= lam < pi and -pi/2 <= phi <= pi/2.
e - Eccentricity of the ellipsoid.
north_square, south_square - (Optional; defaults = 0, 0) Integers between 0 and 3 indicating positions of north polar and south polar squares, respectively. See rhealpix_sphere() docstring for a diagram.
EXAMPLES:
>>> from numpy import arcsin >>> print(tuple(x if type(x) is int else x.tolist() for x in my_round(rhealpix_ellipsoid(0, arcsin(2.0/3)), 15))) (0, 0.785398163397448)
- rhealpixdggs.pj_rhealpix.rhealpix_ellipsoid_inverse(x: float, y: float, e: float = 0, north_square: int = 0, south_square: int = 0, region: str = 'none') tuple[float, float] ¶
Compute the inverse of rhealpix_ellipsoid().
EXAMPLES:
>>> p = (0, pi/4) >>> q = rhealpix_ellipsoid(*p) >>> print(tuple(x.tolist() for x in my_round(rhealpix_ellipsoid_inverse(*q), 15))) (0.0, 0.785398163397448) >>> print(my_round(p, 15)) (0, 0.785398163397448)
- rhealpixdggs.pj_rhealpix.rhealpix_sphere(lam: float, phi: float, north_square: int = 0, south_square: int = 0, region: str = 'none') tuple[float, float] ¶
Compute the signature functions of the rHEALPix map projection of the unit sphere. The north polar square is put in position north_square, and the south polar square is put in position south_square.
INPUT:
lam, phi -Geographic longitude-latitude coordinates in radians. Assume -pi <= lam < pi and -pi/2 <= phi <= pi/2.
north_square, south_square - (Optional; defaults = 0, 0) Integers between 0 and 3 indicating positions of north polar and south polar squares, respectively.
EXAMPLES:
>>> print(tuple(x.tolist() for x in my_round(rhealpix_sphere(0, pi/4), 15))) (-1.619978633413937, 2.307012183573304)
NOTE:
The polar squares are labeled 0, 1, 2, 3 from east to west like this:
east west *---*---*---*---* | 0 | 1 | 2 | 3 | *---*---*---*---* | | | | | *---*---*---*---* | 0 | 1 | 2 | 3 | *---*---*---*---*
- rhealpixdggs.pj_rhealpix.rhealpix_sphere_inverse(x: float, y: float, north_square: int = 0, south_square: int = 0, region: str = 'none') tuple[float, float] ¶
Compute the inverse of rhealpix_sphere().
EXAMPLES:
>>> p = (0, pi/4) >>> q = rhealpix_sphere(*p) >>> print(tuple(x.tolist() for x in my_round(rhealpix_sphere_inverse(*q), 15))) (0.0, 0.785398163397448) >>> print(my_round(p, 15)) (0, 0.785398163397448)
- rhealpixdggs.pj_rhealpix.rhealpix_vertices(north_square: int = 0, south_square: int = 0) list[tuple[float, float]] ¶
Return a list of the planar vertices of the rHEALPix projection of the unit sphere.
- rhealpixdggs.pj_rhealpix.triangle(x: float, y: float, north_square: int = 0, south_square: int = 0, inverse: bool = False) tuple[int, str] ¶
Return the number of the polar triangle and region that (x, y) lies in. If inverse = False, then assume (x,y) lies in the image of the HEALPix projection of the unit sphere. If inverse = True, then assume (x,y) lies in the image of the (north_square, south_square)-rHEALPix projection of the unit sphere.
INPUT:
x, y - Coordinates in the HEALPix or rHEALPix (if inverse = True) projection of the unit sphere.
north_square, south_square - Integers between 0 and 3 indicating the positions of the north_square pole square and south_square pole square respectively. See rhealpix_sphere() docstring for a diagram.
inverse - (Optional; default = False) Boolean. If False, then compute forward function. If True, then compute inverse function.
OUTPUT:
The pair (triangle_number, region). Here region equals ‘north_polar’ (polar), ‘south_polar’ (polar), or ‘equatorial’, indicating where (x, y) lies. If region = ‘equatorial’, then triangle_number = None. Suppose now that region != ‘equatorial’. If inverse = False, then triangle_number is the number (0, 1, 2, or 3) of the HEALPix polar triangle Z that (x, y) lies in. If inverse = True, then triangle_number is the number (0, 1, 2, or 3) of the HEALPix polar triangle that (x, y) will get moved into.
EXAMPLES:
>>> triangle(-pi/4, pi/4 + 0.1) (1, 'north_polar') >>> triangle(-3*pi/4 + 0.1, pi/2, inverse=True) (1, 'north_polar')
NOTES:
In the HEALPix projection, the polar triangles are labeled 0–3 from east to west like this:
* * * * * 0 * * 1 * * 2 * * 3 * *-------*-------*-------*-------* | | | | | | | | | | | | | | | *-------*-------*-------*-------* * 0 * * 1 * * 2 * * 3 * * * * *
In the rHEALPix projection these polar triangles get rearranged into a square with the triangles numbered north_square and south_square remaining fixed. For example, if north_square = 1 and south_square = 3, then the triangles get rearranged this way:
North polar square: *-------* | * 3 * | | 0 * 2 | | * 1 * | ----*-------*---- South polar square: ----*-------*---- | * 3 * | | 2 * 0 | | * 1 * | *-------*