source/geometry/extrudedshape.h

   1 #ifndef MSP_GEOMETRY_EXTRUDEDSHAPE_H_
   2 #define MSP_GEOMETRY_EXTRUDEDSHAPE_H_
   3
   4 #include <algorithm>
   5 #include <cmath>
   6 #include "shape.h"
   7
   8 namespace Msp {
   9 namespace Geometry {
  10
  11 /**
  12 A shape embedded in space of dimension higher by one and extruded towards the
  13 highest dimension.  As an example, extruding a circle creates a cylinder.  The
  14 base shape's orientation is not changed.
  15 */
  16 template<typename T, unsigned D>
  17 class ExtrudedShape: public Shape<T, D>
  18 {
  19 private:
  20         Shape<T, D-1> *base;
  21         T length;
  22
  23 public:
  24         ExtrudedShape(const Shape<T, D-1> &, T);
  25         ExtrudedShape(const ExtrudedShape &);
  26         ExtrudedShape &operator=(const ExtrudedShape &);
  27         virtual ~ExtrudedShape();
  28
  29         virtual ExtrudedShape *clone() const;
  30
  31         const Shape<T, D-1> &get_base() const { return *base; }
  32         T get_length() const { return length; }
  33
  34         virtual HyperBox<T, D> get_axis_aligned_bounding_box() const;
  35         virtual bool contains(const LinAl::Vector<T, D> &) const;
  36         virtual bool check_intersection(const Ray<T, D> &) const;
  37         virtual unsigned get_max_ray_intersections() const;
  38         virtual unsigned get_intersections(const Ray<T, D> &, SurfacePoint<T, D> *, unsigned) const;
  39 };
  40
  41 template<typename T, unsigned D>
  42 inline ExtrudedShape<T, D>::ExtrudedShape(const Shape<T, D-1> &b, T l):
  43         length(l)
  44 {
  45         if(l<=0)
  46                 throw std::invalid_argument("ExtrudedShape::ExtrudedShape");
  47
  48         base = b.clone();
  49 }
  50
  51 template<typename T, unsigned D>
  52 inline ExtrudedShape<T, D>::ExtrudedShape(const ExtrudedShape<T, D> &other):
  53         base(other.base.clone()),
  54         length(other.length)
  55 { }
  56
  57 template<typename T, unsigned D>
  58 inline ExtrudedShape<T, D> &ExtrudedShape<T, D>::operator=(const ExtrudedShape<T, D> &other)
  59 {
  60         delete base;
  61         base = other.base.clone();
  62         length = other.length;
  63 }
  64
  65 template<typename T, unsigned D>
  66 inline ExtrudedShape<T, D>::~ExtrudedShape()
  67 {
  68         delete base;
  69 }
  70
  71 template<typename T, unsigned D>
  72 inline ExtrudedShape<T, D> *ExtrudedShape<T, D>::clone() const
  73 {
  74         return new ExtrudedShape<T, D>(*base, length);
  75 }
  76
  77 template<typename T, unsigned D>
  78 inline HyperBox<T, D> ExtrudedShape<T, D>::get_axis_aligned_bounding_box() const
  79 {
  80         HyperBox<T, D-1> base_bbox = base->get_axis_aligned_bounding_box();
  81         return HyperBox<T, D>(LinAl::Vector<T, D>(base_bbox.get_dimensions(), length));
  82 }
  83
  84 template<typename T, unsigned D>
  85 inline bool ExtrudedShape<T, D>::contains(const LinAl::Vector<T, D> &point) const
  86 {
  87         using std::abs;
  88
  89         if(abs(point[D-1])>length/T(2))
  90                 return false;
  91
  92         return base->contains(LinAl::Vector<T, D-1>(point));
  93 }
  94
  95 template<typename T, unsigned D>
  96 inline bool ExtrudedShape<T, D>::check_intersection(const Ray<T, D> &ray) const
  97 {
  98         return get_intersections(ray, 0, 1);
  99 }
 100
 101 template<typename T, unsigned D>
 102 inline unsigned ExtrudedShape<T, D>::get_max_ray_intersections() const
 103 {
 104         return std::max(base->get_max_ray_intersections(), 2U);
 105 }
 106
 107 template<typename T, unsigned D>
 108 inline unsigned ExtrudedShape<T, D>::get_intersections(const Ray<T, D> &ray, SurfacePoint<T, D> *points, unsigned size) const
 109 {
 110         using std::abs;
 111         using std::sqrt;
 112         using std::swap;
 113
 114         unsigned n = 0;
 115         T half_length = length/T(2);
 116         const LinAl::Vector<T, D> &ray_start = ray.get_start();
 117         const LinAl::Vector<T, D> &ray_direction = ray.get_direction();
 118         LinAl::Vector<T, D-1> base_dir(ray_direction);
 119
 120         /* If the ray does not degenerate to a point in the base space, it could
 121         intersect the base shape. */
 122         if(inner_product(base_dir, base_dir)!=T(0))
 123         {
 124                 T offset = T();
 125                 T limit = T();
 126                 if(ray.get_direction()[D-1]!=T(0))
 127                 {
 128                         offset = (half_length-ray_start[D-1])/ray_direction[D-1];
 129                         limit = (-half_length-ray_start[D-1])/ray_direction[D-1];
 130                         if(offset>limit)
 131                                 swap(offset, limit);
 132                         if(offset<T(0))
 133                                 offset = T(0);
 134                 }
 135                 T distortion = base_dir.norm();
 136                 Ray<T, D-1> base_ray(LinAl::Vector<T, D-1>(ray_start+ray_direction*offset),
 137                         base_dir, (limit-offset)*distortion);
 138
 139                 SurfacePoint<T, D-1> *base_points = 0;
 140                 if(points)
 141                         /* Shamelessly reuse the provided storage.  Align to the end of the array
 142                         so processing can start from the first (nearest) point. */
 143                         base_points = reinterpret_cast<SurfacePoint<T, D-1> *>(points+size)-size;
 144
 145                 unsigned count = base->get_intersections(base_ray, base_points, size);
 146                 for(unsigned i=0; i<count; ++i)
 147                 {
 148                         if(points)
 149                         {
 150                                 T x = offset+base_points[i].distance/distortion;
 151                                 points[n].position = ray_start+ray_direction*x;
 152                                 points[n].normal = LinAl::Vector<T, D>(base_points[i].normal, T(0));
 153                                 points[n].distance = x;
 154                         }
 155
 156                         ++n;
 157                         if(n==size)
 158                                 return n;
 159                 }
 160         }
 161
 162         /* If the ray is not parallel to the base space, it may pass through the
 163         caps. */
 164         if(ray_direction[D-1])
 165         {
 166                 for(int i=-1; i<=1; i+=2)
 167                 {
 168                         T x = (half_length*i-ray_start[D-1])/ray_direction[D-1];
 169                         if(!ray.check_limits(x))
 170                                 continue;
 171
 172                         LinAl::Vector<T, D> p = ray_start+ray_direction*x;
 173                         if(base->contains(LinAl::Vector<T, D-1>(p)) && n<size)
 174                         {
 175                                 if(points)
 176                                 {
 177                                         points[n].position = p;
 178                                         points[n].normal = LinAl::Vector<T, D>();
 179                                         points[n].normal[D-1] = i;
 180                                         points[n].distance = x;
 181
 182                                         if(n==1 && x<points[0].distance)
 183                                                 swap(points[0], points[1]);
 184                                 }
 185
 186                                 ++n;
 187                                 if(n==size)
 188                                         return n;
 189                         }
 190                 }
 191         }
 192
 193         return n;
 194 }
 195
 196 } // namespace Geometry
 197 } // namespace Msp
 198
 199 #endif