import machine
import time

def get_location_and_time():
    # Assuming UART connection to the GPS module at pins GP0 (TX) and GP1 (RX)
    uart = machine.UART(0, baudrate=9600)

    latitude = None
    longitude = None
    datetime = None
    date = None

    # Give the GPS some time to get a fix
    while True:
        time.sleep(1)
        if uart.any():
            line = uart.readline()
            try:
                # Convert bytes to string and strip any trailing whitespace
                line = line.decode('ascii').strip()
                # NMEA sentence for date and time data starts with '$GPRMC'
                if line.startswith('$GPRMC'):
                    data = line.split(',')
                    # Make sure we have a valid fix
                    if data[2] == 'A':  # Data[2] is the status A=active, V=Void.
                        # Parse latitude and longitude with direction
                        latitude = convert_to_degrees(data[3], data[4])
                        longitude = convert_to_degrees(data[5], data[6])
                        # Parse time (hhmmss.ss)
                        time_str = data[1]
                        hour = int(time_str[0:2])
                        minute = int(time_str[2:4])
                        second = int(time_str[4:6])
                        # Parse date (ddmmyy)
                        date_str = data[9]
                        day = int(date_str[0:2])
                        month = int(date_str[2:4])
                        year = int(date_str[4:6]) + 2000  # Adjust for Y2K (assumes a 2000 baseline)
                        datetime = (year, month, day, hour, minute, second)
                        break
            except ValueError:
                # Handle invalid data or incomplete parsing
                pass

    return latitude, longitude, datetime

def convert_to_degrees(raw_value, direction):
    # Split the raw_value into degrees and minutes
    decimal_point_index = raw_value.find('.')  # Find the decimal point position
    degrees = int(raw_value[:decimal_point_index-2])  # Degrees are everything before the minutes
    minutes = float(raw_value[decimal_point_index-2:])  # Minutes are after the degrees
    
    # Convert to decimal degrees
    decimal_degrees = degrees + (minutes/60)
    
    # Check the direction for latitude (N/S) and longitude (E/W)
    if direction == 'S' or direction == 'W':
        decimal_degrees = -decimal_degrees  # South and West are negative
    
    return decimal_degrees

if __name__=="__main__" :
    latitude, longitude, current_time = get_location_and_time()
    print("Latitude:", latitude)
    print("Longitude:", longitude)
    print("Current Time:", current_time)

from gps import get_location_and_time
from stepper import Stepper
from sunrise import improved_sunrise_sunset
import machine
import utime

# Initialize the stepper motor
stepper = Stepper(2, 3, 4, 5, power=False)

# Function to get the number of steps to move based on the current time
def get_steps_to_move(current_time, sunrise_time, sunset_time, total_steps_per_cycle):
    # Calculate the current time in minutes since midnight
    current_minutes = current_time[4] * 60 + current_time[5]
    
    # Calculate minutes for sunrise and sunset, handling the case where they may be on different UTC days
    sunrise_minutes = sunrise_time[0] * 60 + sunrise_time[1]
    sunset_minutes = sunset_time[0] * 60 + sunset_time[1]
    
    # If sunset is earlier in the day than sunrise, it means sunset occurs on the next day
    if sunset_minutes < sunrise_minutes:
        sunset_minutes += 24 * 60
    
    # Check if current time is past midnight but sunset was before midnight
    if current_minutes < sunrise_minutes and sunset_minutes > (24 * 60):
        current_minutes += 24 * 60
    
    # Calculate the ratio of the current time's position within the current day or night cycle
    if current_minutes < sunrise_minutes:
        # Before sunrise (counting from previous sunset to this sunrise)
        cycle_ratio = (current_minutes + (24 * 60 - sunset_minutes)) / (sunrise_minutes + (24 * 60 - sunset_minutes))
    elif current_minutes < sunset_minutes:
        # After sunrise and before sunset (daytime)
        cycle_ratio = (current_minutes - sunrise_minutes) / (sunset_minutes - sunrise_minutes)
    else:
        # After sunset (nighttime)
        cycle_ratio = (current_minutes - sunset_minutes) / ((24 * 60 - sunset_minutes) + sunrise_minutes)
    
    # Convert the ratio to steps
    return int(cycle_ratio * total_steps_per_cycle)

# Initialize the RTC with GPS time
rtc = machine.RTC()
latitude, longitude, current_utc_datetime = get_location_and_time()
if latitude is not None and longitude is not None and current_utc_datetime is not None:
    rtc.datetime((current_utc_datetime[0], current_utc_datetime[1], current_utc_datetime[2], 0, current_utc_datetime[3], current_utc_datetime[4], current_utc_datetime[5], 0))
    sunrise_time, sunset_time = improved_sunrise_sunset(latitude, longitude, current_utc_datetime[0], current_utc_datetime[1], current_utc_datetime[2])

# Define the total steps in a 12-hour cycle (half a day), assuming 4096 steps per revolution and 12 revolutions per cycle
total_steps_per_half_day = 12 * 4096

# Main loop
last_sync_day = -1
while True:
    # Synchronize RTC with GPS once per day
    current_time = rtc.datetime()
    if current_time[2] != last_sync_day:
        latitude, longitude, current_utc_datetime = get_location_and_time()
        if latitude is not None and longitude is not None and current_utc_datetime is not None:
            rtc.datetime((current_utc_datetime[0], current_utc_datetime[1], current_utc_datetime[2], 0, current_utc_datetime[3], current_utc_datetime[4], current_utc_datetime[5], 0))
            sunrise_time, sunset_time = improved_sunrise_sunset(latitude, longitude, current_utc_datetime[0], current_utc_datetime[1], current_utc_datetime[2])
            last_sync_day = current_time[2]
    
    # Update the stepper motor position
    steps_to_move = get_steps_to_move(current_time, sunrise_time, sunset_time, total_steps_per_half_day)
    stepper.step(steps_to_move - stepper.current_position)
    stepper.current_position = steps_to_move
    
    # Sleep the stepper motor to save power
    stepper.sleep()
    
    # Sleep for a short time before next loop iteration
    utime.sleep(1)

from machine import Pin
from time import sleep

class Stepper:
    def __init__(self, pin1, pin2, pin3, pin4, power):
        self.IN1 = Pin(pin1,Pin.OUT)
        self.IN2 = Pin(pin2,Pin.OUT)
        self.IN3 = Pin(pin3,Pin.OUT)
        self.IN4 = Pin(pin4,Pin.OUT)
        self.pins = [self.IN1, self.IN2, self.IN3, self.IN4]
        if power:
            self.sequence = [[1, 1, 0, 0], [0, 1, 1, 0], [0, 0, 1, 1], [1, 0, 0, 1]]
        else:
            self.sequence = [[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1]]
    
    def step(self, steps):
        direction = 1 if steps > 0 else -1
        for _ in range(abs(steps)):
            for step in self.sequence[::direction]:
                for i, pin in enumerate(self.pins):
                    pin.value(step[i])
                    sleep(0.001)
                    
    def sleep(self):
        for i in range(len(self.pins)):
            self.pins[i].value(0)

import math

# Define the improved sunrise and sunset function
def improved_sunrise_sunset(latitude, longitude, year, month, day):
    # Constants
    DEG_TO_RAD = math.pi / 180.0
    RAD_TO_DEG = 180.0 / math.pi
    
    # Calculate Julian Date
    # The day fraction (time of the day) is set to 0.5, which corresponds to 12:00 UT
    if month <= 2:
        year -= 1
        month += 12
    A = math.floor(year / 100)
    B = 2 - A + math.floor(A / 4)
    julian_date = math.floor(365.25 * (year + 4716)) + math.floor(30.6001 * (month + 1)) + day + B - 1524.5 + 0.5
    
    # Calculate the Julian Century
    julian_century = (julian_date - 2451545.0) / 36525.0
    
    # Calculate the Geom Mean Long Sun (deg)
    geom_mean_long_sun = (280.46646 + julian_century * (36000.76983 + julian_century * 0.0003032)) % 360
    
    # Calculate the Geom Mean Anom Sun (deg)
    geom_mean_anom_sun = 357.52911 + julian_century * (35999.05029 - 0.0001537 * julian_century)
    
    # Calculate Eccentricity of Earth's Orbit
    eccent_earth_orbit = 0.016708634 - julian_century * (0.000042037 + 0.0000001267 * julian_century)
    
    # Calculate Sun Eq of Ctr
    sun_eq_of_ctr = math.sin(DEG_TO_RAD * geom_mean_anom_sun) * (1.914602 - julian_century * (0.004817 + 0.000014 * julian_century)) + math.sin(DEG_TO_RAD * (2 * geom_mean_anom_sun)) * (0.019993 - 0.000101 * julian_century) + math.sin(DEG_TO_RAD * (3 * geom_mean_anom_sun)) * 0.000289
    
    # Calculate Sun True Long (deg)
    sun_true_long = geom_mean_long_sun + sun_eq_of_ctr
    
    # Calculate Sun True Anom (deg)
    sun_true_anom = geom_mean_anom_sun + sun_eq_of_ctr
    
    # Calculate Sun Rad Vector (AUs)
    sun_rad_vector = (1.000001018 * (1 - eccent_earth_orbit * eccent_earth_orbit)) / (1 + eccent_earth_orbit * math.cos(DEG_TO_RAD * sun_true_anom))
    
    # Calculate Sun App Long (deg)
    sun_app_long = sun_true_long - 0.00569 - 0.00478 * math.sin(DEG_TO_RAD * (125.04 - 1934.136 * julian_century))
    
    # Calculate Mean Obliq Ecliptic (deg)
    mean_obliq_ecliptic = 23 + (26 + ((21.448 - julian_century * (46.815 + julian_century * (0.00059 - julian_century * 0.001813)))) / 60) / 60
    
    # Calculate Obliq Corr (deg)
    obliq_corr = mean_obliq_ecliptic + 0.00256 * math.cos(DEG_TO_RAD * (125.04 - 1934.136 * julian_century))
    
    # Calculate Sun Declination (deg)
    sun_declination = RAD_TO_DEG * math.asin(math.sin(DEG_TO_RAD * obliq_corr) * math.sin(DEG_TO_RAD * sun_app_long))
    
    # Calculate var y
    var_y = math.tan(DEG_TO_RAD * (obliq_corr / 2)) * math.tan(DEG_TO_RAD * (obliq_corr / 2))
    
    # Calculate Eq of Time (minutes)
    eq_of_time = 4 * RAD_TO_DEG * (var_y * math.sin(2 * DEG_TO_RAD * geom_mean_long_sun) - 2 * eccent_earth_orbit * math.sin(DEG_TO_RAD * geom_mean_anom_sun) + 4 * eccent_earth_orbit * var_y * math.sin(DEG_TO_RAD * geom_mean_anom_sun) * math.cos(2 * DEG_TO_RAD * geom_mean_long_sun) - 0.5 * var_y * var_y * math.sin(4 * DEG_TO_RAD * geom_mean_long_sun) - 1.25 * eccent_earth_orbit * eccent_earth_orbit * math.sin(DEG_TO_RAD * (2 * geom_mean_anom_sun)))
    
    # Calculate Solar Noon (LST)
    solar_noon_LST = (720 - 4 * longitude - eq_of_time + (5 * 60)) / 1440  # (5 * 60) is the time zone offset for EST (UTC-5)
    
    # Calculate Sunrise and Sunset (LST)
    # Solar zenith angle is 90.833 degrees (the solar limb is at the horizon at sunrise and sunset)
    solar_zenith_angle = 90.833
    cos_h = (math.cos(DEG_TO_RAD * solar_zenith_angle) - (math.sin(DEG_TO_RAD * latitude) * math.sin(DEG_TO_RAD * sun_declination))) / (math.cos(DEG_TO_RAD * latitude) * math.cos(DEG_TO_RAD * sun_declination))
    if cos_h > 1:
        return None, None  # Sun never rises on this location (polar night)
    elif cos_h < -1:
        return None, None  # Sun never sets on this location (midnight sun)
    else:
        h = RAD_TO_DEG * math.acos(cos_h)  # Hour angle in degrees
        sunrise_LST = solar_noon_LST - h / 360
        sunset_LST = solar_noon_LST + h / 360
        
        # Convert LST to UTC (for EST time zone)
        sunrise_UTC = sunrise_LST * 1440 - (5 * 60)
        sunset_UTC = sunset_LST * 1440 - (5 * 60)
        
        # Convert decimal time to hours, minutes, and seconds
        def decimal_time_to_hms(decimal_time):
            hours = int(decimal_time // 60)
            minutes = int((decimal_time % 60))
            seconds = int((decimal_time * 60) % 60)
            return (hours, minutes, seconds)
        
        return decimal_time_to_hms(sunrise_UTC), decimal_time_to_hms(sunset_UTC)

# Test the improved function with the provided parameters
if __name__=="__main__" :
    improved_sunrise, improved_sunset = sunrise_sunset(42.73002, -71.48132, 2023, 11, 8)
    print(improved_sunrise, improved_sunset)