Commit ae1dc3b3 authored by Edwin Carlinet's avatar Edwin Carlinet
Browse files

Merge branch 'development/cdt' into 'development/ranges'

Migrate chamfer distance transform.

See merge request !77
parents 775acc37 08a1ca2c
Pipeline #13830 failed with stages
in 94 minutes and 41 seconds
......@@ -6,6 +6,7 @@
#include <mln/core/se/rect2d.hpp>
#include <mln/core/se/mask2d.hpp>
#include <mln/core/neighborhood/c4.hpp>
#include <mln/core/neighborhood/c8.hpp>
#include <mln/io/experimental/imread.hpp>
#include <mln/morpho/experimental/closing.hpp>
......@@ -18,6 +19,7 @@
#include <mln/morpho/experimental/watershed.hpp>
#include <mln/labeling/experimental/local_extrema.hpp>
#include <mln/labeling/experimental/chamfer_distance_transform.hpp>
#include <benchmark/benchmark.h>
......@@ -216,6 +218,20 @@ BENCHMARK_F(BMMorpho, minima)(benchmark::State& st)
this->run(st, f);
}
BENCHMARK_F(BMMorpho, cdt_2_3)(benchmark::State& st)
{
// mln::se::experimental::wmask2d weights = {{3, 2, 3}, //
// {2, 0, 2},
// {3, 2, 3}};
auto f = [](const image_t& input, image_t&) {
auto tmp = mln::transform(input, [](uint8_t x) { return x < 128; });
auto out = mln::labeling::experimental::chamfer_distance_transform<int16_t>(tmp, mln::experimental::c8);
return out;
};
this->run(st, f);
}
BENCHMARK_F(BMMorpho, watershed)(benchmark::State& st)
{
......
......@@ -18,6 +18,7 @@
#include <mln/morpho/watershed.hpp>
#include <mln/labeling/local_extrema.hpp>
#include <mln/transform/chamfer_distance_transform.hpp>
// [legacy]
#include <mln/core/image/image2d.hpp>
......@@ -163,6 +164,18 @@ BENCHMARK_F(BMMorpho, minima)(benchmark::State& st)
this->run(st, f);
}
BENCHMARK_F(BMMorpho, cdt_2_3)(benchmark::State& st)
{
auto f = [](const image_t& input, image_t&) {
auto tmp = mln::transform(input, [](uint8_t x) { return x < 128; });
auto out = mln::transforms::chamfer_distance_transform<int16_t>(tmp, mln::c8);
return out;
};
this->run(st, f);
}
BENCHMARK_F(BMMorpho, watershed)(benchmark::State& st)
{
auto f = [](const image_t& input, image_t&) {
......
......@@ -86,24 +86,29 @@ class BMReferenceLinear : public benchmark::Fixture
virtual void SetUp(const benchmark::State&) override
{
mln::experimental::image2d<mln::rgb8> ima;
mln::io::experimental::imread(filename, ima);
if (!m_initialized)
{
mln::experimental::image2d<mln::rgb8> ima;
mln::io::experimental::imread(filename, ima);
mln::resize(m_exp_input, ima);
mln::resize(m_exp_output, ima);
mln::transform(ima, m_exp_input, [](mln::rgb8 x) -> uint8_t { return x[0]; });
mln::resize(m_exp_input, ima);
mln::resize(m_exp_output, ima);
mln::transform(ima, m_exp_input, [](mln::rgb8 x) -> uint8_t { return x[0]; });
m_lut.resize(256);
for (int i = 0; i < 256; i++)
m_lut[i] = i;
m_lut.resize(256);
for (int i = 0; i < 256; i++)
m_lut[i] = i;
m_exp_input.to(m_input, true);
mln::resize(m_output, m_input);
m_initialized = true;
}
m_height = m_exp_input.height();
m_width = m_exp_input.width();
m_stride = m_exp_input.strides();
m_exp_input.to(m_input, true);
mln::resize(m_output, m_input);
}
protected:
......@@ -144,13 +149,21 @@ protected:
int m_height;
std::ptrdiff_t m_stride;
mln::image2d<uint8_t> m_input;
mln::image2d<uint8_t> m_output;
mln::experimental::image2d<uint8_t> m_exp_input;
mln::experimental::image2d<uint8_t> m_exp_output;
std::vector<uint8_t> m_lut;
static bool m_initialized;
static mln::image2d<uint8_t> m_input;
static mln::image2d<uint8_t> m_output;
static mln::experimental::image2d<uint8_t> m_exp_input;
static mln::experimental::image2d<uint8_t> m_exp_output;
static std::vector<uint8_t> m_lut;
};
bool BMReferenceLinear::m_initialized = false;
mln::image2d<uint8_t> BMReferenceLinear::m_input;
mln::image2d<uint8_t> BMReferenceLinear::m_output;
mln::experimental::image2d<uint8_t> BMReferenceLinear::m_exp_input;
mln::experimental::image2d<uint8_t> BMReferenceLinear::m_exp_output;
std::vector<uint8_t> BMReferenceLinear::m_lut;
BENCHMARK_F(BMReferenceLinear, Mult_C)(benchmark::State& st)
{
......
......@@ -15,19 +15,21 @@ class BMReferenceNeighborhood : public benchmark::Fixture
virtual void SetUp(const benchmark::State&) override
{
mln::experimental::image2d<mln::rgb8> in;
mln::io::experimental::imread(filename, in);
mln::resize(m_exp_input, in);
mln::resize(m_exp_output, in);
mln::transform(in, m_exp_input, [](mln::rgb8 x) -> uint8_t { return x[0]; });
m_exp_input.extension().fill(0);
m_exp_input.to(m_input, true);
mln::resize(m_output, m_input);
m_input.extension().fill(0);
if (!m_initialized)
{
mln::experimental::image2d<mln::rgb8> in;
mln::io::experimental::imread(filename, in);
mln::resize(m_exp_input, in);
mln::resize(m_exp_output, in);
mln::transform(in, m_exp_input, [](mln::rgb8 x) -> uint8_t { return x[0]; });
m_exp_input.extension().fill(0);
m_exp_input.to(m_input, true);
mln::resize(m_output, m_input);
m_input.extension().fill(0);
m_initialized = true;
}
m_h = m_exp_input.height();
m_w = m_exp_input.width();
m_stride = m_exp_input.stride();
......@@ -58,12 +60,19 @@ private:
std::ptrdiff_t m_stride;
mln::uint8* m_ibuffer;
mln::uint8* m_obuffer;
mln::image2d<mln::uint8> m_input;
mln::image2d<mln::uint8> m_output;
mln::experimental::image2d<mln::uint8> m_exp_input;
mln::experimental::image2d<mln::uint8> m_exp_output;
static bool m_initialized;
static mln::image2d<mln::uint8> m_input;
static mln::image2d<mln::uint8> m_output;
static mln::experimental::image2d<mln::uint8> m_exp_input;
static mln::experimental::image2d<mln::uint8> m_exp_output;
};
bool BMReferenceNeighborhood::m_initialized = false;
mln::image2d<mln::uint8> BMReferenceNeighborhood::m_input;
mln::image2d<mln::uint8> BMReferenceNeighborhood::m_output;
mln::experimental::image2d<mln::uint8> BMReferenceNeighborhood::m_exp_input;
mln::experimental::image2d<mln::uint8> BMReferenceNeighborhood::m_exp_output;
BENCHMARK_F(BMReferenceNeighborhood, Sum_C)(benchmark::State& st)
{
runit(st, Sum_C);
......
......@@ -5,3 +5,5 @@ Labeling Module
:maxdepth: 1
labeling/local_extrema
labeling/cdt
Chamfer Distance Transform
==========================
Include :file:`<mln/labeling/chamfer_distance_transform.hpp>`
.. cpp:namespace:: mln::labeling
.. cpp:function:: template <typename DistanceType = int, Image I, WeightedNeighborhood N> \
image_ch_value_t<InputImage, DistanceType> \
chamfer_distance_transform(I input, N nbh, bool background_is_object = false)
Perform the distance transform [1]_ using the given Neighborhood and its weights. Each foreground pixel
gets the distance to the closest background pixel.
Note that the sum saturates to the maximum value that can be taken by `DistanceType`.
:tparam DistanceType: The type used for the accumulation
:param input: Input image
:param nbh: Neighborhood (weighted)
:param background_is_object (optional): Whether to consider the domain "outside" is the foreground (object) or the background.
:return: A distance image
:exception: N/A
Notes
-----
Complexity
----------
Linear in pixels (2 passes).
Example
-------
Example with the chessboard distance::
#include <mln/labeling/chamfer_distance_transform.hpp>
#include <mln/core/neighborhood/c8.hpp>
auto out = mln::labeling::chamfer_distance_transform<int16_t>(input, mln::c8, true);
Example with the 5-7-11 weights::
#include <mln/labeling/chamfer_distance_transform.hpp>
#include <mln/core/se/mask2d.hpp>
auto weights = mln::se::experimental::wmask2d({{+0, 11, +0, 11, +0}, //
{11, +7, +5, +7, 11},
{+0, +5, +0, +5, +0},
{11, +7, +5, +7, 11},
{+0, 11, +0, 11, +0}});
auto out = mln::labeling::chamfer_distance_transform<int16_t>(input, weights, true);
.. list-table::
* - .. figure:: /images/F.png
Original image
- .. figure:: /images/F-4.png
CDT with :cpp:ref:`c4` (Manhattan distance)
- .. figure:: /images/F-8.png
CDT with :cpp:ref:`c8` (Chessboard distance)
* - .. figure:: /images/F-2-3.png
CDT with 2-3 :cpp:ref:`wmask2d` (2-3 distance)
- .. figure:: /images/F-5-7-11.png
CDT with 5-7-11 :cpp:ref:`wmask2d` (5-7-11 distance)
-
References
----------
.. [1] Borgefors, G. (1986). Distance transformations in digital images. Computer vision, graphics, and image processing, 34(3), 344-371.
\ No newline at end of file
......@@ -27,6 +27,12 @@ add_image("erosion-cli;dilation;square;21" "${PYLENE_IMAGE_DIR}/lena.pgm" morpho
add_image("erosion-cli;opening;square;21" "${PYLENE_IMAGE_DIR}/lena.pgm" morpho_opening_1.png)
add_image("erosion-cli;closing;square;21" "${PYLENE_IMAGE_DIR}/lena.pgm" morpho_closing_1.png)
add_image("erosion-cli;median;square;21" "${PYLENE_IMAGE_DIR}/lena.pgm" morpho_median_1.png)
add_image("cdt;1" "${DOCUMENTATION_IMAGE_DIR}/F.png" F-4.png)
add_image("cdt;2" "${DOCUMENTATION_IMAGE_DIR}/F.png" F-8.png)
add_image("cdt;3" "${DOCUMENTATION_IMAGE_DIR}/F.png" F-2-3.png)
add_image("cdt;4" "${DOCUMENTATION_IMAGE_DIR}/F.png" F-5-7-11.png)
add_image("staff_lines" "${DOCUMENTATION_IMAGE_DIR}/staff_lines.pbm"
morpho_hitormiss_1.png morpho_hitormiss_2.png staff_lines_markers.png morpho_reconstruction_1.png morpho_reconstruction_2.png)
add_image("reconstruction"
......@@ -49,3 +55,4 @@ add_executable(staff_lines staff_lines.cpp)
add_executable(component_tree_1 component_tree_1.cpp)
add_executable(blobs_watershed blobs_watershed.cpp)
add_executable(reconstruction reconstruction.cpp)
add_executable(cdt cdt.cpp)
#include <mln/core/algorithm/transform.hpp>
#include <mln/core/image/image2d.hpp>
#include <mln/core/image/experimental/ndimage.hpp>
#include <mln/core/image/view/operators.hpp>
#include <mln/core/neighb2d.hpp>
#include <mln/core/neighborhood/dyn_wneighborhood.hpp>
#include <mln/core/neighborhood/c8.hpp>
#include <mln/core/se/disc.hpp>
#include <mln/data/stretch.hpp>
#include <mln/io/imread.hpp>
#include <mln/io/imsave.hpp>
#include <mln/morpho/structural/dilate.hpp>
#include <mln/morpho/watershed.hpp>
#include <mln/transform/chamfer_distance_transform.hpp>
class c23_t : public mln::dyn_wneighborhood_base<std::array<mln::point2d, 8>, std::array<int, 8>,
mln::constant_neighborhood_tag, c23_t>
{
public:
static constexpr std::array<mln::point2d, 8> dpoints = {
{{-1, -1},
{-1, +0},
{-1, 1}, //
{+0, -1},
{+0, 1}, //
{+1, -1},
{+1, +0},
{+1, 1}} //
};
static constexpr std::array<int, 8> weights = {{3, 2, 3, 2, 2, 3, 2, 3}};
};
#include <mln/core/se/mask2d.hpp>
#include <mln/core/colors.hpp>
#include <mln/data/experimental/stretch.hpp>
#include <mln/io/experimental/imread.hpp>
#include <mln/io/experimental/imsave.hpp>
#include <mln/labeling/experimental/chamfer_distance_transform.hpp>
#include <mln/morpho/experimental/dilation.hpp>
#include <mln/morpho/experimental/watershed.hpp>
#include <cstdint>
constexpr std::array<mln::point2d, 8> c23_t::dpoints;
constexpr std::array<int, 8> c23_t::weights;
static const std::vector<mln::rgb8> regions_table = {
{0, 0, 0}, {255, 255, 255}, {0, 127, 255}, {127, 255, 0}, {255, 0, 127}, {148, 148, 148},
......@@ -123,29 +107,33 @@ int main(int argc, char** argv)
std::abort();
}
mln::image2d<mln::uint8> input;
mln::io::imread(argv[1], input);
mln::experimental::image2d<uint8_t> input;
mln::io::experimental::imread(argv[1], input);
// BEGIN
// (1) Compute the distance transform
auto d = mln::transforms::chamfer_distance_transform<mln::uint8>(input, c23_t());
mln::se::experimental::wmask2d weights = {{3, 2, 3}, //
{2, 0, 2},
{3, 2, 3}};
auto d = mln::labeling::experimental::chamfer_distance_transform<uint8_t>(input, weights);
// Remove non-meaninfull extrema
d = mln::morpho::structural::dilate(d, mln::se::disc(5));
d = mln::morpho::experimental::dilation(d, mln::experimental::se::disc(5));
// (2) Inverse the distance
auto dinv = mln::transform(d, [](mln::uint8 x) -> mln::uint8 { return mln::value_traits<mln::uint8>::max() - x; });
auto dinv = mln::transform(d, [](uint8_t x) -> uint8_t { return UINT8_MAX - x; });
// (3) Run the watershed segmentation
int nlabel;
auto ws = mln::morpho::experimental::watershed<mln::int16>(dinv, mln::c8, nlabel);
auto ws = mln::morpho::experimental::watershed<int16_t>(dinv, mln::experimental::c8, nlabel);
// (4) Labelize input
auto output = mln::view::ifelse(ws == 0, 1, mln::view::ifelse(input, ws, 0));
// END
auto d_stretched = mln::data::stretch<float>(d);
auto d_stretched = mln::data::experimental::stretch<float>(d);
mln::io::imsave(mln::transform(d_stretched, heat_lut), argv[2]);
mln::io::imsave(mln::transform(output, regions_lut), argv[3]);
mln::io::experimental::imsave(mln::transform(d_stretched, heat_lut), argv[2]);
mln::io::experimental::imsave(mln::transform(output, regions_lut), argv[3]);
}
#include <mln/core/image/experimental/ndimage.hpp>
#include <mln/core/se/mask2d.hpp>
#include <mln/core/colors.hpp>
#include <mln/data/experimental/stretch.hpp>
#include <mln/io/experimental/imread.hpp>
#include <mln/io/experimental/imsave.hpp>
#include <mln/labeling/experimental/chamfer_distance_transform.hpp>
#include <cstdint>
mln::rgb8 heat_lut(float x)
{
assert(0 <= x && x <= 1);
float x0 = 1.f / 4.f;
float x1 = 2.f / 4.f;
float x2 = 3.f / 4.f;
if (x < x0)
{
auto g = static_cast<mln::uint8>(x / x0 * 255);
return mln::rgb8{0, g, 255};
}
else if (x < x1)
{
auto b = static_cast<mln::uint8>((x1 - x) / x0 * 255);
return mln::rgb8{0, 255, b};
}
else if (x < x2)
{
auto r = static_cast<mln::uint8>((x - x1) / x0 * 255);
return mln::rgb8{r, 255, 0};
}
else
{
auto b = static_cast<mln::uint8>((1.f - x) / x0 * 255);
return mln::rgb8{255, b, 0};
}
}
int main(int argc, char** argv)
{
if (argc < 3)
{
std::cout << "Usage:" << argv[0] << " distance_type input output\n"
"chamfer distance transform with disttype:\n"
"1: 4-connexion\n"
"2: 8-connexion\n"
"3: 2-3 weights\n"
"4: 5-7-11 weights\n";
std::abort();
}
mln::experimental::image2d<bool> input;
mln::io::experimental::imread(argv[2], input);
// BEGIN
// (1) Compute the distance transform
mln::se::experimental::wmask2d<int> weights;
switch (std::atoi(argv[1]))
{
case 1:
weights = mln::se::experimental::wmask2d({{0, 1, 0}, //
{1, 0, 1},
{0, 1, 0}});
break;
case 2:
weights = mln::se::experimental::wmask2d({{1, 1, 1}, //
{1, 0, 1},
{1, 1, 1}});
break;
case 3:
weights = mln::se::experimental::wmask2d({{3, 2, 3}, //
{2, 0, 2},
{3, 2, 3}});
break;
case 4:
weights = mln::se::experimental::wmask2d({{+0, 11, +0, 11, +0}, //
{11, +7, +5, +7, 11},
{+0, +5, +0, +5, +0},
{11, +7, +5, +7, 11},
{+0, 11, +0, 11, +0}});
break;
default:
std::cerr << "Bad type.\n";
std::abort();
}
auto d = mln::labeling::experimental::chamfer_distance_transform<int16_t>(input, weights, true);
auto d_stretched = mln::data::experimental::stretch<float>(d);
mln::io::experimental::imsave(mln::transform(d_stretched, heat_lut), argv[3]);
}
......@@ -47,6 +47,7 @@ namespace mln
void transform(InputImage in, OutputImage out, UnaryFunction f)
{
static_assert(mln::is_a<InputImage, experimental::Image>());
static_assert(mln::is_a<OutputImage, experimental::Image>());
static_assert(::ranges::Invocable<UnaryFunction, image_reference_t<InputImage>>());
static_assert(std::is_convertible_v<std::invoke_result_t<UnaryFunction, image_reference_t<InputImage>>,
image_value_t<OutputImage>>,
......
......@@ -18,7 +18,6 @@ namespace mln::concepts
// clang-format off
#ifdef PYLENE_CONCEPT_TS_ENABLED
template <typename SE, typename P>
concept Neighborhood =
StructuringElement<SE, P> &&
......@@ -30,10 +29,9 @@ namespace mln::concepts
requires detail::RangeValueTypeConvertibleTo<decltype(se.before(p)), P>;
requires detail::RangeValueTypeConvertibleTo<decltype(se.after(p)), P>;
requires detail::RangeValueTypeConvertibleTo<decltype(se.before(px)), mln::archetypes::PixelT<P>>;
requires detail::RangeValueTypeConvertibleTo<decltype(se.after(px)), mln::archetypes::PixelT<P>>;
requires concepts::Pixel<::ranges::range_value_t<decltype(se.before(px))>>;
requires concepts::Pixel<::ranges::range_value_t<decltype(se.after(px))>>;
};
#endif
// clang-format on
......
......@@ -57,7 +57,7 @@ namespace mln::concepts
{ cse.offsets() } -> stl::ForwardRange&&;
requires detail::RangeValueTypeConvertibleTo<decltype(se(p)), P>;
requires detail::RangeValueTypeConvertibleTo<decltype(se(px)), mln::archetypes::PixelT<P>>;
requires concepts::Pixel<::ranges::range_value_t<decltype(se(px))>>;
requires detail::RangeValueTypeConvertibleTo<decltype(cse.offsets()), P>;
};
......
......@@ -4,42 +4,57 @@
#include <mln/core/image/image.hpp>
namespace mln
namespace mln::extension
{
namespace extension
{
template <typename I>
const I& fill(const Image<I>& ima, mln_value(I) v);
template <typename I>
const I& fill(const Image<I>& ima, mln_value(I) v);
template <typename I>
I&& fill(Image<I>&& ima, mln_value(I) v);
template <typename I>
I&& fill(Image<I>&& ima, mln_value(I) v);
template <class InputImage>
bool try_fill(const InputImage& f, image_value_t<InputImage> v);
/******************************************/
/**** Implementation ****/
/******************************************/
/******************************************/
/**** Implementation ****/
/******************************************/
template <typename I>
const I& fill(const Image<I>& ima, mln_value(I) v)
{
static_assert(image_has_extension<I>::value, "Image must have an extension.");
static_assert(extension_traits<typename image_extension_type<I>::type>::support_fill::value,
"Image extension must support filling.");
template <typename I>
const I& fill(const Image<I>& ima, mln_value(I) v)
{
static_assert(image_has_extension<I>::value, "Image must have an extension.");
static_assert(extension_traits<typename image_extension_type<I>::type>::support_fill::value,
"Image extension must support filling.");
exact(ima).extension().fill(v);
return exact(ima);
}
exact(ima).extension().fill(v);
return exact(ima);
}
</