home *** CD-ROM | disk | FTP | other *** search
/ Nebula 2 / Nebula Two.iso / Apps / Astro / Moon / Source / julian.c < prev    next >
Encoding:
C/C++ Source or Header  |  1995-06-12  |  3.1 KB  |  139 lines

  1. /* julian.c
  2.  * Part of the Moon application for the NeXT computer.
  3.  * Authors:  John Walker of Autodesk
  4.  *           Geoffrey S. Knauth (NeXT port)
  5.  * Date:     January 4, 1992
  6.  *
  7.  * This code was placed in the public domain by John Walker.
  8.  */
  9.  
  10. #import "all.h"
  11.  
  12. /*
  13.  * Convert internal GMT date to Julian day.
  14.  */
  15. static long jdate(t)
  16. struct tm *t;
  17. {
  18.     long c, m, y;
  19.     
  20.     y = t->tm_year + 1900;
  21.     m = t->tm_mon + 1;
  22.     if (m > 2)
  23.     m = m - 3;
  24.     else {
  25.     m = m + 9;
  26.     y--;
  27.     }
  28.     c = y / 100L;            /* Compute century */
  29.     y -= 100L * c;
  30.     return t->tm_mday + (c * 146097L) / 4 + (y * 1461L) / 4 +
  31.     (m * 153L + 2) / 5 + 1721119L;
  32. }
  33.  
  34. /*
  35.  * Convert internal GMT date and time to astronomical Julian time
  36.  * (i.e. Julian date plus day fraction, expressed as a double).
  37.  */
  38. double jtime(t)
  39. struct tm *t;
  40. {
  41.     return (jdate(t) - 0.5) + 
  42.     (t->tm_sec + 60 * (t->tm_min + 60 * t->tm_hour)) / 86400.0;
  43. }
  44.  
  45. /*
  46.  * Convert Julian date to year, month, day, which are returned via integer
  47.  * pointers to integers.
  48.  */
  49. void jyear(td, yy, mm, dd)
  50. double td;
  51. int *yy, *mm, *dd;
  52. {
  53.     double j, d, y, m;
  54.     
  55.     td += 0.5;                /* Astronomical to civil */
  56.     j = floor(td);
  57.     j = j - 1721119.0;
  58.     y = floor(((4 * j) - 1) / 146097.0);
  59.     j = (j * 4.0) - (1.0 + (146097.0 * y));
  60.     d = floor(j / 4.0);
  61.     j = floor(((4.0 * d) + 3.0) / 1461.0);
  62.     d = ((4.0 * d) + 3.0) - (1461.0 * j);
  63.     d = floor((d + 4.0) / 4.0);
  64.     m = floor(((5.0 * d) - 3) / 153.0);
  65.     d = (5.0 * d) - (3.0 + (153.0 * m));
  66.     d = floor((d + 5.0) / 5.0);
  67.     y = (100.0 * y) + j;
  68.     if (m < 10.0)
  69.     m = m + 3;
  70.     else {
  71.     m = m - 9;
  72.     y = y + 1;
  73.     }
  74.     *yy = y;
  75.     *mm = m;
  76.     *dd = d;
  77. }
  78.  
  79. /*
  80.  * Convert Julian time to hour, minutes, and seconds.
  81.  */
  82. void jhms(j, h, m, s)
  83. double j;
  84. int *h, *m, *s;
  85. {
  86.     long ij;
  87.     
  88.     j += 0.5;                /* Astronomical to civil */
  89.  
  90.   /* I added the rint(), because I noticed that local time was sometimes
  91.    * one second behind the value I expected.  --gsk 1/4/92
  92.    */
  93.     ij = rint((j - floor(j)) * 86400.0);
  94.  
  95.     *h = ij / 3600L;
  96.     *m = (ij / 60L) % 60L;
  97.     *s = ij % 60L;
  98. }
  99.  
  100. /*********************************************************************
  101.  * The following routines were added to the original Sun moontool code
  102.  * so that the NeXT version could offer the user time travel.
  103.  *********************************************************************
  104.  */
  105.  
  106. /*
  107.  * Convert {year, month, day} to astronomical Julian date.
  108.  * This is a long, because there is no fraction part (yet).
  109.  */
  110. static long ymdToJdate(int year, int month, int day)
  111. {
  112.     long c, m, y;
  113.     
  114.     y = year;
  115.     m = month;
  116.     if (m > 2)
  117.     m = m - 3;
  118.     else {
  119.     m = m + 9;
  120.     y--;
  121.     }
  122.     c = y / 100L;            /* Compute century */
  123.     y -= 100L * c;
  124.     return day + (c * 146097L) / 4 + (y * 1461L) / 4 +
  125.     (m * 153L + 2) / 5 + 1721119L;
  126. }
  127.  
  128. /*
  129.  * Convert {year, month, day, hour, minute, second} to astronomical Julian
  130.  * date, which will include a fraction, which is why this is a double.
  131.  */
  132. double ymdhmsToJtime
  133.     (int year, int month, int day,
  134.      int hour, int minute, int second)
  135. {
  136.     return (ymdToJdate(year, month, day) - 0.5) + 
  137.     (second + 60 * (minute + 60 * hour)) / 86400.0;
  138. }
  139.