LSSTApplications  16.0-1-gce273f5+17,16.0-1-gf77f410+12,16.0-10-g5a5abec+8,16.0-10-gc1446dd+12,16.0-12-g1dc09ba+6,16.0-12-g569485f,16.0-12-ga22ed6e+1,16.0-13-g4c33ca5+12,16.0-13-gb122224+3,16.0-13-gd9b1b71+12,16.0-14-g22e2ff2,16.0-14-g71e547a+8,16.0-17-g0bdc215+4,16.0-17-g6a7bfb3b+12,16.0-2-g0febb12+14,16.0-2-g839ba83+50,16.0-2-g9d5294e+39,16.0-20-ga7ad2685,16.0-3-g404ea43+9,16.0-3-gbc759ec+10,16.0-3-gcfd6c53+37,16.0-4-g03cf288+28,16.0-4-g13a27c5+14,16.0-4-g5f3a788+13,16.0-4-g8a0f11a+34,16.0-4-ga3eb747+3,16.0-45-g4805a823c,16.0-5-g1991253+12,16.0-5-g1e9226d+1,16.0-5-g865efd9+12,16.0-5-gb3f8a4b+44,16.0-5-gd0f1235+6,16.0-6-gf0acd13+31,16.0-6-gf9cb114+13,16.0-7-g6043bfc,16.0-7-ga8e1655+8,16.0-8-g23bbf3f+3,16.0-8-g4dec96c+25,16.0-8-gfd407c0+2,master-g965b868a3d+1,master-gdc6be1965f+1,w.2018.39
LSSTDataManagementBasePackage
Angle.h
Go to the documentation of this file.
1 /*
2  * Developed for the LSST Data Management System.
3  * This product includes software developed by the LSST Project
4  * (https://www.lsst.org).
5  * See the COPYRIGHT file at the top-level directory of this distribution
6  * for details of code ownership.
7  *
8  * This program is free software: you can redistribute it and/or modify
9  * it under the terms of the GNU General Public License as published by
10  * the Free Software Foundation, either version 3 of the License, or
11  * (at your option) any later version.
12  *
13  * This program is distributed in the hope that it will be useful,
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  * GNU General Public License for more details.
17  *
18  * You should have received a copy of the GNU General Public License
19  * along with this program. If not, see <https://www.gnu.org/licenses/>.
20  */
21 #ifndef LSST_GEOM_ANGLE_H
22 #define LSST_GEOM_ANGLE_H
23 
24 #include <cmath>
25 #include <iostream>
26 #include <type_traits>
27 
28 #include "boost/math/constants/constants.hpp"
29 
30 namespace lsst {
31 namespace geom {
32 
33 /*
34  * None of C99, C++98, and C++11 define M_PI, so we'll do it ourselves
35  */
36 #pragma clang diagnostic push
37 #pragma clang diagnostic ignored "-Wunused-variable"
38 double constexpr PI = boost::math::constants::pi<double>();
40 double constexpr TWOPI = boost::math::constants::pi<double>() * 2.0;
41 double constexpr HALFPI = boost::math::constants::pi<double>() * 0.5;
42 double constexpr ONE_OVER_PI = 1.0 / boost::math::constants::pi<double>();
43 // sqrt is not a constexpr on OS X
44 double const SQRTPI = sqrt(boost::math::constants::pi<double>());
45 double const INVSQRTPI = 1.0 / sqrt(boost::math::constants::pi<double>());
46 double constexpr ROOT2 = boost::math::constants::root_two<double>(); // sqrt(2)
47 #pragma clang diagnostic pop
48 
49 // These shouldn't be necessary if the Angle class is used, but sometimes you just need
50 // them. Better to define them once here than have *180/PI throughout the code...
51 inline constexpr double degToRad(double x) noexcept { return x * PI / 180.; }
52 inline constexpr double radToDeg(double x) noexcept { return x * 180. / PI; }
53 inline constexpr double radToArcsec(double x) noexcept { return x * 3600. * 180. / PI; }
54 inline constexpr double radToMas(double x) noexcept { return x * 1000. * 3600. * 180. / PI; }
55 inline constexpr double arcsecToRad(double x) noexcept { return (x / 3600.) * PI / 180.; }
56 inline constexpr double masToRad(double x) noexcept { return (x / (1000. * 3600.)) * PI / 180.; }
57 
58 class Angle;
70 class AngleUnit final {
71  friend class Angle;
72  template <typename T>
73  friend constexpr Angle operator*(T lhs, AngleUnit rhs) noexcept;
74 
75 public:
83  // Current implementation does not throw exceptions, but somebody may want
84  // to add input validation later
85  explicit constexpr AngleUnit(double val) : _val(val) {}
86 
95  constexpr bool operator==(AngleUnit const& rhs) const noexcept;
96 
97 private:
98  double _val;
99 };
100 
101 inline constexpr bool AngleUnit::operator==(AngleUnit const& rhs) const noexcept {
102  return (_val == rhs._val);
103 }
104 
105 AngleUnit constexpr radians = AngleUnit(1.0);
106 AngleUnit constexpr degrees = AngleUnit(PI / 180.0);
107 AngleUnit constexpr hours = AngleUnit(PI * 15.0 / 180.0);
108 AngleUnit constexpr arcminutes = AngleUnit(PI / 60 / 180.0);
109 AngleUnit constexpr arcseconds = AngleUnit(PI / 180.0 / 3600.0);
110 // Note: if we use PI / 180.0 / 3.6e6 then 60*60*180*1000*lsst.geom.milliarcseconds
111 // does not test exactly equal to 180*lsst.geom.degrees
113  AngleUnit(PI / (180.0 * 3.6e6));
114 
124 class Angle final {
125  friend class AngleUnit;
126 
127 public:
133  explicit constexpr Angle(double val, AngleUnit units = radians) noexcept : _val(val * units._val) {}
134 
136  constexpr Angle() noexcept : _val(0) {}
137 
139  constexpr Angle(Angle const& other) noexcept = default;
140 
142  constexpr Angle(Angle&& other) noexcept = default;
143 
145  Angle& operator=(Angle const& other) noexcept = default;
146 
148  Angle& operator=(Angle&& other) noexcept = default;
149 
150  ~Angle() = default;
151 
153  constexpr operator double() const noexcept { return _val; }
154 
160  constexpr double asAngularUnits(AngleUnit const& units) const noexcept { return _val / units._val; }
161 
163  constexpr double asRadians() const noexcept { return asAngularUnits(radians); }
164 
166  constexpr double asDegrees() const noexcept { return asAngularUnits(degrees); }
167 
169  constexpr double asHours() const noexcept { return asAngularUnits(hours); }
170 
172  constexpr double asArcminutes() const noexcept { return asAngularUnits(arcminutes); }
173 
175  constexpr double asArcseconds() const noexcept { return asAngularUnits(arcseconds); }
176 
178  constexpr double asMilliarcseconds() const noexcept { return asAngularUnits(milliarcseconds); }
179 
192  Angle wrap() const noexcept;
193 
206  Angle wrapCtr() const noexcept;
207 
220  Angle wrapNear(Angle const& refAng) const noexcept;
221 
230  Angle separation(Angle const& other) const noexcept;
231 
232 #define ANGLE_OPUP_TYPE(OP, TYPE) \
233  Angle& operator OP(TYPE const& d) noexcept { \
234  _val OP d; \
235  return *this; \
236  }
237 
239  ANGLE_OPUP_TYPE(*=, double)
243 
244  ANGLE_OPUP_TYPE(+=, double)
246  ANGLE_OPUP_TYPE(+=, int)
248 
249  ANGLE_OPUP_TYPE(-=, double)
251  ANGLE_OPUP_TYPE(-=, int)
253 
254 #undef ANGLE_OPUP_TYPE
255 
256 #define ANGLE_COMP(OP) \
257  constexpr bool operator OP(const Angle& rhs) const noexcept { return _val OP rhs._val; }
258 
260  ANGLE_COMP(==)
264 
265  ANGLE_COMP(<=)
271 
272 #undef ANGLE_COMP
273 
274 private:
275  double _val;
276 };
277 
278 /*
279  * Operators for Angles.
280  */
281 #define ANGLE_OP(OP) \
282  inline constexpr Angle operator OP(Angle a, Angle d) noexcept { \
283  return Angle(static_cast<double>(a) OP static_cast<double>(d)); \
284  }
285 
286 // We need both int and double versions to avoid ambiguous overloading due to
287 // implicit conversion of Angle to double
288 #define ANGLE_OP_TYPE(OP, TYPE) \
289  inline constexpr Angle operator OP(Angle a, TYPE d) noexcept { \
290  return Angle(static_cast<double>(a) OP d); \
291  } \
292  \
293  inline constexpr Angle operator OP(TYPE d, Angle a) noexcept { \
294  return Angle(d OP static_cast<double>(a)); \
295  }
296 
318 
323 ANGLE_OP_TYPE(*, double)
326 
327 #undef ANGLE_OP
328 #undef ANGLE_OP_TYPE
329 
335 inline constexpr Angle operator-(Angle angle) { return Angle(-static_cast<double>(angle)); }
336 
337 // Apparently @relatesalso doesn't work with grouping
343 inline constexpr Angle operator/(Angle a, int d) noexcept { return Angle(static_cast<double>(a) / d); }
344 
350 inline constexpr Angle operator/(Angle a, double d) noexcept { return Angle(static_cast<double>(a) / d); }
351 
352 // Division is different. Don't allow division by an Angle
353 template <typename T>
354 constexpr double operator/(T const lhs, Angle rhs) noexcept = delete;
355 
371 
373 template <typename T>
374 inline constexpr bool isAngle(T) noexcept {
376 };
377 
388 template <typename T>
389 inline constexpr Angle operator*(T lhs, AngleUnit rhs) noexcept {
390  static_assert(std::is_arithmetic<T>::value,
391  "Only numeric types may be multiplied by an AngleUnit to create an Angle!");
392  return Angle(lhs * rhs._val);
393 }
394 
395 // Inline method definitions, placed last in order to benefit from Angle's full API
396 
397 inline Angle Angle::wrap() const noexcept {
398  double wrapped = std::fmod(_val, TWOPI);
399  // wrapped is in the range (-TWOPI, TWOPI)
400  if (wrapped < 0.0) wrapped += TWOPI;
401  // if wrapped is small enough, adding 2 pi gives 2 pi
402  if (wrapped >= TWOPI) wrapped = 0.0;
403  return wrapped * radians;
404 }
405 
406 inline Angle Angle::wrapCtr() const noexcept {
407  double wrapped = std::fmod(_val, TWOPI);
408  // wrapped is in the range [-TWOPI, TWOPI]
409  if (wrapped < -PI) {
410  wrapped += TWOPI;
411  if (wrapped >= PI) {
412  // handle roundoff error, however unlikely
413  wrapped = -PI;
414  }
415  } else if (wrapped >= PI) {
416  wrapped -= TWOPI;
417  if (wrapped < -PI) {
418  // handle roundoff error, however unlikely
419  wrapped = -PI;
420  }
421  }
422  return wrapped * radians;
423 }
424 
425 inline Angle Angle::wrapNear(Angle const& refAng) const noexcept {
426  // compute (this - refAng).wrapCtr() + refAng
427  // which is correct except for roundoff error at the edges
428  double const refAngRad = refAng.asRadians();
429  double wrapped = (*this - refAng).wrapCtr().asRadians() + refAngRad;
430 
431  // roundoff can cause slightly out-of-range values; fix those as bast we can;
432  // note that both conditionals are wanted, since either or both may execute
433  // (though if/else could be used if the lower limit was squishy for radians)
434  if (wrapped - refAngRad >= PI) {
435  wrapped -= TWOPI;
436  }
437  if (wrapped - refAngRad < -PI) {
438  wrapped += TWOPI;
439  }
440  return wrapped * radians;
441 }
442 
443 inline Angle Angle::separation(Angle const& other) const noexcept { return (*this - other).wrapCtr(); }
444 
445 } // namespace geom
446 } // namespace lsst
447 
448 #endif // if !defined(LSST_GEOM_ANGLE_H)
constexpr double asAngularUnits(AngleUnit const &units) const noexcept
Return an Angle&#39;s value in the specified units.
Definition: Angle.h:160
AngleUnit constexpr arcminutes
constant with units of arcminutes
Definition: Angle.h:108
constexpr double asRadians() const noexcept
Return an Angle&#39;s value in radians.
Definition: Angle.h:163
friend class Angle
Definition: Angle.h:71
constexpr double asDegrees() const noexcept
Return an Angle&#39;s value in degrees.
Definition: Angle.h:166
#define ANGLE_OPUP_TYPE(OP, TYPE)
Definition: Angle.h:232
double constexpr TWOPI
Definition: Angle.h:40
AngleUnit constexpr hours
constant with units of hours
Definition: Angle.h:107
constexpr double radToDeg(double x) noexcept
Definition: Angle.h:52
T fmod(T... args)
table::Key< int > a
Angle separation(Angle const &other) const noexcept
The signed difference between two Angles.
Definition: Angle.h:443
#define ANGLE_OP_TYPE(OP, TYPE)
Definition: Angle.h:288
ImageT val
Definition: CR.cc:146
def wrap(ctrl)
Definition: wrap.py:248
AngleUnit constexpr radians
constant with units of radians
Definition: Angle.h:105
A class representing an angle.
Definition: Angle.h:124
Angle wrapCtr() const noexcept
Wrap this angle to the range [-π, π).
Definition: Angle.h:406
constexpr double radToArcsec(double x) noexcept
Definition: Angle.h:53
double constexpr ONE_OVER_PI
Definition: Angle.h:42
constexpr Angle operator-(Angle a, Angle d) noexcept
Difference of two angles.
Definition: Angle.h:308
AngleUnit constexpr degrees
constant with units of degrees
Definition: Angle.h:106
A base class for image defects.
Definition: cameraGeom.dox:3
constexpr double asHours() const noexcept
Return an Angle&#39;s value in hours.
Definition: Angle.h:169
AngleUnit constexpr arcseconds
constant with units of arcseconds
Definition: Angle.h:109
constexpr Angle() noexcept
Construct the zero angle.
Definition: Angle.h:136
constexpr double degToRad(double x) noexcept
Definition: Angle.h:51
std::ostream & operator<<(std::ostream &os, lsst::geom::AffineTransform const &transform)
solver_t * s
constexpr double asArcseconds() const noexcept
Return an Angle&#39;s value in arcseconds.
Definition: Angle.h:175
double x
constexpr Angle(double val, AngleUnit units=radians) noexcept
Construct an Angle with the specified value (interpreted in the given units).
Definition: Angle.h:133
#define ANGLE_COMP(OP)
Definition: Angle.h:256
Angle wrapNear(Angle const &refAng) const noexcept
Wrap this angle to a value x such that -π ≤ x - refAng ≤ π, approximately.
Definition: Angle.h:425
constexpr double asMilliarcseconds() const noexcept
Return an Angle&#39;s value in milliarcseconds.
Definition: Angle.h:178
double const INVSQRTPI
Definition: Angle.h:45
constexpr Angle operator/(Angle a, int d) noexcept
Ratio of an angle and a scalar.
Definition: Angle.h:343
constexpr bool operator==(AngleUnit const &rhs) const noexcept
Test if two units are the same.
Definition: Angle.h:101
double constexpr HALFPI
Definition: Angle.h:41
friend constexpr Angle operator*(T lhs, AngleUnit rhs) noexcept
Use AngleUnit to convert a POD (e.g. int, double) to an Angle; e.g. 180*degrees.
Definition: Angle.h:389
constexpr double radToMas(double x) noexcept
Definition: Angle.h:54
double const SQRTPI
Definition: Angle.h:44
ItemVariant const * other
Definition: Schema.cc:55
constexpr double asArcminutes() const noexcept
Return an Angle&#39;s value in arcminutes.
Definition: Angle.h:172
constexpr AngleUnit(double val)
Define a new angle unit.
Definition: Angle.h:85
Angle wrap() const noexcept
Wrap this angle to the range [0, 2π).
Definition: Angle.h:397
AngleUnit constexpr milliarcseconds
constant with units of milliarcseconds
Definition: Angle.h:112
constexpr bool isAngle(T) noexcept
Allow a user to check if they have an angle.
Definition: Angle.h:374
constexpr double masToRad(double x) noexcept
Definition: Angle.h:56
constexpr double arcsecToRad(double x) noexcept
Definition: Angle.h:55
STL class.
#define ANGLE_OP(OP)
Definition: Angle.h:281
double constexpr ROOT2
Definition: Angle.h:46
double constexpr PI
The ratio of a circle&#39;s circumference to diameter.
Definition: Angle.h:39
A class used to convert scalar POD types such as double to Angle.
Definition: Angle.h:70