Gracefully handle midnight sun/polar night

The sun trajectory calculation previously emitted garbage if midnight
sun or polar night phenomena was in effect. To fix this, the calculation
is overhauled to report if the computation failed, and if so, which
specific phenomena was the cause.

Timing recalulation uses this information to change out of the normal
state that uses the four sun position timings, into a static state where
wlsunset simply waits a day at a time for a normal sun trajectory
pattern to be restored. In this state, a fixedc high/low temperature is
set as appropriate for the phenomena.

A transition phase replicating the previous day's sunrise is used to
smooth out entry into a midnight sun, which would otherwise cause a
jarring instant-daylight setting.
master
Kenny Levinsen 2020-10-17 13:23:48 +02:00
parent 7bd167e9d3
commit e14e8efcf3
3 changed files with 253 additions and 78 deletions

View File

@ -4,51 +4,71 @@
#include <time.h>
#include "color_math.h"
#define SOLAR_HORIZON 90.833
#define SOLAR_START_TWILIGHT 6.0
#define SOLAR_END_TWILIGHT -3.0
static int is_leap(int year) {
return (year % 4 == 0 && year % 100 != 0) || year % 400 == 0;
}
static double SOLAR_START_TWILIGHT = RADIANS(90.833 + 6.0);
static double SOLAR_END_TWILIGHT = RADIANS(90.833 - 3.0);
static int days_in_year(int year) {
return is_leap(year) ? 366 : 365;
int leap = (year % 4 == 0 && year % 100 != 0) || year % 400 == 0;
return leap ? 366 : 365;
}
static double radians(double degrees) {
return degrees * M_PI / 180.0;
static double date_orbit_angle(struct tm *tm) {
return 2 * M_PI / (double)days_in_year(tm->tm_year + 1900) * tm->tm_yday;
}
static double degrees(double radians) {
return radians * 180.0 / M_PI;
}
void calc_sun(struct tm *tm, double longitude, double latitude, struct sun *sun) {
static double equation_of_time(double orbit_angle) {
// https://www.esrl.noaa.gov/gmd/grad/solcalc/solareqns.PDF
double year_rad = 2 * M_PI / days_in_year(tm->tm_year) * (tm->tm_yday - 1);
double eqtime = 229.18 * (0.000075 +
0.001868 * cos(year_rad) -
0.032077 * sin(year_rad) -
0.014615 * cos(2*year_rad) -
0.040849 * sin(2*year_rad));
double decl = 0.006918 -
0.399912 * cos(year_rad) +
0.070257 * sin(year_rad) -
0.006758 * cos(2*year_rad) +
0.000907 * sin(2*year_rad) -
0.002697 * cos(3*year_rad) +
0.00148 * sin(3*year_rad);
double ha1 = degrees(acos(
cos(radians(SOLAR_HORIZON + SOLAR_START_TWILIGHT)) / (cos(radians(latitude)) * cos(decl)) -
tan(radians(latitude)) * tan(decl)));
double ha2 = degrees(acos(
cos(radians(SOLAR_HORIZON + SOLAR_END_TWILIGHT)) / (cos(radians(latitude)) * cos(decl)) -
tan(radians(latitude)) * tan(decl)));
sun->dawn = (720 - 4 * (longitude + fabs(ha1)) - eqtime) * 60;
sun->sunrise = (720 - 4 * (longitude + fabs(ha2)) - eqtime) * 60;
sun->sunset = (720 - 4 * (longitude - fabs(ha2)) - eqtime) * 60;
sun->dusk = (720 - 4 * (longitude - fabs(ha1)) - eqtime) * 60;
return 4 * (0.000075 +
0.001868 * cos(orbit_angle) -
0.032077 * sin(orbit_angle) -
0.014615 * cos(2*orbit_angle) -
0.040849 * sin(2*orbit_angle));
}
static double sun_declination(double orbit_angle) {
// https://www.esrl.noaa.gov/gmd/grad/solcalc/solareqns.PDF
return 0.006918 -
0.399912 * cos(orbit_angle) +
0.070257 * sin(orbit_angle) -
0.006758 * cos(2*orbit_angle) +
0.000907 * sin(2*orbit_angle) -
0.002697 * cos(3*orbit_angle) +
0.00148 * sin(3*orbit_angle);
}
static double sun_hour_angle(double latitude, double declination, double target_sun) {
// https://www.esrl.noaa.gov/gmd/grad/solcalc/solareqns.PDF
return acos(cos(target_sun) /
cos(latitude) * cos(declination) -
tan(latitude) * tan(declination));
}
static time_t hour_angle_to_time(double longitude, double eqtime, double hour_angle) {
// https://www.esrl.noaa.gov/gmd/grad/solcalc/solareqns.PDF
return isnan(hour_angle) ? -1 :
DEGREES((4.0 * M_PI - 4 * (longitude + hour_angle) - eqtime) * 60);
}
static enum sun_condition condition(double latitude_rad, double sun_declination) {
int sign_lat = signbit(latitude_rad) == 0;
int sign_decl = signbit(sun_declination) == 0;
return sign_lat == sign_decl ? MIDNIGHT_SUN : POLAR_NIGHT;
}
enum sun_condition calc_sun(struct tm *tm, double longitude, double latitude, struct sun *sun) {
double orbit_angle = date_orbit_angle(tm);
double decl = sun_declination(orbit_angle);
double eqtime = equation_of_time(orbit_angle);
double ha_twilight = sun_hour_angle(latitude, decl, SOLAR_START_TWILIGHT);
double ha_daylight = sun_hour_angle(latitude, decl, SOLAR_END_TWILIGHT);
sun->dawn = hour_angle_to_time(longitude, eqtime, fabs(ha_twilight));
sun->dusk = hour_angle_to_time(longitude, eqtime, -fabs(ha_twilight));
sun->sunrise = hour_angle_to_time(longitude, eqtime, fabs(ha_daylight));
sun->sunset = hour_angle_to_time(longitude, eqtime, -fabs(ha_daylight));
return isnan(ha_twilight) || isnan(ha_daylight) ? condition(latitude, decl) : NORMAL;
}
static int illuminant_d(int temp, double *x, double *y) {

View File

@ -1,8 +1,20 @@
#ifndef _COLOR_MATH_H
#define _COLOR_MATH_H
#include "math.h"
#include "time.h"
// These are macros so they can be applied to constants
#define DEGREES(rad) ((rad) * 180.0 / M_PI)
#define RADIANS(deg) ((deg) * M_PI / 180.0)
enum sun_condition {
NORMAL,
MIDNIGHT_SUN,
POLAR_NIGHT,
SUN_CONDITION_LAST
};
struct sun {
time_t dawn;
time_t sunrise;
@ -10,8 +22,8 @@ struct sun {
time_t dusk;
};
void calc_sun(struct tm *tm, double longitude, double latitude, struct sun *sun);
double clamp(double value);
enum sun_condition calc_sun(struct tm *tm, double longitude, double latitude, struct sun *sun);
void calc_whitepoint(int temp, double *rw, double *gw, double *bw);
double clamp(double value);
#endif

219
main.c
View File

@ -71,12 +71,23 @@ struct config {
double latitude;
};
enum state {
STATE_INITIAL,
STATE_NORMAL,
STATE_TRANSITION,
STATE_STATIC,
};
struct context {
struct config config;
struct sun sun;
enum state state;
enum sun_condition condition;
time_t dawn_step_time;
time_t dusk_step_time;
time_t calc_day;
bool new_output;
struct wl_list outputs;
@ -96,6 +107,14 @@ struct output {
uint16_t *table;
};
static time_t round_day(time_t now) {
return now - (now % 86400);
}
static time_t tomorrow(time_t now) {
return round_day(now) + 86400;
}
static struct zwlr_gamma_control_manager_v1 *gamma_control_manager = NULL;
static int create_anonymous_file(off_t size) {
@ -267,56 +286,123 @@ static void set_temperature(struct wl_list *outputs, int temp, int gamma) {
}
}
static const char *sun_condition_str[] = {
"normal",
"midnight sun",
"polar night",
"invalid",
NULL,
};
static void print_trajectory(struct context *ctx, time_t day) {
fprintf(stderr, "calculated sun trajectory:");
struct tm tm;
if (ctx->sun.dawn >= day) {
localtime_r(&ctx->sun.dawn, &tm);
fprintf(stderr, " dawn %02d:%02d,", tm.tm_hour, tm.tm_min);
} else {
fprintf(stderr, " dawn N/A,");
}
if (ctx->sun.sunrise >= day) {
localtime_r(&ctx->sun.sunrise, &tm);
fprintf(stderr, " sunrise %02d:%02d,", tm.tm_hour, tm.tm_min);
} else {
fprintf(stderr, " sunrise N/A,");
}
if (ctx->sun.sunset >= day) {
localtime_r(&ctx->sun.sunset, &tm);
fprintf(stderr, " sunset %02d:%02d,", tm.tm_hour, tm.tm_min);
} else {
fprintf(stderr, " sunset N/A,");
}
if (ctx->sun.dusk >= day) {
localtime_r(&ctx->sun.dusk, &tm);
fprintf(stderr, " dusk %02d:%02d,", tm.tm_hour, tm.tm_min);
} else {
fprintf(stderr, " dusk N/A,");
}
fprintf(stderr, " condition: %s\n", sun_condition_str[ctx->condition]);
}
static int anim_kelvin_step = 25;
static void recalc_stops(struct context *ctx, time_t now) {
if (now < ctx->sun.dusk) {
time_t day = round_day(now);
time_t last_day = ctx->calc_day;
if (day == last_day) {
return;
}
ctx->calc_day = day;
time_t day = now - (now % 86400);
if (day < ctx->sun.dusk) {
day += 86400;
}
struct sun sun;
struct tm tm = { 0 };
gmtime_r(&day, &tm);
calc_sun(&tm, ctx->config.longitude, ctx->config.latitude, &ctx->sun);
if (ctx->config.duration != -1) {
ctx->sun.dawn = ctx->sun.sunrise - ctx->config.duration;
ctx->sun.dusk = ctx->sun.sunset + ctx->config.duration;
}
enum sun_condition cond = calc_sun(&tm, ctx->config.longitude, ctx->config.latitude, &sun);
// TODO: Cap on eternal days?
ctx->sun.dawn += day;
ctx->sun.sunrise += day;
ctx->sun.sunset += day;
ctx->sun.dusk += day;
assert(ctx->sun.dusk > now);
switch (cond) {
case NORMAL:
if (ctx->condition == MIDNIGHT_SUN) {
// Yesterday had no sunset, so remove our sunrise.
sun.dawn = -1;
sun.sunrise = -1;
}
ctx->state = STATE_NORMAL;
break;
case MIDNIGHT_SUN:
if (ctx->condition == POLAR_NIGHT) {
fprintf(stderr, "warning: direct polar night to midnight sun transition\n");
}
if (ctx->state != STATE_NORMAL) {
ctx->state = STATE_STATIC;
break;
}
// Borrow yesterday's sunrise to animate into the midnight sun
sun.dawn = sun.dawn != -1 ? sun.dawn :
ctx->sun.dawn - last_day;
sun.sunrise = sun.sunrise != -1 ? sun.sunrise :
ctx->sun.sunrise - last_day;
ctx->state = STATE_TRANSITION;
break;
case POLAR_NIGHT:
if (ctx->condition == MIDNIGHT_SUN) {
fprintf(stderr, "warning: direct midnight sun to polar night transition\n");
}
ctx->state = STATE_STATIC;
break;
default:
abort();
}
ctx->condition = cond;
ctx->sun.dawn = sun.dawn + day;
ctx->sun.sunrise = sun.sunrise + day;
ctx->sun.sunset = sun.sunset + day;
ctx->sun.dusk = sun.dusk + day;
int temp_diff = ctx->config.high_temp - ctx->config.low_temp;
ctx->dawn_step_time = (ctx->sun.sunrise - ctx->sun.dawn) * anim_kelvin_step / temp_diff;
ctx->dusk_step_time = (ctx->sun.dusk - ctx->sun.sunset) * anim_kelvin_step / temp_diff;
struct tm dawn, sunrise, sunset, dusk;
localtime_r(&ctx->sun.dawn, &dawn);
localtime_r(&ctx->sun.sunrise, &sunrise);
localtime_r(&ctx->sun.sunset, &sunset);
localtime_r(&ctx->sun.dusk, &dusk);
fprintf(stderr, "calculated new sun trajectory: dawn %02d:%02d, sunrise %02d:%02d, sunset %02d:%02d, dusk %02d:%02d\n",
dawn.tm_hour, dawn.tm_min,
sunrise.tm_hour, sunrise.tm_min,
sunset.tm_hour, sunset.tm_min,
dusk.tm_hour, dusk.tm_min);
print_trajectory(ctx, day);
}
static int interpolate_temperature(time_t now, time_t start, time_t stop, int temp_start, int temp_stop) {
if (start == stop) {
return stop;
}
double time_pos = clamp((double)(now - start) / (double)(stop - start));
int temp_pos = (double)(temp_stop - temp_start) * time_pos;
return temp_start + temp_pos;
}
static int get_temperature(const struct context *ctx, time_t now) {
static int get_temperature_normal(const struct context *ctx, time_t now) {
if (now < ctx->sun.dawn) {
return ctx->config.low_temp;
} else if (now < ctx->sun.sunrise) {
@ -330,18 +416,73 @@ static int get_temperature(const struct context *ctx, time_t now) {
}
}
static void update_timer(struct context *ctx, timer_t timer, time_t now) {
assert(now < ctx->sun.dusk);
static int get_temperature_transition(const struct context *ctx, time_t now) {
switch (ctx->condition) {
case MIDNIGHT_SUN:
if (now < ctx->sun.sunrise) {
return get_temperature_normal(ctx, now);
}
return ctx->config.high_temp;
default:
abort();
}
}
time_t deadline;
static int get_temperature(const struct context *ctx, time_t now) {
switch (ctx->state) {
case STATE_NORMAL:
return get_temperature_normal(ctx, now);
case STATE_TRANSITION:
return get_temperature_transition(ctx, now);
case STATE_STATIC:
return ctx->condition == MIDNIGHT_SUN ? ctx->config.high_temp : ctx->config.low_temp;
default:
abort();
}
}
static time_t get_deadline_normal(const struct context *ctx, time_t now) {
if (now < ctx->sun.dawn) {
deadline = ctx->sun.dawn;
return ctx->sun.dawn;
} else if (now < ctx->sun.sunrise) {
deadline = now + ctx->dawn_step_time;
return now + ctx->dawn_step_time;
} else if (now < ctx->sun.sunset) {
deadline = ctx->sun.sunset;
return ctx->sun.sunset;
} else if (now < ctx->sun.dusk) {
deadline = now + ctx->dusk_step_time;
return now + ctx->dusk_step_time;
} else {
return tomorrow(now);
}
}
static time_t get_deadline_transition(const struct context *ctx, time_t now) {
switch (ctx->condition) {
case MIDNIGHT_SUN:
if (now < ctx->sun.sunrise) {
return get_deadline_normal(ctx, now);
}
// fallthrough
case POLAR_NIGHT:
return tomorrow(now);
default:
abort();
}
}
static void update_timer(const struct context *ctx, timer_t timer, time_t now) {
time_t deadline;
switch (ctx->state) {
case STATE_NORMAL:
deadline = get_deadline_normal(ctx, now);
break;
case STATE_TRANSITION:
deadline = get_deadline_transition(ctx, now);
break;
case STATE_STATIC:
deadline = tomorrow(now);
break;
default:
abort();
}
assert(deadline > now);
@ -416,6 +557,8 @@ int main(int argc, char *argv[]) {
// Initialize defaults
struct context ctx = {
.sun = { 0 },
.condition = SUN_CONDITION_LAST,
.state = STATE_INITIAL,
.config = {
.high_temp = 6500,
.low_temp = 4000,
@ -455,10 +598,10 @@ int main(int argc, char *argv[]) {
ctx.config.low_temp = strtol(optarg, NULL, 10);
break;
case 'l':
ctx.config.latitude = strtod(optarg, NULL);
ctx.config.latitude = RADIANS(strtod(optarg, NULL));
break;
case 'L':
ctx.config.longitude = strtod(optarg, NULL);
ctx.config.longitude = RADIANS(strtod(optarg, NULL));
break;
case 'd':
fprintf(stderr, "using animation duration override\n");
@ -505,10 +648,10 @@ int main(int argc, char *argv[]) {
time_t now = get_time_sec(NULL);
recalc_stops(&ctx, now);
update_timer(&ctx, ctx.timer, now);
int temp = get_temperature(&ctx, now);
set_temperature(&ctx.outputs, temp, gamma);
update_timer(&ctx, ctx.timer, now);
int old_temp = temp;
while (display_dispatch(display, -1) != -1) {