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 BoundingBox<T, D> get_axis_aligned_bounding_box(unsigned = 0) const;
  35         virtual bool contains(const LinAl::Vector<T, D> &) const;
  36         virtual unsigned get_max_ray_intersections() const;
  37         virtual unsigned get_intersections(const Ray<T, D> &, SurfacePoint<T, D> *, unsigned) const;
  38         virtual Coverage get_coverage(const BoundingBox<T, D> &) 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<=T(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 BoundingBox<T, D> ExtrudedShape<T, D>::get_axis_aligned_bounding_box(unsigned detail) const
  79 {
  80         BoundingBox<T, D-1> base_bbox = base->get_axis_aligned_bounding_box(detail);
  81         T half_length = length/T(2);
  82         return BoundingBox<T, D>(compose(base_bbox.get_minimum_point(), -half_length),
  83                 compose(base_bbox.get_maximum_point(), half_length));
  84 }
  85
  86 template<typename T, unsigned D>
  87 inline bool ExtrudedShape<T, D>::contains(const LinAl::Vector<T, D> &point) const
  88 {
  89         using std::abs;
  90
  91         if(abs(point[D-1])>length/T(2))
  92                 return false;
  93
  94         return base->contains(point.template slice<D-1>(0));
  95 }
  96
  97 template<typename T, unsigned D>
  98 inline unsigned ExtrudedShape<T, D>::get_max_ray_intersections() const
  99 {
 100         return std::max(base->get_max_ray_intersections(), 2U);
 101 }
 102
 103 template<typename T, unsigned D>
 104 inline unsigned ExtrudedShape<T, D>::get_intersections(const Ray<T, D> &ray, SurfacePoint<T, D> *points, unsigned size) const
 105 {
 106         using std::swap;
 107
 108         unsigned n = 0;
 109         T half_length = length/T(2);
 110         const LinAl::Vector<T, D> &ray_start = ray.get_start();
 111         const LinAl::Vector<T, D> &ray_direction = ray.get_direction();
 112         LinAl::Vector<T, D-1> base_dir = ray_direction.template slice<D-1>(0);
 113
 114         /* If the ray does not degenerate to a point in the base space, it could
 115         intersect the base shape. */
 116         if(inner_product(base_dir, base_dir)!=T(0))
 117         {
 118                 T offset = T();
 119                 T limit = T();
 120                 if(ray_direction[D-1]!=T(0))
 121                 {
 122                         offset = (half_length-ray_start[D-1])/ray_direction[D-1];
 123                         limit = (-half_length-ray_start[D-1])/ray_direction[D-1];
 124                         if(offset>limit)
 125                                 swap(offset, limit);
 126                         if(offset<T(0))
 127                                 offset = T(0);
 128                 }
 129
 130                 if(limit>=offset)
 131                 {
 132                         T distortion = base_dir.norm();
 133                         Ray<T, D-1> base_ray((ray_start+ray_direction*offset).template slice<D-1>(0),
 134                                 base_dir, (limit-offset)*distortion);
 135
 136                         SurfacePoint<T, D-1> *base_points = 0;
 137                         if(points)
 138                                 /* Shamelessly reuse the provided storage.  Align to the end of the
 139                                 array so processing can start from the first (nearest) point. */
 140                                 base_points = reinterpret_cast<SurfacePoint<T, D-1> *>(points+size)-size;
 141
 142                         unsigned count = base->get_intersections(base_ray, base_points, size);
 143                         for(unsigned i=0; (n<size && i<count); ++i)
 144                         {
 145                                 if(points)
 146                                 {
 147                                         T x = offset+base_points[i].distance/distortion;
 148                                         points[n].position = ray_start+ray_direction*x;
 149                                         points[n].normal = compose(base_points[i].normal, T(0));
 150                                         points[n].distance = x;
 151                                         points[n].entry = base_points[i].entry;
 152                                 }
 153
 154                                 ++n;
 155                         }
 156                 }
 157         }
 158
 159         /* If the ray is not parallel to the base space, it may pass through the
 160         caps. */
 161         if(n<size && ray_direction[D-1]!=T(0))
 162         {
 163                 for(int i=-1; (n<size && i<=1); i+=2)
 164                 {
 165                         T x = (half_length*i-ray_start[D-1])/ray_direction[D-1];
 166                         if(!ray.check_limits(x))
 167                                 continue;
 168
 169                         LinAl::Vector<T, D> p = ray_start+ray_direction*x;
 170                         if(base->contains(p.template slice<D-1>(0)))
 171                         {
 172                                 if(points)
 173                                 {
 174                                         points[n].position = p;
 175                                         points[n].normal = LinAl::Vector<T, D>();
 176                                         points[n].normal[D-1] = i;
 177                                         points[n].distance = x;
 178                                         points[n].entry = (T(i)*ray_direction[D-1]<T(0));
 179                                 }
 180
 181                                 ++n;
 182                         }
 183                 }
 184
 185                 sort_points(points, n);
 186         }
 187
 188         return n;
 189 }
 190
 191 template<typename T, unsigned D>
 192 inline Coverage ExtrudedShape<T, D>::get_coverage(const BoundingBox<T, D> &bbox) const
 193 {
 194         T half_length = length/T(2);
 195         if(bbox.get_maximum_coordinate(D-1)<-half_length || bbox.get_minimum_coordinate(D-1)>half_length)
 196                 return NO_COVERAGE;
 197
 198         BoundingBox<T, D-1> base_bbox(bbox.get_minimum_point().template slice<D-1>(0), bbox.get_maximum_point().template slice<D-1>(0));
 199         Coverage coverage = base->get_coverage(base_bbox);
 200         if(coverage==NO_COVERAGE)
 201                 return NO_COVERAGE;
 202
 203         if(bbox.get_minimum_coordinate(D-1)<-half_length || bbox.get_maximum_coordinate(D-1)>half_length)
 204                 return PARTIAL_COVERAGE;
 205         else
 206                 return coverage;
 207 }
 208
 209 } // namespace Geometry
 210 } // namespace Msp
 211
 212 #endif