ephcom.c

//
//   Gregorian calendar.
//
    else // Explanatory Supplement to Astronomical Almanac, p. 604
    {
        L = j + 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;
        D = L - ( 2447 * J ) / 80;
        L = J / 11;
        M = J + 2 - 12 * L;
        Y = 100 * ( N - 49 ) + I + L;
    }

    idate[0] = Y;
    idate[1] = M;
    idate[2] = D;
    idate[3] = ihour;
    idate[4] = imin;
    idate[5] = isec;
}

//! Convert calendar date and time to JD.  From pp. 604, 606 in the
//! Explanatory Supplement to the Astronomical Almanac.
//!
//! @param idate [IN ONLY] integer array of 6 values which contains
//! the integer year, month, day, hour, minute, and second.
//! @param calendar_type [IN ONLY]Integer value which controls the
//! kind of calendar used for the transformation: -1=Julian;
//! 0=automatic; 1=Gregorian.  If automatic, use Julian calendar for
//! dates before 15 October 1582.
//! @returns double-precision Julian Day number corresponding to idate.
//!
double
ephcom_cal2jd( int idate[6], int calendar_type )
{
    double tjd;
    int    jd;

//
//   Convert hours, minutes, seconds to fractional JD.
//
    tjd = ( idate[3] + ( idate[4] + idate[5] / 60.0 ) / 60.0 ) / 24.0 - 0.5;
//
//   Julian calendar.  Explanatory Supplement to Astronomical Alamanac, p. 606.
//   If automatic, use Julian calendar for dates before 15 October 1582.
//
    if ( calendar_type == -1 ||
         ( calendar_type == 0 &&
           ( idate[0] < 1582 ||                                // Before 1582
             ( idate[0] == 1582 &&
               ( idate[1] < 10 ||                              // Before October 1582
                 ( idate[1] == 10 && idate[2] < 15 ) ) ) ) ) ) // Before 15 October 1582
    {
        jd = 367 * idate[0] -
             ( 7 * ( idate[0] + 5001 + ( idate[1] - 9 ) / 7 ) ) / 4 +
             ( 275 * idate[1] ) / 9 +
             idate[2] + 1729777;
    }
//
//   Gregorian calendar.
//
    else // Explanatory Supplement to Astronomical Almanac, p. 604
    {
        jd = ( 1461 * ( idate[0] + 4800 + ( idate[1] - 14 ) / 12 ) ) / 4 +
             ( 367 * ( idate[1] - 2 - 12 * ( ( idate[1] - 14 ) / 12 ) ) ) / 12 -
             ( 3 * ( ( idate[0] + 4900 + ( idate[1] - 14 ) / 12 ) / 100 ) ) / 4 +
             idate[2] - 32075;
    }
//
//   Return value is whole JD number plus fractional JD number.
//
    tjd += (double) jd;

    return ( tjd );
}

//! If Julian date is close to an integer + 0.5, return that exact
//! "half" value.  Use of this routine makes the code more robust for
//! an early version of de422 (since corrected) which had numerical
//! noise in its half days.
//!
//! @param time [IN ONLY]Value of the Julian date which should be
//! exactly half integral.
//! @returns exactly half integral Julian date closest to time if
//! input time close to half-integral.  If input time is not close
//! to half-integral the routine returns the unmodified input time.
//!
static double
ephcom_exact_time( double time )
{
    double exact_half_time = ( time >= 0. ) ? (double) ( (int) time + 0.5 ) : (double) ( (int) time - 0.5 );
    if ( IF_SAME_DATE( time, exact_half_time ) )
        return exact_half_time;
    else
        return time;
}

//! Split time into an integer part and a positive remainder in the
//! semi-open range [0., 1.).
//!
//! @param time [IN ONLY]time value to be split.
//! @param integral_time [OUT ONLY]Pointer to a double value which
//! upon return will contain the integral part of the split time.
//! @returns remainder of the split time in the semi-open range [0., 1.).
//!

static double
ephcom_split( double time, double * integral_time )
{
    double retval = modf( time, integral_time );
    if ( retval < 0. )
    {
        retval         += 1.;
        *integral_time -= 1.;
    }
    return retval;
}