|
| 1 | +!function(exports, Math) { |
| 2 | + exports.sun = function sun(lat, lon, time, recalc) { |
| 3 | + time = new Date(time) |
| 4 | + |
| 5 | + var tmp |
| 6 | + , acos = Math.acos |
| 7 | + , cos = Math.cos |
| 8 | + , sin = Math.sin |
| 9 | + , tan = Math.tan |
| 10 | + , PI = Math.PI |
| 11 | + , rad = PI / 180 |
| 12 | + , dayMs = 86400000 |
| 13 | + , startOfDay = Math.floor(time / dayMs) * dayMs |
| 14 | + |
| 15 | + // Julian Century |
| 16 | + , JC = (time / dayMs - 10957.5) / 36525 |
| 17 | + , cosLat = cos(rad * lat) |
| 18 | + , sinLat = sin(rad * lat) |
| 19 | + |
| 20 | + // Mean Longitude of the Sun |
| 21 | + , lonMean = rad * (280.46646 + JC * (36000.76983 + JC*0.0003032)) |
| 22 | + // Mean Anomaly of the Sun |
| 23 | + , anomMean = rad * (357.52911 + JC * (35999.05029 - 0.0001537 * JC)) |
| 24 | + // Mean Obliquity of the Ecliptic (degrees) |
| 25 | + , obliqMean = 23 + (26 + (21.448 - JC*(46.8150 + JC*(0.00059 - JC*0.001813)))/60)/60 |
| 26 | + // Corrected Obliquity |
| 27 | + , obliqCorr = rad * (obliqMean + 0.00256 * cos(tmp = rad * (125.04 - 1934.136 * JC))) |
| 28 | + |
| 29 | + // Sun's Equation of Center |
| 30 | + , sunEqCtr = sin(anomMean) * (1.914602 - JC * (0.004817 + 0.000014 * JC)) + sin(2 * anomMean) * (0.019993 - 0.000101 * JC) + sin(3 * anomMean) * 0.000289 |
| 31 | + |
| 32 | + // Apparent Longitude of the Sun |
| 33 | + , lonApp = lonMean + rad * (sunEqCtr - 0.00569 - 0.00478 * sin(tmp)) |
| 34 | + |
| 35 | + // Declination of the Sun |
| 36 | + , declRad = Math.asin(sin(obliqCorr) * sin(lonApp)) |
| 37 | + , sinDecl = sin(declRad) |
| 38 | + , cosDecl = cos(declRad) |
| 39 | + |
| 40 | + // Eccentricity of Earth's Orbit |
| 41 | + , earthEccent = 0.016708634 - JC * (0.000042037 + 0.0000001267 * JC) |
| 42 | + , varY = (tmp = tan(obliqCorr / 2)) * tmp |
| 43 | + |
| 44 | + // Equation of Time (minutes) |
| 45 | + , eqTime = 4 * ( |
| 46 | + varY * sin(2 * lonMean) - |
| 47 | + 2 * earthEccent * sin(anomMean) + |
| 48 | + 4 * earthEccent * varY * sin(anomMean) * cos(2 * lonMean) - |
| 49 | + 0.5 * varY * varY * sin(4 * lonMean) - |
| 50 | + 1.25 * earthEccent * earthEccent * sin(2 * anomMean) |
| 51 | + ) / rad |
| 52 | + |
| 53 | + // Hour Angle at Sunrise (degrees) |
| 54 | + , hourAngleRise = acos((cos(rad * 90.833) - sinLat * sinDecl) / (cosLat * cosDecl)) / rad |
| 55 | + |
| 56 | + // True Solar Time (minutes) |
| 57 | + , sunTime = (time / 60000 + eqTime + 4 * lon) % 1440 |
| 58 | + , hourAngle = sunTime < 0 ? sunTime / 4 + 180 : sunTime / 4 - 180 |
| 59 | + |
| 60 | + // Solar Zenith/Elevation |
| 61 | + , cosZenith = sinLat * sinDecl + cosLat * cosDecl * cos(rad * hourAngle) |
| 62 | + , elTrue = 90 - acos(cosZenith) / rad |
| 63 | + // Atmospheric Refraction Correction (degrees) |
| 64 | + , atmCorr = ( |
| 65 | + elTrue > 85 ? 0 : (tmp = tan(rad * elTrue), |
| 66 | + elTrue > 5 ? 58.1 / tmp - 0.07 / (tmp*tmp*tmp) + 0.000086 / (tmp*tmp*tmp*tmp*tmp) : |
| 67 | + elTrue > -0.575 ? 1735 + elTrue * (-518.2 + elTrue * (103.4 + elTrue * (-12.79 + elTrue * 0.711))) : |
| 68 | + -20.774 / tmp |
| 69 | + ) / 3600 |
| 70 | + ) |
| 71 | + |
| 72 | + // Solar Noon (minutes) |
| 73 | + , noon = 720 - 4 * lon - eqTime |
| 74 | + // Sunrise & Sunset Time (minutes) |
| 75 | + , rise = noon - hourAngleRise * 4 |
| 76 | + , set = noon + hourAngleRise * 4 |
| 77 | + |
| 78 | + if (recalc) { |
| 79 | + tmp = Math.round((recalc === "rise" ? rise : recalc === "set" ? set : noon) * 60) * 1000 |
| 80 | + time.setTime(startOfDay + tmp) |
| 81 | + } |
| 82 | + |
| 83 | + return { |
| 84 | + time, |
| 85 | + decl: declRad / rad, |
| 86 | + dur: 8 * hourAngleRise, // Total Sunlight Duration (minutes) |
| 87 | + el: elTrue + atmCorr, // Solar Elevation (Corrected for Refraction) |
| 88 | + elTrue, |
| 89 | + eqTime, |
| 90 | + hourAngle, |
| 91 | + get azimuth() { |
| 92 | + return (tmp = acos((sinLat * cosZenith - sinDecl) / (cosLat * Math.sqrt(1 - cosZenith * cosZenith))), hourAngle > 0 ? tmp + PI : 3 * PI - tmp) * 180 / PI % 360 |
| 93 | + }, |
| 94 | + get dist() { |
| 95 | + return 1.000001018 * (1 - earthEccent * earthEccent) / (1 + earthEccent * cos(anomMean + rad * sunEqCtr)) // Distance to Sun (Astronomical Units) |
| 96 | + }, |
| 97 | + get rightAsc() { |
| 98 | + return Math.atan2(cos(obliqCorr) * sin(lonApp), cos(lonApp)) / rad // Right Ascension of the Sun (degrees) |
| 99 | + }, |
| 100 | + get noon() { |
| 101 | + return findNext("noon", noon, -1) |
| 102 | + }, |
| 103 | + get rise() { |
| 104 | + return findNext("rise", rise, 1) |
| 105 | + }, |
| 106 | + get set() { |
| 107 | + return findNext("set", set, -1) |
| 108 | + }, |
| 109 | + } |
| 110 | + |
| 111 | + function findNext(key, minutes, inc) { |
| 112 | + if (minutes === minutes) return sun(lat, lon, startOfDay + minutes * 60000, key) |
| 113 | + var day = startOfDay, i = 366, d = (time - Date.UTC(time.getUTCFullYear(), 0, 0)) / dayMs |
| 114 | + for (inc *= (lat > 66.4 && d > 79 && d < 267 || lat < -66.4 && (d < 83 || d > 263)) ? -dayMs : dayMs; i--; ) |
| 115 | + if ((d = sun(lat, lon, day += inc)).dur === d.dur) return d[key] |
| 116 | + return null |
| 117 | + } |
| 118 | + } |
| 119 | +}(this, Math) // jshint ignore:line |
0 commit comments