PCL
LanczosInterpolation.h
Go to the documentation of this file.
1 // ____ ______ __
2 // / __ \ / ____// /
3 // / /_/ // / / /
4 // / ____// /___ / /___ PixInsight Class Library
5 // /_/ \____//_____/ PCL 2.1.19
6 // ----------------------------------------------------------------------------
7 // pcl/LanczosInterpolation.h - Released 2019-11-07T10:59:34Z
8 // ----------------------------------------------------------------------------
9 // This file is part of the PixInsight Class Library (PCL).
10 // PCL is a multiplatform C++ framework for development of PixInsight modules.
11 //
12 // Copyright (c) 2003-2019 Pleiades Astrophoto S.L. All Rights Reserved.
13 //
14 // Redistribution and use in both source and binary forms, with or without
15 // modification, is permitted provided that the following conditions are met:
16 //
17 // 1. All redistributions of source code must retain the above copyright
18 // notice, this list of conditions and the following disclaimer.
19 //
20 // 2. All redistributions in binary form must reproduce the above copyright
21 // notice, this list of conditions and the following disclaimer in the
22 // documentation and/or other materials provided with the distribution.
23 //
24 // 3. Neither the names "PixInsight" and "Pleiades Astrophoto", nor the names
25 // of their contributors, may be used to endorse or promote products derived
26 // from this software without specific prior written permission. For written
27 // permission, please contact info@pixinsight.com.
28 //
29 // 4. All products derived from this software, in any form whatsoever, must
30 // reproduce the following acknowledgment in the end-user documentation
31 // and/or other materials provided with the product:
32 //
33 // "This product is based on software from the PixInsight project, developed
34 // by Pleiades Astrophoto and its contributors (http://pixinsight.com/)."
35 //
36 // Alternatively, if that is where third-party acknowledgments normally
37 // appear, this acknowledgment must be reproduced in the product itself.
38 //
39 // THIS SOFTWARE IS PROVIDED BY PLEIADES ASTROPHOTO AND ITS CONTRIBUTORS
40 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
41 // TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
42 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL PLEIADES ASTROPHOTO OR ITS
43 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
44 // EXEMPLARY OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, BUSINESS
45 // INTERRUPTION; PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; AND LOSS OF USE,
46 // DATA OR PROFITS) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
47 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
48 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
49 // POSSIBILITY OF SUCH DAMAGE.
50 // ----------------------------------------------------------------------------
51 
52 #ifndef __PCL_LanczosInterpolation_h
53 #define __PCL_LanczosInterpolation_h
54 
56 
57 #include <pcl/Defs.h>
58 #include <pcl/Diagnostics.h>
59 
61 #include <pcl/Math.h>
62 #include <pcl/Utility.h>
63 #include <pcl/Vector.h>
64 
65 namespace pcl
66 {
67 
68 // ----------------------------------------------------------------------------
69 
70 #define m_width this->m_width
71 #define m_height this->m_height
72 #define m_fillBorder this->m_fillBorder
73 #define m_fillValue this->m_fillValue
74 #define m_data this->m_data
75 
76 // ----------------------------------------------------------------------------
77 
78 /*
79  * Default clamping threshold for Lanczos interpolation. This value has been
80  * selected as the best trade-off for a large set of test linear images.
81  */
82 #ifndef __PCL_LANCZOS_CLAMPING_THRESHOLD
83 #define __PCL_LANCZOS_CLAMPING_THRESHOLD 0.3F
84 #endif
85 
86 // ----------------------------------------------------------------------------
87 
88 /*
89  * Floating point and integer LUT-based interpolations for 3rd, 4th and 5th
90  * order Lanczos functions. LUTs are initialized automatically on-demand by
91  * thread-safe internal routines.
92  *
93  * Real Lanczos LUTs are accurate to +/- 1e-7 DN
94  * Integer Lanczos LUTs are accurate to +/- 1 16-bit DN
95  */
96 #define __PCL_LANCZOS_LUT_REAL_RESOLUTION 4096
97 const double** PCL_FUNC PCL_InitializeLanczosRealLUT( int );
98 #define __PCL_LANCZOS_LUT_INT_RESOLUTION 65535
99 const float* PCL_FUNC PCL_InitializeLanczosIntLUT( int );
100 
101 // ----------------------------------------------------------------------------
102 
103 #define PCL_LANCZOS_ACC() \
104  if ( s < 0 ) \
105  sn -= s, wn -= L; \
106  else \
107  sp += s, wp += L;
108 
148 template <typename T>
150 {
151 private:
152 
153  struct Default
154  {
155  template <typename _T> static bool UseLUT( _T* ) { return false; }
156  static bool UseLUT( uint8* ) { return true; }
157  static bool UseLUT( int8* ) { return true; }
158  static bool UseLUT( uint16* ) { return true; }
159  static bool UseLUT( int16* ) { return true; }
160  static bool UseLUT( float* ) { return true; }
161  };
162 
163 public:
164 
197  LanczosInterpolation( int n = 3, float clamp = __PCL_LANCZOS_CLAMPING_THRESHOLD, bool useLUT = Default::UseLUT( (T*)0 ) ) :
198  m_n( Max( 1, n ) ),
199  m_lut( useLUT ? PCL_InitializeLanczosRealLUT( m_n ) : nullptr ),
200  m_clamp( clamp >= 0 ), m_clampTh( Range( clamp, 0.0F, 1.0F ) ),
201  m_clampThInv( 1 - m_clampTh ),
202  m_Lx( 2*m_n )
203  {
204  PCL_PRECONDITION( n >= 1 )
205  PCL_PRECONDITION( clamp < 0 || 0 <= clamp && clamp <= 1 )
206  }
207 
211  LanczosInterpolation( const LanczosInterpolation& ) = default;
212 
217  {
218  }
219 
226  double operator()( double x, double y ) const override
227  {
228  PCL_PRECONDITION( m_data != nullptr )
229  PCL_PRECONDITION( m_width > 0 && m_height > 0 )
230  PCL_PRECONDITION( x >= 0 && x < m_width )
231  PCL_PRECONDITION( y >= 0 && y < m_height )
232  PCL_CHECK( m_n >= 1 )
233 
234  // Central grid coordinates
235  int x0 = Range( TruncInt( x ), 0, m_width-1 );
236  int y0 = Range( TruncInt( y ), 0, m_height-1 );
237 
238  double sp = 0; // positive filter values
239  double sn = 0; // negative filter values
240  double wp = 0; // positive filter weight
241  double wn = 0; // negative filter weight
242  int i; // row index
243 
244  if ( m_lut != nullptr )
245  {
246  // Discrete interpolation increments
247  int dx = TruncInt( __PCL_LANCZOS_LUT_REAL_RESOLUTION*(x - x0) );
248  int dy = TruncInt( __PCL_LANCZOS_LUT_REAL_RESOLUTION*(y - y0) );
249 
250  // Precalculate horizontal filter values
251  for ( int j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
252  m_Lx[k] = m_lut[k][dx];
253 
254  int k; // LUT node index
255 
256  // Clipped rows at top
257  for ( i = -m_n + 1, k = 0; i <= m_n; ++i, ++k )
258  {
259  int y = y0 + i;
260  if ( y >= 0 )
261  break;
262  if ( m_fillBorder )
263  FillRow( sp, sn, wp, wn, m_lut[k][dy] );
264  else
265  InterpolateRow( sp, sn, wp, wn, m_data - 2*int64( y )*m_width, x0, m_lut[k][dy] );
266  }
267 
268  // Unclipped rows
269  for ( ; i <= m_n; ++i, ++k )
270  {
271  int y = y0 + i;
272  if ( y == m_height )
273  break;
274  InterpolateRow( sp, sn, wp, wn, m_data + int64( y )*m_width, x0, m_lut[k][dy] );
275  }
276 
277  // Clipped rows at bottom
278  for ( ; i <= m_n; ++i, ++k )
279  {
280  if ( m_fillBorder )
281  FillRow( sp, sn, wp, wn, m_lut[k][dy] );
282  else
283  InterpolateRow( sp, sn, wp, wn, m_data + int64( 2*m_height - 2 - y0 - i )*m_width, x0, m_lut[k][dy] );
284  }
285  }
286  else
287  {
288  // Interpolation increments
289  double dx = x - x0;
290  double dy = y - y0;
291 
292  // Precalculate horizontal filter values
293  for ( int j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
294  m_Lx[k] = Lanczos( j - dx );
295 
296  // Clipped rows at top
297  for ( i = -m_n + 1; i <= m_n; ++i )
298  {
299  int y = y0 + i;
300  if ( y >= 0 )
301  break;
302  if ( m_fillBorder )
303  FillRow( sp, sn, wp, wn, Lanczos( i - dy ) );
304  else
305  InterpolateRow( sp, sn, wp, wn, m_data - 2*int64( y )*m_width, x0, Lanczos( i - dy ) );
306  }
307 
308  // Unclipped rows
309  for ( ; i <= m_n; ++i )
310  {
311  int y = y0 + i;
312  if ( y == m_height )
313  break;
314  InterpolateRow( sp, sn, wp, wn, m_data + int64( y )*m_width, x0, Lanczos( i - dy ) );
315  }
316 
317  // Clipped rows at bottom
318  for ( ; i <= m_n; ++i )
319  {
320  if ( m_fillBorder )
321  FillRow( sp, sn, wp, wn, Lanczos( i - dy ) );
322  else
323  InterpolateRow( sp, sn, wp, wn, m_data + int64( 2*m_height - 2 - y0 - i )*m_width, x0, Lanczos( i - dy ) );
324  }
325  }
326 
327  // Clamping
328  if ( m_clamp )
329  {
330  // Empty data?
331  if ( sp == 0 )
332  return 0;
333 
334  // Clamping ratio: s-/s+
335  double r = sn/sp;
336 
337  // Clamp for s- >= s+
338  if ( r >= 1 )
339  return sp/wp;
340 
341  // Clamp for c < s- < s+
342  if ( r > m_clampTh )
343  {
344  r = (r - m_clampTh)/m_clampThInv;
345  double c = 1 - r*r;
346  sn *= c, wn *= c;
347  }
348  }
349 
350  // Weighted convolution
351  return (sp - sn)/(wp - wn);
352  }
353 
360  bool IsClampingEnabled() const
361  {
362  return m_clamp;
363  }
364 
370  void EnableClamping( bool enable = true )
371  {
372  m_clamp = enable;
373  }
374 
380  void DisableClamping( bool disable = true )
381  {
382  EnableClamping( !disable );
383  }
384 
393  float ClampingThreshold() const
394  {
395  return m_clampTh;
396  }
397 
418  void SetClampingThreshold( float clamp )
419  {
420  PCL_PRECONDITION( 0 <= clamp && clamp <= 1 )
421  m_clampTh = Range( clamp, 0.0F, 1.0F );
422  }
423 
424 private:
425 
426  int m_n; // filter order
427  const double** m_lut; // precomputed function values
428  bool m_clamp : 1; // clamping enabled ?
429  double m_clampTh; // clamping threshold in [0,1]
430  double m_clampThInv; // 1 - m_clampTh
431  mutable DVector m_Lx; // precalculated row of function values
432 
433  /*
434  * Sinc function for x > 0
435  */
436  static double Sinc( double x )
437  {
438  x *= Const<double>::pi();
439  return (x > 1.0e-07) ? Sin( x )/x : 1.0;
440  }
441 
442  /*
443  * Evaluate Lanczos function at x.
444  */
445  double Lanczos( double x ) const
446  {
447  if ( x < 0 )
448  x = -x;
449  if ( x < m_n )
450  return Sinc( x ) * Sinc( x/m_n );
451  return 0;
452  }
453 
454  /*
455  * Interpolate a row of pixels.
456  * Can be either an unclipped row or a mirrored border row.
457  */
458  void InterpolateRow( double& sp, double& sn, double& wp, double& wn, const T* f, int x0, double Ly ) const
459  {
460  int j, k;
461 
462  // Clipped pixels at the left border
463  for ( j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
464  {
465  int x = x0 + j;
466  if ( x >= 0 )
467  break;
468  double L = m_Lx[k] * Ly;
469  double s = (m_fillBorder ? m_fillValue : double( f[-x] )) * L;
470  PCL_LANCZOS_ACC()
471  }
472 
473  // Unclipped pixels
474  for ( ; j <= m_n; ++j, ++k )
475  {
476  int x = x0 + j;
477  if ( x == m_width )
478  break;
479  double L = m_Lx[k] * Ly;
480  double s = f[x] * L;
481  PCL_LANCZOS_ACC()
482  }
483 
484  // Clipped pixels at the right border
485  for ( ; j <= m_n; ++j, ++k )
486  {
487  int x = x0 + j;
488  double L = m_Lx[k] * Ly;
489  double s = (m_fillBorder ? m_fillValue : double( f[2*m_width - 2 - x] )) * L;
490  PCL_LANCZOS_ACC()
491  }
492  }
493 
494  /*
495  * Interpolate a clipped pixel row with border filling.
496  */
497  void FillRow( double& sp, double& sn, double& wp, double& wn, double Ly ) const
498  {
499  for ( int j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
500  {
501  double L = m_Lx[k] * Ly;
502  double s = m_fillValue * L;
503  PCL_LANCZOS_ACC()
504  }
505  }
506 };
507 
508 // ----------------------------------------------------------------------------
509 
529 template <typename T, int m_n>
530 class PCL_CLASS LanczosLUTInterpolationBase : public BidimensionalInterpolation<T>
531 {
532 public:
533 
546  LanczosLUTInterpolationBase( float clamp ) :
548  m_lut( PCL_InitializeLanczosIntLUT( m_n ) ),
549  m_clamp( clamp >= 0 ), m_clampTh( Range( clamp, 0.0F, 1.0F ) ),
550  m_clampThInv( 1 - m_clampTh ),
551  m_Lx( 2*m_n ), m_Ly( 2*m_n )
552  {
553  PCL_PRECONDITION( m_n >= 1 )
554  PCL_PRECONDITION( clamp < 0 || 0 <= clamp && clamp <= 1 )
555  PCL_CHECK( m_lut != nullptr )
556  }
557 
561  LanczosLUTInterpolationBase( const LanczosLUTInterpolationBase& ) = default;
562 
566  virtual ~LanczosLUTInterpolationBase()
567  {
568  }
569 
576  double operator()( double x, double y ) const override
577  {
578  PCL_PRECONDITION( m_data != nullptr )
579  PCL_PRECONDITION( m_width > 0 && m_height > 0 )
580  PCL_PRECONDITION( x >= 0 && x < m_width )
581  PCL_PRECONDITION( y >= 0 && y < m_height )
582 
583  // Central grid coordinates
584  int x0 = Range( TruncInt( x ), 0, m_width-1 );
585  int y0 = Range( TruncInt( y ), 0, m_height-1 );
586 
587  // Precalculate function values
588  int dx = RoundInt( (x - x0)*__PCL_LANCZOS_LUT_INT_RESOLUTION );
589  int dy = RoundInt( (y - y0)*__PCL_LANCZOS_LUT_INT_RESOLUTION );
590  for ( int j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
591  {
592  int d0 = j*__PCL_LANCZOS_LUT_INT_RESOLUTION;
593  m_Lx[k] = m_lut[Abs( d0 - dx )];
594  m_Ly[k] = m_lut[Abs( d0 - dy )];
595  }
596 
597  double sp = 0; // positive filter values
598  double sn = 0; // negative filter values
599  double wp = 0; // positive filter weight
600  double wn = 0; // negative filter weight
601  int i, k; // row and coefficient indices
602 
603  // Clipped rows at top
604  for ( i = -m_n + 1, k = 0; i <= m_n; ++i, ++k )
605  {
606  int y = y0 + i;
607  if ( y >= 0 )
608  break;
609  if ( m_fillBorder )
610  FillRow( sp, sn, wp, wn, m_Ly[k] );
611  else
612  InterpolateRow( sp, sn, wp, wn, m_data - 2*int64( y )*m_width, x0, m_Ly[k] );
613  }
614 
615  // Unclipped rows
616  for ( ; i <= m_n; ++i, ++k )
617  {
618  int y = y0 + i;
619  if ( y == m_height )
620  break;
621  InterpolateRow( sp, sn, wp, wn, m_data + int64( y )*m_width, x0, m_Ly[k] );
622  }
623 
624  // Clipped rows at bottom
625  for ( ; i <= m_n; ++i, ++k )
626  {
627  if ( m_fillBorder )
628  FillRow( sp, sn, wp, wn, m_Ly[k] );
629  else
630  InterpolateRow( sp, sn, wp, wn, m_data + int64( 2*m_height - 2 - y0 - i )*m_width, x0, m_Ly[k] );
631  }
632 
633  // Clamping
634  if ( m_clamp )
635  {
636  // Empty data?
637  if ( sp == 0 )
638  return 0;
639 
640  // Clamping ratio: s-/s+
641  double r = sn/sp;
642 
643  // Clamp for s- >= s+
644  if ( r >= 1 )
645  return sp/wp;
646 
647  // Clamp for c < s- < s+
648  if ( r > m_clampTh )
649  {
650  r = (r - m_clampTh)/m_clampThInv;
651  double c = 1 - r*r;
652  sn *= c, wn *= c;
653  }
654  }
655 
656  // Weighted convolution
657  return (sp - sn)/(wp - wn);
658  }
659 
666  bool IsClampingEnabled() const
667  {
668  return m_clamp;
669  }
670 
676  void EnableClamping( bool enable = true )
677  {
678  m_clamp = enable;
679  }
680 
686  void DisableClamping( bool disable = true )
687  {
688  EnableClamping( !disable );
689  }
690 
699  float ClampingThreshold() const
700  {
701  return m_clampTh;
702  }
703 
724  void SetClampingThreshold( float clamp )
725  {
726  PCL_PRECONDITION( 0 <= clamp && clamp <= 1 )
727  m_clampTh = Range( clamp, 0.0F, 1.0F );
728  }
729 
730 private:
731 
732  const float* m_lut; // filter LUT
733  bool m_clamp : 1; // clamping enabled ?
734  double m_clampTh; // clamping threshold in [0,1]
735  double m_clampThInv; // 1 - m_clampTh
736  mutable FVector m_Lx, m_Ly; // precalculated function values
737 
738  /*
739  * Interpolate a row of pixels.
740  * Can be either an unclipped row or a mirrored border row.
741  */
742  void InterpolateRow( double& sp, double& sn, double& wp, double& wn, const T* f, int x0, float Ly ) const
743  {
744  int j, k;
745 
746  // Clipped pixels at the left border
747  for ( j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
748  {
749  int x = x0 + j;
750  if ( x >= 0 )
751  break;
752  double L = m_Lx[k] * Ly;
753  double s = (m_fillBorder ? m_fillValue : double( f[-x] )) * L;
754  PCL_LANCZOS_ACC()
755  }
756 
757  // Unclipped pixels
758  for ( ; j <= m_n; ++j, ++k )
759  {
760  int x = x0 + j;
761  if ( x == m_width )
762  break;
763  double L = m_Lx[k] * Ly;
764  double s = f[x] * L;
765  PCL_LANCZOS_ACC()
766  }
767 
768  // Clipped pixels at the right border
769  for ( ; j <= m_n; ++j, ++k )
770  {
771  int x = x0 + j;
772  double L = m_Lx[k] * Ly;
773  double s = (m_fillBorder ? m_fillValue : double( f[2*m_width - 2 - x] )) * L;
774  PCL_LANCZOS_ACC()
775  }
776  }
777 
778  /*
779  * Interpolate a clipped pixel row with border filling.
780  */
781  void FillRow( double& sp, double& sn, double& wp, double& wn, float Ly ) const
782  {
783  for ( int j = -m_n + 1, k = 0; j <= m_n; ++j, ++k )
784  {
785  double L = m_Lx[k] * Ly;
786  double s = m_fillValue * L;
787  PCL_LANCZOS_ACC()
788  }
789  }
790 };
791 
792 // ----------------------------------------------------------------------------
793 
811 template <typename T>
812 class PCL_CLASS Lanczos3LUTInterpolation : public LanczosLUTInterpolationBase<T,3>
813 {
814 public:
815 
828  Lanczos3LUTInterpolation( float clamp = __PCL_LANCZOS_CLAMPING_THRESHOLD ) :
829  LanczosLUTInterpolationBase<T,3>( clamp )
830  {
831  PCL_PRECONDITION( 0 <= clamp && clamp <= 1 )
832  }
833 
838 
843  {
844  }
845 };
846 
847 // ----------------------------------------------------------------------------
848 
866 template <typename T>
867 class PCL_CLASS Lanczos4LUTInterpolation : public LanczosLUTInterpolationBase<T,4>
868 {
869 public:
870 
883  Lanczos4LUTInterpolation( float clamp = __PCL_LANCZOS_CLAMPING_THRESHOLD ) :
884  LanczosLUTInterpolationBase<T,4>( clamp )
885  {
886  PCL_PRECONDITION( 0 <= clamp && clamp <= 1 )
887  }
888 
893 
898  {
899  }
900 };
901 
902 // ----------------------------------------------------------------------------
903 
921 template <typename T>
922 class PCL_CLASS Lanczos5LUTInterpolation : public LanczosLUTInterpolationBase<T,5>
923 {
924 public:
925 
938  Lanczos5LUTInterpolation( float clamp = __PCL_LANCZOS_CLAMPING_THRESHOLD ) :
939  LanczosLUTInterpolationBase<T,5>( clamp )
940  {
941  PCL_PRECONDITION( 0 <= clamp && clamp <= 1 )
942  }
943 
948 
953  {
954  }
955 };
956 
957 // ----------------------------------------------------------------------------
958 
959 #undef PCL_LANCZOS_ACC
960 
961 #undef m_width
962 #undef m_height
963 #undef m_fillBorder
964 #undef m_fillValue
965 #undef m_data
966 
967 // ----------------------------------------------------------------------------
968 
969 } // pcl
970 
971 #endif // __PCL_LanczosInterpolation_h
972 
973 // ----------------------------------------------------------------------------
974 // EOF pcl/LanczosInterpolation.h - Released 2019-11-07T10:59:34Z
unsigned char uint8
Definition: Defs.h:576
Lanczos5LUTInterpolation(float clamp=__PCL_LANCZOS_CLAMPING_THRESHOLD)
static constexpr T pi()
Definition: Constants.h:94
int RoundInt(T x)
Definition: Math.h:1397
Two dimensional LUT-based 3rd-order Lanczos interpolation algorithm.
signed short int16
Definition: Defs.h:582
PCL root namespace.
Definition: AbstractImage.h:76
A generic interface to two-dimensional interpolation algorithms.
Two dimensional Lanczos interpolation algorithm.
Lanczos4LUTInterpolation(float clamp=__PCL_LANCZOS_CLAMPING_THRESHOLD)
constexpr const T & Range(const T &x, const T &a, const T &b)
Definition: Utility.h:190
T Abs(const Complex< T > &c)
Definition: Complex.h:420
Two dimensional LUT-based 4th-order Lanczos interpolation algorithm.
signed char int8
Definition: Defs.h:570
constexpr const T & Max(const T &a, const T &b)
Definition: Utility.h:119
64-bit floating point real vector.
unsigned short uint16
Definition: Defs.h:588
int TruncInt(T x)
Definition: Math.h:1035
Two dimensional LUT-based 5th-order Lanczos interpolation algorithm.
void EnableClamping(bool enable=true)
LanczosInterpolation(int n=3, float clamp=__PCL_LANCZOS_CLAMPING_THRESHOLD, bool useLUT=Default::UseLUT((T *) 0))
double operator()(double x, double y) const override
Complex< T > Sin(const Complex< T > &c)
Definition: Complex.h:786
void DisableClamping(bool disable=true)
void SetClampingThreshold(float clamp)
32-bit floating point real vector.
signed long long int64
Definition: Defs.h:610
Lanczos3LUTInterpolation(float clamp=__PCL_LANCZOS_CLAMPING_THRESHOLD)