pingswept / pysolar

Pysolar is a collection of Python libraries for simulating the irradiation of any point on earth by the sun. It includes code for extremely precise ephemeris calculations.

This URL has Read+Write access

pysolar / testsolar.py
100644 142 lines (110 sloc) 7.826 kb
raw blame history 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142

#!/usr/bin/python
 
#    Library for calculating location of the sun
 
#    Copyright 2007 Brandon Stafford
#
#    This file is part of Pysolar.
#
#    Pysolar is free software; you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation; either version 3 of the License, or
#    (at your option) any later version.
#
#    Pysolar is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License along
#    with Pysolar. If not, see <http://www.gnu.org/licenses/>.
 
import solar
import constants
import datetime
import unittest
 
class testSolar(unittest.TestCase):
 
	def setUp(self):
		self.d = datetime.datetime(2003, 10, 17, 19, 30, 30)
		self.longitude = -105.1786
		self.latitude = 39.742476
		self.pressure = 820.0 # millibars
		self.elevation = 1830.14 # meters
		self.temperature = 11.0 # degrees Celsius
		self.slope = 30.0 # degrees
		self.slope_orientation = -10.0 # degrees east from south
		self.jd = solar.GetJulianDay(self.d)
		self.jc = solar.GetJulianCentury(self.jd)
		self.jde = solar.GetJulianEphemerisDay(self.jd, 67.0) #65.284)
		self.jce = solar.GetJulianEphemerisCentury(self.jde)
		self.jme = solar.GetJulianEphemerisMillenium(self.jce)
		self.geocentric_longitude = solar.GetGeocentricLongitude(self.jme)
		self.geocentric_latitude = solar.GetGeocentricLatitude(self.jme)
		self.nutation = solar.GetNutation(self.jde)
		self.radius_vector = solar.GetRadiusVector(self.jme)
		self.true_ecliptic_obliquity = solar.GetTrueEclipticObliquity(self.jme, self.nutation)
		self.aberration_correction = solar.GetAberrationCorrection(self.radius_vector)
		self.apparent_sun_longitude = solar.GetApparentSunLongitude(self.geocentric_longitude, self.nutation, self.aberration_correction)
		self.apparent_sidereal_time = solar.GetApparentSiderealTime(self.jd, self.jme, self.nutation)
		self.geocentric_sun_right_ascension = solar.GetGeocentricSunRightAscension(self.apparent_sun_longitude, self.true_ecliptic_obliquity, self.geocentric_latitude)
		self.geocentric_sun_declination = solar.GetGeocentricSunDeclination(self.apparent_sun_longitude, self.true_ecliptic_obliquity, self.geocentric_latitude)
		self.local_hour_angle = solar.GetLocalHourAngle(318.5119, self.longitude, self.geocentric_sun_right_ascension) #self.apparent_sidereal_time only correct to 5 sig figs, so override
		self.equatorial_horizontal_parallax = solar.GetEquatorialHorizontalParallax(self.radius_vector)
		self.projected_radial_distance = solar.GetProjectedRadialDistance(self.elevation, self.latitude)
		self.projected_axial_distance = solar.GetProjectedAxialDistance(self.elevation, self.latitude)
		self.topocentric_sun_right_ascension = solar.GetTopocentricSunRightAscension(self.projected_radial_distance,
		    self.equatorial_horizontal_parallax, self.local_hour_angle, self.projected_axial_distance, self.apparent_sun_longitude, self.true_ecliptic_obliquity, self.geocentric_latitude)
		self.parallax_sun_right_ascension = solar.GetParallaxSunRightAscension(self.projected_radial_distance, self.equatorial_horizontal_parallax, self.local_hour_angle, self.geocentric_sun_declination, self.projected_axial_distance)
		self.topocentric_sun_declination = solar.GetTopocentricSunDeclination(self.geocentric_sun_declination, self.projected_axial_distance, self.equatorial_horizontal_parallax, self.parallax_sun_right_ascension, self.local_hour_angle)
		self.topocentric_local_hour_angle = solar.GetTopocentricLocalHourAngle(self.local_hour_angle, self.parallax_sun_right_ascension)
		self.topocentric_zenith_angle = solar.GetTopocentricZenithAngle(self.latitude, self.topocentric_sun_declination, self.topocentric_local_hour_angle, self.pressure, self.temperature)
		self.topocentric_azimuth_angle = solar.GetTopocentricAzimuthAngle(self.topocentric_local_hour_angle, self.latitude, self.topocentric_sun_declination)
		self.incidence_angle = solar.GetIncidenceAngle(self.topocentric_zenith_angle, self.slope, self.slope_orientation, self.topocentric_azimuth_angle)
 
	def testGetJulianDay(self):
		self.assertAlmostEqual(2452930.312847, self.jd, 6) # value from Reda and Andreas (2005)
 
	def testGetJulianEphemerisDay(self):
		self.assertAlmostEqual(2452930.3136, self.jde, 4) # value not validated
 
	def testGetJulianCentury(self):
		self.assertAlmostEqual(0.03792779869191517, self.jc, 12) # value not validated
 
	def testGetJulianEphemerisMillenium(self):
		self.assertAlmostEqual(0.0037927819922933584, self.jme, 12) # value not validated
 
	def testGetGeocentricLongitude(self):
		self.assertAlmostEqual(204.0182635175, self.geocentric_longitude, 4) # value from Reda and Andreas (2005)
 
	def testGetGeocentricLatitude(self):
		self.assertAlmostEqual(0.0001011219, self.geocentric_latitude, 9) # value from Reda and Andreas (2005)
 
	def testGetNutation(self):
		self.assertAlmostEqual(0.00166657, self.nutation['obliquity'], 8) # value from Reda and Andreas (2005)
		self.assertAlmostEqual(-0.00399840, self.nutation['longitude'], 8) # value from Reda and Andreas (2005)
 
	def testGetRadiusVector(self):
		self.assertAlmostEqual(0.9965421031, self.radius_vector, 7) # value from Reda and Andreas (2005)
 
	def testGetTrueEclipticObliquity(self):
		self.assertAlmostEqual(23.440465, self.true_ecliptic_obliquity, 6) # value from Reda and Andreas (2005)
 
	def testGetAberrationCorrection(self):
		self.assertAlmostEqual(-0.0057113603, self.aberration_correction, 9) # value not validated
 
	def testGetApparentSunLongitude(self):
		self.assertAlmostEqual(204.0085537528, self.apparent_sun_longitude, 4) # value from Reda and Andreas (2005)
 
	def testGetApparentSiderealTime(self):
		self.assertAlmostEqual(318.5119, self.apparent_sidereal_time, 2) # value derived from Reda and Andreas (2005)
 
	def testGetGeocentricSunRightAscension(self):
		self.assertAlmostEqual(202.22741, self.geocentric_sun_right_ascension, 4) # value from Reda and Andreas (2005)
 
	def testGetGeocentricSunDeclination(self):
		self.assertAlmostEqual(-9.31434, self.geocentric_sun_declination, 4) # value from Reda and Andreas (2005)
 
	def testGetLocalHourAngle(self):
		self.assertAlmostEqual(11.105900, self.local_hour_angle, 4) # value from Reda and Andreas (2005)
 
	def testGetProjectedRadialDistance(self):
		self.assertAlmostEqual(0.7702006, self.projected_radial_distance, 6) # value not validated
 
	def testGetTopocentricSunRightAscension(self):
		self.assertAlmostEqual(202.22741, self.topocentric_sun_right_ascension, 3) # value from Reda and Andreas (2005)
 
	def testGetParallaxSunRightAscension(self):
		self.assertAlmostEqual(-0.00036598821845849395, self.parallax_sun_right_ascension, 12) # value not validated
 
	def testGetTopocentricSunDeclination(self):
		self.assertAlmostEqual(-9.316179, self.topocentric_sun_declination, 3) # value from Reda and Andreas (2005)
 
	def testGetTopocentricLocalHourAngle(self):
		self.assertAlmostEqual(11.10629, self.topocentric_local_hour_angle, 4) # value from Reda and Andreas (2005)
 
	def testGetTopocentricZenithAngle(self):
		self.assertAlmostEqual(50.11162, self.topocentric_zenith_angle, 3) # value from Reda and Andreas (2005)
 
	def testGetTopocentricAzimuthAngle(self):
		self.assertAlmostEqual(194.34024, self.topocentric_azimuth_angle, 3) # value from Reda and Andreas (2005)
 
	def testGetIncidenceAngle(self):
		self.assertAlmostEqual(25.18700, self.incidence_angle, 3) # value from Reda and Andreas (2005)
 
suite = unittest.TestLoader().loadTestsFromTestCase(testSolar)
unittest.TextTestRunner(verbosity=2).run(suite)
 
# if __name__ == "__main__":
#	unittest.main()