1 # SPDX-License-Identifier: GPL-2.0-or-later
4 from bpy
.app
.handlers
import persistent
5 from mathutils
import Euler
7 from math
import degrees
, radians
, pi
9 from .geo
import parse_position
12 ############################################################################
14 # SunClass is used for storing intermediate sun calculations.
16 ############################################################################
52 use_daylight_savings
= False
58 def sun_update(self
, context
):
62 def parse_coordinates(self
, context
):
63 error_message
= "ERROR: Could not parse coordinates"
64 sun_props
= context
.scene
.sun_pos_properties
66 if sun_props
.co_parser
:
67 parsed_co
= parse_position(sun_props
.co_parser
)
69 if parsed_co
is not None and len(parsed_co
) == 2:
70 sun_props
.latitude
, sun_props
.longitude
= parsed_co
71 elif sun_props
.co_parser
!= error_message
:
72 sun_props
.co_parser
= error_message
75 if sun_props
.co_parser
not in {'', error_message
}:
76 sun_props
.co_parser
= ''
79 def sun_handler(scene
):
80 update_time(bpy
.context
)
84 ############################################################################
86 # move_sun() will cycle through all the selected objects
87 # and call set_sun_position and set_sun_rotations
88 # to place them in the sky.
90 ############################################################################
93 def move_sun(context
):
94 addon_prefs
= context
.preferences
.addons
[__package__
].preferences
95 sun_props
= context
.scene
.sun_pos_properties
97 if sun_props
.usage_mode
== "HDR":
98 nt
= context
.scene
.world
.node_tree
.nodes
99 env_tex
= nt
.get(sun_props
.hdr_texture
)
101 if sun
.bind_to_sun
!= sun_props
.bind_to_sun
:
102 # bind_to_sun was just toggled
103 sun
.bind_to_sun
= sun_props
.bind_to_sun
104 sun
.bind
.az_start_sun
= sun_props
.hdr_azimuth
106 sun
.bind
.az_start_env
= env_tex
.texture_mapping
.rotation
.z
108 if env_tex
and sun_props
.bind_to_sun
:
109 az
= sun_props
.hdr_azimuth
- sun
.bind
.az_start_sun
+ sun
.bind
.az_start_env
110 env_tex
.texture_mapping
.rotation
.z
= az
112 if sun_props
.sun_object
:
113 sun
.theta
= math
.pi
/ 2 - sun_props
.hdr_elevation
114 sun
.phi
= -sun_props
.hdr_azimuth
116 obj
= sun_props
.sun_object
117 set_sun_position(obj
, sun_props
.sun_distance
)
118 rotation_euler
= Euler((sun_props
.hdr_elevation
- pi
/2,
119 0, -sun_props
.hdr_azimuth
))
121 set_sun_rotations(obj
, rotation_euler
)
124 local_time
= sun_props
.time
125 zone
= -sun_props
.UTC_zone
126 sun
.use_daylight_savings
= sun_props
.use_daylight_savings
127 if sun
.use_daylight_savings
:
130 north_offset
= degrees(sun_props
.north_offset
)
132 if addon_prefs
.show_rise_set
:
133 calc_sunrise_sunset(rise
=True)
134 calc_sunrise_sunset(rise
=False)
136 get_sun_position(local_time
, sun_props
.latitude
, sun_props
.longitude
,
137 north_offset
, zone
, sun_props
.month
, sun_props
.day
, sun_props
.year
,
138 sun_props
.sun_distance
)
140 if sun_props
.sky_texture
:
141 sky_node
= bpy
.context
.scene
.world
.node_tree
.nodes
.get(sun_props
.sky_texture
)
142 if sky_node
is not None and sky_node
.type == "TEX_SKY":
143 locX
= math
.sin(sun
.phi
) * math
.sin(-sun
.theta
)
144 locY
= math
.sin(sun
.theta
) * math
.cos(sun
.phi
)
145 locZ
= math
.cos(sun
.theta
)
146 sky_node
.texture_mapping
.rotation
.z
= 0.0
147 sky_node
.sun_direction
= locX
, locY
, locZ
148 sky_node
.sun_elevation
= math
.radians(sun
.elevation
)
149 sky_node
.sun_rotation
= math
.radians(sun
.az_north
)
152 if (sun_props
.sun_object
is not None
153 and sun_props
.sun_object
.name
in context
.view_layer
.objects
):
154 obj
= sun_props
.sun_object
155 set_sun_position(obj
, sun_props
.sun_distance
)
156 rotation_euler
= Euler((math
.radians(sun
.elevation
- 90), 0,
157 math
.radians(-sun
.az_north
)))
158 set_sun_rotations(obj
, rotation_euler
)
161 if sun_props
.object_collection
is not None:
162 sun_objects
= sun_props
.object_collection
.objects
163 object_count
= len(sun_objects
)
164 if sun_props
.object_collection_type
== 'DIURNAL':
167 time_increment
= sun_props
.time_spread
/ (object_count
- 1)
168 local_time
= local_time
+ time_increment
* (object_count
- 1)
170 time_increment
= sun_props
.time_spread
172 for obj
in sun_objects
:
173 get_sun_position(local_time
, sun_props
.latitude
,
174 sun_props
.longitude
, north_offset
, zone
,
175 sun_props
.month
, sun_props
.day
,
176 sun_props
.year
, sun_props
.sun_distance
)
177 set_sun_position(obj
, sun_props
.sun_distance
)
178 local_time
-= time_increment
179 obj
.rotation_euler
= (
180 (math
.radians(sun
.elevation
- 90), 0,
181 math
.radians(-sun
.az_north
)))
184 day_increment
= 365 / object_count
185 day
= sun_props
.day_of_year
+ day_increment
* (object_count
- 1)
186 for obj
in sun_objects
:
187 dt
= (datetime
.date(sun_props
.year
, 1, 1) +
188 datetime
.timedelta(day
- 1))
189 get_sun_position(local_time
, sun_props
.latitude
,
190 sun_props
.longitude
, north_offset
, zone
,
191 dt
.month
, dt
.day
, sun_props
.year
,
192 sun_props
.sun_distance
)
193 set_sun_position(obj
, sun_props
.sun_distance
)
195 obj
.rotation_euler
= (
196 (math
.radians(sun
.elevation
- 90), 0,
197 math
.radians(-sun
.az_north
)))
199 def update_time(context
):
200 sun_props
= context
.scene
.sun_pos_properties
202 if sun_props
.use_day_of_year
:
203 dt
= (datetime
.date(sun_props
.year
, 1, 1) +
204 datetime
.timedelta(sun_props
.day_of_year
- 1))
207 sun
.day_of_year
= sun_props
.day_of_year
208 if sun_props
.day
!= dt
.day
:
209 sun_props
.day
= dt
.day
210 if sun_props
.month
!= dt
.month
:
211 sun_props
.month
= dt
.month
213 dt
= datetime
.date(sun_props
.year
, sun_props
.month
, sun_props
.day
)
214 day_of_year
= dt
.timetuple().tm_yday
215 if sun_props
.day_of_year
!= day_of_year
:
216 sun_props
.day_of_year
= day_of_year
217 sun
.day
= sun_props
.day
218 sun
.month
= sun_props
.month
219 sun
.day_of_year
= day_of_year
220 sun
.year
= sun_props
.year
221 sun
.longitude
= sun_props
.longitude
222 sun
.latitude
= sun_props
.latitude
223 sun
.UTC_zone
= sun_props
.UTC_zone
226 def format_time(the_time
, daylight_savings
, longitude
, UTC_zone
=None):
227 if UTC_zone
is not None:
235 mm
= (the_time
- int(the_time
)) * 60
236 ss
= int((mm
- int(mm
)) * 60)
238 return ("%02i:%02i:%02i" % (hh
, mm
, ss
))
241 def format_hms(the_time
):
242 hh
= str(int(the_time
))
243 min = (the_time
- int(the_time
)) * 60
244 sec
= int((min - int(min)) * 60)
245 mm
= "0" + str(int(min)) if min < 10 else str(int(min))
246 ss
= "0" + str(sec
) if sec
< 10 else str(sec
)
248 return (hh
+ ":" + mm
+ ":" + ss
)
251 def format_lat_long(lat_long
, is_latitude
):
252 hh
= str(abs(int(lat_long
)))
253 min = abs((lat_long
- int(lat_long
)) * 60)
254 sec
= abs(int((min - int(min)) * 60))
255 mm
= "0" + str(int(min)) if min < 10 else str(int(min))
256 ss
= "0" + str(sec
) if sec
< 10 else str(sec
)
261 coord_tag
= " N" if lat_long
> 0 else " S"
263 coord_tag
= " E" if lat_long
> 0 else " W"
265 return hh
+ "° " + mm
+ "' " + ss
+ '"' + coord_tag
270 ############################################################################
272 # Calculate the actual position of the sun based on input parameters.
274 # The sun positioning algorithms below are based on the National Oceanic
275 # and Atmospheric Administration's (NOAA) Solar Position Calculator
276 # which rely on calculations of Jean Meeus' book "Astronomical Algorithms."
277 # Use of NOAA data and products are in the public domain and may be used
278 # freely by the public as outlined in their policies at
279 # www.nws.noaa.gov/disclaimer.php
281 # The calculations of this script can be verified with those of NOAA's
282 # using the Azimuth and Solar Elevation displayed in the SunPos_Panel.
283 # NOAA's web site is:
284 # http://www.esrl.noaa.gov/gmd/grad/solcalc
285 ############################################################################
288 def get_sun_position(local_time
, latitude
, longitude
, north_offset
,
289 utc_zone
, month
, day
, year
, distance
):
291 addon_prefs
= bpy
.context
.preferences
.addons
[__package__
].preferences
292 sun_props
= bpy
.context
.scene
.sun_pos_properties
294 longitude
*= -1 # for internal calculations
295 utc_time
= local_time
+ utc_zone
# Set Greenwich Meridian Time
297 if latitude
> 89.93: # Latitude 90 and -90 gives
298 latitude
= radians(89.93) # erroneous results so nudge it
299 elif latitude
< -89.93:
300 latitude
= radians(-89.93)
302 latitude
= radians(latitude
)
304 t
= julian_time_from_y2k(utc_time
, year
, month
, day
)
306 e
= radians(obliquity_correction(t
))
307 L
= apparent_longitude_of_sun(t
)
308 solar_dec
= sun_declination(e
, L
)
309 eqtime
= calc_equation_of_time(t
)
311 time_correction
= (eqtime
- 4 * longitude
) + 60 * utc_zone
312 true_solar_time
= ((utc_time
- utc_zone
) * 60.0 + time_correction
) % 1440
314 hour_angle
= true_solar_time
/ 4.0 - 180.0
315 if hour_angle
< -180.0:
318 csz
= (math
.sin(latitude
) * math
.sin(solar_dec
) +
319 math
.cos(latitude
) * math
.cos(solar_dec
) *
320 math
.cos(radians(hour_angle
)))
326 zenith
= math
.acos(csz
)
328 az_denom
= math
.cos(latitude
) * math
.sin(zenith
)
330 if abs(az_denom
) > 0.001:
331 az_rad
= ((math
.sin(latitude
) *
332 math
.cos(zenith
)) - math
.sin(solar_dec
)) / az_denom
333 if abs(az_rad
) > 1.0:
334 az_rad
= -1.0 if (az_rad
< 0.0) else 1.0
335 azimuth
= 180.0 - degrees(math
.acos(az_rad
))
339 azimuth
= 180.0 if (latitude
> 0.0) else 0.0
342 azimuth
= azimuth
+ 360.0
344 exoatm_elevation
= 90.0 - degrees(zenith
)
346 if sun_props
.use_refraction
:
347 if exoatm_elevation
> 85.0:
348 refraction_correction
= 0.0
350 te
= math
.tan(radians(exoatm_elevation
))
351 if exoatm_elevation
> 5.0:
352 refraction_correction
= (
353 58.1 / te
- 0.07 / (te
** 3) + 0.000086 / (te
** 5))
354 elif (exoatm_elevation
> -0.575):
355 s1
= (-12.79 + exoatm_elevation
* 0.711)
356 s2
= (103.4 + exoatm_elevation
* (s1
))
357 s3
= (-518.2 + exoatm_elevation
* (s2
))
358 refraction_correction
= 1735.0 + exoatm_elevation
* (s3
)
360 refraction_correction
= -20.774 / te
362 refraction_correction
= refraction_correction
/ 3600
363 solar_elevation
= 90.0 - (degrees(zenith
) - refraction_correction
)
366 solar_elevation
= 90.0 - degrees(zenith
)
368 solar_azimuth
= azimuth
369 solar_azimuth
+= north_offset
371 sun
.az_north
= solar_azimuth
373 sun
.theta
= math
.pi
/ 2 - radians(solar_elevation
)
374 sun
.phi
= radians(solar_azimuth
) * -1
375 sun
.azimuth
= azimuth
376 sun
.elevation
= solar_elevation
379 def set_sun_position(obj
, distance
):
380 locX
= math
.sin(sun
.phi
) * math
.sin(-sun
.theta
) * distance
381 locY
= math
.sin(sun
.theta
) * math
.cos(sun
.phi
) * distance
382 locZ
= math
.cos(sun
.theta
) * distance
384 #----------------------------------------------
385 # Update selected object in viewport
386 #----------------------------------------------
387 obj
.location
= locX
, locY
, locZ
390 def set_sun_rotations(obj
, rotation_euler
):
391 rotation_quaternion
= rotation_euler
.to_quaternion()
392 obj
.rotation_quaternion
= rotation_quaternion
394 if obj
.rotation_mode
in {'XZY', 'YXZ', 'YZX', 'ZXY','ZYX'}:
395 obj
.rotation_euler
= rotation_quaternion
.to_euler(obj
.rotation_mode
)
397 obj
.rotation_euler
= rotation_euler
399 rotation_axis_angle
= obj
.rotation_quaternion
.to_axis_angle()
400 obj
.rotation_axis_angle
= (rotation_axis_angle
[1],
401 *rotation_axis_angle
[0])
404 def calc_sunrise_set_UTC(rise
, jd
, latitude
, longitude
):
405 t
= calc_time_julian_cent(jd
)
406 eq_time
= calc_equation_of_time(t
)
407 solar_dec
= calc_sun_declination(t
)
408 hour_angle
= calc_hour_angle_sunrise(latitude
, solar_dec
)
410 hour_angle
= -hour_angle
411 delta
= longitude
+ degrees(hour_angle
)
412 time_UTC
= 720 - (4.0 * delta
) - eq_time
416 def calc_sun_declination(t
):
417 e
= radians(obliquity_correction(t
))
418 L
= apparent_longitude_of_sun(t
)
419 solar_dec
= sun_declination(e
, L
)
423 def calc_hour_angle_sunrise(lat
, solar_dec
):
424 lat_rad
= radians(lat
)
425 HAarg
= (math
.cos(radians(90.833)) /
426 (math
.cos(lat_rad
) * math
.cos(solar_dec
))
427 - math
.tan(lat_rad
) * math
.tan(solar_dec
))
432 HA
= math
.acos(HAarg
)
436 def calc_solar_noon(jd
, longitude
, timezone
, dst
):
437 t
= calc_time_julian_cent(jd
- longitude
/ 360.0)
438 eq_time
= calc_equation_of_time(t
)
439 noon_offset
= 720.0 - (longitude
* 4.0) - eq_time
440 newt
= calc_time_julian_cent(jd
+ noon_offset
/ 1440.0)
441 eq_time
= calc_equation_of_time(newt
)
443 nv
= 780.0 if dst
else 720.0
444 noon_local
= (nv
- (longitude
* 4.0) - eq_time
+ (timezone
* 60.0)) % 1440
445 sun
.solar_noon
.time
= noon_local
/ 60.0
448 def calc_sunrise_sunset(rise
):
451 jd
= get_julian_day(sun
.year
, sun
.month
, sun
.day
)
452 time_UTC
= calc_sunrise_set_UTC(rise
, jd
, sun
.latitude
, sun
.longitude
)
453 new_time_UTC
= calc_sunrise_set_UTC(rise
, jd
+ time_UTC
/ 1440.0,
454 sun
.latitude
, sun
.longitude
)
455 time_local
= new_time_UTC
+ (-zone
* 60.0)
456 tl
= time_local
/ 60.0
457 get_sun_position(tl
, sun
.latitude
, sun
.longitude
, 0.0,
458 zone
, sun
.month
, sun
.day
, sun
.year
,
460 if sun
.use_daylight_savings
:
462 tl
= time_local
/ 60.0
465 sun
.sunrise
.time
= tl
466 sun
.sunrise
.azimuth
= sun
.azimuth
467 sun
.sunrise
.elevation
= sun
.elevation
468 calc_solar_noon(jd
, sun
.longitude
, -zone
, sun
.use_daylight_savings
)
469 get_sun_position(sun
.solar_noon
.time
, sun
.latitude
, sun
.longitude
,
470 0.0, zone
, sun
.month
, sun
.day
, sun
.year
,
472 sun
.solar_noon
.elevation
= sun
.elevation
475 sun
.sunset
.azimuth
= sun
.azimuth
476 sun
.sunset
.elevation
= sun
.elevation
478 ##########################################################################
479 ## Get the elapsed julian time since 1/1/2000 12:00 gmt
480 ## Y2k epoch (1/1/2000 12:00 gmt) is Julian day 2451545.0
481 ##########################################################################
484 def julian_time_from_y2k(utc_time
, year
, month
, day
):
485 century
= 36525.0 # Days in Julian Century
486 epoch
= 2451545.0 # Julian Day for 1/1/2000 12:00 gmt
487 jd
= get_julian_day(year
, month
, day
)
488 return ((jd
+ (utc_time
/ 24)) - epoch
) / century
491 def get_julian_day(year
, month
, day
):
495 A
= math
.floor(year
/ 100)
496 B
= 2 - A
+ math
.floor(A
/ 4.0)
497 jd
= (math
.floor((365.25 * (year
+ 4716.0))) +
498 math
.floor(30.6001 * (month
+ 1)) + day
+ B
- 1524.5)
502 def calc_time_julian_cent(jd
):
503 t
= (jd
- 2451545.0) / 36525.0
507 def sun_declination(e
, L
):
508 return (math
.asin(math
.sin(e
) * math
.sin(L
)))
511 def calc_equation_of_time(t
):
512 epsilon
= obliquity_correction(t
)
513 ml
= radians(mean_longitude_sun(t
))
514 e
= eccentricity_earth_orbit(t
)
515 m
= radians(mean_anomaly_sun(t
))
516 y
= math
.tan(radians(epsilon
) / 2.0)
518 sin2ml
= math
.sin(2.0 * ml
)
519 cos2ml
= math
.cos(2.0 * ml
)
520 sin4ml
= math
.sin(4.0 * ml
)
522 sin2m
= math
.sin(2.0 * m
)
523 etime
= (y
* sin2ml
- 2.0 * e
* sinm
+ 4.0 * e
* y
*
524 sinm
* cos2ml
- 0.5 * y
** 2 * sin4ml
- 1.25 * e
** 2 * sin2m
)
525 return (degrees(etime
) * 4)
528 def obliquity_correction(t
):
529 ec
= obliquity_of_ecliptic(t
)
530 omega
= 125.04 - 1934.136 * t
531 return (ec
+ 0.00256 * math
.cos(radians(omega
)))
534 def obliquity_of_ecliptic(t
):
535 return ((23.0 + 26.0 / 60 + (21.4480 - 46.8150) / 3600 * t
-
536 (0.00059 / 3600) * t
** 2 + (0.001813 / 3600) * t
** 3))
539 def true_longitude_of_sun(t
):
540 return (mean_longitude_sun(t
) + equation_of_sun_center(t
))
543 def calc_sun_apparent_long(t
):
544 o
= true_longitude_of_sun(t
)
545 omega
= 125.04 - 1934.136 * t
546 lamb
= o
- 0.00569 - 0.00478 * math
.sin(radians(omega
))
550 def apparent_longitude_of_sun(t
):
551 return (radians(true_longitude_of_sun(t
) - 0.00569 - 0.00478 *
552 math
.sin(radians(125.04 - 1934.136 * t
))))
555 def mean_longitude_sun(t
):
556 return (280.46646 + 36000.76983 * t
+ 0.0003032 * t
** 2) % 360
559 def equation_of_sun_center(t
):
560 m
= radians(mean_anomaly_sun(t
))
561 c
= ((1.914602 - 0.004817 * t
- 0.000014 * t
** 2) * math
.sin(m
) +
562 (0.019993 - 0.000101 * t
) * math
.sin(m
* 2) +
563 0.000289 * math
.sin(m
* 3))
567 def mean_anomaly_sun(t
):
568 return (357.52911 + t
* (35999.05029 - 0.0001537 * t
))
571 def eccentricity_earth_orbit(t
):
572 return (0.016708634 - 0.000042037 * t
- 0.0000001267 * t
** 2)