PCL
BilinearInterpolation.h
Go to the documentation of this file.
1 // ____ ______ __
2 // / __ \ / ____// /
3 // / /_/ // / / /
4 // / ____// /___ / /___ PixInsight Class Library
5 // /_/ \____//_____/ PCL 2.7.0
6 // ----------------------------------------------------------------------------
7 // pcl/BilinearInterpolation.h - Released 2024-06-18T15:48:54Z
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-2024 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 (https://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_BilinearInterpolation_h
53 #define __PCL_BilinearInterpolation_h
54 
56 
57 #include <pcl/Defs.h>
58 #include <pcl/Diagnostics.h>
59 
61 #include <pcl/Utility.h>
62 
63 namespace pcl
64 {
65 
66 // ----------------------------------------------------------------------------
67 
68 #define m_width this->m_width
69 #define m_height this->m_height
70 #define m_fillBorder this->m_fillBorder
71 #define m_fillValue this->m_fillValue
72 #define m_data this->m_data
73 
74 // ----------------------------------------------------------------------------
75 
87 template <typename T>
89 {
90 public:
91 
95  BilinearInterpolation() = default;
96 
101 
106  {
107  }
108 
114  double operator()( double x, double y ) const override
115  {
116  PCL_PRECONDITION( m_data != nullptr )
117  PCL_PRECONDITION( m_width > 0 && m_height > 0 )
118 
119  int j0 = pcl::Range( TruncInt( x ), 0, m_width-1 );
120  int i0 = pcl::Range( TruncInt( y ), 0, m_height-1 );
121 
122  int j1 = j0 + 1;
123  int i1 = i0 + 1;
124 
125  double p00, p10, p01, p11;
126  const T* fp = m_data + (int64( i0 )*m_width + j0);
127 
128  p00 = *fp;
129  p10 = (j1 < m_width) ? fp[1] : (m_fillBorder ? m_fillValue : *fp);
130 
131  if ( i1 < m_height )
132  fp += m_width;
133  else if ( m_fillBorder )
134  {
135  p01 = p11 = m_fillValue;
136  goto __1;
137  }
138 
139  p01 = *fp;
140  p11 = (j1 < m_width) ? fp[1] : (m_fillBorder ? m_fillValue : *fp);
141 
142 __1:
143  double dx = x - j0, dx1 = 1 - dx;
144  double dy = y - i0, dy1 = 1 - dy;
145  return p00*dx1*dy1 + p10*dx*dy1 + p01*dx1*dy + p11*dx*dy;
146  }
147 };
148 
149 // ----------------------------------------------------------------------------
150 
151 #undef m_width
152 #undef m_height
153 #undef m_fillBorder
154 #undef m_fillValue
155 #undef m_data
156 
157 // ----------------------------------------------------------------------------
158 
159 } // pcl
160 
161 #endif // __PCL_BilinearInterpolation_h
162 
163 // ----------------------------------------------------------------------------
164 // EOF pcl/BilinearInterpolation.h - Released 2024-06-18T15:48:54Z
A generic interface to two-dimensional interpolation algorithms.
Bilinear interpolation algorithm.
double operator()(double x, double y) const override
BilinearInterpolation(const BilinearInterpolation &)=default
int TruncInt(T x) noexcept
Definition: Math.h:1132
signed long long int64
Definition: Defs.h:676
constexpr const T & Range(const T &x, const T &a, const T &b) noexcept
Definition: Utility.h:190
PCL root namespace.
Definition: AbstractImage.h:77