Hourly solar analemmas (ignoring daylight savings time) as seen from San Francisco in 2014.
| <!DOCTYPE html> | |
| <meta charset="utf-8"> | |
| <style> | |
| path { | |
| fill: none; | |
| stroke-linecap: round; | |
| stroke-linejoin: round; | |
| } | |
| text { | |
| font: 10px sans-serif; | |
| } | |
| .horizon { | |
| stroke: #000; | |
| stroke-width: 1.5px; | |
| } | |
| .graticule { | |
| stroke: #000; | |
| stroke-opacity: .15; | |
| } | |
| .analemmas path { | |
| stroke: #f00; | |
| stroke-width: 2px; | |
| } | |
| .ticks line { | |
| stroke: #000; | |
| } | |
| .ticks text { | |
| text-anchor: middle; | |
| } | |
| .ticks--azimuth text:nth-of-type(9n + 1) { | |
| font-weight: bold; | |
| font-size: 14px; | |
| } | |
| </style> | |
| <body> | |
| <script src="https://cdnjs.cloudflare.com/ajax/libs/d3/3.5.5/d3.min.js"></script> | |
| <script src="solar-calculator.js"></script> | |
| <script> | |
| var π = Math.PI, | |
| τ = 2 * π, | |
| radians = π / 180, | |
| degrees = 180 / π; | |
| var width = 960, | |
| height = 960, | |
| scale = width * .45; | |
| var solar = solarCalculator([-122, 37]), | |
| start = d3.time.year.utc.floor(new Date), | |
| end = d3.time.year.utc.offset(start, 1); | |
| var projection = d3.geo.projection(flippedStereographic) | |
| .scale(scale) | |
| .clipAngle(130) | |
| .rotate([0, -90]) | |
| .translate([width / 2 + .5, height / 2 + .5]) | |
| .precision(.1); | |
| var path = d3.geo.path() | |
| .projection(projection); | |
| var svg = d3.select("body").append("svg") | |
| .attr("width", width) | |
| .attr("height", height); | |
| svg.append("path") | |
| .datum(d3.geo.circle().origin([0, 90]).angle(90)) | |
| .attr("class", "horizon") | |
| .attr("d", path); | |
| svg.append("path") | |
| .datum(d3.geo.graticule()) | |
| .attr("class", "graticule") | |
| .attr("d", path); | |
| var ticksAzimuth = svg.append("g") | |
| .attr("class", "ticks ticks--azimuth"); | |
| ticksAzimuth.selectAll("line") | |
| .data(d3.range(360)) | |
| .enter().append("line") | |
| .each(function(d) { | |
| var p0 = projection([d, 0]), | |
| p1 = projection([d, d % 10 ? -1 : -2]); | |
| d3.select(this) | |
| .attr("x1", p0[0]) | |
| .attr("y1", p0[1]) | |
| .attr("x2", p1[0]) | |
| .attr("y2", p1[1]); | |
| }); | |
| ticksAzimuth.selectAll("text") | |
| .data(d3.range(0, 360, 10)) | |
| .enter().append("text") | |
| .each(function(d) { | |
| var p = projection([d, -4]); | |
| d3.select(this) | |
| .attr("x", p[0]) | |
| .attr("y", p[1]); | |
| }) | |
| .attr("dy", ".35em") | |
| .text(function(d) { return d === 0 ? "N" : d === 90 ? "E" : d === 180 ? "S" : d === 270 ? "W" : d + "°"; }); | |
| svg.append("g") | |
| .attr("class", "ticks ticks--elevation") | |
| .selectAll("text") | |
| .data(d3.range(10, 91, 10)) | |
| .enter().append("text") | |
| .each(function(d) { | |
| var p = projection([0, d]); | |
| d3.select(this) | |
| .attr("x", p[0]) | |
| .attr("y", p[1]); | |
| }) | |
| .attr("dy", ".35em") | |
| .text(function(d) { return d + "°"; }); | |
| svg.insert("g", ".sphere") | |
| .attr("class", "analemmas") | |
| .selectAll("path") | |
| .data(d3.range(24)) | |
| .enter().append("path") | |
| .datum(function(h) { return {type: "LineString", coordinates: d3.time.days.utc(start, end).map(function(d) { return solar.position(d3.time.hour.utc.offset(d, h)); })}; }) | |
| .attr("d", path); | |
| d3.select(self.frameElement).style("height", height + "px"); | |
| function flippedStereographic(λ, φ) { | |
| var cosλ = Math.cos(λ), | |
| cosφ = Math.cos(φ), | |
| k = 1 / (1 + cosλ * cosφ); | |
| return [ | |
| k * cosφ * Math.sin(λ), | |
| -k * Math.sin(φ) | |
| ]; | |
| } | |
| </script> |
| // Equations based on NOAA’s Solar Calculator; all angles in radians. | |
| // http://www.esrl.noaa.gov/gmd/grad/solcalc/ | |
| (function() { | |
| var J2000 = Date.UTC(2000, 0, 1, 12), | |
| π = Math.PI, | |
| τ = 2 * π, | |
| radians = π / 180, | |
| degrees = 180 / π; | |
| solarCalculator = function(location) { | |
| var longitude = location[0], | |
| minutesOffset = 720 - longitude * 4, | |
| λ = location[0] * radians, | |
| φ = location[1] * radians, | |
| cosφ = Math.cos(φ), | |
| sinφ = Math.sin(φ); | |
| function position(date) { | |
| var centuries = (date - J2000) / (864e5 * 36525), | |
| θ = solarDeclination(centuries), | |
| cosθ = Math.cos(θ), | |
| sinθ = Math.sin(θ), | |
| azimuth = ((date - d3.time.day.utc.floor(date)) / 864e5 * τ + equationOfTime(centuries) + λ) % τ - π, | |
| zenith = Math.acos(Math.max(-1, Math.min(1, sinφ * sinθ + cosφ * cosθ * Math.cos(azimuth)))), | |
| azimuthDenominator = cosφ * Math.sin(zenith); | |
| if (azimuth < -π) azimuth += τ; | |
| if (Math.abs(azimuthDenominator) > 1e-6) azimuth = (azimuth > 0 ? -1 : 1) * (π - Math.acos(Math.max(-1, Math.min(1, (sinφ * Math.cos(zenith) - sinθ) / azimuthDenominator)))); | |
| if (azimuth < 0) azimuth += τ; | |
| // Correct for atmospheric refraction. | |
| var atmosphere = 90 - zenith * degrees; | |
| if (atmosphere <= 85) { | |
| var te = Math.tan(atmosphere * radians); | |
| zenith -= (atmosphere > 5 ? 58.1 / te - .07 / (te * te * te) + .000086 / (te * te * te * te * te) | |
| : atmosphere > -.575 ? 1735 + atmosphere * (-518.2 + atmosphere * (103.4 + atmosphere * (-12.79 + atmosphere * .711))) | |
| : -20.774 / te) / 3600 * radians; | |
| } | |
| // Note: if zenith > 108°, it’s dark. | |
| return [azimuth * degrees, 90 - zenith * degrees]; | |
| } | |
| function noon(date) { | |
| var centuries = (d3.time.day.utc.floor(date) - J2000) / (864e5 * 36525), | |
| minutes = (minutesOffset - (equationOfTime(centuries + (minutesOffset - (equationOfTime(centuries - longitude / (360 * 365.25 * 100)) * degrees * 4)) / (1440 * 365.25 * 100)) * degrees * 4) - date.getTimezoneOffset()) % 1440; | |
| if (minutes < 0) minutes += 1440; | |
| return new Date(+d3.time.day.floor(date) + minutes * 60 * 1000); | |
| } | |
| return { | |
| position: position, | |
| noon: noon | |
| }; | |
| }; | |
| function equationOfTime(centuries) { | |
| var e = eccentricityEarthOrbit(centuries), | |
| m = solarGeometricMeanAnomaly(centuries), | |
| l = solarGeometricMeanLongitude(centuries), | |
| y = Math.tan(obliquityCorrection(centuries) / 2); | |
| y *= y; | |
| return y * Math.sin(2 * l) | |
| - 2 * e * Math.sin(m) | |
| + 4 * e * y * Math.sin(m) * Math.cos(2 * l) | |
| - 0.5 * y * y * Math.sin(4 * l) | |
| - 1.25 * e * e * Math.sin(2 * m); | |
| } | |
| function solarDeclination(centuries) { | |
| return Math.asin(Math.sin(obliquityCorrection(centuries)) * Math.sin(solarApparentLongitude(centuries))); | |
| } | |
| function solarApparentLongitude(centuries) { | |
| return solarTrueLongitude(centuries) - (0.00569 + 0.00478 * Math.sin((125.04 - 1934.136 * centuries) * radians)) * radians; | |
| } | |
| function solarTrueLongitude(centuries) { | |
| return solarGeometricMeanLongitude(centuries) + solarEquationOfCenter(centuries); | |
| } | |
| function solarGeometricMeanAnomaly(centuries) { | |
| return (357.52911 + centuries * (35999.05029 - 0.0001537 * centuries)) * radians; | |
| } | |
| function solarGeometricMeanLongitude(centuries) { | |
| var l = (280.46646 + centuries * (36000.76983 + centuries * 0.0003032)) % 360; | |
| return (l < 0 ? l + 360 : l) / 180 * π; | |
| } | |
| function solarEquationOfCenter(centuries) { | |
| var m = solarGeometricMeanAnomaly(centuries); | |
| return (Math.sin(m) * (1.914602 - centuries * (0.004817 + 0.000014 * centuries)) | |
| + Math.sin(m + m) * (0.019993 - 0.000101 * centuries) | |
| + Math.sin(m + m + m) * 0.000289) * radians; | |
| } | |
| function obliquityCorrection(centuries) { | |
| return meanObliquityOfEcliptic(centuries) + 0.00256 * Math.cos((125.04 - 1934.136 * centuries) * radians) * radians; | |
| } | |
| function meanObliquityOfEcliptic(centuries) { | |
| return (23 + (26 + (21.448 - centuries * (46.8150 + centuries * (0.00059 - centuries * 0.001813))) / 60) / 60) * radians; | |
| } | |
| function eccentricityEarthOrbit(centuries) { | |
| return 0.016708634 - centuries * (0.000042037 + 0.0000001267 * centuries); | |
| } | |
| })(); |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment