sun_rise_formula <- function(day_number, declination, latitude, s){
	# altitude correction needed; otherwise centre of the sun disk.
	# - 16 arc minutes for semi-diameter of the sun disk
	# horizontal refraction: - 34 arc minutes
	# total - 50 arc minutes =  - 50 / 60 deg 
	altitude_rad <- - 50 / 60 * pi / 180 + s * pi / 180
	latitude_rad <- latitude * pi / 180
	declination_rad <- declination * pi / 180
	w <- (sin(altitude_rad) - (sin(latitude_rad) * sin(declination_rad))) / 
		(cos(latitude_rad) * cos(declination_rad))
	hour_angle_rad <- acos(w)
	hour_angle <- hour_angle_rad * 180 / pi
	rise_set_seconds <- hour_angle * 240
	result <- rise_set_seconds
}