webdial

Incomplete JS port of moondial.
git clone https://code.literati.org/webdial.git
Log | Files | Refs

sun.js (6221B)


      1 function Sun() {
      2     var _DtoR = Math.PI / 180.0;
      3     function d_to_r(d) {
      4         return d * _DtoR;
      5     }
      6 
      7     function r_to_d(r) {
      8         return r / _DtoR;
      9     }
     10 
     11     function dms_to_d(deg, min, sec) {
     12         var result = Math.abs(deg) + Math.abs(min) / 60.0 + Math.abs(sec) /
     13             3600.0;
     14         if (deg < 0 || min < 0 || sec < 0) {
     15             result = -result;
     16         }
     17         return result;
     18     }
     19 
     20     function jd_to_jcent(jd) {
     21         return (jd - 2451545.0) / 36525.0;
     22     }
     23 
     24     function polynomial(terms, x) {
     25         var i = terms.length - 1;
     26         var result = terms[i];
     27         i--;
     28         while (i >= 0) {
     29             result = result * x + terms[i];
     30             i--;
     31         }
     32         return result;
     33     }
     34 
     35     var PI2 = Math.PI * 2.0;
     36     function modpi2(x) {
     37         var y = x % PI2;
     38         if (y < 0.0) {
     39             y += PI2;
     40         }
     41         return y;
     42     }
     43 
     44     /**
     45      * Convert ecliptic to equitorial coordinates. 
     46      *
     47      * [Meeus-1998: equations 13.3, 13.4]
     48      *
     49      * Parameters:
     50      *   longitude : ecliptic longitude in radians
     51      *   latitude : ecliptic latitude in radians
     52      *   obliquity : obliquity of the ecliptic in radians
     53      *
     54      * Returns:
     55      *   Right accension in radians
     56      *   Declination in radians
     57      */
     58     this.ecl_to_equ = function(lon, lat, obl) {
     59         var cose = Math.cos(obl);
     60         var sine = Math.sin(obl);
     61         var sinl = Math.sin(lon);
     62         var ra = modpi2(Math.atan2(sinl * cose - Math.tan(lat) * sine,
     63                                    Math.cos(lon)));
     64         var dec = Math.asin(Math.sin(lat) * cose + Math.cos(lat) * sine * sinl);
     65         return {"ra": ra, "dec": dec};
     66     };
     67 
     68     this.equ_to_geo = function(ra, dec, st) {
     69         var lon = r_to_d(ra - st);
     70         if (lon > 180.0) {
     71             lon -= 360.0;
     72         }
     73         return {longitude: lon, latitude: r_to_d(dec)};
     74     };
     75 
     76     this.sun_rst_altitude = -0.0145438286569;
     77 
     78     this.terminator = function(lat, lon, alt) {
     79         lat = lat / 180.0 * Math.PI;
     80         var obl = lat - Math.PI / 2;
     81         var points = [];
     82         var deg, H, equ, x, y, px, py;
     83         for (deg = 0; deg < 360; deg += 2) {
     84             H = deg * Math.PI / 180.0;
     85             equ = this.ecl_to_equ(H, alt, obl);
     86             x = equ.ra * 180.0 / Math.PI;
     87             y = (0.5 - equ.dec / Math.PI) * 180;
     88             if (lat < alt) {
     89                 if (x > px) {
     90                     points.push([0,0]);
     91                     points.push([359,0]);
     92                 }
     93             } else if (lat > -alt && x < px) {
     94                 points.push([359,179]);
     95                 points.push([0,179]);
     96             }
     97             points.push([x,y]);
     98             px = x;
     99             py = y;
    100         }
    101         return points
    102     };
    103 
    104     this.cal_to_jd = function(date) {
    105         var yr = date.getUTCFullYear();
    106         var mo = date.getUTCMonth() + 1; // WTF???
    107         var day = date.getUTCDate() + date.getUTCHours() / 24.0
    108             + date.getUTCMinutes() / 1440.0
    109             + date.getUTCSeconds() / 86400.0;
    110         if (mo <= 2) {
    111             yr -= 1;
    112             mo += 12;
    113         }
    114         var A = parseInt(yr / 100);
    115         var B = 2 - A + parseInt(A / 4)
    116 
    117         return parseInt(365.25 * (yr + 4716)) + parseInt(30.6001 * (mo + 1)) +
    118             day + B - 1524.5
    119     };
    120 
    121     this.sidereal_time_greenwich = function(jd) {
    122         var T = jd_to_jcent(jd);
    123         var T2 = T * T;
    124         var T3 = T2 * T;
    125         var theta0 = 280.46061837 + 360.98564736629 * (jd - 2451545.0)  +
    126             0.000387933 * T2 - T3 / 38710000;
    127         return modpi2(d_to_r(theta0));
    128     };
    129 
    130     var _el0 = [d_to_r(dms_to_d(23, 26,  21.448)),
    131                 d_to_r(dms_to_d( 0,  0, -46.8150)),
    132                 d_to_r(dms_to_d( 0,  0,  -0.00059)),
    133                 d_to_r(dms_to_d( 0,  0,   0.001813))];
    134 
    135     this.obliquity = function(jd) {
    136         return polynomial(_el0, jd_to_jcent(jd));
    137     }
    138 
    139     // Constant terms
    140     var _kL0 = [d_to_r(280.46646),  d_to_r(36000.76983),  d_to_r( 0.0003032)];
    141     var _kM  = [d_to_r(357.52911),  d_to_r(35999.05029),  d_to_r(-0.0001537)];
    142     var _kC  = [d_to_r( 1.914602),  d_to_r(  -0.004817),  d_to_r( -0.000014)];
    143 
    144     var _ck3 = d_to_r( 0.019993);
    145     var _ck4 = d_to_r(-0.000101);
    146     var _ck5 = d_to_r( 0.000289);
    147 
    148     /**
    149      * Return geometric longitude and radius vector.
    150      *
    151      * Low precision. The longitude is accurate to 0.01 degree.
    152      * The latitude should be presumed to be 0.0. [Meeus-1998: equations 25.2
    153      * through 25.5]
    154      *
    155      * Parameters:
    156      *   jd : Julian Day in dynamical time
    157      *
    158      * Returns:
    159      *   longitude in radians
    160      *   radius in au
    161      */
    162     this.longitude_radius_low = function(jd) {
    163         var T = jd_to_jcent(jd);
    164         var L0 = polynomial(_kL0, T);
    165         var M = polynomial(_kM, T);
    166         var e = polynomial((0.016708634, -0.000042037, -0.0000001267), T);
    167         var C = polynomial(_kC, T) * Math.sin(M)
    168             + (_ck3 - _ck4 * T) * Math.sin(2 * M)
    169             + _ck5 * Math.sin(3 * M);
    170         var L = modpi2(L0 + C);
    171         var v = M + C;
    172         var R = 1.000001018 * (1 - e * e) / (1 + e * Math.cos(v));
    173         return {longitude: L, radius: R}
    174     };
    175 
    176     // Constant terms
    177     var _lk0 = d_to_r(125.04);
    178     var _lk1 = d_to_r(1934.136);
    179     var _lk2 = d_to_r(0.00569);
    180     var _lk3 = d_to_r(0.00478);
    181 
    182     /**
    183      * Correct the geometric longitude for nutation and aberration.
    184      *
    185      * Low precision. [Meeus-1998: pg 164]
    186      *
    187      * Parameters:
    188      *   jd : Julian Day in dynamical time
    189      *   L : longitude in radians
    190      *
    191      * Returns:
    192      *   corrected longitude in radians
    193      */
    194     this.apparent_longitude_low = function(jd, L) {
    195         var T = jd_to_jcent(jd);
    196         var omega = _lk0 - _lk1 * T;
    197         return modpi2(L - _lk2 - _lk3 * Math.sin(omega));
    198     };
    199 
    200     // Constant terms
    201     var _lk4 = d_to_r(dms_to_d(0, 0, 20.4898));
    202 
    203     /**
    204      * Correct for aberration; low precision, but good enough for most uses. 
    205      *
    206      * [Meeus-1998: pg 164]
    207      *
    208      * Parameters:
    209      *   R : radius in au
    210      *
    211      * Returns:
    212      *   correction in radians
    213      */
    214     this.aberration_low = function(R) {
    215         return -_lk4 / R
    216     };
    217 }