GteIntpTrilinear3.h
Go to the documentation of this file.
1 // David Eberly, Geometric Tools, Redmond WA 98052
2 // Copyright (c) 1998-2017
3 // Distributed under the Boost Software License, Version 1.0.
4 // http://www.boost.org/LICENSE_1_0.txt
5 // http://www.geometrictools.com/License/Boost/LICENSE_1_0.txt
6 // File Version: 3.0.0 (2016/06/19)
7 
8 #pragma once
9 
10 #include <LowLevel/GteLogger.h>
11 
12 // The interpolator is for uniformly spaced(x,y z)-values. The input samples
13 // must be stored in lexicographical order to represent f(x,y,z); that is,
14 // F[c + xBound*(r + yBound*s)] corresponds to f(x,y,z), where c is the index
15 // corresponding to x, r is the index corresponding to y, and s is the index
16 // corresponding to z.
17 
18 namespace gte
19 {
20 
21 template <typename Real>
23 {
24 public:
25  // Construction.
26  IntpTrilinear3(int xBound, int yBound, int zBound, Real xMin,
27  Real xSpacing, Real yMin, Real ySpacing, Real zMin, Real zSpacing,
28  Real const* F);
29 
30  // Member access.
31  inline int GetXBound() const;
32  inline int GetYBound() const;
33  inline int GetZBound() const;
34  inline int GetQuantity() const;
35  inline Real const* GetF() const;
36  inline Real GetXMin() const;
37  inline Real GetXMax() const;
38  inline Real GetXSpacing() const;
39  inline Real GetYMin() const;
40  inline Real GetYMax() const;
41  inline Real GetYSpacing() const;
42  inline Real GetZMin() const;
43  inline Real GetZMax() const;
44  inline Real GetZSpacing() const;
45 
46  // Evaluate the function and its derivatives. The functions clamp the
47  // inputs to xmin <= x <= xmax, ymin <= y <= ymax, and zmin <= z <= zmax.
48  // The first operator is for function evaluation. The second operator is
49  // for function or derivative evaluations. The xOrder argument is the
50  // order of the x-derivative, the yOrder argument is the order of the
51  // y-derivative, and the zOrder argument is the order of the z-derivative.
52  // All orders are zero to get the function value itself.
53  Real operator()(Real x, Real y, Real z) const;
54  Real operator()(int xOrder, int yOrder, int zOrder, Real x, Real y,
55  Real z) const;
56 
57 private:
62  Real const* mF;
63  Real mBlend[2][2];
64 };
65 
66 
67 template <typename Real>
68 IntpTrilinear3<Real>::IntpTrilinear3(int xBound, int yBound, int zBound,
69  Real xMin, Real xSpacing, Real yMin, Real ySpacing, Real zMin,
70  Real zSpacing, Real const* F)
71  :
72  mXBound(xBound),
73  mYBound(yBound),
74  mZBound(zBound),
75  mQuantity(xBound * yBound * zBound),
76  mXMin(xMin),
77  mXSpacing(xSpacing),
78  mYMin(yMin),
79  mYSpacing(ySpacing),
80  mZMin(zMin),
81  mZSpacing(zSpacing),
82  mF(F)
83 {
84  // At least a 2x2x2 block of data points are needed to construct the
85  // trilinear interpolation.
86  LogAssert(mXBound >= 2 && mYBound >= 2 && mZBound >= 2 && mF,
87  "Invalid input.");
88  LogAssert(mXSpacing > (Real)0 && mYSpacing > (Real)0 &&
89  mZSpacing > (Real)0, "Invalid input.");
90 
91  mXMax = mXMin + mXSpacing * static_cast<Real>(mXBound - 1);
92  mInvXSpacing = ((Real)1) / mXSpacing;
93  mYMax = mYMin + mYSpacing * static_cast<Real>(mYBound - 1);
94  mInvYSpacing = ((Real)1) / mYSpacing;
95  mZMax = mZMin + mZSpacing * static_cast<Real>(mZBound - 1);
96  mInvZSpacing = ((Real)1) / mZSpacing;
97 
98  mBlend[0][0] = (Real)1;
99  mBlend[0][1] = (Real)-1;
100  mBlend[1][0] = (Real)0;
101  mBlend[1][1] = (Real)1;
102 }
103 
104 template <typename Real> inline
106 {
107  return mXBound;
108 }
109 
110 template <typename Real> inline
112 {
113  return mYBound;
114 }
115 
116 template <typename Real> inline
118 {
119  return mZBound;
120 }
121 
122 template <typename Real> inline
124 {
125  return mQuantity;
126 }
127 
128 template <typename Real> inline
129 Real const* IntpTrilinear3<Real>::GetF() const
130 {
131  return mF;
132 }
133 
134 template <typename Real> inline
136 {
137  return mXMin;
138 }
139 
140 template <typename Real> inline
142 {
143  return mXMax;
144 }
145 
146 template <typename Real> inline
148 {
149  return mXSpacing;
150 }
151 
152 template <typename Real> inline
154 {
155  return mYMin;
156 }
157 
158 template <typename Real> inline
160 {
161  return mYMax;
162 }
163 
164 template <typename Real> inline
166 {
167  return mYSpacing;
168 }
169 
170 template <typename Real> inline
172 {
173  return mZMin;
174 }
175 
176 template <typename Real> inline
178 {
179  return mZMax;
180 }
181 
182 template <typename Real> inline
184 {
185  return mZSpacing;
186 }
187 
188 template <typename Real>
189 Real IntpTrilinear3<Real>::operator()(Real x, Real y, Real z) const
190 {
191  // Compute x-index and clamp to image.
192  Real xIndex = (x - mXMin) * mInvXSpacing;
193  int ix = static_cast<int>(xIndex);
194  if (ix < 0)
195  {
196  ix = 0;
197  }
198  else if (ix >= mXBound)
199  {
200  ix = mXBound - 1;
201  }
202 
203  // Compute y-index and clamp to image.
204  Real yIndex = (y - mYMin) * mInvYSpacing;
205  int iy = static_cast<int>(yIndex);
206  if (iy < 0)
207  {
208  iy = 0;
209  }
210  else if (iy >= mYBound)
211  {
212  iy = mYBound - 1;
213  }
214 
215  // Compute z-index and clamp to image.
216  Real zIndex = (z - mZMin) * mInvZSpacing;
217  int iz = static_cast<int>(zIndex);
218  if (iz < 0)
219  {
220  iz = 0;
221  }
222  else if (iz >= mZBound)
223  {
224  iz = mZBound - 1;
225  }
226 
227  Real U[2];
228  U[0] = (Real)1;
229  U[1] = xIndex - ix;
230 
231  Real V[2];
232  V[0] = (Real)1;
233  V[1] = yIndex - iy;
234 
235  Real W[2];
236  W[0] = (Real)1;
237  W[1] = zIndex - iz;
238 
239  // Compute P = M*U, Q = M*V, R = M*W.
240  Real P[2], Q[2], R[2];
241  for (int row = 0; row < 2; ++row)
242  {
243  P[row] = (Real)0;
244  Q[row] = (Real)0;
245  R[row] = (Real)0;
246  for (int col = 0; col < 2; ++col)
247  {
248  P[row] += mBlend[row][col] * U[col];
249  Q[row] += mBlend[row][col] * V[col];
250  R[row] += mBlend[row][col] * W[col];
251  }
252  }
253 
254  // compute the tensor product (M*U)(M*V)(M*W)*D where D is the 2x2x2
255  // subimage containing (x,y,z)
256  Real result = (Real)0;
257  for (int slice = 0; slice < 2; ++slice)
258  {
259  int zClamp = iz + slice;
260  if (zClamp >= mZBound)
261  {
262  zClamp = mZBound - 1;
263  }
264 
265  for (int row = 0; row < 2; ++row)
266  {
267  int yClamp = iy + row;
268  if (yClamp >= mYBound)
269  {
270  yClamp = mYBound - 1;
271  }
272 
273  for (int col = 0; col < 2; ++col)
274  {
275  int xClamp = ix + col;
276  if (xClamp >= mXBound)
277  {
278  xClamp = mXBound - 1;
279  }
280 
281  result += P[col] * Q[row] * R[slice] *
282  mF[xClamp + mXBound * (yClamp + mYBound * zClamp)];
283  }
284  }
285  }
286 
287  return result;
288 }
289 
290 template <typename Real>
291 Real IntpTrilinear3<Real>::operator()(int xOrder, int yOrder, int zOrder,
292  Real x, Real y, Real z) const
293 {
294  // Compute x-index and clamp to image.
295  Real xIndex = (x - mXMin) * mInvXSpacing;
296  int ix = static_cast<int>(xIndex);
297  if (ix < 0)
298  {
299  ix = 0;
300  }
301  else if (ix >= mXBound)
302  {
303  ix = mXBound - 1;
304  }
305 
306  // Compute y-index and clamp to image.
307  Real yIndex = (y - mYMin) * mInvYSpacing;
308  int iy = static_cast<int>(yIndex);
309  if (iy < 0)
310  {
311  iy = 0;
312  }
313  else if (iy >= mYBound)
314  {
315  iy = mYBound - 1;
316  }
317 
318  // Compute z-index and clamp to image.
319  Real zIndex = (z - mZMin) * mInvZSpacing;
320  int iz = static_cast<int>(zIndex);
321  if (iz < 0)
322  {
323  iz = 0;
324  }
325  else if (iz >= mZBound)
326  {
327  iz = mZBound - 1;
328  }
329 
330  Real U[2], dx, xMult;
331  switch (xOrder)
332  {
333  case 0:
334  dx = xIndex - ix;
335  U[0] = (Real)1;
336  U[1] = dx;
337  xMult = (Real)1;
338  break;
339  case 1:
340  dx = xIndex - ix;
341  U[0] = (Real)0;
342  U[1] = (Real)1;
343  xMult = mInvXSpacing;
344  break;
345  default:
346  return (Real)0;
347  }
348 
349  Real V[2], dy, yMult;
350  switch (yOrder)
351  {
352  case 0:
353  dy = yIndex - iy;
354  V[0] = (Real)1;
355  V[1] = dy;
356  yMult = (Real)1;
357  break;
358  case 1:
359  dy = yIndex - iy;
360  V[0] = (Real)0;
361  V[1] = (Real)1;
362  yMult = mInvYSpacing;
363  break;
364  default:
365  return (Real)0;
366  }
367 
368  Real W[2], dz, zMult;
369  switch (zOrder)
370  {
371  case 0:
372  dz = zIndex - iz;
373  W[0] = (Real)1;
374  W[1] = dz;
375  zMult = (Real)1;
376  break;
377  case 1:
378  dz = zIndex - iz;
379  W[0] = (Real)0;
380  W[1] = (Real)1;
381  zMult = mInvZSpacing;
382  break;
383  default:
384  return (Real)0;
385  }
386 
387  // Compute P = M*U, Q = M*V, and R = M*W.
388  Real P[2], Q[2], R[2];
389  for (int row = 0; row < 2; ++row)
390  {
391  P[row] = (Real)0;
392  Q[row] = (Real)0;
393  R[row] = (Real)0;
394  for (int col = 0; col < 2; ++col)
395  {
396  P[row] += mBlend[row][col] * U[col];
397  Q[row] += mBlend[row][col] * V[col];
398  R[row] += mBlend[row][col] * W[col];
399  }
400  }
401 
402  // Compute the tensor product (M*U)(M*V)(M*W)*D where D is the 2x2x2
403  // subimage containing (x,y,z).
404  Real result = (Real)0;
405  for (int slice = 0; slice < 2; ++slice)
406  {
407  int zClamp = iz + slice;
408  if (zClamp >= mZBound)
409  {
410  zClamp = mZBound - 1;
411  }
412 
413  for (int row = 0; row < 2; ++row)
414  {
415  int yClamp = iy + row;
416  if (yClamp >= mYBound)
417  {
418  yClamp = mYBound - 1;
419  }
420 
421  for (int col = 0; col < 2; ++col)
422  {
423  int xClamp = ix + col;
424  if (xClamp >= mXBound)
425  {
426  xClamp = mXBound - 1;
427  }
428 
429  result += P[col] * Q[row] * R[slice] *
430  mF[xClamp + mXBound * (yClamp + mYBound * zClamp)];
431  }
432  }
433  }
434  result *= xMult * yMult * zMult;
435 
436  return result;
437 }
438 
439 
440 }
#define LogAssert(condition, message)
Definition: GteLogger.h:86
GLint GLenum GLint x
Definition: glcorearb.h:404
GLenum GLenum GLsizei void * row
Definition: glext.h:2725
GLdouble GLdouble GLdouble z
Definition: glcorearb.h:843
Real GetXSpacing() const
Real const * GetF() const
Real GetYSpacing() const
Real operator()(Real x, Real y, Real z) const
GLuint64EXT * result
Definition: glext.h:10003
GLint y
Definition: glcorearb.h:98
Real GetZSpacing() const
IntpTrilinear3(int xBound, int yBound, int zBound, Real xMin, Real xSpacing, Real yMin, Real ySpacing, Real zMin, Real zSpacing, Real const *F)


geometric_tools_engine
Author(s): Yijiang Huang
autogenerated on Thu Jul 18 2019 04:00:00