LSSTApplications  20.0.0
LSSTDataManagementBasePackage
Namespaces | Classes | Typedefs | Functions | Variables
lsst::sphgeom Namespace Reference

Namespaces

 detail
 
 python
 
 version
 

Classes

class  Angle
 Angle represents an angle in radians. More...
 
class  AngleInterval
 AngleInterval represents closed intervals of arbitrary angles. More...
 
class  BigInteger
 BigInteger is an arbitrary precision signed integer class. More...
 
class  Box
 Box represents a rectangle in spherical coordinate space that contains its boundary. More...
 
class  Box3d
 Box3d represents a box in ℝ³. More...
 
class  Chunker
 Chunker subdivides the unit sphere into longitude-latitude boxes. More...
 
class  Circle
 Circle is a circular region on the unit sphere that contains its boundary. More...
 
class  ConvexPolygon
 ConvexPolygon is a closed convex polygon on the unit sphere. More...
 
class  Ellipse
 Ellipse is an elliptical region on the sphere. More...
 
class  HtmPixelization
 HtmPixelization provides HTM indexing of points and regions. More...
 
class  Interval
 Interval represents a closed interval of the real numbers by its upper and lower bounds. More...
 
class  Interval1d
 Interval1d represents closed intervals of ℝ. More...
 
class  LonLat
 LonLat represents a spherical coordinate (longitude/latitude angle) pair. More...
 
class  Matrix3d
 A 3x3 matrix with real entries stored in double precision. More...
 
class  Mq3cPixelization
 Mq3cPixelization provides modified Q3C indexing of points and regions. More...
 
class  NormalizedAngle
 NormalizedAngle is an angle that lies in the range [0, 2π), with one exception - a NormalizedAngle can be NaN. More...
 
class  NormalizedAngleInterval
 NormalizedAngleInterval represents closed intervals of normalized angles, i.e. More...
 
class  Pixelization
 A Pixelization (or partitioning) of the sphere is a mapping between points on the sphere and a set of pixels (a.k.a. More...
 
class  Q3cPixelization
 Q3cPixelization provides Q3C indexing of points and regions. More...
 
class  RangeSet
 A RangeSet is a set of unsigned 64 bit integers. More...
 
class  Region
 Region is a minimal interface for 2-dimensional regions on the unit sphere. More...
 
struct  SubChunks
 SubChunks represents a set of sub-chunks of a particular chunk. More...
 
class  UnitVector3d
 UnitVector3d is a unit vector in ℝ³ with components stored in double precision. More...
 
class  Vector3d
 Vector3d is a vector in ℝ³ with components stored in double precision. More...
 

Typedefs

using Relationship = std::bitset< 3 >
 Relationship describes how two sets are related. More...
 

Functions

Angle operator* (double a, Angle const &b)
 
std::ostreamoperator<< (std::ostream &, Angle const &)
 
double sin (Angle const &a)
 
double cos (Angle const &a)
 
double tan (Angle const &a)
 
Angle abs (Angle const &a)
 
std::ostreamoperator<< (std::ostream &, AngleInterval const &)
 
std::ostreamoperator<< (std::ostream &, Box const &)
 
std::ostreamoperator<< (std::ostream &, Box3d const &)
 
std::ostreamoperator<< (std::ostream &, Circle const &)
 
void encodeDouble (double item, std::vector< uint8_t > &buffer)
 encode appends an IEEE double in little-endian byte order to the end of buffer. More...
 
double decodeDouble (uint8_t const *buffer)
 decode extracts an IEEE double from the 8 byte little-endian byte sequence in buffer. More...
 
std::ostreamoperator<< (std::ostream &, ConvexPolygon const &)
 
uint64_t mortonIndex (uint32_t x, uint32_t y)
 mortonIndex interleaves the bits of x and y. More...
 
std::tuple< uint32_t, uint32_t > mortonIndexInverse (uint64_t z)
 mortonIndexInverse separates the even and odd bits of z. More...
 
uint64_t mortonToHilbert (uint64_t z, int m)
 mortonToHilbert converts the 2m-bit Morton index z to the corresponding Hilbert index. More...
 
uint64_t hilbertToMorton (uint64_t h, int m)
 hilbertToMorton converts the 2m-bit Hilbert index h to the corresponding Morton index. More...
 
uint64_t hilbertIndex (uint32_t x, uint32_t y, int m)
 hilbertIndex returns the index of (x, y) in a 2-D Hilbert curve. More...
 
std::tuple< uint32_t, uint32_t > hilbertIndexInverse (uint64_t h, int m)
 hilbertIndexInverse returns the point (x, y) with Hilbert index h, where x and y are m bit integers. More...
 
std::ostreamoperator<< (std::ostream &, Ellipse const &)
 
std::ostreamoperator<< (std::ostream &, Interval1d const &)
 
std::ostreamoperator<< (std::ostream &, LonLat const &)
 
std::ostreamoperator<< (std::ostream &, Matrix3d const &)
 
NormalizedAngle const & abs (NormalizedAngle const &a)
 
std::ostreamoperator<< (std::ostream &, NormalizedAngleInterval const &)
 
int orientationExact (Vector3d const &a, Vector3d const &b, Vector3d const &c)
 orientationExact computes and returns the orientations of 3 vectors a, b and c, which need not be normalized but are assumed to have finite components. More...
 
int orientation (UnitVector3d const &a, UnitVector3d const &b, UnitVector3d const &c)
 orientation computes and returns the orientations of 3 unit vectors a, b and c. More...
 
int orientationX (UnitVector3d const &b, UnitVector3d const &c)
 orientationX(b, c) is equivalent to orientation(UnitVector3d::X(), b, c). More...
 
int orientationY (UnitVector3d const &b, UnitVector3d const &c)
 orientationY(b, c) is equivalent to orientation(UnitVector3d::Y(), b, c). More...
 
int orientationZ (UnitVector3d const &b, UnitVector3d const &c)
 orientationZ(b, c) is equivalent to orientation(UnitVector3d::Z(), b, c). More...
 
void swap (RangeSet &a, RangeSet &b)
 
std::ostreamoperator<< (std::ostream &, RangeSet const &)
 
Relationship invert (Relationship r)
 Given the relationship between two sets A and B (i.e. More...
 
std::ostreamoperator<< (std::ostream &, UnitVector3d const &)
 
double getMinSquaredChordLength (Vector3d const &v, Vector3d const &a, Vector3d const &b, Vector3d const &n)
 Let p be the unit vector closest to v that lies on the plane with normal n in the direction of the cross product of a and b. More...
 
double getMaxSquaredChordLength (Vector3d const &v, Vector3d const &a, Vector3d const &b, Vector3d const &n)
 Let p be the unit vector furthest from v that lies on the plane with normal n in the direction of the cross product of a and b. More...
 
Angle getMinAngleToCircle (Angle x, Angle c)
 getMinAngleToCircle returns the minimum angular separation between a point at latitude x and the points on the circle of constant latitude c. More...
 
Angle getMaxAngleToCircle (Angle x, Angle c)
 getMaxAngleToCircle returns the maximum angular separation between a point at latitude x and the points on the circle of constant latitude c. More...
 
Vector3d getWeightedCentroid (UnitVector3d const &v0, UnitVector3d const &v1, UnitVector3d const &v2)
 getWeightedCentroid returns the center of mass of the given spherical triangle (assuming a uniform mass distribution over the triangle surface), weighted by the triangle area. More...
 
Vector3d operator* (double s, Vector3d const &v)
 
std::ostreamoperator<< (std::ostream &, Vector3d const &)
 

Variables

constexpr double PI = 3.1415926535897932384626433832795
 
constexpr double ONE_OVER_PI = 0.318309886183790671537767526745
 
constexpr double RAD_PER_DEG = 0.0174532925199432957692369076849
 
constexpr double DEG_PER_RAD = 57.2957795130823208767981548141
 
constexpr double MAX_ASIN_ERROR = 1.5e-8
 
constexpr double MAX_SQUARED_CHORD_LENGTH_ERROR = 2.5e-15
 
constexpr double EPSILON = 1.1102230246251565e-16
 

Typedef Documentation

◆ Relationship

Relationship describes how two sets are related.

Definition at line 35 of file Relationship.h.

Function Documentation

◆ abs() [1/2]

Angle lsst::sphgeom::abs ( Angle const &  a)
inline

Definition at line 106 of file Angle.h.

106 { return Angle(std::fabs(a.asRadians())); }

◆ abs() [2/2]

NormalizedAngle const& lsst::sphgeom::abs ( NormalizedAngle const &  a)
inline

Definition at line 150 of file NormalizedAngle.h.

150 { return a; }

◆ cos()

double lsst::sphgeom::cos ( Angle const &  a)
inline

Definition at line 103 of file Angle.h.

103 { return std::cos(a.asRadians()); }

◆ decodeDouble()

double lsst::sphgeom::decodeDouble ( uint8_t const *  buffer)
inline

decode extracts an IEEE double from the 8 byte little-endian byte sequence in buffer.

Definition at line 59 of file codec.h.

59  {
60 #if defined(__x86_64__)
61  // x86-64 is little endian and supports unaligned loads.
62  return *reinterpret_cast<double const *>(buffer);
63 #else
64  union { uint64_t u; double d };
65  u = static_cast<uint64_t>(buffer[0]) +
66  (static_cast<uint64_t>(buffer[1]) << 8) +
67  (static_cast<uint64_t>(buffer[2]) << 16) +
68  (static_cast<uint64_t>(buffer[3]) << 24) +
69  (static_cast<uint64_t>(buffer[4]) << 32) +
70  (static_cast<uint64_t>(buffer[5]) << 40) +
71  (static_cast<uint64_t>(buffer[6]) << 48) +
72  (static_cast<uint64_t>(buffer[7]) << 56);
73  return d;
74 #endif
75 }

◆ encodeDouble()

void lsst::sphgeom::encodeDouble ( double  item,
std::vector< uint8_t > &  buffer 
)
inline

encode appends an IEEE double in little-endian byte order to the end of buffer.

Definition at line 38 of file codec.h.

38  {
39 #if defined(__x86_64__)
40  // x86-64 is little endian.
41  auto ptr = reinterpret_cast<uint8_t const *>(&item);
42  buffer.insert(buffer.end(), ptr, ptr + 8);
43 #else
44  union { uint64_t u; double d };
45  d = item;
46  buffer.push_back(static_cast<uint8_t>(value));
47  buffer.push_back(static_cast<uint8_t>(value >> 8));
48  buffer.push_back(static_cast<uint8_t>(value >> 16));
49  buffer.push_back(static_cast<uint8_t>(value >> 24));
50  buffer.push_back(static_cast<uint8_t>(value >> 32));
51  buffer.push_back(static_cast<uint8_t>(value >> 40));
52  buffer.push_back(static_cast<uint8_t>(value >> 48));
53  buffer.push_back(static_cast<uint8_t>(value >> 56));
54 #endif
55 }

◆ getMaxAngleToCircle()

Angle lsst::sphgeom::getMaxAngleToCircle ( Angle  x,
Angle  c 
)
inline

getMaxAngleToCircle returns the maximum angular separation between a point at latitude x and the points on the circle of constant latitude c.

Definition at line 68 of file utils.h.

68  {
70  if (abs(x) <= abs(c)) {
71  return a + Angle(PI) - 2.0 * abs(c);
72  }
73  if (a < abs(x)) {
74  return Angle(PI) - 2.0 * abs(c) - a;
75  }
76  return Angle(PI) + 2.0 * abs(c) - a;
77 }

◆ getMaxSquaredChordLength()

double lsst::sphgeom::getMaxSquaredChordLength ( Vector3d const &  v,
Vector3d const &  a,
Vector3d const &  b,
Vector3d const &  n 
)

Let p be the unit vector furthest from v that lies on the plane with normal n in the direction of the cross product of a and b.

If p is in the interior of the great circle segment from a to b, then this helper function returns the squared chord length between p and v. Otherwise it returns 0 - the minimum squared chord length between any pair of points on the sphere.

Definition at line 58 of file utils.cc.

62 {
63  Vector3d vxn = v.cross(n);
64  if (vxn.dot(a) < 0.0 && vxn.dot(b) > 0.0) {
65  // v is in the lune defined by the half great circle passing through
66  // n and -a and the half great circle passing through n and -b, so p
67  // is in the interior of the great circle segment from a to b. The
68  // angle θ between p and v satisfies ‖v‖ ‖n‖ sin θ = |v·n|,
69  // and ‖v‖ ‖n‖ cos θ = -‖v × n‖. The desired squared chord length is
70  // 4 sin²(θ/2).
71  double s = std::fabs(v.dot(n));
72  double c = - vxn.getNorm();
73  double d = std::sin(0.5 * std::atan2(s, c));
74  return 4.0 * d * d;
75  }
76  return 0.0;
77 }

◆ getMinAngleToCircle()

Angle lsst::sphgeom::getMinAngleToCircle ( Angle  x,
Angle  c 
)
inline

getMinAngleToCircle returns the minimum angular separation between a point at latitude x and the points on the circle of constant latitude c.

Definition at line 62 of file utils.h.

62  {
63  return abs(x - c);
64 }

◆ getMinSquaredChordLength()

double lsst::sphgeom::getMinSquaredChordLength ( Vector3d const &  v,
Vector3d const &  a,
Vector3d const &  b,
Vector3d const &  n 
)

Let p be the unit vector closest to v that lies on the plane with normal n in the direction of the cross product of a and b.

If p is in the interior of the great circle segment from a to b, then this function returns the squared chord length between p and v. Otherwise it returns 4 - the maximum squared chord length between any pair of points on the unit sphere.

Definition at line 36 of file utils.cc.

40 {
41  Vector3d vxn = v.cross(n);
42  if (vxn.dot(a) > 0.0 && vxn.dot(b) < 0.0) {
43  // v is in the lune defined by the half great circle passing through
44  // n and a and the half great circle passing through n and b, so p
45  // is in the interior of the great circle segment from a to b. The
46  // angle θ between p and v satisfies ‖v‖ ‖n‖ sin θ = |v·n|,
47  // and ‖v‖ ‖n‖ cos θ = ‖v × n‖. The desired squared chord length is
48  // 4 sin²(θ/2).
49  double s = std::fabs(v.dot(n));
50  double c = vxn.getNorm();
51  double theta = (c == 0.0) ? 0.5 * PI : std::atan(s / c);
52  double d = std::sin(0.5 * theta);
53  return 4.0 * d * d;
54  }
55  return 4.0;
56 }

◆ getWeightedCentroid()

Vector3d lsst::sphgeom::getWeightedCentroid ( UnitVector3d const &  v0,
UnitVector3d const &  v1,
UnitVector3d const &  v2 
)

getWeightedCentroid returns the center of mass of the given spherical triangle (assuming a uniform mass distribution over the triangle surface), weighted by the triangle area.

Definition at line 79 of file utils.cc.

82 {
83  // For the details, see:
84  //
85  // The centroid and inertia tensor for a spherical triangle
86  // John E. Brock
87  // 1974, Naval Postgraduate School, Monterey Calif.
88  //
89  // https://openlibrary.org/books/OL25493734M/The_centroid_and_inertia_tensor_for_a_spherical_triangle
90 
91  Vector3d x01 = v0.robustCross(v1); // twice the cross product of v0 and v1
92  Vector3d x12 = v1.robustCross(v2);
93  Vector3d x20 = v2.robustCross(v0);
94  double s01 = 0.5 * x01.normalize(); // sine of the angle between v0 and v1
95  double s12 = 0.5 * x12.normalize();
96  double s20 = 0.5 * x20.normalize();
97  double c01 = v0.dot(v1); // cosine of the angle between v0 and v1
98  double c12 = v1.dot(v2);
99  double c20 = v2.dot(v0);
100  double a0 = (s12 == 0.0 && c12 == 0.0) ? 0.0 : std::atan2(s12, c12);
101  double a1 = (s20 == 0.0 && c20 == 0.0) ? 0.0 : std::atan2(s20, c20);
102  double a2 = (s01 == 0.0 && c01 == 0.0) ? 0.0 : std::atan2(s01, c01);
103  return 0.5 * (x01 * a2 + x12 * a0 + x20 * a1);
104 }

◆ hilbertIndex()

uint64_t lsst::sphgeom::hilbertIndex ( uint32_t  x,
uint32_t  y,
int  m 
)
inline

hilbertIndex returns the index of (x, y) in a 2-D Hilbert curve.

Only the m least significant bits of x and y are used. Computing the Hilbert index of a point has been measured to take 4 to 15 times as long as computing its Morton index on an Intel Core i7-3820QM CPU. With Xcode 7.3 and -O3, latency is ~19ns per call at a CPU frequency of 3.5 GHz.

Definition at line 349 of file curve.h.

349  {
350  return mortonToHilbert(mortonIndex(x, y), m);
351 }

◆ hilbertIndexInverse()

std::tuple<uint32_t, uint32_t> lsst::sphgeom::hilbertIndexInverse ( uint64_t  h,
int  m 
)
inline

hilbertIndexInverse returns the point (x, y) with Hilbert index h, where x and y are m bit integers.

Definition at line 361 of file curve.h.

361  {
363 }

◆ hilbertToMorton()

uint64_t lsst::sphgeom::hilbertToMorton ( uint64_t  h,
int  m 
)
inline

hilbertToMorton converts the 2m-bit Hilbert index h to the corresponding Morton index.

Definition at line 290 of file curve.h.

290  {
291  alignas(64) static uint8_t const HILBERT_INVERSE_LUT_3[256] = {
292  0x40, 0x02, 0x03, 0xc1, 0x04, 0x45, 0x47, 0x86,
293  0x0c, 0x4d, 0x4f, 0x8e, 0xcb, 0x89, 0x88, 0x4a,
294  0x20, 0x61, 0x63, 0xa2, 0x68, 0x2a, 0x2b, 0xe9,
295  0x6c, 0x2e, 0x2f, 0xed, 0xa7, 0xe6, 0xe4, 0x25,
296  0x30, 0x71, 0x73, 0xb2, 0x78, 0x3a, 0x3b, 0xf9,
297  0x7c, 0x3e, 0x3f, 0xfd, 0xb7, 0xf6, 0xf4, 0x35,
298  0xdf, 0x9d, 0x9c, 0x5e, 0x9b, 0xda, 0xd8, 0x19,
299  0x93, 0xd2, 0xd0, 0x11, 0x54, 0x16, 0x17, 0xd5,
300  0x00, 0x41, 0x43, 0x82, 0x48, 0x0a, 0x0b, 0xc9,
301  0x4c, 0x0e, 0x0f, 0xcd, 0x87, 0xc6, 0xc4, 0x05,
302  0x50, 0x12, 0x13, 0xd1, 0x14, 0x55, 0x57, 0x96,
303  0x1c, 0x5d, 0x5f, 0x9e, 0xdb, 0x99, 0x98, 0x5a,
304  0x70, 0x32, 0x33, 0xf1, 0x34, 0x75, 0x77, 0xb6,
305  0x3c, 0x7d, 0x7f, 0xbe, 0xfb, 0xb9, 0xb8, 0x7a,
306  0xaf, 0xee, 0xec, 0x2d, 0xe7, 0xa5, 0xa4, 0x66,
307  0xe3, 0xa1, 0xa0, 0x62, 0x28, 0x69, 0x6b, 0xaa,
308  0xff, 0xbd, 0xbc, 0x7e, 0xbb, 0xfa, 0xf8, 0x39,
309  0xb3, 0xf2, 0xf0, 0x31, 0x74, 0x36, 0x37, 0xf5,
310  0x9f, 0xde, 0xdc, 0x1d, 0xd7, 0x95, 0x94, 0x56,
311  0xd3, 0x91, 0x90, 0x52, 0x18, 0x59, 0x5b, 0x9a,
312  0x8f, 0xce, 0xcc, 0x0d, 0xc7, 0x85, 0x84, 0x46,
313  0xc3, 0x81, 0x80, 0x42, 0x08, 0x49, 0x4b, 0x8a,
314  0x60, 0x22, 0x23, 0xe1, 0x24, 0x65, 0x67, 0xa6,
315  0x2c, 0x6d, 0x6f, 0xae, 0xeb, 0xa9, 0xa8, 0x6a,
316  0xbf, 0xfe, 0xfc, 0x3d, 0xf7, 0xb5, 0xb4, 0x76,
317  0xf3, 0xb1, 0xb0, 0x72, 0x38, 0x79, 0x7b, 0xba,
318  0xef, 0xad, 0xac, 0x6e, 0xab, 0xea, 0xe8, 0x29,
319  0xa3, 0xe2, 0xe0, 0x21, 0x64, 0x26, 0x27, 0xe5,
320  0xcf, 0x8d, 0x8c, 0x4e, 0x8b, 0xca, 0xc8, 0x09,
321  0x83, 0xc2, 0xc0, 0x01, 0x44, 0x06, 0x07, 0xc5,
322  0x10, 0x51, 0x53, 0x92, 0x58, 0x1a, 0x1b, 0xd9,
323  0x5c, 0x1e, 0x1f, 0xdd, 0x97, 0xd6, 0xd4, 0x15
324  };
325  uint64_t z = 0;
326  uint64_t i = 0;
327  for (m = 2 * m; m >= 6;) {
328  m -= 6;
329  uint8_t j = HILBERT_INVERSE_LUT_3[i | ((h >> m) & 0x3f)];
330  z = (z << 6) | (j & 0x3f);
331  i = j & 0xc0;
332  }
333  if (m != 0) {
334  // m = 2 or 4
335  int r = 6 - m;
336  uint8_t j = HILBERT_INVERSE_LUT_3[i | ((h << r) & 0x3f)];
337  z = (z << m) | ((j & 0x3f) >> r);
338  }
339  return z;
340 }

◆ invert()

Relationship lsst::sphgeom::invert ( Relationship  r)
inline

Given the relationship between two sets A and B (i.e.

the output of A.relate(B)), invert returns the relationship between B and A (B.relate(A)).

Definition at line 55 of file Relationship.h.

55  {
56  // If A is disjoint from B, then B is disjoint from A. But if A contains B
57  // then B is within A, so the corresponding bits must be swapped.
58  return (r & DISJOINT) | ((r & CONTAINS) << 1) | ((r & WITHIN) >> 1);
59 }

◆ log2() [1/2]

uint8_t lsst::sphgeom::log2 ( uint32_t  x)
inline

log2 returns the index of the most significant 1 bit in x. If x is 0, the return value is 0.

A beautiful algorithm to find this index is presented in:

Using de Bruijn Sequences to Index a 1 in a Computer Word C. E. Leiserson, H. Prokop, and K. H. Randall. http://supertech.csail.mit.edu/papers/debruijn.pdf

Definition at line 125 of file curve.h.

125  {
126  // See https://graphics.stanford.edu/~seander/bithacks.html#IntegerLogDeBruijn
127  alignas(32) static uint8_t const PERFECT_HASH_TABLE[32] = {
128  0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30,
129  8, 12, 20, 28, 15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31
130  };
131  uint32_t const DE_BRUIJN_SEQUENCE = UINT32_C(0x07c4acdd);
132  x |= (x >> 1);
133  x |= (x >> 2);
134  x |= (x >> 4);
135  x |= (x >> 8);
136  x |= (x >> 16);
137  return PERFECT_HASH_TABLE[(DE_BRUIJN_SEQUENCE * x) >> 27];
138 }

◆ log2() [2/2]

uint8_t lsst::sphgeom::log2 ( uint64_t  x)
inline

log2 returns the index of the most significant 1 bit in x. If x is 0, the return value is 0.

A beautiful algorithm to find this index is presented in:

Using de Bruijn Sequences to Index a 1 in a Computer Word C. E. Leiserson, H. Prokop, and K. H. Randall. http://supertech.csail.mit.edu/papers/debruijn.pdf

Definition at line 98 of file curve.h.

98  {
99  alignas(64) static uint8_t const PERFECT_HASH_TABLE[64] = {
100  0, 1, 2, 7, 3, 13, 8, 19, 4, 25, 14, 28, 9, 34, 20, 40,
101  5, 17, 26, 38, 15, 46, 29, 48, 10, 31, 35, 54, 21, 50, 41, 57,
102  63, 6, 12, 18, 24, 27, 33, 39, 16, 37, 45, 47, 30, 53, 49, 56,
103  62, 11, 23, 32, 36, 44, 52, 55, 61, 22, 43, 51, 60, 42, 59, 58
104  };
105  uint64_t const DE_BRUIJN_SEQUENCE = UINT64_C(0x0218a392cd3d5dbf);
106  // First ensure that all bits below the MSB are set.
107  x |= (x >> 1);
108  x |= (x >> 2);
109  x |= (x >> 4);
110  x |= (x >> 8);
111  x |= (x >> 16);
112  x |= (x >> 32);
113  // Then, subtract them away.
114  x = x - (x >> 1);
115  // Multiplication by x is now a shift by the index i of the MSB.
116  //
117  // By definition, the value of the upper 6 bits of a 64-bit De Bruijn
118  // sequence left shifted by i is different for every value of i in
119  // [0, 64). It can therefore be used as an an index into a lookup table
120  // that recovers i. In other words, (DE_BRUIJN_SEQUENCE * x) >> 58 is a
121  // minimal perfect hash function for 64 bit powers of 2.
122  return PERFECT_HASH_TABLE[(DE_BRUIJN_SEQUENCE * x) >> 58];
123 }

◆ mortonIndex()

uint64_t lsst::sphgeom::mortonIndex ( uint32_t  x,
uint32_t  y 
)
inline

mortonIndex interleaves the bits of x and y.

The 32 even bits of the return value will be the bits of x, and the 32 odd bits those of y. This is the z-value of (x,y) defined by the Morton order function. See https://en.wikipedia.org/wiki/Z-order_curve for more information.

Definition at line 148 of file curve.h.

148  {
149  // This is just a 64-bit extension of:
150  // http://graphics.stanford.edu/~seander/bithacks.html#InterleaveBMN
151  uint64_t b = y;
152  uint64_t a = x;
153  b = (b | (b << 16)) & UINT64_C(0x0000ffff0000ffff);
154  a = (a | (a << 16)) & UINT64_C(0x0000ffff0000ffff);
155  b = (b | (b << 8)) & UINT64_C(0x00ff00ff00ff00ff);
156  a = (a | (a << 8)) & UINT64_C(0x00ff00ff00ff00ff);
157  b = (b | (b << 4)) & UINT64_C(0x0f0f0f0f0f0f0f0f);
158  a = (a | (a << 4)) & UINT64_C(0x0f0f0f0f0f0f0f0f);
159  b = (b | (b << 2)) & UINT64_C(0x3333333333333333);
160  a = (a | (a << 2)) & UINT64_C(0x3333333333333333);
161  b = (b | (b << 1)) & UINT64_C(0x5555555555555555);
162  a = (a | (a << 1)) & UINT64_C(0x5555555555555555);
163  return a | (b << 1);
164  }

◆ mortonIndexInverse()

std::tuple<uint32_t, uint32_t> lsst::sphgeom::mortonIndexInverse ( uint64_t  z)
inline

mortonIndexInverse separates the even and odd bits of z.

The 32 even bits of z are returned in the first element of the result tuple, and the 32 odd bits in the second. This is the inverse of mortonIndex().

Definition at line 195 of file curve.h.

195  {
196  uint64_t x = z & UINT64_C(0x5555555555555555);
197  uint64_t y = (z >> 1) & UINT64_C(0x5555555555555555);
198  x = (x | (x >> 1)) & UINT64_C(0x3333333333333333);
199  y = (y | (y >> 1)) & UINT64_C(0x3333333333333333);
200  x = (x | (x >> 2)) & UINT64_C(0x0f0f0f0f0f0f0f0f);
201  y = (y | (y >> 2)) & UINT64_C(0x0f0f0f0f0f0f0f0f);
202  x = (x | (x >> 4)) & UINT64_C(0x00ff00ff00ff00ff);
203  y = (y | (y >> 4)) & UINT64_C(0x00ff00ff00ff00ff);
204  x = (x | (x >> 8)) & UINT64_C(0x0000ffff0000ffff);
205  y = (y | (y >> 8)) & UINT64_C(0x0000ffff0000ffff);
206  return std::make_tuple(static_cast<uint32_t>(x | (x >> 16)),
207  static_cast<uint32_t>(y | (y >> 16)));
208  }

◆ mortonToHilbert()

uint64_t lsst::sphgeom::mortonToHilbert ( uint64_t  z,
int  m 
)
inline

mortonToHilbert converts the 2m-bit Morton index z to the corresponding Hilbert index.

Definition at line 236 of file curve.h.

236  {
237  alignas(64) static uint8_t const HILBERT_LUT_3[256] = {
238  0x40, 0xc3, 0x01, 0x02, 0x04, 0x45, 0x87, 0x46,
239  0x8e, 0x8d, 0x4f, 0xcc, 0x08, 0x49, 0x8b, 0x4a,
240  0xfa, 0x3b, 0xf9, 0xb8, 0x7c, 0xff, 0x3d, 0x3e,
241  0xf6, 0x37, 0xf5, 0xb4, 0xb2, 0xb1, 0x73, 0xf0,
242  0x10, 0x51, 0x93, 0x52, 0xde, 0x1f, 0xdd, 0x9c,
243  0x54, 0xd7, 0x15, 0x16, 0x58, 0xdb, 0x19, 0x1a,
244  0x20, 0x61, 0xa3, 0x62, 0xee, 0x2f, 0xed, 0xac,
245  0x64, 0xe7, 0x25, 0x26, 0x68, 0xeb, 0x29, 0x2a,
246  0x00, 0x41, 0x83, 0x42, 0xce, 0x0f, 0xcd, 0x8c,
247  0x44, 0xc7, 0x05, 0x06, 0x48, 0xcb, 0x09, 0x0a,
248  0x50, 0xd3, 0x11, 0x12, 0x14, 0x55, 0x97, 0x56,
249  0x9e, 0x9d, 0x5f, 0xdc, 0x18, 0x59, 0x9b, 0x5a,
250  0xba, 0xb9, 0x7b, 0xf8, 0xb6, 0xb5, 0x77, 0xf4,
251  0x3c, 0x7d, 0xbf, 0x7e, 0xf2, 0x33, 0xf1, 0xb0,
252  0x60, 0xe3, 0x21, 0x22, 0x24, 0x65, 0xa7, 0x66,
253  0xae, 0xad, 0x6f, 0xec, 0x28, 0x69, 0xab, 0x6a,
254  0xaa, 0xa9, 0x6b, 0xe8, 0xa6, 0xa5, 0x67, 0xe4,
255  0x2c, 0x6d, 0xaf, 0x6e, 0xe2, 0x23, 0xe1, 0xa0,
256  0x9a, 0x99, 0x5b, 0xd8, 0x96, 0x95, 0x57, 0xd4,
257  0x1c, 0x5d, 0x9f, 0x5e, 0xd2, 0x13, 0xd1, 0x90,
258  0x70, 0xf3, 0x31, 0x32, 0x34, 0x75, 0xb7, 0x76,
259  0xbe, 0xbd, 0x7f, 0xfc, 0x38, 0x79, 0xbb, 0x7a,
260  0xca, 0x0b, 0xc9, 0x88, 0x4c, 0xcf, 0x0d, 0x0e,
261  0xc6, 0x07, 0xc5, 0x84, 0x82, 0x81, 0x43, 0xc0,
262  0xea, 0x2b, 0xe9, 0xa8, 0x6c, 0xef, 0x2d, 0x2e,
263  0xe6, 0x27, 0xe5, 0xa4, 0xa2, 0xa1, 0x63, 0xe0,
264  0x30, 0x71, 0xb3, 0x72, 0xfe, 0x3f, 0xfd, 0xbc,
265  0x74, 0xf7, 0x35, 0x36, 0x78, 0xfb, 0x39, 0x3a,
266  0xda, 0x1b, 0xd9, 0x98, 0x5c, 0xdf, 0x1d, 0x1e,
267  0xd6, 0x17, 0xd5, 0x94, 0x92, 0x91, 0x53, 0xd0,
268  0x8a, 0x89, 0x4b, 0xc8, 0x86, 0x85, 0x47, 0xc4,
269  0x0c, 0x4d, 0x8f, 0x4e, 0xc2, 0x03, 0xc1, 0x80
270  };
271  uint64_t h = 0;
272  uint64_t i = 0;
273  for (m = 2 * m; m >= 6;) {
274  m -= 6;
275  uint8_t j = HILBERT_LUT_3[i | ((z >> m) & 0x3f)];
276  h = (h << 6) | (j & 0x3f);
277  i = j & 0xc0;
278  }
279  if (m != 0) {
280  // m = 2 or 4
281  int r = 6 - m;
282  uint8_t j = HILBERT_LUT_3[i | ((z << r) & 0x3f)];
283  h = (h << m) | ((j & 0x3f) >> r);
284  }
285  return h;
286 }

◆ operator*() [1/2]

Angle lsst::sphgeom::operator* ( double  a,
Angle const &  b 
)
inline

Definition at line 98 of file Angle.h.

98 { return b * a; }

◆ operator*() [2/2]

Vector3d lsst::sphgeom::operator* ( double  s,
Vector3d const &  v 
)
inline

Definition at line 165 of file Vector3d.h.

165 { return v * s; }

◆ operator<<() [1/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Angle const &  a 
)

Definition at line 34 of file Angle.cc.

34  {
35  char buf[32];
36  std::snprintf(buf, sizeof(buf), "%.17g", a.asRadians());
37  return os << buf;
38 }

◆ operator<<() [2/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
AngleInterval const &  i 
)

Definition at line 34 of file AngleInterval.cc.

34  {
35  return os << '[' << i.getA() << ", " << i.getB() << ']';
36 }

◆ operator<<() [3/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Box const &  b 
)

Definition at line 475 of file Box.cc.

475  {
476  return os << "{\"Box\": [" << b.getLon() << ", " << b.getLat() << "]}";
477 }

◆ operator<<() [4/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Box3d const &  b 
)

Definition at line 34 of file Box3d.cc.

34  {
35  return os << "{\"Box3d\": [" << b.x() << ", " << b.y() << ", " << b.z() << "]}";
36 }

◆ operator<<() [5/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Circle const &  c 
)

Definition at line 354 of file Circle.cc.

354  {
355  char tail[32];
356  std::snprintf(tail, sizeof(tail), ", %.17g]}", c.getSquaredChordLength());
357  return os << "{\"Circle\": [" << c.getCenter() << tail;
358 }

◆ operator<<() [6/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
ConvexPolygon const &  p 
)

Definition at line 382 of file ConvexPolygon.cc.

382  {
383  typedef std::vector<UnitVector3d>::const_iterator VertexIterator;
384  VertexIterator v = p.getVertices().begin();
385  VertexIterator const end = p.getVertices().end();
386  os << "{\"ConvexPolygon\": [" << *v;
387  for (++v; v != end; ++v) { os << ", " << *v; }
388  os << "]}";
389  return os;
390 }

◆ operator<<() [7/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Ellipse const &  e 
)

Definition at line 388 of file Ellipse.cc.

388  {
389  os << "{\"Ellipse\": [" << e.getTransformMatrix() << ", "
390  << e.getAlpha() << ", " << e.getBeta() << "]}";
391  return os;
392 }

◆ operator<<() [8/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Interval1d const &  i 
)

Definition at line 35 of file Interval1d.cc.

35  {
36  char buf[64];
37  std::snprintf(buf, sizeof(buf), "[%.17g, %.17g]", i.getA(), i.getB());
38  return os << buf;
39 }

◆ operator<<() [9/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
LonLat const &  p 
)

Definition at line 82 of file LonLat.cc.

82  {
83  return os << '[' << p.getLon() << ", " << p.getLat() << ']';
84 }

◆ operator<<() [10/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Matrix3d const &  m 
)

Definition at line 35 of file Matrix3d.cc.

35  {
36  return os << '[' << m.getRow(0) << ", " << m.getRow(1) << ", " << m.getRow(2) << ']';
37 }

◆ operator<<() [11/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
NormalizedAngleInterval const &  i 
)

Definition at line 284 of file NormalizedAngleInterval.cc.

286 {
287  return os << '[' << i.getA() << ", " << i.getB() << ']';
288 }

◆ operator<<() [12/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
RangeSet const &  s 
)

Definition at line 616 of file RangeSet.cc.

616  {
617  os << "{\"RangeSet\": [";
618  bool first = true;
619  for (auto const & t: s) {
620  if (!first) {
621  os << ", ";
622  }
623  first = false;
624  os << '[' << std::get<0>(t) << ", " << std::get<1>(t) << ']';
625  }
626  os << "]}";
627  return os;
628 }

◆ operator<<() [13/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
UnitVector3d const &  v 
)

Definition at line 73 of file UnitVector3d.cc.

73  {
74  return os << static_cast<Vector3d const &>(v);
75 }

◆ operator<<() [14/14]

std::ostream & lsst::sphgeom::operator<< ( std::ostream os,
Vector3d const &  v 
)

Definition at line 133 of file Vector3d.cc.

133  {
134  char buf[128];
135  std::snprintf(buf, sizeof(buf), "[%.17g, %.17g, %.17g]",
136  v.x(), v.y(), v.z());
137  return os << buf;
138 }

◆ orientation()

int lsst::sphgeom::orientation ( UnitVector3d const &  a,
UnitVector3d const &  b,
UnitVector3d const &  c 
)

orientation computes and returns the orientations of 3 unit vectors a, b and c.

The return value is +1 if the vectors are in counter-clockwise orientation, 0 if they are coplanar, colinear or identical, and -1 if they are in clockwise orientation.

This is equivalent to computing the sign of the scalar triple product a · (b x c), which is the sign of the determinant of the 3x3 matrix with a, b and c as columns/rows.

The implementation proceeds by first computing a double precision approximation, and then falling back to arbitrary precision arithmetic when necessary. Consequently, the result is exact.

Definition at line 135 of file orientation.cc.

138 {
139  // This constant is a little more than 5ε, where ε = 2^-53. When multiplied
140  // by the permanent of |M|, it gives an error bound on the determinant of
141  // M. Here, M is a 3x3 matrix and |M| denotes the matrix obtained by
142  // taking the absolute value of each of its components. The derivation of
143  // this proceeds in the same manner as the derivation of the error bounds
144  // in section 4.3 of:
145  //
146  // Adaptive Precision Floating-Point Arithmetic
147  // and Fast Robust Geometric Predicates,
148  // Jonathan Richard Shewchuk,
149  // Discrete & Computational Geometry 18(3):305–363, October 1997.
150  //
151  // available online at http://www.cs.berkeley.edu/~jrs/papers/robustr.pdf
152  static double const relativeError = 5.6e-16;
153  // Because all 3 unit vectors are normalized, the maximum absolute value of
154  // any vector component, cross product component or dot product term in
155  // the calculation is very close to 1. The permanent of |M| must therefore
156  // be below 3 + c, where c is some small multiple of ε. This constant, a
157  // little larger than 3 * 5ε, is an upper bound on the absolute error in
158  // the determinant calculation.
159  static double const maxAbsoluteError = 1.7e-15;
160  // This constant accounts for floating point underflow (assuming hardware
161  // without gradual underflow, just to be conservative) in the computation
162  // of det(M). It is a little more than 14 * 2^-1022.
163  static double const minAbsoluteError = 4.0e-307;
164 
165  double bycz = b.y() * c.z();
166  double bzcy = b.z() * c.y();
167  double bzcx = b.z() * c.x();
168  double bxcz = b.x() * c.z();
169  double bxcy = b.x() * c.y();
170  double bycx = b.y() * c.x();
171  double determinant = a.x() * (bycz - bzcy) +
172  a.y() * (bzcx - bxcz) +
173  a.z() * (bxcy - bycx);
174  if (determinant > maxAbsoluteError) {
175  return 1;
176  } else if (determinant < -maxAbsoluteError) {
177  return -1;
178  }
179  // Expend some more effort on what is hopefully a tighter error bound
180  // before falling back on arbitrary precision arithmetic.
181  double permanent = std::fabs(a.x()) * (std::fabs(bycz) + std::fabs(bzcy)) +
182  std::fabs(a.y()) * (std::fabs(bzcx) + std::fabs(bxcz)) +
183  std::fabs(a.z()) * (std::fabs(bxcy) + std::fabs(bycx));
184  double maxError = relativeError * permanent + minAbsoluteError;
185  if (determinant > maxError) {
186  return 1;
187  } else if (determinant < -maxError) {
188  return -1;
189  }
190  // Avoid the slow path when any two inputs are identical or antipodal.
191  if (a == b || b == c || a == c || a == -b || b == -c || a == -c) {
192  return 0;
193  }
194  return orientationExact(a, b, c);
195 }

◆ orientationExact()

int lsst::sphgeom::orientationExact ( Vector3d const &  a,
Vector3d const &  b,
Vector3d const &  c 
)

orientationExact computes and returns the orientations of 3 vectors a, b and c, which need not be normalized but are assumed to have finite components.

The return value is +1 if the vectors a, b, and c are in counter-clockwise orientation, 0 if they are coplanar, colinear, or identical, and -1 if they are in clockwise orientation. The implementation uses arbitrary precision arithmetic to avoid floating point rounding error, underflow and overflow.

Definition at line 76 of file orientation.cc.

79 {
80  // Product mantissa storage buffers.
81  uint32_t mantissaBuffers[6][6];
82  // Product mantissas.
83  BigInteger mantissas[6] = {
84  BigInteger(mantissaBuffers[0], 6),
85  BigInteger(mantissaBuffers[1], 6),
86  BigInteger(mantissaBuffers[2], 6),
87  BigInteger(mantissaBuffers[3], 6),
88  BigInteger(mantissaBuffers[4], 6),
89  BigInteger(mantissaBuffers[5], 6)
90  };
91  BigFloat products[6] = {
92  BigFloat(&mantissas[0]),
93  BigFloat(&mantissas[1]),
94  BigFloat(&mantissas[2]),
95  BigFloat(&mantissas[3]),
96  BigFloat(&mantissas[4]),
97  BigFloat(&mantissas[5])
98  };
99  // An accumulator and its storage.
100  uint32_t accumulatorBuffer[512];
101  BigInteger accumulator(accumulatorBuffer,
102  sizeof(accumulatorBuffer) / sizeof(uint32_t));
103  // Compute the products in the determinant. Performing all multiplication
104  // up front means that each product mantissa occupies at most 3*53 bits.
105  computeProduct(products[0], a.x(), b.y(), c.z());
106  computeProduct(products[1], a.x(), b.z(), c.y());
107  computeProduct(products[2], a.y(), b.z(), c.x());
108  computeProduct(products[3], a.y(), b.x(), c.z());
109  computeProduct(products[4], a.z(), b.x(), c.y());
110  computeProduct(products[5], a.z(), b.y(), c.x());
111  mantissas[1].negate();
112  mantissas[3].negate();
113  mantissas[5].negate();
114  // Sort the array of products in descending exponent order.
115  std::sort(products, products + 6, [](BigFloat const & a, BigFloat const & b) {
116  return a.exponent > b.exponent;
117  });
118  // First, initialize the accumulator to the product with the highest
119  // exponent, then add the remaining products. Prior to each addition, we
120  // must shift the accumulated value so that its radix point lines up with
121  // the the radix point of the product to add.
122  //
123  // More precisely, at each step we have an accumulated value A·2ʲ and a
124  // product P·2ᵏ, and we update the accumulator to equal (A·2ʲ⁻ᵏ + P)·2ᵏ.
125  // Because the products were sorted beforehand, j ≥ k and 2ʲ⁻ᵏ is an
126  // integer.
127  accumulator = *products[0].mantissa;
128  for (int i = 1; i < 6; ++i) {
129  accumulator.multiplyPow2(products[i - 1].exponent - products[i].exponent);
130  accumulator.add(*products[i].mantissa);
131  }
132  return accumulator.getSign();
133 }

◆ orientationX()

int lsst::sphgeom::orientationX ( UnitVector3d const &  b,
UnitVector3d const &  c 
)

orientationX(b, c) is equivalent to orientation(UnitVector3d::X(), b, c).

Definition at line 227 of file orientation.cc.

227  {
228  int o = _orientationXYZ(b.y() * c.z(), b.z() * c.y());
229  return (o != 0) ? o : orientationExact(UnitVector3d::X(), b, c);
230 }

◆ orientationY()

int lsst::sphgeom::orientationY ( UnitVector3d const &  b,
UnitVector3d const &  c 
)

orientationY(b, c) is equivalent to orientation(UnitVector3d::Y(), b, c).

Definition at line 232 of file orientation.cc.

232  {
233  int o = _orientationXYZ(b.z() * c.x(), b.x() * c.z());
234  return (o != 0) ? o : orientationExact(UnitVector3d::Y(), b, c);
235 }

◆ orientationZ()

int lsst::sphgeom::orientationZ ( UnitVector3d const &  b,
UnitVector3d const &  c 
)

orientationZ(b, c) is equivalent to orientation(UnitVector3d::Z(), b, c).

Definition at line 237 of file orientation.cc.

237  {
238  int o = _orientationXYZ(b.x() * c.y(), b.y() * c.x());
239  return (o != 0) ? o : orientationExact(UnitVector3d::Z(), b, c);
240 }

◆ sin()

double lsst::sphgeom::sin ( Angle const &  a)
inline

Definition at line 102 of file Angle.h.

102 { return std::sin(a.asRadians()); }

◆ swap()

void lsst::sphgeom::swap ( RangeSet a,
RangeSet b 
)
inline

Definition at line 608 of file RangeSet.h.

608  {
609  a.swap(b);
610 }

◆ tan()

double lsst::sphgeom::tan ( Angle const &  a)
inline

Definition at line 104 of file Angle.h.

104 { return std::tan(a.asRadians()); }

Variable Documentation

◆ DEG_PER_RAD

constexpr double lsst::sphgeom::DEG_PER_RAD = 57.2957795130823208767981548141
constexpr

Definition at line 39 of file constants.h.

◆ EPSILON

constexpr double lsst::sphgeom::EPSILON = 1.1102230246251565e-16
constexpr

Definition at line 54 of file constants.h.

◆ MAX_ASIN_ERROR

constexpr double lsst::sphgeom::MAX_ASIN_ERROR = 1.5e-8
constexpr

Definition at line 45 of file constants.h.

◆ MAX_SQUARED_CHORD_LENGTH_ERROR

constexpr double lsst::sphgeom::MAX_SQUARED_CHORD_LENGTH_ERROR = 2.5e-15
constexpr

Definition at line 50 of file constants.h.

◆ ONE_OVER_PI

constexpr double lsst::sphgeom::ONE_OVER_PI = 0.318309886183790671537767526745
constexpr

Definition at line 37 of file constants.h.

◆ PI

constexpr double lsst::sphgeom::PI = 3.1415926535897932384626433832795
constexpr

Definition at line 36 of file constants.h.

◆ RAD_PER_DEG

constexpr double lsst::sphgeom::RAD_PER_DEG = 0.0174532925199432957692369076849
constexpr

Definition at line 38 of file constants.h.

y
int y
Definition: SpanSet.cc:49
std::tan
T tan(T... args)
std::make_tuple
T make_tuple(T... args)
lsst::geom::PI
constexpr double PI
The ratio of a circle's circumference to diameter.
Definition: Angle.h:39
std::fabs
T fabs(T... args)
std::atan2
T atan2(T... args)
std::cos
T cos(T... args)
std::vector
STL class.
lsst::sphgeom::abs
Angle abs(Angle const &a)
Definition: Angle.h:106
lsst::afw::table._match.first
first
Definition: _match.py:76
lsst::afw::image::X
@ X
Definition: ImageUtils.h:36
std::sort
T sort(T... args)
end
int end
Definition: BoundedField.cc:105
std::vector::push_back
T push_back(T... args)
lsst::afw::table::Angle
lsst::geom::Angle Angle
Definition: misc.h:33
std::snprintf
T snprintf(T... args)
z
double z
Definition: Match.cc:44
x
double x
Definition: ChebyshevBoundedField.cc:277
std::atan
T atan(T... args)
lsst::sphgeom::getMinAngleToCircle
Angle getMinAngleToCircle(Angle x, Angle c)
getMinAngleToCircle returns the minimum angular separation between a point at latitude x and the poin...
Definition: utils.h:62
lsst::afw::image::Y
@ Y
Definition: ImageUtils.h:36
ptr
uint64_t * ptr
Definition: RangeSet.cc:88
b
table::Key< int > b
Definition: TransmissionCurve.cc:467
lsst::sphgeom::mortonIndexInverse
std::tuple< uint32_t, uint32_t > mortonIndexInverse(uint64_t z)
mortonIndexInverse separates the even and odd bits of z.
Definition: curve.h:195
lsst::sphgeom::mortonToHilbert
uint64_t mortonToHilbert(uint64_t z, int m)
mortonToHilbert converts the 2m-bit Morton index z to the corresponding Hilbert index.
Definition: curve.h:236
std::sin
T sin(T... args)
lsst::sphgeom::hilbertToMorton
uint64_t hilbertToMorton(uint64_t h, int m)
hilbertToMorton converts the 2m-bit Hilbert index h to the corresponding Morton index.
Definition: curve.h:290
os
std::ostream * os
Definition: Schema.cc:746
exponent
int exponent
Definition: orientation.cc:41
a
table::Key< int > a
Definition: TransmissionCurve.cc:466
std::vector::begin
T begin(T... args)
lsst::sphgeom::orientationExact
int orientationExact(Vector3d const &a, Vector3d const &b, Vector3d const &c)
orientationExact computes and returns the orientations of 3 vectors a, b and c, which need not be nor...
Definition: orientation.cc:76
std::vector::insert
T insert(T... args)
std::vector::end
T end(T... args)
lsst::sphgeom::mortonIndex
uint64_t mortonIndex(uint32_t x, uint32_t y)
mortonIndex interleaves the bits of x and y.
Definition: curve.h:148
mantissa
BigInteger * mantissa
Definition: orientation.cc:40
m
int m
Definition: SpanSet.cc:49