| [bcf653] | 1 | /*
 | 
|---|
 | 2 |  * Project: MoleCuilder
 | 
|---|
 | 3 |  * Description: creates and alters molecular systems
 | 
|---|
| [0aa122] | 4 |  * Copyright (C)  2010-2012 University of Bonn. All rights reserved.
 | 
|---|
| [bcf653] | 5 |  * Please see the LICENSE file or "Copyright notice" in builder.cpp for details.
 | 
|---|
 | 6 |  */
 | 
|---|
 | 7 | 
 | 
|---|
| [997784] | 8 | /*
 | 
|---|
 | 9 |  * Shape.cpp
 | 
|---|
 | 10 |  *
 | 
|---|
 | 11 |  *  Created on: Jun 18, 2010
 | 
|---|
 | 12 |  *      Author: crueger
 | 
|---|
 | 13 |  */
 | 
|---|
 | 14 | 
 | 
|---|
| [bf3817] | 15 | // include config.h
 | 
|---|
 | 16 | #ifdef HAVE_CONFIG_H
 | 
|---|
 | 17 | #include <config.h>
 | 
|---|
 | 18 | #endif
 | 
|---|
 | 19 | 
 | 
|---|
| [ad011c] | 20 | #include "CodePatterns/MemDebug.hpp"
 | 
|---|
| [997784] | 21 | 
 | 
|---|
| [ad011c] | 22 | #include "CodePatterns/Assert.hpp"
 | 
|---|
| [6c438f] | 23 | #include "LinearAlgebra/Vector.hpp"
 | 
|---|
 | 24 | 
 | 
|---|
| [b94634] | 25 | #include "Shapes/Shape.hpp"
 | 
|---|
 | 26 | #include "Shapes/Shape_impl.hpp"
 | 
|---|
 | 27 | #include "Shapes/ShapeExceptions.hpp"
 | 
|---|
| [5de9da] | 28 | 
 | 
|---|
| [cfda65] | 29 | #include <string>
 | 
|---|
 | 30 | 
 | 
|---|
 | 31 | 
 | 
|---|
| [997784] | 32 | Shape::Shape(const Shape& src) :
 | 
|---|
 | 33 |   impl(src.getImpl())
 | 
|---|
 | 34 | {}
 | 
|---|
 | 35 | 
 | 
|---|
 | 36 | Shape::~Shape(){}
 | 
|---|
 | 37 | 
 | 
|---|
| [205d9b] | 38 | bool Shape::isInside(const Vector &point) const{
 | 
|---|
| [997784] | 39 |   return impl->isInside(point);
 | 
|---|
 | 40 | }
 | 
|---|
 | 41 | 
 | 
|---|
| [5de9da] | 42 | bool Shape::isOnSurface(const Vector &point) const{
 | 
|---|
 | 43 |   return impl->isOnSurface(point);
 | 
|---|
 | 44 | }
 | 
|---|
 | 45 | 
 | 
|---|
 | 46 | Vector Shape::getNormal(const Vector &point) const throw (NotOnSurfaceException){
 | 
|---|
 | 47 |   return impl->getNormal(point);
 | 
|---|
 | 48 | }
 | 
|---|
 | 49 | 
 | 
|---|
| [c6f395] | 50 | LineSegmentSet Shape::getLineIntersections(const Line &line){
 | 
|---|
 | 51 |   return impl->getLineIntersections(line);
 | 
|---|
 | 52 | }
 | 
|---|
 | 53 | 
 | 
|---|
| [9c1c89] | 54 | std::vector<Vector> Shape::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [c5186e] | 55 |   return impl->getHomogeneousPointsOnSurface(N);
 | 
|---|
 | 56 | }
 | 
|---|
 | 57 | 
 | 
|---|
| [997784] | 58 | Shape::Shape(Shape::impl_ptr _impl) :
 | 
|---|
 | 59 |     impl(_impl)
 | 
|---|
 | 60 | {}
 | 
|---|
 | 61 | 
 | 
|---|
 | 62 | Shape &Shape::operator=(const Shape& rhs){
 | 
|---|
 | 63 |   if(&rhs!=this){
 | 
|---|
 | 64 |     impl=rhs.getImpl();
 | 
|---|
 | 65 |   }
 | 
|---|
 | 66 |   return *this;
 | 
|---|
 | 67 | }
 | 
|---|
 | 68 | 
 | 
|---|
| [cfda65] | 69 | std::string Shape::toString() const{
 | 
|---|
 | 70 |   return impl->toString();
 | 
|---|
 | 71 | }
 | 
|---|
 | 72 | 
 | 
|---|
| [997784] | 73 | Shape::impl_ptr Shape::getImpl() const{
 | 
|---|
 | 74 |   return impl;
 | 
|---|
 | 75 | }
 | 
|---|
 | 76 | 
 | 
|---|
| [e09b70] | 77 | // allows arbitrary friendship, but only if implementation is known
 | 
|---|
 | 78 | Shape::impl_ptr getShapeImpl(const Shape &shape){
 | 
|---|
 | 79 |   return shape.getImpl();
 | 
|---|
 | 80 | }
 | 
|---|
 | 81 | 
 | 
|---|
| [997784] | 82 | /***************************** Some simple Shapes ***************************/
 | 
|---|
 | 83 | 
 | 
|---|
 | 84 | Shape Everywhere(){
 | 
|---|
 | 85 |   static Shape::impl_ptr impl = Shape::impl_ptr(new Everywhere_impl());
 | 
|---|
 | 86 |   return Shape(impl);
 | 
|---|
 | 87 | }
 | 
|---|
 | 88 | 
 | 
|---|
 | 89 | Shape Nowhere(){
 | 
|---|
 | 90 |   static Shape::impl_ptr impl = Shape::impl_ptr(new Nowhere_impl());
 | 
|---|
 | 91 |   return Shape(impl);
 | 
|---|
 | 92 | }
 | 
|---|
 | 93 | 
 | 
|---|
 | 94 | /****************************** Operators ***********************************/
 | 
|---|
 | 95 | 
 | 
|---|
 | 96 | // AND
 | 
|---|
 | 97 | 
 | 
|---|
 | 98 | AndShape_impl::AndShape_impl(const Shape::impl_ptr &_lhs, const Shape::impl_ptr &_rhs) :
 | 
|---|
 | 99 |   lhs(_lhs),rhs(_rhs)
 | 
|---|
 | 100 | {}
 | 
|---|
 | 101 | 
 | 
|---|
 | 102 | AndShape_impl::~AndShape_impl(){}
 | 
|---|
 | 103 | 
 | 
|---|
 | 104 | bool AndShape_impl::isInside(const Vector &point){
 | 
|---|
 | 105 |   return lhs->isInside(point) && rhs->isInside(point);
 | 
|---|
 | 106 | }
 | 
|---|
 | 107 | 
 | 
|---|
| [5de9da] | 108 | bool AndShape_impl::isOnSurface(const Vector &point){
 | 
|---|
 | 109 |   // check the number of surfaces that this point is on
 | 
|---|
 | 110 |   int surfaces =0;
 | 
|---|
 | 111 |   surfaces += lhs->isOnSurface(point);
 | 
|---|
 | 112 |   surfaces += rhs->isOnSurface(point);
 | 
|---|
 | 113 | 
 | 
|---|
 | 114 |   switch(surfaces){
 | 
|---|
 | 115 |     case 0:
 | 
|---|
 | 116 |       return false;
 | 
|---|
 | 117 |       // no break necessary
 | 
|---|
 | 118 |     case 1:
 | 
|---|
 | 119 |       // if it is inside for the object where it does not lie on
 | 
|---|
 | 120 |       // the surface the whole point lies inside
 | 
|---|
| [ef84ca] | 121 |       return (lhs->isOnSurface(point) && rhs->isInside(point)) ||
 | 
|---|
 | 122 |              (rhs->isOnSurface(point) && lhs->isInside(point));
 | 
|---|
| [5de9da] | 123 |       // no break necessary
 | 
|---|
 | 124 |     case 2:
 | 
|---|
 | 125 |       {
 | 
|---|
 | 126 |         // it lies on both Shapes... could be an edge or an inner point
 | 
|---|
 | 127 |         // test the direction of the normals
 | 
|---|
 | 128 |         Vector direction=lhs->getNormal(point)+rhs->getNormal(point);
 | 
|---|
 | 129 |         // if the directions are opposite we lie on the inside
 | 
|---|
 | 130 |         return !direction.IsZero();
 | 
|---|
 | 131 |       }
 | 
|---|
 | 132 |       // no break necessary
 | 
|---|
 | 133 |     default:
 | 
|---|
 | 134 |       // if this happens there is something wrong
 | 
|---|
 | 135 |       ASSERT(0,"Default case should have never been used");
 | 
|---|
 | 136 |   }
 | 
|---|
 | 137 |   return false; // never reached
 | 
|---|
 | 138 | }
 | 
|---|
 | 139 | 
 | 
|---|
 | 140 | Vector AndShape_impl::getNormal(const Vector &point) throw (NotOnSurfaceException){
 | 
|---|
 | 141 |   Vector res;
 | 
|---|
 | 142 |   if(!isOnSurface(point)){
 | 
|---|
| [b94634] | 143 |     throw NotOnSurfaceException() << ShapeVector(&point);
 | 
|---|
| [5de9da] | 144 |   }
 | 
|---|
 | 145 |   res += lhs->isOnSurface(point)?lhs->getNormal(point):zeroVec;
 | 
|---|
 | 146 |   res += rhs->isOnSurface(point)?rhs->getNormal(point):zeroVec;
 | 
|---|
 | 147 |   res.Normalize();
 | 
|---|
 | 148 |   return res;
 | 
|---|
 | 149 | }
 | 
|---|
 | 150 | 
 | 
|---|
| [c6f395] | 151 | LineSegmentSet AndShape_impl::getLineIntersections(const Line &line){
 | 
|---|
 | 152 |   return intersect(lhs->getLineIntersections(line),rhs->getLineIntersections(line));
 | 
|---|
 | 153 | }
 | 
|---|
 | 154 | 
 | 
|---|
| [cfda65] | 155 | string AndShape_impl::toString(){
 | 
|---|
 | 156 |   return string("(") + lhs->toString() + string("&&") + rhs->toString() + string(")");
 | 
|---|
 | 157 | }
 | 
|---|
 | 158 | 
 | 
|---|
| [9c1c89] | 159 | std::vector<Vector> AndShape_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [23bade] | 160 |   std::vector<Vector> PointsOnSurface_lhs = lhs->getHomogeneousPointsOnSurface(N);
 | 
|---|
 | 161 |   std::vector<Vector> PointsOnSurface_rhs = rhs->getHomogeneousPointsOnSurface(N);
 | 
|---|
| [c5186e] | 162 |   std::vector<Vector> PointsOnSurface;
 | 
|---|
| [23bade] | 163 | 
 | 
|---|
 | 164 |   for (std::vector<Vector>::const_iterator iter = PointsOnSurface_lhs.begin(); iter != PointsOnSurface_lhs.end(); ++iter) {
 | 
|---|
 | 165 |     if (rhs->isInside(*iter))
 | 
|---|
 | 166 |       PointsOnSurface.push_back(*iter);
 | 
|---|
 | 167 |   }
 | 
|---|
 | 168 |   for (std::vector<Vector>::const_iterator iter = PointsOnSurface_rhs.begin(); iter != PointsOnSurface_rhs.end(); ++iter) {
 | 
|---|
 | 169 |     if (lhs->isInside(*iter))
 | 
|---|
 | 170 |       PointsOnSurface.push_back(*iter);
 | 
|---|
 | 171 |   }
 | 
|---|
 | 172 | 
 | 
|---|
| [c5186e] | 173 |   return PointsOnSurface;
 | 
|---|
 | 174 | }
 | 
|---|
 | 175 | 
 | 
|---|
 | 176 | 
 | 
|---|
| [997784] | 177 | Shape operator&&(const Shape &lhs,const Shape &rhs){
 | 
|---|
| [e09b70] | 178 |   Shape::impl_ptr newImpl = Shape::impl_ptr(new AndShape_impl(getShapeImpl(lhs),getShapeImpl(rhs)));
 | 
|---|
| [997784] | 179 |   return Shape(newImpl);
 | 
|---|
 | 180 | }
 | 
|---|
 | 181 | 
 | 
|---|
 | 182 | // OR
 | 
|---|
 | 183 | 
 | 
|---|
 | 184 | OrShape_impl::OrShape_impl(const Shape::impl_ptr &_lhs, const Shape::impl_ptr &_rhs) :
 | 
|---|
 | 185 |   lhs(_lhs),rhs(_rhs)
 | 
|---|
 | 186 | {}
 | 
|---|
 | 187 | 
 | 
|---|
 | 188 | OrShape_impl::~OrShape_impl(){}
 | 
|---|
 | 189 | 
 | 
|---|
 | 190 | bool OrShape_impl::isInside(const Vector &point){
 | 
|---|
 | 191 |   return rhs->isInside(point) || lhs->isInside(point);
 | 
|---|
 | 192 | }
 | 
|---|
 | 193 | 
 | 
|---|
| [5de9da] | 194 | bool OrShape_impl::isOnSurface(const Vector &point){
 | 
|---|
 | 195 |   // check the number of surfaces that this point is on
 | 
|---|
 | 196 |   int surfaces =0;
 | 
|---|
 | 197 |   surfaces += lhs->isOnSurface(point);
 | 
|---|
 | 198 |   surfaces += rhs->isOnSurface(point);
 | 
|---|
 | 199 | 
 | 
|---|
 | 200 |   switch(surfaces){
 | 
|---|
 | 201 |     case 0:
 | 
|---|
 | 202 |       return false;
 | 
|---|
 | 203 |       // no break necessary
 | 
|---|
 | 204 |     case 1:
 | 
|---|
 | 205 |       // if it is inside for the object where it does not lie on
 | 
|---|
 | 206 |       // the surface the whole point lies inside
 | 
|---|
| [ef84ca] | 207 |       return (lhs->isOnSurface(point) && !rhs->isInside(point)) ||
 | 
|---|
 | 208 |              (rhs->isOnSurface(point) && !lhs->isInside(point));
 | 
|---|
| [5de9da] | 209 |       // no break necessary
 | 
|---|
 | 210 |     case 2:
 | 
|---|
 | 211 |       {
 | 
|---|
 | 212 |         // it lies on both Shapes... could be an edge or an inner point
 | 
|---|
 | 213 |         // test the direction of the normals
 | 
|---|
 | 214 |         Vector direction=lhs->getNormal(point)+rhs->getNormal(point);
 | 
|---|
 | 215 |         // if the directions are opposite we lie on the inside
 | 
|---|
 | 216 |         return !direction.IsZero();
 | 
|---|
 | 217 |       }
 | 
|---|
 | 218 |       // no break necessary
 | 
|---|
 | 219 |     default:
 | 
|---|
 | 220 |       // if this happens there is something wrong
 | 
|---|
 | 221 |       ASSERT(0,"Default case should have never been used");
 | 
|---|
 | 222 |   }
 | 
|---|
 | 223 |   return false; // never reached
 | 
|---|
 | 224 | }
 | 
|---|
 | 225 | 
 | 
|---|
 | 226 | Vector OrShape_impl::getNormal(const Vector &point) throw (NotOnSurfaceException){
 | 
|---|
 | 227 |   Vector res;
 | 
|---|
 | 228 |   if(!isOnSurface(point)){
 | 
|---|
| [b94634] | 229 |     throw NotOnSurfaceException() << ShapeVector(&point);
 | 
|---|
| [5de9da] | 230 |   }
 | 
|---|
 | 231 |   res += lhs->isOnSurface(point)?lhs->getNormal(point):zeroVec;
 | 
|---|
 | 232 |   res += rhs->isOnSurface(point)?rhs->getNormal(point):zeroVec;
 | 
|---|
 | 233 |   res.Normalize();
 | 
|---|
 | 234 |   return res;
 | 
|---|
 | 235 | }
 | 
|---|
 | 236 | 
 | 
|---|
| [c6f395] | 237 | LineSegmentSet OrShape_impl::getLineIntersections(const Line &line){
 | 
|---|
 | 238 |   return merge(lhs->getLineIntersections(line),rhs->getLineIntersections(line));
 | 
|---|
 | 239 | }
 | 
|---|
 | 240 | 
 | 
|---|
| [cfda65] | 241 | string OrShape_impl::toString(){
 | 
|---|
 | 242 |   return string("(") + lhs->toString() + string("||") + rhs->toString() + string(")");
 | 
|---|
 | 243 | }
 | 
|---|
 | 244 | 
 | 
|---|
| [9c1c89] | 245 | std::vector<Vector> OrShape_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [c5186e] | 246 |   std::vector<Vector> PointsOnSurface_lhs = lhs->getHomogeneousPointsOnSurface(N);
 | 
|---|
 | 247 |   std::vector<Vector> PointsOnSurface_rhs = rhs->getHomogeneousPointsOnSurface(N);
 | 
|---|
 | 248 |   std::vector<Vector> PointsOnSurface;
 | 
|---|
 | 249 | 
 | 
|---|
 | 250 |   for (std::vector<Vector>::const_iterator iter = PointsOnSurface_lhs.begin(); iter != PointsOnSurface_lhs.end(); ++iter) {
 | 
|---|
 | 251 |     if (!rhs->isInside(*iter))
 | 
|---|
 | 252 |       PointsOnSurface.push_back(*iter);
 | 
|---|
 | 253 |   }
 | 
|---|
 | 254 |   for (std::vector<Vector>::const_iterator iter = PointsOnSurface_rhs.begin(); iter != PointsOnSurface_rhs.end(); ++iter) {
 | 
|---|
 | 255 |     if (!lhs->isInside(*iter))
 | 
|---|
 | 256 |       PointsOnSurface.push_back(*iter);
 | 
|---|
 | 257 |   }
 | 
|---|
 | 258 | 
 | 
|---|
 | 259 |   return PointsOnSurface;
 | 
|---|
 | 260 | }
 | 
|---|
 | 261 | 
 | 
|---|
| [997784] | 262 | Shape operator||(const Shape &lhs,const Shape &rhs){
 | 
|---|
| [e09b70] | 263 |   Shape::impl_ptr newImpl = Shape::impl_ptr(new OrShape_impl(getShapeImpl(lhs),getShapeImpl(rhs)));
 | 
|---|
| [997784] | 264 |   return Shape(newImpl);
 | 
|---|
 | 265 | }
 | 
|---|
 | 266 | 
 | 
|---|
 | 267 | // NOT
 | 
|---|
 | 268 | 
 | 
|---|
 | 269 | NotShape_impl::NotShape_impl(const Shape::impl_ptr &_arg) :
 | 
|---|
 | 270 |   arg(_arg)
 | 
|---|
 | 271 | {}
 | 
|---|
 | 272 | 
 | 
|---|
 | 273 | NotShape_impl::~NotShape_impl(){}
 | 
|---|
 | 274 | 
 | 
|---|
 | 275 | bool NotShape_impl::isInside(const Vector &point){
 | 
|---|
 | 276 |   return !arg->isInside(point);
 | 
|---|
 | 277 | }
 | 
|---|
 | 278 | 
 | 
|---|
| [5de9da] | 279 | bool NotShape_impl::isOnSurface(const Vector &point){
 | 
|---|
 | 280 |   return arg->isOnSurface(point);
 | 
|---|
 | 281 | }
 | 
|---|
 | 282 | 
 | 
|---|
 | 283 | Vector NotShape_impl::getNormal(const Vector &point) throw(NotOnSurfaceException){
 | 
|---|
| [b94634] | 284 |   return -1.*arg->getNormal(point);
 | 
|---|
| [5de9da] | 285 | }
 | 
|---|
 | 286 | 
 | 
|---|
| [c6f395] | 287 | LineSegmentSet NotShape_impl::getLineIntersections(const Line &line){
 | 
|---|
 | 288 |   return invert(arg->getLineIntersections(line));
 | 
|---|
 | 289 | }
 | 
|---|
 | 290 | 
 | 
|---|
| [cfda65] | 291 | string NotShape_impl::toString(){
 | 
|---|
 | 292 |   return string("!") + arg->toString();
 | 
|---|
 | 293 | }
 | 
|---|
 | 294 | 
 | 
|---|
| [9c1c89] | 295 | std::vector<Vector> NotShape_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [c5186e] | 296 |   // surfaces are the same, only normal direction is different
 | 
|---|
 | 297 |   return arg->getHomogeneousPointsOnSurface(N);
 | 
|---|
 | 298 | }
 | 
|---|
 | 299 | 
 | 
|---|
| [997784] | 300 | Shape operator!(const Shape &arg){
 | 
|---|
| [e09b70] | 301 |   Shape::impl_ptr newImpl = Shape::impl_ptr(new NotShape_impl(getShapeImpl(arg)));
 | 
|---|
| [997784] | 302 |   return Shape(newImpl);
 | 
|---|
 | 303 | }
 | 
|---|
| [cfda65] | 304 | 
 | 
|---|
 | 305 | /**************** global operations *********************************/
 | 
|---|
 | 306 | ostream &operator<<(ostream &ost,const Shape &shape){
 | 
|---|
 | 307 |   ost << shape.toString();
 | 
|---|
 | 308 |   return ost;
 | 
|---|
 | 309 | }
 | 
|---|