Index: trunk/psLib/src/astro/psTime.c
===================================================================
--- trunk/psLib/src/astro/psTime.c	(revision 5684)
+++ 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;