| 123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464 | //// Copyright 2019 Miral Shah <miralshah2211@gmail.com>//// Use, modification and distribution are subject to the Boost Software License,// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at// http://www.boost.org/LICENSE_1_0.txt)//#ifndef BOOST_GIL_IMAGE_PROCESSING_THRESHOLD_HPP#define BOOST_GIL_IMAGE_PROCESSING_THRESHOLD_HPP#include <limits>#include <array>#include <type_traits>#include <cstddef>#include <algorithm>#include <vector>#include <cmath>#include <boost/assert.hpp>#include <boost/gil/image.hpp>#include <boost/gil/extension/numeric/kernel.hpp>#include <boost/gil/extension/numeric/convolve.hpp>#include <boost/gil/image_processing/numeric.hpp>namespace boost { namespace gil {namespace detail {template<    typename SourceChannelT,    typename ResultChannelT,    typename SrcView,    typename DstView,    typename Operator>void threshold_impl(SrcView const& src_view, DstView const& dst_view, Operator const& threshold_op){    gil_function_requires<ImageViewConcept<SrcView>>();    gil_function_requires<MutableImageViewConcept<DstView>>();    static_assert(color_spaces_are_compatible    <        typename color_space_type<SrcView>::type,        typename color_space_type<DstView>::type    >::value, "Source and destination views must have pixels with the same color space");    //iterate over the image checking each pixel value for the threshold    for (std::ptrdiff_t y = 0; y < src_view.height(); y++)    {        typename SrcView::x_iterator src_it = src_view.row_begin(y);        typename DstView::x_iterator dst_it = dst_view.row_begin(y);        for (std::ptrdiff_t x = 0; x < src_view.width(); x++)        {            static_transform(src_it[x], dst_it[x], threshold_op);        }    }}} //namespace boost::gil::detail/// \addtogroup ImageProcessing/// @{////// \brief Direction of image segmentation./// The direction specifies which pixels are considered as corresponding to object/// and which pixels correspond to background.enum class threshold_direction{    regular, ///< Consider values greater than threshold value    inverse  ///< Consider values less than or equal to threshold value};/// \ingroup ImageProcessing/// \brief Method of optimal threshold value calculation.enum class threshold_optimal_value{    otsu        ///< \todo TODO};/// \ingroup ImageProcessing/// \brief TODOenum class threshold_truncate_mode{    threshold,  ///< \todo TODO    zero        ///< \todo TODO};enum class threshold_adaptive_method{    mean,    gaussian};/// \ingroup ImageProcessing/// \brief Applies fixed threshold to each pixel of image view./// Performs image binarization by thresholding channel value of each/// pixel of given image view./// \param src_view - TODO/// \param dst_view - TODO/// \param threshold_value - TODO/// \param max_value - TODO/// \param threshold_direction - if regular, values greater than threshold_value are/// set to max_value else set to 0; if inverse, values greater than threshold_value are/// set to 0 else set to max_value.template <typename SrcView, typename DstView>void threshold_binary(    SrcView const& src_view,    DstView const& dst_view,    typename channel_type<DstView>::type threshold_value,    typename channel_type<DstView>::type max_value,    threshold_direction direction = threshold_direction::regular){    //deciding output channel type and creating functor    using source_channel_t = typename channel_type<SrcView>::type;    using result_channel_t = typename channel_type<DstView>::type;    if (direction == threshold_direction::regular)    {        detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,            [threshold_value, max_value](source_channel_t px) -> result_channel_t {                return px > threshold_value ? max_value : 0;            });    }    else    {        detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,            [threshold_value, max_value](source_channel_t px) -> result_channel_t {                return px > threshold_value ? 0 : max_value;            });    }}/// \ingroup ImageProcessing/// \brief Applies fixed threshold to each pixel of image view./// Performs image binarization by thresholding channel value of each/// pixel of given image view./// This variant of threshold_binary automatically deduces maximum value for each channel/// of pixel based on channel type./// If direction is regular, values greater than threshold_value will be set to maximum/// numeric limit of channel else 0./// If direction is inverse, values greater than threshold_value will be set to 0 else maximum/// numeric limit of channel.template <typename SrcView, typename DstView>void threshold_binary(    SrcView const& src_view,    DstView const& dst_view,    typename channel_type<DstView>::type threshold_value,    threshold_direction direction = threshold_direction::regular){    //deciding output channel type and creating functor    using result_channel_t = typename channel_type<DstView>::type;    result_channel_t max_value = (std::numeric_limits<result_channel_t>::max)();    threshold_binary(src_view, dst_view, threshold_value, max_value, direction);}/// \ingroup ImageProcessing/// \brief Applies truncating threshold to each pixel of image view./// Takes an image view and performs truncating threshold operation on each chennel./// If mode is threshold and direction is regular:/// values greater than threshold_value will be set to threshold_value else no change/// If mode is threshold and direction is inverse:/// values less than or equal to threshold_value will be set to threshold_value else no change/// If mode is zero and direction is regular:/// values less than or equal to threshold_value will be set to 0 else no change/// If mode is zero and direction is inverse:/// values more than threshold_value will be set to 0 else no changetemplate <typename SrcView, typename DstView>void threshold_truncate(    SrcView const& src_view,    DstView const& dst_view,    typename channel_type<DstView>::type threshold_value,    threshold_truncate_mode mode = threshold_truncate_mode::threshold,    threshold_direction direction = threshold_direction::regular){    //deciding output channel type and creating functor    using source_channel_t = typename channel_type<SrcView>::type;    using result_channel_t = typename channel_type<DstView>::type;    std::function<result_channel_t(source_channel_t)> threshold_logic;    if (mode == threshold_truncate_mode::threshold)    {        if (direction == threshold_direction::regular)        {            detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,                [threshold_value](source_channel_t px) -> result_channel_t {                    return px > threshold_value ? threshold_value : px;                });        }        else        {            detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,                [threshold_value](source_channel_t px) -> result_channel_t {                    return px > threshold_value ? px : threshold_value;                });        }    }    else    {        if (direction == threshold_direction::regular)        {            detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,                [threshold_value](source_channel_t px) -> result_channel_t {                    return px > threshold_value ? px : 0;                });        }        else        {            detail::threshold_impl<source_channel_t, result_channel_t>(src_view, dst_view,                [threshold_value](source_channel_t px) -> result_channel_t {                    return px > threshold_value ? 0 : px;                });        }    }}namespace detail{template <typename SrcView, typename DstView>void otsu_impl(SrcView const& src_view, DstView const& dst_view, threshold_direction direction){    //deciding output channel type and creating functor    using source_channel_t = typename channel_type<SrcView>::type;    std::array<std::size_t, 256> histogram{};    //initial value of min is set to maximum possible value to compare histogram data    //initial value of max is set to minimum possible value to compare histogram data    auto min = (std::numeric_limits<source_channel_t>::max)(),        max = (std::numeric_limits<source_channel_t>::min)();    if (sizeof(source_channel_t) > 1 || std::is_signed<source_channel_t>::value)    {        //iterate over the image to find the min and max pixel values        for (std::ptrdiff_t y = 0; y < src_view.height(); y++)        {            typename SrcView::x_iterator src_it = src_view.row_begin(y);            for (std::ptrdiff_t x = 0; x < src_view.width(); x++)            {                if (src_it[x] < min) min = src_it[x];                if (src_it[x] > min) min = src_it[x];            }        }        //making histogram        for (std::ptrdiff_t y = 0; y < src_view.height(); y++)        {            typename SrcView::x_iterator src_it = src_view.row_begin(y);            for (std::ptrdiff_t x = 0; x < src_view.width(); x++)            {                histogram[((src_it[x] - min) * 255) / (max - min)]++;            }        }    }    else    {        //making histogram        for (std::ptrdiff_t y = 0; y < src_view.height(); y++)        {            typename SrcView::x_iterator src_it = src_view.row_begin(y);            for (std::ptrdiff_t x = 0; x < src_view.width(); x++)            {                histogram[src_it[x]]++;            }        }    }    //histData = histogram data    //sum = total (background + foreground)    //sumB = sum background    //wB = weight background    //wf = weight foreground    //varMax = tracking the maximum known value of between class variance    //mB = mu background    //mF = mu foreground    //varBeetween = between class variance    //http://www.labbookpages.co.uk/software/imgProc/otsuThreshold.html    //https://www.ipol.im/pub/art/2016/158/    std::ptrdiff_t total_pixel = src_view.height() * src_view.width();    std::ptrdiff_t sum_total = 0, sum_back = 0;    std::size_t weight_back = 0, weight_fore = 0, threshold = 0;    double var_max = 0, mean_back, mean_fore, var_intra_class;    for (std::size_t t = 0; t < 256; t++)    {        sum_total += t * histogram[t];    }    for (int t = 0; t < 256; t++)    {        weight_back += histogram[t];               // Weight Background        if (weight_back == 0) continue;        weight_fore = total_pixel - weight_back;          // Weight Foreground        if (weight_fore == 0) break;        sum_back += t * histogram[t];        mean_back = sum_back / weight_back;            // Mean Background        mean_fore = (sum_total - sum_back) / weight_fore;    // Mean Foreground        // Calculate Between Class Variance        var_intra_class = weight_back * weight_fore * (mean_back - mean_fore) * (mean_back - mean_fore);        // Check if new maximum found        if (var_intra_class > var_max) {            var_max = var_intra_class;            threshold = t;        }    }    if (sizeof(source_channel_t) > 1 && std::is_unsigned<source_channel_t>::value)    {        threshold_binary(src_view, dst_view, (threshold * (max - min) / 255) + min, direction);    }    else {        threshold_binary(src_view, dst_view, threshold, direction);    }}} //namespace detailtemplate <typename SrcView, typename DstView>void threshold_optimal(    SrcView const& src_view,    DstView const& dst_view,    threshold_optimal_value mode = threshold_optimal_value::otsu,    threshold_direction direction = threshold_direction::regular){    if (mode == threshold_optimal_value::otsu)    {        for (std::size_t i = 0; i < src_view.num_channels(); i++)        {            detail::otsu_impl                (nth_channel_view(src_view, i), nth_channel_view(dst_view, i), direction);        }    }}namespace detail {template<    typename SourceChannelT,    typename ResultChannelT,    typename SrcView,    typename DstView,    typename Operator>void adaptive_impl(    SrcView const& src_view,    SrcView const& convolved_view,    DstView const& dst_view,    Operator const& threshold_op){    //template argument validation    gil_function_requires<ImageViewConcept<SrcView>>();    gil_function_requires<MutableImageViewConcept<DstView>>();    static_assert(color_spaces_are_compatible    <        typename color_space_type<SrcView>::type,        typename color_space_type<DstView>::type    >::value, "Source and destination views must have pixels with the same color space");    //iterate over the image checking each pixel value for the threshold    for (std::ptrdiff_t y = 0; y < src_view.height(); y++)    {        typename SrcView::x_iterator src_it = src_view.row_begin(y);        typename SrcView::x_iterator convolved_it = convolved_view.row_begin(y);        typename DstView::x_iterator dst_it = dst_view.row_begin(y);        for (std::ptrdiff_t x = 0; x < src_view.width(); x++)        {            static_transform(src_it[x], convolved_it[x], dst_it[x], threshold_op);        }    }}} //namespace boost::gil::detailtemplate <typename SrcView, typename DstView>void threshold_adaptive(    SrcView const& src_view,    DstView const& dst_view,    typename channel_type<DstView>::type max_value,    std::size_t kernel_size,    threshold_adaptive_method method = threshold_adaptive_method::mean,    threshold_direction direction = threshold_direction::regular,    typename channel_type<DstView>::type constant = 0){    BOOST_ASSERT_MSG((kernel_size % 2 != 0), "Kernel size must be an odd number");    typedef typename channel_type<SrcView>::type source_channel_t;    typedef typename channel_type<DstView>::type result_channel_t;    image<typename SrcView::value_type> temp_img(src_view.width(), src_view.height());    typename image<typename SrcView::value_type>::view_t temp_view = view(temp_img);    SrcView temp_conv(temp_view);    if (method == threshold_adaptive_method::mean)    {        std::vector<float> mean_kernel_values(kernel_size, 1.0f/kernel_size);        kernel_1d<float> kernel(mean_kernel_values.begin(), kernel_size, kernel_size/2);        detail::convolve_1d        <            pixel<float, typename SrcView::value_type::layout_t>        >(src_view, kernel, temp_view);    }    else if (method == threshold_adaptive_method::gaussian)    {        detail::kernel_2d<float> kernel = generate_gaussian_kernel(kernel_size, 1.0);        convolve_2d(src_view, kernel, temp_view);    }    if (direction == threshold_direction::regular)    {        detail::adaptive_impl<source_channel_t, result_channel_t>(src_view, temp_conv, dst_view,            [max_value, constant](source_channel_t px, source_channel_t threshold) -> result_channel_t        { return px > (threshold - constant) ? max_value : 0; });    }    else    {        detail::adaptive_impl<source_channel_t, result_channel_t>(src_view, temp_conv, dst_view,            [max_value, constant](source_channel_t px, source_channel_t threshold) -> result_channel_t        { return px > (threshold - constant) ? 0 : max_value; });    }}template <typename SrcView, typename DstView>void threshold_adaptive(    SrcView const& src_view,    DstView const& dst_view,    std::size_t kernel_size,    threshold_adaptive_method method = threshold_adaptive_method::mean,    threshold_direction direction = threshold_direction::regular,    int constant = 0){    //deciding output channel type and creating functor    typedef typename channel_type<DstView>::type result_channel_t;    result_channel_t max_value = (std::numeric_limits<result_channel_t>::max)();    threshold_adaptive(src_view, dst_view, max_value, kernel_size, method, direction, constant);}/// @}}} //namespace boost::gil#endif //BOOST_GIL_IMAGE_PROCESSING_THRESHOLD_HPP
 |