Building a declination table

from kanon.units import Sexagesimal
from kanon.units.precision import set_precision, get_context, TruncatureMode
from kanon.tables import HTable
from kanon.tables.htable import join_multiple
import numpy as np
import math

# For these calculations we will round after every operations.

get_context().mutate(tmode=TruncatureMode.ROUND)

OBLIQUITY = "23;51,20"

# We construct a table of sine from 0 to 90.
x = list(Sexagesimal.range(91))
y = [round(Sexagesimal.from_float(np.sin(math.radians(n)), 3)) for n in x]
sin_table = HTable([x, y], names=("Arg", "Val"), index="Arg")

sin_table

# We apply a multiplication on all the sine values to get obl_table = Sin(x) * obl

obl = sin_table.get(Sexagesimal(OBLIQUITY))

obl_table = sin_table.apply("Val", lambda x: x * obl)
obl_table

# Arcsin is sin_table taken from Val to Arg instead of Arg to Val.

arcsin_table = sin_table.copy(set_index="Val")

# For example, if we want arcsin(0.3) :

arcsin_table.get(0.3)

# Finally, we apply arcsin_table.get on every obl_table values to find decl_table = arcsin(obl_table)

decl_table = obl_table.apply("Val", arcsin_table.get)
decl_table

# Let's make the same calculation with an interpolated sine table.

# First, extract every 12th row until the 60th then every 6th.

sin_table_grid = sin_table[
    [i for i in range(91) if i <= 60 and i % 12 == 0 or i > 60 and i % 6 == 0]
]

sin_table_grid

# We can populate it with its interpolation method (linear by default).

sin_table_pop_linear = sin_table_grid.populate(list(Sexagesimal.range(91)), method="interpolate")

sin_table_pop_linear

# You can use other types of interpolation, for example a quadratic interpolation.

from kanon.tables.interpolations import quadratic_interpolation

sin_table_grid.interpolate = quadratic_interpolation

# This interpolation method needs a higher precision to avoid weird results.
with set_precision(pmode=5):
    sin_table_pop_quadratic = sin_table_grid.populate(list(Sexagesimal.range(91)), method="interpolate")

sin_table_pop_quadratic

# Another type of interpolation can come from filling a interpolation grid made from populate
# With this you can use methods interpolating the whole table at once, like euclidean distribution interpolations

sin_table_pop_euclidean = sin_table_grid.populate(list(Sexagesimal.range(91)))
sin_table_pop_euclidean = sin_table_pop_euclidean.fill("distributed_convex")

sin_table_pop_euclidean

# Comparing the 2nd difference between the true sine table and the interpolated ones.

true_2diff = sin_table.diff(n=2, new_name="True Sine 2-diff")
interpolated_linear_2diff = sin_table_pop_linear.diff(n=2, new_name="Linear Sine 2-diff")
interpolated_quadratic_2diff = sin_table_pop_quadratic.diff(n=2, new_name="Quadratic Sine 2-diff")
interpolated_euclidean_2diff = sin_table_pop_euclidean.diff(n=2, new_name="Distributed Sine 2-diff")

# Now we have 4 tables with the same arguments. We can compare them side to side by joining them.

join_multiple(true_2diff, interpolated_linear_2diff, interpolated_quadratic_2diff, interpolated_euclidean_2diff)

# With this newly built table we are able to construct a new declination table.

obl_linear = sin_table_pop_linear.get(Sexagesimal(OBLIQUITY))

obl_table_linear = sin_table_pop_linear.apply("Val", lambda x: x * obl)
arcsin_table_linear = sin_table_pop_linear.copy(set_index="Val")

obl_table_linear.apply("Val", arcsin_table_linear.get)