/*** File libwcs/dateutil.c *** October 19, 2012 *** By Jessica Mink, jmink@cfa.harvard.edu *** Harvard-Smithsonian Center for Astrophysics *** Copyright (C) 1999-2012 *** Smithsonian Astrophysical Observatory, Cambridge, MA, USA This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version. This library 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 Lesser General Public License for more details. You should have received a copy of the GNU Lesser General Public License along with this library; if not, write to the Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA Correspondence concerning WCSTools should be addressed as follows: Internet email: jmink@cfa.harvard.edu Postal address: Jessica Mink Smithsonian Astrophysical Observatory 60 Garden St. Cambridge, MA 02138 USA */ /* Date and time conversion routines using the following conventions: ang = Angle in fractional degrees deg = Angle in degrees as dd:mm:ss.ss doy = 2 floating point numbers: year and day, including fraction, of year *** First day of year is 1, not zero. dt = 2 floating point numbers: yyyy.mmdd, hh.mmssssss ep = fractional year, often epoch of a position including proper motion epb = Besselian epoch = 365.242198781-day years based on 1900.0 epj = Julian epoch = 365.25-day years based on 2000.0 fd = FITS date string which may be any of the following: yyyy.ffff (fractional year) dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) hr = Sexigesimal hours as hh:mm:dd.ss jd = Julian Date lt = Local time mjd = modified Julian Date = JD - 2400000.5 ofd = FITS date string (dd/mm/yy before 2000, else no return) time = use fd2* with no date to convert time as hh:mm:ss.ss to sec, day, year ts = UT seconds since 1950-01-01T00:00 (used for ephemeris computations) tsi = local seconds since 1980-01-01T00:00 (used by IRAF as a time tag) tsu = UT seconds since 1970-01-01T00:00 (used as Unix system time) tsd = UT seconds of current day ut = Universal Time (UTC) et = Ephemeris Time (or TDB or TT) = TAI + 32.184 seconds tai = International Atomic Time (Temps Atomique International) = ET - 32.184 seconds gps = GPS time = TAI - 19 seconds mst = Mean Greenwich Sidereal Time gst = Greenwich Sidereal Time (includes nutation) lst = Local Sidereal Time (includes nutation) (longitude must be set) hjd = Heliocentric Julian Date mhjd = modified Heliocentric Julian Date = HJD - 2400000.5 * ang2hr (angle) * Convert angle in decimal floating point degrees to hours as hh:mm:ss.ss * ang2deg (angle) * Convert angle in decimal floating point degrees to degrees as dd:mm:ss.ss * deg2ang (angle as dd:mm:ss.ss) * Convert angle in degrees as dd:mm:ss.ss to decimal floating point degrees * ang2hr (angle) * Convert angle in hours as hh:mm:ss.ss to decimal floating point degrees * * doy2dt (year, doy, date, time) * Convert year and day of year to date as yyyy.ddmm and time as hh.mmsss * doy2ep, doy2epb, doy2epj (date, time) * Convert year and day of year to fractional year * doy2fd (year, doy) * Convert year and day of year to FITS date string * doy2mjd (year, doy) * Convert year and day of year to modified Julian date * * dt2doy (date, time, year, doy) * Convert date as yyyy.ddmm and time as hh.mmsss to year and day of year * dt2ep, dt2epb, dt2epj (date, time) * Convert date as yyyy.ddmm and time as hh.mmsss to fractional year * dt2fd (date, time) * Convert date as yyyy.ddmm and time as hh.mmsss to FITS date string * dt2i (date,time,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert yyyy.mmdd hh.mmssss to year month day hours minutes seconds * dt2jd (date,time) * Convert date as yyyy.ddmm and time as hh.mmsss to Julian date * dt2mjd (date,time) * Convert date as yyyy.ddmm and time as hh.mmsss to modified Julian date * dt2ts (date,time) * Convert date (yyyy.ddmm) and time (hh.mmsss) to seconds since 1950-01-01 * dt2tsi (date,time) * Convert date (yyyy.ddmm) and time (hh.mmsss) to seconds since 1980-01-01 * dt2tsu (date,time) * Convert date (yyyy.ddmm) and time (hh.mmsss) to seconds since 1970-01-01 * * ep2dt, epb2dt, epj2dt (epoch,date, time) * Convert fractional year to date as yyyy.ddmm and time as hh.mmsss * ep2fd, epb2fd, epj2fd (epoch) * Convert epoch to FITS ISO date string * ep2i, epb2i, epj2i (epoch,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert fractional year to year month day hours minutes seconds * ep2jd, epb2jd, epj2jd (epoch) * Convert fractional year as used in epoch to Julian date * ep2mjd, epb2mjd, epj2mjd (epoch) * Convert fractional year as used in epoch to modified Julian date * ep2ts, epb2ts, epj2ts (epoch) * Convert fractional year to seconds since 1950.0 * * et2fd (string) * Convert from ET (or TDT or TT) in FITS format to UT in FITS format * fd2et (string) * Convert from UT in FITS format to ET (or TDT or TT) in FITS format * jd2jed (dj) * Convert from Julian Date to Julian Ephemeris Date * jed2jd (dj) * Convert from Julian Ephemeris Date to Julian Date * dt2et (date, time) * Convert date (yyyy.ddmm) and time (hh.mmsss) to ephemeris time * edt2dt (date, time) * Convert ephemeris date (yyyy.ddmm) and time (hh.mmsss) to UT * dt2tai (date, time) * Convert date (yyyy.ddmm) and time (hh.mmsss) to TAI date and time * tai2dt (date, time) * Convert TAI date (yyyy.ddmm) and time (hh.mmsss) to UT * ts2ets (tsec) * Convert from UT in seconds since 1950-01-01 to ET in same format * ets2ts (tsec) * Convert from ET in seconds since 1950-01-01 to UT in same format * * fd2ep, fd2epb, fd2epj (string) * Convert FITS date string to fractional year * Convert time alone to fraction of Besselian year * fd2doy (string, year, doy) * Convert FITS standard date string to year and day of year * fd2dt (string, date, time) * Convert FITS date string to date as yyyy.ddmm and time as hh.mmsss * Convert time alone to hh.mmssss with date set to 0.0 * fd2i (string,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert FITS standard date string to year month day hours min sec * Convert time alone to hours min sec, year month day are zero * fd2jd (string) * Convert FITS standard date string to Julian date * Convert time alone to fraction of day * fd2mjd (string) * Convert FITS standard date string to modified Julian date * fd2ts (string) * Convert FITS standard date string to seconds since 1950.0 * Convert time alone to seconds of day * fd2fd (string) * Convert FITS standard date string to ISO FITS date string * fd2of (string) * Convert FITS standard date string to old-format FITS date and time * fd2ofd (string) * Convert FITS standard date string to old-format FITS date string * fd2oft (string) * Convert time part of FITS standard date string to FITS date string * * jd2doy (dj, year, doy) * Convert Julian date to year and day of year * jd2dt (dj,date,time) * Convert Julian date to date as yyyy.mmdd and time as hh.mmssss * jd2ep, jd2epb, jd2epj (dj) * Convert Julian date to fractional year as used in epoch * jd2fd (dj) * Convert Julian date to FITS ISO date string * jd2i (dj,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert Julian date to year month day hours min sec * jd2mjd (dj) * Convert Julian date to modified Julian date * jd2ts (dj) * Convert Julian day to seconds since 1950.0 * * lt2dt() * Return local time as yyyy.mmdd and time as hh.mmssss * lt2fd() * Return local time as FITS ISO date string * lt2tsi() * Return local time as IRAF seconds since 1980-01-01 00:00 * lt2tsu() * Return local time as Unix seconds since 1970-01-01 00:00 * lt2ts() * Return local time as Unix seconds since 1950-01-01 00:00 * * mjd2doy (dj,year,doy) * Convert modified Julian date to date as year and day of year * mjd2dt (dj,date,time) * Convert modified Julian date to date as yyyy.mmdd and time as hh.mmssss * mjd2ep, mjd2epb, mjd2epj (dj) * Convert modified Julian date to fractional year as used in epoch * mjd2fd (dj) * Convert modified Julian date to FITS ISO date string * mjd2i (dj,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert modified Julian date to year month day hours min sec * mjd2jd (dj) * Convert modified Julian date to Julian date * mjd2ts (dj) * Convert modified Julian day to seconds since 1950.0 * * ts2dt (tsec,date,time) * Convert seconds since 1950.0 to date as yyyy.ddmm and time as hh.mmsss * ts2ep, ts2epb, ts2epj (tsec) * Convert seconds since 1950.0 to fractional year * ts2fd (tsec) * Convert seconds since 1950.0 to FITS standard date string * ts2i (tsec,iyr,imon,iday,ihr,imn,sec, ndsec) * Convert sec since 1950.0 to year month day hours minutes seconds * ts2jd (tsec) * Convert seconds since 1950.0 to Julian date * ts2mjd (tsec) * Convert seconds since 1950.0 to modified Julian date * tsi2fd (tsec) * Convert seconds since 1980-01-01 to FITS standard date string * tsi2dt (tsec,date,time) * Convert seconds since 1980-01-01 to date as yyyy.ddmm, time as hh.mmsss * tsu2fd (tsec) * Convert seconds since 1970-01-01 to FITS standard date string * tsu2tsi (tsec) * Convert UT seconds since 1970-01-01 to local seconds since 1980-01-01 * tsu2dt (tsec,date,time) * Convert seconds since 1970-01-01 to date as yyyy.ddmm, time as hh.mmsss * * tsd2fd (tsec) * Convert seconds since start of day to FITS time, hh:mm:ss.ss * tsd2dt (tsec) * Convert seconds since start of day to hh.mmssss * * fd2gst (string) * convert from FITS date Greenwich Sidereal Time * dt2gst (date, time) * convert from UT as yyyy.mmdd hh.mmssss to Greenwich Sidereal Time * ts2gst (tsec) * Calculate Greenwich Sidereal Time given Universal Time * in seconds since 1951-01-01T0:00:00 * fd2mst (string) * convert from FITS UT date to Mean Sidereal Time * dt2gmt (date, time) * convert from UT as yyyy.mmdd hh.mmssss to Mean Sidereal Time * ts2mst (tsec) * Calculate Mean Sidereal Time given Universal Time * in seconds since 1951-01-01T0:00:00 * jd2mst (string) * convert from Julian Date to Mean Sidereal Time * mst2fd (string) * convert to current UT in FITS format given Greenwich Mean Sidereal Time * mst2jd (dj) * convert to current UT as Julian Date given Greenwich Mean Sidereal Time * jd2lst (dj) * Calculate Local Sidereal Time from Julian Date * ts2lst (tsec) * Calculate Local Sidereal Time given UT in seconds since 1951-01-01T0:00 * fd2lst (string) * Calculate Local Sidereal Time given Universal Time as FITS ISO date * lst2jd (dj, lst) * Calculate Julian Date given current Julian date and Local Sidereal Time * lst2fd (string, lst) * Calculate Julian Date given current UT date and Local Sidereal Time * gst2fd (string) * Calculate current UT given UT date and Greenwich Sidereal Time * gst2jd (dj) * Calculate current UT given UT date and Greenwich Sidereal Time as JD * * compnut (dj, dpsi, deps, eps0) * Compute the longitude and obliquity components of nutation and * mean obliquity from the IAU 1980 theory * * utdt (dj) * Compute difference between UT and dynamical time (ET-UT) * ut2dt (year, doy) * Current Universal Time to year and day of year * ut2dt (date, time) * Current Universal Time to date (yyyy.mmdd) and time (hh.mmsss) * ut2ep(), ut2epb(), ut2epj() * Current Universal Time to fractional year, Besselian, Julian epoch * ut2fd() * Current Universal Time to FITS ISO date string * ut2jd() * Current Universal Time to Julian Date * ut2mjd() * Current Universal Time to Modified Julian Date * ut2tsi() * Current Universal Time to IRAF seconds since 1980-01-01T00:00 * ut2tsu() * Current Universal Time to Unix seconds since 1970-01-01T00:00 * ut2ts() * Current Universal Time to seconds since 1950-01-01T00:00 * isdate (string) * Return 1 if string is a FITS date (old or ISO) * * Internally-used subroutines * * fixdate (iyr, imon, iday, ihr, imn, sec, ndsec) * Round seconds and make sure date and time numbers are within limits * caldays (year, month) * Calculate days in month 1-12 given year (Gregorian calendar only * dint (dnum) * Return integer part of floating point number * dmod (dnum) * Return Mod of floating point number */ #include #include #include #include #include #include #include "wcs.h" #include "fitsfile.h" static double suntl(); static void fixdate(); static int caldays(); static double dint(); static double dmod(); static double longitude = 0.0; /* longitude of observatory in degrees (+=west) */ void setlongitude (longitude0) double longitude0; { longitude = longitude0; return; } static int ndec = 3; void setdatedec (nd) int nd; { ndec = nd; return; } /* ANG2HR -- Convert angle in fraction degrees to hours as hh:mm:ss.ss */ void ang2hr (angle, lstr, string) double angle; /* Angle in fractional degrees */ int lstr; /* Maximum number of characters in string */ char *string; /* Character string (hh:mm:ss.ss returned) */ { angle = angle / 15.0; dec2str (string, lstr, angle, ndec); return; } /* ANG2DEG -- Convert angle in fraction degrees to degrees as dd:mm:ss.ss */ void ang2deg (angle, lstr, string) double angle; /* Angle in fractional degrees */ int lstr; /* Maximum number of characters in string */ char *string; /* Character string (dd:mm:ss.ss returned) */ { dec2str (string, lstr, angle, ndec); return; } /* DEG2ANG -- Convert angle in degrees as dd:mm:ss.ss to fractional degrees */ double deg2ang (angle) char *angle; /* Angle as dd:mm:ss.ss */ { double deg; deg = str2dec (angle); return (deg); } /* HR2ANG -- Convert angle in hours as hh:mm:ss.ss to fractional degrees */ double hr2ang (angle) char *angle; /* Angle in sexigesimal hours (hh:mm:ss.sss) */ { double deg; deg = str2dec (angle); deg = deg * 15.0; return (deg); } /* DT2FD-- convert vigesimal date and time to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * dt2fd (date, time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { int iyr,imon,iday,ihr,imn; double sec; int nf; char *string; char tstring[32], dstring[32]; char outform[64]; dt2i (date, time, &iyr,&imon,&iday,&ihr,&imn,&sec, ndec); /* Convert to ISO date format */ string = (char *) calloc (32, sizeof (char)); /* Make time string */ if (time != 0.0 || ndec > 0) { if (ndec == 0) nf = 2; else nf = 3 + ndec; if (ndec > 0) { sprintf (outform, "%%02d:%%02d:%%0%d.%df", nf, ndec); sprintf (tstring, outform, ihr, imn, sec); } else { sprintf (outform, "%%02d:%%02d:%%0%dd", nf); sprintf (tstring, outform, ihr, imn, (int)(sec+0.5)); } } /* Make date string */ if (date != 0.0) sprintf (dstring, "%4d-%02d-%02d", iyr, imon, iday); /* Make FITS (ISO) date string */ if (date == 0.0) strcpy (string, tstring); else if (time == 0.0 && ndec < 1) strcpy (string, dstring); else sprintf (string, "%sT%s", dstring, tstring); return (string); } /* DT2JD-- convert from date as yyyy.mmdd and time as hh.mmsss to Julian Date * Return fractional days if date is zero */ double dt2jd (date,time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Julian date (returned) */ double tsec; /* seconds since 1950.0 */ tsec = dt2ts (date, time); if (date == 0.0) dj = tsec / 86400.0; else dj = ts2jd (tsec); return (dj); } /* DT2MJD-- convert from date yyyy.mmdd time hh.mmsss to modified Julian Date * Return fractional days if date is zero */ double dt2mjd (date,time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Modified Julian date (returned) */ double tsec; /* seconds since 1950.0 */ tsec = dt2ts (date, time); if (date == 0.0) dj = tsec / 86400.0; else dj = ts2jd (tsec); return (dj - 2400000.5); } /* HJD2JD-- convert Heliocentric Julian Date to (geocentric) Julian date */ double hjd2jd (dj, ra, dec, sys) double dj; /* Heliocentric Julian date */ double ra; /* Right ascension (degrees) */ double dec; /* Declination (degrees) */ int sys; /* J2000, B1950, GALACTIC, ECLIPTIC */ { double lt; /* Light travel difference to the Sun (days) */ lt = suntl (dj, ra, dec, sys); /* Return Heliocentric Julian Date */ return (dj - lt); } /* JD2HJD-- convert (geocentric) Julian date to Heliocentric Julian Date */ double jd2hjd (dj, ra, dec, sys) double dj; /* Julian date (geocentric) */ double ra; /* Right ascension (degrees) */ double dec; /* Declination (degrees) */ int sys; /* J2000, B1950, GALACTIC, ECLIPTIC */ { double lt; /* Light travel difference to the Sun (days) */ lt = suntl (dj, ra, dec, sys); /* Return Heliocentric Julian Date */ return (dj + lt); } /* MHJD2MJD-- convert modified Heliocentric Julian Date to modified geocentric Julian date */ double mhjd2mjd (mhjd, ra, dec, sys) double mhjd; /* Modified Heliocentric Julian date */ double ra; /* Right ascension (degrees) */ double dec; /* Declination (degrees) */ int sys; /* J2000, B1950, GALACTIC, ECLIPTIC */ { double lt; /* Light travel difference to the Sun (days) */ double hjd; /* Heliocentric Julian date */ hjd = mjd2jd (mhjd); lt = suntl (hjd, ra, dec, sys); /* Return Heliocentric Julian Date */ return (jd2mjd (hjd - lt)); } /* MJD2MHJD-- convert modified geocentric Julian date tp modified Heliocentric Julian Date */ double mjd2mhjd (mjd, ra, dec, sys) double mjd; /* Julian date (geocentric) */ double ra; /* Right ascension (degrees) */ double dec; /* Declination (degrees) */ int sys; /* J2000, B1950, GALACTIC, ECLIPTIC */ { double lt; /* Light travel difference to the Sun (days) */ double dj; /* Julian date (geocentric) */ dj = mjd2jd (mjd); lt = suntl (dj, ra, dec, sys); /* Return Heliocentric Julian Date */ return (jd2mjd (dj + lt)); } /* SUNTL-- compute light travel time to heliocentric correction in days */ /* Translated into C from IRAF SPP noao.astutils.asttools.asthjd.x */ static double suntl (dj, ra, dec, sys) double dj; /* Julian date (geocentric) */ double ra; /* Right ascension (degrees) */ double dec; /* Declination (degrees) */ int sys; /* J2000, B1950, GALACTIC, ECLIPTIC */ { double t; /* Number of Julian centuries since J1900 */ double manom; /* Mean anomaly of the Earth's orbit (degrees) */ double lperi; /* Mean longitude of perihelion (degrees) */ double oblq; /* Mean obliquity of the ecliptic (degrees) */ double eccen; /* Eccentricity of the Earth's orbit (dimensionless) */ double eccen2, eccen3; double tanom; /* True anomaly (approximate formula) (radians) */ double slong; /* True longitude of the Sun from the Earth (radians) */ double rs; /* Distance to the sun (AU) */ double lt; /* Light travel difference to the Sun (days) */ double l; /* Longitude of star in orbital plane of Earth (radians) */ double b; /* Latitude of star in orbital plane of Earth (radians) */ double epoch; /* Epoch of obervation */ double rs1,rs2; t = (dj - 2415020.0) / 36525.0; /* Compute earth orbital parameters */ manom = 358.47583 + (t * (35999.04975 - t * (0.000150 + t * 0.000003))); lperi = 101.22083 + (t * (1.7191733 + t * (0.000453 + t * 0.000003))); oblq = 23.452294 - (t * (0.0130125 + t * (0.00000164 - t * 0.000000503))); eccen = 0.01675104 - (t * (0.00004180 + t * 0.000000126)); eccen2 = eccen * eccen; eccen3 = eccen * eccen2; /* Convert to principle angles */ manom = manom - (360.0 * (dint) (manom / 360.0)); lperi = lperi - (360.0 * (dint) (lperi / 360.0)); /* Convert to radians */ manom = degrad (manom); lperi = degrad (lperi); oblq = degrad (oblq); /* True anomaly */ tanom = manom + (2 * eccen - 0.25 * eccen3) * sin (manom) + 1.25 * eccen2 * sin (2 * manom) + 13./12. * eccen3 * sin (3 * manom); /* Distance to the Sun */ rs1 = 1.0 - eccen2; rs2 = 1.0 + (eccen * cos (tanom)); rs = rs1 / rs2; /* True longitude of the Sun seen from the Earth */ slong = lperi + tanom + PI; /* Longitude and latitude of star in orbital plane of the Earth */ epoch = jd2ep (dj); wcscon (sys, WCS_ECLIPTIC, 0.0, 0.0, &ra, &dec, epoch); l = degrad (ra); b = degrad (dec); /* Light travel difference to the Sun */ lt = -0.005770 * rs * cos (b) * cos (l - slong); /* Return light travel difference */ return (lt); } /* JD2DT-- convert Julian date to date as yyyy.mmdd and time as hh.mmssss */ void jd2dt (dj,date,time) double dj; /* Julian date */ double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { int iyr,imon,iday,ihr,imn; double sec; /* Convert Julian Date to date and time */ jd2i (dj, &iyr, &imon, &iday, &ihr, &imn, &sec, 4); /* Convert date to yyyy.mmdd */ if (iyr < 0) { *date = (double) (-iyr) + 0.01 * (double) imon + 0.0001 * (double) iday; *date = -(*date); } else *date = (double) iyr + 0.01 * (double) imon + 0.0001 * (double) iday; /* Convert time to hh.mmssssss */ *time = (double) ihr + 0.01 * (double) imn + 0.0001 * sec; return; } /* JD2I-- convert Julian date to date as year, month, and day, and time hours, minutes, and seconds */ /* after Fliegel and Van Flander, CACM 11, 657 (1968) */ void jd2i (dj, iyr, imon, iday, ihr, imn, sec, ndsec) double dj; /* Julian date */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double tsec; double frac, dts, ts, sday; int jd, l, n, i, j; tsec = jd2ts (dj); /* ts2i (tsec, iyr, imon, iday, ihr, imn, sec, ndsec); */ /* Round seconds to 0 - 4 decimal places */ if (tsec < 0.0) dts = -0.5; else dts = 0.5; if (ndsec < 1) ts = dint (tsec + dts); else if (ndsec < 2) ts = dint (tsec * 10.0 + dts) / 10.0; else if (ndsec < 3) ts = dint (tsec * 100.0 + dts) / 100.0; else if (ndsec < 4) ts = dint (tsec * 1000.0 + dts) / 1000.0; else ts = dint (tsec * 10000.0 + dts) / 10000.0; /* Convert back to Julian Date */ dj = ts2jd (ts); /* Compute time from fraction of a day */ frac = dmod (dj, 1.0); if (frac < 0.5) { jd = (int) (dj - frac); sday = (frac + 0.5) * 86400.0; } else { jd = (int) (dj - frac) + 1; sday = (frac - 0.5) * 86400.0; } *ihr = (int) (sday / 3600.0); sday = sday - (double) (*ihr * 3600); *imn = (int) (sday / 60.0); *sec = sday - (double) (*imn * 60); /* Compute day, month, year */ l = jd + 68569; n = (4 * l) / 146097; l = l - (146097 * n + 3) / 4; i = (4000 * (l + 1)) / 1461001; l = l - (1461 * i) / 4 + 31; j = (80 * l) / 2447; *iday = l - (2447 * j) / 80; l = j / 11; *imon = j + 2 - (12 * l); *iyr = 100 * (n - 49) + i + l; return; } /* JD2MJD-- convert Julian Date to Modified Julian Date */ double jd2mjd (dj) double dj; /* Julian Date */ { return (dj - 2400000.5); } /* JD2EP-- convert Julian date to fractional year as used in epoch */ double jd2ep (dj) double dj; /* Julian date */ { double date, time; jd2dt (dj, &date, &time); return (dt2ep (date, time)); } /* JD2EPB-- convert Julian date to Besselian epoch */ double jd2epb (dj) double dj; /* Julian date */ { return (1900.0 + (dj - 2415020.31352) / 365.242198781); } /* JD2EPJ-- convert Julian date to Julian epoch */ double jd2epj (dj) double dj; /* Julian date */ { return (2000.0 + (dj - 2451545.0) / 365.25); } /* LT2DT-- Return local time as yyyy.mmdd and time as hh.mmssss */ void lt2dt(date, time) double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { time_t tsec; struct timeval tp; struct timezone tzp; struct tm *ts; gettimeofday (&tp,&tzp); tsec = tp.tv_sec; ts = localtime (&tsec); if (ts->tm_year < 1000) *date = (double) (ts->tm_year + 1900); else *date = (double) ts->tm_year; *date = *date + (0.01 * (double) (ts->tm_mon + 1)); *date = *date + (0.0001 * (double) ts->tm_mday); *time = (double) ts->tm_hour; *time = *time + (0.01 * (double) ts->tm_min); *time = *time + (0.0001 * (double) ts->tm_sec); return; } /* LT2FD-- Return current local time as FITS ISO date string */ char * lt2fd() { time_t tsec; struct tm *ts; struct timeval tp; struct timezone tzp; int month, day, year, hour, minute, second; char *isotime; gettimeofday (&tp,&tzp); tsec = tp.tv_sec; ts = localtime (&tsec); year = ts->tm_year; if (year < 1000) year = year + 1900; month = ts->tm_mon + 1; day = ts->tm_mday; hour = ts->tm_hour; minute = ts->tm_min; second = ts->tm_sec; isotime = (char *) calloc (32, sizeof (char)); sprintf (isotime, "%04d-%02d-%02dT%02d:%02d:%02d", year, month, day, hour, minute, second); return (isotime); } /* LT2TSI-- Return local time as IRAF seconds since 1980-01-01 00:00 */ int lt2tsi() { return ((int)(lt2ts() - 946684800.0)); } /* LT2TSU-- Return local time as Unix seconds since 1970-01-01 00:00 */ time_t lt2tsu() { return ((time_t)(lt2ts() - 631152000.0)); } /* LT2TS-- Return local time as Unix seconds since 1950-01-01 00:00 */ double lt2ts() { double tsec; char *datestring; datestring = lt2fd(); tsec = fd2ts (datestring); free (datestring); return (tsec); } /* MJD2DT-- convert Modified Julian Date to date (yyyy.mmdd) time (hh.mmssss) */ void mjd2dt (dj,date,time) double dj; /* Modified Julian Date */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double tsec; tsec = jd2ts (dj + 2400000.5); ts2dt (tsec, date, time); return; } /* MJD2I-- convert Modified Julian Date to date as year, month, day and time as hours, minutes, seconds */ void mjd2i (dj, iyr, imon, iday, ihr, imn, sec, ndsec) double dj; /* Modified Julian Date */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double tsec; tsec = jd2ts (dj + 2400000.5); ts2i (tsec, iyr, imon, iday, ihr, imn, sec, ndsec); return; } /* MJD2DOY-- convert Modified Julian Date to Year,Day-of-Year */ void mjd2doy (dj, year, doy) double dj; /* Modified Julian Date */ int *year; /* Year (returned) */ double *doy; /* Day of year with fraction (returned) */ { jd2doy (dj + 2400000.5, year, doy); return; } /* MJD2JD-- convert Modified Julian Date to Julian Date */ double mjd2jd (dj) double dj; /* Modified Julian Date */ { return (dj + 2400000.5); } /* MJD2EP-- convert Modified Julian Date to fractional year */ double mjd2ep (dj) double dj; /* Modified Julian Date */ { double date, time; jd2dt (dj + 2400000.5, &date, &time); return (dt2ep (date, time)); } /* MJD2EPB-- convert Modified Julian Date to Besselian epoch */ double mjd2epb (dj) double dj; /* Modified Julian Date */ { return (1900.0 + (dj - 15019.81352) / 365.242198781); } /* MJD2EPJ-- convert Modified Julian Date to Julian epoch */ double mjd2epj (dj) double dj; /* Modified Julian Date */ { return (2000.0 + (dj - 51544.5) / 365.25); } /* MJD2FD-- convert modified Julian date to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * mjd2fd (dj) double dj; /* Modified Julian date */ { return (jd2fd (dj + 2400000.5)); } /* MJD2TS-- convert modified Julian date to seconds since 1950.0 */ double mjd2ts (dj) double dj; /* Modified Julian date */ { return ((dj - 33282.0) * 86400.0); } /* EP2FD-- convert fractional year to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * ep2fd (epoch) double epoch; /* Date as fractional year */ { double tsec; /* seconds since 1950.0 (returned) */ tsec = ep2ts (epoch); return (ts2fd (tsec)); } /* EPB2FD-- convert Besselian epoch to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * epb2fd (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { double dj; /* Julian Date */ dj = epb2jd (epoch); return (jd2fd (dj)); } /* EPJ2FD-- convert Julian epoch to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * epj2fd (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { double dj; /* Julian Date */ dj = epj2jd (epoch); return (jd2fd (dj)); } /* EP2TS-- convert fractional year to seconds since 1950.0 */ double ep2ts (epoch) double epoch; /* Date as fractional year */ { double dj; dj = ep2jd (epoch); return ((dj - 2433282.5) * 86400.0); } /* EPB2TS-- convert Besselian epoch to seconds since 1950.0 */ double epb2ts (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { double dj; dj = epb2jd (epoch); return ((dj - 2433282.5) * 86400.0); } /* EPJ2TS-- convert Julian epoch to seconds since 1950.0 */ double epj2ts (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { double dj; dj = epj2jd (epoch); return ((dj - 2433282.5) * 86400.0); } /* EPB2EP-- convert Besselian epoch to fractional years */ double epb2ep (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { double dj; dj = epb2jd (epoch); return (jd2ep (dj)); } /* EP2EPB-- convert fractional year to Besselian epoch */ double ep2epb (epoch) double epoch; /* Fractional year */ { double dj; dj = ep2jd (epoch); return (jd2epb (dj)); } /* EPJ2EP-- convert Julian epoch to fractional year */ double epj2ep (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { double dj; dj = epj2jd (epoch); return (jd2ep (dj)); } /* EP2EPJ-- convert fractional year to Julian epoch */ double ep2epj (epoch) double epoch; /* Fractional year */ { double dj; dj = ep2jd (epoch); return (jd2epj (dj)); } /* EP2I-- convert fractional year to year month day hours min sec */ void ep2i (epoch, iyr, imon, iday, ihr, imn, sec, ndsec) double epoch; /* Date as fractional year */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double date, time; ep2dt (epoch, &date, &time); dt2i (date, time, iyr,imon,iday,ihr,imn,sec, ndsec); return; } /* EPB2I-- convert Besselian epoch to year month day hours min sec */ void epb2i (epoch, iyr, imon, iday, ihr, imn, sec, ndsec) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double date, time; epb2dt (epoch, &date, &time); dt2i (date, time, iyr,imon,iday,ihr,imn,sec, ndsec); return; } /* EPJ2I-- convert Julian epoch to year month day hours min sec */ void epj2i (epoch, iyr, imon, iday, ihr, imn, sec, ndsec) double epoch; /* Julian epoch (fractional 365.25-day years) */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double date, time; epj2dt (epoch, &date, &time); dt2i (date, time, iyr,imon,iday,ihr,imn,sec, ndsec); return; } /* EP2JD-- convert fractional year as used in epoch to Julian date */ double ep2jd (epoch) double epoch; /* Date as fractional year */ { double dj; /* Julian date (returned)*/ double date, time; ep2dt (epoch, &date, &time); dj = dt2jd (date, time); return (dj); } /* EPB2JD-- convert Besselian epoch to Julian Date */ double epb2jd (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { return (2415020.31352 + ((epoch - 1900.0) * 365.242198781)); } /* EPJ2JD-- convert Julian epoch to Julian Date */ double epj2jd (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { return (2451545.0 + ((epoch - 2000.0) * 365.25)); } /* EP2MJD-- convert fractional year as used in epoch to modified Julian date */ double ep2mjd (epoch) double epoch; /* Date as fractional year */ { double dj; /* Julian date (returned)*/ double date, time; ep2dt (epoch, &date, &time); dj = dt2jd (date, time); return (dj - 2400000.5); } /* EPB2MJD-- convert Besselian epoch to modified Julian Date */ double epb2mjd (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { return (15019.81352 + ((epoch - 1900.0) * 365.242198781)); } /* EPJ2MJD-- convert Julian epoch to modified Julian Date */ double epj2mjd (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { return (51544.5 + ((epoch - 2000.0) * 365.25)); } /* EPB2EPJ-- convert Besselian epoch to Julian epoch */ double epb2epj (epoch) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ { double dj; /* Julian date */ dj = epb2jd (epoch); return (jd2epj (dj)); } /* EPJ2EPB-- convert Julian epoch to Besselian epoch */ double epj2epb (epoch) double epoch; /* Julian epoch (fractional 365.25-day years) */ { double dj; /* Julian date */ dj = epj2jd (epoch); return (jd2epb (dj)); } /* JD2FD-- convert Julian date to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * jd2fd (dj) double dj; /* Julian date */ { double tsec; /* seconds since 1950.0 (returned) */ tsec = (dj - 2433282.5) * 86400.0; return (ts2fd (tsec)); } /* JD2TS-- convert Julian date to seconds since 1950.0 */ double jd2ts (dj) double dj; /* Julian date */ { return ((dj - 2433282.5) * 86400.0); } /* JD2TSI-- convert Julian date to IRAF seconds since 1980-01-01T0:00 */ int jd2tsi (dj) double dj; /* Julian date */ { double ts; ts = (dj - 2444239.5) * 86400.0; return ((int) ts); } /* JD2TSU-- convert Julian date to Unix seconds since 1970-01-01T0:00 */ time_t jd2tsu (dj) double dj; /* Julian date */ { return ((time_t)((dj - 2440587.5) * 86400.0)); } /* DT2DOY-- convert yyyy.mmdd hh.mmss to year and day of year */ void dt2doy (date, time, year, doy) double date; /* Date as yyyy.mmdd */ double time; /* Time as hh.mmssxxxx */ int *year; /* Year (returned) */ double *doy; /* Day of year with fraction (returned) */ { double dj; /* Julian date */ double dj0; /* Julian date on January 1 0:00 */ double date0; /* January first of date's year */ double dyear; dyear = floor (date); date0 = dyear + 0.0101; dj0 = dt2jd (date0, 0.0); dj = dt2jd (date, time); *year = (int) (dyear + 0.00000001); *doy = dj - dj0 + 1.0; return; } /* DOY2DT-- convert year and day of year to yyyy.mmdd hh.mmss */ void doy2dt (year, doy, date, time) int year; /* Year */ double doy; /* Day of year with fraction */ double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { double dj; /* Julian date */ double dj0; /* Julian date on January 1 0:00 */ double date0; /* January first of date's year */ date0 = year + 0.0101; dj0 = dt2jd (date0, 0.0); dj = dj0 + doy - 1.0; jd2dt (dj, date, time); return; } /* DOY2EP-- convert year and day of year to fractional year as used in epoch */ double doy2ep (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double date, time; doy2dt (year, doy, &date, &time); return (dt2ep (date, time)); } /* DOY2EPB-- convert year and day of year to Besellian epoch */ double doy2epb (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; dj = doy2jd (year, doy); return (jd2epb (dj)); } /* DOY2EPJ-- convert year and day of year to Julian epoch */ double doy2epj (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; dj = doy2jd (year, doy); return (jd2epj (dj)); } /* DOY2FD-- convert year and day of year to FITS date */ char * doy2fd (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; /* Julian date */ dj = doy2jd (year, doy); return (jd2fd (dj)); } /* DOY2JD-- convert year and day of year to Julian date */ double doy2jd (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj0; /* Julian date */ double date; /* Date as yyyy.mmdd (returned) */ double time; /* Time as hh.mmssxxxx (returned) */ date = (double) year + 0.0101; time = 0.0; dj0 = dt2jd (date, time); return (dj0 + doy - 1.0); } /* DOY2MJD-- convert year and day of year to Julian date */ double doy2mjd (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj0; /* Julian date */ double date; /* Date as yyyy.mmdd (returned) */ double time; /* Time as hh.mmssxxxx (returned) */ date = (double) year + 0.0101; time = 0.0; dj0 = dt2jd (date, time); return (dj0 + doy - 1.0 - 2400000.5); } /* DOY2TSU-- convert from FITS date to Unix seconds since 1970-01-01T0:00 */ time_t doy2tsu (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; dj = doy2jd (year, doy); return ((time_t)jd2ts (dj)); } /* DOY2TSI-- convert from FITS date to IRAF seconds since 1980-01-01T0:00 */ int doy2tsi (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; dj = doy2jd (year, doy); return ((int)jd2tsi (dj)); } /* DOY2TS-- convert year, day of year to seconds since 1950 */ double doy2ts (year, doy) int year; /* Year */ double doy; /* Day of year with fraction */ { double dj; dj = doy2jd (year, doy); return (jd2ts (dj)); } /* FD2DOY-- convert FITS date to year and day of year */ void fd2doy (string, year, doy) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ int *year; /* Year (returned) */ double *doy; /* Day of year with fraction (returned) */ { double dj; /* Julian date */ dj = fd2jd (string); jd2doy (dj, year, doy); return; } /* JD2DOY-- convert Julian date to year and day of year */ void jd2doy (dj, year, doy) double dj; /* Julian date */ int *year; /* Year (returned) */ double *doy; /* Day of year with fraction (returned) */ { double date; /* Date as yyyy.mmdd (returned) */ double time; /* Time as hh.mmssxxxx (returned) */ double dj0; /* Julian date at 0:00 on 1/1 */ double dyear; jd2dt (dj, &date, &time); *year = (int) date; dyear = (double) *year; dj0 = dt2jd (dyear+0.0101, 0.0); *doy = dj - dj0 + 1.0; return; } /* TS2JD-- convert seconds since 1950.0 to Julian date */ double ts2jd (tsec) double tsec; /* seconds since 1950.0 */ { return (2433282.5 + (tsec / 86400.0)); } /* TS2MJD-- convert seconds since 1950.0 to modified Julian date */ double ts2mjd (tsec) double tsec; /* seconds since 1950.0 */ { return (33282.0 + (tsec / 86400.0)); } /* TS2EP-- convert seconds since 1950.0 to fractional year as used in epoch */ double ts2ep (tsec) double tsec; /* Seconds since 1950.0 */ { double date, time; ts2dt (tsec, &date, &time); return (dt2ep (date, time)); } /* TS2EPB-- convert seconds since 1950.0 to Besselian epoch */ double ts2epb (tsec) double tsec; /* Seconds since 1950.0 */ { double dj; /* Julian Date */ dj = ts2jd (tsec); return (jd2epb (dj)); } /* TS2EPB-- convert seconds since 1950.0 to Julian epoch */ double ts2epj (tsec) double tsec; /* Seconds since 1950.0 */ { double dj; /* Julian Date */ dj = ts2jd (tsec); return (jd2epj (dj)); } /* DT2EP-- convert from date, time as yyyy.mmdd hh.mmsss to fractional year */ double dt2ep (date, time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double epoch; /* Date as fractional year (returned) */ double dj, dj0, dj1, date0, time0, date1; dj = dt2jd (date, time); if (date == 0.0) epoch = dj / 365.2422; else { time0 = 0.0; date0 = dint (date) + 0.0101; date1 = dint (date) + 1.0101; dj0 = dt2jd (date0, time0); dj1 = dt2jd (date1, time0); epoch = dint (date) + ((dj - dj0) / (dj1 - dj0)); } return (epoch); } /* DT2EPB-- convert from date, time as yyyy.mmdd hh.mmsss to Besselian epoch */ double dt2epb (date, time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Julian date */ double epoch; /* Date as fractional year (returned) */ dj = dt2jd (date, time); if (date == 0.0) epoch = dj / 365.242198781; else epoch = jd2epb (dj); return (epoch); } /* DT2EPJ-- convert from date, time as yyyy.mmdd hh.mmsss to Julian epoch */ double dt2epj (date, time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Julian date */ double epoch; /* Date as fractional year (returned) */ dj = dt2jd (date, time); if (date == 0.0) epoch = dj / 365.25; else epoch = jd2epj (dj); return (epoch); } /* EP2DT-- convert from fractional year to date, time as yyyy.mmdd hh.mmsss */ void ep2dt (epoch, date, time) double epoch; /* Date as fractional year */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj, dj0, dj1, date0, time0, date1, epochi, epochf; time0 = 0.0; epochi = dint (epoch); epochf = epoch - epochi; date0 = epochi + 0.0101; date1 = epochi + 1.0101; dj0 = dt2jd (date0, time0); dj1 = dt2jd (date1, time0); dj = dj0 + epochf * (dj1 - dj0); jd2dt (dj, date, time); return; } /* EPB2DT-- convert from Besselian epoch to date, time as yyyy.mmdd hh.mmsss */ void epb2dt (epoch, date, time) double epoch; /* Besselian epoch (fractional 365.242198781-day years) */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Julian date */ dj = epb2jd (epoch); jd2dt (dj, date, time); } /* EPJ2DT-- convert from Julian epoch to date, time as yyyy.mmdd hh.mmsss */ void epj2dt (epoch, date, time) double epoch; /* Julian epoch (fractional 365.25-day years) */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double dj; /* Julian date */ dj = epj2jd (epoch); jd2dt (dj, date, time); } /* FD2JD-- convert FITS standard date to Julian date */ double fd2jd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double date, time; fd2dt (string, &date, &time); return (dt2jd (date, time)); } /* FD2MJD-- convert FITS standard date to modified Julian date */ double fd2mjd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { return (fd2jd (string) - 2400000.5); } /* FD2TSU-- convert from FITS date to Unix seconds since 1970-01-01T0:00 */ time_t fd2tsu (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double date, time; fd2dt (string, &date, &time); return (dt2tsu (date, time)); } /* FD2TSI-- convert from FITS date to IRAF seconds since 1980-01-01T0:00 */ int fd2tsi (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double date, time; fd2dt (string, &date, &time); return (dt2tsi (date, time)); } /* FD2TS-- convert FITS standard date to seconds since 1950 */ double fd2ts (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double date, time; fd2dt (string, &date, &time); return (dt2ts (date, time)); } /* FD2FD-- convert any FITS standard date to ISO FITS standard date */ char * fd2fd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double date, time; fd2dt (string, &date, &time); return (dt2fd (date, time)); } /* FD2OF-- convert any FITS standard date to old FITS standard date time */ char * fd2of (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { int iyr,imon,iday,ihr,imn; double sec; fd2i (string,&iyr,&imon,&iday,&ihr,&imn,&sec, 3); /* Convert to old FITS date format */ string = (char *) calloc (32, sizeof (char)); if (iyr < 1900) sprintf (string, "*** date out of range ***"); else if (iyr < 2000) sprintf (string, "%02d/%02d/%02d %02d:%02d:%06.3f", iday, imon, iyr-1900, ihr, imn, sec); else if (iyr < 2900.0) sprintf (string, "%02d/%02d/%3d %02d:%02d:%6.3f", iday, imon, iyr-1900, ihr, imn, sec); else sprintf (string, "*** date out of range ***"); return (string); } /* TAI-UTC from the U.S. Naval Observatory */ /* ftp://maia.usno.navy.mil/ser7/tai-utc.dat */ static double taijd[26]={2441317.5, 2441499.5, 2441683.5, 2442048.5, 2442413.5, 2442778.5, 2443144.5, 2443509.5, 2443874.5, 2444239.5, 2444786.5, 2445151.5, 2445516.5, 2446247.5, 2447161.5, 2447892.5, 2448257.5, 2448804.5, 2449169.5, 2449534.5, 2450083.5, 2450630.5, 2451179.5, 2453736.5, 2454832.5, 2456293.5}; static double taidt[26]={10.0,11.0,12.0,13.0,14.0,15.0,16.0,17.0,18.0,19.0, 20.0,21.0,22.0,23.0,24.0,25.0,26.0,27.0,28.0,29.0,30.0,31.0,32.0, 33.0,34.0,35.0}; static double dttab[173]={13.7,13.4,13.1,12.9,12.7,12.6,12.5,12.5,12.5,12.5, 12.5,12.5,12.5,12.5,12.5,12.5,12.5,12.4,12.3,12.2,12.0,11.7,11.4, 11.1,10.6,10.2, 9.6, 9.1, 8.6, 8.0, 7.5, 7.0, 6.6, 6.3, 6.0, 5.8, 5.7, 5.6, 5.6, 5.6, 5.7, 5.8, 5.9, 6.1, 6.2, 6.3, 6.5, 6.6, 6.8, 6.9, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.7, 7.8, 7.8,7.88,7.82, 7.54, 6.97, 6.40, 6.02, 5.41, 4.10, 2.92, 1.82, 1.61, 0.10,-1.02, -1.28,-2.69,-3.24,-3.64,-4.54,-4.71,-5.11,-5.40,-5.42,-5.20,-5.46, -5.46,-5.79,-5.63,-5.64,-5.80,-5.66,-5.87,-6.01,-6.19,-6.64,-6.44, -6.47,-6.09,-5.76,-4.66,-3.74,-2.72,-1.54,-0.02, 1.24, 2.64, 3.86, 5.37, 6.14, 7.75, 9.13,10.46,11.53,13.36,14.65,16.01,17.20,18.24, 19.06,20.25,20.95,21.16,22.25,22.41,23.03,23.49,23.62,23.86,24.49, 24.34,24.08,24.02,24.00,23.87,23.95,23.86,23.93,23.73,23.92,23.96, 24.02,24.33,24.83,25.30,25.70,26.24,26.77,27.28,27.78,28.25,28.71, 29.15,29.57,29.97,30.36,30.72,31.07,31.35,31.68,32.18,32.68,33.15, 33.59,34.00,34.47,35.03,35.73,36.54,37.43,38.29,39.20,40.18,41.17, 42.23}; /* TAI2FD-- convert from TAI in FITS format to UT in FITS format */ char * tai2fd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj0, dj, tsec, dt; dj0 = fd2jd (string); dt = utdt (dj0); dj = dj0 - (dt / 86400.0); dt = utdt (dj); tsec = fd2ts (string); tsec = tsec - dt + 32.184; return (ts2fd (tsec)); } /* FD2TAI-- convert from UT in FITS format to TAI in FITS format */ char * fd2tai (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj, tsec, dt; dj = fd2jd (string); dt = utdt (dj); tsec = fd2ts (string); tsec = tsec + dt - 32.184; return (ts2fd (tsec)); } /* DT2TAI-- convert from UT as yyyy.mmdd hh.mmssss to TAI in same format */ void dt2tai (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, dt, tsec; dj = dt2jd (*date, *time); dt = utdt (dj); tsec = dt2ts (*date, *time); tsec = tsec + dt - 32.184; ts2dt (tsec, date, time); return; } /* TAI2DT-- convert from TAI as yyyy.mmdd hh.mmssss to UT in same format */ void tai2dt (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, dt, tsec, tsec0; dj = dt2jd (*date, *time); dt = utdt (dj); tsec0 = dt2ts (*date, *time); tsec = tsec0 + dt; dj = ts2jd (tsec); dt = utdt (dj); tsec = tsec0 + dt + 32.184; ts2dt (tsec, date, time); return; } /* ET2FD-- convert from ET (or TDT or TT) in FITS format to UT in FITS format */ char * et2fd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj0, dj, tsec, dt; dj0 = fd2jd (string); dt = utdt (dj0); dj = dj0 - (dt / 86400.0); dt = utdt (dj); tsec = fd2ts (string); tsec = tsec - dt; return (ts2fd (tsec)); } /* FD2ET-- convert from UT in FITS format to ET (or TDT or TT) in FITS format */ char * fd2et (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj, tsec, dt; dj = fd2jd (string); dt = utdt (dj); tsec = fd2ts (string); tsec = tsec + dt; return (ts2fd (tsec)); } /* DT2ET-- convert from UT as yyyy.mmdd hh.mmssss to ET in same format */ void dt2et (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, dt, tsec; dj = dt2jd (*date, *time); dt = utdt (dj); tsec = dt2ts (*date, *time); tsec = tsec + dt; ts2dt (tsec, date, time); return; } /* EDT2DT-- convert from ET as yyyy.mmdd hh.mmssss to UT in same format */ void edt2dt (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, dt, tsec, tsec0; dj = dt2jd (*date, *time); dt = utdt (dj); tsec0 = dt2ts (*date, *time); tsec = tsec0 + dt; dj = ts2jd (tsec); dt = utdt (dj); tsec = tsec0 + dt; ts2dt (tsec, date, time); return; } /* JD2JED-- convert from Julian Date to Julian Ephemeris Date */ double jd2jed (dj) double dj; /* Julian Date */ { double dt; dt = utdt (dj); return (dj + (dt / 86400.0)); } /* JED2JD-- convert from Julian Ephemeris Date to Julian Date */ double jed2jd (dj) double dj; /* Julian Ephemeris Date */ { double dj0, dt; dj0 = dj; dt = utdt (dj); dj = dj0 - (dt / 86400.0); dt = utdt (dj); return (dj - (dt / 86400.0)); } /* TS2ETS-- convert from UT in seconds since 1950-01-01 to ET in same format */ double ts2ets (tsec) double tsec; { double dj, dt; dj = ts2jd (tsec); dt = utdt (dj); return (tsec + dt); } /* ETS2TS-- convert from ET in seconds since 1950-01-01 to UT in same format */ double ets2ts (tsec) double tsec; { double dj, dj0, dt; dj0 = ts2jd (tsec); dt = utdt (dj0); dj = dj0 - (dt / 86400.0); dt = utdt (dj); return (tsec - dt); } /* UTDT-- Compute difference between UT and dynamical time (ET-UT) */ double utdt (dj) double dj; /* Julian Date (UT) */ { double dt, date, time, ts, ts1, ts0, date0, yfrac, diff, cj; int i, iyr, iyear; /* If after 1972-01-01, use tabulated TAI-UT */ if (dj >= 2441317.5) { dt = 0.0; for (i = 22; i > 0; i--) { if (dj >= taijd[i]) dt = taidt[i]; } dt = dt + 32.184; } /* For 1800-01-01 to 1972-01-01, use table of ET-UT from AE */ else if (dj >= 2378496.5) { jd2dt (dj, &date, &time); ts = jd2ts (dj); iyear = (int) date; iyr = iyear - 1800; date0 = (double) iyear + 0.0101; ts0 = dt2ts (date0, 0.0); date0 = (double) (iyear + 1) + 0.0101; ts1 = dt2ts (date0, 0.0); yfrac = (ts - ts0) / (ts1 - ts0); diff = dttab[iyr+1] - dttab[iyr]; dt = dttab[iyr] + (diff * yfrac); } /* Compute back to 1600 using formula from McCarthy and Babcock (1986) */ else if (dj >= 2305447.5) { cj = (dj - 2378496.5) / 36525.0; dt = 5.156 + 13.3066 * (cj - 0.19) * (cj - 0.19); } /* Compute back to 948 using formula from Stephenson and Morrison (1984) */ else if (dj >= 2067309.5) { cj = (dj - 2378496.5) / 36525.0; dt = 25.5 * cj * cj; } /*Compute back to 390 BC using formula from Stephenson and Morrison (1984)*/ else if (dj >= 0.0) { cj = (dj = 2378496.5) / 36525.0; dt = 1360.0 + (320.0 * cj) + (44.3 * cj * cj); } else dt = 0.0; return (dt); } /* FD2OFD-- convert any FITS standard date to old FITS standard date */ char * fd2ofd (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { int iyr,imon,iday,ihr,imn; double sec; fd2i (string,&iyr,&imon,&iday,&ihr,&imn,&sec, 3); /* Convert to old FITS date format */ string = (char *) calloc (32, sizeof (char)); if (iyr < 1900) sprintf (string, "*** date out of range ***"); else if (iyr < 2000) sprintf (string, "%02d/%02d/%02d", iday, imon, iyr-1900); else if (iyr < 2900.0) sprintf (string, "%02d/%02d/%3d", iday, imon, iyr-1900); else sprintf (string, "*** date out of range ***"); return (string); } /* FD2OFT-- convert any FITS standard date to old FITS standard time */ char * fd2oft (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { int iyr,imon,iday,ihr,imn; double sec; fd2i (string,&iyr,&imon,&iday,&ihr,&imn,&sec, 3); /* Convert to old FITS date format */ string = (char *) calloc (32, sizeof (char)); sprintf (string, "%02d:%02d:%06.3f", ihr, imn, sec); return (string); } /* FD2DT-- convert FITS standard date to date, time as yyyy.mmdd hh.mmsss */ void fd2dt (string, date, time) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { int iyr,imon,iday,ihr,imn; double sec; fd2i (string,&iyr,&imon,&iday,&ihr,&imn,&sec, 4); /* Convert date to yyyy.mmdd */ if (iyr < 0) { *date = (double) (-iyr) + 0.01 * (double) imon + 0.0001 * (double) iday; *date = -(*date); } else *date = (double) iyr + 0.01 * (double) imon + 0.0001 * (double) iday; /* Convert time to hh.mmssssss */ *time = (double) ihr + 0.01 * (double) imn + 0.0001 * sec; return; } /* FD2EP-- convert from FITS standard date to fractional year */ double fd2ep (string) char *string; /* FITS date string, which may be: yyyy.ffff (fractional year) dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj; /* Julian date */ dj = fd2jd (string); if (dj < 1.0) return (dj / 365.2422); else return (jd2ep (dj)); } /* FD2EPB-- convert from FITS standard date to Besselian epoch */ double fd2epb (string) char *string; /* FITS date string, which may be: yyyy.ffff (fractional year) dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj; /* Julian date */ dj = fd2jd (string); if (dj < 1.0) return (dj / 365.242198781); else return (jd2epb (dj)); } /* FD2EPJ-- convert from FITS standard date to Julian epoch */ double fd2epj (string) char *string; /* FITS date string, which may be: yyyy.ffff (fractional year) dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj; /* Julian date */ dj = fd2jd (string); if (dj < 1.0) return (dj / 365.25); else return (jd2epj (dj)); } /* DT2TSU-- convert from date and time to Unix seconds since 1970-01-01T0:00 */ time_t dt2tsu (date,time) double date; /* Date as yyyy.mmdd */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { return ((time_t)(dt2ts (date, time) - 631152000.0)); } /* DT2TSI-- convert from date and time to IRAF seconds since 1980-01-01T0:00 */ int dt2tsi (date,time) double date; /* Date as yyyy.mmdd */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { return ((int)(dt2ts (date, time) - 946684800.0)); } /* DT2TS-- convert from date, time as yyyy.mmdd hh.mmsss to sec since 1950.0 */ double dt2ts (date,time) double date; /* Date as yyyy.mmdd yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { double tsec; /* Seconds past 1950.0 (returned) */ double dh,dm,dd; int iy,im,id; /* Calculate the number of full years, months, and days already * elapsed since 0h, March 1, -1 (up to most recent midnight). */ /* convert time of day to elapsed seconds */ /* If time is < 0, it is assumed to be a fractional day */ if (time < 0.0) tsec = time * -86400.0; else { dh = (int) (time + 0.0000000001); dm = (int) (((time - dh) * 100.0) + 0.0000000001); tsec = (time * 10000.0) - (dh * 10000.0) - (dm * 100.0); tsec = (int) (tsec * 100000.0 + 0.0001) / 100000.0; tsec = tsec + (dm * 60.0) + (dh * 3600.0); } /* Calculate the number of full months elapsed since * the current or most recent March */ if (date >= 0.0301) { iy = (int) (date + 0.0000000001); im = (int) (((date - (double) (iy)) * 10000.0) + 0.00000001); id = im % 100; im = (im / 100) + 9; if (im < 12) iy = iy - 1; im = im % 12; id = id - 1; /* starting with March as month 0 and ending with the following * February as month 11, the calculation of the number of days * per month reduces to a simple formula. the following statement * determines the number of whole days elapsed since 3/1/-1 and then * subtracts the 712163 days between then and 1/1/1950. it converts * the result to seconds and adds the accumulated seconds above. */ id = id + ((im+1+im/6+im/11)/2 * 31) + ((im-im/6-im/11)/2 * 30) + (iy / 4) - (iy / 100) + (iy / 400); dd = (double) id + (365.0 * (double) iy) - 712163.0; tsec = tsec + (dd * 86400.0); } return (tsec); } /* TS2DT-- convert seconds since 1950.0 to date, time as yyyy.mmdd hh.mmssss */ void ts2dt (tsec,date,time) double tsec; /* Seconds past 1950.0 */ double *date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double *time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ { int iyr,imon,iday,ihr,imn; double sec; ts2i (tsec,&iyr,&imon,&iday,&ihr,&imn,&sec, 4); /* Convert date to yyyy.mmdd */ if (iyr < 0) { *date = (double) (-iyr) + 0.01 * (double) imon + 0.0001 * (double) iday; *date = -(*date); } else *date = (double) iyr + 0.01 * (double) imon + 0.0001 * (double) iday; /* Convert time to hh.mmssssss */ *time = (double) ihr + 0.01 * (double) imn + 0.0001 * sec; return; } /* TSI2DT-- Convert seconds since 1980-01-01 to date yyyy.ddmm, time hh.mmsss */ void tsi2dt (isec,date,time) int isec; /* Seconds past 1980-01-01 */ double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { ts2dt (tsi2ts (isec), date, time); } /* TSI2FD-- Convert seconds since 1980-01-01 to FITS standard date string */ char * tsi2fd (isec) int isec; /* Seconds past 1980-01-01 */ { return (ts2fd (tsi2ts (isec))); } /* TSI2TS-- Convert seconds since 1980-01-01 to seconds since 1950-01-01 */ double tsi2ts (isec) int isec; /* Seconds past 1980-01-01 */ { return ((double) isec + 946684800.0); } /* TSU2FD-- Convert seconds since 1970-01-01 to FITS standard date string */ char * tsu2fd (isec) time_t isec; /* Seconds past 1970-01-01 */ { return (ts2fd (tsu2ts (isec))); } /* TSU2DT-- Convert seconds since 1970-01-01 to date yyyy.ddmm, time hh.mmsss */ void tsu2dt (isec,date,time) time_t isec; /* Seconds past 1970-01-01 */ double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { ts2dt (tsu2ts (isec), date, time); } /* TSU2TS-- Convert seconds since 1970-01-01 to seconds since 1950-01-01 */ double tsu2ts (isec) time_t isec; /* Seconds past 1970-01-01 */ { return ((double) isec + 631152000.0); } /* TSU2TSI-- UT seconds since 1970-01-01 to local seconds since 1980-01-01 */ int tsu2tsi (isec) time_t isec; /* Seconds past 1970-01-01 */ { double date, time; struct tm *ts; /* Get local time from UT seconds */ ts = localtime (&isec); if (ts->tm_year < 1000) date = (double) (ts->tm_year + 1900); else date = (double) ts->tm_year; date = date + (0.01 * (double) (ts->tm_mon + 1)); date = date + (0.0001 * (double) ts->tm_mday); time = (double) ts->tm_hour; time = time + (0.01 * (double) ts->tm_min); time = time + (0.0001 * (double) ts->tm_sec); return ((int)(dt2ts (date, time) - 631152000.0)); } /* TS2FD-- convert seconds since 1950.0 to FITS date, yyyy-mm-ddThh:mm:ss.ss */ char * ts2fd (tsec) double tsec; /* Seconds past 1950.0 */ { double date, time; ts2dt (tsec, &date, &time); return (dt2fd (date, time)); } /* TSD2FD-- convert seconds since start of day to FITS time, hh:mm:ss.ss */ char * tsd2fd (tsec) double tsec; /* Seconds since start of day */ { double date, time; char *thms, *fdate; int lfd, nbc; ts2dt (tsec, &date, &time); fdate = dt2fd (date, time); thms = (char *) calloc (16, 1); lfd = strlen (fdate); nbc = lfd - 11; strncpy (thms, fdate+11, nbc); return (thms); } /* TSD2DT-- convert seconds since start of day to hh.mmssss */ double tsd2dt (tsec) double tsec; /* Seconds since start of day */ { double date, time; ts2dt (tsec, &date, &time); return (time); } /* DT2I-- convert vigesimal date and time to year month day hours min sec */ void dt2i (date, time, iyr, imon, iday, ihr, imn, sec, ndsec) double date; /* Date as yyyy.mmdd (returned) yyyy = calendar year (e.g. 1973) mm = calendar month (e.g. 04 = april) dd = calendar day (e.g. 15) */ double time; /* Time as hh.mmssxxxx (returned) *if time<0, it is time as -(fraction of a day) hh = hour of day (0 .le. hh .le. 23) nn = minutes (0 .le. nn .le. 59) ss = seconds (0 .le. ss .le. 59) xxxx = tenths of milliseconds (0 .le. xxxx .le. 9999) */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double t,d; t = time; if (date < 0.0) d = -date; else d = date; /* Extract components of time */ *ihr = dint (t + 0.000000001); t = 100.0 * (t - (double) *ihr); *imn = dint (t + 0.0000001); *sec = 100.0 * (t - (double) *imn); /* Extract components of date */ *iyr = dint (d + 0.00001); d = 100.0 * (d - (double) *iyr); if (date < 0.0) *iyr = - *iyr; *imon = dint (d + 0.001); d = 100.0 * (d - (double) *imon); *iday = dint (d + 0.1); /* Make sure date and time are legal */ fixdate (iyr, imon, iday, ihr, imn, sec, ndsec); return; } /* FD2I-- convert from FITS standard date to year, mon, day, hours, min, sec */ void fd2i (string, iyr, imon, iday, ihr, imn, sec, ndsec) char *string; /* FITS date string, which may be: yyyy.ffff (fractional year) dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double tsec, fday, hr, mn; int i; char *sstr, *dstr, *tstr, *cstr, *nval, *fstr; /* Initialize all returned data to zero */ *iyr = 0; *imon = 0; *iday = 0; *ihr = 0; *imn = 0; *sec = 0.0; /* Return if no input string */ if (string == NULL) return; /* Check for various non-numeric characters */ sstr = strchr (string,'/'); dstr = strchr (string,'-'); if (dstr == string) dstr = strchr (string+1, '-'); fstr = strchr (string, '.'); tstr = strchr (string,'T'); if (tstr == NULL) tstr = strchr (string, 'Z'); if (tstr == NULL) tstr = strchr (string, 'S'); if (fstr != NULL && tstr != NULL && fstr > tstr) fstr = NULL; cstr = strchr (string,':'); /* Original FITS date format: dd/mm/yy */ if (sstr > string) { *sstr = '\0'; *iday = (int) atof (string); if (*iday > 31) { *iyr = *iday; if (*iyr >= 0 && *iyr <= 49) *iyr = *iyr + 2000; else if (*iyr < 1000) *iyr = *iyr + 1900; *sstr = '/'; nval = sstr + 1; sstr = strchr (nval,'/'); if (sstr > string) { *sstr = '\0'; *imon = (int) atof (nval); *sstr = '/'; nval = sstr + 1; *iday = (int) atof (nval); } } else { *sstr = '/'; nval = sstr + 1; sstr = strchr (nval,'/'); if (sstr == NULL) sstr = strchr (nval,'-'); if (sstr > string) { *sstr = '\0'; *imon = (int) atof (nval); *sstr = '/'; nval = sstr + 1; *iyr = (int) atof (nval); if (*iyr >= 0 && *iyr <= 49) *iyr = *iyr + 2000; else if (*iyr < 1000) *iyr = *iyr + 1900; } } tstr = strchr (string,'_'); if (tstr == NULL) return; } /* New FITS date format: yyyy-mm-ddThh:mm:ss[.sss] */ else if (dstr > string) { *dstr = '\0'; *iyr = (int) atof (string); *dstr = '-'; nval = dstr + 1; dstr = strchr (nval,'-'); *imon = 1; *iday = 1; /* Decode year, month, and day */ if (dstr > string) { *dstr = '\0'; *imon = (int) atof (nval); *dstr = '-'; nval = dstr + 1; if (tstr > string) *tstr = '\0'; *iday = (int) atof (nval); /* If fraction of a day is present, turn it into a time */ if (fstr != NULL) { fday = atof (fstr); hr = fday * 24.0; *ihr = (int) hr; mn = 60.0 * (hr - (double) *ihr); *imn = (int) mn; *sec = 60.0 * (mn - (double) *imn); } if (tstr > string) *tstr = 'T'; } /* If date is > 31, it is really year in old format */ if (*iday > 31) { i = *iyr; if (*iday < 100) *iyr = *iday + 1900; else *iyr = *iday; *iday = i; } } /* In rare cases, a FITS time is entered as an epoch */ else if (tstr == NULL && cstr == NULL && isnum (string)) { tsec = ep2ts (atof (string)); ts2i (tsec,iyr,imon,iday,ihr,imn,sec, ndsec); return; } /* Extract time, if it is present */ if (tstr > string || cstr > string) { if (tstr > string) nval = tstr + 1; else nval = string; cstr = strchr (nval,':'); if (cstr > string) { *cstr = '\0'; *ihr = (int) atof (nval); *cstr = ':'; nval = cstr + 1; cstr = strchr (nval,':'); if (cstr > string) { *cstr = '\0'; *imn = (int) atof (nval); *cstr = ':'; nval = cstr + 1; *sec = atof (nval); } else *imn = (int) atof (nval); } else *ihr = (int) atof (nval); } else ndsec = -1; /* Make sure date and time are legal */ fixdate (iyr, imon, iday, ihr, imn, sec, ndsec); return; } /* TS2I-- convert sec since 1950.0 to year month day hours minutes seconds */ void ts2i (tsec,iyr,imon,iday,ihr,imn,sec, ndsec) double tsec; /* seconds since 1/1/1950 0:00 */ int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double t,days, ts, dts; int nc,nc4,nly,ny,m,im; /* Round seconds to 0 - 4 decimal places */ ts = tsec + 61530883200.0; if (ts < 0.0) dts = -0.5; else dts = 0.5; if (ndsec < 1) t = dint (ts + dts) * 10000.0; else if (ndsec < 2) t = dint (ts * 10.0 + dts) * 1000.0; else if (ndsec < 3) t = dint (ts * 100.0 + dts) * 100.0; else if (ndsec < 4) t = dint (ts * 1000.0 + dts) * 10.0; else t = dint (ts * 10000.0 + dts); ts = t / 10000.0; /* Time of day (hours, minutes, seconds */ *ihr = (int) (dmod (ts/3600.0, 24.0)); *imn = (int) (dmod (ts/60.0, 60.0)); *sec = dmod (ts, 60.0); /* Number of days since 0 hr 0/0/0000 */ days = dint ((t / 864000000.0) + 0.000001); /* Number of leap centuries (400 years) */ nc4 = (int) ((days / 146097.0) + 0.00001); /* Number of centuries since last /400 */ days = days - (146097.0 * (double) (nc4)); nc = (int) ((days / 36524.0) + 0.000001); if (nc > 3) nc = 3; /* Number of leap years since last century */ days = days - (36524.0 * nc); nly = (int) ((days / 1461.0) + 0.0000000001); /* Number of years since last leap year */ days = days - (1461.0 * (double) nly); ny = (int) ((days / 365.0) + 0.00000001); if (ny > 3) ny = 3; /* Day of month */ days = days - (365.0 * (double) ny); if (days < 0) { m = 0; *iday = 29; } else { *iday = (int) (days + 0.00000001) + 1; for (m = 1; m <= 12; m++) { im = (m + ((m - 1) / 5)) % 2; /* fprintf (stderr,"%d %d %d %d\n", m, im, *iday, nc); */ if (*iday-1 < im+30) break; *iday = *iday - im - 30; } } /* Month */ *imon = ((m+1) % 12) + 1; /* Year */ *iyr = nc4*400 + nc*100 + nly*4 + ny + m/11; /* Make sure date and time are legal */ fixdate (iyr, imon, iday, ihr, imn, sec, ndsec); return; } /* UT2DOY-- Current Universal Time as year, day of year */ void ut2doy (year, doy) int *year; /* Year (returned) */ double *doy; /* Day of year (returned) */ { double date, time; ut2dt (&date, &time); dt2doy (date, time, year, doy); return; } /* UT2DT-- Current Universal Time as date (yyyy.mmdd) and time (hh.mmsss) */ void ut2dt(date, time) double *date; /* Date as yyyy.mmdd (returned) */ double *time; /* Time as hh.mmssxxxx (returned) */ { time_t tsec; struct timeval tp; struct timezone tzp; struct tm *ts; gettimeofday (&tp,&tzp); tsec = tp.tv_sec; ts = gmtime (&tsec); if (ts->tm_year < 1000) *date = (double) (ts->tm_year + 1900); else *date = (double) ts->tm_year; *date = *date + (0.01 * (double) (ts->tm_mon + 1)); *date = *date + (0.0001 * (double) ts->tm_mday); *time = (double) ts->tm_hour; *time = *time + (0.01 * (double) ts->tm_min); *time = *time + (0.0001 * (double) ts->tm_sec); return; } /* UT2EP-- Return current Universal Time as fractional year */ double ut2ep() { return (jd2ep (ut2jd())); } /* UT2EPB-- Return current Universal Time as Besselian epoch */ double ut2epb() { return (jd2epb (ut2jd())); } /* UT2EPJ-- Return current Universal Time as Julian epoch */ double ut2epj() { return (jd2epj (ut2jd())); } /* UT2FD-- Return current Universal Time as FITS ISO date string */ char * ut2fd() { int year, month, day, hour, minute, second; time_t tsec; struct timeval tp; struct timezone tzp; struct tm *ts; char *isotime; gettimeofday (&tp,&tzp); tsec = tp.tv_sec; ts = gmtime (&tsec); year = ts->tm_year; if (year < 1000) year = year + 1900; month = ts->tm_mon + 1; day = ts->tm_mday; hour = ts->tm_hour; minute = ts->tm_min; second = ts->tm_sec; isotime = (char *) calloc (32, sizeof (char)); sprintf (isotime, "%04d-%02d-%02dT%02d:%02d:%02d", year, month, day, hour, minute, second); return (isotime); } /* UT2JD-- Return current Universal Time as Julian Date */ double ut2jd() { return (fd2jd (ut2fd())); } /* UT2MJD-- convert current UT to Modified Julian Date */ double ut2mjd () { return (ut2jd() - 2400000.5); } /* UT2TS-- current Universal Time as IRAF seconds since 1950-01-01T00:00 */ double ut2ts() { double tsec; char *datestring; datestring = ut2fd(); tsec = fd2ts (datestring); free (datestring); return (tsec); } /* UT2TSI-- current Universal Time as IRAF seconds since 1980-01-01T00:00 */ int ut2tsi() { return ((int)(ut2ts() - 946684800.0)); } /* UT2TSU-- current Universal Time as IRAF seconds since 1970-01-01T00:00 */ time_t ut2tsu() { return ((time_t)(ut2ts () - 631152000.0)); } /* FD2GST-- convert from FITS date to Greenwich Sidereal Time */ char * fd2gst (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double dj, gsec, date, time; dj = fd2jd (string); gsec = jd2gst (dj); ts2dt (gsec, &date, &time); date = 0.0; return (dt2fd (date, time)); } /* DT2GST-- convert from UT as yyyy.mmdd hh.mmssss to Greenwich Sidereal Time*/ void dt2gst (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, gsec; dj = dt2ts (*date, *time); gsec = jd2gst (dj); ts2dt (gsec, date, time); *date = 0.0; return; } /* JD2LST - Local Sidereal Time in seconds from Julian Date */ double jd2lst (dj) double dj; /* Julian Date */ { double gst, lst; /* Compute Greenwich Sidereal Time at this epoch */ gst = jd2gst (dj); /* Subtract longitude (degrees to seconds of time) */ lst = gst - (240.0 * longitude); if (lst < 0.0) lst = lst + 86400.0; else if (lst > 86400.0) lst = lst - 86400.0; return (lst); } /* FD2LST - Local Sidereal Time as hh:mm:ss.ss from Universal Time as FITS ISO date */ char * fd2lst (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) */ { double dj, date, time, lst; dj = fd2jd (string); lst = jd2lst (dj); ts2dt (lst, &date, &time); date = 0.0; return (dt2fd (date, time)); } /* DT2LST - Local Sidereal Time as hh.mmssss from Universal Time as yyyy.mmdd hh.mmssss */ void dt2lst (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double dj, lst, date0; dj = dt2jd (*date, *time); lst = jd2lst (dj); date0 = 0.0; ts2dt (lst, &date0, time); return; } /* TS2LST - Local Sidereal Time in seconds of day * from Universal Time in seconds since 1951-01-01T0:00:00 */ double ts2lst (tsec) double tsec; /* time since 1950.0 in UT seconds */ { double gst; /* Greenwich Sidereal Time in seconds since 0:00 */ double lst; /* Local Sidereal Time in seconds since 0:00 */ double gsec, date; /* Greenwich Sidereal Time */ gsec = ts2gst (tsec); date = 0.0; ts2dt (gsec, &date, &gst); lst = gst - (longitude / 15.0); if (lst < 0.0) lst = lst + 86400.0; else if (lst > 86400.0) lst = lst - 86400.0; return (lst); } /* LST2FD - calculate current UT given Local Sidereal Time * plus date in FITS ISO format (yyyy-mm-dd) * Return UT date and time in FITS ISO format */ char * lst2fd (string) char *string; /* UT Date, LST as yyyy-mm-ddShh:mm:ss.ss */ { double sdj, dj; sdj = fd2jd (string); dj = lst2jd (sdj); return (jd2fd (dj)); } /* LST2JD - calculate current Julian Date given Local Sidereal Time * plus current Julian Date (0.5 at 0:00 UT) * Return UT date and time as Julian Date */ double lst2jd (sdj) double sdj; /* Julian Date of desired day at 0:00 UT + sidereal time */ { double gst; /* Greenwich Sidereal Time in seconds since 0:00 */ double lsd; /* Local Sidereal Time in seconds since 0:00 */ double gst0, tsd, dj1, dj0, eqnx; int idj; /* Julian date at 0:00 UT */ idj = (int) sdj; dj0 = (double) idj + 0.5; if (dj0 > sdj) dj0 = dj0 - 1.0; /* Greenwich Sidereal Time at 0:00 UT in seconds */ gst0 = jd2gst (dj0); /* Sidereal seconds since 0:00 */ lsd = (sdj - dj0) * 86400.0; /* Remove longitude for current Greenwich Sidereal Time in seconds */ /* (convert longitude from degrees to seconds of time) */ gst = lsd + (longitude * 240.0); /* Time since 0:00 UT */ tsd = (gst - gst0) / 1.0027379093; /* Julian Date (UT) */ dj1 = dj0 + (tsd / 86400.0); /* Equation of the equinoxes converted to UT seconds */ eqnx = eqeqnx (dj1) / 1.002739093; /* Remove equation of equinoxes */ dj1 = dj1 - (eqnx / 86400.0); if (dj1 < dj0) dj1 = dj1 + 1.0; return (dj1); } /* MST2FD - calculate current UT given Greenwich Mean Sidereal Time * plus date in FITS ISO format (yyyy-mm-ddShh:mm:ss.ss) * Return UT date and time in FITS ISO format */ char * mst2fd (string) char *string; /* UT Date, MST as yyyy-mm-ddShh:mm:ss.ss */ { double sdj, dj; sdj = fd2jd (string); dj = mst2jd (sdj); return (jd2fd (dj)); } /* MST2JD - calculate current UT given Greenwich Mean Sidereal Time * plus date in Julian Date (0:00 UT + Mean Sidereal Time) * Return UT date and time as Julian Date */ double mst2jd (sdj) double sdj; /* UT Date, MST as Julian Date */ { double tsd, djd, st0, dj0, dj; dj0 = (double) ((int) sdj) + 0.5; /* Greenwich Mean Sidereal Time at 0:00 UT in seconds */ st0 = jd2mst (dj0); /* Mean Sidereal Time in seconds */ tsd = (sdj - dj0) * 86400.0; if (tsd < 0.0) tsd = tsd + 86400.0; /* Convert to fraction of a day since 0:00 UT */ djd = ((tsd - st0) / 1.0027379093) / 86400.0; /* Julian Date */ dj = dj0 + djd; if (dj < dj0) dj = dj + (1.0 / 1.0027379093); return (dj); } /* GST2FD - calculate current UT given Greenwich Sidereal Time * plus date in FITS ISO format (yyyy-mm-ddShh:mm:ss.ss) * Return UT date and time in FITS ISO format */ char * gst2fd (string) char *string; /* UT Date, GST as yyyy-mm-ddShh:mm:ss.ss */ { double sdj, dj; sdj = fd2jd (string); dj = gst2jd (sdj); return (jd2fd (dj)); } /* GST2JD - calculate current UT given Greenwich Sidereal Time * plus date as Julian Date (JD at 0:00 UT + sidereal time) * Return UT date and time as Julian Date */ double gst2jd (sdj) double sdj; /* UT Date, GST as Julian Date */ { double dj, tsd, djd, st0, dj0, eqnx; dj0 = (double) ((int) sdj) + 0.5; /* Greenwich Mean Sidereal Time at 0:00 UT in seconds */ st0 = jd2mst (dj0); /* Mean Sidereal Time in seconds */ tsd = (sdj - dj0) * 86400.0; if (tsd < 0.0) tsd = tsd + 86400.0; /* Convert to fraction of a day since 0:00 UT */ djd = ((tsd - st0) / 1.0027379093) / 86400.0; /* Julian Date */ dj = dj0 + djd; /* Equation of the equinoxes (converted to UT seconds) */ eqnx = eqeqnx (dj) / 1.002737909; dj = dj - eqnx / 86400.0; if (dj < dj0) dj = dj + 1.0; return (dj); } /* LST2DT - calculate current UT given Local Sidereal Time as hh.mmsss * plus date as yyyy.mmdd * Return UT time as hh.mmssss */ double lst2dt (date0, time0) double date0; /* UT date as yyyy.mmdd */ double time0; /* LST as hh.mmssss */ { double gst; /* Greenwich Sidereal Time in seconds since 0:00 */ double lst; /* Local Sidereal Time in seconds since 0:00 */ double date1; /* UT date as yyyy.mmdd */ double time1; /* UT as hh.mmssss */ double tsec0, gst0, tsd, tsec; /* Greenwich Sidereal Time at 0:00 UT */ tsec0 = dt2ts (date0, 0.0); gst0 = ts2gst (tsec0); /* Current Greenwich Sidereal Time in seconds */ /* (convert longitude from degrees to seconds of time) */ lst = dt2ts (0.0, time0); gst = lst + (longitude * 240.0); /* Time since 0:00 UT */ tsd = (gst - gst0) / 1.0027379093; /* UT date and time */ tsec = tsec0 + tsd; ts2dt (tsec, &date1, &time1); return (time1); } /* TS2GST - calculate Greenwich Sidereal Time given Universal Time * in seconds since 1951-01-01T0:00:00 * Return sidereal time of day in seconds */ double ts2gst (tsec) double tsec; /* time since 1950.0 in UT seconds */ { double gst; /* Greenwich Sidereal Time in seconds since 0:00 */ double tsd, eqnx, dj; int its; /* Elapsed time as of 0:00 UT */ if (tsec >= 0.0) { its = (int) (tsec + 0.5); tsd = (double) (its % 86400); } else { its = (int) (-tsec + 0.5); tsd = (double) (86400 - (its % 86400)); } /* Mean sidereal time */ gst = ts2mst (tsec); /* Equation of the equinoxes */ dj = ts2jd (tsec); eqnx = eqeqnx (dj); /* Apparent sidereal time at 0:00 ut */ gst = gst + eqnx; /* Current sidereal time */ gst = gst + (tsd * 1.0027379093); gst = dmod (gst,86400.0); return (gst); } /* FD2MST-- convert from FITS date Mean Sidereal Time */ char * fd2mst (string) char *string; /* FITS date string, which may be: fractional year dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { double gsec, date, time, dj; dj = fd2jd (string); gsec = jd2mst (dj); ts2dt (gsec, &date, &time); date = 0.0; return (dt2fd (date, time)); } /* DT2MST-- convert from UT as yyyy.mmdd hh.mmssss to Mean Sidereal Time in the same format */ void dt2mst (date, time) double *date; /* Date as yyyy.mmdd */ double *time; /* Time as hh.mmssxxxx *if time<0, it is time as -(fraction of a day) */ { double date0, gsec, dj; date0 = *date; dj = dt2jd (*date, *time); gsec = jd2mst (dj); ts2dt (gsec, date, time); *date = date0; return; } /* TS2MST - calculate Greenwich Mean Sidereal Time given Universal Time * in seconds since 1951-01-01T0:00:00 */ double ts2mst (tsec) double tsec; /* time since 1950.0 in UT seconds */ { double dj; dj = ts2jd (tsec); return (jd2mst (dj)); } /* JD2MST - Julian Date to Greenwich Mean Sidereal Time using IAU 2000 * Return sideral time in seconds of time * (from USNO NOVAS package * http://aa.usno.navy.mil/software/novas/novas_info.html */ double jd2mst2 (dj) double dj; /* Julian Date */ { double dt, t, t2, t3, mst, st; dt = dj - 2451545.0; t = dt / 36525.0; t2 = t * t; t3 = t2 * t; /* Compute Greenwich Mean Sidereal Time in seconds */ st = (8640184.812866 * t) + (3155760000.0 * t) - (0.0000062 * t3) + (0.093104 * t2) + 67310.54841; mst = dmod (st, 86400.0); if (mst < 0.0) mst = mst + 86400.0; return (mst); } /* MJD2MST - Modified Julian Date to Greenwich Mean Sidereal Time using IAU 2000 * Return sideral time in seconds of time * (from USNO NOVAS package * http://aa.usno.navy.mil/software/novas/novas_info.html */ double mjd2mst (dj) double dj; /* Modified Julian Date */ { double dt, t, t2, t3, mst, st; dt = dj - 51544.5; t = dt / 36525.0; t2 = t * t; t3 = t2 * t; /* Compute Greenwich Mean Sidereal Time in seconds */ st = (8640184.812866 * t) + (3155760000.0 * t) - (0.0000062 * t3) + (0.093104 * t2) + 67310.54841; mst = dmod (st, 86400.0); if (mst < 0.0) mst = mst + 86400.0; return (mst); } /* JD2GST - Julian Date to Greenwich Sideral Time * Return sideral time in seconds of time * (Jean Meeus, Astronomical Algorithms, Willmann-Bell, 1991, pp 83-84) */ double jd2gst (dj) double dj; /* Julian Date */ { double dj0, gmt, gst, tsd, eqnx, ssd, l0; double ts2ss = 1.00273790935; int ijd; /* Julian date at 0:00 UT */ ijd = (int) dj; dj0 = (double) ijd + 0.5; if (dj0 > dj) dj0 = dj0 - 1.0; /* Greenwich mean sidereal time at 0:00 UT in seconds */ l0 = longitude; longitude = 0.0; gmt = jd2mst (dj0); longitude = l0; /* Equation of the equinoxes */ eqnx = eqeqnx (dj); /* Apparent sidereal time at 0:00 ut */ gst = gmt + eqnx; /* UT seconds since 0:00 */ tsd = (dj - dj0) * 86400.0; ssd = tsd * ts2ss; /* Current sidereal time */ gst = gst + ssd; gst = dmod (gst, 86400.0); return (gst); } /* EQEQNX - Compute equation of the equinoxes for apparent sidereal time */ double eqeqnx (dj) double dj; /* Julian Date */ { double dt, edj, dpsi, deps, obl, eqnx; double rad2tsec = 13750.98708; /* Convert UT to Ephemeris Time (TDB or TT)*/ dt = utdt (dj); edj = dj + dt / 86400.0; /* Nutation and obliquity */ compnut (edj, &dpsi, &deps, &obl); /* Correct obliquity for nutation */ obl = obl + deps; /* Equation of the equinoxes in seconds */ eqnx = (dpsi * cos (obl)) * rad2tsec; return (eqnx); } /* JD2MST - Julian Date to Mean Sideral Time * Return sideral time in seconds of time * (Jean Meeus, Astronomical Algorithms, Willmann-Bell, 1991, pp 83-84) */ double jd2mst (dj) double dj; /* Julian Date */ { double dt, t, mst; dt = dj - 2451545.0; t = dt / 36525.0; /* Compute Greenwich mean sidereal time in degrees (Meeus, page 84) */ mst = 280.46061837 + (360.98564736629 * dt) + (0.000387933 * t * t) - (t * t * t / 38710000.0); /* Keep degrees between 0 and 360 */ while (mst > 360.0) mst = mst - 360.0; while (mst < 0.0) mst = mst + 360.0; /* Convert to time in seconds (3600 / 15) */ mst = mst * 240.0; /* Subtract longitude (degrees to seconds of time) */ mst = mst - (240.0 * longitude); if (mst < 0.0) mst = mst + 86400.0; else if (mst > 86400.0) mst = mst - 86400.0; return (mst); } /* COMPNUT - Compute nutation using the IAU 2000b model */ /* Translated from Pat Wallace's Fortran subroutine iau_nut00b (June 26 2007) into C by Jessica Mink on September 5, 2008 */ #define NLS 77 /* number of terms in the luni-solar nutation model */ void compnut (dj, dpsi, deps, eps0) double dj; /* Julian Date */ double *dpsi; /* Nutation in longitude in radians (returned) */ double *deps; /* Nutation in obliquity in radians (returned) */ double *eps0; /* Mean obliquity in radians (returned) */ /* This routine is translated from the International Astronomical Union's * Fortran SOFA (Standards Of Fundamental Astronomy) software collection. * * notes: * * 1) the nutation components in longitude and obliquity are in radians * and with respect to the equinox and ecliptic of date. the * obliquity at j2000 is assumed to be the lieske et al. (1977) value * of 84381.448 arcsec. (the errors that result from using this * routine with the iau 2006 value of 84381.406 arcsec can be * neglected.) * * the nutation model consists only of luni-solar terms, but includes * also a fixed offset which compensates for certain long-period * planetary terms (note 7). * * 2) this routine is an implementation of the iau 2000b abridged * nutation model formally adopted by the iau general assembly in * 2000. the routine computes the mhb_2000_short luni-solar nutation * series (luzum 2001), but without the associated corrections for * the precession rate adjustments and the offset between the gcrs * and j2000 mean poles. * * 3) the full IAU 2000a (mhb2000) nutation model contains nearly 1400 * terms. the IAU 2000b model (mccarthy & luzum 2003) contains only * 77 terms, plus additional simplifications, yet still delivers * results of 1 mas accuracy at present epochs. this combination of * accuracy and size makes the IAU 2000b abridged nutation model * suitable for most practical applications. * * the routine delivers a pole accurate to 1 mas from 1900 to 2100 * (usually better than 1 mas, very occasionally just outside 1 mas). * the full IAU 2000a model, which is implemented in the routine * iau_nut00a (q.v.), delivers considerably greater accuracy at * current epochs; however, to realize this improved accuracy, * corrections for the essentially unpredictable free-core-nutation * (fcn) must also be included. * * 4) the present routine provides classical nutation. the * mhb_2000_short algorithm, from which it is adapted, deals also * with (i) the offsets between the gcrs and mean poles and (ii) the * adjustments in longitude and obliquity due to the changed * precession rates. these additional functions, namely frame bias * and precession adjustments, are supported by the sofa routines * iau_bi00 and iau_pr00. * * 6) the mhb_2000_short algorithm also provides "total" nutations, * comprising the arithmetic sum of the frame bias, precession * adjustments, and nutation (luni-solar + planetary). these total * nutations can be used in combination with an existing IAU 1976 * precession implementation, such as iau_pmat76, to deliver gcrs-to- * true predictions of mas accuracy at current epochs. however, for * symmetry with the iau_nut00a routine (q.v. for the reasons), the * sofa routines do not generate the "total nutations" directly. * should they be required, they could of course easily be generated * by calling iau_bi00, iau_pr00 and the present routine and adding * the results. * * 7) the IAU 2000b model includes "planetary bias" terms that are fixed * in size but compensate for long-period nutations. the amplitudes * quoted in mccarthy & luzum (2003), namely dpsi = -1.5835 mas and * depsilon = +1.6339 mas, are optimized for the "total nutations" * method described in note 6. the luzum (2001) values used in this * sofa implementation, namely -0.135 mas and +0.388 mas, are * optimized for the "rigorous" method, where frame bias, precession * and nutation are applied separately and in that order. during the * interval 1995-2050, the sofa implementation delivers a maximum * error of 1.001 mas (not including fcn). * * References from original Fortran subroutines: * * Hilton, J. et al., 2006, Celest.Mech.Dyn.Astron. 94, 351 * * Lieske, J.H., Lederle, T., Fricke, W., Morando, B., "Expressions * for the precession quantities based upon the IAU 1976 system of * astronomical constants", Astron.Astrophys. 58, 1-2, 1-16. (1977) * * Luzum, B., private communication, 2001 (Fortran code * mhb_2000_short) * * McCarthy, D.D. & Luzum, B.J., "An abridged model of the * precession-nutation of the celestial pole", Cel.Mech.Dyn.Astron. * 85, 37-49 (2003) * * Simon, J.-L., Bretagnon, P., Chapront, J., Chapront-Touze, M., * Francou, G., Laskar, J., Astron.Astrophys. 282, 663-683 (1994) * */ { double as2r = 0.000004848136811095359935899141; /* arcseconds to radians */ double dmas2r = as2r / 1000.0; /* milliarcseconds to radians */ double as2pi = 1296000.0; /* arc seconds in a full circle */ double d2pi = 6.283185307179586476925287; /* 2pi */ double u2r = as2r / 10000000.0; /* units of 0.1 microarcsecond to radians */ double dj0 = 2451545.0; /* reference epoch (j2000), jd */ double djc = 36525.0; /* Days per julian century */ /* Miscellaneous */ double t, el, elp, f, d, om, arg, dp, de, sarg, carg; double dpsils, depsls, dpsipl, depspl; int i, j; int nls = NLS; /* number of terms in the luni-solar nutation model */ /* Fixed offset in lieu of planetary terms (radians) */ double dpplan = - 0.135 * dmas2r; double deplan = + 0.388 * dmas2r; /* Tables of argument and term coefficients */ /* Coefficients for fundamental arguments */ /* Luni-solar argument multipliers: */ /* l l' f d om */ static int nals[5*NLS]= {0, 0, 0, 0, 1, 0, 0, 2, -2, 2, 0, 0, 2, 0, 2, 0, 0, 0, 0, 2, 0, 1, 0, 0, 0, 0, 1, 2, -2, 2, 1, 0, 0, 0, 0, 0, 0, 2, 0, 1, 1, 0, 2, 0, 2, 0, -1, 2, -2, 2, 0, 0, 2, -2, 1, -1, 0, 2, 0, 2, -1, 0, 0, 2, 0, 1, 0, 0, 0, 1, -1, 0, 0, 0, 1, -1, 0, 2, 2, 2, 1, 0, 2, 0, 1, -2, 0, 2, 0, 1, 0, 0, 0, 2, 0, 0, 0, 2, 2, 2, 0, -2, 2, -2, 2, -2, 0, 0, 2, 0, 2, 0, 2, 0, 2, 1, 0, 2, -2, 2, -1, 0, 2, 0, 1, 2, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 1, 0, 0, 1, -1, 0, 0, 2, 1, 0, 2, 2, -2, 2, 0, 0, -2, 2, 0, 1, 0, 0, -2, 1, 0, -1, 0, 0, 1, -1, 0, 2, 2, 1, 0, 2, 0, 0, 0, 1, 0, 2, 2, 2, -2, 0, 2, 0, 0, 0, 1, 2, 0, 2, 0, 0, 2, 2, 1, 0, -1, 2, 0, 2, 0, 0, 0, 2, 1, 1, 0, 2, -2, 1, 2, 0, 2, -2, 2, -2, 0, 0, 2, 1, 2, 0, 2, 0, 1, 0, -1, 2, -2, 1, 0, 0, 0, -2, 1, -1, -1, 0, 2, 0, 2, 0, 0, -2, 1, 1, 0, 0, 2, 0, 0, 1, 2, -2, 1, 1, -1, 0, 0, 0, -2, 0, 2, 0, 2, 3, 0, 2, 0, 2, 0, -1, 0, 2, 0, 1, -1, 2, 0, 2, 0, 0, 0, 1, 0, -1, -1, 2, 2, 2, -1, 0, 2, 0, 0, 0, -1, 2, 2, 2, -2, 0, 0, 0, 1, 1, 1, 2, 0, 2, 2, 0, 0, 0, 1, -1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 2, 0, 0, -1, 0, 2, -2, 1, 1, 0, 0, 0, 2, -1, 0, 0, 1, 0, 0, 0, 2, 1, 2, -1, 0, 2, 4, 2, -1, 1, 0, 1, 1, 0, -2, 2, -2, 1, 1, 0, 2, 2, 1, -2, 0, 2, 2, 2, -1, 0, 0, 0, 2, 1, 1, 2, -2, 2}; /* Luni-solar nutation coefficients, in 1e-7 arcsec */ /* longitude (sin, t*sin, cos), obliquity (cos, t*cos, sin) */ static double cls[6*NLS]= {-172064161.0, -174666.0, 33386.0, 92052331.0, 9086.0, 15377.0, -13170906.0, -1675.0, -13696.0, 5730336.0, -3015.0, -4587.0, -2276413.0, -234.0, 2796.0, 978459.0, -485.0, 1374.0, 2074554.0, 207.0, -698.0, -897492.0, 470.0, -291.0, 1475877.0, -3633.0, 11817.0, 73871.0, -184.0, -1924.0, -516821.0, 1226.0, -524.0, 224386.0, -677.0, -174.0, 711159.0, 73.0, -872.0, -6750.0, 0.0, 358.0, -387298.0, -367.0, 380.0, 200728.0, 18.0, 318.0, -301461.0, -36.0, 816.0, 129025.0, -63.0, 367.0, 215829.0, -494.0, 111.0, -95929.0, 299.0, 132.0, 128227.0, 137.0, 181.0, -68982.0, -9.0, 39.0, 123457.0, 11.0, 19.0, -53311.0, 32.0, -4.0, 156994.0, 10.0, -168.0, -1235.0, 0.0, 82.0, 63110.0, 63.0, 27.0, -33228.0, 0.0, -9.0, -57976.0, -63.0, -189.0, 31429.0, 0.0, -75.0, -59641.0, -11.0, 149.0, 25543.0, -11.0, 66.0, -51613.0, -42.0, 129.0, 26366.0, 0.0, 78.0, 45893.0, 50.0, 31.0, -24236.0, -10.0, 20.0, 63384.0, 11.0, -150.0, -1220.0, 0.0, 29.0, -38571.0, -1.0, 158.0, 16452.0, -11.0, 68.0, 32481.0, 0.0, 0.0, -13870.0, 0.0, 0.0, -47722.0, 0.0, -18.0, 477.0, 0.0, -25.0, -31046.0, -1.0, 131.0, 13238.0, -11.0, 59.0, 28593.0, 0.0, -1.0, -12338.0, 10.0, -3.0, 20441.0, 21.0, 10.0, -10758.0, 0.0, -3.0, 29243.0, 0.0, -74.0, -609.0, 0.0, 13.0, 25887.0, 0.0, -66.0, -550.0, 0.0, 11.0, -14053.0, -25.0, 79.0, 8551.0, -2.0, -45.0, 15164.0, 10.0, 11.0, -8001.0, 0.0, -1.0, -15794.0, 72.0, -16.0, 6850.0, -42.0, -5.0, 21783.0, 0.0, 13.0, -167.0, 0.0, 13.0, -12873.0, -10.0, -37.0, 6953.0, 0.0, -14.0, -12654.0, 11.0, 63.0, 6415.0, 0.0, 26.0, -10204.0, 0.0, 25.0, 5222.0, 0.0, 15.0, 16707.0, -85.0, -10.0, 168.0, -1.0, 10.0, -7691.0, 0.0, 44.0, 3268.0, 0.0, 19.0, -11024.0, 0.0, -14.0, 104.0, 0.0, 2.0, 7566.0, -21.0, -11.0, -3250.0, 0.0, -5.0, -6637.0, -11.0, 25.0, 3353.0, 0.0, 14.0, -7141.0, 21.0, 8.0, 3070.0, 0.0, 4.0, -6302.0, -11.0, 2.0, 3272.0, 0.0, 4.0, 5800.0, 10.0, 2.0, -3045.0, 0.0, -1.0, 6443.0, 0.0, -7.0, -2768.0, 0.0, -4.0, -5774.0, -11.0, -15.0, 3041.0, 0.0, -5.0, -5350.0, 0.0, 21.0, 2695.0, 0.0, 12.0, -4752.0, -11.0, -3.0, 2719.0, 0.0, -3.0, -4940.0, -11.0, -21.0, 2720.0, 0.0, -9.0, 7350.0, 0.0, -8.0, -51.0, 0.0, 4.0, 4065.0, 0.0, 6.0, -2206.0, 0.0, 1.0, 6579.0, 0.0, -24.0, -199.0, 0.0, 2.0, 3579.0, 0.0, 5.0, -1900.0, 0.0, 1.0, 4725.0, 0.0, -6.0, -41.0, 0.0, 3.0, -3075.0, 0.0, -2.0, 1313.0, 0.0, -1.0, -2904.0, 0.0, 15.0, 1233.0, 0.0, 7.0, 4348.0, 0.0, -10.0, -81.0, 0.0, 2.0, -2878.0, 0.0, 8.0, 1232.0, 0.0, 4.0, -4230.0, 0.0, 5.0, -20.0, 0.0, -2.0, -2819.0, 0.0, 7.0, 1207.0, 0.0, 3.0, -4056.0, 0.0, 5.0, 40.0, 0.0, -2.0, -2647.0, 0.0, 11.0, 1129.0, 0.0, 5.0, -2294.0, 0.0, -10.0, 1266.0, 0.0, -4.0, 2481.0, 0.0, -7.0, -1062.0, 0.0, -3.0, 2179.0, 0.0, -2.0, -1129.0, 0.0, -2.0, 3276.0, 0.0, 1.0, -9.0, 0.0, 0.0, -3389.0, 0.0, 5.0, 35.0, 0.0, -2.0, 3339.0, 0.0, -13.0, -107.0, 0.0, 1.0, -1987.0, 0.0, -6.0, 1073.0, 0.0, -2.0, -1981.0, 0.0, 0.0, 854.0, 0.0, 0.0, 4026.0, 0.0, -353.0, -553.0, 0.0, -139.0, 1660.0, 0.0, -5.0, -710.0, 0.0, -2.0, -1521.0, 0.0, 9.0, 647.0, 0.0, 4.0, 1314.0, 0.0, 0.0, -700.0, 0.0, 0.0, -1283.0, 0.0, 0.0, 672.0, 0.0, 0.0, -1331.0, 0.0, 8.0, 663.0, 0.0, 4.0, 1383.0, 0.0, -2.0, -594.0, 0.0, -2.0, 1405.0, 0.0, 4.0, -610.0, 0.0, 2.0, 1290.0, 0.0, 0.0, -556.0, 0.0, 0.0}; /* Interval between fundamental epoch J2000.0 and given date (JC) */ t = (dj - dj0) / djc; /* Luni-solar nutation */ /* Fundamental (delaunay) arguments from Simon et al. (1994) */ /* Mean anomaly of the moon */ el = fmod (485868.249036 + (1717915923.2178 * t), as2pi) * as2r; /* Mean anomaly of the sun */ elp = fmod (1287104.79305 + (129596581.0481 * t), as2pi) * as2r; /* Mean argument of the latitude of the moon */ f = fmod (335779.526232 + (1739527262.8478 * t), as2pi) * as2r; /* Mean elongation of the moon from the sun */ d = fmod (1072260.70369 + (1602961601.2090 * t), as2pi ) * as2r; /* Mean longitude of the ascending node of the moon */ om = fmod (450160.398036 - (6962890.5431 * t), as2pi ) * as2r; /* Initialize the nutation values */ dp = 0.0; de = 0.0; /* Summation of luni-solar nutation series (in reverse order) */ for (i = nls; i > 0; i=i-1) { j = i - 1; /* Argument and functions */ arg = fmod ( (double) (nals[5*j]) * el + (double) (nals[1+5*j]) * elp + (double) (nals[2+5*j]) * f + (double) (nals[3+5*j]) * d + (double) (nals[4+5*j]) * om, d2pi); sarg = sin (arg); carg = cos (arg); /* Terms */ dp = dp + (cls[6*j] + cls[1+6*j] * t) * sarg + cls[2+6*j] * carg; de = de + (cls[3+6*j] + cls[4+6*j] * t) * carg + cls[5+6*j] * sarg; } /* Convert from 0.1 microarcsec units to radians */ dpsils = dp * u2r; depsls = de * u2r; /* In lieu of planetary nutation */ /* Fixed offset to correct for missing terms in truncated series */ dpsipl = dpplan; depspl = deplan; /* Results */ /* Add luni-solar and planetary components */ *dpsi = dpsils + dpsipl; *deps = depsls + depspl; /* Mean Obliquity in radians (IAU 2006, Hilton, et al.) */ *eps0 = ( 84381.406 + ( -46.836769 + ( -0.0001831 + ( 0.00200340 + ( -0.000000576 + ( -0.0000000434 ) * t ) * t ) * t ) * t ) * t ) * as2r; } /* ISDATE - Return 1 if string is an old or ISO FITS standard date */ int isdate (string) char *string; /* Possible FITS date string, which may be: dd/mm/yy (FITS standard before 2000) dd-mm-yy (nonstandard FITS use before 2000) yyyy-mm-dd (FITS standard after 1999) yyyy-mm-ddThh:mm:ss.ss (FITS standard after 1999) */ { int iyr = 0; /* year (returned) */ int imon = 0; /* month (returned) */ int iday = 0; /* day (returned) */ int i; char *sstr, *dstr, *tstr, *nval; /* Translate string from ASCII to binary */ if (string == NULL) return (0); sstr = strchr (string,'/'); dstr = strchr (string,'-'); if (dstr == string) dstr = strchr (string+1,'-'); tstr = strchr (string,'T'); /* Original FITS date format: dd/mm/yy */ if (sstr > string) { *sstr = '\0'; iday = (int) atof (string); *sstr = '/'; nval = sstr + 1; sstr = strchr (nval,'/'); if (sstr == NULL) sstr = strchr (nval,'-'); if (sstr > string) { *sstr = '\0'; imon = (int) atof (nval); *sstr = '/'; nval = sstr + 1; iyr = (int) atof (nval); if (iyr < 1000) iyr = iyr + 1900; } if (imon > 0 && iday > 0) return (1); else return (0); } /* New FITS date format: yyyy-mm-ddThh:mm:ss[.sss] */ else if (dstr > string) { *dstr = '\0'; iyr = (int) atof (string); nval = dstr + 1; *dstr = '-'; dstr = strchr (nval,'-'); imon = 0; iday = 0; /* Decode year, month, and day */ if (dstr > string) { *dstr = '\0'; imon = (int) atof (nval); *dstr = '-'; nval = dstr + 1; if (tstr > string) *tstr = '\0'; iday = (int) atof (nval); if (tstr > string) *tstr = 'T'; } /* If day is > 31, it is really year in old format */ if (iday > 31) { i = iyr; if (iday < 100) iyr = iday + 1900; else iyr = iday; iday = i; } if (imon > 0 && iday > 0) return (1); else return (0); } /* If FITS date is entered as an epoch, return 0 anyway */ else return (0); } /* Round seconds and make sure date and time numbers are within limits */ static void fixdate (iyr, imon, iday, ihr, imn, sec, ndsec) int *iyr; /* year (returned) */ int *imon; /* month (returned) */ int *iday; /* day (returned) */ int *ihr; /* hours (returned) */ int *imn; /* minutes (returned) */ double *sec; /* seconds (returned) */ int ndsec; /* Number of decimal places in seconds (0=int) */ { double days; /* Round seconds to 0 - 4 decimal places (no rounding if <0, >4) */ if (ndsec == 0) *sec = dint (*sec + 0.5); else if (ndsec < 2) *sec = dint (*sec * 10.0 + 0.5) / 10.0; else if (ndsec < 3) *sec = dint (*sec * 100.0 + 0.5) / 100.0; else if (ndsec < 4) *sec = dint (*sec * 1000.0 + 0.5) / 1000.0; else if (ndsec < 5) *sec = dint (*sec * 10000.0 + 0.5) / 10000.0; /* Adjust minutes and hours */ if (*sec > 60.0) { *sec = *sec - 60.0; *imn = *imn + 1; } if (*imn > 60) { *imn = *imn - 60; *ihr = *ihr + 1; } /* Return if no date */ if (*iyr == 0 && *imon == 0 && *iday == 0) return; /* Adjust date */ if (*ihr > 23) { *ihr = *ihr - 24; *iday = *iday + 1; } days = caldays (*iyr, *imon); if (*iday > days) { *iday = *iday - days; *imon = *imon + 1; } if (*iday < 1) { *imon = *imon - 1; if (*imon < 1) { *imon = *imon + 12; *iyr = *iyr - 1; } days = caldays (*iyr, *imon); *iday = *iday + days; } if (*imon < 1) { *imon = *imon + 12; *iyr = *iyr - 1; days = caldays (*iyr, *imon); if (*iday > days) { *iday = *iday - days; *imon = *imon + 1; } } if (*imon > 12) { *imon = *imon - 12; *iyr = *iyr + 1; } return; } /* Calculate days in month 1-12 given year (Gregorian calendar only) */ static int caldays (year, month) int year; /* 4-digit year */ int month; /* Month (1=January, 2=February, etc.) */ { if (month < 1) { month = month + 12; year = year + 1; } if (month > 12) { month = month - 12; year = year + 1; } switch (month) { case 1: return (31); case 2: if (year%400 == 0) return (29); else if (year%100 == 0) return (28); else if (year%4 == 0) return (29); else return (28); case 3: return (31); case 4: return (30); case 5: return (31); case 6: return (30); case 7: return (31); case 8: return (31); case 9: return (30); case 10: return (31); case 11: return (30); case 12: return (31); default: return (0); } } static double dint (dnum) double dnum; { double dn; if (dnum < 0.0) dn = -floor (-dnum); else dn = floor (dnum); return (dn); } static double dmod (dnum, dm) double dnum, dm; { double dnumx, dnumi, dnumf; if (dnum < 0.0) dnumx = -dnum; else dnumx = dnum; dnumi = dint (dnumx / dm); if (dnum < 0.0) dnumf = dnum + (dnumi * dm); else if (dnum > 0.0) dnumf = dnum - (dnumi * dm); else dnumf = 0.0; return (dnumf); } /* Jul 1 1999 New file, based on iolib/jcon.f and iolib/vcon.f and hgetdate() * Oct 21 1999 Fix declarations after lint * Oct 27 1999 Fix bug to return epoch if fractional year input * Dec 9 1999 Fix bug in ts2jd() found by Pete Ratzlaff (SAO) * Dec 17 1999 Add all unimplemented conversions * Dec 20 1999 Add isdate(); leave date, time strings unchanged in fd2i() * Dec 20 1999 Make all fd2*() subroutines deal with time alone * * Jan 3 2000 In old FITS format, year 100 is assumed to be 2000 * Jan 11 2000 Fix epoch to date conversion so .0 is 0:00, not 12:00 * Jan 21 2000 Add separate Besselian and Julian epoch computations * Jan 28 2000 Add Modified Julian Date conversions * Mar 2 2000 Implement decimal places for FITS date string * Mar 14 2000 Fix bug in dealing with 2000-02-29 in ts2i() * Mar 22 2000 Add lt2* and ut2* to get current time as local and UT * Mar 24 2000 Fix calloc() calls * Mar 24 2000 Add tsi2* and tsu2* to convert IRAF and Unix seconds * May 1 2000 In old FITS format, all years < 1000 get 1900 added to them * Aug 1 2000 Make ep2jd and jd2ep consistently starting at 1/1 0:00 * * Jan 11 2001 Print all messages to stderr * May 21 2001 Add day of year conversions * May 25 2001 Allow fraction of day in FITS date instead of time * * Apr 8 2002 Change all long declaration to time_t * May 13 2002 Fix bugs found by lint * Jul 5 2002 Fix bug in fixdate() so fractional seconds come out * Jul 8 2002 Fix rounding bug in t2i() * Jul 8 2002 Try Fliegel and Van Flandern's algorithm for JD to UT date * Jul 8 2002 If first character of string is -, check for other -'s in isdate * Sep 10 2002 Add ET/TDT/TT conversion from UT subroutines * Sep 10 2002 Add sidereal time conversions * * Jan 30 2003 Fix typo in ts2gst() * Mar 7 2003 Add conversions for heliocentric julian dates * May 20 2003 Declare nd in setdatedec() * Jul 18 2003 Add code to parse Las Campanas dates * * Mar 24 2004 If ndec > 0, add UT to FITS date even if it is 0:00:00 * * Oct 14 2005 Add tsd2fd() and tsd2dt() * * May 3 2006 Drop declaration of unused variables * Jun 20 2006 Initialized uninitialized variables * Aug 2 2006 Add local sidereal time * Sep 13 2006 Add more local sidereal time subroutines * Oct 2 2006 Add UT to old FITS date conversions * Oct 6 2006 Add eqeqnx() to compute equation of the equinoxes * * Jan 8 2007 Remove unused variables * * Sep 5 2008 Replace nutation with IAU 2006 model translated from SOFA * Sep 9 2008 Add ang2hr(), ang2deg(), hr2ang(), deg2ang() * Sep 10 2008 Add longitude to mean standard time (default = Greenwich) * Oct 8 2008 Clean up sidereal time computations * * Sep 24 2009 Add end to comment "Coefficients for fundamental arguments" * * Jan 11 2012 Add TAI, TT, GPS time * Oct 19 2012 Unused l0 dropped from jd2lst(); ts2ss from jd2mst() */