LSST Applications g0265f82a02+c6dfa2ddaf,g1162b98a3f+ffe7eabc7e,g2079a07aa2+1b2e822518,g2bbee38e9b+c6dfa2ddaf,g337abbeb29+c6dfa2ddaf,g36da64cc00+ea84795170,g3ddfee87b4+955a963fd8,g50ff169b8f+2eb0e556e8,g52b1c1532d+90ebb246c7,g555ede804d+955a963fd8,g591dd9f2cf+bac198a2cb,g5ec818987f+420292cfeb,g858d7b2824+d6c9a0a3b8,g876c692160+aabc49a3c3,g8a8a8dda67+90ebb246c7,g8cdfe0ae6a+4fd9e222a8,g99cad8db69+e6cd765486,g9ddcbc5298+a1346535a5,ga1e77700b3+df8f93165b,ga8c6da7877+acd47f83f4,gae46bcf261+c6dfa2ddaf,gb0e22166c9+8634eb87fb,gb3f2274832+12c8382528,gba4ed39666+1ac82b564f,gbb8dafda3b+0574160a1f,gbeb006f7da+dea2fbb49f,gc28159a63d+c6dfa2ddaf,gc86a011abf+d6c9a0a3b8,gcf0d15dbbd+955a963fd8,gdaeeff99f8+1cafcb7cd4,gdc0c513512+d6c9a0a3b8,ge79ae78c31+c6dfa2ddaf,geb67518f79+ba1859f325,gee10cc3b42+90ebb246c7,gf1cff7945b+d6c9a0a3b8,w.2024.13
LSST Data Management Base Package
Loading...
Searching...
No Matches
Functions
lsst.pipe.tasks.associationUtils Namespace Reference

Functions

 toIndex (nside, ra, dec)
 
 toRaDec (nside, index)
 
 eq2xyz (ra, dec)
 
 eq2xyzVec (ra, dec)
 
 convert_spherical (ra, dec)
 
 convert_spherical_array (array)
 
 query_disc (nside, ra, dec, max_rad, min_rad=0)
 

Detailed Description

Utilities for interfacing with hpgeom. Originally implemented in
http://github.com/LSSTDESC/dia_pipe and then translated to hpgeom.

Function Documentation

◆ convert_spherical()

lsst.pipe.tasks.associationUtils.convert_spherical ( ra,
dec )
Convert from ra,dec to spherical coordinates.

Used in query_disc.

Parameters
----------
ra : `float`
    RA in radians.
dec : `float`
    Declination in radians

Definition at line 124 of file associationUtils.py.

124def convert_spherical(ra, dec):
125 """Convert from ra,dec to spherical coordinates.
126
127 Used in query_disc.
128
129 Parameters
130 ----------
131 ra : `float`
132 RA in radians.
133 dec : `float`
134 Declination in radians
135 """
136 return np.dstack([np.cos(dec*np.pi/180)*np.cos(ra*np.pi/180),
137 np.cos(dec*np.pi/180)*np.sin(ra*np.pi/180),
138 np.sin(dec*np.pi/180)])[0]
139
140

◆ convert_spherical_array()

lsst.pipe.tasks.associationUtils.convert_spherical_array ( array)
Convert from and a array ra,dec to spherical coordinates.

Used in query_disc

Parameters
----------
array : `numpy.ndarray`, (N, 2)
    (N, 2) Array of RA, DEC values.

Returns
-------
vecs : `numpy.ndarray`, (N, 3)
    Vectors on the unit sphere

Definition at line 141 of file associationUtils.py.

141def convert_spherical_array(array):
142 """Convert from and a array ra,dec to spherical coordinates.
143
144 Used in query_disc
145
146 Parameters
147 ----------
148 array : `numpy.ndarray`, (N, 2)
149 (N, 2) Array of RA, DEC values.
150
151 Returns
152 -------
153 vecs : `numpy.ndarray`, (N, 3)
154 Vectors on the unit sphere
155 """
156 ra = array[:, 0]
157 dec = array[:, 1]
158 return convert_spherical(ra, dec)
159
160

◆ eq2xyz()

lsst.pipe.tasks.associationUtils.eq2xyz ( ra,
dec )
Convert from equatorial ra,dec in degrees to x,y,z on unit sphere.

Parameters
----------
ra : `float`
    RA in degrees.
dec : `float`
    Declination in degrees

Returns
-------
xyz : `numpy.ndarray`, (3,)
    Float xyz positions on the unit sphere.

Definition at line 72 of file associationUtils.py.

72def eq2xyz(ra, dec):
73 """Convert from equatorial ra,dec in degrees to x,y,z on unit sphere.
74
75 Parameters
76 ----------
77 ra : `float`
78 RA in degrees.
79 dec : `float`
80 Declination in degrees
81
82 Returns
83 -------
84 xyz : `numpy.ndarray`, (3,)
85 Float xyz positions on the unit sphere.
86 """
87 phi = np.deg2rad(ra)
88 theta = np.pi/2 - np.deg2rad(dec)
89 sintheta = np.sin(theta)
90 x = sintheta*np.cos(phi)
91 y = sintheta*np.sin(phi)
92 z = np.cos(theta)
93 return np.array([x, y, z])
94
95

◆ eq2xyzVec()

lsst.pipe.tasks.associationUtils.eq2xyzVec ( ra,
dec )
Convert equatorial ra,dec in degrees to x,y,z on the unit sphere
parameters

Vectorized version of ``eq2xyz``

Parameters
----------
ra : array_like, (N,)
    Array of RA in degrees.
dec : array_like, (N,)
    Declination in degrees

Returns
-------
vec : `numpy.ndarray`, (N,3)
    Array of unitsphere 3-vectors.

Definition at line 96 of file associationUtils.py.

96def eq2xyzVec(ra, dec):
97 """Convert equatorial ra,dec in degrees to x,y,z on the unit sphere
98 parameters
99
100 Vectorized version of ``eq2xyz``
101
102 Parameters
103 ----------
104 ra : array_like, (N,)
105 Array of RA in degrees.
106 dec : array_like, (N,)
107 Declination in degrees
108
109 Returns
110 -------
111 vec : `numpy.ndarray`, (N,3)
112 Array of unitsphere 3-vectors.
113 """
114 ra = np.array(ra, dtype='f8', ndmin=1, copy=False)
115 dec = np.array(dec, dtype='f8', ndmin=1, copy=False)
116 if ra.size != dec.size:
117 raise ValueError("ra,dec not same size: %s,%s" % (ra.size, dec.size))
118
119 vec = eq2xyz(ra, dec)
120
121 return vec
122
123

◆ query_disc()

lsst.pipe.tasks.associationUtils.query_disc ( nside,
ra,
dec,
max_rad,
min_rad = 0 )
Get the list of healpix indices within max_rad, min_rad given in radians
around ra,dec given in degrees

Parameters
----------
nside : `int`
    Resolution of the healpixels to search/return.
ra : `float`
    RA in degrees.
dec : `float`
    Declination in degrees
max_rad : `float`
    Max distance in radians to search nearby healpixels.
min_rad : `float`, optional
    Minimum distance in radians to search healpixels. Default = 0.

Definition at line 161 of file associationUtils.py.

161def query_disc(nside, ra, dec, max_rad, min_rad=0):
162 """Get the list of healpix indices within max_rad, min_rad given in radians
163 around ra,dec given in degrees
164
165 Parameters
166 ----------
167 nside : `int`
168 Resolution of the healpixels to search/return.
169 ra : `float`
170 RA in degrees.
171 dec : `float`
172 Declination in degrees
173 max_rad : `float`
174 Max distance in radians to search nearby healpixels.
175 min_rad : `float`, optional
176 Minimum distance in radians to search healpixels. Default = 0.
177 """
178 ra = np.atleast_1d(ra)
179 dec = np.atleast_1d(dec)
180
181 max_rad_deg = np.rad2deg(max_rad)
182
183 pixels = np.unique(
184 [hpg.query_circle(nside, a, b, max_rad_deg, nest=False)
185 for (a, b) in zip(ra, dec)])
186
187 if min_rad > 0 and len(pixels) > 0:
188 vec0 = convert_spherical(ra, dec)
189 min_rad2 = min_rad**2
190 vecs = convert_spherical_array(toRaDec(nside, pixels))
191 dsq = np.sum((vecs - vec0)**2, axis=1)
192 match = dsq > min_rad2
193 pixels = pixels[match]
194
195 return pixels

◆ toIndex()

lsst.pipe.tasks.associationUtils.toIndex ( nside,
ra,
dec )
Return healpix index given ra, dec in degrees

Parameters
----------
nside : `int`
    Power of 2 nside healpix resolution.
ra : `float`
    RA in degrees.
dec : `float`
    Declination in degrees

Returns
-------
index : `int`
    Unique healpix pixel ID containing point RA, DEC at resolution nside.

Definition at line 33 of file associationUtils.py.

33def toIndex(nside, ra, dec):
34 """Return healpix index given ra, dec in degrees
35
36 Parameters
37 ----------
38 nside : `int`
39 Power of 2 nside healpix resolution.
40 ra : `float`
41 RA in degrees.
42 dec : `float`
43 Declination in degrees
44
45 Returns
46 -------
47 index : `int`
48 Unique healpix pixel ID containing point RA, DEC at resolution nside.
49 """
50 return hpg.angle_to_pixel(nside, ra, dec, nest=False)
51
52

◆ toRaDec()

lsst.pipe.tasks.associationUtils.toRaDec ( nside,
index )
Convert from healpix index to ra,dec in degrees

Parameters
----------
nside : `int`
    Resolution of healpixel "grid".
index : `int`
    Index of the healpix pixel we want to find the location of.

Returns
-------
pos : `numpy.ndarray`, (2,)
    RA and DEC of healpix pixel location in degrees.

Definition at line 53 of file associationUtils.py.

53def toRaDec(nside, index):
54 """Convert from healpix index to ra,dec in degrees
55
56 Parameters
57 ----------
58 nside : `int`
59 Resolution of healpixel "grid".
60 index : `int`
61 Index of the healpix pixel we want to find the location of.
62
63 Returns
64 -------
65 pos : `numpy.ndarray`, (2,)
66 RA and DEC of healpix pixel location in degrees.
67 """
68 ra, dec = hpg.pixel_to_angle(nside, index, nest=False)
69 return np.dstack((ra, dec))[0]
70
71