Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Bicubic sampler #588

Open
wants to merge 3 commits into
base: develop
Choose a base branch
from
Open

Bicubic sampler #588

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
793cd33
Bicubic sampler with bugs
Scramjet911 Mar 13, 2021
7979ec6
Cleaned bicubic sampler with bugs
Scramjet911 Mar 13, 2021
b81d3cc
Bicubic Sampler and static_for_each with 5 sources
Scramjet911 Mar 30, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
16 changes: 16 additions & 0 deletions include/boost/gil/color_base_algorithm.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -363,6 +363,13 @@ struct element_recursion
op2(semantic_at_c<N-1>(p1), semantic_at_c<N-1>(p2), semantic_at_c<N-1>(p3));
return op2;
}
//static_for_each with 5 sources
template <typename P1,typename P2,typename P3,typename P4,typename P5,typename Op>
static Op static_for_each(const P1& p1, const P2& p2, const P3& p3, const P4& p4, P5& p5, Op op) {
Op op2(element_recursion<N-1>::static_for_each(p1,p2,p3,p4,p5,op));
op2(semantic_at_c<N-1>(p1), semantic_at_c<N-1>(p2), semantic_at_c<N-1>(p3), semantic_at_c<N-1>(p4), semantic_at_c<N-1>(p5));
return op2;
}
//static_transform with one source
template <typename P1,typename Dst,typename Op>
static Op static_transform(P1& src, Dst& dst, Op op) {
Expand Down Expand Up @@ -401,6 +408,8 @@ struct element_recursion
semantic_at_c<N-1>(dst)=op2(semantic_at_c<N-1>(src1), semantic_at_c<N-1>(src2));
return op2;
}


};

// Termination condition of the compile-time recursion for element operations on a color base
Expand All @@ -426,12 +435,16 @@ template<> struct element_recursion<0> {
//static_for_each with three sources
template <typename P1,typename P2,typename P3,typename Op>
static Op static_for_each(const P1&,const P2&,const P3&,Op op){return op;}
//static_for_each with 5 sources
template <typename P1,typename P2,typename P3,typename P4,typename P5,typename Op>
static Op static_for_each(const P1&,const P2&,const P3&,const P4&,const P5&,Op op){return op;}
//static_transform with one source
template <typename P1,typename Dst,typename Op>
static Op static_transform(const P1&,const Dst&,Op op){return op;}
//static_transform with two sources
template <typename P1,typename P2,typename Dst,typename Op>
static Op static_transform(const P1&,const P2&,const Dst&,Op op){return op;}

};

// std::min and std::max don't have the mutable overloads...
Expand Down Expand Up @@ -706,6 +719,9 @@ Op static_for_each(const P1& p1,const P2& p2,P3& p3,Op op) { return detail::elem
template <typename P1,typename P2,typename P3,typename Op>
BOOST_FORCEINLINE
Op static_for_each(const P1& p1,const P2& p2,const P3& p3,Op op) { return detail::element_recursion<size<P1>::value>::static_for_each(p1,p2,p3,op); }
template <typename P1,typename P2,typename P3,typename P4,typename P5,typename Op>
BOOST_FORCEINLINE
Op static_for_each(const P1& p1,const P2& p2,const P3& p3,const P4& p4,P5& p5,Op op) { return detail::element_recursion<size<P1>::value>::static_for_each(p1,p2,p3,p4,p5,op); }
///\}

} } // namespace boost::gil
Expand Down
315 changes: 314 additions & 1 deletion include/boost/gil/extension/numeric/sampler.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@

namespace boost { namespace gil {

// Nearest-neighbor and bilinear image samplers.
// Nearest-neighbor, bilinear and bicubic image samplers.
// NOTE: The code is for example use only. It is not optimized for performance

///////////////////////////////////////////////////////////////////////////
Expand Down Expand Up @@ -55,11 +55,18 @@ struct cast_channel_fn {
}
};



template <typename SrcPixel, typename DstPixel>
void cast_pixel(const SrcPixel& src, DstPixel& dst) {
static_for_each(src,dst,cast_channel_fn());
}

// template <typename SrcPixel, typename DstPixel>
// void cast_pixel_bicubic(const SrcPixel& src, DstPixel& dst) {
// static_for_each(src, dst, cubic_interpolate_check());
// }

namespace detail {

template <typename Weight>
Expand All @@ -85,6 +92,56 @@ struct add_dst_mul_src {
// dst);
}
};

struct cubic_interpolate_check {
template <typename SrcChannel, typename DstChannel>
void operator()(const SrcChannel& , DstChannel& dst) {
// using dst_value_t = typename channel_traits<DstChannel>::value_type;
using src_value_t = typename channel_traits<SrcChannel>::value_type;

if(dst > channel_traits<SrcChannel>::max_value()){
dst = src_value_t(channel_traits<SrcChannel>::max_value());
}
else if(dst < channel_traits<SrcChannel>::min_value()){
dst = src_value_t(channel_traits<SrcChannel>::min_value());
}
else{
dst = src_value_t(dst);
}
}
};

template <typename Weight>
struct cubic_interpolate_channel {
Weight _w1;
Weight _w2;
Weight _w3;
Weight _w4;
cubic_interpolate_channel(Weight w1, Weight w2, Weight w3, Weight w4) : _w1(w1), _w2(w2), _w3(w3), _w4(w4) {}

template <typename SrcChannel, typename DstChannel>
void operator()(const SrcChannel& src1, const SrcChannel& src2, const SrcChannel& src3, const SrcChannel& src4, DstChannel& dst) const {
dst = DstChannel(src1*_w1 + src2*_w2 + src3*_w3 + src4*_w4);
}
};

template <typename SrcView, typename SrcP, typename Weight, typename DstP>
struct cubic_interpolate
{
void operator()(const SrcP& p0, const SrcP& p1, const SrcP& p2, const SrcP& p3, const Weight x, DstP& dst) const
{
Weight w0 = ((-0.5 * x) + (x*x) - (0.5 * x*x*x));
Weight w1 = (1 - (2.5 * x*x) + (1.5 * x*x*x));
Weight w2 = ((0.5 * x) + (2 * x*x) - (1.5 * x*x*x));
Weight w3 = ((-0.5 * x*x) + (0.5 * x*x*x));

static_for_each(p0, p1, p2, p3, dst, cubic_interpolate_channel<Weight>(w0, w1, w2, w3));

// Overflow check. Should optimally be after all cubic calculations
static_for_each(p0, dst, cubic_interpolate_check());
}
};

} // namespace detail

/// \brief A sampler that sets the destination pixel as the bilinear interpolation of the four closest pixels from the source.
Expand Down Expand Up @@ -183,6 +240,262 @@ bool sample(bilinear_sampler, SrcView const& src, point<F> const& p, DstP& resul
return true;
}

/// \brief A sampler that sets the destination pixel as the bicubic interpolation of the 16 closest pixels from the source.
/// If outside the bounds, it doesn't change the destination. Boundary pixels are the cubic interpolation of 4 adjacent
/// pixels.
/// \ingroup ImageAlgorithms
struct bicubic_sampler {};

template <typename DstP, typename SrcView, typename F>
bool sample(bicubic_sampler, SrcView const& src, point<F> const& p, DstP& result)
{
using SrcP = typename SrcView::value_type;
using cubic_interpolate_current = typename detail::cubic_interpolate<SrcView, SrcP, F, pixel<F, devicen_layout_t<num_channels<SrcView>::value> > >;

point_t p0(ifloor(p.x), ifloor(p.y)); // the closest integer coordinate top left from p
point<F> frac(p.x-p0.x, p.y-p0.y);

if (p0.x < -1 || p0.y < -1 || p0.x>=src.width() || p0.y>=src.height())
{
return false;
}

pixel<F,devicen_layout_t<num_channels<SrcView>::value> > mp0(0), mp1(0), mp2(0), mp3(0), mpres(0); // suboptimal

typename SrcView::xy_locator loc=src.xy_at(p0.x,p0.y);

if(p0.y == -1)
{
// Topmost row pixels
++loc.y();
if(p0.x == -1)
{
// Top left corner
detail::add_dst_mul_src<SrcP,F,pixel<F,devicen_layout_t<num_channels<SrcView>::value> > >()(loc.x()[1], 1, mpres);
}
else if(p0.x == 0)
{
// Top left after corner pixel
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mpres);
}
else if(p0.x+1 == src.width())
{
// Top right pixel
detail::add_dst_mul_src<SrcP,F,pixel<F,devicen_layout_t<num_channels<SrcView>::value> > >()(*loc, 1, mpres);
}
else if(p0.x+2 == src.width())
{
// Top right pixel before corner pixel
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mpres);
}
else
{
// All other top row pixels
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mpres);
}
}
else if(p0.y == 0)
{
// Row after topmost row, 2nd row
if(p0.x == -1)
{
// Left corner pixel of 2nd row
++loc.x();
cubic_interpolate_current()(*loc, *loc, *(loc.y()+1), *(loc.y()+2), frac.y, mpres);
}
else if(p0.x == 0)
{
// Pixel after left corner pixel of 2nd row
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
mp0 = mp1;
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else if(p0.x+1 == src.width())
{
// Right corner pixel of 2nd row
cubic_interpolate_current()(*loc, *loc, *(loc.y()+1), *(loc.y()+2), frac.y, mpres);
}
else if(p0.x+2 == src.width())
{
// Pixel before right corner pixel of 2nd row
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp1);
mp0 = mp1;
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else
{
// All other pixels of 2nd row
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
mp0 = mp1;
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}

}
else if(p0.y+1 == src.height())
{
// Bottom most row
if(p0.x == -1)
{
// Left corner pixel of last row
detail::add_dst_mul_src<SrcP,F,pixel<F,devicen_layout_t<num_channels<SrcView>::value> > >()(loc.x()[1], 1, mpres);
}
else if(p0.x == 0)
{
// Pixel after left corner of last row
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mpres);
}
else if(p0.x+1 == src.width())
{
// Right corner of last row
detail::add_dst_mul_src<SrcP,F,pixel<F,devicen_layout_t<num_channels<SrcView>::value> > >()(*loc, 1, mpres);
}
else if(p0.x+2 == src.width())
{
// Pixel before right corner of last row
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mpres);
}
else
{
// All other pixels of last row
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mpres);
}
}
else if(p0.y+2 == src.height())
{
// Row before bottom most row
if(p0.x == -1)
{
// Leftmost corner of 2nd last row
++loc.x();
cubic_interpolate_current()(*(loc.y()-1), *loc, *(loc.y()+1), *(loc.y()+1), frac.y, mpres);
}
else if(p0.x == 0)
{
// Pixel after left corner of second last row
--loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
mp3 = mp2;

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else if(p0.x+1 == src.width())
{
// Right corner of second last row
cubic_interpolate_current()(*(loc.y()-1), *loc, *(loc.y()+1), *(loc.y()+1), frac.y, mpres);
}
else if(p0.x+2 == src.width())
{
// Pixel before right corner of second last row
--loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp2);
mp3 = mp2;

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else
{
// All other pixels in second last row
--loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
mp3 = mp2;

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
}
else
{
// All other rows
if(p0.x == -1)
{
// First column
++loc.x();
cubic_interpolate_current()(*(loc.y()-1), *loc, *(loc.y()+1), *(loc.y()+2), frac.y, mpres);
}
else if(p0.x == 0)
{
// 2nd column
--loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*loc, *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else if(p0.x+1 == src.width())
{
// Last column
cubic_interpolate_current()(*(loc.y()-1), *loc, *(loc.y()+1), *(loc.y()+2), frac.y, mpres);
}
else if(p0.x+2 == src.width())
{
// 2nd last column
--loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+1), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);
}
else
{
// All general cases, center pixels
--loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp0);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp1);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp2);
++loc.y();
cubic_interpolate_current()(*(loc.x()-1), *loc, *(loc.x()+1), *(loc.x()+2), frac.x, mp3);

cubic_interpolate_current()(mp0, mp1, mp2, mp3, frac.y, mpres);

}
}

SrcP src_result;

cast_pixel(mpres, src_result);

color_convert(src_result, result);
return true;
}

}} // namespace boost::gil

#endif
Loading