TerrainGen coverage

terraingen/
Coverage:
low: ≥ 0%
medium: ≥ 75.0%
high: ≥ 90.0%
Lines:
213 / 215
99.1%
Branches:
241 / 264
91.3%
Functions:
25 / 25
100.0%
layers/craters/generation/population.cpp
Line Branch Exec Source
1 #include "terraingen/layers/craters/generation/population.h"
2
3 #include "terraingen/layers/craters/generation/spacing.h"
4 #include "terraingen/random/philox.h"
5 #include "terraingen/random/samplers.h"
6 #include "terraingen/util/checks.h"
7 #include "terraingen/util/log.h"
8
9 #include <algorithm>
10 #include <cmath>
11 #include <memory>
12 #include <stdexcept>
13 #include <string>
14 #include <utility>
15
16 namespace terraingen {
17
18 namespace {
19
20 constexpr std::uint32_t kCraterStream = 0x43524154u;
21 constexpr std::uint32_t kCraterAgeStream = 0x43524147u;
22 constexpr std::uint32_t kCraterIdStream = 0x43524944u;
23 constexpr std::uint32_t kCraterKeepStream = 0x43524B50u;
24
25 3092554 std::uint64_t CraterId(std::uint32_t seed, const CraterBand& band, const CellIndex& cell,
26 std::uint32_t index) {
27 const U32x4 block =
28 3092554 Philox4x32({static_cast<std::uint32_t>(cell.i), static_cast<std::uint32_t>(cell.j),
29 3092554 static_cast<std::uint32_t>(band.exponent), index},
30 {seed, kCraterIdStream});
31 3092554 return (static_cast<std::uint64_t>(block.x) << 32) | block.y;
32 }
33
34 364680 void CheckOptions(const CraterGenerationOptions& options) {
35 364680 ValidateCraterSpacing(options.spacing);
36
1/1
✓ Branch 2 taken 364677 times.
 364679 RequirePositive(options.dMin, "dMin");
37
2/2
✓ Branch 0 taken 1 time.
✓ Branch 1 taken 364676 times.
 364677 if (!(options.dMax > options.dMin)) {
38
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("dMax must exceed dMin");
39 }
40
1/1
✓ Branch 2 taken 364675 times.
 364676 RequireNonNegative(options.marginRadii, "marginRadii");
41
1/1
✓ Branch 2 taken 364674 times.
 364675 RequirePositive(options.cellSizeFactor, "cellSizeFactor");
42
2/2
✓ Branch 0 taken 1 time.
✓ Branch 1 taken 364673 times.
 364674 if (!(options.maxExpectedCraters >= 0.0)) {
43
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("maxExpectedCraters must not be negative");
44 }
45
1/1
✓ Branch 2 taken 364672 times.
 364673 RequireNonNegative(options.minResolvedDiameter, "minResolvedDiameter");
46
2/2
✓ Branch 0 taken 1 time.
✓ Branch 1 taken 364671 times.
 364672 if (!(options.maxBandDiameter > 0.0)) {
47
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("maxBandDiameter must be positive");
48 }
49 364671 }
50
51 /// One band of a generation: the cells drawn and the craters they are expected to hold.
52 struct BandPlan {
53 CraterBand band;
54 std::vector<CellIndex> cells;
55 double expected = 0.0;
56 };
57
58 /// What a generation keeps and drops, for its log line.
59 struct GenerationCounts {
60 std::size_t bands = 0;
61 double expected = 0.0;
62 std::size_t removed = 0; ///< By spacing.
63 std::size_t capped = 0; ///< Over `maxCraters`.
64 };
65
66 /// One call of `GenerateCraters`, in steps.
67 class Generation {
68 public:
69 364671 Generation(const SurfaceCells& cells, const Bounds2d& region,
70 const SizeFrequencyDistribution& density, const CraterAgeModel& ages,
71 std::uint32_t seed, const CraterGenerationOptions& options)
72 364671 : cells_(cells), region_(region), density_(density), ages_(ages), seed_(seed),
73 364671 options_(options), surface_(options.surface) {
74
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 364670 times.
 364671 if (density.Kind() != SfdKind::Density) {
75
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("Craters are generated from a Density distribution");
76 }
77 364670 dMax_ = ResolveDMax(options.dMax, density);
78
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 364669 times.
 364670 if (std::isinf(dMax_)) {
79
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("An unbounded density needs a finite dMax");
80 }
81
2/2
✓ Branch 0 taken 1 time.
✓ Branch 1 taken 364668 times.
 364669 if (!(dMax_ > options.dMin)) {
82
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("dMax must exceed dMin");
83 }
84 364668 }
85
86 364668 CraterList Run(GenerationCounts& counts) const {
87
1/1
✓ Branch 1 taken 364666 times.
 364668 const std::vector<BandPlan> plan = PlanBands();
88 364666 counts.bands = plan.size();
89
2/2
✓ Branch 4 taken 1142969 times.
✓ Branch 5 taken 364666 times.
 1507635 for (const BandPlan& entry : plan) {
90 1142969 counts.expected += entry.expected;
91 }
92
2/2
✓ Branch 0 taken 1 time.
✓ Branch 1 taken 364665 times.
 364666 if (counts.expected > options_.maxExpectedCraters) {
93
2/2
✓ Branch 1 taken 1 time.
✓ Branch 4 taken 1 time.
 2 throw std::invalid_argument("Expected crater count " + std::to_string(counts.expected) +
94
2/2
✓ Branch 2 taken 1 time.
✓ Branch 5 taken 1 time.
 3 " exceeds maxExpectedCraters");
95 }
96
1/1
✓ Branch 1 taken 364665 times.
 364665 CraterList craters = DrawPlannedCells(plan, counts.expected);
97 364665 const std::size_t generated = craters.Size();
98
1/1
✓ Branch 2 taken 364662 times.
 364668 AppendManualCraters(options_.manual, craters);
99
2/2
✓ Branch 1 taken 37 times.
✓ Branch 2 taken 364625 times.
 364662 if (options_.spacing.Enabled()) {
100
1/1
✓ Branch 1 taken 37 times.
 37 craters = RemoveSpacingConflicts(craters, counts.removed);
101
2/2
✓ Branch 1 taken 2950 times.
✓ Branch 2 taken 361675 times.
 364625 } else if (options_.manual.Size() > 0) {
102
1/1
✓ Branch 1 taken 2950 times.
 2950 craters = CropToRegion(craters, generated);
103 }
104
1/1
✓ Branch 1 taken 364662 times.
 729324 return KeepLargest(craters, counts.capped);
105 364669 }
106
107 364662 void Log(std::size_t kept, const GenerationCounts& counts, const Stopwatch& stopwatch) const {
108 364662 const bool spaced = options_.spacing.Enabled();
109 364662 const std::size_t manual = options_.manual.Size();
110
1/1
✓ Branch 1 taken 364662 times.
 364662 LogDebug(options_.logName, "Generated ", kept, " ", options_.logName, " (", options_.dMin,
111 364662 "-", dMax_, " m, ", counts.bands, " bands, ",
112 364662 static_cast<long long>(std::llround(counts.expected)), " expected",
113
12/18
✓ Branch 0 taken 37 times.
✓ Branch 1 taken 364625 times.
✓ Branch 3 taken 37 times.
✓ Branch 6 taken 37 times.
✓ Branch 9 taken 37 times.
✓ Branch 12 taken 364625 times.
✓ Branch 15 taken 364625 times.
✓ Branch 16 taken 37 times.
✓ Branch 17 taken 37 times.
✓ Branch 18 taken 364625 times.
✓ Branch 20 taken 37 times.
✓ Branch 21 taken 364625 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
 1458574 spaced ? ", " + std::to_string(counts.removed) + " removed by spacing" : "",
114
12/18
✓ Branch 0 taken 2965 times.
✓ Branch 1 taken 361697 times.
✓ Branch 3 taken 2965 times.
✓ Branch 6 taken 2965 times.
✓ Branch 9 taken 2965 times.
✓ Branch 12 taken 361697 times.
✓ Branch 15 taken 361697 times.
✓ Branch 16 taken 2965 times.
✓ Branch 17 taken 2965 times.
✓ Branch 18 taken 361697 times.
✓ Branch 20 taken 2965 times.
✓ Branch 21 taken 361697 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
 1452718 manual > 0 ? ", " + std::to_string(manual) + " manual supplied" : "",
115 364662 counts.capped > 0
116
12/18
✓ Branch 0 taken 6 times.
✓ Branch 1 taken 364656 times.
✓ Branch 3 taken 6 times.
✓ Branch 6 taken 6 times.
✓ Branch 9 taken 6 times.
✓ Branch 12 taken 364656 times.
✓ Branch 15 taken 364656 times.
✓ Branch 16 taken 6 times.
✓ Branch 17 taken 6 times.
✓ Branch 18 taken 364656 times.
✓ Branch 20 taken 6 times.
✓ Branch 21 taken 364656 times.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
 1458636 ? ", " + std::to_string(counts.capped) + " over maxCraters dropped"
117 : "",
118
2/2
✓ Branch 1 taken 364662 times.
✓ Branch 4 taken 364662 times.
 729324 ") in ", FormatDuration(stopwatch.Milliseconds()));
119 364662 }
120
121 private:
122 /// How far beyond the region a crater `diameter` across is drawn: its margin, and with spacing
123 /// as far again as the craters it can conflict with reach.
124 4075650 double CandidateMargin(double diameter) const {
125 4075650 const CraterSpacing& spacing = options_.spacing;
126
2/2
✓ Branch 1 taken 4070811 times.
✓ Branch 2 taken 4839 times.
 4075650 if (!spacing.Enabled()) {
127 4070811 return 0.5 * options_.marginRadii * diameter;
128 }
129 4839 const double neighbour = std::fmin(spacing.diameterRatio * diameter, dMax_);
130 4839 return SpacingReach(diameter, dMax_, spacing) + 0.5 * options_.marginRadii * neighbour;
131 }
132
133 /// Whether `crater` is centred in the region or within its margin of it.
134 679625 bool InRegion(const Crater& crater, double margin) const {
135
1/1
✓ Branch 1 taken 679625 times.
 679625 return surface_.Contains(surface_.Expand(region_, margin), {crater.x, crater.y});
136 }
137
138 364668 std::vector<BandPlan> PlanBands() const {
139 364668 const CraterSpacing& spacing = options_.spacing;
140 const double skipBelow =
141
2/2
✓ Branch 1 taken 37 times.
✓ Branch 2 taken 364631 times.
 364668 options_.minResolvedDiameter / (spacing.Enabled() ? spacing.diameterRatio : 1.0);
142 364668 std::vector<BandPlan> plan;
143
3/3
✓ Branch 1 taken 364668 times.
✓ Branch 8 taken 1190665 times.
✓ Branch 9 taken 364666 times.
 1555331 for (const CraterBand& band : CraterBands(options_.dMin, dMax_, options_.cellSizeFactor)) {
144
4/4
✓ Branch 1 taken 1190497 times.
✓ Branch 2 taken 168 times.
✓ Branch 3 taken 1190459 times.
✓ Branch 4 taken 38 times.
 2381124 const bool tooLarge = !spacing.Enabled() && options_.maxCraters == 0 &&
145
2/2
✓ Branch 0 taken 28 times.
✓ Branch 1 taken 1190431 times.
 1190459 band.low > options_.maxBandDiameter;
146
4/4
✓ Branch 0 taken 1143003 times.
✓ Branch 1 taken 47662 times.
✓ Branch 2 taken 28 times.
✓ Branch 3 taken 1142975 times.
 1190665 if (band.high <= skipBelow || tooLarge) {
147 47690 continue;
148 }
149
1/1
✓ Branch 1 taken 1142975 times.
 1142975 const double margin = CandidateMargin(band.high);
150 1142975 BandPlan entry{band,
151
2/2
✓ Branch 1 taken 1142975 times.
✓ Branch 4 taken 1142974 times.
 1142975 cells_.CellsCovering(surface_.Expand(region_, margin), band.cellSize)};
152
1/1
✓ Branch 1 taken 1142973 times.
 1142974 const double countPerArea = density_.ExpectedCount(1.0, band.low, band.high);
153
2/2
✓ Branch 4 taken 3080380 times.
✓ Branch 5 taken 1142973 times.
 4223353 for (const CellIndex& cell : entry.cells) {
154
1/1
✓ Branch 1 taken 3080380 times.
 3080380 entry.expected += countPerArea * cells_.CellArea(cell, band.cellSize);
155 }
156
1/1
✓ Branch 1 taken 1142973 times.
 1142973 LogDebug(options_.logName, "Band [", band.low, ", ", band.high,
157 1142973 ") m: ", entry.cells.size(), " cells of ", band.cellSize, " m, ",
158 1142973 entry.expected, " ", options_.logName, " expected");
159
1/1
✓ Branch 2 taken 1142973 times.
 1142973 plan.push_back(std::move(entry));
160 1507642 }
161 364666 return plan;
162 2 }
163
164 364665 CraterList DrawPlannedCells(const std::vector<BandPlan>& plan, double expected) const {
165 364665 CraterList craters;
166
1/1
✓ Branch 1 taken 364665 times.
 364665 craters.Reserve(static_cast<std::size_t>(expected));
167
2/2
✓ Branch 5 taken 1142958 times.
✓ Branch 6 taken 364665 times.
 1507623 for (const BandPlan& entry : plan) {
168
1/1
✓ Branch 1 taken 1142958 times.
 1142958 const BandRegion region = BandRegionOf(entry.band);
169
2/2
✓ Branch 5 taken 2866993 times.
✓ Branch 6 taken 1142958 times.
 4009951 for (const CellIndex& cell : entry.cells) {
170
1/1
✓ Branch 1 taken 2866993 times.
 2866993 const std::shared_ptr<const CraterList> drawn = DrawCell(entry.band, cell);
171
2/2
✓ Branch 2 taken 112124646 times.
✓ Branch 3 taken 2866993 times.
 114991639 for (std::size_t n = 0; n < drawn->Size(); ++n) {
172
1/1
✓ Branch 2 taken 112124646 times.
 112124646 const Crater crater = drawn->At(n);
173
8/8
✓ Branch 1 taken 112124646 times.
✓ Branch 3 taken 8932640 times.
✓ Branch 4 taken 103192006 times.
✓ Branch 6 taken 8932640 times.
✓ Branch 8 taken 8897595 times.
✓ Branch 9 taken 35045 times.
✓ Branch 10 taken 8897595 times.
✓ Branch 11 taken 103227051 times.
 112124646 if (InBandRegion(region, crater) && Kept(crater)) {
174
1/1
✓ Branch 1 taken 8897595 times.
 8897595 craters.Append(crater);
175 }
176 }
177 2866993 }
178 }
179 364665 return craters;
180 }
181
182 /// The region expanded by the candidate margins of a band's smallest and largest craters,
183 /// which bound every margin in between since `Expand` and `CandidateMargin` only grow.
184 struct BandRegion {
185 CraterBand band;
186 Bounds2d inner;
187 Bounds2d outer;
188 };
189
190 1142958 BandRegion BandRegionOf(const CraterBand& band) const {
191 1142958 return {band, surface_.Expand(region_, CandidateMargin(band.low)),
192 2285916 surface_.Expand(region_, CandidateMargin(band.high))};
193 }
194
195 /// `InRegion` with the crater's candidate margin, expanding the region only for craters of the
196 /// band between its inner and outer regions.
197 112124646 bool InBandRegion(const BandRegion& region, const Crater& crater) const {
198
3/4
✓ Branch 0 taken 112124646 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 40 times.
✓ Branch 3 taken 112124606 times.
 112124646 if (!(crater.diameter >= region.band.low && crater.diameter <= region.band.high)) {
199
2/2
✓ Branch 1 taken 40 times.
✓ Branch 4 taken 40 times.
 40 return InRegion(crater, CandidateMargin(crater.diameter));
200 }
201 112124606 const Vec2d centre{crater.x, crater.y};
202
2/2
✓ Branch 1 taken 102751342 times.
✓ Branch 2 taken 9373264 times.
 112124606 if (!surface_.Contains(region.outer, centre)) {
203 102751342 return false;
204 }
205
2/2
✓ Branch 1 taken 646719 times.
✓ Branch 2 taken 8726545 times.
 10019983 return surface_.Contains(region.inner, centre) ||
206
4/4
✓ Branch 1 taken 646719 times.
✓ Branch 4 taken 646719 times.
✓ Branch 6 taken 206075 times.
✓ Branch 7 taken 440644 times.
 10019983 InRegion(crater, CandidateMargin(crater.diameter));
207 }
208
209 /// Whether a generated crater stays: below the catalog's completeness limit at its centre, and
210 /// within the keep share there by a uniform keyed to its ID.
211 8932640 bool Kept(const Crater& crater) const {
212
4/4
✓ Branch 0 taken 37031 times.
✓ Branch 1 taken 8895609 times.
✓ Branch 2 taken 1000 times.
✓ Branch 3 taken 8931640 times.
 8969671 if (options_.completeness &&
213
3/3
✓ Branch 3 taken 37031 times.
✓ Branch 5 taken 1000 times.
✓ Branch 6 taken 36031 times.
 8969671 crater.diameter >= options_.completeness->ValueAt(crater.x, crater.y)) {
214 1000 return false;
215 }
216
2/2
✓ Branch 1 taken 8384135 times.
✓ Branch 2 taken 547505 times.
 8931640 if (!options_.keepShare) {
217 8384135 return true;
218 }
219 547505 const IdWords id = SplitId(crater.id);
220
1/1
✓ Branch 1 taken 547505 times.
 547505 RandomStream random(seed_, kCraterKeepStream, id.low, id.high);
221
2/2
✓ Branch 1 taken 547505 times.
✓ Branch 4 taken 547505 times.
 547505 return random.NextDouble() < options_.keepShare(crater);
222 }
223
224 2866993 std::shared_ptr<const CraterList> DrawCell(const CraterBand& band,
225 const CellIndex& cell) const {
226
2/2
✓ Branch 1 taken 2764900 times.
✓ Branch 2 taken 102093 times.
 2866993 if (options_.cellCraters) {
227 2764900 return options_.cellCraters(band, cell);
228 }
229 return std::make_shared<const CraterList>(
230
2/2
✓ Branch 1 taken 102093 times.
✓ Branch 4 taken 102093 times.
 102093 GenerateCellCraters(cells_, band, cell, density_, ages_, seed_));
231 }
232
233 /// The craters within their margin of the region, the first `generated` of which already are.
234 2950 CraterList CropToRegion(const CraterList& craters, std::size_t generated) const {
235 2950 return FilterCraters(craters, [&](std::size_t n, const Crater& crater) {
236
4/4
✓ Branch 0 taken 2965 times.
✓ Branch 1 taken 7280332 times.
✓ Branch 3 taken 897 times.
✓ Branch 4 taken 2068 times.
 7283297 return n < generated || InRegion(crater, 0.5 * options_.marginRadii * crater.diameter);
237
1/1
✓ Branch 1 taken 2950 times.
 5900 });
238 }
239
240 /// The craters within their margin of the region that spacing keeps, counting those removed.
241 37 CraterList RemoveSpacingConflicts(const CraterList& craters, std::size_t& removed) const {
242
1/1
✓ Branch 1 taken 37 times.
 37 const std::vector<bool> conflicts = SpacingConflicts(craters, options_.spacing, surface_);
243 ✗ return FilterCraters(craters, [&](std::size_t n, const Crater& crater) {
244
2/2
✓ Branch 1 taken 4764 times.
✓ Branch 2 taken 25137 times.
 29901 if (!InRegion(crater, 0.5 * options_.marginRadii * crater.diameter)) {
245 4764 return false;
246 }
247
2/2
✓ Branch 1 taken 8848 times.
✓ Branch 2 taken 16289 times.
 25137 if (conflicts[n]) {
248 8848 ++removed;
249 8848 return false;
250 }
251 16289 return true;
252
1/1
✓ Branch 1 taken 37 times.
 74 });
253 37 }
254
255 /// The manual craters and the `maxCraters` largest generated ones (ties to the smaller ID),
256 /// counting the generated craters dropped.
257 364662 CraterList KeepLargest(const CraterList& craters, std::size_t& capped) const {
258 364662 const std::size_t limit = options_.maxCraters;
259 364662 std::vector<std::size_t> generated;
260
2/2
✓ Branch 1 taken 8884876 times.
✓ Branch 2 taken 364662 times.
 9249538 for (std::size_t n = 0; n < craters.Size(); ++n) {
261
2/2
✓ Branch 2 taken 8883962 times.
✓ Branch 3 taken 914 times.
 8884876 if (craters.Manual()[n] == 0) {
262
1/1
✓ Branch 1 taken 8883962 times.
 8883962 generated.push_back(n);
263 }
264 }
265
6/6
✓ Branch 0 taken 7 times.
✓ Branch 1 taken 364655 times.
✓ Branch 3 taken 1 time.
✓ Branch 4 taken 6 times.
✓ Branch 5 taken 364656 times.
✓ Branch 6 taken 6 times.
 364662 if (limit == 0 || generated.size() <= limit) {
266
1/1
✓ Branch 1 taken 364656 times.
 364656 return craters;
267 }
268 6 const auto& diameters = craters.Diameters();
269 6 const auto& ids = craters.Ids();
270
1/1
✓ Branch 5 taken 6 times.
 6 std::nth_element(generated.begin(), generated.begin() + static_cast<std::ptrdiff_t>(limit),
271 32818 generated.end(), [&](std::size_t a, std::size_t b) {
272
1/2
✓ Branch 2 taken 32818 times.
✗ Branch 3 not taken.
 32818 if (diameters[a] != diameters[b]) {
273 32818 return diameters[a] > diameters[b];
274 }
275 ✗ return ids[a] < ids[b];
276 });
277
1/1
✓ Branch 2 taken 6 times.
 6 std::vector<bool> keep(craters.Size(), false);
278
2/2
✓ Branch 0 taken 660 times.
✓ Branch 1 taken 6 times.
 666 for (std::size_t k = 0; k < limit; ++k) {
279
1/1
✓ Branch 2 taken 660 times.
 660 keep[generated[k]] = true;
280 }
281 6 capped = generated.size() - limit;
282 return FilterCraters(
283
6/6
✓ Branch 1 taken 6 times.
✓ Branch 4 taken 14546 times.
✓ Branch 7 taken 13886 times.
✓ Branch 8 taken 660 times.
✓ Branch 9 taken 1 time.
✓ Branch 10 taken 13885 times.
 14552 craters, [&](std::size_t n, const Crater& crater) { return keep[n] || crater.manual; });
284 364662 }
285
286 const SurfaceCells& cells_;
287 const Bounds2d& region_;
288 const SizeFrequencyDistribution& density_;
289 const CraterAgeModel& ages_;
290 std::uint32_t seed_;
291 const CraterGenerationOptions& options_;
292 const Surface& surface_;
293 double dMax_ = 0.0;
294 };
295
296 } // namespace
297
298 364670 std::vector<CraterBand> CraterBands(double dMin, double dMax, double cellSizeFactor) {
299 364670 std::vector<CraterBand> bands;
300
2/2
✓ Branch 0 taken 1190674 times.
✓ Branch 1 taken 364670 times.
 1555344 for (int k = std::ilogb(dMin); std::ldexp(1.0, k) < dMax; ++k) {
301 1190674 const double low = std::fmax(dMin, std::ldexp(1.0, k));
302 1190674 const double high = std::fmin(dMax, std::ldexp(1.0, k + 1));
303
1/1
✓ Branch 1 taken 1190674 times.
 1190674 bands.push_back({k, low, high, cellSizeFactor * std::ldexp(1.0, k + 1)});
304 }
305 364670 return bands;
306 }
307
308 208356 CraterList GenerateCellCraters(const SurfaceCells& cells, const CraterBand& band,
309 const CellIndex& cell, const SizeFrequencyDistribution& density,
310 const CraterAgeModel& ages, std::uint32_t seed) {
311 208356 const auto k = static_cast<std::uint32_t>(band.exponent);
312 208356 const auto i = static_cast<std::uint32_t>(cell.i);
313 208356 const auto j = static_cast<std::uint32_t>(cell.j);
314
1/1
✓ Branch 1 taken 208356 times.
 208356 RandomStream geometry(seed, kCraterStream, i, j, k);
315
1/1
✓ Branch 1 taken 208356 times.
 208356 RandomStream ageStream(seed, kCraterAgeStream, i, j, k);
316
1/1
✓ Branch 1 taken 208356 times.
 208356 const double countPerArea = density.ExpectedCount(1.0, band.low, band.high);
317
1/1
✓ Branch 1 taken 208356 times.
 208356 const SurfaceCells::Sampling sampling = cells.CellSampling(cell, band.cellSize);
318
1/1
✓ Branch 1 taken 208356 times.
 208356 const std::uint64_t count = SamplePoisson(geometry, countPerArea * sampling.area);
319 208356 CraterList craters;
320
1/1
✓ Branch 1 taken 208356 times.
 208356 craters.Reserve(static_cast<std::size_t>(count));
321
2/2
✓ Branch 0 taken 3179938 times.
✓ Branch 1 taken 208356 times.
 3388294 for (std::uint64_t index = 0; index < count; ++index) {
322
1/1
✓ Branch 1 taken 3179938 times.
 3179938 const double diameter = density.SampleDiameter(geometry, band.low, band.high);
323
1/1
✓ Branch 1 taken 3179938 times.
 3179938 const double u = geometry.NextDouble();
324
1/1
✓ Branch 1 taken 3179938 times.
 3179938 const double v = geometry.NextDouble();
325
1/1
✓ Branch 1 taken 3179938 times.
 3179938 const Vec2d centre = cells.SamplePoint(u, v, cell, band.cellSize);
326
4/4
✓ Branch 0 taken 160202 times.
✓ Branch 1 taken 3019736 times.
✓ Branch 2 taken 87384 times.
✓ Branch 3 taken 3092554 times.
 3340140 if (sampling.bound > 0.0 &&
327
4/4
✓ Branch 1 taken 160202 times.
✓ Branch 4 taken 160202 times.
✓ Branch 6 taken 87384 times.
✓ Branch 7 taken 72818 times.
 160202 geometry.NextDouble() * sampling.bound >= cells.Weight(centre)) {
328 87384 continue;
329 }
330
1/1
✓ Branch 1 taken 3092554 times.
 3092554 const double ageGyr = ages.SampleAgeGyr(ageStream, diameter);
331
1/1
✓ Branch 1 taken 3092554 times.
 3092554 craters.Append({centre.x, centre.y, diameter, ageGyr,
332 3092554 CraterId(seed, band, cell, static_cast<std::uint32_t>(index))});
333 }
334 416712 return craters;
335 }
336
337 583992 double ResolveDMax(double dMax, const SizeFrequencyDistribution& density) {
338
2/2
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 583987 times.
 583992 return std::isinf(dMax) ? density.ValidRange().max : dMax;
339 }
340
341 364680 CraterList GenerateCraters(const SurfaceCells& cells, const Bounds2d& region,
342 const SizeFrequencyDistribution& density, const CraterAgeModel& ages,
343 std::uint32_t seed, const CraterGenerationOptions& options) {
344
1/1
✓ Branch 1 taken 364671 times.
 364680 CheckOptions(options);
345 364671 const Stopwatch stopwatch;
346
1/1
✓ Branch 1 taken 364668 times.
 364671 const Generation generation(cells, region, density, ages, seed, options);
347 364668 GenerationCounts counts;
348
1/1
✓ Branch 1 taken 364662 times.
 364668 CraterList craters = generation.Run(counts);
349
1/1
✓ Branch 2 taken 364662 times.
 364662 generation.Log(craters.Size(), counts, stopwatch);
350 729324 return craters;
351 }
352
353 6 double MinDiameterForCount(const SizeFrequencyDistribution& density, double area, double count,
354 double dMin, double dMax) {
355
6/6
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 1 time.
✓ Branch 4 taken 1 time.
✓ Branch 5 taken 4 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 4 times.
 6 if (!(IsPositive(count) && IsPositive(area))) {
356
1/1
✓ Branch 2 taken 2 times.
 2 throw std::invalid_argument("Count and area must be positive and finite");
357 }
358
6/8
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 3 times.
✓ Branch 4 taken 1 time.
✗ Branch 6 not taken.
✓ Branch 7 taken 3 times.
✓ Branch 8 taken 1 time.
✓ Branch 9 taken 3 times.
 4 if (!(IsPositive(dMin) && dMax > dMin && std::isfinite(dMax))) {
359
1/1
✓ Branch 2 taken 1 time.
 1 throw std::invalid_argument("Diameter interval must satisfy 0 < dMin < dMax");
360 }
361
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 2 times.
 3 if (density.ExpectedCount(area, dMin, dMax) <= count) {
362 1 return dMin;
363 }
364 2 double low = std::log(dMin);
365 2 double high = std::log(dMax);
366
3/4
✓ Branch 0 taken 88 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 86 times.
✓ Branch 3 taken 2 times.
 88 for (int i = 0; i < 200 && high - low > 1e-12; ++i) {
367 86 const double middle = 0.5 * (low + high);
368
2/2
✓ Branch 1 taken 40 times.
✓ Branch 2 taken 46 times.
 86 if (density.ExpectedCount(area, std::exp(middle), dMax) > count) {
369 40 low = middle;
370 } else {
371 46 high = middle;
372 }
373 }
374 2 return std::exp(high);
375 }
376
377 } // namespace terraingen
378