Index: /trunk/psLib/src/astro/psEarthOrientation.c =================================================================== --- /trunk/psLib/src/astro/psEarthOrientation.c (revision 6035) +++ /trunk/psLib/src/astro/psEarthOrientation.c (revision 6036) @@ -8,6 +8,6 @@ * @author Robert Daniel DeSonia, MHPCC * - * @version $Revision: 1.27 $ $Name: not supported by cvs2svn $ - * @date $Date: 2006-01-18 00:41:29 $ + * @version $Revision: 1.28 $ $Name: not supported by cvs2svn $ + * @date $Date: 2006-01-18 20:59:31 $ * * Copyright 2005 Maui High Performance Computing Center, University of Hawaii @@ -846,7 +846,7 @@ if (tidalCorr != NULL && tidalCorr->s != 0.0) { int nsec = in->nsec + (int)(tidalCorr->s * 1e9); - if (nsec < 0.0) { + if (nsec < 0) { in->sec += -1; - in->nsec = (int)(1e9) - nsec; + in->nsec = (int)(1e9) + nsec; } else { in->nsec = nsec; @@ -1009,6 +1009,4 @@ // Calculate number of Julian centuries since 2000 - //XXX: NOT SURE IF THIS IS CORRECT FOR THIS SITUATION - // double RJD = ( MJD - MJD_2000 ) / JULIAN_CENTURY; double RJD = MJD; Index: /trunk/psLib/src/astro/psTime.c =================================================================== --- /trunk/psLib/src/astro/psTime.c (revision 6035) +++ /trunk/psLib/src/astro/psTime.c (revision 6036) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.76 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-12-05 22:00:48 $ + * @version $Revision: 1.77 $ $Name: not supported by cvs2svn $ + * @date $Date: 2006-01-18 20:59:31 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -96,5 +96,6 @@ /** Removes leading and trailing whitespace and # characters from a string. The cleaned string is a new null * terminated copy of the original input string. */ -static char *cleanString(char *inString, int sLen) +static char *cleanString(char *inString, + int sLen) { char *ptrB = NULL; @@ -126,5 +127,7 @@ /** Returns cleaned token based on delimiter, but not including delimiter. Also changes the pointer location * the beginning of the string. Tokens are newly allocated null terminated strings. */ -static char* getToken(char **inString, char *delimiter, psParseErrorType *status) +static char* getToken(char **inString, + char *delimiter, + psParseErrorType *status) { char *cleanToken = NULL; @@ -169,6 +172,9 @@ // Searches time tables in priority order and performs interpolation if input index value is within a table. // If the index value is out of range, the status is set accordingly. -psF64 p_psTimeSearchTables(psF64 index, psU64 column, char *metadataTableNames[], - psU32 nTables, psLookupStatusType* status) +psF64 p_psTimeSearchTables(psF64 index, + psU64 column, + char *metadataTableNames[], + psU32 nTables, + psLookupStatusType* status) { char* tableName = NULL; @@ -660,5 +666,6 @@ } -psTime* psTimeConvert(psTime *time, psTimeType type) +psTime* psTimeConvert(psTime *time, + psTimeType type) { // Error checks @@ -745,5 +752,6 @@ } -double psTimeToLMST(psTime *time, double longitude) +double psTimeToLMST(psTime *time, + double longitude) { psF64 jdTdtDays = 0.0; @@ -818,5 +826,6 @@ } -double psTimeGetUT1Delta(const psTime *time, psTimeBulletin bulletin) +double psTimeGetUT1Delta(const psTime *time, + psTimeBulletin bulletin) { psU32 nTables = 2; @@ -899,4 +908,99 @@ return result; +} + +static double DMOD(double x, double y) +{ + double value = x - y * trunc(x/y); + return value; +} + +psTime *psTime_TideUT1Corr(const psTime *time) +{ + PS_ASSERT_PTR_NON_NULL(time, NULL); + psTime *out = NULL; + + // Convert psTime to MJD + double MJD = psTimeToMJD(time); + if (MJD == NAN) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, + "Time conversion to MJD failed. Invalid input time.\n"); + return NULL; + } + + // Calculate number of Julian centuries since 2000 + double RJD = MJD; + + //Formula comes from fortran reference + //DMOD in fortran ref. = double remainder -> x - y * trunc(x/y) + double T, L, LPRIME, CAPF, CAPD, OMEGA, THETA, CORZ; + double ARG1, ARG2, ARG3, ARG4, ARG5, ARG6, ARG7, ARG8; + double T2, T3, T4; + T = (RJD - 51544.5) / 36525.0; + T2 = T*T; + T3 = T*T*T; + T4 = T*T*T*T; + L = -0.0002447 * T4 + 0.051635 * T3 + 31.8792 * T2 + 1717915923.2178 * T + 485868.249036; + L = DMOD(L, 1296000.0); + LPRIME = -0.00001149 * T4 - 0.000136 * T3 - 0.5532 * T2 + 129596581.0481 * T + 1287104.79305; + LPRIME = DMOD(LPRIME, 1296000.0); + CAPF = 0.00000417 * T4 - 0.001037 * T3 - 12.7512 * T2 + 1739527262.8478 * T + 335779.526232; + CAPF = DMOD(CAPF, 1296000.0); + CAPD = -0.00003169 * T4 + 0.006593 * T3 - 6.3706 * T2 + 1602961601.209 * T + 1072260.70369; + CAPD = DMOD(CAPD, 1296000.0); + OMEGA = -0.00005939 * T4 + 0.007702 * T3 + 7.4722 * T2 - 6962890.2665 * T + 450160.398036; + OMEGA = DMOD(OMEGA, 1296000.0); + THETA = (67310.54841 + (876600.0 * 3600.0 + 8640184.812866) * T + 0.093104 * T2 - + 6.2e-6 * T3) * 15.0 + 648000.0; + ARG7 = DMOD((-L - 2.0 * CAPF - 2.0 * OMEGA + THETA) * M_PI / 648000.0, 2.0 * M_PI) + - M_PI / 2.0; + ARG1 = DMOD((-2.0 * CAPF - 2.0 * OMEGA + THETA) * M_PI / 648000.0, 2.0 * M_PI) - M_PI / 2.0; + ARG2 = DMOD((-2.0 * CAPF + 2.0 * CAPD - 2.0 * OMEGA + THETA) * M_PI / 648000.0, 2.0 * M_PI) + - M_PI / 2.0; + ARG3 = DMOD(THETA * M_PI / 648000.0, 2.0 * M_PI) - M_PI / 2.0; + ARG4 = DMOD((-L - 2.0 * CAPF - 2.0 * OMEGA + 2.0 * THETA) * M_PI / 648000.0, 2.0 * M_PI); + ARG5 = DMOD((-2.0 * CAPF - 2.0 * OMEGA + 2.0 * THETA) * M_PI / 648000.0, 2.0 * M_PI); + ARG6 = DMOD((-2.0 * CAPF + 2.0 * CAPD - 2.0 * OMEGA + 2.0 * THETA) * M_PI / 648000.0, + 2.0 * M_PI); + ARG8 = DMOD((2.0 * THETA) * M_PI / 648000.0, 2.0 * M_PI); + CORZ = 0.0245 * sin(ARG7) + 0.0503 * cos(ARG7) + +0.1210 * sin(ARG1) + 0.1605 * cos(ARG1) + +0.0286 * sin(ARG2) + 0.0516 * cos(ARG2) + +0.0864 * sin(ARG3) + 0.1771 * cos(ARG3) + -0.0380 * sin(ARG4) - 0.0154 * cos(ARG4) + -0.1617 * sin(ARG5) - 0.0720 * cos(ARG5) + -0.0759 * sin(ARG6) - 0.0004 * cos(ARG6) + -0.0196 * sin(ARG8) - 0.0038 * cos(ARG8); + CORZ = CORZ * 0.1e-3; + + double timeCheck = (double)(time->sec) + (double)(1e-9*time->nsec); + if ( (timeCheck + CORZ) < 0.0 ) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, + "Invalid time for Tide Correction.\n"); + return NULL; + } + out = psTimeAlloc(time->type); + *out = *time; + // out->sec = time->sec; + // out->nsec = time->nsec; + // out->leapsecond = time->leapsecond; + if (out->type != PS_TIME_UT1) { + out = psTimeConvert(out, PS_TIME_UT1); + } + if (fabs(CORZ) > 1.0) { + int sec = (int)CORZ; + out->sec += sec; + int nsec = (int)((CORZ - sec)*1e9); + out->nsec += nsec; + } else { + int nsec = out->nsec + (int)(CORZ * 1e9); + if (nsec < 0) { + out->sec += -1; + out->nsec = (int)(1e9) + nsec; + } else { + out->nsec = nsec; + } + } + return out; } @@ -1090,5 +1194,6 @@ } -long psTimeLeapSecondDelta(const psTime *time1, const psTime *time2) +long psTimeLeapSecondDelta(const psTime *time1, + const psTime *time2) { psS64 diff = 0; @@ -1397,5 +1502,6 @@ } -psTime* psTimeFromTT(psS64 sec, psU32 nsec) +psTime* psTimeFromTT(psS64 sec, + psU32 nsec) { psTime* outTime = NULL; @@ -1415,5 +1521,7 @@ } -psTime* psTimeFromUTC(psS64 sec, psU32 nsec, bool leapsecond) +psTime* psTimeFromUTC(psS64 sec, + psU32 nsec, + bool leapsecond) { psTime* outTime = NULL; @@ -1519,5 +1627,6 @@ } -psTime* psTimeMath(const psTime *time, double delta) +psTime* psTimeMath(const psTime *time, + double delta) { psF64 sec = 0.0; @@ -1560,5 +1669,6 @@ } -double psTimeDelta(const psTime *time1, const psTime *time2) +double psTimeDelta(const psTime *time1, + const psTime *time2) { psF64 out = 0.0; Index: /trunk/psLib/src/astro/psTime.h =================================================================== --- /trunk/psLib/src/astro/psTime.h (revision 6035) +++ /trunk/psLib/src/astro/psTime.h (revision 6036) @@ -11,6 +11,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.43 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-11-12 03:37:34 $ + * @version $Revision: 1.44 $ $Name: not supported by cvs2svn $ + * @date $Date: 2006-01-18 20:59:31 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -76,5 +76,5 @@ * Reads config and data files associated with various time conversions. * - * @return bool: True for success, false for failure. + * @return bool: True for success, false for failure. */ psBool p_psTimeInit( @@ -94,5 +94,5 @@ * Frees time data to be held in memory until the end of successful program execution. * - * @return void: void. + * @return void: void. */ psBool p_psTimeFinalize(void); @@ -110,5 +110,5 @@ * of the struct are set to zero. * - * @return psTime*: Struct with empty time. + * @return psTime*: Struct with empty time. */ psTime* psTimeAlloc( @@ -133,5 +133,5 @@ * (PS_TIME_TAI or PS_TIME_UTC) in the argument. * - * @return psTime*: Struct with current time. + * @return psTime*: Struct with current time. */ psTime* psTimeGetNow( @@ -143,5 +143,5 @@ * Converts psTime to UTC, TAI, UT1, or TT time based on the psTimeType argument. * - * @return psTime*: Pointer to psTime. + * @return psTime*: Pointer to psTime. */ psTime* psTimeConvert( @@ -155,5 +155,5 @@ * in UTC format, then it is converted. * - * @return double: LST Time. + * @return double: LST Time. */ double psTimeToLMST( @@ -166,5 +166,5 @@ * This function is necessary to for various SLALIB functions. * - * @return double: Time difference. + * @return double: Time difference. */ double psTimeGetUT1Delta( @@ -173,9 +173,19 @@ ); +/** Provides tidal corrections to UT1-UTC. + * + * Uses the Ray model of Simon et al. + * + * @return psTime*: The corrected time in UT1. + */ +psTime *psTime_TideUT1Corr( + const psTime *time ///< psTime to be corrected. +); + /** Determine TAI - UTC from table lookup. * * This function is necessary to for various psTime functions. * - * @return psF64: Time difference. + * @return psF64: Time difference. */ psF64 p_psTimeGetTAIDelta( @@ -187,5 +197,5 @@ * Determines the orientation of the polar axis at the given time. * - * @return psSphere*: Spherical coordinates of Earth's polar axias. + * @return psSphere*: Spherical coordinates of Earth's polar axias. */ psSphere* p_psTimeGetPoleCoords( @@ -197,5 +207,5 @@ * Calculates the number of leapseconds between two times. * - * @return long: leapseconds added between given times + * @return long: leapseconds added between given times */ long psTimeLeapSecondDelta( @@ -208,5 +218,5 @@ * Determines if the specified UTC time is a valid leapsecond. * - * @return bool: valid leap second + * @return bool: valid leap second */ bool psTimeIsLeapSecond( @@ -219,5 +229,5 @@ * subtract leapseconds. * - * @return double: Julian Date (JD) time. + * @return double: Julian Date (JD) time. */ double psTimeToJD( @@ -229,5 +239,5 @@ * add or subtract leapseconds. * - * @return double: Modified Julian Days (MJD) time. + * @return double: Modified Julian Days (MJD) time. */ double psTimeToMJD( @@ -240,5 +250,5 @@ * This function does not add or subtract leapseconds. * - * @return psString: Pointer null terminated array of chars in ISO time. + * @return psString: Pointer null terminated array of chars in ISO time. */ psString psTimeToISO( @@ -250,5 +260,5 @@ * Converts psTime to timeval time. This function does not add or subtract leapseconds. * - * @return timeval*: timeval struct time. + * @return timeval*: timeval struct time. */ struct timeval* psTimeToTimeval( Index: /trunk/psLib/test/astro/tst_psTime_04.c =================================================================== --- /trunk/psLib/test/astro/tst_psTime_04.c (revision 6035) +++ /trunk/psLib/test/astro/tst_psTime_04.c (revision 6036) @@ -11,10 +11,13 @@ * Test F - Attempt to read data file with typo in number * Test G - Free data + * Test H - Attempt to use all timer functions + * Test I - Tidal Corrections to UT1-UTC * * @author Ross Harman, MHPCC * @author Eric Van Alst, MHPCC + * @author David Robbins, MHPCC * - * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-01 02:22:20 $ + * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ + * @date $Date: 2006-01-18 20:59:32 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -33,4 +36,5 @@ static psS32 testTimeInit6(void); static psS32 testTimer1(void); +static psS32 testTideUT1Corr(void); testDescription tests[] = { @@ -42,4 +46,5 @@ {testTimeInit6,-6,"p_psTimeInit",0,false}, {testTimer1,-7,"psTimer",0,false}, + {testTideUT1Corr,-8,"psTime_TideUT1Corr",0,false}, {NULL} }; @@ -136,2 +141,36 @@ } +psS32 testTideUT1Corr(void) +{ + psTime *empty = NULL; + psTime *tide = NULL; + psTime *noTide = psTimeAlloc(PS_TIME_UTC); + noTide->sec = 1049160600; + noTide->nsec = 0; + noTide->leapsecond = false; + + empty = psTime_TideUT1Corr(tide); + if (empty != NULL) { + psError(PS_ERR_BAD_PARAMETER_NULL, false, + "psTime_TideUT1Corr failed to return NULL for NULL input time.\n"); + return 1; + } + + empty = psTime_TideUT1Corr(noTide); + if (empty->sec != 1049160599 || empty->nsec != 656982272) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, + "psTime_TideUT1Corr failed to return correct values.\n"); + printf("\nsec = %ld, nsec = %u\n", empty->sec, empty->nsec); + return 2; + } + + if (!p_psTimeFinalize()) { + return 3; + } + + psFree(empty); + psFree(noTide); + + return 0; +} +