| 1 | /*
|
|---|
| 2 | * BaseShapes_impl.cpp
|
|---|
| 3 | *
|
|---|
| 4 | * Created on: Jun 18, 2010
|
|---|
| 5 | * Author: crueger
|
|---|
| 6 | */
|
|---|
| 7 |
|
|---|
| 8 | #include "Shapes/BaseShapes.hpp"
|
|---|
| 9 | #include "Shapes/BaseShapes_impl.hpp"
|
|---|
| 10 | #include "Shapes/ShapeOps.hpp"
|
|---|
| 11 |
|
|---|
| 12 | #include "vector.hpp"
|
|---|
| 13 | #include "Helpers/Assert.hpp"
|
|---|
| 14 |
|
|---|
| 15 | #include "Line.hpp"
|
|---|
| 16 | #include "Plane.hpp"
|
|---|
| 17 | #include "LineSegment.hpp"
|
|---|
| 18 | #include "LineSegmentSet.hpp"
|
|---|
| 19 |
|
|---|
| 20 | #include <cmath>
|
|---|
| 21 | #include <algorithm>
|
|---|
| 22 |
|
|---|
| 23 | bool Sphere_impl::isInside(const Vector &point){
|
|---|
| 24 | return point.NormSquared()<=1;
|
|---|
| 25 | }
|
|---|
| 26 |
|
|---|
| 27 | bool Sphere_impl::isOnSurface(const Vector &point){
|
|---|
| 28 | return fabs(point.NormSquared()-1)<MYEPSILON;
|
|---|
| 29 | }
|
|---|
| 30 |
|
|---|
| 31 | Vector Sphere_impl::getNormal(const Vector &point) throw(NotOnSurfaceException){
|
|---|
| 32 | if(!isOnSurface(point)){
|
|---|
| 33 | throw NotOnSurfaceException(__FILE__,__LINE__);
|
|---|
| 34 | }
|
|---|
| 35 | return point;
|
|---|
| 36 | }
|
|---|
| 37 |
|
|---|
| 38 | LineSegmentSet Sphere_impl::getLineIntersections(const Line &line){
|
|---|
| 39 | LineSegmentSet res(line);
|
|---|
| 40 | std::vector<Vector> intersections = line.getSphereIntersections();
|
|---|
| 41 | if(intersections.size()==2){
|
|---|
| 42 | res.insert(LineSegment(intersections[0],intersections[1]));
|
|---|
| 43 | }
|
|---|
| 44 | return res;
|
|---|
| 45 | }
|
|---|
| 46 |
|
|---|
| 47 | string Sphere_impl::toString(){
|
|---|
| 48 | return "Sphere()";
|
|---|
| 49 | }
|
|---|
| 50 |
|
|---|
| 51 | Shape Sphere(){
|
|---|
| 52 | Shape::impl_ptr impl = Shape::impl_ptr(new Sphere_impl());
|
|---|
| 53 | return Shape(impl);
|
|---|
| 54 | }
|
|---|
| 55 |
|
|---|
| 56 | Shape Sphere(const Vector ¢er,double radius){
|
|---|
| 57 | return translate(resize(Sphere(),radius),center);
|
|---|
| 58 | }
|
|---|
| 59 |
|
|---|
| 60 | Shape Ellipsoid(const Vector ¢er, const Vector &radius){
|
|---|
| 61 | return translate(stretch(Sphere(),radius),center);
|
|---|
| 62 | }
|
|---|
| 63 |
|
|---|
| 64 | bool Cuboid_impl::isInside(const Vector &point){
|
|---|
| 65 | return fabs(point[0])<=1 && fabs(point[1])<=1 && fabs(point[2])<=1;
|
|---|
| 66 | }
|
|---|
| 67 |
|
|---|
| 68 | bool Cuboid_impl::isOnSurface(const Vector &point){
|
|---|
| 69 | bool retVal = isInside(point);
|
|---|
| 70 | // test all borders of the cuboid
|
|---|
| 71 | // double fabs
|
|---|
| 72 | retVal = retVal &&
|
|---|
| 73 | ((fabs(fabs(point[0])-1) < MYEPSILON) ||
|
|---|
| 74 | (fabs(fabs(point[1])-1) < MYEPSILON) ||
|
|---|
| 75 | (fabs(fabs(point[2])-1) < MYEPSILON));
|
|---|
| 76 | return retVal;
|
|---|
| 77 | }
|
|---|
| 78 |
|
|---|
| 79 | Vector Cuboid_impl::getNormal(const Vector &point) throw(NotOnSurfaceException){
|
|---|
| 80 | if(!isOnSurface(point)){
|
|---|
| 81 | throw NotOnSurfaceException(__FILE__,__LINE__);
|
|---|
| 82 | }
|
|---|
| 83 | Vector res;
|
|---|
| 84 | // figure out on which sides the Vector lies (maximum 3, when it is in a corner)
|
|---|
| 85 | for(int i=NDIM;i--;){
|
|---|
| 86 | if(fabs(fabs(point[i])-1)<MYEPSILON){
|
|---|
| 87 | // add the scaled (-1/+1) Vector to the set of surface vectors
|
|---|
| 88 | res[i] = point[i];
|
|---|
| 89 | }
|
|---|
| 90 | }
|
|---|
| 91 | ASSERT(res.NormSquared()>=1 && res.NormSquared()<=3,"To many or to few sides found for this Vector");
|
|---|
| 92 |
|
|---|
| 93 | res.Normalize();
|
|---|
| 94 | return res;
|
|---|
| 95 | }
|
|---|
| 96 |
|
|---|
| 97 | LineSegmentSet Cuboid_impl::getLineIntersections(const Line &line){
|
|---|
| 98 | LineSegmentSet res(line);
|
|---|
| 99 | // get the intersection on each of the six faces
|
|---|
| 100 | vector<Vector> intersections;
|
|---|
| 101 | intersections.resize(2);
|
|---|
| 102 | int c=0;
|
|---|
| 103 | int x[2]={-1,+1};
|
|---|
| 104 | for(int i=NDIM;i--;){
|
|---|
| 105 | for(int p=0;p<2;++p){
|
|---|
| 106 | if(c==2) goto end; // I know this sucks, but breaking two loops is stupid
|
|---|
| 107 | Vector base;
|
|---|
| 108 | base[i]=x[p];
|
|---|
| 109 | // base now points to the surface and is normal to it at the same time
|
|---|
| 110 | Plane p(base,base);
|
|---|
| 111 | Vector intersection = p.GetIntersection(line);
|
|---|
| 112 | if(isInside(intersection)){
|
|---|
| 113 | // if we have a point on the edge it might already be contained
|
|---|
| 114 | if(c==1 && intersections[0]==intersection)
|
|---|
| 115 | continue;
|
|---|
| 116 | intersections[c++]=intersection;
|
|---|
| 117 | }
|
|---|
| 118 | }
|
|---|
| 119 | }
|
|---|
| 120 | end:
|
|---|
| 121 | if(c==2){
|
|---|
| 122 | res.insert(LineSegment(intersections[0],intersections[1]));
|
|---|
| 123 | }
|
|---|
| 124 | return res;
|
|---|
| 125 | }
|
|---|
| 126 |
|
|---|
| 127 | string Cuboid_impl::toString(){
|
|---|
| 128 | return "Cuboid()";
|
|---|
| 129 | }
|
|---|
| 130 |
|
|---|
| 131 | Shape Cuboid(){
|
|---|
| 132 | Shape::impl_ptr impl = Shape::impl_ptr(new Cuboid_impl());
|
|---|
| 133 | return Shape(impl);
|
|---|
| 134 | }
|
|---|
| 135 |
|
|---|
| 136 | Shape Cuboid(const Vector &corner1, const Vector &corner2){
|
|---|
| 137 | // make sure the two edges are upper left front and lower right back
|
|---|
| 138 | Vector sortedC1;
|
|---|
| 139 | Vector sortedC2;
|
|---|
| 140 | for(int i=NDIM;i--;){
|
|---|
| 141 | sortedC1[i] = min(corner1[i],corner2[i]);
|
|---|
| 142 | sortedC2[i] = max(corner1[i],corner2[i]);
|
|---|
| 143 | ASSERT(corner1[i]!=corner2[i],"Given points for cuboid edges did not define a valid space");
|
|---|
| 144 | }
|
|---|
| 145 | // get the middle point
|
|---|
| 146 | Vector middle = (1./2.)*(sortedC1+sortedC2);
|
|---|
| 147 | Vector factors = sortedC2-middle;
|
|---|
| 148 | return translate(stretch(Cuboid(),factors),middle);
|
|---|
| 149 | }
|
|---|