CppEphem
CECoordinates.h
Go to the documentation of this file.
1 /***************************************************************************
2  * CECoordinates.h: CppEphem *
3  * ----------------------------------------------------------------------- *
4  * Copyright © 2016-2019 JCardenzana *
5  * ----------------------------------------------------------------------- *
6  * *
7  * This program is free software: you can redistribute it and/or modify *
8  * it under the terms of the GNU General Public License as published by *
9  * the Free Software Foundation, either version 3 of the License, or *
10  * (at your option) any later version. *
11  * *
12  * This program is distributed in the hope that it will be useful, *
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15  * GNU General Public License for more details. *
16  * *
17  * You should have received a copy of the GNU General Public License *
18  * along with this program. If not, see <http://www.gnu.org/licenses/>. *
19  * *
20  ***************************************************************************/
21 
22 #ifndef CECoordinates_h
23 #define CECoordinates_h
24 
25 #include <string>
26 #include <vector>
27 
28 // CppEphem HEADERS
29 #include "CEAngle.h"
30 #include "CEDate.h"
31 #include "CENamespace.h"
32 #include "CEException.h"
33 
34 // SOFA HEADER
35 #include "sofa.h"
36 
37 class CEObserver ;
38 
40 enum class CECoordinateType
41 {
42  CIRS=0,
43  ICRS=1,
44  GALACTIC=2,
45  OBSERVED=3
46 };
47 
48 // Initiate the class that holds the coordinate information
49 class CECoordinates {
50 
51  friend bool operator==(const CECoordinates& lhs, const CECoordinates& rhs);
52  friend bool operator!=(const CECoordinates& lhs, const CECoordinates& rhs);
53 
54 public:
55 
56  /****** CONSTRUCTORS ******/
57  CECoordinates() ;
58  CECoordinates(const CEAngle& xcoord,
59  const CEAngle& ycoord,
60  const CECoordinateType& coord_type=CECoordinateType::ICRS) ;
61  CECoordinates(const std::vector<double>& xcoord,
62  const std::vector<double>& ycoord,
63  const CECoordinateType& coord_type=CECoordinateType::ICRS);
64  CECoordinates(const CECoordinateType& coord_type) ;
65  CECoordinates(const CECoordinates& other) ;
66  virtual ~CECoordinates() ;
67 
68  CECoordinates& operator=(const CECoordinates& other);
69 
70  /*********************************************************
71  * Angular separation between two coordinate positions
72  *********************************************************/
73  virtual CEAngle AngularSeparation(const CECoordinates& coords) const;
74  static CEAngle AngularSeparation(const CECoordinates& coords1,
75  const CECoordinates& coords2);
76  static CEAngle AngularSeparation(const CEAngle& xcoord_first,
77  const CEAngle& ycoord_first,
78  const CEAngle& xcoord_second,
79  const CEAngle& ycoord_second);
80 
81  /**********************************************************
82  * Methods for accessing the coordinate information
83  **********************************************************/
84 
85  virtual CEAngle XCoord(const double& jd=CppEphem::julian_date_J2000()) const;
86  virtual CEAngle YCoord(const double& jd=CppEphem::julian_date_J2000()) const;
87 
88  // The following methods will be depricated in the future
89  virtual double XCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const;
90  virtual double XCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const;
91  virtual double YCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const;
92  virtual double YCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const;
93 
94  // Convert an angle into Hours:minutes:seconds format
95  static std::vector<double> GetHMS(const double& angle,
96  const CEAngleType& angle_type = CEAngleType::DEGREES);
97  static std::vector<double> GetDMS(const double& angle,
98  const CEAngleType& angle_type = CEAngleType::DEGREES);
99  static double HMSToAngle(const std::vector<double>& angle,
100  const CEAngleType& return_type = CEAngleType::DEGREES);
101  static double DMSToAngle(const std::vector<double>& angle,
102  const CEAngleType& return_type = CEAngleType::DEGREES);
103 
104  // Return coordinate system
105  CECoordinateType GetCoordSystem(void) const;
106 
107  /**********************************************************
108  * Methods for converting between coordinate types
109  **********************************************************/
110  // Note that whenever a date is required, the default will be set to the
111  // start of the J2000 epoch (January 1, 2000 at 12:00 GMT). This corresponds
112  // to the Julian Date of 2451545.0.
113 
114  // Convert from CIRS to other coordinates
115  static void CIRS2ICRS(double input_ra, double input_dec, double *return_ra, double *return_dec,
116  const CEDate& date=CEDate(),
117  const CEAngleType& angle_type=CEAngleType::RADIANS);
118  static void CIRS2Galactic(double ra, double dec, double *glon, double *glat,
119  const CEDate& date=CEDate(),
120  const CEAngleType& angle_type=CEAngleType::RADIANS);
121  static int CIRS2Observed(double ra,
122  double dec,
123  double *az,
124  double *zen,
125  const CEDate& date,
126  const CEObserver& observer,
127  const CEAngleType& angle_type=CEAngleType::RADIANS,
128  double wavelength_um=0.5,
129  double *observed_ra=nullptr,
130  double *observed_dec=nullptr,
131  double *hour_angle=nullptr);
132 
133  // Convert from ICRS to other coordinates
134  static void ICRS2CIRS(double input_ra,
135  double input_dec,
136  double *return_ra,
137  double *return_dec,
138  const CEDate& date=CEDate(),
139  const CEAngleType& angle_type=CEAngleType::RADIANS);
140  static void ICRS2Galactic(double ra,
141  double dec,
142  double *glon,
143  double *glat,
144  const CEAngleType& angle_type=CEAngleType::RADIANS);
145  static int ICRS2Observed(double ra,
146  double dec,
147  double *az,
148  double *zen,
149  const CEDate& date,
150  const CEObserver& observer,
151  const CEAngleType& angle_type=CEAngleType::RADIANS,
152  double wavelength_um=0.5,
153  double *observed_ra=nullptr,
154  double *observed_dec=nullptr,
155  double *hour_angle=nullptr);
156 
157  // Convert from GALACTIC to other coordinates
158  static void Galactic2CIRS(double glon, double glat, double *ra, double *dec,
159  const CEDate& date=CEDate(),
160  const CEAngleType& angle_type=CEAngleType::RADIANS);
161  static void Galactic2ICRS(double glon, double glat, double *ra, double *dec,
162  const CEAngleType& angle_type=CEAngleType::RADIANS);
163  static int Galactic2Observed(double glon,
164  double glat,
165  double *az,
166  double *zen,
167  const CEDate& date,
168  const CEObserver& observer,
169  const CEAngleType& angle_type=CEAngleType::RADIANS,
170  double wavelength_um=0.5,
171  double *observed_glon=nullptr,
172  double *observed_glat=nullptr,
173  double *hour_angle=nullptr);
174  // Convert from OBSERVED to other coordinates
175  static int Observed2CIRS(double az,
176  double zen,
177  double *ra,
178  double *dec,
179  const CEDate& date,
180  const CEObserver& observer,
181  const CEAngleType& angle_type=CEAngleType::RADIANS);
182  static int Observed2ICRS(double az,
183  double zen,
184  double *ra,
185  double *dec,
186  const CEDate& date,
187  const CEObserver& observer,
188  const CEAngleType& angle_type=CEAngleType::RADIANS);
189  static int Observed2Galactic(double az,
190  double zen,
191  double *glon,
192  double *glat,
193  const CEDate& date,
194  const CEObserver& observer,
195  const CEAngleType& angle_type=CEAngleType::RADIANS);
196 
197 
198  // The following are provided to allow converting to OBSERVED
199  // azimuth, zenith without the need to use the CEObserver class.
200  // These are the methods that are ultimately called by the above
201  // methods that do use a CEObserver object
202 
203  /* Raw methods for converting CIRS <-> Observed (observer specific) coordinates. */
204  static int CIRS2Observed(double ra,
205  double dec,
206  double *az,
207  double *zen,
208  double julian_date,
209  double longitude,
210  double latitude,
211  double elevation_m=0.0,
212  double pressure_hPa=-1.0,
213  double temperature_celsius=-1000,
214  double relative_humidity=0.0,
215  double dut1=0.0,
216  double xp=0.0,
217  double yp=0.0,
218  double wavelength_um=0.5,
219  double *observed_ra=nullptr,
220  double *observed_dec=nullptr,
221  double *hour_angle=nullptr);
222  static int Observed2CIRS(double az,
223  double zen,
224  double *ra,
225  double *dec,
226  double julian_date,
227  double longitude,
228  double latitude,
229  double elevation_m=0.0,
230  double pressure_hPa=-1.0,
231  double temperature_celsius=-1000,
232  double relative_humidity=0.0,
233  double dut1=0.0,
234  double xp=0.0,
235  double yp=0.0,
236  double wavelength_um=0.5);
237 
238  /* Raw methods for converting ICRS <-> Observed (observer specific) coordinates. */
239  static int ICRS2Observed(double ra,
240  double dec,
241  double *az,
242  double *zen,
243  double julian_date,
244  double longitude,
245  double latitude,
246  double elevation_m=0.0,
247  double pressure_hPa=-1.0,
248  double temperature_celsius=-1000,
249  double relative_humidity=0.0,
250  double dut1=0.0,
251  double xp=0.0, double yp=0.0,
252  double wavelength_um=0.5,
253  double *observed_ra=nullptr,
254  double *observed_dec=nullptr,
255  double *hour_angle=nullptr);
256  static int Observed2ICRS(double az,
257  double zen,
258  double *ra,
259  double *dec,
260  double julian_date,
261  double longitude,
262  double latitude,
263  double elevation_m=0.0,
264  double pressure_hPa=-1.0,
265  double temperature_celsius=-1000,
266  double relative_humidity=0.0,
267  double dut1=0.0,
268  double xp=0.0, double yp=0.0,
269  double wavelength_um=0.5);
270 
271  /* Raw methods for converting Galactic <-> Observed (observer specific) coordinates. */
272  static int Galactic2Observed(double glon,
273  double glat,
274  double *az,
275  double *zen,
276  double julian_date,
277  double longitude,
278  double latitude,
279  double elevation_m=0.0,
280  double pressure_hPa=-1.0,
281  double temperature_celsius=-1000,
282  double relative_humidity=0.0,
283  double dut1=0.0,
284  double xp=0.0, double yp=0.0,
285  double wavelength_um=0.50,
286  double *observed_glon=nullptr,
287  double *observed_glat=nullptr,
288  double *hour_angle=nullptr);
289  static int Observed2Galactic(double az,
290  double zen,
291  double *glon,
292  double *glat,
293  double julian_date,
294  double longitude,
295  double latitude,
296  double elevation_m=0.0,
297  double pressure_hPa=-1.0,
298  double temperature_celsius=-1000,
299  double relative_humidity=0.0,
300  double dut1=0.0,
301  double xp=0.0, double yp=0.0,
302  double wavelength_um=0.50);
303 
304  virtual CECoordinates GetObservedCoords(const double& julian_date,
305  const double& longitude,
306  const double& latitude,
307  const double& elevation_m=0.0,
308  const double& pressure_hPa=-1.0,
309  const double& temperature_celsius=-1000,
310  const double& relative_humidity=0.0,
311  const double& dut1=0.0,
312  const double& xp=0.0, const double& yp=0.0,
313  const double& wavelength_um=0.5) const;
314  virtual CECoordinates GetObservedCoords(const CEDate& date,
315  const CEObserver& observer) const;
316 
317  /*********************************************************
318  * More generic methods for converting between coordinate types
319  *********************************************************/
320  CECoordinates ConvertTo(CECoordinateType output_coord_type,
321  const CEObserver& observer,
322  const CEDate& date=CEDate::CurrentJD());
323 
324  CECoordinates ConvertTo(CECoordinateType output_coord_type,
325  double jd=CEDate::CurrentJD(),
326  double longitude=0.0,
327  double latitude=0.0,
328  double elevation_m=0.0,
329  double pressure_hPa=-1.0,
330  double temperature_celsius=-1000,
331  double relative_humidity=0.0,
332  double dut1=0.0,
333 
334  double xp=0.0, double yp=0.0,
335  double wavelength_um=0.5);
336 
338  double longitude=0.0,
339  double latitude=0.0,
340  double elevation_m=0.0,
341  double pressure_hPa=-1.0,
342  double temperature_celsius=-1000,
343  double relative_humidity=0.0,
344  double dut1=0.0,
345  double xp=0.0, double yp=0.0,
346  double wavelength_um=0.5);
348  double longitude=0.0,
349  double latitude=0.0,
350  double elevation_m=0.0,
351  double pressure_hPa=-1.0,
352  double temperature_celsius=-1000,
353  double relative_humidity=0.0,
354  double dut1=0.0,
355  double xp=0.0, double yp=0.0,
356  double wavelength_um=0.5);
358  double longitude=0.0,
359  double latitude=0.0,
360  double elevation_m=0.0,
361  double pressure_hPa=-1.0,
362  double temperature_celsius=-1000,
363  double relative_humidity=0.0,
364  double dut1=0.0,
365  double xp=0.0, double yp=0.0,
366  double wavelength_um=0.5);
368  double longitude=0.0,
369  double latitude=0.0,
370  double elevation_m=0.0,
371  double pressure_hPa=-1.0,
372  double temperature_celsius=-1000,
373  double relative_humidity=0.0,
374  double dut1=0.0,
375  double xp=0.0, double yp=0.0,
376  double wavelength_um=0.5);
377 
378  /*********************************************************
379  * Methods for setting the coordinates of this object
380  *********************************************************/
381  virtual void SetCoordinates(const CEAngle& xcoord,
382  const CEAngle& ycoord,
383  const CECoordinateType& coord_type = CECoordinateType::ICRS);
384  virtual void SetCoordinates(const CECoordinates& coords);
385 
386  // Support methods
387  std::string print(void) const;
388 
389 protected:
390  // Coordinate variables
391  mutable CEAngle xcoord_; //<! X coordinate
392  mutable CEAngle ycoord_; //<! Y coordinate
393  CECoordinateType coord_type_; //<! Coordinate system to which 'xcoord_' and 'ycoord_' belong.
394 
395 private:
396  /*********************************************************
397  * Private methods
398  *********************************************************/
399  void copy_members(const CECoordinates& other);
400  void free_members(void);
401  void init_members(void);
402 };
403 
404 
405 /**********************************************************************/
411 inline
412 CEAngle CECoordinates::XCoord(const double& jd) const
413 {
414  return xcoord_;
415 }
416 
417 
418 /**********************************************************************/
424 inline
425 CEAngle CECoordinates::YCoord(const double& jd) const
426 {
427  return ycoord_;
428 }
429 
430 
431 /**********************************************************************/
437 inline
438 double CECoordinates::XCoordinate_Rad(double jd) const
439 {
440  return XCoord(jd).Rad();
441 }
442 
443 
444 /**********************************************************************/
450 inline
451 double CECoordinates::XCoordinate_Deg(double jd) const
452 {
453  return XCoord(jd).Deg();
454 }
455 
456 
457 /**********************************************************************/
463 inline
464 double CECoordinates::YCoordinate_Rad(double jd) const
465 {
466  return YCoord(jd).Rad();
467 }
468 
469 /**********************************************************************/
475 inline
476 double CECoordinates::YCoordinate_Deg(double jd) const
477 {
478  return YCoord(jd).Deg();
479 }
480 
481 
482 /**********************************************************************/
486 inline
488 {
489  return coord_type_;
490 }
491 
492 #endif /* CECoordinates_h */
CECoordinates::DMSToAngle
static double DMSToAngle(const std::vector< double > &angle, const CEAngleType &return_type=CEAngleType::DEGREES)
Convert a given angle vector from {degrees, minutes, seconds} to an angle.
Definition: CECoordinates.cpp:1662
CECoordinates::Observed2Galactic
static int Observed2Galactic(double az, double zen, double *glon, double *glat, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS)
Convert Observed -> Galactic coordinates.
Definition: CECoordinates.cpp:653
CECoordinates::AngularSeparation
virtual CEAngle AngularSeparation(const CECoordinates &coords) const
Get the angular separation between the coordinates represented by this object and another coordinate ...
Definition: CECoordinates.cpp:1169
CECoordinates::operator=
CECoordinates & operator=(const CECoordinates &other)
Overloaded '=' (assignment) operator.
Definition: CECoordinates.cpp:122
CENamespace.h
CECoordinates::ycoord_
CEAngle ycoord_
Definition: CECoordinates.h:386
CECoordinateType::GALACTIC
Galacitc longitude, latitude.
CEDate
Definition: CEDate.h:43
CppEphem::yp
double yp(const double &mjd)
Polar motion (x) for a given modified julian date (radians)
Definition: CENamespace.cpp:160
CECoordinates::CIRS2Galactic
static void CIRS2Galactic(double ra, double dec, double *glon, double *glat, const CEDate &date=CEDate(), const CEAngleType &angle_type=CEAngleType::RADIANS)
CIRS -> Galactic coordinate conversion.
Definition: CECoordinates.cpp:181
CECoordinates::ICRS2Galactic
static void ICRS2Galactic(double ra, double dec, double *glon, double *glat, const CEAngleType &angle_type=CEAngleType::RADIANS)
ICRS -> Galactic coordinate conversion (uses the SOFA 'iauIcrs2g' function)
Definition: CECoordinates.cpp:336
CECoordinates::XCoordinate_Deg
virtual double XCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const
Returns x coordinate at given julian date in degrees.
Definition: CECoordinates.h:450
CECoordinates::ConvertToICRS
CECoordinates ConvertToICRS(double jd=CEDate::CurrentJD(), double longitude=0.0, double latitude=0.0, double elevation_m=0.0, double pressure_hPa=-1.0, double temperature_celsius=-1000, double relative_humidity=0.0, double dut1=0.0, double xp=0.0, double yp=0.0, double wavelength_um=0.5)
Convert these coordinates to ICRS coordinates.
Definition: CECoordinates.cpp:1418
CECoordinates::GetObservedCoords
virtual CECoordinates GetObservedCoords(const double &julian_date, const double &longitude, const double &latitude, const double &elevation_m=0.0, const double &pressure_hPa=-1.0, const double &temperature_celsius=-1000, const double &relative_humidity=0.0, const double &dut1=0.0, const double &xp=0.0, const double &yp=0.0, const double &wavelength_um=0.5) const
Return the local sky coordinates of this object as a CECoordinates object.
Definition: CECoordinates.cpp:1087
CECoordinates::CIRS2ICRS
static void CIRS2ICRS(double input_ra, double input_dec, double *return_ra, double *return_dec, const CEDate &date=CEDate(), const CEAngleType &angle_type=CEAngleType::RADIANS)
CIRS -> ICRS coordinate conversion.
Definition: CECoordinates.cpp:144
CECoordinates::Observed2CIRS
static int Observed2CIRS(double az, double zen, double *ra, double *dec, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS)
Convert Observed -> CIRS coordinates.
Definition: CECoordinates.cpp:566
CEAngleType
CEAngleType
Definition: CEAngle.h:30
CECoordinates::YCoord
virtual CEAngle YCoord(const double &jd=CppEphem::julian_date_J2000()) const
Return y coordinate at given Julian date.
Definition: CECoordinates.h:424
CECoordinates::XCoordinate_Rad
virtual double XCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const
Return x coordinate at given Julian date in radians.
Definition: CECoordinates.h:437
CECoordinates::free_members
void free_members(void)
Cleanup data members that need to be freed or cleared.
Definition: CECoordinates.cpp:1730
CECoordinates::ConvertToGalactic
CECoordinates ConvertToGalactic(double jd=CEDate::CurrentJD(), double longitude=0.0, double latitude=0.0, double elevation_m=0.0, double pressure_hPa=-1.0, double temperature_celsius=-1000, double relative_humidity=0.0, double dut1=0.0, double xp=0.0, double yp=0.0, double wavelength_um=0.5)
Convert these coordinates to GALACTIC coordinates.
Definition: CECoordinates.cpp:1471
CECoordinates::YCoordinate_Deg
virtual double YCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const
Returns y coordinate at given Julian date in degrees.
Definition: CECoordinates.h:475
CECoordinateType::ICRS
RA, Dec (referenced at the barycenter of the solarsystem)
CEDate.h
CECoordinates::ConvertTo
CECoordinates ConvertTo(CECoordinateType output_coord_type, const CEObserver &observer, const CEDate &date=CEDate::CurrentJD())
Convert these coordinates to another coordinate system NOTE: If this object is not OBSERVED coordinat...
Definition: CECoordinates.cpp:1261
CECoordinates::Observed2ICRS
static int Observed2ICRS(double az, double zen, double *ra, double *dec, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS)
Convert Observed -> ICRS coordinates.
Definition: CECoordinates.cpp:608
CECoordinates::copy_members
void copy_members(const CECoordinates &other)
Copy data members from another CECoordinates object.
Definition: CECoordinates.cpp:1719
CECoordinates::xcoord_
CEAngle xcoord_
Definition: CECoordinates.h:385
CECoordinates::GetDMS
static std::vector< double > GetDMS(const double &angle, const CEAngleType &angle_type=CEAngleType::DEGREES)
Convert a given angle into degrees, arcminutes, arcseconds.
Definition: CECoordinates.cpp:1594
CEException.h
CECoordinates::operator==
friend bool operator==(const CECoordinates &lhs, const CECoordinates &rhs)
Compare two coordinate objects.
Definition: CECoordinates.cpp:1751
CEAngle
Definition: CEAngle.h:38
CECoordinates::ConvertToCIRS
CECoordinates ConvertToCIRS(double jd=CEDate::CurrentJD(), double longitude=0.0, double latitude=0.0, double elevation_m=0.0, double pressure_hPa=-1.0, double temperature_celsius=-1000, double relative_humidity=0.0, double dut1=0.0, double xp=0.0, double yp=0.0, double wavelength_um=0.5)
Convert these coordinates to CIRS coordinates.
Definition: CECoordinates.cpp:1358
CECoordinates::CECoordinates
CECoordinates()
Default constructor.
Definition: CECoordinates.cpp:38
CECoordinates::GetCoordSystem
CECoordinateType GetCoordSystem(void) const
Return coordinate system.
Definition: CECoordinates.h:486
CppEphem::julian_date_J2000
double julian_date_J2000()
Julian Date corresponding to J2000.
Definition: CENamespace.h:67
CEAngle::Rad
static CEAngle Rad(const double &angle)
Return angle constructed from a radians angle.
Definition: CEAngle.cpp:340
CEAngle.h
CEAngleType::DEGREES
CECoordinates::ConvertToObserved
CECoordinates ConvertToObserved(double jd=CEDate::CurrentJD(), double longitude=0.0, double latitude=0.0, double elevation_m=0.0, double pressure_hPa=-1.0, double temperature_celsius=-1000, double relative_humidity=0.0, double dut1=0.0, double xp=0.0, double yp=0.0, double wavelength_um=0.5)
Convert these coordinates to observed coordinates.
Definition: CECoordinates.cpp:1522
CECoordinates::Galactic2CIRS
static void Galactic2CIRS(double glon, double glat, double *ra, double *dec, const CEDate &date=CEDate(), const CEAngleType &angle_type=CEAngleType::RADIANS)
Galactic -> CIRS coordinate conversion.
Definition: CECoordinates.cpp:441
CECoordinates::CIRS2Observed
static int CIRS2Observed(double ra, double dec, double *az, double *zen, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS, double wavelength_um=0.5, double *observed_ra=nullptr, double *observed_dec=nullptr, double *hour_angle=nullptr)
CIRS -> Observed (or observer specific) coordinate conversion This function takes in verious observat...
Definition: CECoordinates.cpp:224
CECoordinates::Galactic2Observed
static int Galactic2Observed(double glon, double glat, double *az, double *zen, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS, double wavelength_um=0.5, double *observed_glon=nullptr, double *observed_glat=nullptr, double *hour_angle=nullptr)
Galactic -> Observed (i.e.
Definition: CECoordinates.cpp:508
CECoordinates::HMSToAngle
static double HMSToAngle(const std::vector< double > &angle, const CEAngleType &return_type=CEAngleType::DEGREES)
Convert from {hours, minutes, seconds} to an angle.
Definition: CECoordinates.cpp:1636
CECoordinates::GetHMS
static std::vector< double > GetHMS(const double &angle, const CEAngleType &angle_type=CEAngleType::DEGREES)
Convert a given angle into hours, minutes, seconds.
Definition: CECoordinates.cpp:1613
CECoordinates
Definition: CECoordinates.h:48
CECoordinates::coord_type_
CECoordinateType coord_type_
Definition: CECoordinates.h:387
CEDate::CurrentJD
static double CurrentJD()
Static method for getting the current Julian date.
Definition: CEDate.cpp:632
CEAngleType::RADIANS
CECoordinates::init_members
void init_members(void)
Set initial values and allocate memory for data members.
Definition: CECoordinates.cpp:1739
CEObserver
Definition: CEObserver.h:28
CECoordinates::XCoord
virtual CEAngle XCoord(const double &jd=CppEphem::julian_date_J2000()) const
Return x coordinate at given Julian date.
Definition: CECoordinates.h:411
CECoordinateType::CIRS
RA, Dec (referenced at the center of the Earth)
CECoordinateType
CECoordinateType
The following enum specifies what coordinates this object represents.
Definition: CECoordinates.h:39
CECoordinates::ICRS2CIRS
static void ICRS2CIRS(double input_ra, double input_dec, double *return_ra, double *return_dec, const CEDate &date=CEDate(), const CEAngleType &angle_type=CEAngleType::RADIANS)
ICRS -> CIRS coordinate conversion.
Definition: CECoordinates.cpp:287
CECoordinates::SetCoordinates
virtual void SetCoordinates(const CEAngle &xcoord, const CEAngle &ycoord, const CECoordinateType &coord_type=CECoordinateType::ICRS)
Set the coordinates of this object.
Definition: CECoordinates.cpp:1682
CppEphem::dut1
double dut1(const double &mjd)
Return dut1 based on a given modified julian date (seconds)
Definition: CENamespace.cpp:111
CECoordinates::operator!=
friend bool operator!=(const CECoordinates &lhs, const CECoordinates &rhs)
Compare two coordinate objects.
Definition: CECoordinates.cpp:1778
CECoordinates::YCoordinate_Rad
virtual double YCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const
Returns y coordinate at given Julian date in radians.
Definition: CECoordinates.h:463
CECoordinates::print
std::string print(void) const
Generate a message string that specifies the information about this coordinate.
Definition: CECoordinates.cpp:1705
CECoordinateType::OBSERVED
Azimuth, Zenith (requires additional observer information)
CEAngle::Deg
static CEAngle Deg(const double &angle)
Return angle (radians) constructed from a degree angle.
Definition: CEAngle.cpp:317
CECoordinates::ICRS2Observed
static int ICRS2Observed(double ra, double dec, double *az, double *zen, const CEDate &date, const CEObserver &observer, const CEAngleType &angle_type=CEAngleType::RADIANS, double wavelength_um=0.5, double *observed_ra=nullptr, double *observed_dec=nullptr, double *hour_angle=nullptr)
ICRS -> Observed coordinate conversion.
Definition: CECoordinates.cpp:373
CECoordinates::Galactic2ICRS
static void Galactic2ICRS(double glon, double glat, double *ra, double *dec, const CEAngleType &angle_type=CEAngleType::RADIANS)
Galactic -> ICRS coordinate conversion.
Definition: CECoordinates.cpp:474
CppEphem::xp
double xp(const double &mjd)
Polar motion (x) for a given modified julian date (radians)
Definition: CENamespace.cpp:148
CECoordinates::~CECoordinates
virtual ~CECoordinates()
Destructor.
Definition: CECoordinates.cpp:112