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