import math import pandas as pd from datetime import datetime from pathlib import Path class CitySolarAngles: def __init__(self, file_name, location_latitude, location_longitude, tilt_angle, surface_azimuth_angle, standard_meridian=-75): self.file_name = file_name self.location_latitude = location_latitude self.location_longitude = location_longitude self.location_latitude_rad = math.radians(location_latitude) self.surface_azimuth_angle = surface_azimuth_angle self.surface_azimuth_rad = math.radians(surface_azimuth_angle) self.tilt_angle = tilt_angle self.tilt_angle_rad = math.radians(tilt_angle) self.standard_meridian = standard_meridian self.longitude_correction = (location_longitude - standard_meridian) * 4 self.timezone = 'Etc/GMT+5' self.eot = [] self.ast = [] self.hour_angles = [] self.declinations = [] self.solar_altitudes = [] self.solar_azimuths = [] self.zeniths = [] self.incidents = [] self.beam_tilted = [] self.factor = [] self.times = pd.date_range(start='2023-01-01', end='2024-01-01', freq='H', tz=self.timezone) self.df = pd.DataFrame(index=self.times) self.day_of_year = self.df.index.dayofyear def solar_time(self, datetime_val, day_of_year): b = (day_of_year - 81) * 2 * math.pi / 364 eot = 9.87 * math.sin(2 * b) - 7.53 * math.cos(b) - 1.5 * math.sin(b) self.eot.append(eot) # Calculate Local Solar Time (LST) lst_hour = datetime_val.hour lst_minute = datetime_val.minute lst_second = datetime_val.second lst = lst_hour + lst_minute / 60 + lst_second / 3600 # Calculate Apparent Solar Time (AST) in decimal hours ast_decimal = lst + eot / 60 + self.longitude_correction / 60 ast_hours = int(ast_decimal) ast_minutes = round((ast_decimal - ast_hours) * 60) # Ensure ast_minutes is within valid range if ast_minutes == 60: ast_hours += 1 ast_minutes = 0 elif ast_minutes < 0: ast_minutes = 0 ast_time = datetime(year=datetime_val.year, month=datetime_val.month, day=datetime_val.day, hour=ast_hours, minute=ast_minutes) self.ast.append(ast_time) return ast_time def declination_angle(self, day_of_year): declination = 23.45 * math.sin(math.radians(360 / 365 * (284 + day_of_year))) declination_radian = math.radians(declination) self.declinations.append(declination) return declination_radian def hour_angle(self, ast_time): hour_angle = ((ast_time.hour * 60 + ast_time.minute) - 720) / 4 hour_angle_radian = math.radians(hour_angle) self.hour_angles.append(hour_angle) return hour_angle_radian def solar_altitude(self, declination_radian, hour_angle_radian): solar_altitude_radians = math.asin(math.cos(self.location_latitude_rad) * math.cos(declination_radian) * math.cos(hour_angle_radian) + math.sin(self.location_latitude_rad) * math.sin(declination_radian)) solar_altitude = math.degrees(solar_altitude_radians) self.solar_altitudes.append(solar_altitude) return solar_altitude_radians def zenith(self, solar_altitude_radians): solar_altitude = math.degrees(solar_altitude_radians) zenith_degree = 90 - solar_altitude zenith_radian = math.radians(zenith_degree) self.zeniths.append(zenith_degree) return zenith_radian def solar_azimuth_analytical(self, hourangle, declination, zenith): numer = (math.cos(zenith) * math.sin(self.location_latitude_rad) - math.sin(declination)) denom = (math.sin(zenith) * math.cos(self.location_latitude_rad)) if math.isclose(denom, 0.0, abs_tol=1e-8): cos_azi = 1.0 else: cos_azi = numer / denom cos_azi = max(-1.0, min(1.0, cos_azi)) sign_ha = math.copysign(1, hourangle) solar_azimuth_radians = sign_ha * math.acos(cos_azi) + math.pi solar_azimuth_degrees = math.degrees(solar_azimuth_radians) self.solar_azimuths.append(solar_azimuth_degrees) return solar_azimuth_radians def incident_angle(self, solar_altitude_radians, solar_azimuth_radians): incident_radian = math.acos(math.cos(solar_altitude_radians) * math.cos(abs(solar_azimuth_radians - self.surface_azimuth_rad)) * math.sin(self.tilt_angle_rad) + math.sin(solar_altitude_radians) * math.cos(self.tilt_angle_rad)) incident_angle_degrees = math.degrees(incident_radian) self.incidents.append(incident_angle_degrees) return incident_radian @property def calculate(self) -> pd.DataFrame: for i in range(len(self.times)): datetime_val = self.times[i] day_of_year = self.day_of_year[i] declination_radians = self.declination_angle(day_of_year) ast_time = self.solar_time(datetime_val, day_of_year) hour_angle_radians = self.hour_angle(ast_time) solar_altitude_radians = self.solar_altitude(declination_radians, hour_angle_radians) zenith_radians = self.zenith(solar_altitude_radians) solar_azimuth_radians = self.solar_azimuth_analytical(hour_angle_radians, declination_radians, zenith_radians) incident_angle_radian = self.incident_angle(solar_altitude_radians, solar_azimuth_radians) self.df['DateTime'] = self.times self.df['AST'] = self.ast self.df['hour angle'] = self.hour_angles self.df['eot'] = self.eot self.df['declination angle'] = self.declinations self.df['solar altitude'] = self.solar_altitudes self.df['zenith'] = self.zeniths self.df['solar azimuth'] = self.solar_azimuths self.df['incident angle'] = self.incidents return self.df