27#ifndef VIGRA_EXT_INTERPOLATORS_H
28#define VIGRA_EXT_INTERPOLATORS_H
37#include <vigra/accessor.hxx>
38#include <vigra/diff2d.hxx>
46static double sinc (
double x );
47static double cubic01 (
double x );
48static double cubic12 (
double x );
50static double sinc(
double x )
62static const double A(-0.75);
67 return ((
A + 2.0 )*x - (
A + 3.0 ))*x*x +1.0;
73 return ((
A * x - 5.0 *
A ) * x + 8.0 *
A ) * x - 4.0 *
A;
103 w[1] = (x >= 0.5) ? 1 : 0;
104 w[0] = (x < 0.5) ? 1 : 0;
108 oss <<
" return (i == 0.0) ? float(f < 0.5) : float(f >= 0.5);" << endl;
126 oss <<
" return abs(i + f - 1.0);" << endl;
146 oss <<
" float A = " <<
A <<
";" << endl
147 <<
" float c = abs(i - 1.0);" << endl
148 <<
" float m = (i > 1.0) ? -1.0 : 1.0;" << endl
149 <<
" float p = c + m * f;" << endl
150 <<
" if (i == 1.0 || i == 2.0) {" << endl
151 <<
" return (( A + 2.0 )*p - ( A + 3.0 ))*p*p + 1.0;" << endl
152 <<
" } else {" << endl
153 <<
" return (( A * p - 5.0 * A ) * p + 8.0 * A ) * p - 4.0 * A;" << endl
167 w[3] = ( ( 1.0/3.0 * x - 1.0/5.0 ) * x - 2.0/15.0 ) * x;
168 w[2] = ( ( 6.0/5.0 - x ) * x + 4.0/5.0 ) * x;
169 w[1] = ( ( x - 9.0/5.0 ) * x - 1.0/5.0 ) * x + 1.0;
170 w[0] = ( ( -1.0/3.0 * x + 4.0/5.0 ) * x - 7.0/15.0 ) * x;
174 oss <<
" return (i > 1.0) ? (i == 3.0) ? (( ( 1.0/3.0 * f - 1.0/5.0 ) * f - 2.0/15.0 ) * f)" << endl
175 <<
" : (( ( 6.0/5.0 - f ) * f + 4.0/5.0 ) * f)" << endl
176 <<
" : (i == 1.0) ? (( ( f - 9.0/5.0 ) * f - 1.0/5.0 ) * f + 1.0)" << endl
177 <<
" : (( ( -1.0/3.0 * f + 4.0/5.0 ) * f - 7.0/15.0 ) * f);" << endl;
195 w[5] = ( ( - 1.0/11.0 * x + 12.0/ 209.0 ) * x + 7.0/ 209.0 ) * x;
196 w[4] = ( ( 6.0/11.0 * x - 72.0/ 209.0 ) * x - 42.0/ 209.0 ) * x;
197 w[3] = ( ( - 13.0/11.0 * x + 288.0/ 209.0 ) * x + 168.0/ 209.0 ) * x;
198 w[2] = ( ( 13.0/11.0 * x - 453.0/ 209.0 ) * x - 3.0/ 209.0 ) * x + 1.0;
199 w[1] = ( ( - 6.0/11.0 * x + 270.0/ 209.0 ) * x - 156.0/ 209.0 ) * x;
200 w[0] = ( ( 1.0/11.0 * x - 45.0/ 209.0 ) * x + 26.0/ 209.0 ) * x;
204 oss <<
" return (i > 3.0) ? (i == 5.0) ? (( ( - 1.0/11.0 * f + 12.0/ 209.0 ) * f + 7.0/ 209.0 ) * f)" << endl
205 <<
" : (( ( 6.0/11.0 * f - 72.0/ 209.0 ) * f - 42.0/ 209.0 ) * f)" << endl
206 <<
" : (i > 1.0) ? (i == 3.0) ? (( ( - 13.0/11.0 * f + 288.0/ 209.0 ) * f + 168.0/ 209.0 ) * f)" << endl
207 <<
" : (( ( 13.0/11.0 * f - 453.0/ 209.0 ) * f - 3.0/ 209.0 ) * f + 1.0)" << endl
208 <<
" : (i == 1.0) ? (( ( - 6.0/11.0 * f + 270.0/ 209.0 ) * f - 156.0/ 209.0 ) * f)" << endl
209 <<
" : (( ( 1.0/11.0 * f - 45.0/ 209.0 ) * f + 26.0/ 209.0 ) * f);" << endl;
223 w[7] = (( 1.0/41.0 * x - 45.0/2911.0) * x - 26.0/2911.0) * x;
224 w[6] = ((- 6.0/41.0 * x + 270.0/2911.0) * x + 156.0/2911.0) * x;
225 w[5] = (( 24.0/41.0 * x - 1080.0/2911.0) * x - 624.0/2911.0) * x;
226 w[4] = ((-49.0/41.0 * x + 4050.0/2911.0) * x + 2340.0/2911.0) * x;
227 w[3] = (( 49.0/41.0 * x - 6387.0/2911.0) * x - 3.0/2911.0) * x + 1.0;
228 w[2] = ((-24.0/41.0 * x + 4032.0/2911.0) * x - 2328.0/2911.0) * x;
229 w[1] = (( 6.0/41.0 * x - 1008.0/2911.0) * x + 582.0/2911.0) * x;
230 w[0] = ((- 1.0/41.0 * x + 168.0/2911.0) * x - 97.0/2911.0) * x;
234 oss <<
" return (i > 3.0) ? (i > 5.0) ? (i == 7.0) ? ((( 1.0/41.0 * f - 45.0/2911.0) * f - 26.0/2911.0) * f)" << endl
235 <<
" : (((- 6.0/41.0 * f + 270.0/2911.0) * f + 156.0/2911.0) * f)" << endl
236 <<
" : (i == 5.0) ? ((( 24.0/41.0 * f - 1080.0/2911.0) * f - 624.0/2911.0) * f)" << endl
237 <<
" : (((-49.0/41.0 * f + 4050.0/2911.0) * f + 2340.0/2911.0) * f)" << endl
238 <<
" : (i > 1.0) ? (i == 3.0) ? ((( 49.0/41.0 * f - 6387.0/2911.0) * f - 3.0/2911.0) * f + 1.0)" << endl
239 <<
" : (((-24.0/41.0 * f + 4032.0/2911.0) * f - 2328.0/2911.0) * f)" << endl
240 <<
" : (i == 1.0) ? ((( 6.0/41.0 * f - 1008.0/2911.0) * f + 582.0/2911.0) * f)" << endl
241 <<
" : (((- 1.0/41.0 * f + 168.0/2911.0) * f - 97.0/2911.0) * f);" << endl;
272 oss <<
" float c = (i < " << (
size/2.0) <<
") ? 1.0 : -1.0;" << endl
273 <<
" float x = c * (" << (
size/2 - 1.0) <<
" - i + f);" << endl
274 <<
" vec2 xpi = vec2(x, x / " << (
size/2.0) <<
") * " <<
M_PI <<
";" << endl
275 <<
" vec2 xsin = sin(xpi);" << endl
276 <<
" vec2 result = vec2(1.0, 1.0);" << endl
277 <<
" if (xpi.x != 0.0) result.x = xsin.x / xpi.x;" << endl
278 <<
" if (xpi.y != 0.0) result.y = xsin.y / xpi.y;" << endl
279 <<
" return result.x * result.y;" << endl;
297 typedef typename vigra::NumericTraits<PixelType>::RealPromote
RealPixelType;
362 if (
srcx > INTERPOLATOR::size/2 &&
srcx <
m_w -INTERPOLATOR::size/2 &&
363 srcy > INTERPOLATOR::size/2 &&
srcy <
m_h - INTERPOLATOR::size/2)
368 double wx[INTERPOLATOR::size];
369 double wy[INTERPOLATOR::size];
375 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
377 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++) {
385 for (
int kx = 0;
kx < INTERPOLATOR::size;
kx++) {
417 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
426 double w[INTERPOLATOR::size];
435 vigra::Diff2D offset(
srcx - INTERPOLATOR::size/2 + 1,
436 srcy - INTERPOLATOR::size/2 + 1);
438 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++, ++(
ys.y)) {
439 p = vigra::NumericTraits<RealPixelType>::zero();
440 typename SrcImageIterator::row_iterator
xs(
ys.rowIterator());
442 for (
int kx = 0;
kx < INTERPOLATOR::size;
kx++, ++
xs) {
450 p = vigra::NumericTraits<RealPixelType>::zero();
451 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++) {
455 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
478 typedef typename MaskAccessor::value_type
MaskType;
480 typedef typename vigra::NumericTraits<PixelType>::RealPromote
RealPixelType;
496 std::pair<MaskIterator, MaskAccessor> mask,
568 double w[INTERPOLATOR::size];
570 double weightsX[INTERPOLATOR::size];
578 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++) {
592 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
595 for (
int kx = 0;
kx < INTERPOLATOR::size;
kx++) {
626 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
628 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++) {
637 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
676 if (
srcx > INTERPOLATOR::size/2 &&
srcx <
m_w -INTERPOLATOR::size/2 &&
677 srcy > INTERPOLATOR::size/2 &&
srcy <
m_h - INTERPOLATOR::size/2)
682 double wx[INTERPOLATOR::size];
683 double wy[INTERPOLATOR::size];
691 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
694 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++) {
702 for (
int kx = 0;
kx < INTERPOLATOR::size;
kx++) {
742 mask = vigra::detail::RequiresExplicitCast<MaskType>::cast(m);
743 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
753 double wx[INTERPOLATOR::size];
754 double wy[INTERPOLATOR::size];
760 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
763 vigra::Diff2D offset(
srcx - INTERPOLATOR::size/2 + 1,
764 srcy - INTERPOLATOR::size/2 + 1);
767 for (
int ky = 0;
ky < INTERPOLATOR::size;
ky++, ++(
ys.y), ++(
yms.y)) {
769 typename SrcImageIterator::row_iterator
xs(
ys.rowIterator());
770 typename MaskIterator::row_iterator
xms(
yms.rowIterator());
771 for (
int kx = 0;
kx < INTERPOLATOR::size;
kx++, ++
xs, ++
xms) {
794 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
795 mask = vigra::detail::RequiresExplicitCast<MaskType>::cast(m);
856template <
class ACCESSOR,
class VALUETYPE,
class INTERPOLATOR>
886 template <
class ITERATOR>
893 double wx[INTERPOLATOR::size];
894 double wy[INTERPOLATOR::size];
897 typename vigra::NumericTraits<value_type>::RealPromote
898 ret (vigra::NumericTraits<value_type>::zero());
906 for(
int y = 0; y <
inter_y.size; ++y, ++
ys.y) {
908 for(
int x = 0; x <
inter_x.size; x++, ++
xs.x) {
912 return vigra::detail::RequiresExplicitCast<value_type>::cast(
ret);
936 template <
class ITERATOR,
class ALPHAITERATOR,
class ALPHAACCESSOR>
944 double wx[INTERPOLATOR::size];
945 double wy[INTERPOLATOR::size];
948 typename vigra::NumericTraits<value_type>::RealPromote
949 ret (vigra::NumericTraits<value_type>::zero());
959 for(
int y = 0; y <
inter_y.size; ++y, ++
ys.y, ++
ays.y) {
962 for(
int x = 0; x <
inter_x.size; x++, ++
xs.x, ++
axs.x) {
963 if (alpha.second(
axs) <= 0 ) {
969 result = vigra::detail::RequiresExplicitCast<value_type>::cast(
ret);
"wrapper" for efficient interpolation access to an image
bool operator()(double x, double y, PixelType &result, MaskType &mask) const
Interpolate without mask, but return dummy alpha value nevertheless.
ImageInterpolator(SrcImageIterator src_upperleft, SrcImageIterator src_lowerright, SrcAccessor sa, INTERPOLATOR &inter, bool warparound)
Construct interpolator for specific image.
bool operator()(double x, double y, PixelType &result) const
Interpolate without mask.
ImageInterpolator(vigra::triple< SrcImageIterator, SrcImageIterator, SrcAccessor > const &src, INTERPOLATOR &inter, bool warparound)
Construct interpolator for an given image.
vigra::NumericTraits< PixelType >::RealPromote RealPixelType
SrcAccessor::value_type PixelType
bool interpolateNoMaskInside(int srcx, int srcy, double dx, double dy, PixelType &result) const
Interpolate without boundary check and mask.
void emitGLSL(std::ostringstream &oss) const
"wrapper" for efficient interpolation access to an image
bool operator()(double x, double y, PixelType &result, MaskType &mask) const
Interpolate the data item at a non-integer position x, y.
MaskAccessor::value_type MaskType
bool interpolateInside(int srcx, int srcy, double dx, double dy, PixelType &result, MaskType &mask) const
Interpolate without boundary check.
ImageMaskInterpolator(SrcImageIterator src_upperleft, SrcImageIterator src_lowerright, SrcAccessor sa, MaskIterator mask_upperleft, MaskAccessor ma, INTERPOLATOR &inter, bool warparound)
Construct interpolator for specific image.
ImageMaskInterpolator(vigra::triple< SrcImageIterator, SrcImageIterator, SrcAccessor > const &src, std::pair< MaskIterator, MaskAccessor > mask, INTERPOLATOR &inter, bool warparound)
Construct interpolator for an given image.
vigra::NumericTraits< PixelType >::RealPromote RealPixelType
SrcAccessor::value_type PixelType
interpolation at non-integer positions.
VALUETYPE value_type
the iterators' pixel type
InterpolatingAccessor(ACCESSOR a, INTERPOLATOR inter)
init from given accessor
value_type operator()(ITERATOR const &i, float x, float y) const
Interpolate the data item at a non-integer position x, y.
bool operator()(ITERATOR const &i, std::pair< ALPHAITERATOR, ALPHAACCESSOR > const &alpha, float x, float y, value_type &result) const
Interpolate the data item at a non-integer position x, y.
misc math function & classes used by other parts of the program
static double sinc(double x)
static const double A(-0.75)
static double cubic01(double x)
Interpolator
enum with all interpolation methods
@ INTERP_NEAREST_NEIGHBOUR
static double cubic12(double x)
simple bilinear interpolation
void calc_coeff(double x, double *w) const
void emitGLSL(std::ostringstream &oss) const
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
initialize weights for given x
several classes to calculate interpolator weights,
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
sinc interpolation, with variable width
void calc_coeff(double x, double *w) const
initialize weights for given offset x
void emitGLSL(std::ostringstream &oss) const
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
initialize weights for given x
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
calculate weights for given offset x.
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
initialize weights for given offset x
std::vector< deghosting::BImagePtr > threshold(const std::vector< deghosting::FImagePtr > &inputImages, const double threshold, const uint16_t flags)
Threshold function used for creating alpha masks for images.