From 8531377d9c23143db77c25037cc597c8b98b3295 Mon Sep 17 00:00:00 2001 From: Nemo Date: Fri, 22 May 2020 19:56:49 +0530 Subject: [PATCH] Initial commit --- .editorconfig | 9 +++ .gitignore | 9 +++ .travis.yml | 6 ++ LICENSE | 21 ++++++ README.md | 44 ++++++++++++ shard.yml | 9 +++ spec/spec_helper.cr | 3 + spec/suntime_spec.cr | 56 ++++++++++++++++ src/suntime.cr | 155 +++++++++++++++++++++++++++++++++++++++++++ 9 files changed, 312 insertions(+) create mode 100644 .editorconfig create mode 100644 .gitignore create mode 100644 .travis.yml create mode 100644 LICENSE create mode 100644 README.md create mode 100644 shard.yml create mode 100644 spec/spec_helper.cr create mode 100644 spec/suntime_spec.cr create mode 100644 src/suntime.cr diff --git a/.editorconfig b/.editorconfig new file mode 100644 index 0000000..163eb75 --- /dev/null +++ b/.editorconfig @@ -0,0 +1,9 @@ +root = true + +[*.cr] +charset = utf-8 +end_of_line = lf +insert_final_newline = true +indent_style = space +indent_size = 2 +trim_trailing_whitespace = true diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..0bbd4a9 --- /dev/null +++ b/.gitignore @@ -0,0 +1,9 @@ +/docs/ +/lib/ +/bin/ +/.shards/ +*.dwarf + +# Libraries don't need dependency lock +# Dependencies will be locked in applications that use them +/shard.lock diff --git a/.travis.yml b/.travis.yml new file mode 100644 index 0000000..765f0e9 --- /dev/null +++ b/.travis.yml @@ -0,0 +1,6 @@ +language: crystal + +# Uncomment the following if you'd like Travis to run specs and check code formatting +# script: +# - crystal spec +# - crystal tool format --check diff --git a/LICENSE b/LICENSE new file mode 100644 index 0000000..8dce3a4 --- /dev/null +++ b/LICENSE @@ -0,0 +1,21 @@ +The MIT License (MIT) + +Copyright (c) 2020 Nemo + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in +all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +THE SOFTWARE. diff --git a/README.md b/README.md new file mode 100644 index 0000000..3f9962e --- /dev/null +++ b/README.md @@ -0,0 +1,44 @@ +# suntime + +Crystal library for calculating sunrise and sunset times. Uses the algorithm from . + +## Installation + +1. Add the dependency to your `shard.yml`: + + ```yaml + dependencies: + suntime: + github: captn3m0/suntime + ``` + +2. Run `shards install` + +## Usage + +```crystal +require "suntime" + +# Time is optional, local time is used otherwise +Suntime.new(lat,long, time) +``` + +## Development + +TODO: Write development instructions here + +## Contributing + +1. Fork it () +2. Create your feature branch (`git checkout -b my-new-feature`) +3. Commit your changes (`git commit -am 'Add some feature'`) +4. Push to the branch (`git push origin my-new-feature`) +5. Create a new Pull Request + +## Contributors + +- [Nemo](https://github.com/captn3m0) - creator and maintainer + +## License + +Licensed under the [MIT License](https://nemo.mit-license.org/). See LICENSE file for details. diff --git a/shard.yml b/shard.yml new file mode 100644 index 0000000..6b1ccf5 --- /dev/null +++ b/shard.yml @@ -0,0 +1,9 @@ +name: suntime +version: 0.1.0 + +authors: + - Nemo + +crystal: 0.34.0 + +license: MIT diff --git a/spec/spec_helper.cr b/spec/spec_helper.cr new file mode 100644 index 0000000..355a50c --- /dev/null +++ b/spec/spec_helper.cr @@ -0,0 +1,3 @@ +require "spec" +require "time" +require "../src/suntime" diff --git a/spec/suntime_spec.cr b/spec/suntime_spec.cr new file mode 100644 index 0000000..a2d29ac --- /dev/null +++ b/spec/suntime_spec.cr @@ -0,0 +1,56 @@ +require "./spec_helper" + +describe Suntime do + it "should return correct timestamp for example" do + lat = 40.9 + lng = -74.3 + + l = Time::Location.load("America/New_York") + t = Suntime::Suntime.new(lat, lng, Time.local(1990, 6, 25, 0, 0, 0, location: l)) + + t.day_of_year.should eq(176) + t.approx_time.should eq(176.45638888888888) + t.approx_time(false).should eq(176.95638888888888) + t.sun_mean_anomaly(176.456).should eq(170.6260336) + t.sun_true_longitude(170.6260336).should eq(93.56567911851744) + t.sun_right_ascension(93.56567911851744).should eq(6.2589844084367705) + t.sun_local_hour_angle(93.56567911851744).should eq(-0.39570523787054585) + t.calculate_h(-0.39570523787054585).should eq(16.44600232000648) + t.local_mean_time(16.44600232000648, 6.2589844084367705, 176.456).should eq(5.441396301776585) + + a = t.sunrise + a.year.should eq(1990) + a.month.should eq(6) + a.day.should eq(25) + a.hour.should eq(5) + a.minute.should eq(26) + a.second.should eq(29) + a.offset.should eq(-4 * 60 * 60) + a.location.should eq l + end + + it "should work for bangalore" do + # +5.5 + bangalore_tz = Time::Location.load("Asia/Kolkata") + t = Suntime::Suntime.new(12.955800, 77.620979, Time.local(2020, 5, 23, 0, 0, 0, location: bangalore_tz)) + a = t.sunrise + a.year.should eq(2020) + a.month.should eq(5) + a.day.should eq(23) + a.hour.should eq(5) + a.minute.should eq(52) + a.second.should eq(41) + a.offset.should eq(5.5 * 60 * 60) + a.location.should eq bangalore_tz + + b = t.sunset + b.year.should eq(2020) + b.month.should eq(5) + b.day.should eq(23) + b.hour.should eq(18) + b.minute.should eq(40) + b.second.should eq(0) + b.offset.should eq(5.5 * 60 * 60) + b.location.should eq bangalore_tz + end +end diff --git a/src/suntime.cr b/src/suntime.cr new file mode 100644 index 0000000..ca7bc83 --- /dev/null +++ b/src/suntime.cr @@ -0,0 +1,155 @@ +# TODO: Write documentation for `Suntime` +module Suntime + VERSION = "0.1.0" + + # Currently based on http://www.edwilliams.org/sunrise_sunset_algorithm.htm + # TODO: Switch to https://www.esrl.noaa.gov/gmd/grad/solcalc/solareqns.PDF + # to use radians directly + class Suntime + @lng : Float64 + @lat : Float64 + @t : Time + + def initialize(lat : Float64, lng : Float64, t : Time | Nil = nil) + @t = Time.local + @t = t unless t.nil? + @lng = lng + @lat = lat + end + + def calculate_sun_time(sunrise = true) + t = approx_time(sunrise) + m = sun_mean_anomaly(t) + l = sun_true_longitude(m) + ra = sun_right_ascension(l) + cosH = sun_local_hour_angle(l) + h = calculate_h(cosH, sunrise) + fractional_time_to_proper_time local_mean_time(h, ra, t) + end + + def fractional_time_to_proper_time(t_in_fractional_hours) + minutes = t_in_fractional_hours * 60 + seconds = (minutes * 60).to_i + new_time = @t.at_beginning_of_day + new_time.shift(seconds: seconds) + end + + def sunrise + calculate_sun_time(true) + end + + def sunset + calculate_sun_time(false) + end + + # 1. first calculate the day of the year + def day_of_year + n1 = (275 * @t.month / 9).floor + n2 = ((@t.month + 9) / 12).floor + n3 = (1 + ((@t.year - 4 * (@t.year / 4).floor + 2) / 3).floor) + return n1 - (n2 * n3) + @t.day - 30 + end + + # 2. convert the longitude to hour value and calculate an approximate time + # pass false for sunset + def approx_time(sunrise = true) + lngHour = @lng / 15 + if sunrise + day_of_year + ((6 - lngHour) / 24) + else + day_of_year + ((18 - lngHour) / 24) + end + end + + # 3. calculate the Sun's mean anomaly + def sun_mean_anomaly(t : Float64) : Float64 + (0.9856 * t) - 3.289 + end + + def put_in_range(number, lower, upper, adjuster) + if number > upper + number -= adjuster + elsif number < lower + number += adjuster + else + number + end + end + + # 4. calculate the Sun's true longitude + def sun_true_longitude(m) + m_in_radians = (Math::PI/180) * m + l = m + (1.916 * Math.sin(m_in_radians)) + (0.020 * Math.sin(2 * m_in_radians)) + 282.634 + + # NOTE: L potentially needs to be adjusted into the range [0,360) by adding/subtracting 360 + put_in_range(l, 0, 360, 360) + end + + # 5. calculate the Sun's right ascension + def sun_right_ascension(l) + l_in_radians = (Math::PI/180) * l + ra = (180/Math::PI) * Math.atan(0.91764 * Math.tan(l_in_radians)) + + # 5b. right ascension value needs to be in the same quadrant as L + l_quadrant = (l/90).floor * 90 + ra_quadrant = (ra/90).floor * 90 + ra = ra + (l_quadrant - ra_quadrant) + + # NOTE: RA potentially needs to be adjusted into the range [0,360) by adding/subtracting 360 + ra = put_in_range(ra, 0, 360, 360) + + # 5c. right ascension value needs to be converted into hours + ra/15 + end + + # 6. calculate the Sun's declination + # 7a. calculate the Sun's local hour angle + def sun_local_hour_angle(l, zenith = :official) + case zenith + when :official + z = 1.58534074 + when :civil + z = 1.67552 + when :nautical + z = 1.78024 + when :astronomical + z = 1.88496 + else + z = 1.58534074 + end + + l_in_radians = (Math::PI/180) * l + sinDec = 0.39782 * Math.sin(l_in_radians) + cosDec = Math.cos(Math.asin(sinDec)) + + lat_in_radians = (Math::PI/180) * @lat + + cosH = (Math.cos(z) - (sinDec * Math.sin(lat_in_radians))) / (cosDec * Math.cos(lat_in_radians)) + + raise "the sun never rises on this location (on the specified date)" if cosH > 1 + raise "the sun never sets on this location (on the specified date)" if cosH < -1 + cosH + end + + # 7b. finish calculating H and convert into hours + def calculate_h(cosH, sunrise = true) + h = 0 + if sunrise + h = 360 - (Math.acos(cosH) * (180/Math::PI)) + else + h = Math.acos(cosH) * (180/Math::PI) + end + # convert it to hours + h/15 + end + + # 8. calculate local mean time of rising/setting + def local_mean_time(h, ra, t) + t = h + ra - (0.06571 * t) - 6.622 + lngHour = @lng / 15 + t = (t - lngHour) + t = put_in_range(t, 0, 24, 24) + t + (@t.offset / 3600) + end + end +end