distribution/piecewise_power_law.cpp

Line	Branch	Exec	Source
1			#include "terraingen/distribution/piecewise_power_law.h"
2
3			#include "terraingen/distribution/combined.h"
4			#include "terraingen/util/checks.h"
5
6			#include <cmath>
7			#include <stdexcept>
8			#include <utility>
9
10			namespace terraingen {
11
12		45	PiecewisePowerLawSfd::PiecewisePowerLawSfd(std::vector<PowerLawSegment> segments, SfdKind kind,
13			DiameterRange range, bool allowDiscontinuities,
14		45	double tolerance)
15		45	: SizeFrequencyDistribution(kind, range), segments_(std::move(segments)) {
16	2/2 ✓ Branch 1 taken 1 time. ✓ Branch 2 taken 44 times.	45	if (segments_.empty()) {
17	1/1 ✓ Branch 2 taken 1 time.	1	throw std::invalid_argument("A piecewise power law needs at least one segment");
18			}
19	2/2 ✓ Branch 1 taken 1 time. ✓ Branch 2 taken 43 times.	44	if (segments_.front().start != range.min) {
20	1/1 ✓ Branch 2 taken 1 time.	1	throw std::invalid_argument("The first segment must start at the range minimum");
21			}
22	2/2 ✓ Branch 1 taken 89 times. ✓ Branch 2 taken 35 times.	124	for (std::size_t i = 0; i < segments_.size(); ++i) {
23		89	const PowerLawSegment& segment = segments_[i];
24	6/6 ✓ Branch 1 taken 88 times. ✓ Branch 2 taken 1 time. ✓ Branch 4 taken 1 time. ✓ Branch 5 taken 87 times. ✓ Branch 6 taken 2 times. ✓ Branch 7 taken 87 times.	89	if (!IsPositive(segment.coefficient) \|\| !IsPositive(segment.exponent)) {
25	1/1 ✓ Branch 2 taken 2 times.	2	throw std::invalid_argument("Segment coefficients and exponents must be positive");
26			}
27	2/2 ✓ Branch 0 taken 41 times. ✓ Branch 1 taken 46 times.	87	if (i == 0) {
28		41	continue;
29			}
30	6/6 ✓ Branch 1 taken 45 times. ✓ Branch 2 taken 1 time. ✓ Branch 4 taken 44 times. ✓ Branch 5 taken 1 time. ✓ Branch 6 taken 3 times. ✓ Branch 7 taken 43 times.	90	if (!(std::isfinite(segment.start) && segment.start > segments_[i - 1].start &&
31	2/2 ✓ Branch 0 taken 1 time. ✓ Branch 1 taken 43 times.	44	segment.start < range.max)) {
32	1/1 ✓ Branch 2 taken 3 times.	3	throw std::invalid_argument("Segment starts must increase within the range");
33			}
34	1/1 ✓ Branch 3 taken 40 times.	43	detail::CheckJoin(segments_[i - 1].Cumulative(segment.start),
35		43	segment.Cumulative(segment.start), segment.start, tolerance,
36			allowDiscontinuities);
37			}
38		55	}
39
40		388763	const PowerLawSegment& PiecewisePowerLawSfd::SegmentAt(double diameter) const {
41		388763	return segments_[detail::IndexOfStart(segments_, diameter)];
42			}
43
44		388741	double PiecewisePowerLawSfd::CumulativeAt(double diameter) const {
45		388741	return SegmentAt(diameter).Cumulative(diameter);
46			}
47
48		22	double PiecewisePowerLawSfd::IncrementalAt(double diameter) const {
49		22	return SegmentAt(diameter).Incremental(diameter);
50			}
51
52		89285	double PiecewisePowerLawSfd::InvertCumulative(double target, double dMin, double dMax) const {
53		89285	std::size_t index = 0;
54	2/2 ✓ Branch 1 taken 91969 times. ✓ Branch 2 taken 3722 times.	95691	while (index + 1 < segments_.size()) {
55		91969	const PowerLawSegment& segment = segments_[index];
56		91969	const double high = std::fmin(dMax, segments_[index + 1].start);
57	6/6 ✓ Branch 0 taken 91298 times. ✓ Branch 1 taken 671 times. ✓ Branch 3 taken 85563 times. ✓ Branch 4 taken 5735 times. ✓ Branch 5 taken 85563 times. ✓ Branch 6 taken 6406 times.	91969	if (std::fmax(dMin, segment.start) < high && target >= segment.Cumulative(high)) {
58		85563	break;
59			}
60		6406	++index;
61			}
62		89285	const PowerLawSegment& segment = segments_[index];
63		89285	const double low = std::fmax(dMin, segment.start);
64	2/2 ✓ Branch 1 taken 930 times. ✓ Branch 2 taken 88355 times.	89285	if (target > segment.Cumulative(low)) {
65		930	return low;
66			}
67		88355	return segment.Invert(target);
68			}
69
70			} // namespace terraingen
71

TerrainGen coverage