TerrainGen coverage

terraingen/
Coverage:
low: ≥ 0%
medium: ≥ 75.0%
high: ≥ 90.0%
Lines:
360 / 365
98.6%
Branches:
186 / 198
93.9%
Functions:
28 / 28
100.0%
layers/kernels/craters.h
Line Branch Exec Source
1 /// @file
2 /// @brief The kernel of crater shaping (`ApplyCraters`); a kernel header for C++ and Slang.
3 ///
4 /// The host prepares each crater (`CraterRecord`), fits its datum and gathers what varies in
5 /// length, outlines and diffusion tables, into one array of floats (`CraterData`); the kernel
6 /// shapes one sample from its offset to a crater's centre, in single precision.
7 #pragma once
8
9 #include "terraingen/noise/fractal.h"
10 #include "terraingen/util/portable.h"
11
12 #ifndef __SLANG__
13 #include "terraingen/gpu/program.h"
14 #endif
15
16 namespace terraingen {
17
18 /// The fresh profile of simple craters.
19 enum class CraterProfileKind : uint32_t {
20 FassettThomson = 0, ///< `FassettThomsonProfile`.
21 MartinParkesDunstan = 1, ///< `MartinParkesDunstanProfile`.
22 };
23
24 /// A crater profile as the kernel evaluates it: a simple profile, complex craters from a
25 /// transition diameter (`ComplexCraterProfile`), and diffusion of the fresh shape by tables
26 /// (`DiffusedProfile`).
27 struct CraterProfileKernel {
28 CraterProfileKind simple = CraterProfileKind::FassettThomson;
29 uint32_t rim = 0; ///< `FassettThomsonRim` as a number.
30 uint32_t complex = 0; ///< 1 when complex craters apply.
31 uint32_t diffused = 0; ///< 1 when the profile is diffused by tables.
32 float transitionDiameter = 0.0f;
33 float floorCoefficient = 0.0f;
34 float floorExponent = 0.0f;
35 float peakHeightCoefficient = 0.0f;
36 float peakHeightExponent = 0.0f;
37 float peakWidthCoefficient = 0.0f;
38 float peakWidthExponent = 0.0f;
39 uint32_t peakFromVolume = 0; ///< 1 when the peak's base follows from its volume.
40 float peakFitMaxDiameter = 0.0f;
41 float ejectaExponent = 3.0f;
42 float complexReachRadii = 3.0f;
43 float diffusedExtent = 1.5f; ///< Half the diffused domain's reach, in diameters.
44 uint32_t diffusedCells = 0;
45 uint32_t agesOffset = 0; ///< Where the tabulated diffusion ages start in `CraterData`.
46 uint32_t ageCount = 0;
47 };
48
49 /// How craters are shaped, as `CraterApplyOptions` and the profile give it.
50 struct CraterShapeKernel {
51 CraterProfileKernel profile;
52 float reachRadii = 1.5f; ///< `CraterReachRadii`.
53 float profileReachRadii = 1.5f; ///< The profile's own reach.
54 float blendStart = 0.8f;
55 float blendEnd = 1.2f;
56 uint32_t blendShape = 0;
57 float interiorReplacement = 1.0f;
58 float mantling = 0.0f;
59 float ejectaReachRadii = 3.0f;
60 float rimSigma = 0.0f;
61 float rimInner = 0.8f;
62 float rimOuter = 1.2f;
63 float rimStrength = 1.0f;
64 uint32_t rimShape = 0;
65 float noiseAmplitude = 0.0f;
66 float noiseWavelength = 0.5f;
67 int noiseOctaves = 4;
68 float noiseGain = 0.5f;
69 float noiseRimAmplitude = 0.0f;
70 float noiseRimWavelength = 0.1f;
71 float noiseRimInner = 0.8f;
72 float noiseRimOuter = 1.2f;
73 float wallStreaks = 0.0f;
74 float wallStreakWidth = 0.15f;
75 float ejecta = 0.0f;
76 float ejectaLineation = 0.15f;
77 float minFreshness = 0.5f;
78 };
79
80 /// One crater ready to shape.
81 struct CraterRecord {
82 float x = 0.0f; ///< Centre east of the dispatch's origin (m); read on the GPU only.
83 float y = 0.0f; ///< Centre north of the dispatch's origin (m); read on the GPU only.
84 float diameter = 0.0f;
85 float resolution = 1.0f;
86 float freshDepth = 0.0f;
87 float freshRim = 0.0f;
88 float currentDepth = 0.0f;
89 float currentRim = 0.0f;
90 float diffusion = 0.0f;
91 float ellipticity = 1.0f;
92 float heading = 0.0f;
93 uint32_t firstHarmonic = 2;
94 uint32_t harmonicCount = 0;
95 uint32_t harmonicOffset = 0; ///< Amplitude and phase pairs in `CraterData`.
96 uint32_t tableCount = 0;
97 uint32_t tableOffset = 0; ///< The outline's deviations in `CraterData`.
98 uint32_t diffusionOffset = 0; ///< Its diffusion table in `CraterData`.
99 uint32_t hasDatum = 0;
100 float datumA = 0.0f;
101 float datumB = 0.0f;
102 float datumC = 0.0f;
103 uint32_t noiseKey = 0;
104 uint32_t textureKey = 0;
105 uint32_t markCount = 0;
106 /// The marks' inner, outer and fade radii, then their knots, in `CraterData`.
107 uint32_t markOffset = 0;
108 float shiftX = 0.0f; ///< `CraterSlopeShape::shift` east (m).
109 float shiftY = 0.0f; ///< `CraterSlopeShape::shift` north (m).
110 float upslopeX = 0.0f; ///< `CraterSlopeShape::upslope` east.
111 float upslopeY = 0.0f; ///< `CraterSlopeShape::upslope` north.
112 };
113
114 #ifdef __SLANG__
115 /// Where the variable-length values of craters are read from: a buffer the dispatch binds.
116 typedef StructuredBuffer<float> CraterData;
117
118 inline float CraterValue(CraterData data, uint32_t index) {
119 return data[index];
120 }
121
122 /// Crater records, floats and unsigned integers a dispatch binds, indexed as arrays.
123 typedef StructuredBuffer<CraterRecord> CraterRecords;
124 typedef StructuredBuffer<float> KernelFloats;
125 typedef StructuredBuffer<uint> KernelUints;
126 #else
127 /// Crater records, floats and unsigned integers, indexed as arrays.
128 using CraterRecords = const CraterRecord*;
129 using KernelFloats = const float*;
130 using KernelUints = const uint32_t*;
131
132 /// Where the variable-length values of craters are read from.
133 struct CraterData {
134 const float* values = nullptr;
135 };
136
137 1032219415 inline float CraterValue(CraterData data, uint32_t index) {
138 1032219415 return data.values[index];
139 }
140 #endif
141
142 /// 1 - `t` shaped by `BlendShape` number `shape` (linear, smoothstep, smootherstep, Gaussian), for
143 /// `t` in [0, 1].
144 9954953 inline float KernelFalloffShape(float t, uint32_t shape) {
145
2/2
✓ Branch 0 taken 12287 times.
✓ Branch 1 taken 9942666 times.
 9954953 if (shape == 0u) {
146 12287 return 1.0f - t;
147 }
148
2/2
✓ Branch 0 taken 252260 times.
✓ Branch 1 taken 9690406 times.
 9942666 if (shape == 1u) {
149 252260 return 1.0f - t * t * (3.0f - 2.0f * t);
150 }
151
2/2
✓ Branch 0 taken 9685006 times.
✓ Branch 1 taken 5400 times.
 9690406 if (shape == 2u) {
152 9685006 return 1.0f - t * t * t * (t * (t * 6.0f - 15.0f) + 10.0f);
153 }
154 5400 const float tail = exp(-4.5f);
155 5400 const float z = 3.0f * t;
156 5400 return (exp(-0.5f * z * z) - tail) / (1.0f - tail);
157 }
158
159 /// 1 up to `start`, 0 from `end`, shaped by `shape` between (`Falloff`).
160 42986002 inline float KernelFalloff(float x, float start, float end, uint32_t shape) {
161
2/2
✓ Branch 0 taken 2527202 times.
✓ Branch 1 taken 40458800 times.
 42986002 if (x <= start) {
162 2527202 return 1.0f;
163 }
164
2/2
✓ Branch 0 taken 30503847 times.
✓ Branch 1 taken 9954953 times.
 40458800 if (x >= end) {
165 30503847 return 0.0f;
166 }
167 9954953 return KernelFalloffShape((x - start) / (end - start), shape);
168 }
169
170 15258947 inline float KernelSmoothstep(float t) {
171
1/1
✓ Branch 1 taken 15258947 times.
 15258947 const float u = clamp(t, 0.0f, 1.0f);
172 15258947 return u * u * (3.0f - 2.0f * u);
173 }
174
175 /// The direction of (`dx`, `dy`) in radians from east, 0 at the origin: fast maths on some GPUs
176 /// leaves `atan2(0, 0)` undefined.
177 27762947 inline float KernelAngle(float dy, float dx) {
178
4/4
✓ Branch 0 taken 29233 times.
✓ Branch 1 taken 27733714 times.
✓ Branch 2 taken 260 times.
✓ Branch 3 taken 28973 times.
 27762947 if (dx == 0.0f && dy == 0.0f) {
179 260 return 0.0f;
180 }
181 27762687 return atan2(dy, dx);
182 }
183
184 /// `EjectaThicknessShare`.
185 8312658 inline float KernelEjectaShare(float radii, float reachRadii) {
186
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 8312658 times.
 8312658 if (radii <= 1.0f) {
187 ✗ return 1.0f;
188 }
189
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 8312658 times.
 8312658 if (radii >= reachRadii) {
190 ✗ return 0.0f;
191 }
192 8312658 const float atReach = 1.0f / (reachRadii * reachRadii * reachRadii);
193 8312658 return (1.0f / (radii * radii * radii) - atReach) / (1.0f - atReach);
194 }
195
196 /// `FassettThomsonProfile::Height` over the diameter.
197 9780118 inline float FassettThomsonShape(float x, uint32_t rim, float freshRim, float currentDepth) {
198
2/2
✓ Branch 0 taken 373323 times.
✓ Branch 1 taken 9406795 times.
 9780118 if (x > 1.5f) {
199 373323 return 0.0f;
200 }
201 9406795 const float interior =
202 9406795 -0.228809953f + x * (0.227533882f + x * (0.083116795f - x * 0.039499407f));
203 9406795 const float exterior = 0.188253307f + x * (-0.187050452f + x * (0.01844746f + x * 0.01505647f));
204 9406795 float h = 0.0f;
205
2/2
✓ Branch 0 taken 78943 times.
✓ Branch 1 taken 9327852 times.
 9406795 if (rim == 2u) {
206
2/2
✓ Branch 0 taken 37087 times.
✓ Branch 1 taken 41856 times.
 78943 if (x <= 0.98f) {
207 37087 h = interior;
208
2/2
✓ Branch 0 taken 3037 times.
✓ Branch 1 taken 38819 times.
 41856 } else if (x <= 1.02f) {
209 3037 h = freshRim;
210 } else {
211 38819 h = exterior;
212 }
213 } else {
214
2/2
✓ Branch 0 taken 9248909 times.
✓ Branch 1 taken 78943 times.
 9327852 const float edge = rim == 0u ? 0.98f : 1.0f;
215
2/2
✓ Branch 0 taken 176031 times.
✓ Branch 1 taken 9151821 times.
 9327852 if (x <= 0.2f) {
216 176031 h = -0.181f;
217
2/2
✓ Branch 0 taken 3881806 times.
✓ Branch 1 taken 5270015 times.
 9151821 } else if (x <= edge) {
218 3881806 h = interior;
219 } else {
220 5270015 h = exterior;
221 }
222 }
223 9406795 return max(h, freshRim - currentDepth);
224 }
225
226 /// `MartinParkesDunstanProfile::Height` (m); the host checks the state.
227 17287777 inline float MartinParkesDunstanHeight(float distance, float diameter, KERNEL_IN(CraterRecord) c) {
228 17287777 const float depth = c.currentDepth * diameter;
229 17287777 const float rim = c.currentRim * diameter;
230 17287777 const float freshRim = c.freshRim * diameter;
231 17287777 const float radius = 0.5f * diameter;
232 17287777 const float x = distance / radius - 1.0f;
233
2/2
✓ Branch 0 taken 656489 times.
✓ Branch 1 taken 16631288 times.
 17287777 if (x >= 1.0f) {
234 656489 return 0.0f;
235 }
236 16631288 const float k = freshRim - rim + depth;
237 16631288 const float alpha = (rim - freshRim) / k;
238
2/2
✓ Branch 0 taken 3520799 times.
✓ Branch 1 taken 13110489 times.
 16631288 if (x <= alpha) {
239 3520799 return (k * x + 2.0f * k) * x + freshRim;
240 }
241
2/2
✓ Branch 0 taken 1210563 times.
✓ Branch 1 taken 11899926 times.
 13110489 if (x <= 0.0f) {
242 1210563 return k * (alpha + 1.0f) / alpha * x * x + rim;
243 }
244 11899926 const float ejectaScale = 0.14f * pow(radius, 0.74f);
245 11899926 const float offset = rim - ejectaScale;
246 11899926 const float gradient = -ejectaScale / 8.0f - offset;
247 11899926 const float beta = min(1.0f, 1.5f * (freshRim - rim) / k);
248 11899926 const float atBeta =
249 11899926 ejectaScale / ((beta + 1.0f) * (beta + 1.0f) * (beta + 1.0f)) + offset + gradient * beta;
250
2/2
✓ Branch 0 taken 2232261 times.
✓ Branch 1 taken 9667665 times.
 11899926 if (x <= beta) {
251 2232261 const float a = (rim - atBeta) / (2.0f * beta * beta * beta);
252 2232261 return (a * x - 3.0f * a * beta) * x * x + rim;
253 }
254 9667665 return ejectaScale / ((x + 1.0f) * (x + 1.0f) * (x + 1.0f)) + offset + gradient * x;
255 }
256
257 /// The simple profile's height (m) in the fresh state when `fresh` is 1.
258 27067895 inline float SimpleProfileHeight(KERNEL_IN(CraterProfileKernel) p, float distance, float diameter,
259 KERNEL_IN(CraterRecord) c) {
260
2/2
✓ Branch 0 taken 17287777 times.
✓ Branch 1 taken 9780118 times.
 27067895 if (p.simple == CraterProfileKind::MartinParkesDunstan) {
261 17287777 return MartinParkesDunstanHeight(distance, diameter, c);
262 }
263 9780118 return FassettThomsonShape(distance / (0.5f * diameter), p.rim, c.freshRim, c.currentDepth) *
264 9780118 diameter;
265 }
266
267 /// y (m) of a power law `coefficient` D^`exponent` in kilometres, or m^3 with `volume`.
268 388314 inline float KilometreLaw(float coefficient, float exponent, float diameter, bool volume) {
269 388314 const float kilometres = coefficient * pow(diameter / 1000.0f, exponent);
270
2/2
✓ Branch 0 taken 121776 times.
✓ Branch 1 taken 266538 times.
 388314 return volume ? kilometres * 1e9f : kilometres * 1000.0f;
271 }
272
273 /// `ComplexCraterProfile::Height` of a complex crater (m).
274 269221 inline float ComplexProfileHeight(KERNEL_IN(CraterProfileKernel) p, float distance, float diameter,
275 KERNEL_IN(CraterRecord) c) {
276 269221 const float radius = 0.5f * diameter;
277 269221 const float x = distance / radius;
278 269221 const float depth = c.currentDepth * diameter;
279 269221 const float rim = c.currentRim * diameter;
280
2/2
✓ Branch 0 taken 139783 times.
✓ Branch 1 taken 129438 times.
 269221 if (x >= 1.0f) {
281
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 139783 times.
 139783 if (x >= p.complexReachRadii) {
282 ✗ return 0.0f;
283 }
284 139783 const float edge = pow(p.complexReachRadii, -p.ejectaExponent);
285 139783 return rim * (pow(x, -p.ejectaExponent) - edge) / (1.0f - edge);
286 }
287 const float floorShare =
288 129438 min(0.9f, KilometreLaw(p.floorCoefficient, p.floorExponent, diameter, false) / diameter);
289 129438 float height = rim - depth;
290
2/2
✓ Branch 0 taken 28259 times.
✓ Branch 1 taken 101179 times.
 129438 if (x > floorShare) {
291 28259 const float t = (x - floorShare) / (1.0f - floorShare);
292 28259 height += depth * t * t;
293 }
294 const float fitted =
295
1/2
✓ Branch 0 taken 129438 times.
✗ Branch 1 not taken.
 129438 p.peakFitMaxDiameter > 0.0f ? min(diameter, p.peakFitMaxDiameter) : diameter;
296 const float peakHeight =
297 129438 KilometreLaw(p.peakHeightCoefficient, p.peakHeightExponent, fitted, false);
298 129438 float base = 0.0f;
299
2/2
✓ Branch 0 taken 121776 times.
✓ Branch 1 taken 7662 times.
 129438 if (p.peakFromVolume == 1u) {
300 const float volume =
301 121776 KilometreLaw(p.peakWidthCoefficient, p.peakWidthExponent, fitted, true);
302 121776 base = sqrt(3.0f * volume / (3.14159265358979f * peakHeight));
303 } else {
304 7662 base = 0.5f * KilometreLaw(p.peakWidthCoefficient, p.peakWidthExponent, fitted, false);
305 }
306 129438 base *= diameter / fitted;
307
2/2
✓ Branch 0 taken 94160 times.
✓ Branch 1 taken 35278 times.
 129438 if (distance < base) {
308
1/2
✓ Branch 0 taken 94160 times.
✗ Branch 1 not taken.
 94160 const float kept = c.freshDepth > 0.0f ? c.currentDepth / c.freshDepth : 1.0f;
309 94160 const float u = 1.0f - (distance / base) * (distance / base);
310 94160 height += kept * peakHeight * u * u;
311 }
312 129438 return height;
313 }
314
315 /// The fresh profile's height (m): simple or complex.
316 27337116 inline float FreshProfileHeight(KERNEL_IN(CraterProfileKernel) p, float distance, float diameter,
317 KERNEL_IN(CraterRecord) c) {
318
4/4
✓ Branch 0 taken 437438 times.
✓ Branch 1 taken 26899678 times.
✓ Branch 2 taken 269221 times.
✓ Branch 3 taken 168217 times.
 27337116 if (p.complex == 1u && diameter >= p.transitionDiameter) {
319 269221 return ComplexProfileHeight(p, distance, diameter, c);
320 }
321 27067895 return SimpleProfileHeight(p, distance, diameter, c);
322 }
323
324 /// `DiffusedProfile::Height` (m), reading the crater's table from `data`.
325 35484607 inline float DiffusedProfileHeight(KERNEL_IN(CraterProfileKernel) p, CraterData data,
326 float distance, float diameter, KERNEL_IN(CraterRecord) c) {
327
2/2
✓ Branch 0 taken 1435576 times.
✓ Branch 1 taken 34049031 times.
 35484607 if (!(c.diffusion > 0.0f)) {
328 1435576 CraterRecord fresh = c;
329 1435576 fresh.currentDepth = c.freshDepth;
330 1435576 fresh.currentRim = c.freshRim;
331
1/1
✓ Branch 1 taken 1435576 times.
 1435576 return FreshProfileHeight(p, distance, diameter, fresh);
332 }
333 34049031 const float x = distance / diameter;
334
2/2
✓ Branch 0 taken 3 times.
✓ Branch 1 taken 34049028 times.
 34049031 if (x >= p.diffusedExtent) {
335 3 return 0.0f;
336 }
337 34049028 const float cell = x / p.diffusedExtent * float(p.diffusedCells);
338 34049028 const uint32_t i = min(uint32_t(cell), p.diffusedCells - 1u);
339 34049028 const float t = cell - float(i);
340 34049028 uint32_t row = 0u;
341 34049028 float weight = 0.0f;
342 34049028 const uint32_t last = p.ageCount - 1u;
343 34049028 const float first = CraterValue(data, p.agesOffset + 1u);
344
2/2
✓ Branch 1 taken 238401 times.
✓ Branch 2 taken 33810627 times.
 34049028 if (c.diffusion >= CraterValue(data, p.agesOffset + last)) {
345 238401 row = last;
346
2/2
✓ Branch 0 taken 52782 times.
✓ Branch 1 taken 33757845 times.
 33810627 } else if (c.diffusion < first) {
347 52782 weight = c.diffusion / first;
348 } else {
349 33757845 uint32_t low = 1u;
350 33757845 uint32_t high = last;
351
2/2
✓ Branch 0 taken 249828272 times.
✓ Branch 1 taken 33757845 times.
 283586117 while (high - low > 1u) {
352 249828272 const uint32_t middle = (low + high) / 2u;
353
2/2
✓ Branch 1 taken 140830406 times.
✓ Branch 2 taken 108997866 times.
 249828272 if (CraterValue(data, p.agesOffset + middle) <= c.diffusion) {
354 140830406 low = middle;
355 } else {
356 108997866 high = middle;
357 }
358 }
359 33757845 row = low;
360 33757845 const float below = CraterValue(data, p.agesOffset + low);
361 33757845 const float above = CraterValue(data, p.agesOffset + high);
362 33757845 weight = log(c.diffusion / below) / log(above / below);
363 }
364 34049028 const uint32_t width = p.diffusedCells + 1u;
365 34049028 const uint32_t at = c.diffusionOffset + row * width + i;
366 34049028 float value = (1.0f - t) * CraterValue(data, at) + t * CraterValue(data, at + 1u);
367
2/2
✓ Branch 0 taken 33810627 times.
✓ Branch 1 taken 238401 times.
 34049028 if (weight > 0.0f) {
368 const float next =
369 33810627 (1.0f - t) * CraterValue(data, at + width) + t * CraterValue(data, at + width + 1u);
370 33810627 value += weight * (next - value);
371 }
372 34049028 return diameter * value;
373 }
374
375 /// The profile's height (m) at `distance`.
376 61386147 inline float KernelProfileHeight(KERNEL_IN(CraterShapeKernel) k, CraterData data, float distance,
377 float diameter, KERNEL_IN(CraterRecord) c) {
378
2/2
✓ Branch 0 taken 35484607 times.
✓ Branch 1 taken 25901540 times.
 61386147 if (k.profile.diffused == 1u) {
379 35484607 return DiffusedProfileHeight(k.profile, data, distance, diameter, c);
380 }
381 25901540 return FreshProfileHeight(k.profile, distance, diameter, c);
382 }
383
384 /// The profile's height with its rim rounded (`RimSmoothing`), by 9-point Gauss-Hermite.
385 52926992 inline float KernelRoundedHeight(KERNEL_IN(CraterShapeKernel) k, CraterData data, float distance,
386 float diameter, KERNEL_IN(CraterRecord) c) {
387
1/1
✓ Branch 1 taken 52926992 times.
 52926992 const float height = KernelProfileHeight(k, data, distance, diameter, c);
388
2/2
✓ Branch 0 taken 49987287 times.
✓ Branch 1 taken 2939705 times.
 52926992 if (k.rimSigma == 0.0f) {
389 49987287 return height;
390 }
391 2939705 const float radius = 0.5f * diameter;
392 2939705 const float radii = distance / radius;
393 const float weight =
394 radii <= 1.0f
395
2/2
✓ Branch 0 taken 639780 times.
✓ Branch 1 taken 2299925 times.
 5239630 ? k.rimStrength * KernelFalloff(2.0f - radii, 1.0f, 2.0f - k.rimInner, k.rimShape)
396 2299925 : k.rimStrength * KernelFalloff(radii, 1.0f, k.rimOuter, k.rimShape);
397
2/2
✓ Branch 0 taken 2392928 times.
✓ Branch 1 taken 546777 times.
 2939705 if (weight == 0.0f) {
398 2392928 return height;
399 }
400 546777 const float nodes[5] = {0.0f, 0.723551018752838f, 1.46855328921667f, 2.26658058453184f,
401 3.19099320178153f};
402 546777 const float weights[5] = {0.720235215606051f, 0.432651559002556f, 0.0884745273943766f,
403 0.00494362427553695f, 0.0000396069772632645f};
404 546777 const float step = 1.41421356237310f * k.rimSigma * radius;
405 546777 float smoothed = weights[0] * height;
406
2/2
✓ Branch 0 taken 2187108 times.
✓ Branch 1 taken 546777 times.
 2733885 for (int i = 1; i < 5; ++i) {
407 2187108 smoothed +=
408
1/1
✓ Branch 2 taken 2187108 times.
 2187108 weights[i] * KernelProfileHeight(k, data, abs(distance - step * nodes[i]), diameter, c);
409 2187108 smoothed +=
410
1/1
✓ Branch 2 taken 2187108 times.
 2187108 weights[i] * KernelProfileHeight(k, data, abs(distance + step * nodes[i]), diameter, c);
411 }
412 546777 smoothed /= 1.77245385090552f;
413 546777 return height + weight * (smoothed - height);
414 }
415
416 /// Fractal parameters with every field set, as kernels need: Slang does not apply a struct's
417 /// default values to a local.
418 11486975 inline FractalParams KernelFractalParams(uint32_t seed, float wavelength, int octaves, float gain) {
419 11486975 FractalParams params;
420 11486975 params.seed = seed;
421 11486975 params.wavelength = wavelength;
422 11486975 params.amplitude = 1.0f;
423 11486975 params.octaves = octaves;
424 11486975 params.lacunarity = 2.0f;
425 11486975 params.gain = gain;
426 11486975 params.basis = NoiseBasis::Perlin;
427 11486975 params.type = FractalType::Fbm;
428 11486975 params.offset = 1.0f;
429 11486975 return params;
430 }
431
432 /// The crater noise (`CraterNoise`) at `dx`, `dy` metres from the centre and `radii` radii out.
433 52926992 inline float KernelNoiseHeight(KERNEL_IN(CraterShapeKernel) k, KERNEL_IN(CraterRecord) c, float dx,
434 float dy, float radii) {
435
4/4
✓ Branch 0 taken 46308072 times.
✓ Branch 1 taken 6618920 times.
✓ Branch 2 taken 46308012 times.
✓ Branch 3 taken 60 times.
 52926992 if (k.noiseAmplitude == 0.0f && k.noiseRimAmplitude == 0.0f) {
436 46308012 return 0.0f;
437 }
438 6618980 const float radius = 0.5f * c.diameter;
439 const float footprint =
440
2/2
✓ Branch 0 taken 843737 times.
✓ Branch 1 taken 5775243 times.
 6618980 radii <= 1.0f ? 1.0f : KernelFalloff(radii, 1.0f, k.profileReachRadii, 2u);
441
2/2
✓ Branch 0 taken 843737 times.
✓ Branch 1 taken 5775243 times.
 6618980 const float rim = radii <= 1.0f ? KernelFalloff(2.0f - radii, 1.0f, 2.0f - k.noiseRimInner, 2u)
442 5775243 : KernelFalloff(radii, 1.0f, k.noiseRimOuter, 2u);
443 6618980 float height = 0.0f;
444
4/4
✓ Branch 0 taken 6618920 times.
✓ Branch 1 taken 60 times.
✓ Branch 2 taken 6605267 times.
✓ Branch 3 taken 13653 times.
 6618980 if (k.noiseAmplitude > 0.0f && footprint > 0.0f) {
445 6605267 const FractalParams params = KernelFractalParams(c.noiseKey, k.noiseWavelength * radius,
446 6605267 k.noiseOctaves, k.noiseGain);
447
1/1
✓ Branch 1 taken 6605267 times.
 6605267 height += k.noiseAmplitude * footprint * Fractal2(dx, dy, params);
448 }
449
4/4
✓ Branch 0 taken 4939266 times.
✓ Branch 1 taken 1679714 times.
✓ Branch 2 taken 453527 times.
✓ Branch 3 taken 4485739 times.
 6618980 if (k.noiseRimAmplitude > 0.0f && rim > 0.0f) {
450 453527 const FractalParams params = KernelFractalParams(
451 453527 c.noiseKey + 1u, k.noiseRimWavelength * radius, k.noiseOctaves, k.noiseGain);
452
1/1
✓ Branch 1 taken 453527 times.
 453527 height += k.noiseRimAmplitude * rim * Fractal2(dx, dy, params);
453 }
454 6618980 return c.diameter * height;
455 }
456
457 /// Fractal noise streaked radially (`RadialStreaks`).
458 4428181 inline float KernelRadialStreaks(float radii, float angle, float across, float along, int octaves,
459 uint32_t seed) {
460 4428181 const FractalParams params = KernelFractalParams(seed, 1.0f, octaves, 0.5f);
461 4428181 const float ring = 1.0f / across;
462
1/1
✓ Branch 3 taken 4428181 times.
 8856362 return Fractal3Float(ring * cos(angle), ring * sin(angle), radii / along, params);
463 }
464
465 /// The wall gullies and ejecta texture (`CraterTexture`).
466 52926992 inline float KernelTextureHeight(KERNEL_IN(CraterShapeKernel) k, KERNEL_IN(CraterRecord) c,
467 float dx, float dy, float radii) {
468
4/4
✓ Branch 0 taken 46220349 times.
✓ Branch 1 taken 6706643 times.
✓ Branch 2 taken 46207651 times.
✓ Branch 3 taken 12698 times.
 52926992 if (k.wallStreaks == 0.0f && k.ejecta == 0.0f) {
469 46207651 return 0.0f;
470 }
471
1/2
✓ Branch 0 taken 6719341 times.
✗ Branch 1 not taken.
 6719341 const float freshness = c.freshDepth > 0.0f ? c.currentDepth / c.freshDepth : 1.0f;
472 6719341 const float fade = KernelSmoothstep((freshness - k.minFreshness) / (1.0f - k.minFreshness));
473
2/2
✓ Branch 0 taken 2165927 times.
✓ Branch 1 taken 4553414 times.
 6719341 if (!(fade > 0.0f)) {
474 2165927 return 0.0f;
475 }
476 4553414 const float angle = KernelAngle(dy, dx);
477 4553414 float height = 0.0f;
478
6/6
✓ Branch 0 taken 4540716 times.
✓ Branch 1 taken 12698 times.
✓ Branch 2 taken 4416974 times.
✓ Branch 3 taken 123742 times.
✓ Branch 4 taken 473417 times.
✓ Branch 5 taken 3943557 times.
 4553414 if (k.wallStreaks > 0.0f && radii > 0.25f && radii < 1.0f) {
479 const float band =
480 473417 KernelSmoothstep((radii - 0.25f) / 0.35f) * KernelSmoothstep((1.0f - radii) / 0.1f);
481 946834 height += k.wallStreaks * c.diameter * band *
482 473417 KernelRadialStreaks(radii, angle, k.wallStreakWidth, 0.5f, 2, c.textureKey);
483 }
484
5/6
✓ Branch 0 taken 4553198 times.
✓ Branch 1 taken 216 times.
✓ Branch 2 taken 3954764 times.
✓ Branch 3 taken 598434 times.
✓ Branch 4 taken 3954764 times.
✗ Branch 5 not taken.
 4553414 if (k.ejecta > 0.0f && radii > 1.0f && radii < k.ejectaReachRadii) {
485 3954764 const float radius = 0.5f * c.diameter;
486 const float thickness =
487 3954764 0.14f * pow(radius, 0.74f) * KernelEjectaShare(radii, k.ejectaReachRadii);
488 3954764 const float rise = KernelSmoothstep((radii - 1.0f) / 0.1f);
489 7909528 height += k.ejecta * thickness * rise *
490 3954764 KernelRadialStreaks(radii, angle, k.ejectaLineation, 2.0f * k.ejectaLineation, 3,
491 3954764 c.textureKey + 1u);
492 }
493 4553414 return fade * height;
494 }
495
496 /// How far the outline reaches toward (`dx`, `dy`), over the round radius (`CraterOutline::Scale`).
497 23015362 inline float KernelOutlineScale(CraterData data, KERNEL_IN(CraterRecord) c, float dx, float dy) {
498 23015362 const float theta = KernelAngle(dy, dx);
499 23015362 float scale = 1.0f;
500
2/2
✓ Branch 0 taken 2247557 times.
✓ Branch 1 taken 20767805 times.
 23015362 if (c.ellipticity != 1.0f) {
501 2247557 const float along = sqrt(c.ellipticity);
502 2247557 const float across = 1.0f / along;
503 2247557 const float u = cos(theta - c.heading) / along;
504 2247557 const float v = sin(theta - c.heading) / across;
505 2247557 scale = 1.0f / sqrt(u * u + v * v);
506 }
507 23015362 float wobble = 1.0f;
508
2/2
✓ Branch 0 taken 2000 times.
✓ Branch 1 taken 23013362 times.
 23015362 if (c.tableCount > 0u) {
509 2000 const float turns = theta / 6.28318530717959f;
510 2000 const float position = (turns - floor(turns)) * float(c.tableCount);
511 2000 const uint32_t index = min(uint32_t(position), c.tableCount - 1u);
512 2000 const float t = position - float(index);
513 2000 wobble += (1.0f - t) * CraterValue(data, c.tableOffset + index) +
514 2000 t * CraterValue(data, c.tableOffset + (index + 1u) % c.tableCount);
515 }
516
2/2
✓ Branch 0 taken 252016567 times.
✓ Branch 1 taken 23015362 times.
 275031929 for (uint32_t n = 0u; n < c.harmonicCount; ++n) {
517 252016567 const float order = float(c.firstHarmonic + n);
518 504033134 wobble += CraterValue(data, c.harmonicOffset + 2u * n) *
519 252016567 cos(order * theta + CraterValue(data, c.harmonicOffset + 2u * n + 1u));
520 }
521 23015362 return scale * wobble;
522 }
523
524 /// The shift in radii crater `c`'s wall marks add at `radii` toward (`dx`, `dy`)
525 /// (`CraterOutline::MarkShift`).
526 2081423 inline float KernelMarkShift(CraterData data, KERNEL_IN(CraterRecord) c, float dx, float dy,
527 float radii) {
528 2081423 const float inner = CraterValue(data, c.markOffset);
529 2081423 const float outer = CraterValue(data, c.markOffset + 1u);
530 2081423 const float fade = CraterValue(data, c.markOffset + 2u);
531
4/4
✓ Branch 0 taken 1982482 times.
✓ Branch 1 taken 98941 times.
✓ Branch 2 taken 1788311 times.
✓ Branch 3 taken 194171 times.
 2081423 if (!(radii > inner - fade) || !(radii < outer + fade)) {
532 1887252 return 0.0f;
533 }
534 194171 const float window = KernelSmoothstep((radii - inner + fade) / fade) *
535 194171 KernelSmoothstep((outer + fade - radii) / fade);
536 194171 const float turns = KernelAngle(dy, dx) / 6.28318530717959f;
537 194171 const uint32_t size = c.markCount;
538 194171 const float position = (turns - floor(turns)) * float(size);
539 194171 const uint32_t index = min(uint32_t(position), size - 1u);
540 194171 const float t = position - float(index);
541 194171 const uint32_t knots = c.markOffset + 3u;
542 194171 const float p0 = CraterValue(data, knots + (index + size - 1u) % size);
543 194171 const float p1 = CraterValue(data, knots + index);
544 194171 const float p2 = CraterValue(data, knots + (index + 1u) % size);
545 194171 const float p3 = CraterValue(data, knots + (index + 2u) % size);
546 194171 const float shift =
547 194171 0.5f * (2.0f * p1 + (p2 - p0) * t + (2.0f * p0 - 5.0f * p1 + 4.0f * p2 - p3) * t * t +
548 194171 (3.0f * (p1 - p2) + p3 - p0) * t * t * t);
549 194171 return window * shift;
550 }
551
552 /// The distance (m) crater `c`'s round profile is read at, `dx` east and `dy` north of its centre
553 /// (`CraterProfileDistance`).
554 81839007 inline float KernelProfileDistance(CraterData data, KERNEL_IN(CraterRecord) c, float dx, float dy) {
555 81839007 float distance = sqrt(dx * dx + dy * dy);
556
6/8
✓ Branch 0 taken 79591450 times.
✓ Branch 1 taken 2247557 times.
✓ Branch 2 taken 58823645 times.
✓ Branch 3 taken 20767805 times.
✓ Branch 4 taken 58823645 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 58823645 times.
 81839007 if (c.ellipticity != 1.0f || c.harmonicCount > 0u || c.tableCount > 0u || c.markCount > 0u) {
557 23015362 distance /= KernelOutlineScale(data, c, dx, dy);
558 }
559
2/2
✓ Branch 0 taken 2081423 times.
✓ Branch 1 taken 79757584 times.
 81839007 if (c.markCount > 0u) {
560 2081423 const float radius = 0.5f * c.diameter;
561 2081423 distance += radius * KernelMarkShift(data, c, dx, dy, distance / radius);
562 }
563 81839007 return distance;
564 }
565
566 /// How much of crater `c`'s floor shift applies at `dx` east and `dy` north of its centre
567 /// (`SlopeShiftedOffset`): all within 0.3 radii, fading along a smoothstep to none at the rim.
568 3249666 inline float KernelSlopeShiftWeight(KERNEL_IN(CraterRecord) c, float dx, float dy) {
569 3249666 const float radii = sqrt(dx * dx + dy * dy) / (0.5f * c.diameter);
570 3249666 return 1.0f - KernelSmoothstep((radii - 0.3f) / 0.7f);
571 }
572
573 /// The share of crater `c`'s raised relief kept at `dx` east and `dy` north of its centre
574 /// (`SlopeRaisedShare`).
575 1517783 inline float KernelSlopeRaisedShare(KERNEL_IN(CraterRecord) c, float dx, float dy) {
576 1517783 const float distance = sqrt(dx * dx + dy * dy);
577
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1517783 times.
 1517783 if (!(distance > 0.0f)) {
578 ✗ return 1.0f;
579 }
580 1517783 const float removed = sqrt(c.upslopeX * c.upslopeX + c.upslopeY * c.upslopeY);
581 1517783 const float cosine = max(0.0f, (c.upslopeX * dx + c.upslopeY * dy) / (removed * distance));
582 1517783 return 1.0f - removed * cosine * cosine;
583 }
584
585 /// `height` after crater `c`, at `dx` east and `dy` north of its centre in metres
586 /// (`ApplyCraterTo`); unchanged beyond its reach.
587 73007976 inline float ApplyCraterKernel(KERNEL_IN(CraterShapeKernel) k, CraterData data,
588 KERNEL_IN(CraterRecord) c, float height, float dx, float dy) {
589 73007976 const float diameter = c.diameter;
590 73007976 const float radius = 0.5f * diameter;
591
4/4
✓ Branch 0 taken 70240593 times.
✓ Branch 1 taken 2767383 times.
✓ Branch 2 taken 157188 times.
✓ Branch 3 taken 70083405 times.
 73007976 const bool shifted = c.shiftX != 0.0f || c.shiftY != 0.0f;
592
2/2
✓ Branch 0 taken 2924571 times.
✓ Branch 1 taken 70083405 times.
 73007976 const float shift = shifted ? KernelSlopeShiftWeight(c, dx, dy) : 0.0f;
593 const float distance =
594 73007976 KernelProfileDistance(data, c, dx - shift * c.shiftX, dy - shift * c.shiftY);
595
2/2
✓ Branch 0 taken 20080984 times.
✓ Branch 1 taken 52926992 times.
 73007976 if (distance >= k.reachRadii * radius) {
596 20080984 return height;
597 }
598 52926992 const float radii = distance / radius;
599 52926992 float replaced = 0.0f;
600
2/2
✓ Branch 0 taken 27652074 times.
✓ Branch 1 taken 25274918 times.
 52926992 if (c.hasDatum == 1u) {
601 float weight =
602 27652074 k.interiorReplacement * KernelFalloff(radii, k.blendStart, k.blendEnd, k.blendShape);
603
4/4
✓ Branch 0 taken 4919239 times.
✓ Branch 1 taken 22732835 times.
✓ Branch 2 taken 4357894 times.
✓ Branch 3 taken 561345 times.
 27652074 if (k.mantling > 0.0f && radii > 1.0f) {
604 4357894 weight = max(weight, k.mantling * KernelEjectaShare(radii, k.ejectaReachRadii));
605 }
606 27652074 replaced = weight * (c.datumA + c.datumB * dx + c.datumC * dy - height);
607 }
608 52926992 float bare = KernelRoundedHeight(k, data, distance, diameter, c);
609
6/6
✓ Branch 0 taken 40305377 times.
✓ Branch 1 taken 12621615 times.
✓ Branch 2 taken 38992703 times.
✓ Branch 3 taken 1312674 times.
✓ Branch 4 taken 91742 times.
✓ Branch 5 taken 38900961 times.
 52926992 if (bare > 0.0f && (c.upslopeX != 0.0f || c.upslopeY != 0.0f)) {
610 1404416 bare *= KernelSlopeRaisedShare(c, dx, dy);
611 }
612 const float shaped =
613 52926992 bare + KernelNoiseHeight(k, c, dx, dy, radii) + KernelTextureHeight(k, c, dx, dy, radii);
614 52926992 return height + c.resolution * shaped + c.resolution * replaced;
615 }
616
617 /// What fitting crater datums takes: the craters' shapes, the fit footprint and the lattice the
618 /// ground is sampled on (`FitCraterDatum`), and the size bands that index the craters.
619 struct CraterDatumKernel {
620 CraterShapeKernel shape;
621 float fitRadii = 1.0f; ///< The fit footprint, in crater radii.
622 uint32_t latticeSide = 7u; ///< Lattice points per side.
623 uint32_t bandCount = 0u; ///< Size bands of the crater index.
624 uint32_t unused = 0u;
625 };
626
627 /// The relief (m) crater `c` adds at `dx` east and `dy` north of its centre without its datum,
628 /// as `CraterReliefIndex` sums it: its resolution times the bare profile, within the profile's
629 /// reach. `bound` (m) is at least that reach stretched by the crater's outline.
630 23008264 inline float CraterRelief(KERNEL_IN(CraterShapeKernel) k, CraterData data,
631 KERNEL_IN(CraterRecord) c, float bound, float dx, float dy) {
632 23008264 const float squared = dx * dx + dy * dy;
633
2/2
✓ Branch 0 taken 14177233 times.
✓ Branch 1 taken 8831031 times.
 23008264 if (!(squared < bound * bound)) {
634 14177233 return 0.0f;
635 }
636 const float shift =
637
4/4
✓ Branch 0 taken 8505973 times.
✓ Branch 1 taken 325058 times.
✓ Branch 2 taken 37 times.
✓ Branch 3 taken 8505936 times.
 8831031 c.shiftX != 0.0f || c.shiftY != 0.0f ? KernelSlopeShiftWeight(c, dx, dy) : 0.0f;
638 const float distance =
639 8831031 KernelProfileDistance(data, c, dx - shift * c.shiftX, dy - shift * c.shiftY);
640
2/2
✓ Branch 0 taken 4746092 times.
✓ Branch 1 taken 4084939 times.
 8831031 if (!(distance < k.profileReachRadii * 0.5f * c.diameter)) {
641 4746092 return 0.0f;
642 }
643 4084939 float height = KernelProfileHeight(k, data, distance, c.diameter, c);
644
5/6
✓ Branch 0 taken 3638676 times.
✓ Branch 1 taken 446263 times.
✓ Branch 2 taken 3525309 times.
✓ Branch 3 taken 113367 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 3525309 times.
 4084939 if (height > 0.0f && (c.upslopeX != 0.0f || c.upslopeY != 0.0f)) {
645 113367 height *= KernelSlopeRaisedShare(c, dx, dy);
646 }
647 4084939 return c.resolution * height;
648 }
649
650 /// The offset (m) of lattice point `index` of `side` across a footprint of `fitRadius`.
651 2664906 inline float KernelDatumOffset(uint32_t index, uint32_t side, float fitRadius) {
652 2664906 return fitRadius * ((2.0f * float(index) + 1.0f) / float(side) - 1.0f);
653 }
654
655 /// Relative float error allowed for, in cell lookups, at a position's magnitude.
656 static const float kDatumCellSlack = 4.8e-7f;
657
658 /// The relief (m) of the craters before crater `k` of `records` at lattice point `point` of `k`'s
659 /// datum footprint, numbered row by row; 0 outside the footprint.
660 ///
661 /// Records hold positions as the floats nearest them and `lows` (x, y pairs) what those drop.
662 /// Craters are found through size bands: per band, `bandFloats` holds its cell size, which is at
663 /// least the reach of its craters, and its first cell's column and row, counted from the
664 /// surface's origin; `bandCells` its cells across and down and its first cell; per cell, `cells`
665 /// holds its first entry in `entries` and its entry count, entries being record indices in
666 /// ascending order. Cells and order are absolute, so the sum depends only on the craters reaching
667 /// the point.
668 1332453 inline float KernelDatumRelief(KERNEL_IN(CraterDatumKernel) kernel, CraterData data,
669 CraterRecords records, KernelFloats lows, KernelFloats bandFloats,
670 KernelUints bandCells, KernelUints cells, KernelUints entries,
671 uint32_t k, uint32_t point) {
672 1332453 const uint32_t side = kernel.latticeSide;
673 1332453 const float fitRadius = kernel.fitRadii * 0.5f * records[k].diameter;
674 1332453 const float dx = KernelDatumOffset(point % side, side, fitRadius);
675 1332453 const float dy = KernelDatumOffset(point / side, side, fitRadius);
676
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1332453 times.
 1332453 if (dx * dx + dy * dy >= fitRadius * fitRadius) {
677 ✗ return 0.0f;
678 }
679 1332453 const float lowX = lows[2u * k] + dx;
680 1332453 const float lowY = lows[2u * k + 1u] + dy;
681 1332453 const float px = records[k].x + lowX;
682 1332453 const float py = records[k].y + lowY;
683 1332453 const float slackX = kDatumCellSlack * abs(px);
684 1332453 const float slackY = kDatumCellSlack * abs(py);
685 1332453 float relief = 0.0f;
686
2/2
✓ Branch 0 taken 7315183 times.
✓ Branch 1 taken 1332453 times.
 8647636 for (uint32_t b = 0u; b < kernel.bandCount; ++b) {
687 7315183 const float size = bandFloats[3u * b];
688 7315183 const int32_t across = int32_t(bandCells[3u * b]);
689 7315183 const int32_t down = int32_t(bandCells[3u * b + 1u]);
690 const int32_t i0 =
691 7315183 max(int32_t(floor((px - size - slackX) / size) - bandFloats[3u * b + 1u]), 0);
692 const int32_t i1 =
693 7315183 min(int32_t(floor((px + size + slackX) / size) - bandFloats[3u * b + 1u]), across - 1);
694 const int32_t j0 =
695 7315183 max(int32_t(floor((py - size - slackY) / size) - bandFloats[3u * b + 2u]), 0);
696 const int32_t j1 =
697 7315183 min(int32_t(floor((py + size + slackY) / size) - bandFloats[3u * b + 2u]), down - 1);
698
2/2
✓ Branch 0 taken 20717011 times.
✓ Branch 1 taken 7315183 times.
 28032194 for (int32_t j = j0; j <= j1; ++j) {
699
2/2
✓ Branch 0 taken 60117104 times.
✓ Branch 1 taken 20717011 times.
 80834115 for (int32_t i = i0; i <= i1; ++i) {
700 60117104 const uint32_t cell = bandCells[3u * b + 2u] + uint32_t(j * across + i);
701 60117104 const uint32_t first = cells[2u * cell];
702 60117104 const uint32_t end = first + cells[2u * cell + 1u];
703
2/2
✓ Branch 0 taken 30681252 times.
✓ Branch 1 taken 52444116 times.
 83125368 for (uint32_t e = first; e < end; ++e) {
704 30681252 const uint32_t n = entries[e];
705
2/2
✓ Branch 0 taken 7672988 times.
✓ Branch 1 taken 23008264 times.
 30681252 if (n >= k) {
706 7672988 break;
707 }
708 23008264 relief +=
709 23008264 CraterRelief(kernel.shape, data, records[n], size,
710 23008264 (records[k].x - records[n].x) + (lowX - lows[2u * n]),
711 23008264 (records[k].y - records[n].y) + (lowY - lows[2u * n + 1u]));
712 }
713 }
714 }
715 }
716 1332453 return relief;
717 }
718
719 #ifndef __SLANG__
720 static_assert(sizeof(CraterProfileKernel) == 76,
721 "CraterProfileKernel must match its std430 layout");
722 static_assert(sizeof(CraterShapeKernel) == 180, "CraterShapeKernel must match its std430 layout");
723 static_assert(sizeof(CraterRecord) == 116, "CraterRecord must match its std430 layout");
724
725 /// The compiled forms of the crater kernel.
726 template <> KernelProgram CompiledKernel<CraterShapeKernel>();
727 static_assert(sizeof(CraterDatumKernel) == 196, "CraterDatumKernel must match its std430 layout");
728
729 /// The compiled forms of the crater datum kernel.
730 template <> KernelProgram CompiledKernel<CraterDatumKernel>();
731 #endif
732
733 } // namespace terraingen
734