Solar Compass? なんて呼べばいいのかなコレ。

日時計的なもの、とも言える。

Solar Compass を kivyLauncher で、またの名を、「True North を知る」の別解

動機

聞くな。ほんと個人的興味だよ。

だけではなんなので…。

  1. これとスマホの組み合わせは相性よかろ?
  2. これこれこれ。個人的興味から「実際の磁偏差ってどんだけなの?」を知りたい。

お手本

コンパスに連動してないものでは2つ:

地図背景はワタシはいらんの。これらの、sun rise / sun set の時刻・方位、現在時刻での方位と高度(角度ね)をリアルタイムで、なおかつコンパスに連動したものが欲しい、と。

本質以外の準備

今回の本質ではない前提3つ:

本題: noaa_solar_calc.py

日出日没計算、やってみようを書いた時点では紹介してないが、先の NOAA Solar Calculator の javascript はかなり直裁的に書かれてる。つまり移植が容易。

一部意図的に機能を削ったが、たぶんこんな具合:

noaa_solar_calc.py (一気に移植して真剣なテストしてないので信頼しすぎるべからず)
  1 import math
  2 from math import radians, degrees, floor
  3 
  4 def _sind(d):
  5     return math.sin(radians(d))
  6 
  7 def _cosd(d):
  8     return math.cos(radians(d))
  9 
 10 def _tand(d):
 11     return math.tan(radians(d))
 12 
 13 
 14 def calcTimeJulianCent(jd):
 15     return (jd - 2451545.0) / 36525.0
 16 
 17 def calcJDFromJulianCent(t):
 18     return t * 36525.0 + 2451545.0
 19 
 20 def isLeap(yr):
 21     return ((yr % 4 == 0 and yr % 100 != 0) or yr % 400 == 0)
 22 
 23 def calcDoyFromJD(jd):
 24     z = floor(jd + 0.5)
 25     f = (jd + 0.5) - z
 26     if (z < 2299161):
 27         A = z
 28     else:
 29         alpha = floor((z - 1867216.25) / 36524.25)
 30         A = z + 1 + alpha - floor(alpha / 4)
 31 
 32     B = A + 1524
 33     C = floor((B - 122.1) / 365.25)
 34     D = floor(365.25 * C)
 35     E = floor((B - D) / 30.6001)
 36     day = B - D - floor(30.6001 * E) + f
 37     month = E - 1 if (E < 14) else E - 13
 38     year = C - 4716 if (month > 2) else C - 4715
 39 
 40     k = 1 if isLeap(year) else 2
 41     doy = floor((275 * month) / 9) - k * floor((month + 9) / 12) + day - 30
 42     return doy
 43 
 44 def calcGeomMeanLongSun(t):
 45     L0 = 280.46646 + t * (36000.76983 + t * (0.0003032))
 46     while (L0 > 360.0):
 47         L0 -= 360.0
 48     while (L0 < 0.0):
 49         L0 += 360.0
 50     return L0  # in degrees
 51 
 52 def calcGeomMeanAnomalySun(t):
 53     M = 357.52911 + t * (35999.05029 - 0.0001537 * t)
 54     return M  # in degrees
 55 
 56 def calcEccentricityEarthOrbit(t):
 57     e = 0.016708634 - t * (0.000042037 + 0.0000001267 * t)
 58     return e  # unitless
 59 
 60 def calcSunEqOfCenter(t):
 61     m = calcGeomMeanAnomalySun(t)
 62     mrad = radians(m)
 63     sinm = math.sin(mrad)
 64     sin2m = math.sin(mrad + mrad)
 65     sin3m = math.sin(mrad + mrad + mrad)
 66     C = sinm * (1.914602 - t * (0.004817 + 0.000014 * t)) + sin2m * (0.019993 - 0.000101 * t) + sin3m * 0.000289
 67     return C  # in degrees
 68 
 69 def calcSunTrueLong(t):
 70     l0 = calcGeomMeanLongSun(t)
 71     c = calcSunEqOfCenter(t)
 72     O = l0 + c
 73     return O  # in degrees
 74 
 75 def calcSunTrueAnomaly(t):
 76     m = calcGeomMeanAnomalySun(t)
 77     c = calcSunEqOfCenter(t)
 78     v = m + c
 79     return v  # in degrees
 80 
 81 def calcSunRadVector(t):
 82     v = calcSunTrueAnomaly(t)
 83     e = calcEccentricityEarthOrbit(t)
 84     R = (1.000001018 * (1 - e * e)) / (1 + e * _cosd(v))
 85     return R  # in AUs
 86 
 87 def calcSunApparentLong(t):
 88     o = calcSunTrueLong(t)
 89     omega = 125.04 - 1934.136 * t
 90     lambda_ = o - 0.00569 - 0.00478 * _sind(omega)
 91     return lambda_  # in degrees
 92 
 93 def calcMeanObliquityOfEcliptic(t):
 94     seconds = 21.448 - t * (46.8150 + t * (0.00059 - t * (0.001813)))
 95     e0 = 23.0 + (26.0 + (seconds / 60.0)) / 60.0
 96     return e0  # in degrees
 97 
 98 def calcObliquityCorrection(t):
 99     e0 = calcMeanObliquityOfEcliptic(t)
100     omega = 125.04 - 1934.136 * t
101     e = e0 + 0.00256 * _cosd(omega)
102     return e  # in degrees
103 
104 def calcSunRtAscension(t):
105     e = calcObliquityCorrection(t)
106     lambda_ = calcSunApparentLong(t)
107     tananum = (_cosd(e) * _sind(lambda_))
108     tanadenom = (_cosd(lambda_))
109     alpha = degrees(math.atan2(tananum, tanadenom))
110     return alpha  # in degrees
111 
112 def calcSunDeclination(t):
113     e = calcObliquityCorrection(t)
114     lambda_ = calcSunApparentLong(t)
115     sint = _sind(e) * _sind(lambda_)
116     theta = degrees(math.asin(sint))
117     return theta  # in degrees
118 
119 def calcEquationOfTime(t):
120     epsilon = calcObliquityCorrection(t)
121     l0 = calcGeomMeanLongSun(t)
122     e = calcEccentricityEarthOrbit(t)
123     m = calcGeomMeanAnomalySun(t)
124 
125     y = math.tan(radians(epsilon) / 2.0)
126     y *= y
127 
128     sin2l0 = math.sin(2.0 * radians(l0))
129     sinm = _sind(m)
130     cos2l0 = math.cos(2.0 * radians(l0))
131     sin4l0 = math.sin(4.0 * radians(l0))
132     sin2m = math.sin(2.0 * radians(m))
133 
134     Etime = y * sin2l0 - 2.0 * e * sinm + 4.0 * e * y * sinm * cos2l0 - 0.5 * y * y * sin4l0 - 1.25 * e * e * sin2m
135     return degrees(Etime) * 4.0  # in minutes of time
136 
137 def calcHourAngleSunrise(lat, solarDec):
138     HAarg = (_cosd(90.833) / (_cosd(lat) * _cosd(solarDec)) - _tand(lat) * _tand(solarDec))
139     HA = math.acos(HAarg)
140     return HA  # in radians (for sunset, use -HA)
141 
142 _monthList = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]
143 
144 def getJD(year, month, day):
145     if ((isLeap(year)) and (month == 2)):
146         day = min(day, 29)
147     else:
148         day = min(day, _monthList[month - 1])
149     if (month <= 2):
150         year -= 1
151         month += 12
152     A = floor(year / 100)
153     B = 2 - A + floor(A / 4)
154     JD = floor(365.25 * (year + 4716)) + floor(30.6001 * (month + 1)) + day + B - 1524.5
155     return JD
156 
157 def getTimeLocal(hr, mn, sc, dst=False):
158     if dst:
159         hr -= 1
160     mins = hr * 60 + mn + sc / 60.0
161     return mins
162 
163 def calcAzEl(T, localtime, longitude, latitude, zone):
164     eqTime = calcEquationOfTime(T)
165     theta = calcSunDeclination(T)
166     #
167     solarTimeFix = eqTime + 4.0 * longitude - 60.0 * zone
168     earthRadVec = calcSunRadVector(T)
169     trueSolarTime = localtime + solarTimeFix
170     while (trueSolarTime > 1440):
171         trueSolarTime -= 1440
172     hourAngle = trueSolarTime / 4.0 - 180.0
173     if (hourAngle < -180):
174         hourAngle += 360.0
175     csz = _sind(latitude) * _sind(theta) + _cosd(latitude) * _cosd(theta) * _cosd(hourAngle)
176     if (csz > 1.0):
177         csz = 1.0
178     elif (csz < -1.0):
179         csz = -1.0
180     zenith = degrees(math.acos(csz))
181     azDenom = (_cosd(latitude) * _sind(zenith))
182     if (abs(azDenom) > 0.001):
183         azRad = ((_sind(latitude) * _cosd(zenith)) - _sind(theta)) / azDenom
184         if (abs(azRad) > 1.0):
185             if (azRad < 0):
186                 azRad = -1.0
187             else:
188                 azRad = 1.0
189         azimuth = 180.0 - degrees(math.acos(azRad))
190         if (hourAngle > 0.0):
191             azimuth = -azimuth
192     else:
193         if (latitude > 0.0):
194             azimuth = 180.0
195         else:
196             azimuth = 0.0
197     if (azimuth < 0.0):
198         azimuth += 360.0
199     exoatmElevation = 90.0 - zenith
200 
201     # Atmospheric Refraction correction
202 
203     if (exoatmElevation > 85.0):
204         refractionCorrection = 0.0
205     else:
206         te = _tand(exoatmElevation)
207         if (exoatmElevation > 5.0):
208             refractionCorrection = 58.1 / te - 0.07 / (te*te*te) + 0.000086 / (te*te*te*te*te)
209         elif (exoatmElevation > -0.575):
210             refractionCorrection = 1735.0 + exoatmElevation * (
211                 -518.2 + exoatmElevation * (103.4 + exoatmElevation * (-12.79 + exoatmElevation * 0.711)))
212         else:
213             refractionCorrection = -20.774 / te
214         refractionCorrection = refractionCorrection / 3600.0
215 
216     solarZen = zenith - refractionCorrection
217 
218     return {
219         'LocalTime': localtime,
220         'Azimuth': azimuth,  # clockwise from TN in decrees
221         'Elevation': 90.0 - solarZen,  # solar angle in degrees
222         'Equation of Time': eqTime,  # in minutes
223         'Solar Declination': theta,  # in degrees
224         }
225 
226 def calcSolNoon(jd, longitude, timezone, dst):
227     tnoon = calcTimeJulianCent(jd - longitude / 360.0)
228     eqTime = calcEquationOfTime(tnoon)
229     solNoonOffset = 720.0 - (longitude * 4) - eqTime  # in minutes
230     newt = calcTimeJulianCent(jd + solNoonOffset / 1440.0)
231     eqTime = calcEquationOfTime(newt)
232     solNoonLocal = 720 - (longitude * 4) - eqTime + (timezone*60.0)  # in minutes
233     if (dst):
234         solNoonLocal += 60.0
235     while (solNoonLocal < 0.0):
236         solNoonLocal += 1440.0
237     while (solNoonLocal >= 1440.0):
238         solNoonLocal -= 1440.0
239     return solNoonLocal
240 
241 def calcSunriseSetUTC(rise, JD, longitude, latitude):
242     t = calcTimeJulianCent(JD);
243     eqTime = calcEquationOfTime(t);
244     solarDec = calcSunDeclination(t);
245     hourAngle = calcHourAngleSunrise(latitude, solarDec);
246     if (not rise):
247         hourAngle = -hourAngle;
248     delta = longitude + degrees(hourAngle);
249     timeUTC = 720 - (4.0 * delta) - eqTime  # in minutes
250     return timeUTC
251 
252 def calcSunriseSet(rise, JD, longitude, latitude, timezone, dst):
253     """rise = 1 for sunrise, 0 for sunset"""
254 
255     timeUTC = calcSunriseSetUTC(rise, JD, longitude, latitude)
256     newTimeUTC = calcSunriseSetUTC(rise, JD + timeUTC / 1440.0, longitude, latitude)
257     # if newTimeUTC is NaN, it means that no sunrise/set found,
258     # but I'm not interesting those case, so I'd omitted.
259     # See http://www.esrl.noaa.gov/gmd/grad/solcalc/main.js.
260 
261     timeLocal = newTimeUTC + (timezone * 60.0)
262 
263     riseT = calcTimeJulianCent(JD + newTimeUTC / 1440.0)
264     riseAz = calcAzEl(riseT, timeLocal, longitude, latitude, timezone)['Azimuth']
265 
266     timeLocal += 60.0 if dst else 0.0
267     if not ((timeLocal >= 0.0) and (timeLocal < 1440.0)):
268         jday = JD
269         increment = 1 if (timeLocal < 0) else -1
270         while ((timeLocal < 0.0) or (timeLocal >= 1440.0)):
271             timeLocal += increment * 1440.0
272             jday -= increment
273     return {
274         'LocalTime': timeLocal,
275         'Azimuth': riseAz,
276         }
277 
278 def calculate(lon, lat, year, month, day, hr, mn, sc, tz, dst=False):
279     jday = getJD(year, month, day)
280     tl = getTimeLocal(hr, mn, sc, dst)
281     total = jday + tl / 1440.0 - tz / 24.0
282     T = calcTimeJulianCent(total)
283 
284     # --- followings are now testing...
285     #solnoonlt = calcSolNoon(jday, lon, tz, dst)
286     #at_solnoon = calcAzEl(
287     #    calcTimeJulianCent(jday + solnoonlt / 1440.0 - tz / 24.0),
288     #    solnoonlt, lon, lat, tz)
289     #
290 
291     return {
292         "sun rise": calcSunriseSet(1, jday, lon, lat, tz, dst),
293         #"solar noon": at_solnoon,
294         "sun set": calcSunriseSet(0, jday, lon, lat, tz, dst),
295         "current": calcAzEl(T, tl, lon, lat, tz),
296         }

日出日没計算、やってみようでちょろっと説明してるんだけれど、このモジュールの calculate の戻りの「LocalTime」は 00:00 からの積算分。ユリウス通日なんぞ経由しちゃうんで、なんとなく得体の知れないものなんじゃないかと恐れ戦いてしまいそうなんだけれど、冷静に順を追って理解すれば単にグレゴリウス暦現地時間の深夜零時からの積算分であることはわかる、はず。色んな同業モジュールが「スマートな日付時刻オブジェクト」に一所懸命読み替えてたりするんだけれど、別に積算分で誰が困るんだ、と思うのでそのまま返してる。これは本家のノリとは違う。

ハイライトした部分が削った機能で、ここは本物の NOAA Solar Calculator では「sun rise / sun set が計算出来ない」、つまり百夜・百昼時は「次に sun rise / sun set が来るのはいつぞ?」を計算している。ちょっと javascript と Python の振る舞いの違いを考えるのが面倒だったこともあるし、今のところは興味なかったんで。

あと solar noon 時の計算、本家は時刻しか計算してなくて、このときの azimuth, elevation も計算したい、って思ったのだけど、まだ合ってるかどうか未検証なのでコメントアウトしてある。

astronomical twilight、voyage twilight、citizen twilight も計算出来ると面白いのだけどなぁ、とは思うんだけれど、NOAA Solar Calculator でやってないので、移植版でもやってない。これはやれそうな気はする。

あと NOAA Solar Calculator は現在高度の考慮がない。これを入れるのはちょっと理論的なことを正確に理解してないと書けないのでキビシイかも。(ただしすげー簡単な可能性もある。)

全てを組み合わせる

android.txt
1 title=Simple Solar Compass
2 author=hhsprings
3 orientation=vertical
sensor_data.py
 1 # -*- coding: utf-8 -*-
 2 import time
 3 from datetime import datetime
 4 from kivy.clock import Clock, mainthread
 5 from plyer import gps
 6 from kivy.logger import Logger
 7 
 8 #------------------------------------------------------
 9 #from plyer import compass
10 from my_compass import AndroidCompass
11 compass = AndroidCompass()
12 #------------------------------------------------------
13 
14 class SensorData(object):
15     def __init__(self, **kwargs):
16         self._ready = False
17         self.lon = None
18         self.lat = None
19         self.alt = None
20 
21         compass.enable()
22         try:
23             gps.configure(
24                 on_status=self._on_status,
25                 on_location=self._on_location)
26             gps.start()
27         except NotImplementedError:
28             import traceback
29             traceback.print_exc()
30             #self.gps_status = 'GPS is not implemented for your platform'
31 
32     @property
33     def compass_azimuth(self):
34         # clockwise from north
35         azimuth = compass.orientation[0]
36         return azimuth
37 
38     @property
39     def ready(self):
40         return self.compass_azimuth is not None and self._ready
41 
42     @mainthread
43     def _on_status(self, stype, status):
44         #self.gps_status = 'type={}\n{}'.format(stype, status)
45         pass
46 
47     @mainthread
48     def _on_location(self, **kwargs):
49         Logger.debug("{} _on_location".format(str(datetime.now())))
50 
51         self.lon = float(kwargs['lon'])
52         self.lat = float(kwargs['lat'])
53         self.alt = float(kwargs['altitude'])
54 
55         self._ready = True
main.py (雑ですまんの)
  1 from datetime import datetime
  2 from math import cos, sin, pi, radians
  3 from kivy.app import App
  4 from kivy.uix.widget import Widget
  5 from kivy.graphics import Color, Line, Ellipse
  6 from kivy.uix.label import Label
  7 from kivy.uix.floatlayout import FloatLayout
  8 from kivy.uix.boxlayout import BoxLayout
  9 from kivy.clock import Clock, mainthread
 10 from kivy.lang import Builder
 11 from kivy.properties import StringProperty, NumericProperty
 12 from plyer.platforms.android import activity
 13 from kivy.logger import Logger
 14 from kivy.config import Config
 15 Config.set('kivy', 'log_level', 'debug')
 16 Config.set('kivy', 'log_enable', 1) 
 17 
 18 from sensor_data import SensorData
 19 from noaa_magnetic_declination import get_magnetic_declination
 20 import noaa_solar_calc
 21 from my_accordion import (
 22     AllCollapsableAccordion, AllCollapsableAccordionItem)
 23 
 24 kv = '''
 25 #:import math math
 26 
 27 <SolarCompassRoot>:
 28     orientation: 'vertical'
 29 
 30     top_pane: top_pane
 31     face: face
 32     extra: extra
 33 
 34     BoxLayout:
 35         size_hint_y: None
 36         height: '48dp'
 37         padding: '4dp'
 38         id: top_pane
 39 
 40         now: now
 41         toggle_fixnorth_button: toggle_fixnorth_button
 42 
 43         Label:
 44             id: now
 45 
 46         ToggleButton:
 47             id: toggle_fixnorth_button
 48             text: 'Toggle fix north'
 49             on_press: root.do_toggle_fix_north()
 50 
 51     SolarCompassFace:
 52         id: face
 53         r: (min(root.size[0], root.size[1] - extra.height - top_pane.height)) / 2
 54 
 55     AllCollapsableAccordion:
 56         id: extra
 57         orientation: 'vertical'
 58 
 59         AllCollapsableAccordionItem:
 60             collapse: True
 61             title: 'Magnetic Declination Config'
 62 
 63             BoxLayout:
 64                 orientation: 'vertical'
 65                 Slider:
 66                     id: linear_pos_slider
 67                     orientation: 'horizontal'
 68                     max: 20
 69                     min: -20
 70                     padding: 1
 71                     value: face.magnetic_declination
 72                     on_value: face.magnetic_declination = self.value
 73                     step: 0.1
 74                 BoxLayout:
 75                     orientation: 'horizontal'
 76                     Button:
 77                         text: 'from NOAA'
 78                         on_release: face.get_magnetic_declination_from_NOAA()
 79                     TextInput:
 80                         id: linear_pos_ti
 81                         font_size: 70
 82                         input_filter: 'float'
 83                         multiline: 'False'
 84                         text: str(face.magnetic_declination)
 85                         on_text: face.magnetic_declination = float(self.text)
 86 '''
 87 Builder.load_string(kv)
 88 
 89 _fix_north = False
 90 
 91 class SolarCompassRoot(BoxLayout):
 92     def __init__(self):
 93         super(SolarCompassRoot, self).__init__()
 94         Clock.schedule_interval(self._draw_now, 1.)
 95 
 96     def do_toggle_fix_north(self):
 97         global _fix_north
 98         _fix_north = not _fix_north
 99 
100     def _draw_now(self, df):
101         self.top_pane.now.text = datetime.now().strftime('%Y-%m-%d %H:%M:%S')
102 
103 class SolarCompassFace(Widget):
104     magnetic_declination = NumericProperty(-7.0)
105 
106     def __init__(self, **kwargs):
107         super(SolarCompassFace, self).__init__(**kwargs)
108         Clock.schedule_interval(self._draw, 2.)
109 
110     def get_magnetic_declination_from_NOAA(self):
111         if not _sensor.ready:
112             return
113         self.magnetic_declination = get_magnetic_declination(
114             _sensor.lon, _sensor.lat)['declination']
115 
116     def _deg2pos(self, deg, rat):
117         pos = (self.center_x + rat * self.r * sin(radians(deg)),
118                self.center_y + rat * self.r * cos(radians(deg)))
119         return pos
120 
121     @mainthread
122     def _draw(self, df):
123         if not _sensor.ready:
124             return
125 
126         # --------------------------------------------
127         if not _fix_north:
128             magnorth_deg = 360 - _sensor.compass_azimuth
129             truenorth_deg = magnorth_deg - self.magnetic_declination
130         else:
131             magnorth_deg = self.magnetic_declination
132             truenorth_deg = 0.0
133 
134         magnorth_pos = self._deg2pos(magnorth_deg, 0.9)
135         truenorth_pos = self._deg2pos(truenorth_deg, 0.9)
136 
137         # --------------------------------------------
138         def _lt_to_disp(lt):
139             return "%02d:%02d" % (int(lt / 60), int(lt % 60))
140 
141         now = datetime.now()
142         #
143         nres = noaa_solar_calc.calculate(
144             _sensor.lon, _sensor.lat,
145             now.year, now.month, now.day,
146             now.hour, now.minute, now.second, 9)
147         #
148         sunrise_deg = nres['sun rise']['Azimuth'] + truenorth_deg
149         sunrize_pos = self._deg2pos(sunrise_deg, 0.8)
150         sunset_deg = nres['sun set']['Azimuth'] + truenorth_deg
151         sunset_pos = self._deg2pos(sunset_deg, 0.8)
152 
153         t = now.hour * 60 + now.minute
154         current_solar_deg = nres['current']['Azimuth'] + truenorth_deg
155         current_solar_pos = self._deg2pos(
156             current_solar_deg, cos(radians(nres['current']['Elevation'])))
157         current_solar_pos2 = self._deg2pos(
158             current_solar_deg, 1.0)
159 
160         self.canvas.clear()
161 
162         with self.canvas:
163             # --- TRUE NORTH ---
164             Color(0.2, 0.2, 0.7)
165             Line(points=[
166                     self.center_x, self.center_y,
167                     truenorth_pos[0], truenorth_pos[1]
168                     ], width=5, cap="round")
169 
170             # --- MAGNETIC NORTH ---
171             Color(0.7, 0.2, 0.2)
172             Line(points=[
173                     self.center_x, self.center_y,
174                     magnorth_pos[0], magnorth_pos[1]
175                     ], width=3, cap="round")
176             label_mn = Label(size_hint_y=None,
177                              text='MN',
178                              pos=self._deg2pos(magnorth_deg, 0.8),
179                              color=(0.7, 0.2, 0.2))
180             label_mn.text_size = label_mn.width * 1.5, None
181             label_mn.height = label_mn.texture_size[1]
182 
183             # --- SUNRISE ---
184             Color(0.7, 0.7, 0.2)
185             Line(points=[
186                     self.center_x, self.center_y,
187                     sunrize_pos[0], sunrize_pos[1]
188                     ], width=5, cap="round")
189             label_sr = Label(size_hint_y=None,
190                              text=_lt_to_disp(nres['sun rise']['LocalTime']),
191                              pos=self._deg2pos(sunrise_deg, 0.8))
192             label_sr.text_size = label_sr.width * 1.5, None
193             label_sr.height = label_sr.texture_size[1]
194 
195             # --- SUNSET ---
196             Color(0.7, 0.2, 0.7)
197             Line(points=[
198                     self.center_x, self.center_y,
199                     sunset_pos[0], sunset_pos[1]
200                     ], width=5, cap="round")
201             label_ss = Label(size_hint_y=None,
202                              text=_lt_to_disp(nres['sun set']['LocalTime']),
203                              pos=self._deg2pos(sunset_deg, 0.8))
204             label_ss.text_size = label_ss.width * 1.5, None
205             label_ss.height = label_ss.texture_size[1]
206 
207             # --- CURRENT SOLAR DIRECTION ---
208             Color(0.7, 0.7, 0.2)
209             Line(points=[
210                     self.center_x, self.center_y,
211                     current_solar_pos[0], current_solar_pos[1]
212                     ], width=3, cap="round")
213             Color(0.6, 0.6, 0.1)
214             Line(points=[
215                     current_solar_pos[0], current_solar_pos[1],
216                     current_solar_pos2[0], current_solar_pos2[1]
217                     ], width=2, cap="round")
218             Color(0.7, 0.7, 0.2)
219             ellips_r = 25
220             Ellipse(pos=[current_solar_pos[0] - ellips_r,
221                          current_solar_pos[1] - ellips_r],
222                     size=[ellips_r * 2, ellips_r * 2])
223             label_cu = Label(
224                 size_hint_y=None,
225                 text=(u"%.1f" % nres['current']['Elevation']),
226                 pos=self._deg2pos(
227                     current_solar_deg,
228                     0.5 * cos(radians(nres['current']['Elevation']))))
229             label_cu.text_size = label_cu.width * 1.5, None
230             label_cu.height = label_cu.texture_size[1]
231 
232 class SolarCompassApp(App):
233     def build(self):
234         w = SolarCompassRoot()
235         return w
236 
237 
238 if __name__ == '__main__':
239     _sensor = SensorData()
240     SolarCompassApp().run()

動かすとこんななの

日の出方位と時刻、日没方位と時刻、現在時刻での太陽の方位と高度(角度ね)、それと True North、Magnetic North ね。現在時刻の太陽の方位の黄色い丸は高度を反映していて、つまり黄色丸は高度が低いほど描画の半径が「長」くなり、高いほど「短」くなります。今は夏なんで結構短くなるし、冬はとっても長くなる、はず。また、北海道と沖縄で較べれば当然北海道の方が「長」くなる。

Solar Calculation 部分は、当たり前なんだけどあくまでも地理座標系に沿った計算であって、磁場の影響は考慮してないわけね。なので、日の出方位、日没方位、現在時刻での太陽の方位は True North 基準。

そしてマグネットコンパスの磁北は真北からは必ずずれる、と。

ジャイロセンサーを駆使して正確な True North を指せないかなぁ、と色々探ってはみたんだけれど、簡単ではなさそうなので、マニュアルで磁偏差を調整出来るようにした:

ダサい GUI ですまんの。あたしゃ気にしてない。

完全に手動で調整してもいいし、「from NOAA」ボタンは NOAA の Magnetic declination calculator から情報を持ってくる。

一番上のボタン「Toggle Fix North」は、True North をデバイスの heading に固定しちゃうモード:

こうなるともはやコンパスではない、だったりもするんだけれど、ちょっと面白い使い方も出来てな。つまり、実際の太陽方位にぴったり合わせることで、逆に True North が割り出せてしまうという。(太陽を直視するのはキツいので、電柱とかの影を利用するといい。) まぁ「太陽の方向」を知るセンサーはないから、方位合わせは至極アナログなわけだけどな。

ちなみに、描画に関する些細な問題があるけどひとまず気にしてない。文字(日の出日没時刻と現在の太陽高度)がね、ときどき消えちゃうんですな。これはまだアタシが kivy に慣れてないから、あんましわかってないんだよねー。

何の役に立つの?

知らんよ、アタシの遊びだもん。

なんつーか、「ちょっと科学的なことを身近に置いておくことで、賢くなった気分になれる」てことかもな。

そもそもさ、「スマホを持ち歩く」って、昔じゃ考えられなかった「色んなセンサーを持ち歩く」ことなわけだね。つまりは「アマチュア科学」が以前よりずっとやりやすいってことなのなー。