Index: /trunk/psLib/src/astro/psEarthOrientation.c =================================================================== --- /trunk/psLib/src/astro/psEarthOrientation.c (revision 5506) +++ /trunk/psLib/src/astro/psEarthOrientation.c (revision 5507) @@ -8,6 +8,6 @@ * @author Robert Daniel DeSonia, MHPCC * -* @version $Revision: 1.11 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-11-10 00:13:50 $ +* @version $Revision: 1.12 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-11-12 03:37:34 $ * * Copyright 2005 Maui High Performance Computing Center, University of Hawaii @@ -487,9 +487,137 @@ } - psEarthPole *psEOC_PrecessionCorr(const psTime *time, psTimeBulletin bulletin) { - return NULL; + PS_ASSERT_PTR_NON_NULL(time,NULL); + PS_ASSERT_INT_WITHIN_RANGE(bulletin, PS_IERS_A, PS_IERS_B, NULL); + + int tableColumn; + double mjd; + psLookupStatusType status; + char* tableName[1] = {"dailyTable"}; + double x, y, s; + // psMetadataItem *tableMetadataItem = NULL; + // psVector *temp = NULL; + // double calc = 0.0; + // double t = 0.0; + // psMetadata *timeMetadata = p_psTimeGetMetadata(); + + mjd = psTimeToMJD(time); + + if (bulletin == PS_IERS_A) { + tableColumn = 1; + } else { + tableColumn = 4; + } + + x = p_psTimeSearchTables(mjd, tableColumn, tableName, 1, &status); + if ( status == PS_LOOKUP_ERROR ) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, + "Time table search returned unsuccessful status.\n"); + return NULL; + } else if (status != PS_LOOKUP_SUCCESS) { + psSphere *sphere = p_psTimeGetPoleCoords(time); + x = sphere->r; + psFree(sphere); + } else { + x = SEC_TO_RAD(x); + } + + + /* }else if(status == PS_LOOKUP_PAST_TOP) { + // Date too early for tables. Get default time delta value from metadata, and issue warning. + psLogMsg(__func__,PS_LOG_WARN,PS_ERRORTEXT_psTime_TIME_PREDATES_TABLES,mjd,"xp"); + + // Lookup value from time metadata loaded from pslib.config + tableMetadataItem = psMetadataLookup(timeMetadata, "psLib.time.before.xp"); + if(tableMetadataItem == NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, PS_ERRORTEXT_psTime_LOOKUP_METADATA_FAILED, + "psLib.time.before.xp"); + return NULL; + } + x = tableMetadataItem->data.F64; + + } else if(status == PS_LOOKUP_PAST_BOTTOM) { + tableMetadataItem = psMetadataLookup(timeMetadata, "psLib.time.predict.mjd"); + if(tableMetadataItem == NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, PS_ERRORTEXT_psTime_LOOKUP_METADATA_FAILED, + "psLib.time.predict.mjd"); + return NULL; + } + mjdPred = tableMetadataItem->data.F64; + + // Get xp + tableMetadataItem = psMetadataLookup(timeMetadata, "psLib.time.predict.xp"); + if(tableMetadataItem == NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, PS_ERRORTEXT_psTime_LOOKUP_METADATA_FAILED, "psLib.time.predict.xp"); + return NULL; + } + xp = (psVector*)tableMetadataItem->data.V; + PS_ASSERT_PTR_NON_NULL(xp,NULL); + + // Get yp + tableMetadataItem = psMetadataLookup(timeMetadata, "psLib.time.predict.yp"); + if(tableMetadataItem == NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, PS_ERRORTEXT_psTime_LOOKUP_METADATA_FAILED, "psLib.time.predict.yp"); + return NULL; + } + yp = (psVector*)tableMetadataItem->data.V; + PS_ASSERT_PTR_NON_NULL(yp,NULL); + + // Calculate "a" and "c" constants + a = TWOPI*(mjd - mjdPred)/365.25; + c = TWOPI*(mjd - mjdPred)/435.0; + + // Calculate x and y polar coordinates + x = xp->data.F64[0] + + xp->data.F64[1]*cos(a) + + xp->data.F64[2]*sin(a) + + xp->data.F64[3]*cos(c) + + xp->data.F64[4]*sin(c); + + + + + + } + */ + + + tableColumn++; + y = p_psTimeSearchTables(mjd, tableColumn, tableName, 1, &status); + if ( status == PS_LOOKUP_ERROR ) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, + "Time table search returned unsuccessful status.\n"); + return NULL; + } else if (status != PS_LOOKUP_SUCCESS) { + psSphere *sphere = p_psTimeGetPoleCoords(time); + y = sphere->d; + psFree(sphere); + } else { + y = SEC_TO_RAD(y); + } + + tableColumn++; + s = p_psTimeSearchTables(mjd, tableColumn, tableName, 1, &status); + if ( status == PS_LOOKUP_ERROR ) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, + "Time table search returned unsuccessful status.\n"); + return NULL; + } else if (status != PS_LOOKUP_SUCCESS) { + s = psTimeGetUT1Delta(time, bulletin); + } else { + s = SEC_TO_RAD(s); + } + + //XXX: Do I need to do something for nutation terms w/periods of less than 2 days??? + psEarthPole *out = psEarthPoleAlloc(); + out->x = x; + out->y = y; + out->s = s; + + // psFree(tableMetadataItem); + // psFree(timeMetadata); + return out; } @@ -566,5 +694,4 @@ } - psSphereRot* psSphereRot_TEOtoCEO(const psTime *time) { @@ -589,6 +716,34 @@ psTimeBulletin bulletin) { - - return NULL; + PS_ASSERT_PTR_NON_NULL(time, NULL); + + psEarthPole *out = psEarthPoleAlloc(); + psSphere *in = NULL; + double s; + + //XXX: This may be the wrong idea... ADD says to use 3rd-order polys to interpolate + //polar motion coordinates. + in = p_psTimeGetPoleCoords(time); + if (in == NULL) { + psError(PS_ERR_BAD_PARAMETER_NULL, true, + "p_psTimeGetPoleCoords return NULL psSphere for non-NULL input time.\n"); + return NULL; + } + out->x = in->r; + out->y = in->d; + + s = psTimeGetUT1Delta(time, bulletin); + out->s = s; + + //XXX: Apply polar tide correction here??? + psEarthPole *correction = psEOC_PolarTideCorr(time); + out->x += correction->x; + out->y += correction->y; + out->s += correction->s; + + + psFree(in); + psFree(correction); + return out; } @@ -711,5 +866,5 @@ diag_sum[2] = 1.0 - A[0][0] - A[1][1] + A[2][2]; diag_sum[3] = 1.0 + A[0][0] + A[1][1] + A[2][2]; - + maxi = 0; for (int i = 1; i < 4; ++i) { @@ -718,8 +873,8 @@ } } - + double p = 0.5 * sqrt(diag_sum[maxi]); recip = 1.0 / (4.0 * p); - + if (maxi == 0) { out->q0 = p; Index: /trunk/psLib/src/astro/psTime.c =================================================================== --- /trunk/psLib/src/astro/psTime.c (revision 5506) +++ /trunk/psLib/src/astro/psTime.c (revision 5507) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.74 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-05 03:51:41 $ + * @version $Revision: 1.75 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-12 03:37:34 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii Index: /trunk/psLib/src/astro/psTime.h =================================================================== --- /trunk/psLib/src/astro/psTime.h (revision 5506) +++ /trunk/psLib/src/astro/psTime.h (revision 5507) @@ -11,6 +11,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.42 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-11-03 04:31:22 $ + * @version $Revision: 1.43 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-12 03:37:34 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -375,15 +375,19 @@ char* p_psGetConfigFileName(); -/** - * - * - * +/** Searches the IERS time tables for a specified entry location. + * + * Returns the interpolated double precision (arcsec) value at the specified entry + * location. Inputs to specify are the time index in mjd, the column number + * corresponding to Xp, Yp, or Sp (UT1-UTC) in IERS A or B, the time table names, + * and the number of time tables. + * + * @return psF64: Resulting table entry for specified parameters. */ psF64 p_psTimeSearchTables( - psF64 index, ///< - psU64 column, ///< - char *metadataTableNames[], ///< - psU32 nTables, ///< - psLookupStatusType* status ///< + psF64 index, ///< time index for which to search + psU64 column, ///< column number of specified index + char *metadataTableNames[], ///< names of IERS tables to search + psU32 nTables, ///< number of IERS tables to search + psLookupStatusType* status ///< status of table search ); Index: /trunk/psLib/test/astro/tst_psEarthOrientation.c =================================================================== --- /trunk/psLib/test/astro/tst_psEarthOrientation.c (revision 5506) +++ /trunk/psLib/test/astro/tst_psEarthOrientation.c (revision 5507) @@ -5,6 +5,6 @@ * @author d-Rob, MHPCC * -* @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-11-10 00:13:51 $ +* @version $Revision: 1.10 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-11-12 03:37:34 $ * * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii @@ -139,5 +139,55 @@ psS32 testEOCPrecession(void) { - + psTime *empty = NULL; + psTime *time = psTimeAlloc(PS_TIME_TAI); + time->sec = 1131579114; + time->nsec = 498489000; + time->leapsecond = false; + + + //Tests for Precession Correction function// + //Return NULL for NULL time input + psEarthPole *pcorr = NULL; + pcorr = psEOC_PrecessionCorr(empty, PS_IERS_A); + if (pcorr != NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, + "psEOC_PrecessionCorr failed to return NULL for NULL time input.\n"); + return 4; + } + + //Return NULL for Invalid IERS table + pcorr = psEOC_PrecessionCorr(time, 3); + if (pcorr != NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, + "psEOC_PrecessionCorr failed to return NULL for incorrect IERS table.\n"); + return 5; + } + psFree(pcorr); + + //Check values from IERS table A + pcorr = psEOC_PrecessionCorr(time, PS_IERS_A); + if ( pcorr == NULL ) { + psError(PS_ERR_BAD_PARAMETER_NULL, false, + "psEOC_PrecessionCorr returned NULL for valid inputs.\n"); + return 6; + } else { + printf("Precession Correction output (IERSA) = x,y,s = %.8g, %.8g, %.8g\n", + pcorr->x, pcorr->y, pcorr->s); + } + psFree(pcorr); + + //Check values from IERS table B + pcorr = psEOC_PrecessionCorr(time, PS_IERS_B); + if ( pcorr == NULL ) { + psError(PS_ERR_BAD_PARAMETER_NULL, false, + "psEOC_PrecessionCorr returned NULL for valid inputs.\n"); + return 7; + } else { + printf("Precession Correction output (IERSB) = x,y,s = %.8g, %.8g, %.8g\n", + pcorr->x, pcorr->y, pcorr->s); + } + psFree(pcorr); + + psFree(time); return 0; } @@ -167,7 +217,27 @@ } + psEarthPole *polarMotion = NULL; + polarMotion = psEOC_GetPolarMotion(empty, PS_IERS_B); + if (polarMotion != NULL) { + psError(PS_ERR_BAD_PARAMETER_VALUE, false, + "psEOC_GetPolarMotion failed to return NULL for NULL input time.\n"); + return 3; + } + polarMotion = psEOC_GetPolarMotion(in, PS_IERS_B); + if (polarMotion == NULL) { + psError(PS_ERR_BAD_PARAMETER_NULL, false, + "psEOC_GetPolarMotion returned NULL for valid input time.\n"); + return 4; + } + + psSphere *interm = p_psTimeGetPoleCoords(in); + printf(" -- sphere = r=%.8g, d=%.8g\n", interm->r, interm->d); printf(" -- Eop = x=%.8g, y=%.8g, s=%.8g\n", eop->x, eop->y, eop->s); + printf(" -- PolarMotion = x=%.8g, y=%.8g, s=%.8g\n", polarMotion->x, + polarMotion->y, polarMotion->s); psFree(in); psFree(eop); + psFree(polarMotion); + psFree(interm); return 0; }