Merge branch 'sharp/jungle' of gitlab.com:veloren/veloren into sharp/jungle

This commit is contained in:
Joshua Barretto 2019-08-23 00:09:45 +01:00
commit 6136edd631
4 changed files with 261 additions and 210 deletions

View File

@ -1,10 +1,13 @@
use std::env; use std::env;
use std::fs::File; use std::fs::File;
use std::io::Read;
use std::io::Write; use std::io::Write;
use std::path::Path; use std::path::Path;
use std::path::PathBuf;
use std::process::Command; use std::process::Command;
fn main() { fn main() {
// Get the current githash
match Command::new("git") match Command::new("git")
.args(&["rev-parse", "--short", "HEAD"]) .args(&["rev-parse", "--short", "HEAD"])
.output() .output()
@ -26,4 +29,33 @@ fn main() {
}, },
Err(e) => panic!("failed to retrieve current git commit hash: {}", e), Err(e) => panic!("failed to retrieve current git commit hash: {}", e),
} }
// Check if git-lfs is working
let asset_path: PathBuf = ["..", "assets", "voxygen", "background", "bg_main.png"]
.iter()
.collect();
let asset_file = match File::open(&asset_path) {
Ok(file) => file,
Err(e) => panic!(
"failed to open asset file {}: {}",
asset_path.to_str().unwrap(),
e
),
};
const LFS_MARKER: &[u8] = b"version https://git-lfs.github.com/spec/";
let mut buffer = Vec::new();
let bytes_read = asset_file
.take(LFS_MARKER.len() as u64)
.read_to_end(&mut buffer)
.expect("failed to read asset file");
if bytes_read == LFS_MARKER.len() && buffer == LFS_MARKER {
panic!(
"\n\nGit Large File Storage (git-lfs) has not been set up correctly.\n\
Most common reasons:\n\
\t- git-lfs was not installed before cloning this repository\n\
\t- this repository was not cloned from the primary gitlab mirror.\n\
\t The github mirror does not support lfs.\n\
See the book at https://book.veloren.net/ for details.\n\n"
);
}
} }

View File

@ -574,14 +574,14 @@ impl Hud {
.color(TEXT_COLOR) .color(TEXT_COLOR)
.set(self.ids.version, ui_widgets); .set(self.ids.version, ui_widgets);
// Ticks per second // Ticks per second
Text::new(&format!("FPS: {:.1}", debug_info.tps)) Text::new(&format!("FPS: {:.0}", debug_info.tps))
.color(TEXT_COLOR) .color(TEXT_COLOR)
.down_from(self.ids.version, 5.0) .down_from(self.ids.version, 5.0)
.font_id(self.fonts.opensans) .font_id(self.fonts.opensans)
.font_size(14) .font_size(14)
.set(self.ids.fps_counter, ui_widgets); .set(self.ids.fps_counter, ui_widgets);
// Ping // Ping
Text::new(&format!("Ping: {:.1}ms", debug_info.ping_ms)) Text::new(&format!("Ping: {:.0}ms", debug_info.ping_ms))
.color(TEXT_COLOR) .color(TEXT_COLOR)
.down_from(self.ids.fps_counter, 5.0) .down_from(self.ids.fps_counter, 5.0)
.font_id(self.fonts.opensans) .font_id(self.fonts.opensans)
@ -589,7 +589,10 @@ impl Hud {
.set(self.ids.ping, ui_widgets); .set(self.ids.ping, ui_widgets);
// Player's position // Player's position
let coordinates_text = match debug_info.coordinates { let coordinates_text = match debug_info.coordinates {
Some(coordinates) => format!("Coordinates: {:.1}", coordinates.0), Some(coordinates) => format!(
"Coordinates: ({:.0}, {:.0}, {:.0})",
coordinates.0.x, coordinates.0.y, coordinates.0.z,
),
None => "Player has no Pos component".to_owned(), None => "Player has no Pos component".to_owned(),
}; };
Text::new(&coordinates_text) Text::new(&coordinates_text)
@ -601,8 +604,11 @@ impl Hud {
// Player's velocity // Player's velocity
let velocity_text = match debug_info.velocity { let velocity_text = match debug_info.velocity {
Some(velocity) => format!( Some(velocity) => format!(
"Velocity: ({:.3}, {:.3}, {:.3})", "Velocity: ({:.1}, {:.1}, {:.1}) [{:.1} u/s]",
velocity.0.x, velocity.0.y, velocity.0.z velocity.0.x,
velocity.0.y,
velocity.0.z,
velocity.0.magnitude()
), ),
None => "Player has no Vel component".to_owned(), None => "Player has no Vel component".to_owned(),
}; };

View File

@ -1,9 +1,13 @@
mod util;
mod location; mod location;
mod settlement; mod settlement;
// Reexports // Reexports
pub use self::location::Location; pub use self::location::Location;
pub use self::settlement::Settlement; pub use self::settlement::Settlement;
use self::util::{
cdf_irwin_hall, InverseCdf, uniform_idx_as_vec2, uniform_noise,
};
use crate::{ use crate::{
all::ForestKind, all::ForestKind,
@ -27,161 +31,6 @@ use vek::*;
pub const WORLD_SIZE: Vec2<usize> = Vec2 { x: 1024, y: 1024 }; pub const WORLD_SIZE: Vec2<usize> = Vec2 { x: 1024, y: 1024 };
/// Computes the cumulative distribution function of the weighted sum of k independent,
/// uniformly distributed random variables between 0 and 1. For each variable i, we use weights[i]
/// as the weight to give samples[i] (the weights should all be positive).
///
/// If the precondition is met, the distribution of the result of calling this function will be
/// uniformly distributed while preserving the same information that was in the original average.
///
/// For N > 33 the function will no longer return correct results since we will overflow u32.
///
/// NOTE:
///
/// Per [1], the problem of determing the CDF of
/// the sum of uniformly distributed random variables over *different* ranges is considerably more
/// complicated than it is for the same-range case. Fortunately, it also provides a reference to
/// [2], which contains a complete derivation of an exact rule for the density function for
/// this case. The CDF is just the integral of the cumulative distribution function [3],
/// which we use to convert this into a CDF formula.
///
/// This allows us to sum weighted, uniform, independent random variables.
///
/// At some point, we should probably contribute this back to stats-rs.
///
/// 1. https://www.r-bloggers.com/sums-of-random-variables/,
/// 2. Sadooghi-Alvandi, S., A. Nematollahi, & R. Habibi, 2009.
/// On the Distribution of the Sum of Independent Uniform Random Variables.
/// Statistical Papers, 50, 171-175.
/// 3. hhttps://en.wikipedia.org/wiki/Cumulative_distribution_function
fn cdf_irwin_hall<const N: usize>(weights: &[f32; N], samples: [f32; N]) -> f32 {
// Let J_k = {(j_1, ... , j_k) : 1 ≤ j_1 < j_2 < ··· < j_k ≤ N }.
//
// Let A_N = Π{k = 1 to n}a_k.
//
// The density function for N ≥ 2 is:
//
// 1/(A_N * (N - 1)!) * (x^(N-1) + Σ{k = 1 to N}((-1)^k *
// Σ{(j_1, ..., j_k) ∈ J_k}(max(0, x - Σ{l = 1 to k}(a_(j_l)))^(N - 1))))
//
// So the cumulative distribution function is its integral, i.e. (I think)
//
// 1/(product{k in A}(k) * N!) * (x^N + sum(k in 1 to N)((-1)^k *
// sum{j in Subsets[A, {k}]}(max(0, x - sum{l in j}(l))^N)))
//
// which is also equivalent to
//
// (letting B_k = { a in Subsets[A, {k}] : sum {l in a} l }, B_(0,1) = 0 and
// H_k = { i : 1 ≤ 1 ≤ N! / (k! * (N - k)!) })
//
// 1/(product{k in A}(k) * N!) * sum(k in 0 to N)((-1)^k *
// sum{l in H_k}(max(0, x - B_(k,l))^N))
//
// We should be able to iterate through the whole power set
// instead, and figure out K by calling count_ones(), so we can compute the result in O(2^N)
// iterations.
let x: f32 = weights
.iter()
.zip(samples.iter())
.map(|(weight, sample)| weight * sample)
.sum();
let mut y = 0.0f32;
for subset in 0u32..(1 << N) {
// Number of set elements
let k = subset.count_ones();
// Add together exactly the set elements to get B_subset
let z = weights
.iter()
.enumerate()
.filter(|(i, _)| subset & (1 << i) as u32 != 0)
.map(|(_, k)| k)
.sum::<f32>();
// Compute max(0, x - B_subset)^N
let z = (x - z).max(0.0).powi(N as i32);
// The parity of k determines whether the sum is negated.
y += if k & 1 == 0 { z } else { -z };
}
// Divide by the product of the weights.
y /= weights.iter().product::<f32>();
// Remember to multiply by 1 / N! at the end.
y / (1..=N as i32).product::<i32>() as f32
}
/// First component of each element of the vector is the computed CDF of the noise function at this
/// index (i.e. its position in a sorted list of value returned by the noise function applied to
/// every chunk in the game). Second component is the cached value of the noise function that
/// generated the index.
type InverseCdf = Box<[(f32, f32); WORLD_SIZE.x * WORLD_SIZE.y]>;
/// Computes the position Vec2 of a SimChunk from an index, where the index was generated by
/// uniform_noise.
fn uniform_idx_as_vec2(idx: usize) -> Vec2<i32> {
Vec2::new((idx / WORLD_SIZE.x) as i32, (idx % WORLD_SIZE.x) as i32)
}
/// Compute inverse cumulative distribution function for arbitrary function f, the hard way. We
/// pre-generate noise values prior to worldgen, then sort them in order to determine the correct
/// position in the sorted order. That lets us use `(index + 1) / (WORLDSIZE.y * WORLDSIZE.x)` as
/// a uniformly distributed (from almost-0 to 1) regularization of the chunks. That is, if we
/// apply the computed "function" F⁻¹(x, y) to (x, y) and get out p, it means that approximately
/// (100 * p)% of chunks have a lower value for F⁻¹ than p. The main purpose of doing this is to
/// make sure we are using the entire range we want, and to allow us to apply the numerous results
/// about distributions on uniform functions to the procedural noise we generate, which lets us
/// much more reliably control the *number* of features in the world while still letting us play
/// with the *shape* of those features, without having arbitrary cutoff points / discontinuities
/// (which tend to produce ugly-looking / unnatural terrain).
///
/// As a concrete example, before doing this it was very hard to tweak humidity so that either most
/// of the world wasn't dry, or most of it wasn't wet, by combining the billow noise function and
/// the computed altitude. This is because the billow noise function has a very unusual
/// distribution that is heavily skewed towards 0. By correcting for this tendency, we can start
/// with uniformly distributed billow noise and altitudes and combine them to get uniformly
/// distributed humidity, while still preserving the existing shapes that the billow noise and
/// altitude functions produce.
///
/// f takes an index, which represents the index corresponding to this chunk in any any SimChunk
/// vector returned by uniform_noise, and (for convenience) the float-translated version of those
/// coordinates.
/// f should return a value with no NaNs. If there is a NaN, it will panic. There are no other
/// conditions on f.
///
/// Returns a vec of (f32, f32) pairs consisting of the percentage of chunks with a value lower than
/// this one, and the actual noise value (we don't need to cache it, but it makes ensuring that
/// subsequent code that needs the noise value actually uses the same one we were using here
/// easier).
fn uniform_noise(f: impl Fn(usize, Vec2<f64>) -> f32) -> InverseCdf {
let mut noise = (0..WORLD_SIZE.x * WORLD_SIZE.y)
.map(|i| {
(
i,
f(
i,
(uniform_idx_as_vec2(i) * TerrainChunkSize::SIZE.map(|e| e as i32))
.map(|e| e as f64),
),
)
})
.collect::<Vec<_>>();
// sort_unstable_by is equivalent to sort_by here since we include the index in the
// comparison. We could leave out the index, but this might make the order not
// reproduce the same way between different versions of Rust (for example).
noise.sort_unstable_by(|f, g| (f.1, f.0).partial_cmp(&(g.1, g.0)).unwrap());
// Construct a vector that associates each chunk position with the 1-indexed
// position of the noise in the sorted vector (divided by the vector length).
// This guarantees a uniform distribution among the samples.
let mut uniform_noise = box [(0.0, 0.0); WORLD_SIZE.x * WORLD_SIZE.y];
let total = (WORLD_SIZE.x * WORLD_SIZE.y) as f32;
for (noise_idx, (chunk_idx, noise_val)) in noise.into_iter().enumerate() {
uniform_noise[chunk_idx] = ((1 + noise_idx) as f32 / total, noise_val);
}
uniform_noise
}
/// Calculates the smallest distance along an axis (x, y) from an edge of /// Calculates the smallest distance along an axis (x, y) from an edge of
/// the world. This value is maximal at WORLD_SIZE / 2 and minimized at the extremes /// the world. This value is maximal at WORLD_SIZE / 2 and minimized at the extremes
/// (0 or WORLD_SIZE on one or more axes). It then divides the quantity by cell_size, /// (0 or WORLD_SIZE on one or more axes). It then divides the quantity by cell_size,
@ -197,13 +46,15 @@ fn map_edge_factor(posi: usize) -> f32 {
.min(1.0) .min(1.0)
} }
/// A structure that holds cached noise values and cumulative distribution functions for the input
/// that led to those values. See the definition of InverseCdf for a description of how to
/// interpret the types of its fields.
struct GenCdf { struct GenCdf {
humid_base: InverseCdf, humid_base: InverseCdf,
temp_base: InverseCdf, temp_base: InverseCdf,
alt_base: InverseCdf, alt_base: InverseCdf,
chaos: InverseCdf, chaos: InverseCdf,
alt_main: InverseCdf, alt: InverseCdf,
alt_pre: InverseCdf,
} }
pub(crate) struct GenCtx { pub(crate) struct GenCtx {
@ -339,6 +190,8 @@ impl WorldSim {
// [0, 1.24] * [0.35, 1.0] = [0, 1.24]. // [0, 1.24] * [0.35, 1.0] = [0, 1.24].
// Sharply decreases (towards 0.35) when temperature is near desert_temp (from below), // Sharply decreases (towards 0.35) when temperature is near desert_temp (from below),
// then saturates just before it actually becomes desert. Otherwise stays at 1. // then saturates just before it actually becomes desert. Otherwise stays at 1.
// Note that this is not the *final* temperature, only the initial noise value for
// temperature.
.mul( .mul(
temp_base[posi] temp_base[posi]
.1 .1
@ -352,36 +205,41 @@ impl WorldSim {
.max(0.1) .max(0.1)
}); });
// This is the extension upwards from the base added to some extra noise from -1 to 1.
// The extra noise is multiplied by alt_main (the mountain part of the extension) clamped to
// be between 0.25 and 1, and made 60% larger (so the extra noise is between -1.6 and 1.6,
// and the final noise is never more than 160% or less than 40% of the original noise,
// depending on altitude).
// Adding this to alt_main thus yields a value between -0.4 (if alt_main = 0 and
// gen_ctx = -1) and 2.6 (if alt_main = 1 and gen_ctx = 1). When the generated small_nz
// value hits -0.625 the value crosses 0, so most of the points are above 0.
//
// Then, we add 1 and divide by 2 to get a value between 0.3 and 1.8.
let alt_main = uniform_noise(|_, wposf| {
// Extension upwards from the base. A positive number from 0 to 1 curved to be maximal
// at 0. Also to be multiplied by CONFIG.mountain_scale.
let alt_main = (gen_ctx.alt_nz.get((wposf.div(2_000.0)).into_array()) as f32)
.abs()
.powf(1.45);
(0.0 + alt_main
+ (gen_ctx.small_nz.get((wposf.div(300.0)).into_array()) as f32)
.mul(alt_main.max(0.25))
.mul(0.2))
.add(1.0)
.mul(0.5)
});
// We ignore sea level because we actually want to be relative to sea level here and want // We ignore sea level because we actually want to be relative to sea level here and want
// things in CONFIG.mountain_scale units, and we are using the version of chaos that doesn't // things in CONFIG.mountain_scale units, but otherwise this is a correct altitude
// know about temperature. Otherwise, this is a correct altitude calculation. // calculation. Note that this is using the "unadjusted" temperature.
let alt_pre = uniform_noise(|posi, _| { let alt = uniform_noise(|posi, wposf| {
(alt_base[posi].1 + alt_main[posi].1.mul(chaos[posi].1.max(0.1))) // This is the extension upwards from the base added to some extra noise from -1 to 1.
// The extra noise is multiplied by alt_main (the mountain part of the extension)
// clamped to [0.25, 1], and made 60% larger (so the extra noise is between [-1.6, 1.6],
// and the final noise is never more than 160% or less than 40% of the original noise,
// depending on altitude).
// Adding this to alt_main thus yields a value between -0.4 (if alt_main = 0 and
// gen_ctx = -1) and 2.6 (if alt_main = 1 and gen_ctx = 1). When the generated small_nz
// value hits -0.625 the value crosses 0, so most of the points are above 0.
//
// Then, we add 1 and divide by 2 to get a value between 0.3 and 1.8.
let alt_main = {
// Extension upwards from the base. A positive number from 0 to 1 curved to be
// maximal at 0. Also to be multiplied by CONFIG.mountain_scale.
let alt_main = (gen_ctx.alt_nz.get((wposf.div(2_000.0)).into_array()) as f32)
.abs()
.powf(1.45);
(0.0 + alt_main
+ (gen_ctx.small_nz.get((wposf.div(300.0)).into_array()) as f32)
.mul(alt_main.max(0.25))
.mul(0.2))
.add(1.0)
.mul(0.5)
};
// Now we can compute the final altitude using chaos.
// We multiply by chaos clamped to [0.1, 1.24] to get a value between 0.03 and 2.232 for
// alt_pre, then multiply by CONFIG.mountain_scale and add to the base and sea level to
// get an adjusted value, then multiply the whole thing by map_edge_factor
// (TODO: compute final bounds).
(alt_base[posi].1 + alt_main.mul(chaos[posi].1))
.mul(map_edge_factor(posi)) .mul(map_edge_factor(posi))
}); });
@ -390,8 +248,7 @@ impl WorldSim {
temp_base, temp_base,
alt_base, alt_base,
chaos, chaos,
alt_main, alt,
alt_pre,
}; };
let mut chunks = Vec::new(); let mut chunks = Vec::new();
@ -656,36 +513,28 @@ impl SimChunk {
let (_, alt_base) = gen_cdf.alt_base[posi]; let (_, alt_base) = gen_cdf.alt_base[posi];
let map_edge_factor = map_edge_factor(posi); let map_edge_factor = map_edge_factor(posi);
let (_, chaos) = gen_cdf.chaos[posi]; let (_, chaos) = gen_cdf.chaos[posi];
let (_, alt_pre) = gen_cdf.alt_main[posi]; let (humid_uniform, _) = gen_cdf.humid_base[posi];
let (humid_base, _) = gen_cdf.humid_base[posi]; let (alt_uniform, alt_pre) = gen_cdf.alt[posi];
let (alt_uniform, _) = gen_cdf.alt_pre[posi]; let (temp_uniform, _) = gen_cdf.temp_base[posi];
// Take the weighted average of our randomly generated base humidity, the scaled // Take the weighted average of our randomly generated base humidity, the scaled
// negative altitude, and other random variable (to add some noise) to yield the // negative altitude, and other random variable (to add some noise) to yield the
// final humidity. Note that we are using the "old" version of chaos here. // final humidity. Note that we are using the "old" version of chaos here.
const HUMID_WEIGHTS: [f32; 2] = [1.0, 1.0]; const HUMID_WEIGHTS: [f32; 2] = [1.0, 1.0];
let humidity = cdf_irwin_hall(&HUMID_WEIGHTS, [humid_base, 1.0 - alt_uniform]); let humidity = cdf_irwin_hall(&HUMID_WEIGHTS, [humid_uniform, 1.0 - alt_uniform]);
let (temp_base, temp_old) = gen_cdf.temp_base[posi];
// We also correlate temperature negatively with altitude using different weighting than we // We also correlate temperature negatively with altitude using different weighting than we
// use for humidity. // use for humidity.
const TEMP_WEIGHTS: [f32; 2] = [2.0, 1.0]; const TEMP_WEIGHTS: [f32; 2] = [2.0, 1.0];
let temp = cdf_irwin_hall(&TEMP_WEIGHTS, [temp_base, 1.0 - alt_uniform]) let temp = cdf_irwin_hall(&TEMP_WEIGHTS, [temp_uniform, 1.0 - alt_uniform])
// Convert to [-1, 1] // Convert to [-1, 1]
.sub(0.5) .sub(0.5)
.mul(2.0); .mul(2.0);
// Now we can recompute altitude using the correct verison of chaos.
// We multiply by chaos clamped to [0.1, 1.24] to get a value between 0.03 and 2.232 for
// alt_pre, then multiply by CONFIG.mountain_scale and add to the base and sea level to get
// an adjusted value, then multiply the whole thing by map_edge_factor (TODO: compute final bounds).
let alt_base = alt_base.mul(CONFIG.mountain_scale); let alt_base = alt_base.mul(CONFIG.mountain_scale);
let alt = CONFIG let alt = CONFIG
.sea_level .sea_level.mul(map_edge_factor)
.add(alt_base) .add(alt_pre.mul(CONFIG.mountain_scale));
.add(alt_pre.mul(chaos).mul(CONFIG.mountain_scale))
.mul(map_edge_factor);
let cliff = gen_ctx.cliff_nz.get((wposf.div(2048.0)).into_array()) as f32 + chaos * 0.2; let cliff = gen_ctx.cliff_nz.get((wposf.div(2048.0)).into_array()) as f32 + chaos * 0.2;
@ -695,8 +544,6 @@ impl SimChunk {
let logistic_2_base = 3.0f32.sqrt().mul(f32::consts::FRAC_2_PI); let logistic_2_base = 3.0f32.sqrt().mul(f32::consts::FRAC_2_PI);
// Assumes μ = 0, σ = 1 // Assumes μ = 0, σ = 1
let logistic_cdf = |x: f32| x.div(logistic_2_base).tanh().mul(0.5).add(0.5); let logistic_cdf = |x: f32| x.div(logistic_2_base).tanh().mul(0.5).add(0.5);
// Weighted logit sum.
let f = |humidity, density| logistic_cdf(logit(humidity) + 0.5 * logit(density));
// No trees in the ocean or with zero humidity (currently) // No trees in the ocean or with zero humidity (currently)
let tree_density = if alt <= CONFIG.sea_level + 5.0 { let tree_density = if alt <= CONFIG.sea_level + 5.0 {
@ -716,6 +563,7 @@ impl SimChunk {
} else if humidity >= 1.0 || tree_density >= 1.0 { } else if humidity >= 1.0 || tree_density >= 1.0 {
1.0 1.0
} else { } else {
// Weighted logit sum.
logistic_cdf(logit(humidity) + 0.5 * logit(tree_density)) logistic_cdf(logit(humidity) + 0.5 * logit(tree_density))
} }
// rescale to (-0.9, 0.9) // rescale to (-0.9, 0.9)

165
world/src/sim/util.rs Normal file
View File

@ -0,0 +1,165 @@
use common::{
terrain::TerrainChunkSize,
vol::VolSize,
};
use super::WORLD_SIZE;
use vek::*;
/// Computes the cumulative distribution function of the weighted sum of k independent,
/// uniformly distributed random variables between 0 and 1. For each variable i, we use weights[i]
/// as the weight to give samples[i] (the weights should all be positive).
///
/// If the precondition is met, the distribution of the result of calling this function will be
/// uniformly distributed while preserving the same information that was in the original average.
///
/// For N > 33 the function will no longer return correct results since we will overflow u32.
///
/// NOTE:
///
/// Per [1], the problem of determing the CDF of
/// the sum of uniformly distributed random variables over *different* ranges is considerably more
/// complicated than it is for the same-range case. Fortunately, it also provides a reference to
/// [2], which contains a complete derivation of an exact rule for the density function for
/// this case. The CDF is just the integral of the cumulative distribution function [3],
/// which we use to convert this into a CDF formula.
///
/// This allows us to sum weighted, uniform, independent random variables.
///
/// At some point, we should probably contribute this back to stats-rs.
///
/// 1. https://www.r-bloggers.com/sums-of-random-variables/,
/// 2. Sadooghi-Alvandi, S., A. Nematollahi, & R. Habibi, 2009.
/// On the Distribution of the Sum of Independent Uniform Random Variables.
/// Statistical Papers, 50, 171-175.
/// 3. hhttps://en.wikipedia.org/wiki/Cumulative_distribution_function
pub fn cdf_irwin_hall<const N: usize>(weights: &[f32; N], samples: [f32; N]) -> f32 {
// Let J_k = {(j_1, ... , j_k) : 1 ≤ j_1 < j_2 < ··· < j_k ≤ N }.
//
// Let A_N = Π{k = 1 to n}a_k.
//
// The density function for N ≥ 2 is:
//
// 1/(A_N * (N - 1)!) * (x^(N-1) + Σ{k = 1 to N}((-1)^k *
// Σ{(j_1, ..., j_k) ∈ J_k}(max(0, x - Σ{l = 1 to k}(a_(j_l)))^(N - 1))))
//
// So the cumulative distribution function is its integral, i.e. (I think)
//
// 1/(product{k in A}(k) * N!) * (x^N + sum(k in 1 to N)((-1)^k *
// sum{j in Subsets[A, {k}]}(max(0, x - sum{l in j}(l))^N)))
//
// which is also equivalent to
//
// (letting B_k = { a in Subsets[A, {k}] : sum {l in a} l }, B_(0,1) = 0 and
// H_k = { i : 1 ≤ 1 ≤ N! / (k! * (N - k)!) })
//
// 1/(product{k in A}(k) * N!) * sum(k in 0 to N)((-1)^k *
// sum{l in H_k}(max(0, x - B_(k,l))^N))
//
// We should be able to iterate through the whole power set
// instead, and figure out K by calling count_ones(), so we can compute the result in O(2^N)
// iterations.
let x: f32 = weights
.iter()
.zip(samples.iter())
.map(|(weight, sample)| weight * sample)
.sum();
let mut y = 0.0f32;
for subset in 0u32..(1 << N) {
// Number of set elements
let k = subset.count_ones();
// Add together exactly the set elements to get B_subset
let z = weights
.iter()
.enumerate()
.filter(|(i, _)| subset & (1 << i) as u32 != 0)
.map(|(_, k)| k)
.sum::<f32>();
// Compute max(0, x - B_subset)^N
let z = (x - z).max(0.0).powi(N as i32);
// The parity of k determines whether the sum is negated.
y += if k & 1 == 0 { z } else { -z };
}
// Divide by the product of the weights.
y /= weights.iter().product::<f32>();
// Remember to multiply by 1 / N! at the end.
y / (1..=N as i32).product::<i32>() as f32
}
/// First component of each element of the vector is the computed CDF of the noise function at this
/// index (i.e. its position in a sorted list of value returned by the noise function applied to
/// every chunk in the game). Second component is the cached value of the noise function that
/// generated the index.
pub type InverseCdf = Box<[(f32, f32); WORLD_SIZE.x * WORLD_SIZE.y]>;
/// Computes the position Vec2 of a SimChunk from an index, where the index was generated by
/// uniform_noise.
pub fn uniform_idx_as_vec2(idx: usize) -> Vec2<i32> {
Vec2::new((idx / WORLD_SIZE.x) as i32, (idx % WORLD_SIZE.x) as i32)
}
/// Compute inverse cumulative distribution function for arbitrary function f, the hard way. We
/// pre-generate noise values prior to worldgen, then sort them in order to determine the correct
/// position in the sorted order. That lets us use `(index + 1) / (WORLDSIZE.y * WORLDSIZE.x)` as
/// a uniformly distributed (from almost-0 to 1) regularization of the chunks. That is, if we
/// apply the computed "function" F⁻¹(x, y) to (x, y) and get out p, it means that approximately
/// (100 * p)% of chunks have a lower value for F⁻¹ than p. The main purpose of doing this is to
/// make sure we are using the entire range we want, and to allow us to apply the numerous results
/// about distributions on uniform functions to the procedural noise we generate, which lets us
/// much more reliably control the *number* of features in the world while still letting us play
/// with the *shape* of those features, without having arbitrary cutoff points / discontinuities
/// (which tend to produce ugly-looking / unnatural terrain).
///
/// As a concrete example, before doing this it was very hard to tweak humidity so that either most
/// of the world wasn't dry, or most of it wasn't wet, by combining the billow noise function and
/// the computed altitude. This is because the billow noise function has a very unusual
/// distribution that is heavily skewed towards 0. By correcting for this tendency, we can start
/// with uniformly distributed billow noise and altitudes and combine them to get uniformly
/// distributed humidity, while still preserving the existing shapes that the billow noise and
/// altitude functions produce.
///
/// f takes an index, which represents the index corresponding to this chunk in any any SimChunk
/// vector returned by uniform_noise, and (for convenience) the float-translated version of those
/// coordinates.
/// f should return a value with no NaNs. If there is a NaN, it will panic. There are no other
/// conditions on f.
///
/// Returns a vec of (f32, f32) pairs consisting of the percentage of chunks with a value lower than
/// this one, and the actual noise value (we don't need to cache it, but it makes ensuring that
/// subsequent code that needs the noise value actually uses the same one we were using here
/// easier).
pub fn uniform_noise(f: impl Fn(usize, Vec2<f64>) -> f32) -> InverseCdf {
let mut noise = (0..WORLD_SIZE.x * WORLD_SIZE.y)
.map(|i| {
(
i,
f(
i,
(uniform_idx_as_vec2(i) * TerrainChunkSize::SIZE.map(|e| e as i32))
.map(|e| e as f64),
),
)
})
.collect::<Vec<_>>();
// sort_unstable_by is equivalent to sort_by here since we include the index in the
// comparison. We could leave out the index, but this might make the order not
// reproduce the same way between different versions of Rust (for example).
noise.sort_unstable_by(|f, g| (f.1, f.0).partial_cmp(&(g.1, g.0)).unwrap());
// Construct a vector that associates each chunk position with the 1-indexed
// position of the noise in the sorted vector (divided by the vector length).
// This guarantees a uniform distribution among the samples.
//
// NOTE: Currently there doesn't seem to be a way to create a large fixed-size
// array on the heap without overflowing the stack unless you use placement box (at least on
// debug mode). So I want to keep using this until a better alternative is made available.
let mut uniform_noise = box [(0.0, 0.0); WORLD_SIZE.x * WORLD_SIZE.y];
let total = (WORLD_SIZE.x * WORLD_SIZE.y) as f32;
for (noise_idx, (chunk_idx, noise_val)) in noise.into_iter().enumerate() {
uniform_noise[chunk_idx] = ((1 + noise_idx) as f32 / total, noise_val);
}
uniform_noise
}