OpenVDB  2.1.0
RayIntersector.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 //
57 
58 
59 #ifndef OPENVDB_TOOLS_RAYINTERSECTOR_HAS_BEEN_INCLUDED
60 #define OPENVDB_TOOLS_RAYINTERSECTOR_HAS_BEEN_INCLUDED
61 
62 #include <openvdb/math/Ray.h>
63 #include <openvdb/math/Stencils.h>
64 #include <openvdb/Grid.h>
65 #include <openvdb/Types.h>
66 #include <openvdb/math/Math.h>
67 #include <boost/utility.hpp>
68 #include <boost/type_traits/is_floating_point.hpp>
69 
70 namespace openvdb {
72 namespace OPENVDB_VERSION_NAME {
73 namespace tools {
74 
75 // Helper class that implements hierarchical Digital Differential Analyzers
76 // specialized for ray intersections with level sets
77 template <typename GridT, int NodeLevel> struct LevelSetHDDA;
78 
79 
82 template <typename GridT, int NodeLevel> struct VolumeHDDA;
83 
84 // Helper class that implements the actual search of the zero-crossing
85 // of the level set along the direction of a ray. This particular
86 // implementation uses iterative linear search.
87 template<typename GridT, int Iterations = 0, typename RealT = double>
89 
90 
92 
93 
111 template<typename GridT,
112  typename SearchImplT = LinearSearchImpl<GridT>,
113  int NodeLevel = GridT::TreeType::RootNodeType::ChildNodeType::LEVEL,
114  typename RayT = math::Ray<Real> >
116 {
117 public:
118  typedef GridT GridType;
119  typedef RayT RayType;
120  typedef typename RayT::RealType RealType;
121  typedef typename RayT::Vec3T Vec3Type;
122  typedef typename GridT::ValueType ValueT;
123  typedef typename GridT::TreeType TreeT;
124 
125  BOOST_STATIC_ASSERT( NodeLevel >= -1 && NodeLevel < int(TreeT::DEPTH)-1);
126  BOOST_STATIC_ASSERT(boost::is_floating_point<ValueT>::value);
127 
130  LevelSetRayIntersector(const GridT& grid): mTester(grid)
131  {
132  if (!grid.hasUniformVoxels() ) {
134  "LevelSetRayIntersector only supports uniform voxels!");
135  }
136  if (grid.getGridClass() != GRID_LEVEL_SET) {
138  "LevelSetRayIntersector only supports level sets!"
139  "\nUse Grid::setGridClass(openvdb::GRID_LEVEL_SET)");
140  }
141  }
142 
145  bool intersectsIS(const RayType& iRay) const
146  {
147  if (!mTester.setIndexRay(iRay)) return false;//missed bbox
148  return LevelSetHDDA<TreeT, NodeLevel>::test(mTester);
149  }
150 
155  bool intersectsIS(const RayType& iRay, Vec3Type& xyz) const
156  {
157  if (!mTester.setIndexRay(iRay)) return false;//missed bbox
158  if (!LevelSetHDDA<TreeT, NodeLevel>::test(mTester)) return false;//missed level set
159  mTester.getIndexPos(xyz);
160  return true;
161  }
162 
165  bool intersectsWS(const RayType& wRay) const
166  {
167  if (!mTester.setWorldRay(wRay)) return false;//missed bbox
168  return LevelSetHDDA<TreeT, NodeLevel>::test(mTester);
169  }
170 
175  bool intersectsWS(const RayType& wRay, Vec3Type& world) const
176  {
177  if (!mTester.setWorldRay(wRay)) return false;//missed bbox
178  if (!LevelSetHDDA<TreeT, NodeLevel>::test(mTester)) return false;//missed level set
179  mTester.getWorldPos(world);
180  return true;
181  }
182 
189  bool intersectsWS(const RayType& wRay, Vec3Type& world, Vec3Type& normal) const
190  {
191  if (!mTester.setWorldRay(wRay)) return false;//missed bbox
192  if (!LevelSetHDDA<TreeT, NodeLevel>::test(mTester)) return false;//missed level set
193  mTester.getWorldPosAndNml(world, normal);
194  return true;
195  }
196 
197 private:
198 
199  mutable SearchImplT mTester;
200 
201 };// LevelSetRayIntersector
202 
203 
205 
206 
233 template<typename GridT,
234  int NodeLevel = GridT::TreeType::RootNodeType::ChildNodeType::LEVEL,
235  typename RayT = math::Ray<Real> >
237 {
238 public:
239  typedef GridT GridType;
240  typedef RayT RayType;
241  typedef typename RayT::RealType RealType;
242  typedef typename GridT::TreeType TreeT;
243 
244  BOOST_STATIC_ASSERT( NodeLevel >= 0 && NodeLevel < int(TreeT::DEPTH)-1);
245 
249  VolumeRayIntersector(const GridT& grid): mGrid(&grid), mAccessor(grid.tree())
250  {
251  if (!grid.hasUniformVoxels() ) {
253  "VolumeRayIntersector only supports uniform voxels!");
254  }
255  if ( grid.empty() ) {
256  OPENVDB_THROW(RuntimeError, "LinearSearchImpl does not supports empty grids");
257  }
258  grid.tree().root().evalActiveBoundingBox(mBBox, /*visit individual voxels*/false);
259 
260  mBBox.max().offset(1);//padding so the bbox of a node becomes (origin,origin + node_dim)
261  }
262 
268  inline bool setIndexRay(const RayT& iRay)
269  {
270  mRay = iRay;
271  const bool hit = mRay.clip(mBBox);
272  if (hit) mTmax = mRay.t1();
273  return hit;
274  }
275 
287  inline bool setWorldRay(const RayT& wRay)
288  {
289  return this->setIndexRay(wRay.worldToIndex(*mGrid));
290  }
291 
301  inline int march(Real& t0, Real& t1)
302  {
303  const int n = mRay.test() ? VolumeHDDA<TreeT, NodeLevel>::test(*this) : 0;
304  if (n>0) {
305  t0 = mRay.t0();
306  t1 = mRay.t1();
307  }
308  mRay.setTimes(mRay.t1() + math::Delta<RealType>::value(), mTmax);
309  return n;
310  }
311 
314  inline Vec3R getIndexPos(Real time) const { return mRay(time); }
315 
318  inline Vec3R getWorldPos(Real time) const { return mGrid->indexToWorld(mRay(time)); }
319 
320 private:
321 
322  inline void setRange(Real t0, Real t1) { mRay.setTimes(t0, t1); }
323 
325  inline const RayT& ray() const { return mRay; }
326 
328  template <typename NodeT>
329  inline bool hasNode(const Coord& ijk)
330  {
331  return mAccessor.template probeConstNode<NodeT>(ijk) != NULL;
332  }
333 
334  bool isValueOn(const Coord& ijk) { return mAccessor.isValueOn(ijk); }
335 
336  template<typename, int> friend struct VolumeHDDA;
337 
338  typedef typename GridT::ConstAccessor AccessorT;
339 
340  const GridT* mGrid;
341  AccessorT mAccessor;
342  RayT mRay;
343  Real mTmax;
344  math::CoordBBox mBBox;
345 
346 };// VolumeRayIntersector
347 
348 
350 
351 
371 template<typename GridT, int Iterations, typename RealT>
372 class LinearSearchImpl
373 {
374 public:
376  typedef typename GridT::ValueType ValueT;
377  typedef typename GridT::ConstAccessor AccessorT;
379  typedef typename StencilT::Vec3Type Vec3T;
380 
382  LinearSearchImpl(const GridT& grid)
383  : mStencil(grid), mThreshold(2*grid.voxelSize()[0])
384  {
385  if ( grid.empty() ) {
386  OPENVDB_THROW(RuntimeError, "LinearSearchImpl does not supports empty grids");
387  }
388  grid.tree().root().evalActiveBoundingBox(mBBox, /*visit individual voxels*/false);
389  }
390 
394  inline bool setIndexRay(const RayT& iRay)
395  {
396  mRay = iRay;
397  return mRay.clip(mBBox);//did it hit the bbox
398  }
399 
403  inline bool setWorldRay(const RayT& wRay)
404  {
405  mRay = wRay.worldToIndex(mStencil.grid());
406  return mRay.clip(mBBox);//did it hit the bbox
407  }
408 
411  inline void getIndexPos(Vec3d& xyz) const { xyz = mRay(mTime); }
412 
415  inline void getWorldPos(Vec3d& xyz) const { xyz = mStencil.grid().indexToWorld(mRay(mTime)); }
416 
420  inline void getWorldPosAndNml(Vec3d& xyz, Vec3d& nml)
421  {
422  this->getIndexPos(xyz);
423  mStencil.moveTo(xyz);
424  nml = mStencil.gradient(xyz);
425  nml.normalize();
426  xyz = mStencil.grid().indexToWorld(xyz);
427  }
428 
429 private:
430 
433  inline void init(RealT t0)
434  {
435  mT[0] = t0;
436  mV[0] = this->interpValue(t0);
437  }
438 
439  inline void setRange(RealT t0, RealT t1) { mRay.setTimes(t0, t1); }
440 
442  inline const RayT& ray() const { return mRay; }
443 
445  template <typename NodeT>
446  inline bool hasNode(const Coord& ijk)
447  {
448  return mStencil.accessor().template probeConstNode<NodeT>(ijk) != NULL;
449  }
450 
456  inline bool operator()(const Coord& ijk, RealT time)
457  {
458  ValueT V;
459  if (mStencil.accessor().probeValue(ijk, V) &&//inside narrow band?
460  math::Abs(V)<mThreshold) {// close to zero-crossing?
461  mT[1] = time;
462  mV[1] = this->interpValue(time);
463  if (math::ZeroCrossing(mV[0], mV[1])) {
464  mTime = this->interpTime();
466  for (int n=0; Iterations>0 && n<Iterations; ++n) {//resolved at compile-time
467  V = this->interpValue(mTime);
468  const int m = math::ZeroCrossing(mV[0], V) ? 1 : 0;
469  mV[m] = V;
470  mT[m] = mTime;
471  mTime = this->interpTime();
472  }
474  return true;
475  }
476  mT[0] = mT[1];
477  mV[0] = mV[1];
478  }
479  return false;
480  }
481 
482  inline RealT interpTime()
483  {
484  assert(math::isApproxLarger(mT[1], mT[0], 1e-6));
485  return mT[0]+(mT[1]-mT[0])*mV[0]/(mV[0]-mV[1]);
486  }
487 
488  inline RealT interpValue(RealT time)
489  {
490  const Vec3R pos = mRay(time);
491  mStencil.moveTo(pos);
492  return mStencil.interpolation(pos);
493  }
494 
495  template<typename, int> friend struct LevelSetHDDA;
496 
497  RayT mRay;
498  StencilT mStencil;
499  RealT mTime;
500  ValueT mV[2];
501  RealT mT[2];
502  ValueT mThreshold;
503  math::CoordBBox mBBox;
504 };// LinearSearchImpl
505 
507 
508 
511 template<typename TreeT, int NodeLevel>
512 struct LevelSetHDDA
513 {
514  typedef typename TreeT::RootNodeType::NodeChainType ChainT;
515  typedef typename boost::mpl::at<ChainT, boost::mpl::int_<NodeLevel> >::type NodeT;
516 
517  template <typename TesterT>
518  static bool test(TesterT& tester)
519  {
521  do {
522  if (tester.template hasNode<NodeT>(dda.voxel())) {
523  tester.setRange(dda.time(), dda.next());
524  if (LevelSetHDDA<TreeT, NodeLevel-1>::test(tester)) return true;
525  }
526  } while(dda.step());
527  return false;
528  }
529 };
530 
533 template<typename TreeT>
534 struct LevelSetHDDA<TreeT, -1>
535 {
536  template <typename TesterT>
537  static bool test(TesterT& tester)
538  {
539  math::DDA<typename TesterT::RayT, 0> dda(tester.ray());
540  tester.init(dda.time());
541  do { if (tester(dda.voxel(), dda.next())) return true; } while(dda.step());
542  return false;
543  }
544 };
545 
546 
548 
549 
552 template <typename TreeT, int NodeLevel>
553 struct VolumeHDDA
554 {
555  typedef typename TreeT::RootNodeType::NodeChainType ChainT;
556  typedef typename boost::mpl::at<ChainT, boost::mpl::int_<NodeLevel> >::type NodeT;
557 
558  template <typename TesterT>
559  static int test(TesterT& tester)
560  {
562  do {
563  if (tester.template hasNode<NodeT>(dda.voxel())) {//child node
564  tester.setRange(dda.time(), dda.next());
565  int hit = VolumeHDDA<TreeT, NodeLevel-1>::test(tester);
566  if (hit > 0) return hit;
567  } else if (tester.isValueOn(dda.voxel())) {//active tile
568  tester.setRange(dda.time(), dda.next());
569  return 1;//active tile
570  }
571  } while (dda.step());
572  return 0;//no hits
573  }
574 };
575 
578 template <typename TreeT>
579 struct VolumeHDDA<TreeT, -1>
580 {
581  template <typename TesterT>
582  static int test(TesterT&) { return 2; }//hit leaf so don't traverse voxels.
583 };
584 
585 } // namespace tools
586 } // namespace OPENVDB_VERSION_NAME
587 } // namespace openvdb
588 
589 #endif // OPENVDB_TOOLS_RAYINTERSECTOR_HAS_BEEN_INCLUDED
590 
591 // Copyright (c) 2012-2013 DreamWorks Animation LLC
592 // All rights reserved. This software is distributed under the
593 // Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ )
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
Implements linear iterative search for a zero-crossing of the level set along along the direction of ...
Definition: RayIntersector.h:88
static bool test(TesterT &tester)
Definition: RayIntersector.h:537
RayT RayType
Definition: RayIntersector.h:240
Definition: RayIntersector.h:82
TreeT::RootNodeType::NodeChainType ChainT
Definition: RayIntersector.h:555
LinearSearchImpl(const GridT &grid)
Constructor from a grid.
Definition: RayIntersector.h:382
math::Vec3< Real > Vec3R
Definition: Types.h:74
LevelSetRayIntersector(const GridT &grid)
Constructor.
Definition: RayIntersector.h:130
Definition: Types.h:137
#define OPENVDB_THROW(exception, message)
Definition: Exceptions.h:97
bool ZeroCrossing(const Type &a, const Type &b)
Return true if the interval [a, b] includes zero, i.e., if either a or b is zero or if they have diff...
Definition: Math.h:647
GridT::TreeType TreeT
Definition: RayIntersector.h:242
GridT GridType
Definition: RayIntersector.h:118
void getWorldPos(Vec3d &xyz) const
Get the intersection point in world space.
Definition: RayIntersector.h:415
int32_t Abs(int32_t i)
Return the absolute value of the given quantity.
Definition: Math.h:253
Vec3R getWorldPos(Real time) const
Return the floating-point world position along the current index ray at the specified time...
Definition: RayIntersector.h:318
Delta for small floating-point offsets.
Definition: Math.h:123
VolumeRayIntersector(const GridT &grid)
Constructor.
Definition: RayIntersector.h:249
bool setIndexRay(const RayT &iRay)
Return false the ray misses the bbox of the grid.
Definition: RayIntersector.h:394
int march(Real &t0, Real &t1)
Return 0 if not hit was detected. A return value of 1 means it hit an active tile, and a return value of 2 means it hit a LeafNode. Only when a hit is detected are t0 and t1 updated with the corresponding entry and exit times along the INDEX ray!
Definition: RayIntersector.h:301
TreeT::RootNodeType::NodeChainType ChainT
Definition: RayIntersector.h:514
RayT::RealType RealType
Definition: RayIntersector.h:241
bool intersectsIS(const RayType &iRay) const
Return true if the index-space ray intersects the level set.
Definition: RayIntersector.h:145
bool intersectsIS(const RayType &iRay, Vec3Type &xyz) const
Return true if the index-space ray intersects the level set.
Definition: RayIntersector.h:155
bool intersectsWS(const RayType &wRay) const
Return true if the world-space ray intersects the level set.
Definition: RayIntersector.h:165
Definition: Exceptions.h:86
Definition: Ray.h:53
GridT GridType
Definition: RayIntersector.h:239
void getIndexPos(Vec3d &xyz) const
Get the intersection point in index space.
Definition: RayIntersector.h:411
Axis-aligned bounding box of signed integer coordinates.
Definition: Coord.h:229
#define OPENVDB_VERSION_NAME
Definition: version.h:45
This class provides the public API for intersecting a ray with a generic (e.g. density) volume...
Definition: RayIntersector.h:236
static int test(TesterT &tester)
Definition: RayIntersector.h:559
StencilT::Vec3Type Vec3T
Definition: RayIntersector.h:379
#define OPENVDB_NO_UNREACHABLE_CODE_WARNING_BEGIN
Definition: Platform.h:96
GridT::TreeType TreeT
Definition: RayIntersector.h:123
boost::mpl::at< ChainT, boost::mpl::int_< NodeLevel > >::type NodeT
Definition: RayIntersector.h:556
RayT RayType
Definition: RayIntersector.h:119
math::BoxStencil< GridT > StencilT
Definition: RayIntersector.h:378
bool setIndexRay(const RayT &iRay)
Return false if the index ray misses the bbox of the grid.
Definition: RayIntersector.h:268
math::Ray< RealT > RayT
Definition: RayIntersector.h:375
Vec3< double > Vec3d
Definition: Vec3.h:605
A Digital Differential Analyzer specialized for OpenVDB grids.
Definition: Ray.h:321
RayT::Vec3T Vec3Type
Definition: RayIntersector.h:121
bool setWorldRay(const RayT &wRay)
Return false the ray misses the bbox of the grid.
Definition: RayIntersector.h:403
double Real
Definition: Types.h:62
void init(const RayT &ray, RealT startTime, RealT maxTime)
Definition: Ray.h:335
#define OPENVDB_NO_UNREACHABLE_CODE_WARNING_END
Definition: Platform.h:97
static bool test(TesterT &tester)
Definition: RayIntersector.h:518
Helper class that implements Hierarchical Digital Differential Analyzers and is specialized for ray i...
Definition: RayIntersector.h:77
GridT::ValueType ValueT
Definition: RayIntersector.h:376
void getWorldPosAndNml(Vec3d &xyz, Vec3d &nml)
Get the intersection point and normal in world space.
Definition: RayIntersector.h:420
GridT::ConstAccessor AccessorT
Definition: RayIntersector.h:377
bool setWorldRay(const RayT &wRay)
Return false if the world ray misses the bbox of the grid.
Definition: RayIntersector.h:287
static int test(TesterT &)
Definition: RayIntersector.h:582
This class provides the public API for intersecting a ray with a narrow-band level set...
Definition: RayIntersector.h:115
#define OPENVDB_USE_VERSION_NAMESPACE
Definition: version.h:67
GridT::ValueType ValueT
Definition: RayIntersector.h:122
bool intersectsWS(const RayType &wRay, Vec3Type &world) const
Return true if the world-space ray intersects the level set.
Definition: RayIntersector.h:175
Vec3R getIndexPos(Real time) const
Return the floating-point index position along the current index ray at the specified time...
Definition: RayIntersector.h:314
bool intersectsWS(const RayType &wRay, Vec3Type &world, Vec3Type &normal) const
Return true if the world-space ray intersects the level set.
Definition: RayIntersector.h:189
bool isApproxLarger(const Type &a, const Type &b, const Type &tolerance)
Return true if a is larger than b to within the given tolerance, i.e., if b - a &lt; tolerance...
Definition: Math.h:351
bool clip(const Vec3T &center, RealT radius)
Return true if this ray intersects the specified sphere.
Definition: Ray.h:212
boost::mpl::at< ChainT, boost::mpl::int_< NodeLevel > >::type NodeT
Definition: RayIntersector.h:515
RayT::RealType RealType
Definition: RayIntersector.h:120