Commit 319395da authored by Purves, Murray's avatar Purves, Murray
Browse files

Added logarithmic option to radixmath/interpolate functions

parent 60af0a0f
......@@ -15,13 +15,14 @@ namespace radix
template <typename T>
std::vector<T> RADIX_PUBLIC interpolateValues(
const std::vector<T> &baseValues, const std::vector<T> &valuesToInterpolate,
const bool circular = false, const T missingValue = T(-9999));
const bool logarithmic = false, const bool circular = false,
const T missingValue = T(-9999));
template <typename T>
std::vector<T> RADIX_PUBLIC interpolateToOtherBaseValues(
const std::vector<T> &baseValues, const std::vector<T> &newBaseValues,
const std::vector<T> &valuesToInterpolate, const bool circular = false,
const T missingValues = T(-9999));
const std::vector<T> &valuesToInterpolate, const bool logarithmic = false,
const bool circular = false, const T missingValues = T(-9999));
} // namespace radix
......
......@@ -22,7 +22,8 @@ namespace radix
template <typename T>
std::vector<T> interpolateValues(const std::vector<T> &baseValues,
const std::vector<T> &valuesToInterpolate,
const bool circular, const T missingValue)
const bool logarithmic, const bool circular,
const T missingValue)
{
std::vector<T> interpolatedValues = valuesToInterpolate;
......@@ -109,10 +110,21 @@ std::vector<T> interpolateValues(const std::vector<T> &baseValues,
}
}
interpolatedValues[j] =
lastGoodValue + ((nextGoodValue - lastGoodValue) *
((baseValues[j] - baseValues[lastGood]) /
(baseValues[nextGood] - baseValues[lastGood])));
if (logarithmic)
{
interpolatedValues[j] =
lastGoodValue +
((nextGoodValue - lastGoodValue) *
((log(baseValues[j]) - log(baseValues[lastGood])) /
(log(baseValues[nextGood]) - log(baseValues[lastGood]))));
}
else
{
interpolatedValues[j] =
lastGoodValue + ((nextGoodValue - lastGoodValue) *
((baseValues[j] - baseValues[lastGood]) /
(baseValues[nextGood] - baseValues[lastGood])));
}
if (circular)
{
......@@ -138,8 +150,8 @@ std::vector<T> interpolateValues(const std::vector<T> &baseValues,
template <typename T>
std::vector<T> interpolateToOtherBaseValues(
const std::vector<T> &baseValues, const std::vector<T> &newBaseValues,
const std::vector<T> &valuesToInterpolate, const bool circular,
const T missingValues)
const std::vector<T> &valuesToInterpolate, const bool logarithmic,
const bool circular, const T missingValues)
{
std::vector<T> finalInterpolatedValues(newBaseValues.size());
......@@ -148,7 +160,7 @@ std::vector<T> interpolateToOtherBaseValues(
// First interpolate across the existing data to ensure there are no missing
radix_line(" Initial interpolation to remove missing values in input data");
std::vector<T> initialInterpolatedValues = interpolateValues(
baseValues, valuesToInterpolate, circular, missingValues);
baseValues, valuesToInterpolate, logarithmic, circular, missingValues);
// If this fails, we have to fail this one too
if (initialInterpolatedValues.size() == 0)
{
......@@ -250,11 +262,23 @@ std::vector<T> interpolateToOtherBaseValues(
}
else
{
finalInterpolatedValues[i] =
lastGoodValue +
((nextGoodValue - lastGoodValue) *
((newBaseValues[i] - baseValues[lowerBaseIndex]) /
(baseValues[higherBaseIndex] - baseValues[lowerBaseIndex])));
if (logarithmic)
{
finalInterpolatedValues[i] =
lastGoodValue +
((nextGoodValue - lastGoodValue) *
((log(newBaseValues[i]) - log(baseValues[lowerBaseIndex])) /
(log(baseValues[higherBaseIndex]) -
log(baseValues[lowerBaseIndex]))));
}
else
{
finalInterpolatedValues[i] =
lastGoodValue +
((nextGoodValue - lastGoodValue) *
((newBaseValues[i] - baseValues[lowerBaseIndex]) /
(baseValues[higherBaseIndex] - baseValues[lowerBaseIndex])));
}
}
if (circular)
......
......@@ -14,10 +14,10 @@ namespace radix
{
template<typename T>
std::vector<T> interpolateValues(std::vector<T> const & baseValues, std::vector<T> const & valuesToInterpolate, const bool circular, const T missingValue);
std::vector<T> interpolateValues(std::vector<T> const & baseValues, std::vector<T> const & valuesToInterpolate, const bool logarithmic, const bool circular, const T missingValue);
template<typename T>
std::vector<T> interpolateToOtherBaseValues(std::vector<T> const & baseValues, std::vector<T> const & newBaseValues, std::vector<T> const & valuesToInterpolate, const bool circular, const T missingValue);
std::vector<T> interpolateToOtherBaseValues(std::vector<T> const & baseValues, std::vector<T> const & newBaseValues, std::vector<T> const & valuesToInterpolate, const bool logarithmic, const bool circular, const T missingValue);
%template(interpolate_values) interpolateValues<double>;
%template(interpolate_to_other_base_values) interpolateToOtherBaseValues<double>;
......
......@@ -6,7 +6,7 @@
using namespace radix;
TEST(Interpolate, InternalInterpolation)
TEST(Interpolate, InternalInterpolationLinear)
{
float missingValue = -9999.f;
float tolerance = 0.0001f;
......@@ -42,7 +42,38 @@ TEST(Interpolate, InternalInterpolation)
std::vector<float> interpValues3{16.f, missingValue, missingValue, 2.f};
std::vector<float> expectValues3{16.f, 8.f, 4.f, 2.f};
std::vector<float> testValues3 =
interpolateValues(baseValues3, interpValues3, false, missingValue);
interpolateValues(baseValues3, interpValues3, false, false, missingValue);
ASSERT_EQ(expectValues3.size(), testValues3.size());
for (size_t i = 0; i < testValues3.size(); ++i)
{
EXPECT_NEAR(expectValues3[i], testValues3[i], expectValues3[i] * tolerance);
}
}
TEST(Interpolate, InternalInterpolationLog)
{
float missingValue = -9999.f;
float tolerance = 0.0001f;
std::vector<float> baseValues1{1.f, 10.f, 100.f, 1000.f,
10000.f, 100000.f, 1000000.f};
std::vector<float> interpValues1{1.f, missingValue, 3.f, missingValue,
6.f, missingValue, 1.f};
std::vector<float> expectValues1{1.f, 2.f, 3.f, 4.5, 6.f, 3.5, 1.f};
std::vector<float> testValues1 =
interpolateValues(baseValues1, interpValues1, true);
ASSERT_EQ(expectValues1.size(), testValues1.size());
for (size_t i = 0; i < testValues1.size(); ++i)
{
EXPECT_NEAR(expectValues1[i], testValues1[i], expectValues1[i] * tolerance);
}
missingValue = 3e12f;
std::vector<float> baseValues3{800.f, 400.f, 200.f, 100.f};
std::vector<float> interpValues3{8.f, missingValue, missingValue, 2.f};
std::vector<float> expectValues3{8.f, 6.f, 4.f, 2.f};
std::vector<float> testValues3 =
interpolateValues(baseValues3, interpValues3, true, false, missingValue);
ASSERT_EQ(expectValues3.size(), testValues3.size());
for (size_t i = 0; i < testValues3.size(); ++i)
{
......@@ -81,7 +112,7 @@ TEST(Interpolate, CircularInterpolation)
350.f, missingValue, 30.f};
std::vector<float> expectValues1{90.f, 180.f, 270.f, 350.f, 10.f, 30.f};
std::vector<float> testValues1 =
interpolateValues(baseValues1, interpValues1, true, missingValue);
interpolateValues(baseValues1, interpValues1, false, true, missingValue);
ASSERT_EQ(expectValues1.size(), testValues1.size());
for (size_t i = 0; i < testValues1.size(); ++i)
{
......@@ -164,6 +195,47 @@ TEST(Interpolate, NewBaseValues)
}
}
TEST(Interpolate, NewBaseValuesLog)
{
float missingValue = -9999.f;
float tolerance = 0.0001f;
std::vector<float> baseValues1{1000.f, 4000.f, 16000.f};
std::vector<float> newBaseValues1{500.f, 2000.f, 8000.f, 32000.f};
std::vector<float> interpValues1{10.5f, 20.5f, 30.5f};
std::vector<float> expectValues1{10.5f, 15.5f, 25.5f, 30.5f};
std::vector<float> testValues1 = interpolateToOtherBaseValues(
baseValues1, newBaseValues1, interpValues1, true);
ASSERT_EQ(expectValues1.size(), testValues1.size());
for (size_t i = 0; i < testValues1.size(); ++i)
{
EXPECT_NEAR(expectValues1[i], testValues1[i], expectValues1[i] * tolerance);
}
std::vector<float> interpValuesMissing{10.5f, missingValue, 30.5f};
std::vector<float> testValuesMissing = interpolateToOtherBaseValues(
baseValues1, newBaseValues1, interpValuesMissing, true);
ASSERT_EQ(expectValues1.size(), testValuesMissing.size());
for (size_t i = 0; i < testValuesMissing.size(); ++i)
{
EXPECT_NEAR(expectValues1[i], testValuesMissing[i],
expectValues1[i] * tolerance);
}
// High-to-low base values
std::vector<float> baseValues2{16000.f, 4000.f, 1000.f};
std::vector<float> newBaseValues2{32000.f, 8000.f, 2000.f, 500.f};
std::vector<float> interpValues2{30.5f, 20.5f, 10.5f};
std::vector<float> expectValues2{30.5f, 25.5f, 15.5f, 10.5f};
std::vector<float> testValues2 = interpolateToOtherBaseValues(
baseValues2, newBaseValues2, interpValues2, true);
ASSERT_EQ(expectValues2.size(), testValues2.size());
for (size_t i = 0; i < testValues2.size(); ++i)
{
EXPECT_NEAR(expectValues2[i], testValues2[i], expectValues2[i] * tolerance);
}
}
TEST(Interpolate, NewBaseValuesDoubles)
{
double missingValue = -9999.0;
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment