Index: trunk/psLib/src/astronomy/psTime.c =================================================================== --- trunk/psLib/src/astronomy/psTime.c (revision 2032) +++ trunk/psLib/src/astronomy/psTime.c (revision 2048) @@ -11,6 +11,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.30 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-10-09 00:17:14 $ + * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-12 01:34:09 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -27,6 +27,15 @@ #include "psAbort.h" #include "psImage.h" +#include "psCoord.h" #include "psAstronomyErrors.h" + +#ifndef SER7_FILE +#define SER7_FILE "../../data/ser7.dat" +#endif + +#ifndef TAIUTC_FILE +#define TAIUTC_FILE "../../data/tai-utc.dat" +#endif /** Sidereal angular conversion from seconds to radians for GMST in seconds (i.e. pi/(180*240)) */ @@ -61,22 +70,88 @@ /** Preprocessor macro to check for null psTime struct */ - #define CHECK_NULL_TIME(TIME,RETURN) \ - if(TIME == NULL) { \ - psError(__func__,"NULL value not allowed"); \ - return RETURN; \ + #define CHECK_NULL_TIME(TIME,RETURN) \ + if(TIME == NULL) { \ + psError(__func__,"NULL value not allowed"); \ + return RETURN; \ } /** Preprocessor macro to allocate psTime struct */ - #define ALLOC_TIME(TIME) \ - TIME = (psTime*)psAlloc(sizeof(psTime)); \ -if (TIME == NULL) { \ - psAbort(__func__, " : Line %d - Failed to allocate memory", __LINE__); \ -} \ -TIME->sec = 0; \ -TIME->usec = 0; \ + #define ALLOC_TIME(TIME) \ + TIME = (psTime*)psAlloc(sizeof(psTime)); \ +if (TIME == NULL) { \ + psAbort(__func__, " : Line %d - Failed to allocate memory", __LINE__); \ +} \ +TIME->sec = 0; \ +TIME->usec = 0; \ TIME->type = PS_TIME_TAI; - -psImage* readTaiUtcFile(char *fileName) +/** Read TAI/UTC data file. + * + * Reads TAI/UTC data file. + * + * @return psImage*: Two-dimensional array with table data. + */ +static psImage* readTaiUtcFile( + const char *fileName ///< Name of file to read. +); + + +/** Read SER7 data file. + * + * Reads SER7 data file. + * + * @return psImage*: Two-dimensional array with table data. + */ +static psImage* readSer7File( + const char *fileName ///< Name of file to read. +); + + +/** Read finals data file. + * + * Reads finals data file. + * + * @return psImage*: Two-dimensional array with table data. + */ +/* not yet used +static psImage* readFinalsFile( + char *fileName ///< Name of file to read. +); +*/ + +/** Read EOPC data file. + * + * Reads EOPC data file. + * + * @return psImage*: Two-dimensional array with table data. + */ +/* not yet used +static psImage* readEopcFile( + char *fileName ///< Name of file to read. +); +*/ + +/** Lookup and interpolate TAI-UTC data + * + * Interpolates TAI-UTC data, given time in MJD format. + * + * @return double: Interpolated value. + */ +static double lookupTaiUtcTable( + const psTime *time ///< time to lookup. +); + +/** Lookup and interpolate polar or UT1-UTC data + * + * Interpolates TAI-UT1 or polar coordinate data, given time in MJD format. + * + * @return double: Interpolated value. + */ +static double lookupSer7Table( + const psTime *time, ///< time to lookup. + int col ///< Column to lookup. +); + +static psImage* readTaiUtcFile(const char *fileName) { char line[LINESIZE]; @@ -113,8 +188,12 @@ *(unsigned int *)&table->numRows = j; + p_psMemSetPersistent(table,true); + p_psMemSetPersistent(table->data.V,true); + p_psMemSetPersistent(table->rawDataBuffer,true); + return table; } -psImage* readSer7File(char *fileName) +static psImage* readSer7File(const char *fileName) { bool beginRecord = false; @@ -160,8 +239,13 @@ *(unsigned int *)&table->numRows = j; + p_psMemSetPersistent(table,true); + p_psMemSetPersistent(table->data.V,true); + p_psMemSetPersistent(table->rawDataBuffer,true); + return table; } -psImage* readEopcFile(char *fileName) +/* Not yet used... +static psImage* readEopcFile(char *fileName) { char line[LINESIZE]; @@ -171,6 +255,6 @@ psImage *table = NULL; FILE *fd = NULL; - - + + fd=fopen(fileName, "r"); if(fd==NULL) { @@ -178,5 +262,5 @@ return NULL ; } - + // Table shouldn't be larger than 2500 rows. All columns are read. table = psImageAlloc(11, maxLines, PS_TYPE_F64); @@ -186,5 +270,5 @@ } } - + while(fgets(line, LINESIZE, fd) != NULL) { if(j>=maxLines) { @@ -197,11 +281,11 @@ ptr, ptr+1, ptr+2, ptr+3, ptr+4, ptr+5, ptr+6, ptr+7, ptr+8, ptr+9, ptr+10); } - + *(unsigned int *)&table->numRows = j; - + return table; } - -psImage* readFinalsFile(char *fileName) + +static psImage* readFinalsFile(char *fileName) { char line[LINESIZE]; @@ -210,6 +294,6 @@ psImage *table = NULL; FILE *fd = NULL; - - + + fd=fopen(fileName, "r"); if(fd==NULL) { @@ -217,5 +301,5 @@ return NULL ; } - + // Table shouldn't be larger than 15000 rows. Select columns are read. table = psImageAlloc(6, maxLines, PS_TYPE_F64); @@ -226,8 +310,8 @@ psFree(table); } - + psF64 *ptr = table->data.F64[j++]; ptr = table->data.F64[j++]; - + // Sometimes there are blank entries in this table, so check for a space to determine if there is no entry if(*(line+7) == ' ') { // MJD @@ -236,5 +320,5 @@ *ptr = atof(line+7); } - + if(*(line+18) == ' ') { // Polar motion X *(ptr+1) = 0.0; @@ -242,5 +326,5 @@ *(ptr+1) = atof(line+18); } - + if(*(line+27) == ' ') { // Polar motion X error *(ptr+2) = 0.0; @@ -248,5 +332,5 @@ *(ptr+2) = atof(line+27); } - + if(*(line+37) == ' ') { // Polar motion Y *(ptr+3) = 0.0; @@ -254,5 +338,5 @@ *(ptr+3) = atof(line+37); } - + if(*(line+46) == ' ') { // Polar motion Y error *(ptr+4) = 0.0; @@ -260,5 +344,5 @@ *(ptr+4) = atof(line+46); } - + if(*(line+58) == ' ') { // UT1-UTC *(ptr+5) = 0.0; @@ -267,12 +351,12 @@ } } - + *(unsigned int *)&table->numRows = j; - + return table; } - - -double lookupTaiUtcTable(psImage *table, psTime *time) +*/ + +static double lookupTaiUtcTable(const psTime *time) { int hiIdx = 0; @@ -286,5 +370,12 @@ double const2 = 0.0; double const3 = 0.0; - + static psImage* table = NULL; + + if (table == NULL) { + table = readTaiUtcFile(TAIUTC_FILE); + if (table == NULL) { + psAbort(__func__,"Failed to read tai-utc.dat file"); + } + } // Number of rows is the number of rows of data read from the data file @@ -337,5 +428,5 @@ } -double lookupSer7Table(psImage *table, psTime *time, int col) +static double lookupSer7Table(const psTime *time, int col) { int hiIdx = 0; @@ -345,5 +436,12 @@ double denom = 0.0; double mjdUtc = 0.0; - + static psImage* table = NULL; + + if (table == NULL) { + table = readSer7File(SER7_FILE); + if (table == NULL) { + psAbort(__func__,"Failed to read the ser7.dat file."); + } + } // Number of rows is the number of rows of data read from the data file @@ -404,5 +502,5 @@ // Add most leapseconds to UTC time to get TAI time if necessary if(type == PS_TIME_TAI) { - time->sec += psGetTAIDelta(time); + time->sec += psTimeGetTAIDelta(time); } @@ -420,5 +518,5 @@ } - delta = psGetTAIDelta(time); + delta = psTimeGetTAIDelta(time); if(type == PS_TIME_UTC) { @@ -461,5 +559,5 @@ taiTime = psMemIncrRefCounter(time); ALLOC_TIME(utcTime); - utcTime->sec = taiTime->sec - psGetTAIDelta(time); + utcTime->sec = taiTime->sec - psTimeGetTAIDelta(time); utcTime->usec = taiTime->usec; utcTime->type = PS_TIME_UTC; @@ -467,5 +565,5 @@ utcTime = psMemIncrRefCounter(time); ALLOC_TIME(taiTime); - taiTime->sec = utcTime->sec + psGetTAIDelta(time); + taiTime->sec = utcTime->sec + psTimeGetTAIDelta(time); taiTime->usec = utcTime->usec; taiTime->type = PS_TIME_TAI; @@ -475,5 +573,5 @@ ALLOC_TIME(ut1UtcDelta); ut1UtcDelta->type = PS_TIME_UTC; - ut1UtcDelta->usec = psGetUT1Delta(utcTime)*1e6; + ut1UtcDelta->usec = psTimeGetUT1Delta(utcTime)*1e6; ut1Time = psTimeAdd(utcTime, ut1UtcDelta); @@ -516,42 +614,30 @@ } -double psGetUT1Delta(psTime *time) -{ - static psImage* iersTable = NULL; - double ut1UtcDeltaUsec = 0.0; - - - if (iersTable == NULL) { - iersTable = readSer7File("../../data/ser7.dat"); - p_psMemSetPersistent(iersTable,true); - p_psMemSetPersistent(iersTable->data.V,true); - p_psMemSetPersistent(iersTable->rawDataBuffer,true); - } - - ut1UtcDeltaUsec = lookupSer7Table(iersTable, time, 6)*1.0e6; - - return ut1UtcDeltaUsec; -} - -double psGetTAIDelta(psTime *time) -{ - double delta = 0.0; - static psImage* taiUtcTable = NULL; - - - if (taiUtcTable == NULL) { - taiUtcTable = readTaiUtcFile("../../data/tai-utc.dat"); - p_psMemSetPersistent(taiUtcTable,true); - p_psMemSetPersistent(taiUtcTable->data.V,true); - p_psMemSetPersistent(taiUtcTable->rawDataBuffer,true); - } - delta = lookupTaiUtcTable(taiUtcTable, time); - - return delta; -} - -psF64 psTimeLeapseconds(const psTime *time1, const psTime *time2) -{ - psF64 diff = 0.0; +double psTimeGetUT1Delta(const psTime *time) +{ + return lookupSer7Table(time, 6); +} + +struct psSphere* psTimeGetPoleCoords(const psTime* time) +{ + + double x = lookupSer7Table((psTime*)time, 4); + double y = lookupSer7Table((psTime*)time, 5); + + struct psSphere* output = psAlloc(sizeof(psSphere)); + output->r = x * M_PI / 648000.0; // x in arcsec converted to radians, i.e., x/60/60*M_PI/180 + output->d = y * M_PI / 648000.0; + + return output; +} + +double psTimeGetTAIDelta(const psTime *time) +{ + return lookupTaiUtcTable(time); +} + +long psTimeLeapseconds(const psTime *time1, const psTime *time2) +{ + long diff = 0; // NULL error checks @@ -559,10 +645,10 @@ CHECK_NULL_TIME(time2,0); - diff = fabs(psGetTAIDelta((psTime*)time1)-psGetTAIDelta((psTime*)time2)); + diff = abs(psTimeGetTAIDelta((psTime*)time1)-psTimeGetTAIDelta((psTime*)time2)); return diff; } -double psTimeToJD(psTime *time) +double psTimeToJD(const psTime *time) { double jd = 0.0; @@ -582,5 +668,5 @@ } -double psTimeToMJD(psTime *time) +double psTimeToMJD(const psTime *time) { double mjd = 0.0; @@ -600,5 +686,5 @@ } -char* psTimeToISO(psTime *time) +char* psTimeToISOTime(const psTime *time) { int ms = 0; @@ -634,5 +720,5 @@ } -struct timeval psTimeToTimeval(psTime *time) +struct timeval psTimeToTimeval(const psTime *time) { struct timeval timevalTime; @@ -648,5 +734,5 @@ } -struct tm* psTimeToTM(psTime *time) +struct tm* psTimeToTM(const psTime *time) { long cent = 0; @@ -711,5 +797,5 @@ } -psTime* psJDToTime(double time) +psTime* psTimeFromJD(double time) { double days = 0.0; @@ -734,5 +820,5 @@ } -psTime* psMJDToTime(double time) +psTime* psTimeFromMJD(double time) { double days = 0.0; @@ -757,5 +843,5 @@ } -psTime* psISOToTime(char *time) +psTime* psTimeFromISOTime(const char *time) { char tempString[MAX_TIME_STRING_LENGTH]; @@ -826,5 +912,5 @@ // Convert tm time to psTime - outTime = psTMToTime(&tmTime); + outTime = psTimeFromTM(&tmTime); outTime->usec = millisecond * 1000; @@ -832,5 +918,5 @@ } -psTime* psTimevalToTime(struct timeval *time) +psTime* psTimeFromTimeval(const struct timeval *time) { psTime *outTime = NULL; @@ -850,5 +936,5 @@ } -psTime* psTMToTime(struct tm* time) +psTime* psTimeFromTM(const struct tm* time) { long year; @@ -905,5 +991,5 @@ } -psTime* psTimeAdd(psTime *tai1, psTime *tai2) +psTime* psTimeAdd(const psTime *tai1, const psTime *tai2) { psTime *outTime = NULL; @@ -929,5 +1015,5 @@ } -psTime* psTimeSub(psTime *tai1, psTime *tai2) +psTime* psTimeSubtract(const psTime *tai1, const psTime *tai2) { psTime *outTime = NULL; @@ -953,5 +1039,5 @@ } -psTime* psTimeDelta(psTime *tai1, psTime *tai2) +psTime* psTimeDelta(const psTime *tai1, const psTime *tai2) { psTime *outTime = NULL;