TerrainGen coverage

terraingen/
Coverage:
low: ≥ 0%
medium: ≥ 75.0%
high: ≥ 90.0%
Lines:
118 / 119
99.2%
Branches:
80 / 85
94.1%
Functions:
13 / 13
100.0%
layers/craters/shaping/datum.cpp
Line Branch Exec Source
1 #include "terraingen/layers/craters/shaping/datum.h"
2
3 #include "terraingen/noise/hash.h"
4
5 #include <algorithm>
6 #include <array>
7 #include <cmath>
8
9 namespace terraingen {
10
11 namespace {
12
13 using Matrix3 = std::array<std::array<double, 3>, 3>;
14
15 147127 double Determinant(const Matrix3& m) {
16 147127 return m[0][0] * (m[1][1] * m[2][2] - m[1][2] * m[2][1]) -
17 147127 m[0][1] * (m[1][0] * m[2][2] - m[1][2] * m[2][0]) +
18 147127 m[0][2] * (m[1][0] * m[2][1] - m[1][1] * m[2][0]);
19 }
20
21 /// Side of the buckets of a band: twice the largest reach of the band's craters.
22 2681 double BucketSize(int band, const CraterProfile& profile, double outlineScale) {
23 2681 return outlineScale * profile.ReachRadii() * std::ldexp(1.0, band + 1);
24 }
25
26 /// The ground under the lattice points inside a crater's fit footprint, fitted by least squares.
27 6 PlaneFit FitFootprint(const Crater& crater, double fitRadii, const Surface& surface,
28 const HeightSampler& ground) {
29 6 const double fitRadius = fitRadii * 0.5 * crater.diameter;
30 6 PlaneFit fit;
31
2/2
✓ Branch 0 taken 42 times.
✓ Branch 1 taken 6 times.
 48 for (int j = 0; j < kDatumLatticeSide; ++j) {
32 42 const double dy = DatumLatticeOffset(j, fitRadius);
33
2/2
✓ Branch 0 taken 294 times.
✓ Branch 1 taken 42 times.
 336 for (int i = 0; i < kDatumLatticeSide; ++i) {
34 294 const double dx = DatumLatticeOffset(i, fitRadius);
35
2/2
✓ Branch 0 taken 72 times.
✓ Branch 1 taken 222 times.
 294 if (dx * dx + dy * dy >= fitRadius * fitRadius) {
36 72 continue;
37 }
38
1/1
✓ Branch 1 taken 222 times.
 222 const Vec2d at = surface.Locate({crater.x, crater.y}, {dx, dy});
39
1/1
✓ Branch 1 taken 222 times.
 222 const float height = ground(at.x, at.y);
40
2/2
✓ Branch 1 taken 185 times.
✓ Branch 2 taken 37 times.
 222 if (!std::isnan(height)) {
41 185 fit.Add(dx, dy, height);
42 }
43 }
44 }
45 6 return fit;
46 }
47
48 } // namespace
49
50 1384339 void PlaneFit::Add(double dx, double dy, double h) {
51 1384339 n_ += 1.0;
52 1384339 sx_ += dx;
53 1384339 sy_ += dy;
54 1384339 sxx_ += dx * dx;
55 1384339 sxy_ += dx * dy;
56 1384339 syy_ += dy * dy;
57 1384339 sh_ += h;
58 1384339 sxh_ += dx * h;
59 1384339 syh_ += dy * h;
60 1384339 }
61
62 36889 Plane PlaneFit::Solve() const {
63 36889 const Matrix3 normal{{{n_, sx_, sy_}, {sx_, sxx_, sxy_}, {sy_, sxy_, syy_}}};
64 36889 const double det = Determinant(normal);
65
4/4
✓ Branch 0 taken 36871 times.
✓ Branch 1 taken 18 times.
✓ Branch 2 taken 125 times.
✓ Branch 3 taken 36746 times.
 36889 if (n_ < 3.0 || !(std::fabs(det) > 1e-9 * n_ * sxx_ * syy_)) {
66 143 return {sh_ / n_, 0.0, 0.0};
67 }
68 36746 const std::array<double, 3> rhs{sh_, sxh_, syh_};
69 36746 std::array<double, 3> solution{};
70
2/2
✓ Branch 0 taken 110238 times.
✓ Branch 1 taken 36746 times.
 146984 for (std::size_t column = 0; column < 3; ++column) {
71 110238 Matrix3 replaced = normal;
72
2/2
✓ Branch 0 taken 330714 times.
✓ Branch 1 taken 110238 times.
 440952 for (std::size_t row = 0; row < 3; ++row) {
73 330714 replaced[row][column] = rhs[row];
74 }
75 110238 solution[column] = Determinant(replaced) / det;
76 }
77 36746 return {solution[0], solution[1], solution[2]};
78 }
79
80 187 std::optional<Plane> FitCraterDatum(const Crater& crater, const CraterApplyOptions& options,
81 const HeightSampler& ground) {
82
2/2
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 185 times.
 187 if (options.datum == CraterDatum::None) {
83 2 return std::nullopt;
84 }
85
2/2
✓ Branch 0 taken 179 times.
✓ Branch 1 taken 6 times.
 185 if (options.datum == CraterDatum::CentreHeight) {
86
1/1
✓ Branch 1 taken 179 times.
 179 const float height = ground(crater.x, crater.y);
87
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 178 times.
 179 if (std::isnan(height)) {
88 1 return std::nullopt;
89 }
90 178 return Plane{static_cast<double>(height), 0.0, 0.0};
91 }
92
1/1
✓ Branch 1 taken 6 times.
 6 const PlaneFit fit = FitFootprint(crater, options.datumFitRadii, options.surface, ground);
93
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 5 times.
 6 if (fit.Empty()) {
94 1 return std::nullopt;
95 }
96
2/2
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 1 time.
 5 return options.datum == CraterDatum::Plane ? fit.Solve() : Plane{fit.Mean(), 0.0, 0.0};
97 }
98
99 10901 std::size_t CraterReliefIndex::KeyHash::operator()(const Key& key) const {
100 const UVec3 hash =
101 10901 Pcg3dHash({static_cast<std::uint32_t>(key.band) ^ static_cast<std::uint32_t>(key.k),
102 10901 static_cast<std::uint32_t>(key.i) ^ static_cast<std::uint32_t>(key.i >> 32),
103 10901 static_cast<std::uint32_t>(key.j) ^ static_cast<std::uint32_t>(key.j >> 32)});
104 10901 return (static_cast<std::size_t>(hash.x) << 32) ^ hash.y;
105 }
106
107 69 CraterReliefIndex::CraterReliefIndex(const std::vector<PreparedCrater>& craters,
108 const CraterProfile& profile, double outlineScale,
109 69 const Surface& surface)
110 69 : craters_(&craters), profile_(profile), outlineScale_(outlineScale), surface_(surface) {
111
2/2
✓ Branch 4 taken 1820 times.
✓ Branch 5 taken 69 times.
 1889 for (const PreparedCrater& prepared : craters) {
112
1/1
✓ Branch 1 taken 1820 times.
 1820 outlineScale_ = std::fmax(outlineScale_, prepared.outline.MaxScale());
113 }
114
2/2
✓ Branch 1 taken 1820 times.
✓ Branch 2 taken 69 times.
 1889 for (std::size_t n = 0; n < craters.size(); ++n) {
115 1820 const Crater& crater = craters[n].crater;
116 1820 const int band = std::ilogb(crater.diameter);
117
1/1
✓ Branch 1 taken 1820 times.
 1820 const double size = BucketSize(band, profile_, outlineScale_);
118
1/1
✓ Branch 1 taken 1820 times.
 1820 const BucketIndex bucket = surface_.Bucket({crater.x, crater.y}, size);
119
2/2
✓ Branch 1 taken 1820 times.
✓ Branch 4 taken 1820 times.
 1820 buckets_[{band, bucket.i, bucket.j, bucket.k}].push_back(n);
120
3/3
✓ Branch 4 taken 1820 times.
✓ Branch 7 taken 262 times.
✓ Branch 8 taken 1558 times.
 1820 if (std::find(bands_.begin(), bands_.end(), band) == bands_.end()) {
121
1/1
✓ Branch 1 taken 262 times.
 262 bands_.push_back(band);
122 }
123 }
124
1/1
✓ Branch 3 taken 69 times.
 69 std::sort(bands_.begin(), bands_.end());
125 69 }
126
127 2108 void CraterReliefIndex::AddBucket(const std::vector<std::size_t>& bucket, double x, double y,
128 std::size_t k, double& relief) const {
129
2/2
✓ Branch 5 taken 2388 times.
✓ Branch 6 taken 1283 times.
 3671 for (const std::size_t n : bucket) {
130
2/2
✓ Branch 0 taken 825 times.
✓ Branch 1 taken 1563 times.
 2388 if (n >= k) {
131 825 break;
132 }
133 1563 const PreparedCrater& prepared = (*craters_)[n];
134 1563 const Crater& crater = prepared.crater;
135 1563 const double radius = 0.5 * crater.diameter;
136
1/1
✓ Branch 1 taken 1563 times.
 1563 const Vec2d offset = surface_.Offset({crater.x, crater.y}, {x, y});
137 1563 const CraterSlopeShape& slope = prepared.slope;
138 3126 const double distance = CraterProfileDistance(
139
1/1
✓ Branch 1 taken 1563 times.
 1563 crater, prepared.outline,
140
1/3
✓ Branch 1 taken 1563 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
 1563 slope.Level() ? offset : SlopeShiftedOffset(slope, offset, radius));
141
3/3
✓ Branch 1 taken 1563 times.
✓ Branch 3 taken 73 times.
✓ Branch 4 taken 1490 times.
 1563 if (distance < profile_.ReachRadii() * radius) {
142
1/1
✓ Branch 1 taken 73 times.
 73 double height = profile_.Height(distance, crater.diameter, prepared.state);
143
4/6
✓ Branch 0 taken 43 times.
✓ Branch 1 taken 30 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 43 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 73 times.
 73 if (height > 0.0 && !slope.Level()) {
144
0/1
✗ Branch 1 not taken.
 ✗ height *= SlopeRaisedShare(slope, offset);
145 }
146 73 relief += prepared.resolution * height;
147 }
148 }
149 2108 }
150
151 257 double CraterReliefIndex::ReliefBefore(double x, double y, std::size_t k) const {
152 257 double relief = 0.0;
153
2/2
✓ Branch 1 taken 74 times.
✓ Branch 2 taken 183 times.
 257 const long long depth = surface_.IsBody() ? 1 : 0;
154
2/2
✓ Branch 5 taken 861 times.
✓ Branch 6 taken 257 times.
 1118 for (const int band : bands_) {
155
1/1
✓ Branch 1 taken 861 times.
 861 const double size = BucketSize(band, profile_, outlineScale_);
156
1/1
✓ Branch 1 taken 861 times.
 861 const BucketIndex centre = surface_.Bucket({x, y}, size);
157
2/2
✓ Branch 0 taken 2583 times.
✓ Branch 1 taken 861 times.
 3444 for (std::int64_t j = centre.j - 1; j <= centre.j + 1; ++j) {
158
2/2
✓ Branch 0 taken 7749 times.
✓ Branch 1 taken 2583 times.
 10332 for (std::int64_t i = centre.i - 1; i <= centre.i + 1; ++i) {
159
2/2
✓ Branch 0 taken 9081 times.
✓ Branch 1 taken 7749 times.
 16830 for (std::int64_t z = centre.k - depth; z <= centre.k + depth; ++z) {
160
1/1
✓ Branch 1 taken 9081 times.
 9081 const auto found = buckets_.find({band, i, j, z});
161
2/2
✓ Branch 2 taken 2108 times.
✓ Branch 3 taken 6973 times.
 9081 if (found != buckets_.end()) {
162
1/1
✓ Branch 2 taken 2108 times.
 2108 AddBucket(found->second, x, y, k, relief);
163 }
164 }
165 }
166 }
167 }
168 257 return relief;
169 }
170
171 181 HeightSampler DatumGround(const HeightSampler& ground, const CraterReliefIndex& index,
172 std::size_t k) {
173 181 return [&ground, &index, k](double x, double y) {
174 253 const float height = ground(x, y);
175
2/2
✓ Branch 1 taken 1 time.
✓ Branch 2 taken 252 times.
 253 if (std::isnan(height)) {
176 1 return height;
177 }
178 252 return static_cast<float>(static_cast<double>(height) + index.ReliefBefore(x, y, k));
179
1/1
✓ Branch 1 taken 181 times.
 362 };
180 }
181
182 244607 double DatumMargin(const CraterApplyOptions& options, const CraterProfile& profile,
183 double maxDiameter) {
184
2/2
✓ Branch 0 taken 239583 times.
✓ Branch 1 taken 5024 times.
 244607 if (options.datum == CraterDatum::None) {
185 239583 return 0.0;
186 }
187 5024 return (options.datumFitRadii + MaxOutlineScale(options.irregularity) * profile.ReachRadii()) *
188 5024 0.5 * maxDiameter;
189 }
190
191 } // namespace terraingen
192