gil/image_processing/hessian.hpp

0001 //
0002 // Copyright 2019 Olzhas Zhumabek <anonymous.from.applecity@gmail.com>
0003 // Copyright 2021 Scramjet911 <36035352+Scramjet911@users.noreply.github.com>
0004 // Use, modification and distribution are subject to the Boost Software License,
0005 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
0006 // http://www.boost.org/LICENSE_1_0.txt)
0007
0008 #ifndef BOOST_GIL_IMAGE_PROCESSING_HESSIAN_HPP
0009 #define BOOST_GIL_IMAGE_PROCESSING_HESSIAN_HPP
0010
0011 #include <boost/gil/image_view.hpp>
0012 #include <boost/gil/typedefs.hpp>
0013 #include <boost/gil/image_processing/kernel.hpp>
0014 #include <stdexcept>
0015
0016 namespace boost { namespace gil {
0017
0018 /// \brief Computes Hessian response
0019 ///
0020 /// Computes Hessian response based on computed entries of Hessian matrix, e.g. second order
0021 /// derivates in x and y, and derivatives in both x, y.
0022 /// d stands for derivative, and x or y stand for derivative direction. For example,
0023 /// ddxx means taking two derivatives (gradients) in horizontal direction.
0024 /// Weights change perception of surroinding pixels.
0025 /// Additional filtering is strongly advised.
0026 template <typename GradientView, typename T, typename Allocator, typename OutputView>
0027 inline void compute_hessian_responses(
0028     GradientView ddxx,
0029     GradientView dxdy,
0030     GradientView ddyy,
0031     const detail::kernel_2d<T, Allocator>& weights,
0032     OutputView dst)
0033 {
0034     if (ddxx.dimensions() != ddyy.dimensions()
0035         || ddyy.dimensions() != dxdy.dimensions()
0036         || dxdy.dimensions() != dst.dimensions()
0037         || weights.center_x() != weights.center_y())
0038     {
0039         throw std::invalid_argument("dimensions of views are not the same"
0040             " or weights don't have equal width and height"
0041             " or weights' dimensions are not odd");
0042     }
0043     // Use pixel type of output, as values will be written to output
0044     using pixel_t = typename std::remove_reference<decltype(std::declval<OutputView>()(0, 0))>::type;
0045
0046     using channel_t = typename std::remove_reference
0047         <
0048             decltype(std::declval<pixel_t>().at(std::integral_constant<int, 0>{}))
0049         >::type;
0050
0051
0052     auto center = weights.center_y();
0053     for (auto y = center; y < dst.height() - center; ++y)
0054     {
0055         for (auto x = center; x < dst.width() - center; ++x)
0056         {
0057             auto ddxx_i = channel_t();
0058             auto ddyy_i = channel_t();
0059             auto dxdy_i = channel_t();
0060             for (typename OutputView::coord_t w_y = 0; w_y < static_cast<std::ptrdiff_t>(weights.size()); ++w_y)
0061             {
0062                 for (typename OutputView::coord_t w_x = 0; w_x < static_cast<std::ptrdiff_t>(weights.size()); ++w_x)
0063                 {
0064                     ddxx_i += ddxx(x + w_x - center, y + w_y - center)
0065                         .at(std::integral_constant<int, 0>{}) * weights.at(w_x, w_y);
0066                     ddyy_i += ddyy(x + w_x - center, y + w_y - center)
0067                         .at(std::integral_constant<int, 0>{}) * weights.at(w_x, w_y);
0068                     dxdy_i += dxdy(x + w_x - center, y + w_y - center)
0069                         .at(std::integral_constant<int, 0>{}) * weights.at(w_x, w_y);
0070                 }
0071             }
0072             auto determinant = ddxx_i * ddyy_i - dxdy_i * dxdy_i;
0073             dst(x, y).at(std::integral_constant<int, 0>{}) = determinant;
0074         }
0075     }
0076 }
0077
0078 }} // namespace boost::gil
0079
0080 #endif