Mercurial > hg > orthanc-stone
changeset 1286:ddb6676bbcbf bugs/2020-02-invisible-slice
Added Histogram functions (new files: ImageToolbox...)
| author | Benjamin Golinvaux <bgo@osimis.io> |
|---|---|
| date | Fri, 14 Feb 2020 14:59:32 +0100 |
| parents | 1bb54005e2bd |
| children | 8e82fdc6200e |
| files | Applications/Samples/CMakeLists.txt Framework/Toolbox/ImageToolbox.cpp Framework/Toolbox/ImageToolbox.h Resources/CMake/OrthancStoneConfiguration.cmake UnitTestsSources/ImageToolboxTests.cpp |
| diffstat | 5 files changed, 644 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- a/Applications/Samples/CMakeLists.txt Tue Feb 11 14:38:36 2020 +0100 +++ b/Applications/Samples/CMakeLists.txt Fri Feb 14 14:59:32 2020 +0100 @@ -251,6 +251,7 @@ add_executable(UnitTests ${GOOGLE_TEST_SOURCES} ${ORTHANC_STONE_ROOT}/UnitTestsSources/GenericToolboxTests.cpp + ${ORTHANC_STONE_ROOT}/UnitTestsSources/ImageToolboxTests.cpp ${ORTHANC_STONE_ROOT}/UnitTestsSources/PixelTestPatternsTests.cpp ${ORTHANC_STONE_ROOT}/UnitTestsSources/TestCommands.cpp ${ORTHANC_STONE_ROOT}/UnitTestsSources/TestMessageBroker.cpp
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/Framework/Toolbox/ImageToolbox.cpp Fri Feb 14 14:59:32 2020 +0100 @@ -0,0 +1,304 @@ +/** + * Stone of Orthanc + * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics + * Department, University Hospital of Liege, Belgium + * Copyright (C) 2017-2020 Osimis S.A., Belgium + * + * This program is free software: you can redistribute it and/or + * modify it under the terms of the GNU Affero General Public License + * as published by the Free Software Foundation, either version 3 of + * the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, but + * WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + **/ + +#include "ImageToolbox.h" + +#include "../StoneException.h" + +#include <Core/Images/ImageProcessing.h> +#include <Core/Images/PixelTraits.h> + +#include <Core/Logging.h> +#include <Core/OrthancException.h> + +#include <boost/static_assert.hpp> +#include <boost/type_traits.hpp> + +#include <vector> + +namespace OrthancStone +{ + namespace + { + using Orthanc::PixelTraits; + using Orthanc::PixelFormat; + using Orthanc::ImageAccessor; + using Orthanc::PixelFormat; + + template<typename PixelFormat Format> + class PixelBinner + { + // "PixelBinner requires an arithmetic (integer or floating-point) pixel format" + typedef typename PixelTraits<Format>::PixelType PixelType; + BOOST_STATIC_ASSERT(boost::is_arithmetic<PixelType>::value); + + public: + PixelBinner(HistogramData& hd, double minValue, double maxValue) + : hd_(hd) + , minValue_(minValue) + , maxValue_(maxValue) + , division_(1.0 / hd_.binSize) + { + ORTHANC_ASSERT(hd_.bins.size() > 0); + ORTHANC_ASSERT(maxValue > minValue); + } + + ORTHANC_FORCE_INLINE void AddPixel(PixelType p) + { + if (p <= minValue_) + { + hd_.bins[0] += 1; + } + else if (p >= maxValue_) + { + hd_.bins.back() += 1; + } + else + { + double distanceFromMin = p - minValue_; + size_t binIndex = static_cast<size_t>( + std::floor(distanceFromMin * division_)); + if (binIndex >= hd_.bins.size()) + binIndex = hd_.bins.size() - 1; + hd_.bins[binIndex] += 1; + } + } + private: + HistogramData& hd_; + double minValue_; + double maxValue_; + double division_; + }; + + template<PixelFormat Format> + struct Histogram + { + typedef typename PixelTraits<Format>::PixelType PixelType; + + static void Apply(const Orthanc::ImageAccessor& img, HistogramData& hd, + double minValue = 0, + double maxValue = 0) + { + ORTHANC_ASSERT(Format == img.GetFormat(), + "Internal error. Wrong template histogram type"); + + const uint8_t* buffer = reinterpret_cast<const uint8_t*>( + img.GetConstBuffer()); + + const size_t pitch = img.GetPitch(); + const size_t bytesPerPix = img.GetBytesPerPixel(); + + const size_t height = img.GetHeight(); + const size_t width = img.GetHeight(); + + if ((minValue == 0) && (maxValue == 0)) + { + throw Orthanc::OrthancException(Orthanc::ErrorCode_NotImplemented); + //ORTHANC_ASSERT(boost::is_integral<PixelType>::value, + // "Min and max values must be supplied for float-based histogram"); + // + //PixelTraits<Format>::SetMinValue(minValue); + //PixelTraits<Format>::SetMaxValue(maxValue); + } + + hd.minValue = minValue; + + // the following code is not really pretty but ensures + size_t numBins = static_cast<size_t>( + std::ceil((maxValue - minValue) / hd.binSize)); + + hd.bins.resize(numBins); + std::fill(hd.bins.begin(), hd.bins.end(), 0); + + PixelBinner<Format> binner(hd, minValue, maxValue); + for (uint32_t y = 0; y < height; ++y) + { + const PixelType* curPix = reinterpret_cast<const PixelType*>( + img.GetConstRow(y)); + + for (uint32_t x = 0; x < width; x++, curPix++) + { + binner.AddPixel(*curPix); + } + } + } + }; + + + template<PixelFormat Format> + struct ComputeMinMax__ + { + typedef typename PixelTraits<Format>::PixelType PixelType; + + static void Apply(const Orthanc::ImageAccessor& img, + PixelType& minValue, PixelType& maxValue) + { + ORTHANC_ASSERT(Format == img.GetFormat(), + "Internal error. Wrong template histogram type"); + + const uint8_t* buffer = reinterpret_cast<const uint8_t*>( + img.GetConstBuffer()); + + const size_t pitch = img.GetPitch(); + const size_t bytesPerPix = img.GetBytesPerPixel(); + + const size_t height = img.GetHeight(); + const size_t width = img.GetHeight(); + + if (height * width == 0) + { + throw Orthanc::OrthancException(Orthanc::ErrorCode_IncompatibleImageFormat); + } + + // min and max are crossed below. Think about it. This is OK :) + PixelTraits<Format>::SetMaxValue(minValue); + PixelTraits<Format>::SetMinValue(maxValue); + + for (uint32_t y = 0; y < height; ++y) + { + const PixelType* curPix = reinterpret_cast<const PixelType*>( + img.GetConstRow(y)); + + for (uint32_t x = 0; x < width; x++, curPix++) + { + if (*curPix <= minValue) + minValue = *curPix; + if (*curPix >= maxValue) + maxValue = *curPix; + } + } + } + }; + + template<PixelFormat Format> + void ComputeMinMax_(const Orthanc::ImageAccessor& img, + double& minValue, double& maxValue) + { + typedef typename PixelTraits<Format>::PixelType PixelType; + PixelType minValuePix = PixelType(); + PixelType maxValuePix = PixelType(); + ComputeMinMax__<Format>::Apply(img, minValuePix, maxValuePix); + minValue = static_cast<double>(minValuePix); + maxValue = static_cast<double>(maxValuePix); + } + + template<PixelFormat Format> + void ComputeHistogram_(const Orthanc::ImageAccessor& img, HistogramData& hd) + { + typedef typename PixelTraits<Format>::PixelType PixelType; + PixelType minValue = PixelType(); + PixelType maxValue = PixelType(); + ComputeMinMax__<Format>::Apply(img, minValue, maxValue); + + // make bins a little bigger to center integer pixel values + Histogram<Format>::Apply(img, hd, + static_cast<double>(minValue) - 0.5, + static_cast<double>(maxValue) + 0.5); + } + } + + void ComputeHistogram(const Orthanc::ImageAccessor& img, + HistogramData& hd, double binSize) + { + using namespace Orthanc; + + hd.binSize = binSize; + + // dynamic/static bridge + switch (img.GetFormat()) + { + case PixelFormat_Grayscale8: + ComputeHistogram_<PixelFormat_Grayscale8> (img, hd); + break; + case PixelFormat_Grayscale16: + ComputeHistogram_<PixelFormat_Grayscale16> (img, hd); + break; + case PixelFormat_SignedGrayscale16: + ComputeHistogram_<PixelFormat_SignedGrayscale16>(img, hd); + break; + case PixelFormat_Float32: + ComputeHistogram_<PixelFormat_Float32> (img, hd); + break; + case PixelFormat_Grayscale32: + ComputeHistogram_<PixelFormat_Grayscale32> (img, hd); + break; + case PixelFormat_Grayscale64: + ComputeHistogram_<PixelFormat_Grayscale64> (img, hd); + break; + default: + throw Orthanc::OrthancException(Orthanc::ErrorCode_IncompatibleImageFormat); + } + } + + void ComputeMinMax(const Orthanc::ImageAccessor& img, + double& minValue, double& maxValue) + { + using namespace Orthanc; + + // dynamic/static bridge + switch (img.GetFormat()) + { + case PixelFormat_Grayscale8: + ComputeMinMax_<PixelFormat_Grayscale8> (img, minValue, maxValue); + break; + case PixelFormat_Grayscale16: + ComputeMinMax_<PixelFormat_Grayscale16> (img, minValue, maxValue); + break; + case PixelFormat_SignedGrayscale16: + ComputeMinMax_<PixelFormat_SignedGrayscale16>(img, minValue, maxValue); + break; + case PixelFormat_Float32: + ComputeMinMax_<PixelFormat_Float32> (img, minValue, maxValue); + break; + case PixelFormat_Grayscale32: + ComputeMinMax_<PixelFormat_Grayscale32> (img, minValue, maxValue); + break; + case PixelFormat_Grayscale64: + ComputeMinMax_<PixelFormat_Grayscale64> (img, minValue, maxValue); + break; + default: + throw Orthanc::OrthancException(Orthanc::ErrorCode_IncompatibleImageFormat); + } + + } + + void DumpHistogramResult(std::string& s, const HistogramData& hd) + { + std::stringstream ss; + ss << "Histogram:\n"; + ss << "==========\n"; + ss << "\n"; + ss << "minValue : " << hd.minValue << "\n"; + ss << "binSize : " << hd.binSize << "\n"; + ss << "bins.size() : " << hd.bins.size() << "\n"; + ss << "bins :\n"; + double curBinStart = hd.minValue; + size_t pixCount = 0; + for (size_t i = 0; i < hd.bins.size(); ++i) + { + ss << "index: " << i << " (from " << curBinStart << " to " + << curBinStart + hd.binSize << ") : " << hd.bins[i] << " pixels\n"; + curBinStart += hd.binSize; + pixCount += hd.bins[i]; + } + ss << "total pix. count: " << pixCount << "\n"; + s = ss.str(); + } +}
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/Framework/Toolbox/ImageToolbox.h Fri Feb 14 14:59:32 2020 +0100 @@ -0,0 +1,78 @@ +/** + * Stone of Orthanc + * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics + * Department, University Hospital of Liege, Belgium + * Copyright (C) 2017-2020 Osimis S.A., Belgium + * + * This program is free software: you can redistribute it and/or + * modify it under the terms of the GNU Affero General Public License + * as published by the Free Software Foundation, either version 3 of + * the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, but + * WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + **/ + + +#pragma once + +#include "../StoneEnumerations.h" +#include "LinearAlgebra.h" + +#include <Core/Images/ImageAccessor.h> + +namespace OrthancStone +{ + + /** + This structure represents the result of an histogram computation + + bins[0] contains the values in [minValue , minValue + binSize [ + bins[1] contains the values in [minValue + binSize, minValue + 2*binSize [ + bins[2] contains the values in [minValue + 2*binSize, minValue + 3*binSize [ + ... + bins[N-1] contains the values in [minValue + (N-1)*binSize, minValue + N*binSize [ + + */ + struct HistogramData + { + std::vector<size_t> bins; + double minValue; + double binSize; + }; + + /** + Dumps the supplied histogram to the supplied strings + */ + void DumpHistogramResult(std::string& s, const HistogramData& hd); + + /** + This will compute the histogram of the supplied image (count the number of + pixels). + + The image must contain arithmetic pixels (that is, having a single component, + integer or float). Compound pixel types like RGB, YUV are not supported and + will cause this function to throw an exception. + + The range of available values will be split in sets of size `binSize`, and + each set will contain the number of pixels in the given bin + (see HistogramResult above). + */ + void ComputeHistogram(const Orthanc::ImageAccessor& img, + HistogramData& hd, double binSize); + + + /** + Computes the min max values in an image + */ + void ComputeMinMax(const Orthanc::ImageAccessor& img, + double& minValue, double& maxValue); + +} + +extern int OrthancStone_Internals_dump_LoadTexture_histogram;
--- a/Resources/CMake/OrthancStoneConfiguration.cmake Tue Feb 11 14:38:36 2020 +0100 +++ b/Resources/CMake/OrthancStoneConfiguration.cmake Fri Feb 14 14:59:32 2020 +0100 @@ -570,6 +570,8 @@ ${ORTHANC_STONE_ROOT}/Framework/Toolbox/GeometryToolbox.h ${ORTHANC_STONE_ROOT}/Framework/Toolbox/ImageGeometry.cpp ${ORTHANC_STONE_ROOT}/Framework/Toolbox/ImageGeometry.h + ${ORTHANC_STONE_ROOT}/Framework/Toolbox/ImageToolbox.cpp + ${ORTHANC_STONE_ROOT}/Framework/Toolbox/ImageToolbox.h ${ORTHANC_STONE_ROOT}/Framework/Toolbox/LinearAlgebra.cpp ${ORTHANC_STONE_ROOT}/Framework/Toolbox/LinearAlgebra.h ${ORTHANC_STONE_ROOT}/Framework/Toolbox/PixelTestPatterns.cpp
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/UnitTestsSources/ImageToolboxTests.cpp Fri Feb 14 14:59:32 2020 +0100 @@ -0,0 +1,259 @@ +/** + * Stone of Orthanc + * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics + * Department, University Hospital of Liege, Belgium + * Copyright (C) 2017-2020 Osimis S.A., Belgium + * + * This program is free software: you can redistribute it and/or + * modify it under the terms of the GNU Affero General Public License + * as published by the Free Software Foundation, either version 3 of + * the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, but + * WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + **/ + +#include <Framework/Toolbox/ImageToolbox.h> + +// #include <boost/chrono.hpp> +// #include <boost/lexical_cast.hpp> + +#include <Core/Images/Image.h> +#include <Core/Images/PixelTraits.h> + +#include "stdint.h" + +#include "gtest/gtest.h" + +#include <cmath> + + +TEST(ImageToolbox, SimpleHisto_Grayscale8_BinSize1) +{ + using OrthancStone::HistogramData; + using OrthancStone::DumpHistogramResult; + using OrthancStone::ComputeHistogram; + + const unsigned int W = 16; + const unsigned int H = 16; + + // 256/17 = 15,... + // 256 % 17 = 1 + // 0 will be 16 times + // 1 will be 15 times + // 2 will be 15 times + // ... + // 16 will be 15 times + + size_t pixCounter = 0; + + std::auto_ptr<Orthanc::Image> image(new Orthanc::Image( + Orthanc::PixelFormat_Grayscale8, W, H, false)); + + for (unsigned int y = 0; y < H; ++y) + { + uint8_t* buffer = reinterpret_cast<uint8_t*>(image->GetRow(y)); + for (unsigned int x = 0; x < W; ++x, ++buffer, ++pixCounter) + { + *buffer = static_cast<uint8_t>(pixCounter % 17); + } + } + + HistogramData hd; + ComputeHistogram(*image, hd, 1); + ASSERT_EQ(-0.5, hd.minValue); + ASSERT_EQ(17u, hd.bins.size()); + ASSERT_EQ(16u, hd.bins[0]); + for (size_t i = 1; i < hd.bins.size(); ++i) + ASSERT_EQ(15u, hd.bins[i]); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale8_BinSize1_FormatString) +{ + using OrthancStone::HistogramData; + using OrthancStone::DumpHistogramResult; + using OrthancStone::ComputeHistogram; + + const unsigned int W = 16; + const unsigned int H = 16; + + // 256/17 = 15,... + // 256 % 17 = 1 + // 0 will be 16 times + // 1 will be 15 times + // 2 will be 15 times + // ... + // 16 will be 15 times + + size_t pixCounter = 0; + + std::auto_ptr<Orthanc::Image> image(new Orthanc::Image( + Orthanc::PixelFormat_Grayscale8, W, H, false)); + + for (unsigned int y = 0; y < H; ++y) + { + uint8_t* buffer = reinterpret_cast<uint8_t*>(image->GetRow(y)); + for (unsigned int x = 0; x < W; ++x, ++buffer, ++pixCounter) + { + *buffer = static_cast<uint8_t>(pixCounter % 17); + } + } + + HistogramData hd; + ComputeHistogram(*image, hd, 1); + + // void DumpHistogramResult(std::string& s, const HistogramData& hd) + std::string s; + DumpHistogramResult(s, hd); + std::cout << s; +} + +template<Orthanc::PixelFormat Format> +void SimpleHisto_T_BinSize1_2() +{ + using OrthancStone::HistogramData; + using OrthancStone::DumpHistogramResult; + using OrthancStone::ComputeHistogram; + + const unsigned int W = 16; + const unsigned int H = 16; + + // 256/17 = 15,... + // 256 % 17 = 1 + // 0 will be 16 times + // 1 will be 15 times + // 2 will be 15 times + // ... + // 16 will be 15 times + + size_t pixCounter = 0; + + std::auto_ptr<Orthanc::Image> image(new Orthanc::Image( + Format, W, H, false)); + + typedef Orthanc::PixelTraits<Format>::PixelType PixelType; + + PixelType pixValue = 0; + + for (unsigned int y = 0; y < H; ++y) + { + PixelType* buffer = reinterpret_cast<PixelType*>(image->GetRow(y)); + for (unsigned int x = 0; x < W; ++x, ++buffer, ++pixCounter) + { + // 0..99 0..99 0..55 + *buffer = pixValue; + pixValue++; + if (pixValue >= 100) + pixValue = 0; + } + } + + HistogramData hd; + ComputeHistogram(*image, hd, 1); + ASSERT_EQ(-0.5, hd.minValue); + ASSERT_EQ(100u, hd.bins.size()); + for (size_t i = 0; i <= 55; ++i) + ASSERT_EQ(3u, hd.bins[i]); + for (size_t i = 56; i <= 99; ++i) + ASSERT_EQ(2u, hd.bins[i]); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale8_BinSize1_2) +{ + SimpleHisto_T_BinSize1_2<Orthanc::PixelFormat_Grayscale8>(); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale16_BinSize1_2) +{ + SimpleHisto_T_BinSize1_2<Orthanc::PixelFormat_Grayscale16>(); +} + +TEST(ImageToolbox, SimpleHisto_SignedGrayscale16_BinSize1_2) +{ + SimpleHisto_T_BinSize1_2<Orthanc::PixelFormat_SignedGrayscale16>(); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale32_BinSize1_2) +{ + SimpleHisto_T_BinSize1_2<Orthanc::PixelFormat_Grayscale32>(); +} + +template<Orthanc::PixelFormat Format> +void SimpleHisto_T_BinSize10_2() +{ + using OrthancStone::HistogramData; + using OrthancStone::DumpHistogramResult; + using OrthancStone::ComputeHistogram; + + const unsigned int W = 16; + const unsigned int H = 16; + + // 256/17 = 15,... + // 256 % 17 = 1 + // 0 will be 16 times + // 1 will be 15 times + // 2 will be 15 times + // ... + // 16 will be 15 times + + size_t pixCounter = 0; + + std::auto_ptr<Orthanc::Image> image(new Orthanc::Image( + Format, W, H, false)); + + typedef Orthanc::PixelTraits<Format>::PixelType PixelType; + + PixelType pixValue = 0; + + for (unsigned int y = 0; y < H; ++y) + { + PixelType* buffer = reinterpret_cast<PixelType*>(image->GetRow(y)); + for (unsigned int x = 0; x < W; ++x, ++buffer, ++pixCounter) + { + // 0..99 0..99 0..55 + *buffer = pixValue; + pixValue++; + if (pixValue >= 100) + pixValue = 0; + } + } + + HistogramData hd; + ComputeHistogram(*image, hd, 10); + ASSERT_EQ(-0.5, hd.minValue); + ASSERT_EQ(10u, hd.bins.size()); + + for (size_t i = 0; i <= 4; ++i) + ASSERT_EQ(30u, hd.bins[i]); + + ASSERT_EQ(26u, hd.bins[5]); + + for (size_t i = 6; i <= 9; ++i) + ASSERT_EQ(20u, hd.bins[i]); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale8_BinSize10_2) +{ + SimpleHisto_T_BinSize10_2<Orthanc::PixelFormat_Grayscale8>(); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale16_BinSize10_2) +{ + SimpleHisto_T_BinSize10_2<Orthanc::PixelFormat_Grayscale16>(); +} + +TEST(ImageToolbox, SimpleHisto_SignedGrayscale16_BinSize10_2) +{ + SimpleHisto_T_BinSize10_2<Orthanc::PixelFormat_SignedGrayscale16>(); +} + +TEST(ImageToolbox, SimpleHisto_Grayscale32_BinSize10_2) +{ + SimpleHisto_T_BinSize10_2<Orthanc::PixelFormat_Grayscale32>(); +} +