Go to the documentation of this file.
22 #ifndef CECoordinates_h
23 #define CECoordinates_h
62 const std::vector<double>& ycoord,
95 static std::vector<double>
GetHMS(
const double& angle,
97 static std::vector<double>
GetDMS(
const double& angle,
99 static double HMSToAngle(
const std::vector<double>& angle,
101 static double DMSToAngle(
const std::vector<double>& angle,
115 static void CIRS2ICRS(
double input_ra,
double input_dec,
double *return_ra,
double *return_dec,
118 static void CIRS2Galactic(
double ra,
double dec,
double *glon,
double *glat,
128 double wavelength_um=0.5,
129 double *observed_ra=
nullptr,
130 double *observed_dec=
nullptr,
131 double *hour_angle=
nullptr);
152 double wavelength_um=0.5,
153 double *observed_ra=
nullptr,
154 double *observed_dec=
nullptr,
155 double *hour_angle=
nullptr);
158 static void Galactic2CIRS(
double glon,
double glat,
double *ra,
double *dec,
161 static void Galactic2ICRS(
double glon,
double glat,
double *ra,
double *dec,
170 double wavelength_um=0.5,
171 double *observed_glon=
nullptr,
172 double *observed_glat=
nullptr,
173 double *hour_angle=
nullptr);
211 double elevation_m=0.0,
212 double pressure_hPa=-1.0,
213 double temperature_celsius=-1000,
214 double relative_humidity=0.0,
218 double wavelength_um=0.5,
219 double *observed_ra=
nullptr,
220 double *observed_dec=
nullptr,
221 double *hour_angle=
nullptr);
229 double elevation_m=0.0,
230 double pressure_hPa=-1.0,
231 double temperature_celsius=-1000,
232 double relative_humidity=0.0,
236 double wavelength_um=0.5);
246 double elevation_m=0.0,
247 double pressure_hPa=-1.0,
248 double temperature_celsius=-1000,
249 double relative_humidity=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);
263 double elevation_m=0.0,
264 double pressure_hPa=-1.0,
265 double temperature_celsius=-1000,
266 double relative_humidity=0.0,
268 double xp=0.0,
double yp=0.0,
269 double wavelength_um=0.5);
279 double elevation_m=0.0,
280 double pressure_hPa=-1.0,
281 double temperature_celsius=-1000,
282 double relative_humidity=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);
296 double elevation_m=0.0,
297 double pressure_hPa=-1.0,
298 double temperature_celsius=-1000,
299 double relative_humidity=0.0,
301 double xp=0.0,
double yp=0.0,
302 double wavelength_um=0.50);
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;
326 double longitude=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,
334 double xp=0.0,
double yp=0.0,
335 double wavelength_um=0.5);
338 double longitude=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,
345 double xp=0.0,
double yp=0.0,
346 double wavelength_um=0.5);
348 double longitude=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,
355 double xp=0.0,
double yp=0.0,
356 double wavelength_um=0.5);
358 double longitude=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,
365 double xp=0.0,
double yp=0.0,
366 double wavelength_um=0.5);
368 double longitude=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,
375 double xp=0.0,
double yp=0.0,
376 double wavelength_um=0.5);
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.
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.
virtual CEAngle AngularSeparation(const CECoordinates &coords) const
Get the angular separation between the coordinates represented by this object and another coordinate ...
CECoordinates & operator=(const CECoordinates &other)
Overloaded '=' (assignment) operator.
Galacitc longitude, latitude.
double yp(const double &mjd)
Polar motion (x) for a given modified julian date (radians)
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.
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)
virtual double XCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const
Returns x coordinate at given julian date in degrees.
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.
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.
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.
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.
virtual CEAngle YCoord(const double &jd=CppEphem::julian_date_J2000()) const
Return y coordinate at given Julian date.
virtual double XCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const
Return x coordinate at given Julian date in radians.
void free_members(void)
Cleanup data members that need to be freed or cleared.
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.
virtual double YCoordinate_Deg(double jd=CppEphem::julian_date_J2000()) const
Returns y coordinate at given Julian date in degrees.
RA, Dec (referenced at the barycenter of the solarsystem)
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...
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.
void copy_members(const CECoordinates &other)
Copy data members from another CECoordinates object.
static std::vector< double > GetDMS(const double &angle, const CEAngleType &angle_type=CEAngleType::DEGREES)
Convert a given angle into degrees, arcminutes, arcseconds.
friend bool operator==(const CECoordinates &lhs, const CECoordinates &rhs)
Compare two coordinate objects.
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.
CECoordinates()
Default constructor.
CECoordinateType GetCoordSystem(void) const
Return coordinate system.
double julian_date_J2000()
Julian Date corresponding to J2000.
static CEAngle Rad(const double &angle)
Return angle constructed from a radians angle.
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.
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.
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...
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.
static double HMSToAngle(const std::vector< double > &angle, const CEAngleType &return_type=CEAngleType::DEGREES)
Convert from {hours, minutes, seconds} to an angle.
static std::vector< double > GetHMS(const double &angle, const CEAngleType &angle_type=CEAngleType::DEGREES)
Convert a given angle into hours, minutes, seconds.
CECoordinateType coord_type_
static double CurrentJD()
Static method for getting the current Julian date.
void init_members(void)
Set initial values and allocate memory for data members.
virtual CEAngle XCoord(const double &jd=CppEphem::julian_date_J2000()) const
Return x coordinate at given Julian date.
RA, Dec (referenced at the center of the Earth)
CECoordinateType
The following enum specifies what coordinates this object represents.
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.
virtual void SetCoordinates(const CEAngle &xcoord, const CEAngle &ycoord, const CECoordinateType &coord_type=CECoordinateType::ICRS)
Set the coordinates of this object.
double dut1(const double &mjd)
Return dut1 based on a given modified julian date (seconds)
friend bool operator!=(const CECoordinates &lhs, const CECoordinates &rhs)
Compare two coordinate objects.
virtual double YCoordinate_Rad(double jd=CppEphem::julian_date_J2000()) const
Returns y coordinate at given Julian date in radians.
std::string print(void) const
Generate a message string that specifies the information about this coordinate.
Azimuth, Zenith (requires additional observer information)
static CEAngle Deg(const double &angle)
Return angle (radians) constructed from a degree angle.
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.
static void Galactic2ICRS(double glon, double glat, double *ra, double *dec, const CEAngleType &angle_type=CEAngleType::RADIANS)
Galactic -> ICRS coordinate conversion.
double xp(const double &mjd)
Polar motion (x) for a given modified julian date (radians)
virtual ~CECoordinates()
Destructor.