Hugin trunk 0.1
Loading...
Searching...
No Matches
Interpolators.h
Go to the documentation of this file.
1// -*- c-basic-offset: 4 -*-
27#ifndef VIGRA_EXT_INTERPOLATORS_H
28#define VIGRA_EXT_INTERPOLATORS_H
29
30#include <iostream>
31#include <iomanip>
32
33#include <math.h>
35#include <algorithm>
36
37#include <vigra/accessor.hxx>
38#include <vigra/diff2d.hxx>
39
40using std::endl;
41
42namespace vigra_ext {
43
44// Some locally needed math functions
45
46static double sinc ( double x );
47static double cubic01 ( double x );
48static double cubic12 ( double x );
49
50static double sinc( double x )
51{
52 x *= M_PI;
53 if(x != 0.0)
54 return(sin(x) / x);
55 return(1.0);
56}
57
58
59// Cubic polynomial with parameter A
60// A = -1: sharpen; A = - 0.5 homogeneous
61// make sure x >= 0
62static const double A(-0.75);
63
64// 0 <= x < 1
65static double cubic01( double x )
66{
67 return (( A + 2.0 )*x - ( A + 3.0 ))*x*x +1.0;
68}
69// 1 <= x < 2
70
71static double cubic12( double x )
72{
73 return (( A * x - 5.0 * A ) * x + 8.0 * A ) * x - 4.0 * A;
74
75}
76
88
89
97{
98 // size of neighbourhood
99 static const int size = 2;
100
101 void calc_coeff(double x, double * w) const
102 {
103 w[1] = (x >= 0.5) ? 1 : 0;
104 w[0] = (x < 0.5) ? 1 : 0;
105 }
106
107 void emitGLSL(std::ostringstream& oss) const {
108 oss << " return (i == 0.0) ? float(f < 0.5) : float(f >= 0.5);" << endl;
109 }
110};
111
112
115{
116 // size of neighbourhood
117 static const int size = 2;
118
119 void calc_coeff(double x, double * w) const
120 {
121 w[1] = x;
122 w[0] = 1.0-x;
123 }
124
125 void emitGLSL(std::ostringstream& oss) const {
126 oss << " return abs(i + f - 1.0);" << endl;
127 }
128};
129
132{
133 // size of neighbourhood
134 static const int size = 4;
135
137 void calc_coeff(double x, double * w) const
138 {
139 w[3] = cubic12( 2.0 - x );
140 w[2] = cubic01( 1.0 - x );
141 w[1] = cubic01( x );
142 w[0] = cubic12( x + 1.0 );
143 }
144
145 void emitGLSL(std::ostringstream& oss) const {
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
154 << " }" << endl;
155 }
156};
157
160{
161 // size of neighbourhood
162 static const int size = 4;
163
165 void calc_coeff(double x, double * w) const
166 {
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;
171 }
172
173 void emitGLSL(std::ostringstream& oss) const {
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;
178 }
179};
180
183{
184 // size of neighbourhood
185 static const int size = 6;
186
193 void calc_coeff(double x, double* w) const
194 {
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;
201 }
202
203 void emitGLSL(std::ostringstream& oss) const {
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;
210 }
211};
212
213
216{
217 // size of neighbourhood
218 static const int size = 8;
219
221 void calc_coeff(double x, double * w) const
222 {
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;
231 }
232
233 void emitGLSL(std::ostringstream& oss) const {
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;
242 }
243};
244
246template <int size_>
248{
249 // size of neighbourhood
250 static const int size = size_;
251
253 void calc_coeff(double x, double * w) const
254 {
255 int idx;
256 double xadd;
257 for( idx = 0, xadd = size / 2 - 1.0 + x;
258 idx < size / 2;
259 xadd-=1.0)
260 {
261 w[idx++] = sinc( xadd ) * sinc( xadd / ( size / 2 ));
262 }
263 for( xadd = 1.0 - x;
264 idx < size;
265 xadd+=1.0)
266 {
267 w[idx++] = sinc( xadd ) * sinc( xadd / ( size / 2 ));
268 }
269 }
270
271 void emitGLSL(std::ostringstream& oss) const {
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;
280 }
281};
282
283
288template <typename SrcImageIterator, typename SrcAccessor,
289 typename INTERPOLATOR>
291{
292public:
293 typedef typename SrcAccessor::value_type PixelType;
294 // dummy mask type to be compatible to algorithms expecting a ImageMaskInterpolator object
295 typedef typename vigra::UInt8 MaskType;
296private:
297 typedef typename vigra::NumericTraits<PixelType>::RealPromote RealPixelType;
298
301 int m_w;
302 int m_h;
304
306
307public:
309 ImageInterpolator(vigra::triple<SrcImageIterator, SrcImageIterator,SrcAccessor> const & src,
311 bool warparound)
312 : m_sIter(src.first),
314 m_w(src.second.x - src.first.x),
315 m_h(src.second.y - src.first.y),
318 {
319 }
320
337
339 bool operator()(double x, double y,
340 PixelType & result, MaskType & mask) const
341 {
342 mask = 255;
343 return operator()(x,y, result);
344 }
345
347 bool operator()(double x, double y,
348 PixelType & result) const
349 {
350
351 // skip all further interpolation if we cannot interpolate anything
352 if (x < -INTERPOLATOR::size/2 || x > m_w + INTERPOLATOR::size/2) return false;
353 if (y < -INTERPOLATOR::size/2 || y > m_h + INTERPOLATOR::size/2) return false;
354
355 double t = floor(x);
356 double dx = x - t;
357 int srcx = int(t);
358 t = floor(y);
359 double dy = y - t;
360 int srcy = int(t);
361
362 if ( srcx > INTERPOLATOR::size/2 && srcx < m_w -INTERPOLATOR::size/2 &&
363 srcy > INTERPOLATOR::size/2 && srcy < m_h - INTERPOLATOR::size/2)
364 {
365 return interpolateNoMaskInside(srcx, srcy, dx, dy, result);
366 }
367
368 double wx[INTERPOLATOR::size];
369 double wy[INTERPOLATOR::size];
370
371 // calculate x interpolation coefficients
372 m_inter.calc_coeff(dx, wx);
373 m_inter.calc_coeff(dy, wy);
374
375 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
376 double weightsum = 0.0;
377 for (int ky = 0; ky < INTERPOLATOR::size; ky++) {
378 int bounded_ky = srcy + 1 + ky - INTERPOLATOR::size/2;
379
380 // Boundary condition: do not replicate top and bottom
382 continue;
383 }
384
385 for (int kx = 0; kx < INTERPOLATOR::size; kx++) {
386 int bounded_kx = srcx + 1 + kx - INTERPOLATOR::size/2;
387
388 if (m_warparound) {
389 // Boundary condition: wrap around the image.
390 if (bounded_kx < 0)
391 bounded_kx += m_w;
392 if (bounded_kx >= m_w)
393 bounded_kx -= m_w;
394 } else {
395 // Boundary condition: replicate first and last column.
396 // if (srcx + kx < 0) bounded_kx -= (srcx + kx);
397 // if (srcx + kx >= src_w) bounded_kx -= (srcx + kx - (src_w - 1));
398 // Boundary condition: do not replicate left and right
399 if (bounded_kx < 0)
400 continue;
401 if (bounded_kx >= m_w)
402 continue;
403 }
404
405 // check mask
406 double f = wx[kx]*wy[ky];
407 p += f * m_sAcc(m_sIter, vigra::Diff2D(bounded_kx, bounded_ky));
408 weightsum += f;
409 }
410 }
411
412 // force a certain weight
413 if (weightsum <= 0.2) return false;
414 // Adjust filter for any ignored transparent pixels.
415 if (weightsum != 1.0) p /= weightsum;
416
417 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
418 return true;
419 }
420
421
423 bool interpolateNoMaskInside(int srcx, int srcy, double dx, double dy,
424 PixelType & result) const
425 {
426 double w[INTERPOLATOR::size];
427 RealPixelType resX[INTERPOLATOR::size];
428
429 // calculate x interpolation coefficients
430 m_inter.calc_coeff(dx, w);
431
433
434 // first pass of separable filter, x pass
435 vigra::Diff2D offset(srcx - INTERPOLATOR::size/2 + 1,
436 srcy - INTERPOLATOR::size/2 + 1);
437 SrcImageIterator ys(m_sIter + offset);
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());
441 //SrcImageIterator xs(ys);
442 for (int kx = 0; kx < INTERPOLATOR::size; kx++, ++xs) {
443 p += w[kx] * m_sAcc(xs);
444 }
445 resX[ky] = p;
446 }
447
448 // y pass.
449 m_inter.calc_coeff(dy, w);
450 p = vigra::NumericTraits<RealPixelType>::zero();
451 for (int ky = 0; ky < INTERPOLATOR::size; ky++) {
452 p += w[ky] * resX[ky];
453 }
454
455 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
456 return true;
457 }
458
459 void emitGLSL(std::ostringstream& oss) const {
460 m_inter.emitGLSL(oss);
461 }
462
463};
464
465
471template <typename SrcImageIterator, typename SrcAccessor,
472 typename MaskIterator, typename MaskAccessor,
473 typename INTERPOLATOR>
475{
476public:
477 typedef typename SrcAccessor::value_type PixelType;
478 typedef typename MaskAccessor::value_type MaskType;
479private:
480 typedef typename vigra::NumericTraits<PixelType>::RealPromote RealPixelType;
481
485 MaskAccessor m_mAcc;
486 int m_w;
487 int m_h;
489
491
492public:
493
495 ImageMaskInterpolator(vigra::triple<SrcImageIterator, SrcImageIterator,SrcAccessor> const & src,
496 std::pair<MaskIterator, MaskAccessor> mask,
498 bool warparound)
499 : m_sIter(src.first),
501 m_mIter(mask.first),
502 m_mAcc(mask.second),
503 m_w(src.second.x - src.first.x),
504 m_h(src.second.y - src.first.y),
507 {
508 }
509
530#if 0
550// bool operator()(float x, float y,
551 // this is slower than the full version, thanks to the normalized interpolation (with masks).
552 bool interpolateSeperable(float x, float y,
553 PixelType & result) const
554 {
555
556 // skip all further interpolation if we cannot interpolate anything
557 if (x < -INTERPOLATOR::size/2 || x > m_w + INTERPOLATOR::size/2) return false;
558 if (y < -INTERPOLATOR::size/2 || y > m_h + INTERPOLATOR::size/2) return false;
559
560 double t = floor(x);
561 double dx = x - t;
562 int srcx = int(t);
563 t = floor(y);
564 double dy = y - t;
565 int srcy = int(t);
566
567
568 double w[INTERPOLATOR::size];
569
570 double weightsX[INTERPOLATOR::size];
571 PixelType resX[INTERPOLATOR::size];
572
573 // calculate x interpolation coefficients
574 m_inter.calc_coeff(dx, w);
575
576 // first pass of separable filter
577
578 for (int ky = 0; ky < INTERPOLATOR::size; ky++) {
579 int bounded_ky = srcy + 1 + ky - INTERPOLATOR::size/2;
580
581 // Boundary condition: replicate top and bottom rows.
582 // if (srcy + ky < 0) bounded_ky -= (srcy + ky);
583 // if (srcy + ky >= src_h) bounded_ky -= (srcy + ky - (src_h - 1));
584
585 // Boundary condition: do not replicate top and bottom
587 weightsX[ky] = 0;
588 resX[ky] = 0;
589 continue;
590 }
591
592 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
593 double weightsum = 0.0;
594
595 for (int kx = 0; kx < INTERPOLATOR::size; kx++) {
596 int bounded_kx = srcx + 1 + kx - INTERPOLATOR::size/2;
597
598 if (m_warparound) {
599 // Boundary condition: wrap around the image.
600 if (bounded_kx < 0)
601 bounded_kx += m_w;
602 if (bounded_kx >= m_w)
603 bounded_kx -= m_w;
604 } else {
605 // Boundary condition: replicate first and last column.
606 // if (srcx + kx < 0) bounded_kx -= (srcx + kx);
607 // if (srcx + kx >= src_w) bounded_kx -= (srcx + kx - (src_w - 1));
608 // Boundary condition: do not replicate left and right
609 if (bounded_kx < 0)
610 continue;
611 if (bounded_kx >= m_w)
612 continue;
613 }
615 // check mask
616 p += w[kx] * m_sIter(bounded_kx, bounded_ky);
617 weightsum += w[kx];
618 }
619 }
621 resX[ky] = p;
622 }
623
624 // y pass.
625 m_inter.calc_coeff(dy, w);
626 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
627 double weightsum = 0.0;
628 for (int ky = 0; ky < INTERPOLATOR::size; ky++) {
629 weightsum += weightsX[ky] * w[ky];
630 p += w[ky] * resX[ky];
631 }
632
633 if (weightsum == 0.0) return false;
634 // Adjust filter for any ignored transparent pixels.
635 if (weightsum != 1.0) p /= weightsum;
636
637 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
638 return true;
639 }
640#endif
641
661 bool operator()(double x, double y,
662 PixelType & result, MaskType & mask) const
663 {
664
665 // skip all further interpolation if we cannot interpolate anything
666 if (x < -INTERPOLATOR::size/2 || x > m_w + INTERPOLATOR::size/2) return false;
667 if (y < -INTERPOLATOR::size/2 || y > m_h + INTERPOLATOR::size/2) return false;
668
669 double t = floor(x);
670 double dx = x - t;
671 int srcx = int(t);
672 t = floor(y);
673 double dy = y - t;
674 int srcy = int(t);
675
676 if ( srcx > INTERPOLATOR::size/2 && srcx < m_w -INTERPOLATOR::size/2 &&
677 srcy > INTERPOLATOR::size/2 && srcy < m_h - INTERPOLATOR::size/2)
678 {
679 return interpolateInside(srcx, srcy, dx, dy, result, mask);
680 }
681
682 double wx[INTERPOLATOR::size];
683 double wy[INTERPOLATOR::size];
684
685 // calculate x interpolation coefficients
686 m_inter.calc_coeff(dx, wx);
687 m_inter.calc_coeff(dy, wy);
688
689 // first pass of separable filter
690
691 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
692 double m = 0;
693 double weightsum = 0.0;
694 for (int ky = 0; ky < INTERPOLATOR::size; ky++) {
695 int bounded_ky = srcy + 1 + ky - INTERPOLATOR::size/2;
696
697 // Boundary condition: do not replicate top and bottom
699 continue;
700 }
701
702 for (int kx = 0; kx < INTERPOLATOR::size; kx++) {
703 int bounded_kx = srcx + 1 + kx - INTERPOLATOR::size/2;
704
705 if (m_warparound) {
706 // Boundary condition: wrap around the image.
707 if (bounded_kx < 0)
708 bounded_kx += m_w;
709 if (bounded_kx >= m_w)
710 bounded_kx -= m_w;
711 } else {
712 // Boundary condition: replicate first and last column.
713 // if (srcx + kx < 0) bounded_kx -= (srcx + kx);
714 // if (srcx + kx >= src_w) bounded_kx -= (srcx + kx - (src_w - 1));
715 // Boundary condition: do not replicate left and right
716 if (bounded_kx < 0)
717 continue;
718 if (bounded_kx >= m_w)
719 continue;
720 }
721
723 if (cmask) {
724 // check mask
725 double f = wx[kx]*wy[ky];
726 // TODO: check if this is good, influences the HDR stitching masks
727 m += f * cmask;
728 p += f * m_sAcc(m_sIter, vigra::Diff2D(bounded_kx, bounded_ky));
729 weightsum += f;
730 }
731 }
732 }
733
734 // force a certain weight
735 if (weightsum <= 0.2) return false;
736 // Adjust filter for any ignored transparent pixels.
737 if (weightsum != 1.0) {
738 p /= weightsum;
739 m /= weightsum;
740 }
741
742 mask = vigra::detail::RequiresExplicitCast<MaskType>::cast(m);
743 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
744 return true;
745 }
746
747
749 bool interpolateInside(int srcx, int srcy, double dx, double dy,
750 PixelType & result, MaskType & mask) const
751 {
752
753 double wx[INTERPOLATOR::size];
754 double wy[INTERPOLATOR::size];
755
756 // calculate x interpolation coefficients
757 m_inter.calc_coeff(dx, wx);
758 m_inter.calc_coeff(dy, wy);
759
760 RealPixelType p(vigra::NumericTraits<RealPixelType>::zero());
761 double weightsum = 0.0;
762 double m = 0.0;
763 vigra::Diff2D offset(srcx - INTERPOLATOR::size/2 + 1,
764 srcy - INTERPOLATOR::size/2 + 1);
765 SrcImageIterator ys(m_sIter + offset);
766 MaskIterator yms(m_mIter + offset);
767 for (int ky = 0; ky < INTERPOLATOR::size; ky++, ++(ys.y), ++(yms.y)) {
768// int bounded_ky = srcy + 1 + ky - INTERPOLATOR::size/2;
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) {
772// int bounded_kx = srcx + 1 + kx - INTERPOLATOR::size/2;
773
774 MaskType cmask = *xms;
775 if (cmask) {
776 // check mask
777 double f = wx[kx]*wy[ky];
778 // TODO: check if this is good, influences the HDR stitching masks
779 m += f * cmask;
780 p += f * m_sAcc(xs);
781 weightsum += f;
782 }
783 }
784 }
785
786 // force a certain weight
787 if (weightsum <= 0.2) return false;
788 // Adjust filter for any ignored transparent pixels.
789 if (weightsum != 1.0) {
790 p /= weightsum;
791 m /= weightsum;
792 }
793
794 result = vigra::detail::RequiresExplicitCast<PixelType>::cast(p);
795 mask = vigra::detail::RequiresExplicitCast<MaskType>::cast(m);
796 return true;
797 }
798
799};
800
801/********************************************************/
802/* */
803/* InterpolatingAccessor */
804/* */
805/********************************************************/
806
856template <class ACCESSOR, class VALUETYPE, class INTERPOLATOR>
858{
859public:
863
869
886 template <class ITERATOR>
887 value_type operator()(ITERATOR const & i, float x, float y) const
888 {
889 int ix = int(x);
890 int iy = int(y);
891 float dx = x - ix;
892 float dy = y - iy;
893 double wx[INTERPOLATOR::size];
894 double wy[INTERPOLATOR::size];
895
896 // promote value_type for multiplication
897 typename vigra::NumericTraits<value_type>::RealPromote
898 ret (vigra::NumericTraits<value_type>::zero());
899
900 // calculate interpolation coefficients
901 inter_x.calc_coeff(dx, wx);
902 inter_y.calc_coeff(dy, wy);
903
904 ITERATOR ys(i + vigra::Diff2D(ix - inter_x.size/2 + 1,
905 iy - inter_y.size/2 + 1));
906 for(int y = 0; y < inter_y.size; ++y, ++ys.y) {
907 ITERATOR xs(ys);
908 for(int x = 0; x < inter_x.size; x++, ++xs.x) {
909 ret += wx[x] * wy[y] * a_(xs);
910 }
911 }
912 return vigra::detail::RequiresExplicitCast<value_type>::cast(ret);
913 }
914
936 template <class ITERATOR, class ALPHAITERATOR, class ALPHAACCESSOR>
937 bool operator()(ITERATOR const & i, std::pair<ALPHAITERATOR, ALPHAACCESSOR> const & alpha,
938 float x, float y, value_type & result) const
939 {
940 int ix = int(x);
941 int iy = int(y);
942 float dx = x - ix;
943 float dy = y - iy;
944 double wx[INTERPOLATOR::size];
945 double wy[INTERPOLATOR::size];
946
947 // promote value_type for multiplication
948 typename vigra::NumericTraits<value_type>::RealPromote
949 ret (vigra::NumericTraits<value_type>::zero());
950
951 // calculate interpolation coefficients
952 inter_x.calc_coeff(dx, wx);
953 inter_y.calc_coeff(dy, wy);
954
955 ITERATOR ys(i + vigra::Diff2D(ix - inter_x.size/2 + 1,
956 iy - inter_y.size/2 + 1));
957 ALPHAITERATOR ays(alpha.first + vigra::Diff2D(ix - inter_x.size/2 + 1,
958 iy - inter_y.size/2 + 1));
959 for(int y = 0; y < inter_y.size; ++y, ++ys.y, ++ays.y) {
960 ITERATOR xs(ys);
962 for(int x = 0; x < inter_x.size; x++, ++xs.x, ++axs.x) {
963 if (alpha.second(axs) <= 0 ) {
964 return false;
965 }
966 ret += wx[x] * wy[y] * a_(xs);
967 }
968 }
969 result = vigra::detail::RequiresExplicitCast<value_type>::cast(ret);
970 return true;
971 }
972
973
974private:
977};
978
979} // namespace
980
981#endif // VIGRA_EXT_INTERPOLATORS_H
#define M_PI
"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
static const int size
void calc_coeff(double x, double *w) const
void emitGLSL(std::ostringstream &oss) const
cubic interpolation
static const int size
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
static const int size
void calc_coeff(double x, double *w) const
initialize weights for given offset x
void emitGLSL(std::ostringstream &oss) const
spline16 interpolation
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
initialize weights for given x
spline36 interpolation
void emitGLSL(std::ostringstream &oss) const
void calc_coeff(double x, double *w) const
calculate weights for given offset x.
spline64 interpolation
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.
Definition threshold.h:41