Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix eternal season affect on weather #54100

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
36 changes: 14 additions & 22 deletions src/calendar.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -92,26 +92,18 @@ static constexpr time_duration angle_to_time( const units::angle a )
return a / 15.0_degrees * 1_hours;
}

// To support the eternal season option we create a strong typedef of timepoint
// which is a solar_effective_time. This converts a regular time to a time
// which would be relevant for sun position calculations. Normally the two
// times are the same, but when eternal seasons are used the effective time is
// always set to the same day, so that the sun position doesn't change from day
// to day.
struct solar_effective_time {
explicit solar_effective_time( const time_point &t_ )
: t( t_ ) {
if( calendar::eternal_season() ) {
const time_point start_midnight =
calendar::start_of_game - time_past_midnight( calendar::start_of_game );
t = start_midnight + time_past_midnight( t_ );
}
season_effective_time::season_effective_time( const time_point &t_ )
: t( t_ )
{
if( calendar::eternal_season() ) {
const time_point start_midnight =
calendar::start_of_game - time_past_midnight( calendar::start_of_game );
t = start_midnight + time_past_midnight( t_ );
}
time_point t;
};
}

static std::pair<units::angle, units::angle> sun_ra_declination(
solar_effective_time t, time_duration timezone )
season_effective_time t, time_duration timezone )
{
// This derivation is mostly from
// https://en.wikipedia.org/wiki/Position_of_the_Sun
Expand Down Expand Up @@ -155,7 +147,7 @@ static std::pair<units::angle, units::angle> sun_ra_declination(
return { right_ascension, declination };
}

static units::angle sidereal_time_at( solar_effective_time t, units::angle longitude,
static units::angle sidereal_time_at( season_effective_time t, units::angle longitude,
time_duration timezone )
{
// Repeat some calculations from sun_ra_declination
Expand All @@ -176,7 +168,7 @@ static units::angle sidereal_time_at( solar_effective_time t, units::angle longi
std::pair<units::angle, units::angle> sun_azimuth_altitude(
time_point ti )
{
const solar_effective_time t = solar_effective_time( ti );
const season_effective_time t = season_effective_time( ti );
const lat_long location = location_boston;
units::angle right_ascension;
units::angle declination;
Expand Down Expand Up @@ -257,7 +249,7 @@ static time_point solar_noon_near( const time_point &t )

static units::angle offset_to_sun_altitude(
const units::angle altitude, const units::angle longitude,
const solar_effective_time approx_time, const bool evening )
const season_effective_time approx_time, const bool evening )
{
units::angle ra;
units::angle declination;
Expand Down Expand Up @@ -286,7 +278,7 @@ static time_point sun_at_altitude( const units::angle altitude, const units::ang
{
const time_point solar_noon = solar_noon_near( t );
units::angle initial_offset =
offset_to_sun_altitude( altitude, longitude, solar_effective_time( solar_noon ), evening );
offset_to_sun_altitude( altitude, longitude, season_effective_time( solar_noon ), evening );
if( !evening ) {
initial_offset -= 360_degrees;
}
Expand All @@ -295,7 +287,7 @@ static time_point sun_at_altitude( const units::angle altitude, const units::ang
// Now we should have the correct time to within a few minutes; iterate to
// get a more precise estimate
units::angle correction_offset =
offset_to_sun_altitude( altitude, longitude, solar_effective_time( initial_approximation ),
offset_to_sun_altitude( altitude, longitude, season_effective_time( initial_approximation ),
evening );
if( correction_offset > 180_degrees ) {
correction_offset -= 360_degrees;
Expand Down
12 changes: 12 additions & 0 deletions src/calendar.h
Original file line number Diff line number Diff line change
Expand Up @@ -623,4 +623,16 @@ enum class weekdays : int {

weekdays day_of_week( const time_point &p );

// To support the eternal season option we create a strong typedef of timepoint
// which is a season_effective_time. This converts a regular time to a time
// which would be relevant for sun position and weather calculations. Normally
// the two times are the same, but when eternal seasons are used the effective
// time is always set to the same day, so that the sun position and weather
// doesn't change from day to day.
struct season_effective_time {
season_effective_time() = default;
explicit season_effective_time( const time_point & );
time_point t;
};

#endif // CATA_SRC_CALENDAR_H
30 changes: 17 additions & 13 deletions src/weather_gen.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -47,34 +47,35 @@ struct weather_gen_common {
season_type season = season_type::SPRING;
};

static weather_gen_common get_common_data( const tripoint &location, const time_point &t,
static weather_gen_common get_common_data( const tripoint &location, const time_point &real_t,
unsigned seed )
{
season_effective_time t( real_t );
weather_gen_common result;
// Integer x position / widening factor of the Perlin function.
result.x = location.x / 2000.0;
// Integer y position / widening factor of the Perlin function.
result.y = location.y / 2000.0;
// Integer turn / widening factor of the Perlin function.
result.z = to_days<double>( t - calendar::turn_zero );
result.z = to_days<double>( real_t - calendar::turn_zero );
// Limit the random seed during noise calculation, a large value flattens the noise generator to zero
// Windows has a rand limit of 32768, other operating systems can have higher limits
result.modSEED = seed % SIMPLEX_NOISE_RANDOM_SEED_LIMIT;
const double year_fraction( time_past_new_year( t ) /
const double year_fraction( time_past_new_year( t.t ) /
calendar::year_length() ); // [0,1)

result.cyf = std::cos( tau * ( year_fraction + .125 ) ); // [-1, 1]
// We add one-eighth to line up `cyf` so that 1 is at
// midwinter and -1 at midsummer. (Cataclsym DDA years
// start when spring starts. Gregorian years start when
// winter starts.)
result.season = season_of_year( t );
result.season = season_of_year( t.t );

return result;
}

static double weather_temperature_from_common_data( const weather_generator &wg,
const weather_gen_common &common, const time_point &t )
const weather_gen_common &common, const season_effective_time &t )
{
const double x( common.x );
const double y( common.y );
Expand All @@ -84,7 +85,7 @@ static double weather_temperature_from_common_data( const weather_generator &wg,
const double seasonality = -common.cyf;
// -1 in midwinter, +1 in midsummer
const season_type season = common.season;
const double dayFraction = time_past_midnight( t ) / 1_days;
const double dayFraction = time_past_midnight( t.t ) / 1_days;
const double dayv = std::cos( tau * ( dayFraction + .5 - coldest_hour / 24 ) );
// -1 at coldest_hour, +1 twelve hours later

Expand All @@ -102,15 +103,17 @@ static double weather_temperature_from_common_data( const weather_generator &wg,
return T * 9 / 5 + 32;
}

double weather_generator::get_weather_temperature( const tripoint &location, const time_point &t,
unsigned seed ) const
double weather_generator::get_weather_temperature(
const tripoint &location, const time_point &real_t, unsigned seed ) const
{
return weather_temperature_from_common_data( *this, get_common_data( location, t, seed ), t );
return weather_temperature_from_common_data( *this, get_common_data( location, real_t, seed ),
season_effective_time( real_t ) );
}
w_point weather_generator::get_weather( const tripoint &location, const time_point &t,
w_point weather_generator::get_weather( const tripoint &location, const time_point &real_t,
unsigned seed ) const
{
const weather_gen_common common = get_common_data( location, t, seed );
season_effective_time t( real_t );
const weather_gen_common common = get_common_data( location, real_t, seed );

const double x( common.x );
const double y( common.y );
Expand Down Expand Up @@ -250,9 +253,10 @@ int weather_generator::get_water_temperature() const
source : http://www.grandriver.ca/index/document.cfm?Sec=2&Sub1=7&sub2=1
**/

season_effective_time t( calendar::turn );
int season_length = to_days<int>( calendar::season_length() );
int day = to_days<int>( time_past_new_year( calendar::turn ) );
int hour = hour_of_day<int>( calendar::turn );
int day = to_days<int>( time_past_new_year( t.t ) );
int hour = hour_of_day<int>( t.t );

int water_temperature = 0;

Expand Down
2 changes: 1 addition & 1 deletion src/weather_gen.h
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ struct w_point {
double windpower = 0;
std::string wind_desc;
int winddirection = 0;
time_point time;
season_effective_time time;
};

class weather_generator
Expand Down
163 changes: 111 additions & 52 deletions tests/weather_test.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -41,83 +41,142 @@ static double proportion_gteq_x( std::vector<double> const &v, double x )
return static_cast<double>( count ) / v.size();
}

TEST_CASE( "weather realism" )
// Check our simulated weather against numbers from real data
// from a few years in a few locations in New England. The numbers
// are based on NOAA's Local Climatological Data (LCD). Analysis code
// can be found at:
// https://gist.github.com/Kodiologist/e2f1e6685e8fd865650f97bb6a67ad07
static constexpr int n_hours = to_hours<int>( 1_days );
static constexpr int n_minutes = to_minutes<int>( 1_days );

struct year_of_weather_data {
explicit year_of_weather_data( int n_days )
: temperature( n_days, std::vector<double>( n_minutes, 0 ) )
, hourly_precip( n_days * n_hours, 0 )
, highs( n_days )
, lows( n_days )
{}

std::vector<std::vector<double>> temperature;
std::vector<double> hourly_precip;
std::vector<double> highs;
std::vector<double> lows;
};

static year_of_weather_data collect_weather_data( unsigned seed )
{
// Try a few randomly selected seeds.
const std::vector<unsigned> seeds = {317'024'741, 870'078'684, 1'192'447'748};

scoped_weather_override null_weather( WEATHER_NULL );
const weather_generator &wgen = get_weather().get_cur_weather_gen();

const time_point begin = calendar::turn_zero;
const time_point end = begin + calendar::year_length();
const int n_days = to_days<int>( end - begin );
const int n_hours = to_hours<int>( 1_days );
const int n_minutes = to_minutes<int>( 1_days );
year_of_weather_data result( n_days );

// Collect generated weather data for a single year.
for( time_point i = begin ; i < end ; i += 1_minutes ) {
w_point w = wgen.get_weather( tripoint_zero, i, seed );
int day = to_days<int>( time_past_new_year( i ) );
int minute = to_minutes<int>( time_past_midnight( i ) );
result.temperature[day][minute] = w.temperature;
int hour = to_hours<int>( time_past_new_year( i ) );
*get_weather().weather_precise = w;
result.hourly_precip[hour] +=
precip_mm_per_hour(
wgen.get_weather_conditions( w )->precip )
/ 60;
}

// Collect daily highs and lows.
for( int do_highs = 0 ; do_highs < 2 ; ++do_highs ) {
std::vector<double> &t = do_highs ? result.highs : result.lows;
std::transform( result.temperature.begin(), result.temperature.end(), t.begin(),
[&]( std::vector<double> const & day ) {
return do_highs
? *std::max_element( day.begin(), day.end() )
: *std::min_element( day.begin(), day.end() );
} );
}

return result;
}

// Try a few randomly selected seeds.
static const std::array<unsigned, 3> seeds = { {317'024'741, 870'078'684, 1'192'447'748} };

TEST_CASE( "weather realism", "[weather]" )
// Check our simulated weather against numbers from real data
// from a few years in a few locations in New England. The numbers
// are based on NOAA's Local Climatological Data (LCD). Analysis code
// can be found at:
// https://gist.github.com/Kodiologist/e2f1e6685e8fd865650f97bb6a67ad07
{
for( unsigned int seed : seeds ) {
std::vector<std::vector<double>> temperature;
temperature.resize( n_days, std::vector<double>( n_minutes, 0 ) );
std::vector<double> hourly_precip;
hourly_precip.resize( n_days * n_hours, 0 );

// Collect generated weather data for a single year.
for( time_point i = begin ; i < end ; i += 1_minutes ) {
w_point w = wgen.get_weather( tripoint_zero, i, seed );
int day = to_days<int>( time_past_new_year( i ) );
int minute = to_minutes<int>( time_past_midnight( i ) );
temperature[day][minute] = w.temperature;
int hour = to_hours<int>( time_past_new_year( i ) );
*get_weather().weather_precise = w;
hourly_precip[hour] +=
precip_mm_per_hour(
wgen.get_weather_conditions( w )->precip )
/ 60;
}
year_of_weather_data data = collect_weather_data( seed );

// Collect daily highs and lows.
std::vector<double> highs( n_days );
std::vector<double> lows( n_days );
for( int do_highs = 0 ; do_highs < 2 ; ++do_highs ) {
std::vector<double> &t = do_highs ? highs : lows;
std::transform( temperature.begin(), temperature.end(), t.begin(),
[&]( std::vector<double> const & day ) {
return do_highs
? *std::max_element( day.begin(), day.end() )
: *std::min_element( day.begin(), day.end() );
} );

// Check the mean absolute difference between the highs or lows
// of adjacent days (Fahrenheit).
const double d = mean_abs_running_diff( t );
CHECK( d >= ( do_highs ? 5.5 : 4 ) );
CHECK( d <= ( do_highs ? 7.5 : 7 ) );
}
// Check the mean absolute difference between the highs or lows
// of adjacent days (Fahrenheit).
const double mad_highs = mean_abs_running_diff( data.highs );
CHECK( mad_highs >= 5.5 );
CHECK( mad_highs <= 7.5 );
const double mad_lows = mean_abs_running_diff( data.lows );
CHECK( mad_lows >= 4 );
CHECK( mad_lows <= 7 );

// Check the daily mean of the range in temperatures (Fahrenheit).
const double mean_of_ranges = mean_pairwise_diffs( highs, lows );
const double mean_of_ranges = mean_pairwise_diffs( data.highs, data.lows );
CHECK( mean_of_ranges >= 14 );
CHECK( mean_of_ranges <= 25 );

// Check that summer and winter temperatures are very different.
size_t half = data.highs.size() / 4;
double summer_low = data.lows[half];
double winter_high = data.highs[0];
{
CAPTURE( summer_low );
CAPTURE( winter_high );
CHECK( summer_low - winter_high >= 10 );
}

// Check the proportion of hours with light precipitation
// or more, counting snow (mm of rain equivalent per hour).
const double at_least_light_precip = proportion_gteq_x(
hourly_precip, 1 );
const double at_least_light_precip = proportion_gteq_x( data.hourly_precip, 1 );
CHECK( at_least_light_precip >= .025 );
CHECK( at_least_light_precip <= .05 );

// Likewise for heavy precipitation.
const double heavy_precip = proportion_gteq_x(
hourly_precip, 2.5 );
const double heavy_precip = proportion_gteq_x( data.hourly_precip, 2.5 );
CHECK( heavy_precip >= .005 );
CHECK( heavy_precip <= .02 );
}
}

TEST_CASE( "eternal_season", "[weather]" )
{
on_out_of_scope restore_eternal_season( []() {
calendar::set_eternal_season( false );
} );
calendar::set_eternal_season( true );
override_option override_eternal_season( "ETERNAL_SEASON", "true" );

for( unsigned int seed : seeds ) {
year_of_weather_data data = collect_weather_data( seed );

// Check that summer and winter temperatures are very similar.
size_t half = data.highs.size() / 4;
double summer_low = data.lows[half];
double winter_high = data.highs[0];
{
CAPTURE( summer_low );
CAPTURE( winter_high );
CHECK( summer_low - winter_high <= -10 );
}

// Check the temperatures still vary from day to day.
const double mad_highs = mean_abs_running_diff( data.highs );
CHECK( mad_highs >= 5.5 );
CHECK( mad_highs <= 7.5 );
const double mad_lows = mean_abs_running_diff( data.lows );
CHECK( mad_lows >= 4 );
CHECK( mad_lows <= 7 );
}
}

TEST_CASE( "local wind chill calculation", "[weather][wind_chill]" )
{
// `get_local_windchill` returns degrees F offset from current temperature,
Expand Down