Computing the true position of the Sun¶
[1]:
# Imports
import numpy as np
from kanon.calendars import Calendar, Date
from kanon.tables import HTable
from kanon.units import Sexagesimal
from kanon.units.precision import set_precision, TruncatureMode, PrecisionMode
import astropy.units as u
[2]:
# 3rd of July 1327
year = 1327
month = 7
day = 3
[3]:
# Day computation
# We're using the Julian A.D. calendar
calendar = Calendar.registry["Julian A.D."]
date = Date(calendar, (year, month, day))
# We need to express the numbers of days from the start of the calendar
# in Sexagesimal representation.
days = Sexagesimal.from_float(date.days_from_epoch(), 0)
days
[3]:
02,14,35,04 ;
[4]:
# Defining a function resolving a position from a mean motion table
def position_from_table(ndays, tab, radix, zodiac_offset=4, width=9):
# Starting with the radix
result = radix
# Adding days
with set_precision(pmode=PrecisionMode.MAX):
for i, v in enumerate(ndays[:]):
result += tab.get(v) >> i + zodiac_offset
# Conversion in degrees modulo 6 zodiacal signs
result %= Sexagesimal(6, 0) * u.degree
return result
[5]:
# Sun mean position
# Reading the table from DISHAS
tab_mean_motion = HTable.read(193, format="dishas")
tab_mean_motion
[5]:
HTable length=60
Days | Entries |
---|---|
d | deg |
Sexagesimal | Sexagesimal |
01 ; | 59,08,19,37,19,13,56 ; |
02 ; | 01,58,16,39,14,38,27,52 ; |
03 ; | 02,57,24,58,51,57,41,48 ; |
04 ; | 03,56,33,18,29,16,55,44 ; |
05 ; | 04,55,41,38,06,36,09,40 ; |
06 ; | 05,54,49,57,43,55,23,36 ; |
07 ; | 06,53,58,17,21,14,37,32 ; |
08 ; | 07,53,06,36,58,33,51,28 ; |
09 ; | 08,52,14,56,35,53,05,24 ; |
... | ... |
51 ; | 50,16,04,40,43,20,50,36 ; |
52 ; | 51,15,13,00,20,40,04,32 ; |
53 ; | 52,14,21,19,57,59,18,28 ; |
54 ; | 53,13,29,39,35,18,32,24 ; |
55 ; | 54,12,37,59,12,37,46,20 ; |
56 ; | 55,11,46,18,49,57,00,16 ; |
57 ; | 56,10,54,38,27,16,14,12 ; |
58 ; | 57,10,02,58,04,35,28,08 ; |
59 ; | 58,09,11,17,41,54,42,04 ; |
01,00 ; | 59,08,19,37,19,13,56,00 ; |
[6]:
# Mean position from days, table, and radix
mean_sun_pos = position_from_table(days, tab_mean_motion, Sexagesimal("4,38;21,0,30,28") * u.degree)
mean_sun_pos
[6]:
$01,47 ; 58,21,32,17,28,35,44 \; \mathrm{{}^{\circ}}$
[7]:
# Fixed stars
tab_fixed_stars = HTable.read(236, format="dishas")
mean_fixed_star_pos = position_from_table(days, tab_fixed_stars, 0)
mean_fixed_star_pos
[7]:
$09 ; 44,44,33,52,35,08,48 \; \mathrm{{}^{\circ}}$
[8]:
# Access and recess position
tab_access_recess = HTable.read(237, format="dishas")
access_recess_pos = position_from_table(days, tab_access_recess, Sexagesimal("5,59;12,34")*u.degree, zodiac_offset=2, width=7)
access_recess_pos
[8]:
$01,07 ; 25,45,56,14,16 \; \mathrm{{}^{\circ}}$
[9]:
# Access and recess equation
tab_eq_access_recess = HTable.read(238, format="dishas")
with set_precision(pmode=PrecisionMode.MAX, tmode=TruncatureMode.ROUND):
eq_access_recess = tab_eq_access_recess.get(access_recess_pos.value)
eq_access_recess
[9]:
$08 ; 18,18,36,54,20 \; \mathrm{{}^{\circ}}$
[10]:
# Sun apogy
with set_precision(pmode=PrecisionMode.MAX, tmode=TruncatureMode.ROUND):
solar_apogee_pos = mean_fixed_star_pos + eq_access_recess + Sexagesimal("1,11;25,23") * u.degree
solar_apogee_pos
[10]:
$01,29 ; 28,26,10,46,55,08,48 \; \mathrm{{}^{\circ}}$
[11]:
# Sun mean argument
with set_precision(pmode=PrecisionMode.MAX, tmode=TruncatureMode.ROUND):
mean_arg_sun = mean_sun_pos + (Sexagesimal(6,0) * u.degree if mean_sun_pos < solar_apogee_pos else 0) - solar_apogee_pos
mean_arg_sun
[11]:
$18 ; 29,55,21,30,33,26,56 \; \mathrm{{}^{\circ}}$
[12]:
# Sun equation
tab_eq_sun = HTable.read(19, format="dishas")
with set_precision(pmode=PrecisionMode.MAX, tmode=TruncatureMode.ROUND):
eq_sun = tab_eq_sun.get(mean_arg_sun.value)
eq_sun
[12]:
$00 ; 39,33,20,10,31,40,48 \; \mathrm{{}^{\circ}}$
[13]:
# Sun true position
with set_precision(pmode=2, tmode=TruncatureMode.ROUND):
true_pos_sun = mean_sun_pos + (Sexagesimal(6,0) * u.degree if mean_sun_pos < eq_sun else 0) - eq_sun
true_pos_sun
[13]:
$01,47 ; 18,49 \; \mathrm{{}^{\circ}}$