138 lines
5.6 KiB
Python
138 lines
5.6 KiB
Python
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 |