Mercurial > hg > orthanc-stone
changeset 144:9b83f30fc1c0 wasm
bilinear and trilinear interpolation
| author | Sebastien Jodogne <s.jodogne@gmail.com> |
|---|---|
| date | Mon, 22 Jan 2018 15:55:03 +0100 |
| parents | 58c545177c1c |
| children | 6941b98cf0fd |
| files | Framework/Toolbox/GeometryToolbox.h UnitTestsSources/UnitTestsMain.cpp |
| diffstat | 2 files changed, 104 insertions(+), 3 deletions(-) [+] |
line wrap: on
line diff
--- a/Framework/Toolbox/GeometryToolbox.h Mon Jan 22 15:16:19 2018 +0100 +++ b/Framework/Toolbox/GeometryToolbox.h Mon Jan 22 15:55:03 2018 +0100 @@ -28,10 +28,11 @@ # include <boost/serialization/array_wrapper.hpp> #endif +#include <Core/DicomFormat/DicomMap.h> + #include <boost/numeric/ublas/matrix.hpp> #include <boost/numeric/ublas/vector.hpp> - -#include <Core/DicomFormat/DicomMap.h> +#include <cassert> namespace OrthancStone { @@ -125,6 +126,93 @@ Matrix CreateRotationMatrixAlongY(double a); - Matrix CreateRotationMatrixAlongZ(double a); + Matrix CreateRotationMatrixAlongZ(double a); + + inline float ComputeBilinearInterpolationInternal(float x, + float y, + float f00, // source(x, y) + float f01, // source(x + 1, y) + float f10, // source(x, y + 1) + float f11); // source(x + 1, y + 1) + + inline float ComputeBilinearInterpolation(float x, + float y, + float f00, // source(x, y) + float f01, // source(x + 1, y) + float f10, // source(x, y + 1) + float f11); // source(x + 1, y + 1) + + inline float ComputeTrilinearInterpolation(float x, + float y, + float z, + float f000, // source(x, y, z) + float f001, // source(x + 1, y, z) + float f010, // source(x, y + 1, z) + float f011, // source(x + 1, y + 1, z) + float f100, // source(x, y, z + 1) + float f101, // source(x + 1, y, z + 1) + float f110, // source(x, y + 1, z + 1) + float f111); // source(x + 1, y + 1, z + 1) }; } + + +float OrthancStone::GeometryToolbox::ComputeBilinearInterpolationInternal(float x, + float y, + float f00, + float f01, + float f10, + float f11) +{ + // This function works on fractional parts + assert(x >= 0 && y >= 0 && x < 1 && y < 1); + + // https://en.wikipedia.org/wiki/Bilinear_interpolation#Unit_square + return f00 * (1 - x) * (1 - y) + f01 * x * (1 - y) + f10 * (1 - x) * y + f11 * x * y; +} + + +float OrthancStone::GeometryToolbox::ComputeBilinearInterpolation(float x, + float y, + float f00, + float f01, + float f10, + float f11) +{ + assert(x >= 0 && y >= 0); + + // Compute the fractional part of (x,y) (this is the "floor()" + // operation on positive integers) + float xx = x - static_cast<float>(static_cast<int>(x)); + float yy = y - static_cast<float>(static_cast<int>(y)); + + return ComputeBilinearInterpolationInternal(xx, yy, f00, f01, f10, f11); +} + + +float OrthancStone::GeometryToolbox::ComputeTrilinearInterpolation(float x, + float y, + float z, + float f000, + float f001, + float f010, + float f011, + float f100, + float f101, + float f110, + float f111) +{ + assert(x >= 0 && y >= 0 && z >= 0); + + float xx = x - static_cast<float>(static_cast<int>(x)); + float yy = y - static_cast<float>(static_cast<int>(y)); + float zz = z - static_cast<float>(static_cast<int>(z)); + + // "In practice, a trilinear interpolation is identical to two + // bilinear interpolation combined with a linear interpolation" + // https://en.wikipedia.org/wiki/Trilinear_interpolation#Method + float a = ComputeBilinearInterpolationInternal(xx, yy, f000, f001, f010, f011); + float b = ComputeBilinearInterpolationInternal(xx, yy, f100, f101, f110, f111); + + return (1 - zz) * a + zz * b; +}
--- a/UnitTestsSources/UnitTestsMain.cpp Mon Jan 22 15:16:19 2018 +0100 +++ b/UnitTestsSources/UnitTestsMain.cpp Mon Jan 22 15:55:03 2018 +0100 @@ -42,6 +42,7 @@ #include <boost/math/special_functions/round.hpp> +#if 0 namespace OrthancStone { class Tata : public OrthancSlicesLoader::ICallback @@ -129,6 +130,18 @@ oracle.Stop(); } +#endif + + +TEST(GeometryToolbox, Interpolation) +{ + // https://en.wikipedia.org/wiki/Bilinear_interpolation#Application_in_image_processing + ASSERT_FLOAT_EQ(146.1f, OrthancStone::GeometryToolbox::ComputeBilinearInterpolation(20.5f, 14.2f, 91, 210, 162, 95)); + + ASSERT_FLOAT_EQ(123.35f, OrthancStone::GeometryToolbox::ComputeTrilinearInterpolation(20.5f, 15.2f, 5.7f, + 91, 210, 162, 95, + 51, 190, 80, 92)); +} int main(int argc, char **argv)