aephem.h

Go to the documentation of this file.

//==============================================================================
// AEPHEM - an astronomical ephemeris and reduction library.
// Copyright 2012 Adam Hincks, Canadian Institute for Theoretical Astrophysics
//
// This file is part of AEPHEM.
//
// AEPHEM is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// AEPHEM is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with AEPHEM.  If not, see <http://www.gnu.org/licenses/>.
//==============================================================================

#ifndef AEPHEM_H
#define AEPHEM_H

#ifdef __cplusplus
  extern "C" {
#endif

#include <math.h>
#include <stdio.h>

// Constants.
#define AE_DTR          1.7453292519943295769e-2 
#define AE_RTD          5.7295779513082320877e1  
#define AE_RTS          2.0626480624709635516e5  
#define AE_STR          4.8481368110953599359e-6 
#define AE_STD          (1.0 / 3600.0)           
#define AE_STDAY        0.00001157407407407407   
#define AE_COSSUN       -0.014543897651582657    
#define AE_COSZEN       -9.8900378587411476e-3   
#define AE_J2000        2451545.0                
#define AE_B1950        2433282.423              
#define AE_J1900        2415020.0                
#define AE_MJD_START    2400000.5
#define AE_AU           1.49597870691e8          
#define AE_CLIGHT       2.99792458e5             
#define AE_S_PER_D      86400.0                  
#define AE_D_PER_S      (1.0 / AE_S_PER_D)       
#define AE_CLIGHT_AUD   (86400.0 * AE_CLIGHT / AE_AU)
#define AE_R_EARTH      6378137.0                
#define AE_R_EARTH_AU   (0.001 * AE_R_EARTH / AE_AU) 
#define AE_E_M_RAT      81.300585                
#define AE_FLAT         298.257222101
#define AE_T_CMB        2.725
#define AE_JPL_HEAD_SIZE (5 * sizeof(double) + 41 * sizeof(int32_t))
#define AE_JPL_CONST_LEN 6      

enum ae_retcode_t {
  AE_RET_BAD_PATH = 1,          
  AE_RET_UNEXPECTED_EOF,        
  AE_RET_READ_ERROR,            
  AE_RET_BAD_JPL_HEADER,        
  AE_RET_BAD_JPL_HEADER_KSIZE,  
  AE_RET_BAD_JPL_HEADER_GROUP,  
  AE_RET_BAD_JPL_HEADER_NUM, 
  AE_RET_BAD_JPL_HEADER_RANGE,  
  AE_RET_BAD_JPL_NOT_ASCII,     
  AE_RET_BAD_JPL_NOT_BIN,       
  AE_RET_BAD_JPL_BLOCK,         
  AE_RET_BAD_JPL_CORRUPT,       
  AE_RET_BAD_JPL_DATE,          
  AE_RET_BAD_JPL_FLOAT,         
  AE_RET_BAD_JPL_CHEBY_TIME,    
  AE_RET_BAD_OBJ_INDEX,         
  AE_RET_OBJ_NOT_FOUND,         
  AE_RET_CFITSIO_NOT_ENABLED,
  AE_RET_CFITSIO_NO_KEYWORD,    
  AE_RET_CFITSIO_BAD_KEYWORD,   
  AE_RET_WCS_UNKNOWN_PROJ,      
  AE_RET_CFITSIO_ERROR,         
  AE_RET_BAD_NATIVE_POLE,       
  AE_RET_NO_INVERSE,            
  AE_RET_UNDER_CONSTRUCTION
};

#define AE_N_SS_BODIES  16
#define AE_N_SS_BODIES_JPL  13
enum ae_ss_bodies_t {
  AE_SS_MERCURY = 0,    
  AE_SS_VENUS,          
  AE_SS_EMBARY,         
  AE_SS_MARS,           
  AE_SS_JUPITER,        
  AE_SS_SATURN,         
  AE_SS_URANUS,         
  AE_SS_NEPTUNE,        
  AE_SS_PLUTO,          
  AE_SS_MOON_EMB,       
  AE_SS_SUN,            
  AE_SS_NUTATION,       
  AE_SS_LIBRATION,      
  AE_SS_EARTH,          
  AE_SS_MOON,           
  AE_SS_SSBARY          
};

enum ae_precess_direction_t {
  AE_FROM_J2000 = -1,   
  AE_TO_J2000 = 1       
};

enum ae_aberration_direction_t {
  AE_ADD_ABERRATION = -1,  
  AE_REMOVE_ABERRATION = 1 
};

enum ae_physical_use_d_or_t {
  AE_PHYSICAL_END = 0,   
  AE_PHYSICAL_D,         
  AE_PHYSICAL_T          
};

#define AE_PLANTBL_N_HARMONIC     18  
struct ae_plantbl_t {
  char maxargs;     
  char max_harmonic[AE_PLANTBL_N_HARMONIC]; 
  char max_power_of_t;  
  char *arg_tbl;        
  void *lon_tbl;        
  void *lat_tbl;        
  void *rad_tbl;        
  double distance;      
  double timescale;     
  double trunclvl;      
};

struct ae_orbit_t {
  char name[32];  
  double epoch;     
  double i;             
  double W;             
  double w;             
  double a;             
  double dm;        
  double ecc;       
  double M;             
  double equinox;       

  // The following used by perterbation formulas.
  struct ae_plantbl_t *ptable;  
  double L;         
  double r;             
  double plat;      
};

struct ae_star_t {
  char name[32];  
  double epoch;     
  double ra;        
  double dec;       
  double px;        
  double mura;      
  double mudec;     
  double v;         
  double mag;       
};

struct ae_jpl_handle_t {
  double start;             
  double end;               
  double step;              
  int n_const;              
  double au;                
  double em_ratio;          
  int ipt[13][3];           
  int ephemeris_version;    
  char **const_name;     
  double *const_val;        

  // The following are computed within the code.
  FILE *fp;                 
  int is_bin;               
  int n_coeff;              
  int swop;                 
  int rec_size;             
  int last_block_num;       
  int max_cheby;            
  long ascii_start;     
  long ascii_len_block;     
  double sun_pv[6];         
  double sun_pv_t;          
  double *block;            
  double *pc;               
  double *vc;               
  double last_x;            
};

struct ae_physical_term_t {
  int time_var;         
  double a;             
  double b;             
  double c;             
  double d;             
};

struct ae_physical_t {
  double r_mean;                
  double r_eq;                  
  double r_pole;                
  double rms_dev_spheroid;   
  double max_elev;              
  double max_depress;           
  double pole_ra;               
  double pole_ra_t;
  struct ae_physical_term_t *ra_sin_term;
  double pole_dec;              
  double pole_dec_t;
  struct ae_physical_term_t *dec_cos_term;
  double w;                     
  double w_d;
  double w_d_sq;
  struct ae_physical_term_t *w_sin_term;
};

// Predefined orbits.
extern struct ae_orbit_t ae_orb_mercury;
extern struct ae_orbit_t ae_orb_venus;
extern struct ae_orbit_t ae_orb_earth;
extern struct ae_orbit_t ae_orb_mars;
extern struct ae_orbit_t ae_orb_jupiter;
extern struct ae_orbit_t ae_orb_saturn;
extern struct ae_orbit_t ae_orb_uranus;
extern struct ae_orbit_t ae_orb_neptune;
extern struct ae_orbit_t ae_orb_pluto;
extern struct ae_orbit_t *ae_orb_planet[];

// Predefined perturbation tables.
extern struct ae_plantbl_t ae_mer404;
extern struct ae_plantbl_t ae_ven404;
extern struct ae_plantbl_t ae_ear404;
extern struct ae_plantbl_t ae_mar404;
extern struct ae_plantbl_t ae_jup404;
extern struct ae_plantbl_t ae_sat404;
extern struct ae_plantbl_t ae_ura404;
extern struct ae_plantbl_t ae_nep404;
extern struct ae_plantbl_t ae_plu404;
extern struct ae_plantbl_t ae_mlat404;
extern struct ae_plantbl_t ae_mlr404;

// Predefined physical parameters.
extern struct ae_physical_t ae_phys_sun;
extern struct ae_physical_t ae_phys_mercury;
extern struct ae_physical_t ae_phys_venus;
extern struct ae_physical_t ae_phys_earth;
extern struct ae_physical_t ae_phys_mars;
extern struct ae_physical_t ae_phys_jupiter;
extern struct ae_physical_t ae_phys_saturn;
extern struct ae_ring_geometry_t ae_saturn_rings;
extern struct ae_physical_t ae_phys_uranus;
extern struct ae_physical_t ae_phys_neptune;
extern struct ae_physical_t ae_phys_pluto;
extern struct ae_physical_t *ae_phys_planet[];
extern struct ae_physical_term_t ae_phys_no_term[];

// Predefined names.
extern const char *ae_constel_name[];
extern const char *ae_ss_name[];

void ae_altaz_to_radec(double last, double glat, double alt, double az,
                       double *ra, double *dec);
void ae_annual_aberration(double v_earth[], double p[], int direction);
double ae_cal_to_jd(long year, int month, double day);
int ae_cat_to_constel_index(const char *cat_name, char *full_name, int len);
int ae_coord_to_constel_index(double ra, double dec, double epoch);
double ae_ctime_to_jd(double t);
double ae_ctime_to_last(double t, double delta_t, double tlong, double nutl, 
                        double eps);
double aes_ctime_to_last(double t, double tlong);
double ae_delta_t(double jd_ut1);
void ae_delta_q(double q0[], double q1[], double *dra, double *ddec);
double ae_disc_semiminor(double a, double c, double lat);
double ae_disc_solid_angle(double a, double b, double dist);
double aes_disc_solid_angle(double jd_ut1, const struct ae_orbit_t *o_orb,
                            const struct ae_orbit_t *q_orb,
                            const struct ae_physical_t *phys);
void ae_diurnal_aberration(double last, double tlat, double trho, int direction,
                           double *ra, double *dec);
void ae_diurnal_parallax(double last, double tlat, double trho, double dist,
                         double *ra, double *dec);
double ae_dut1(double jd_utc);
double ae_epsilon(double jd_tt);
void ae_fk4_to_fk5(struct ae_star_t *star);
void ae_geocentric(double jd_tt, double q[], double e[], double v_e[], 
                   double *ra, double *dec, double *dist);
void ae_geocentric_moon_from_orbit(double jd_tt, const struct ae_orbit_t *o_orb,
                                   double *ra, double *dec, double *dist);
void ae_geocentric_from_cat(double jd_tt, double e[], double v_e[],
                            const struct ae_star_t *star, double *ra,
                            double *dec);
int ae_geocentric_from_jpl(struct ae_jpl_handle_t *jh, double jd_tt, 
                           int obj_num, double *ra, double *dec, double *dist);
void ae_geocentric_from_orbit(double jd_tt, const struct ae_orbit_t *o_orb, 
                              const struct ae_orbit_t *q_orb, double *ra, 
                              double *dec, double *dist);
void ae_geocentric_sun_from_orbit(double jd_tt, const struct ae_orbit_t *o_orb,
                                   double *ra, double *dec, double *dist);
void ae_geocentric_lat(double glat, double height, double *tlat, double *trho);
int ae_read_orbit_from_cat(const char *path, const char *name, 
                           struct ae_orbit_t *orb);
int ae_read_star_from_cat(const char *path, const char *name, 
                          struct ae_star_t *star);
void ae_gmoon(double jd, double rect[], double pol[]);
double ae_gmst(double jd_ut1, double jd_tt);
double ae_gast(double jd_ut1, double jd_tt, double nutl, double eps);
double ae_lmst(double jd_ut1, double jd_tt, double tlong);
double ae_last(double jd_ut1, double jd_tt, double tlong, double nutl,
               double eps);
double aes_last(double jd_ut1, double tlong);
void ae_gplan(double jd_tt, struct ae_plantbl_t *plan, double pobj[]);
double ae_g1plan(double jd_tt, struct ae_plantbl_t *plan);
void ae_g2plan(double jd_tt, struct ae_plantbl_t *plan, double pobj[],
               double *lp_equinox);
void ae_g3plan(double jd_tt, struct ae_plantbl_t *plan, double pobj[], 
               int objnum);
void ae_jd_to_cal(double jd, int *year, int *month, double *day);
double ae_jpl_get_const_val(const struct ae_jpl_handle_t *jh,
                            const char *const_name);
int ae_jpl_init(const char *path, struct ae_jpl_handle_t *jh);
int ae_jpl_init_ascii(const char *path, const char *header_path,
                      int search_dates, struct ae_jpl_handle_t *jh);
int ae_jpl_init_bin(const char *path, struct ae_jpl_handle_t *jh);
int ae_jpl_get_coords(struct ae_jpl_handle_t *jh, double jd_tt, int obj_num,
                      double r[], double v[], char is_planetary);
void ae_jpl_close(struct ae_jpl_handle_t *jh);
void ae_kepler(double jd_tt, const struct ae_orbit_t *orb, double q[]);
double ae_light_t(double p[], double q[], int do_retardation);
int ae_light_t_jpl(struct ae_jpl_handle_t *jh, double jd_tt, double e[], 
                   int obj_num, double q[], double v_q[], char is_planetary);
void ae_light_t_orbit(double jd_tt, double e[], const struct ae_orbit_t *orb,
                      double q[]);
double ae_mjd(double jd);
double ae_mod_2pi(double x);
double ae_mod_360(double x);
double ae_mod_180(double x);
void ae_nutation_lon_ob(double jd_tt, double *nutl, double *nuto);
void ae_nutate(double nutl, double nuto, double eps, double p[], int direction);
void aes_nutate(double jd_tt, double p[], int direction);
void ae_phys_pole(double jd_ut1, double jd_tt, const struct ae_physical_t *phys,
                  double n[], double *w);
int ae_is_retrograde(const struct ae_physical_t *phys);
void ae_polar_to_rect(double ra, double dec, double radius, double rect[]);
void ae_precess(double jd_tt, double r[], int direction);
void ae_radec_to_altaz(double last, double glat, double ra, double dec,
                       double *alt, double *az);
void ae_radec_to_gal(double ra, double dec, double *l, double *b, int fk5);
void ae_gal_to_radec(double ra, double dec, double *l, double *b, int fk5);
void ae_rect_to_polar(const double rect[], double *ra, double *dec, 
                      double *radius);
void ae_rect_to_polar_with_jd(double pp[], double jd, char ofdate, 
                              double polar[]);
double ae_refrac_visible(double alt, void *param);
double ae_refrac_ulich(double alt, void *param);
void ae_relativity(double p[], double q[], double o[]);
const char *ae_retcode_string(int retcode);
double ae_tdb(double jd);
void ae_subobs_point(double j[], double n[], double w, double f, int retrograde,
                     double *lat, double *lon);
void aes_subobs_point(double jd_ut1, const struct ae_orbit_t *o_orb,
                      const struct ae_orbit_t *q_orb,
                      const struct ae_physical_t *phys,
                      double *lat, double *lon, double *dist);
double ae_flattening(const struct ae_physical_t *phys);
void ae_topocentric(double jd_tt, double jd_ut1, double tlat, double glat, 
                    double tlong, double trho, double (*refrac)(double, void *),
                    void *refrac_param, double dist, double *ra, double *dec);
void aes_topocentric(double jd_ut1, double glat, double tlong, double dist,
                     double *ra, double *dec);
void ae_v_orbit(double jd_tt, const struct ae_orbit_t *orb, double v[]);
double ae_zatan2(double x, double y);

// Useful macros.

#define ae_get_d(a) ((double)((int)(a)))

#define ae_get_m(a) ((double)fabs(((int)(((a) - ae_get_d(a)) * 60))))

#define ae_get_s(a) (fabs(((a) - ae_get_d(a) - ((double)((int)(((a) - ae_get_d(a)) * 60))) / 60) * 3600))

#define ae_dms_fmt "%dd %d' %.2f\""

#define ae_dms_sfmt "%4dd %02d' %05.2f\""

#define ae_dms_arg(a) (int)ae_get_d(a), (int)ae_get_m(a), ae_get_s(a)

#define ae_hms_fmt "%dh %dm %.2fs"

#define ae_hms_sfmt "%3dh %02dm %05.2fs"

#define ae_hms_arg(a) (int)ae_get_d(a / 15.0), (int)ae_get_m(a / 15.0), ae_get_s(a / 15.0)

// DEVEL START
struct ae_ring_geometry_t {
  int n;              
  const char **name;        
  const double *tau;        
  const double *a_inner;
  const double *a_outer;
};

struct ae_stokes_t {
  double i;     
  double q;     
  double u;     
  double v;     
};

struct ae_saturn_temp_model_t {
  double t_disc;                
  double *t_ring;               
  double t_background;          
  double *pol_ring;
  double pol_disc;              
  double pol_disc_angle;
};

struct ae_saturn_components_t {
  double n_disc;            
  double *n_ring_disc;
  double *n_ring_back;
  double *exp_tau;          
  double cos_ring_inc;      
};

struct ae_planet_map_t {
  int n_pix;        
  double res;
  struct ae_stokes_t **data;    
};

void ae_map_saturn(double res, const struct ae_physical_t *b,
                   const struct ae_ring_geometry_t *r, double inclination,
                   const struct ae_saturn_temp_model_t *sm,
                   struct ae_planet_map_t *map,
                   struct ae_saturn_components_t *comp);
int ae_map_to_gnuplot(const struct ae_planet_map_t *map, char stokes,
                      const char *path);
void ae_map_t_eff(const struct ae_planet_map_t *map,
                  const struct ae_physical_t *b, double subobs_lat,
                  double background, struct ae_stokes_t *t_eff);
void ae_saturn_t_eff(const struct ae_ring_geometry_t *r,
                     const struct ae_saturn_temp_model_t *sm,
                     const struct ae_saturn_components_t *components,
                     struct ae_stokes_t *t_eff);
void ae_set_stokes(struct ae_stokes_t *s, double i, double q, double u,
                   double v);
void ae_copy_stokes(const struct ae_stokes_t *origin, struct ae_stokes_t *dest);
void ae_saturn_components_alloc(int n_ring,
                                struct ae_saturn_components_t *comp);
void ae_saturn_components_free(struct ae_saturn_components_t *comp);
// DEVEL END

#ifdef __cplusplus
  }
#endif

#endif