| [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 | 
 | 
|---|
| [e09b70] | 8 | /*
 | 
|---|
 | 9 |  * ShapeOps.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"
 | 
|---|
| [bbbad5] | 21 | 
 | 
|---|
| [b94634] | 22 | #include "Shapes/ShapeExceptions.hpp"
 | 
|---|
| [e09b70] | 23 | #include "Shapes/ShapeOps.hpp"
 | 
|---|
 | 24 | #include "Shapes/ShapeOps_impl.hpp"
 | 
|---|
 | 25 | 
 | 
|---|
| [6c438f] | 26 | #include "LinearAlgebra/Vector.hpp"
 | 
|---|
| [ad011c] | 27 | #include "CodePatterns/Assert.hpp"
 | 
|---|
| [394529] | 28 | 
 | 
|---|
| [5de9da] | 29 | /*************** Base case ***********************/
 | 
|---|
 | 30 | 
 | 
|---|
 | 31 | ShapeOpsBase_impl::ShapeOpsBase_impl(const Shape::impl_ptr &_arg) :
 | 
|---|
 | 32 |   arg(_arg){}
 | 
|---|
 | 33 | 
 | 
|---|
 | 34 | ShapeOpsBase_impl::~ShapeOpsBase_impl(){}
 | 
|---|
 | 35 | 
 | 
|---|
 | 36 | bool ShapeOpsBase_impl::isInside(const Vector &point){
 | 
|---|
 | 37 |   return arg->isInside(translateIn(point));
 | 
|---|
 | 38 | }
 | 
|---|
 | 39 | 
 | 
|---|
 | 40 | bool ShapeOpsBase_impl::isOnSurface(const Vector &point){
 | 
|---|
 | 41 |   return arg->isOnSurface(translateIn(point));
 | 
|---|
 | 42 | }
 | 
|---|
 | 43 | 
 | 
|---|
 | 44 | Vector ShapeOpsBase_impl::getNormal(const Vector &point) throw (NotOnSurfaceException){
 | 
|---|
 | 45 |   Vector helper = translateIn(point);
 | 
|---|
 | 46 |   if(!arg->isOnSurface(helper)){
 | 
|---|
| [b94634] | 47 |     throw NotOnSurfaceException() << ShapeVector(&helper);
 | 
|---|
| [5de9da] | 48 |   }
 | 
|---|
| [f12805] | 49 |   Vector res = translateOutNormal(arg->getNormal(helper));
 | 
|---|
 | 50 |   res.Normalize();
 | 
|---|
 | 51 |   return res;
 | 
|---|
| [5de9da] | 52 | }
 | 
|---|
 | 53 | 
 | 
|---|
| [c6f395] | 54 | LineSegmentSet ShapeOpsBase_impl::getLineIntersections(const Line &line){
 | 
|---|
 | 55 |   Line newLine(translateIn(line.getOrigin()),translateIn(line.getDirection()));
 | 
|---|
 | 56 |   LineSegmentSet res(line);
 | 
|---|
 | 57 |   LineSegmentSet helper = getArg()->getLineIntersections(newLine);
 | 
|---|
 | 58 |   for(LineSegmentSet::iterator iter = helper.begin();iter!=helper.end();++iter){
 | 
|---|
 | 59 |     LinePoint lpBegin = iter->getBegin();
 | 
|---|
 | 60 |     LinePoint lpEnd = iter->getBegin();
 | 
|---|
 | 61 |     // translate both linepoints
 | 
|---|
 | 62 |     lpBegin = lpBegin.isNegInfinity()?
 | 
|---|
 | 63 |                 line.negEndpoint():
 | 
|---|
 | 64 |                 line.getLinePoint(translateOutPos(lpBegin.getPoint()));
 | 
|---|
 | 65 |     lpEnd = lpEnd.isPosInfinity()?
 | 
|---|
 | 66 |               line.posEndpoint():
 | 
|---|
 | 67 |               line.getLinePoint(translateOutPos(lpEnd.getPoint()));
 | 
|---|
 | 68 |     res.insert(LineSegment(lpBegin,lpEnd));
 | 
|---|
 | 69 |   }
 | 
|---|
 | 70 |   return res;
 | 
|---|
 | 71 | }
 | 
|---|
 | 72 | 
 | 
|---|
| [6c438f] | 73 | std::vector<Vector> ShapeOpsBase_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
 | 74 |   return getArg()->getHomogeneousPointsOnSurface(N);;
 | 
|---|
 | 75 | }
 | 
|---|
 | 76 | 
 | 
|---|
 | 77 | Shape::impl_ptr ShapeOpsBase_impl::getArg() const{
 | 
|---|
| [cfda65] | 78 |   return arg;
 | 
|---|
 | 79 | }
 | 
|---|
 | 80 | 
 | 
|---|
| [5e588b5] | 81 | /********************* Resize ********************/
 | 
|---|
 | 82 | 
 | 
|---|
| [e09b70] | 83 | Resize_impl::Resize_impl(const Shape::impl_ptr &_arg,double _size) :
 | 
|---|
| [5de9da] | 84 |   ShapeOpsBase_impl(_arg), size(_size)
 | 
|---|
| [394529] | 85 | {
 | 
|---|
 | 86 |   ASSERT(size>0,"Cannot resize a Shape to size zero or below");
 | 
|---|
 | 87 | }
 | 
|---|
| [e09b70] | 88 | 
 | 
|---|
 | 89 | Resize_impl::~Resize_impl(){}
 | 
|---|
 | 90 | 
 | 
|---|
 | 91 | bool Resize_impl::isInside(const Vector& point){
 | 
|---|
| [6c438f] | 92 |   return getArg()->isInside((1/size) * point);
 | 
|---|
 | 93 | }
 | 
|---|
 | 94 | 
 | 
|---|
| [5de9da] | 95 | Vector Resize_impl::translateIn(const Vector& point){
 | 
|---|
 | 96 |   return (1/size) * point;
 | 
|---|
 | 97 | }
 | 
|---|
 | 98 | 
 | 
|---|
 | 99 | Vector Resize_impl::translateOutPos(const Vector& point){
 | 
|---|
 | 100 |   return size * point;
 | 
|---|
 | 101 | }
 | 
|---|
 | 102 | 
 | 
|---|
 | 103 | Vector Resize_impl::translateOutNormal(const Vector& point){
 | 
|---|
 | 104 |   return point;
 | 
|---|
| [e09b70] | 105 | }
 | 
|---|
 | 106 | 
 | 
|---|
| [cfda65] | 107 | string Resize_impl::toString(){
 | 
|---|
 | 108 |   stringstream sstr;
 | 
|---|
 | 109 |   sstr << "resize(" << getArg()->toString() << "," << size << ")";
 | 
|---|
 | 110 |   return sstr.str();
 | 
|---|
 | 111 | }
 | 
|---|
 | 112 | 
 | 
|---|
| [9c1c89] | 113 | std::vector<Vector> Resize_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [6c438f] | 114 |   std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N);
 | 
|---|
| [c5186e] | 115 |   for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) {
 | 
|---|
 | 116 |     *iter *= size;
 | 
|---|
 | 117 |   }
 | 
|---|
 | 118 |   return PointsOnSurface;
 | 
|---|
 | 119 | }
 | 
|---|
 | 120 | 
 | 
|---|
 | 121 | 
 | 
|---|
| [e09b70] | 122 | Shape resize(const Shape &arg,double size){
 | 
|---|
 | 123 |   Shape::impl_ptr impl = Shape::impl_ptr(new Resize_impl(getShapeImpl(arg),size));
 | 
|---|
 | 124 |   return Shape(impl);
 | 
|---|
 | 125 | }
 | 
|---|
 | 126 | 
 | 
|---|
| [5e588b5] | 127 | /*************************** translate *******************/
 | 
|---|
 | 128 | 
 | 
|---|
| [e09b70] | 129 | Translate_impl::Translate_impl(const Shape::impl_ptr &_arg, const Vector &_offset) :
 | 
|---|
| [5de9da] | 130 |   ShapeOpsBase_impl(_arg),offset(_offset)
 | 
|---|
| [e09b70] | 131 | {}
 | 
|---|
 | 132 | 
 | 
|---|
 | 133 | Translate_impl::~Translate_impl(){}
 | 
|---|
 | 134 | 
 | 
|---|
 | 135 | bool Translate_impl::isInside(const Vector& point){
 | 
|---|
| [6c438f] | 136 |   return getArg()->isInside(point-offset);
 | 
|---|
 | 137 | }
 | 
|---|
 | 138 | 
 | 
|---|
| [5de9da] | 139 | Vector Translate_impl::translateIn(const Vector& point){
 | 
|---|
 | 140 |   return point-offset;
 | 
|---|
 | 141 | }
 | 
|---|
 | 142 | 
 | 
|---|
 | 143 | Vector Translate_impl::translateOutPos(const Vector& point){
 | 
|---|
 | 144 |   return point+offset;
 | 
|---|
 | 145 | }
 | 
|---|
 | 146 | 
 | 
|---|
 | 147 | Vector Translate_impl::translateOutNormal(const Vector& point){
 | 
|---|
 | 148 |   return point;
 | 
|---|
| [e09b70] | 149 | }
 | 
|---|
 | 150 | 
 | 
|---|
| [cfda65] | 151 | string Translate_impl::toString(){
 | 
|---|
 | 152 |   stringstream sstr;
 | 
|---|
 | 153 |   sstr << "translate(" << getArg()->toString() << "," << offset << ")";
 | 
|---|
| [79dd0e] | 154 |   return sstr.str();
 | 
|---|
| [cfda65] | 155 | }
 | 
|---|
 | 156 | 
 | 
|---|
| [9c1c89] | 157 | std::vector<Vector> Translate_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [6c438f] | 158 |   std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N);
 | 
|---|
| [c5186e] | 159 |   for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) {
 | 
|---|
 | 160 |     *iter += offset;
 | 
|---|
 | 161 |   }
 | 
|---|
 | 162 |   return PointsOnSurface;
 | 
|---|
 | 163 | }
 | 
|---|
 | 164 | 
 | 
|---|
| [e09b70] | 165 | Shape translate(const Shape &arg, const Vector &offset){
 | 
|---|
 | 166 |   Shape::impl_ptr impl = Shape::impl_ptr(new Translate_impl(getShapeImpl(arg),offset));
 | 
|---|
 | 167 |   return Shape(impl);
 | 
|---|
 | 168 | }
 | 
|---|
| [5e588b5] | 169 | 
 | 
|---|
 | 170 | /*********************** stretch ******************/
 | 
|---|
 | 171 | 
 | 
|---|
 | 172 | Stretch_impl::Stretch_impl(const Shape::impl_ptr &_arg, const Vector &_factors) :
 | 
|---|
| [5de9da] | 173 |   ShapeOpsBase_impl(_arg),factors(_factors)
 | 
|---|
| [5e588b5] | 174 | {
 | 
|---|
 | 175 |   ASSERT(factors[0]>0,"cannot stretch a shape by a negative amount");
 | 
|---|
 | 176 |   ASSERT(factors[1]>0,"cannot stretch a shape by a negative amount");
 | 
|---|
 | 177 |   ASSERT(factors[2]>0,"cannot stretch a shape by a negative amount");
 | 
|---|
 | 178 |   for(int i = NDIM;i--;){
 | 
|---|
 | 179 |     reciFactors[i] = 1/factors[i];
 | 
|---|
 | 180 |   }
 | 
|---|
 | 181 | }
 | 
|---|
 | 182 | 
 | 
|---|
 | 183 | Stretch_impl::~Stretch_impl(){}
 | 
|---|
 | 184 | 
 | 
|---|
 | 185 | bool Stretch_impl::isInside(const Vector& point){
 | 
|---|
 | 186 |   Vector helper=point;
 | 
|---|
 | 187 |   helper.ScaleAll(reciFactors);
 | 
|---|
| [6c438f] | 188 |   return getArg()->isInside(helper);
 | 
|---|
 | 189 | }
 | 
|---|
 | 190 | 
 | 
|---|
| [5de9da] | 191 | Vector Stretch_impl::translateIn(const Vector& point){
 | 
|---|
| [5e588b5] | 192 |   Vector helper=point;
 | 
|---|
 | 193 |   helper.ScaleAll(reciFactors);
 | 
|---|
| [5de9da] | 194 |   return helper;
 | 
|---|
 | 195 | }
 | 
|---|
 | 196 | 
 | 
|---|
 | 197 | Vector Stretch_impl::translateOutPos(const Vector& point){
 | 
|---|
 | 198 |   Vector helper=point;
 | 
|---|
 | 199 |   helper.ScaleAll(factors);
 | 
|---|
 | 200 |   return helper;
 | 
|---|
 | 201 | }
 | 
|---|
 | 202 | 
 | 
|---|
 | 203 | Vector Stretch_impl::translateOutNormal(const Vector& point){
 | 
|---|
 | 204 |   Vector helper=point;
 | 
|---|
 | 205 |   // the normalFactors are derived from appearances of the factors
 | 
|---|
 | 206 |   // with in the vectorproduct
 | 
|---|
 | 207 |   Vector normalFactors;
 | 
|---|
 | 208 |   normalFactors[0]=factors[1]*factors[2];
 | 
|---|
 | 209 |   normalFactors[1]=factors[0]*factors[2];
 | 
|---|
 | 210 |   normalFactors[2]=factors[0]*factors[1];
 | 
|---|
 | 211 |   helper.ScaleAll(normalFactors);
 | 
|---|
 | 212 |   return helper;
 | 
|---|
| [5e588b5] | 213 | }
 | 
|---|
 | 214 | 
 | 
|---|
| [cfda65] | 215 | string Stretch_impl::toString(){
 | 
|---|
 | 216 |   stringstream sstr;
 | 
|---|
 | 217 |   sstr << "stretch(" << getArg()->toString() << "," << factors << ")";
 | 
|---|
 | 218 |   return sstr.str();
 | 
|---|
 | 219 | }
 | 
|---|
 | 220 | 
 | 
|---|
| [9c1c89] | 221 | std::vector<Vector> Stretch_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [6c438f] | 222 |   std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N);
 | 
|---|
| [c5186e] | 223 |   for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) {
 | 
|---|
 | 224 |     (*iter).ScaleAll(reciFactors);
 | 
|---|
 | 225 |   }
 | 
|---|
 | 226 |   return PointsOnSurface;
 | 
|---|
 | 227 | }
 | 
|---|
 | 228 | 
 | 
|---|
| [5e588b5] | 229 | Shape stretch(const Shape &arg, const Vector &factors){
 | 
|---|
 | 230 |   Shape::impl_ptr impl = Shape::impl_ptr(new Stretch_impl(getShapeImpl(arg),factors));
 | 
|---|
 | 231 |   return Shape(impl);
 | 
|---|
 | 232 | }
 | 
|---|
 | 233 | 
 | 
|---|
 | 234 | /************************* transform *****************/
 | 
|---|
 | 235 | 
 | 
|---|
| [cca9ef] | 236 | Transform_impl::Transform_impl(const Shape::impl_ptr &_arg, const RealSpaceMatrix &_transformation) :
 | 
|---|
| [5de9da] | 237 |   ShapeOpsBase_impl(_arg),transformation(_transformation)
 | 
|---|
| [5e588b5] | 238 | {
 | 
|---|
 | 239 |   transformationInv = transformation.invert();
 | 
|---|
 | 240 | }
 | 
|---|
 | 241 | 
 | 
|---|
 | 242 | Transform_impl::~Transform_impl(){}
 | 
|---|
 | 243 | 
 | 
|---|
 | 244 | bool Transform_impl::isInside(const Vector& point){
 | 
|---|
| [6c438f] | 245 |   return getArg()->isInside(transformationInv * point);
 | 
|---|
 | 246 | }
 | 
|---|
 | 247 | 
 | 
|---|
| [5de9da] | 248 | Vector Transform_impl::translateIn(const Vector& point){
 | 
|---|
 | 249 |   return transformationInv * point;
 | 
|---|
 | 250 | }
 | 
|---|
 | 251 | 
 | 
|---|
 | 252 | Vector Transform_impl::translateOutPos(const Vector& point){
 | 
|---|
 | 253 |   return transformation * point;
 | 
|---|
 | 254 | }
 | 
|---|
 | 255 | 
 | 
|---|
 | 256 | Vector Transform_impl::translateOutNormal(const Vector& point){
 | 
|---|
| [cca9ef] | 257 |   RealSpaceMatrix mat = transformation.invert().transpose();
 | 
|---|
| [5de9da] | 258 |   return mat * point;
 | 
|---|
| [5e588b5] | 259 | }
 | 
|---|
 | 260 | 
 | 
|---|
| [cfda65] | 261 | string Transform_impl::toString(){
 | 
|---|
 | 262 |   stringstream sstr;
 | 
|---|
 | 263 |   sstr << "transform(" << getArg()->toString() << "," << transformation << ")";
 | 
|---|
 | 264 |   return sstr.str();
 | 
|---|
 | 265 | }
 | 
|---|
 | 266 | 
 | 
|---|
| [9c1c89] | 267 | std::vector<Vector> Transform_impl::getHomogeneousPointsOnSurface(const size_t N) const {
 | 
|---|
| [6c438f] | 268 |   std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N);
 | 
|---|
| [c5186e] | 269 |   for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) {
 | 
|---|
 | 270 |     *iter = transformation * (*iter);
 | 
|---|
 | 271 |   }
 | 
|---|
 | 272 |   return PointsOnSurface;
 | 
|---|
 | 273 | }
 | 
|---|
 | 274 | 
 | 
|---|
| [cca9ef] | 275 | Shape transform(const Shape &arg, const RealSpaceMatrix &transformation){
 | 
|---|
| [5e588b5] | 276 |   Shape::impl_ptr impl = Shape::impl_ptr(new Transform_impl(getShapeImpl(arg),transformation));
 | 
|---|
 | 277 |   return Shape(impl);
 | 
|---|
 | 278 | }
 | 
|---|