OpenVDB  2.1.0
Ray.h
Go to the documentation of this file.
1 //
3 // Copyright (c) 2012-2013 DreamWorks Animation LLC
4 //
5 // All rights reserved. This software is distributed under the
6 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
7 //
8 // Redistributions of source code must retain the above copyright
9 // and license notice and the following restrictions and disclaimer.
10 //
11 // * Neither the name of DreamWorks Animation nor the names of
12 // its contributors may be used to endorse or promote products derived
13 // from this software without specific prior written permission.
14 //
15 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
16 // "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
17 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
18 // A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
19 // OWNER OR CONTRIBUTORS BE LIABLE FOR ANY INDIRECT, INCIDENTAL,
20 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
21 // LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
22 // DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
23 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
25 // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26 // IN NO EVENT SHALL THE COPYRIGHT HOLDERS' AND CONTRIBUTORS' AGGREGATE
27 // LIABILITY FOR ALL CLAIMS REGARDLESS OF THEIR BASIS EXCEED US$250.00.
28 //
30 //
36 
37 #ifndef OPENVDB_MATH_RAY_HAS_BEEN_INCLUDED
38 #define OPENVDB_MATH_RAY_HAS_BEEN_INCLUDED
39 
40 #include "Math.h"
41 #include "Vec3.h"
42 #include "Transform.h"
43 #include <iostream> // for std::ostream
44 #include <boost/type_traits/is_floating_point.hpp>
45 #include <limits>// for std::numeric_limits<Type>::max()
46 
47 namespace openvdb {
49 namespace OPENVDB_VERSION_NAME {
50 namespace math {
51 
52 template<typename RealT = double>
53 class Ray
54 {
55 public:
56  BOOST_STATIC_ASSERT(boost::is_floating_point<RealT>::value);
57  typedef RealT RealType;
59  typedef Vec3Type Vec3T;
60 
61  Ray(const Vec3Type& eye = Vec3Type(0,0,0),
62  const Vec3Type& direction = Vec3Type(1,0,0),
63  RealT t0 = math::Delta<RealT>::value(),
65  : mEye(eye), mDir(direction), mInvDir(1/mDir), mT0(t0), mT1(t1)
66  {
67  }
68 
69  inline void setEye(const Vec3Type& eye) { mEye = eye; }
70 
71  inline void setDir(const Vec3Type& dir)
72  {
73  mDir = dir;
74  mInvDir = 1/mDir;
75  }
76 
77  inline void setMinTime(RealT t0) { assert(t0>0); mT0 = t0; }
78 
79  inline void setMaxTime(RealT t1) { assert(t1>0); mT1 = t1; }
80 
81  inline void setTimes(RealT t0, RealT t1) { assert(t0>0 && t1>0);mT0 = t0; mT1 = t1; }
82 
83  inline void scaleTimes(RealT scale) { assert(scale>0); mT0 *= scale; mT1 *= scale; }
84 
85  inline void reset(const Vec3Type& eye, const Vec3Type& direction,
86  RealT t0 = 0, RealT t1 = std::numeric_limits<RealT>::max())
87  {
88  this->setEye(eye);
89  this->setDir(direction);
90  this->setTimes(t0, t1);
91  }
92 
93  inline const Vec3T& eye() const {return mEye;}
94 
95  inline const Vec3T& dir() const {return mDir;}
96 
97  inline const Vec3T& invDir() const {return mInvDir;}
98 
99  inline RealT t0() const {return mT0;}
100 
101  inline RealT t1() const {return mT1;}
102 
104  inline Vec3R operator()(RealT time) const { return mEye + mDir * time; }
105 
107  inline Vec3R start() const { return (*this)(mT0); }
108 
110  inline Vec3R end() const { return (*this)(mT1); }
111 
113  inline Vec3R mid() const { return (*this)(0.5*(mT0+mT1)); }
114 
116  inline bool test() const { return (mT0 < mT1); }
117 
119  inline bool test(RealT time) const { return (time>=mT0 && time<=mT1); }
120 
127  template<typename MapType>
128  inline Ray applyMap(const MapType& map) const
129  {
130  assert(map.isLinear());
131  assert(math::isApproxEqual(mDir.length(), RealT(1)));
132  const Vec3T eye = map.applyMap(mEye);
133  const Vec3T dir = map.applyJacobian(mDir);
134  const RealT length = dir.length();
135  return Ray(eye, dir/length, length*mT0, length*mT1);
136  }
137 
144  template<typename MapType>
145  inline Ray applyInverseMap(const MapType& map) const
146  {
147  assert(map.isLinear());
148  assert(math::isApproxEqual(mDir.length(), RealT(1)));
149  const Vec3T eye = map.applyInverseMap(mEye);
150  const Vec3T dir = map.applyInverseJacobian(mDir);
151  const RealT length = dir.length();
152  return Ray(eye, dir/length, length*mT0, length*mT1);
153  }
154 
157  template<typename GridType>
158  inline Ray indexToWorld(const GridType& grid) const
159  {
160  return this->applyMap(*(grid.transform().baseMap()));
161  }
162 
165  template<typename GridType>
166  inline Ray worldToIndex(const GridType& grid) const
167  {
168  return this->applyInverseMap(*(grid.transform().baseMap()));
169  }
170 
178  inline bool intersects(const Vec3T& center, RealT radius, RealT& t0, RealT& t1) const
179  {
180  const Vec3T origin = mEye - center;
181  const RealT A = mDir.lengthSqr();
182  const RealT B = 2 * mDir.dot(origin);
183  const RealT C = origin.lengthSqr() - radius * radius;
184  const RealT D = B * B - 4 * A * C;
185 
186  if (D < 0) return false;
187 
188  const RealT Q = RealT(-0.5)*(B<0 ? (B + Sqrt(D)) : (B - Sqrt(D)));
189 
190  t0 = Q / A;
191  t1 = C / Q;
192 
193  if (t0 > t1) std::swap(t0, t1);
194  if (t0 < mT0) t0 = mT0;
195  if (t1 > mT1) t1 = mT1;
196  return t0 <= t1;
197  }
198 
202  inline bool intersects(const Vec3T& center, RealT radius) const
203  {
204  RealT t0, t1;
205  return this->intersects(center, radius, t0, t1)>0;
206  }
207 
212  inline bool clip(const Vec3T& center, RealT radius)
213  {
214  RealT t0, t1;
215  const bool hit = this->intersects(center, radius, t0, t1);
216  if (hit) {
217  mT0 = t0;
218  mT1 = t1;
219  }
220  return hit;
221  }
222 
230  template<typename BBoxT>
231  inline bool intersects(const BBoxT& bbox, RealT& t0, RealT& t1) const
232  {
233  t0 = mT0;
234  t1 = mT1;
235  for (size_t i = 0; i < 3; ++i) {
236  RealT a = (bbox.min()[i] - mEye[i]) * mInvDir[i];
237  RealT b = (bbox.max()[i] - mEye[i]) * mInvDir[i];
238  if (a > b) std::swap(a, b);
239  if (a > t0) t0 = a;
240  if (b < t1) t1 = b;
241  if (t0 > t1) return false;
242  }
243  return true;
244  }
245 
248  template<typename BBoxT>
249  inline bool intersects(const BBoxT& bbox) const
250  {
251  RealT t0, t1;
252  return this->intersects(bbox, t0, t1);
253  }
254 
258  template<typename BBoxT>
259  inline bool clip(const BBoxT& bbox)
260  {
261  RealT t0, t1;
262  const bool hit = this->intersects(bbox, t0, t1);
263  if (hit) {
264  mT0 = t0;
265  mT1 = t1;
266  }
267  return hit;
268  }
269 
275  inline bool intersects(const Vec3T& normal, RealT distance, RealT& t) const
276  {
277  const RealT cosAngle = mDir.dot(normal);
278  if (math::isApproxZero(cosAngle)) return false;//parallel
279  t = (distance - mEye.dot(normal))/cosAngle;
280  return this->test(t);
281  }
282 
288  inline bool intersects(const Vec3T& normal, const Vec3T& point, RealT& t) const
289  {
290  return this->intersects(normal, point.dot(normal), t);
291  }
292 
293 private:
294  Vec3T mEye, mDir, mInvDir;
295  RealT mT0, mT1;
296 }; // end of Ray class
297 
300 template<typename RealT>
301 inline std::ostream& operator<<(std::ostream& os, const Ray<RealT>& r)
302 {
303  os << "eye=" << r.eye() << " dir=" << r.dir() << " 1/dir="<<r.invDir()
304  << " t0=" << r.t0() << " t1=" << r.t1();
305  return os;
306 }
307 
308 
310 
311 
320 template<typename RayT, Index Log2Dim = 0>
321 class DDA
322 {
323 public:
324  typedef typename RayT::RealType RealType;
325  typedef RealType RealT;
326  typedef typename RayT::Vec3Type Vec3Type;
327  typedef Vec3Type Vec3T;
328 
329  DDA(const RayT& ray) { this->init(ray, ray.t0(), ray.t1()); }
330 
331  DDA(const RayT& ray, RealT startTime) { this->init(ray, startTime, ray.t1()); }
332 
333  DDA(const RayT& ray, RealT startTime, RealT maxTime) { this->init(ray, startTime, maxTime); }
334 
335  inline void init(const RayT& ray, RealT startTime, RealT maxTime)
336  {
337  assert(startTime <= maxTime);
338  static const int DIM = 1 << Log2Dim;
339  mT0 = startTime;
340  mT1 = maxTime;
341  const Vec3T &pos = ray(mT0), &dir = ray.dir(), &inv = ray.invDir();
342  mVoxel = Coord::floor(pos) & (~(DIM-1));
343  for (size_t axis = 0; axis < 3; ++axis) {
344  if (math::isZero(dir[axis])) {//handles dir = +/- 0
345  mStep[axis] = 0;//dummy value
346  mNext[axis] = std::numeric_limits<RealT>::max();//i.e. disabled!
347  mDelta[axis] = std::numeric_limits<RealT>::max();//dummy value
348  } else if (inv[axis] > 0) {
349  mStep[axis] = DIM;
350  mNext[axis] = mT0 + (mVoxel[axis] + DIM - pos[axis]) * inv[axis];
351  mDelta[axis] = mStep[axis] * inv[axis];
352  } else {
353  mStep[axis] = -DIM;
354  mNext[axis] = mT0 + (mVoxel[axis] - pos[axis]) * inv[axis];
355  mDelta[axis] = mStep[axis] * inv[axis];
356  }
357  }
358  }
359 
362  inline bool step()
363  {
364  const size_t stepAxis = math::MinIndex(mNext);
365  mT0 = mNext[stepAxis];
366  mNext[stepAxis] += mDelta[stepAxis];
367  mVoxel[stepAxis] += mStep[stepAxis];
368  return mT0 <= mT1;
369  }
370 
376  inline const Coord& voxel() const { return mVoxel; }
377 
383  inline RealType time() const { return mT0; }
384 
388  inline RealType next() const { return math::Min(mT1, mNext[0], mNext[1], mNext[2]); }
389 
392  void print(std::ostream& os = std::cout) const
393  {
394  os << "Dim=" << (1<<Log2Dim) << " time=" << mT0 << " next()="
395  << this->next() << " voxel=" << mVoxel << " next=" << mNext
396  << " delta=" << mDelta << " step=" << mStep << std::endl;
397  }
398 
399 private:
400  RealT mT0, mT1;
401  Coord mVoxel, mStep;
402  Vec3T mDelta, mNext;
403 }; // class DDA
404 
407 template<typename RayT, Index Log2Dim>
408 inline std::ostream& operator<<(std::ostream& os, const DDA<RayT, Log2Dim>& dda)
409 {
410  os << "Dim=" << (1<<Log2Dim) << " time=" << dda.time()
411  << " next()=" << dda.next() << " voxel=" << dda.voxel();
412  return os;
413 }
414 
415 } // namespace math
416 } // namespace OPENVDB_VERSION_NAME
417 } // namespace openvdb
418 
419 #endif // OPENVDB_MATH_RAY_HAS_BEEN_INCLUDED
420 
421 // Copyright (c) 2012-2013 DreamWorks Animation LLC
422 // All rights reserved. This software is distributed under the
423 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
const Vec3T & dir() const
Definition: Ray.h:95
Vec3R end() const
Return the endpoint of the ray.
Definition: Ray.h:110
bool intersects(const BBoxT &bbox) const
Return true if this ray intersects the specified bounding box.
Definition: Ray.h:249
Ray worldToIndex(const GridType &grid) const
Return a new ray in the index space of the specified grid, assuming the existing ray is represented i...
Definition: Ray.h:166
bool isApproxZero(const Type &x)
Return true if x is equal to zero to within the default floating-point comparison tolerance...
Definition: Math.h:295
const Vec3T & eye() const
Definition: Ray.h:93
RealT t1() const
Definition: Ray.h:101
bool test(RealT time) const
Return true if time is within t0 and t1, both inclusive.
Definition: Ray.h:119
RayT::Vec3Type Vec3Type
Definition: Ray.h:326
General-purpose arithmetic and comparison routines, most of which accept arbitrary value types (or at...
Vec3Type Vec3T
Definition: Ray.h:59
Ray(const Vec3Type &eye=Vec3Type(0, 0, 0), const Vec3Type &direction=Vec3Type(1, 0, 0), RealT t0=math::Delta< RealT >::value(), RealT t1=std::numeric_limits< RealT >::max())
Definition: Ray.h:61
bool intersects(const Vec3T &center, RealT radius) const
Return true if this ray intersects the specified sphere.
Definition: Ray.h:202
void scaleTimes(RealT scale)
Definition: Ray.h:83
Delta for small floating-point offsets.
Definition: Math.h:123
DDA(const RayT &ray, RealT startTime)
Definition: Ray.h:331
bool clip(const BBoxT &bbox)
Return true if this ray intersects the specified bounding box.
Definition: Ray.h:259
MatType scale(const Vec3< typename MatType::value_type > &scaling)
Definition: Mat.h:595
void setMaxTime(RealT t1)
Definition: Ray.h:79
float Sqrt(float x)
Return the square root of a floating-point value.
Definition: Math.h:655
Vec3R start() const
Return the starting point of the ray.
Definition: Ray.h:107
void setEye(const Vec3Type &eye)
Definition: Ray.h:69
RayT::RealType RealType
Definition: Ray.h:324
bool intersects(const Vec3T &center, RealT radius, RealT &t0, RealT &t1) const
Return true if this ray intersects the specified sphere.
Definition: Ray.h:178
Definition: Ray.h:53
RealT RealType
Definition: Ray.h:57
Vec3< Real > Vec3Type
Definition: Ray.h:58
RealType next() const
Return the time (parameterized along the Ray) of the second (i.e. next) hit of a tree node of size 2^...
Definition: Ray.h:388
#define OPENVDB_VERSION_NAME
Definition: version.h:45
Ray indexToWorld(const GridType &grid) const
Return a new ray in world space, assuming the existing ray is represented in the index space of the s...
Definition: Ray.h:158
Signed (x, y, z) 32-bit integer coordinates.
Definition: Coord.h:47
const Vec3T & invDir() const
Definition: Ray.h:97
bool test() const
Return true if t0 is strictly less then t1.
Definition: Ray.h:116
Vec3R operator()(RealT time) const
Return the position along the ray at the specified time.
Definition: Ray.h:104
DDA(const RayT &ray, RealT startTime, RealT maxTime)
Definition: Ray.h:333
DDA(const RayT &ray)
Definition: Ray.h:329
const Coord & voxel() const
Return the index coordinates of the next node or voxel intersected by the ray. If Log2Dim = 0 the ret...
Definition: Ray.h:376
void reset(const Vec3Type &eye, const Vec3Type &direction, RealT t0=0, RealT t1=std::numeric_limits< RealT >::max())
Definition: Ray.h:85
bool step()
Increment the voxel index to next intersected voxel or node and returns true if the step in time does...
Definition: Ray.h:362
T length() const
Length of the vector.
Definition: Vec3.h:208
void setDir(const Vec3Type &dir)
Definition: Ray.h:71
T lengthSqr() const
Definition: Vec3.h:219
size_t MinIndex(const Vec3T &v)
Return the index [0,1,2] of the smallest value in a 3D vector.
Definition: Math.h:803
OPENVDB_API Hermite max(const Hermite &, const Hermite &)
min and max operations done directly on the compressed data.
bool intersects(const BBoxT &bbox, RealT &t0, RealT &t1) const
Return true if the Ray intersects the specified axisaligned bounding box.
Definition: Ray.h:231
A Digital Differential Analyzer specialized for OpenVDB grids.
Definition: Ray.h:321
Ray applyMap(const MapType &map) const
Return a new Ray that is transformed with the specified map.
Definition: Ray.h:128
bool intersects(const Vec3T &normal, RealT distance, RealT &t) const
Return true if the Ray intersects the plane specified by a normal and distance from the origin...
Definition: Ray.h:275
RealType time() const
Return the time (parameterized along the Ray) of the first hit of a tree node of size 2^Log2Dim...
Definition: Ray.h:383
void init(const RayT &ray, RealT startTime, RealT maxTime)
Definition: Ray.h:335
Vec3R mid() const
Return the midpoint of the ray.
Definition: Ray.h:113
const Type & Min(const Type &a, const Type &b)
Return the minimum of two values.
Definition: Math.h:575
bool isZero(const Type &x)
Return true if x is exactly equal to zero.
Definition: Math.h:283
bool isApproxEqual(const Hermite &lhs, const Hermite &rhs)
Definition: Hermite.h:470
RealT t0() const
Definition: Ray.h:99
Ray applyInverseMap(const MapType &map) const
Return a new Ray that is transformed with the inverse of the specified map.
Definition: Ray.h:145
void setTimes(RealT t0, RealT t1)
Definition: Ray.h:81
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:67
RealType RealT
Definition: Ray.h:325
void print(std::ostream &os=std::cout) const
Print information about this DDA for debugging.
Definition: Ray.h:392
T dot(const Vec3< T > &v) const
Dot product.
Definition: Vec3.h:199
Vec3Type Vec3T
Definition: Ray.h:327
bool clip(const Vec3T &center, RealT radius)
Return true if this ray intersects the specified sphere.
Definition: Ray.h:212
void setMinTime(RealT t0)
Definition: Ray.h:77
bool intersects(const Vec3T &normal, const Vec3T &point, RealT &t) const
Return true if the Ray intersects the plane specified by a normal and point.
Definition: Ray.h:288