| 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 |  | 
|---|
| 8 | /* | 
|---|
| 9 | * ShapeOps.cpp | 
|---|
| 10 | * | 
|---|
| 11 | *  Created on: Jun 18, 2010 | 
|---|
| 12 | *      Author: crueger | 
|---|
| 13 | */ | 
|---|
| 14 |  | 
|---|
| 15 | // include config.h | 
|---|
| 16 | #ifdef HAVE_CONFIG_H | 
|---|
| 17 | #include <config.h> | 
|---|
| 18 | #endif | 
|---|
| 19 |  | 
|---|
| 20 | #include "CodePatterns/MemDebug.hpp" | 
|---|
| 21 |  | 
|---|
| 22 | #include "Shapes/ShapeExceptions.hpp" | 
|---|
| 23 | #include "Shapes/ShapeOps.hpp" | 
|---|
| 24 | #include "Shapes/ShapeOps_impl.hpp" | 
|---|
| 25 |  | 
|---|
| 26 | #include "LinearAlgebra/Vector.hpp" | 
|---|
| 27 | #include "CodePatterns/Assert.hpp" | 
|---|
| 28 |  | 
|---|
| 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)){ | 
|---|
| 47 | throw NotOnSurfaceException() << ShapeVector(&helper); | 
|---|
| 48 | } | 
|---|
| 49 | Vector res = translateOutNormal(arg->getNormal(helper)); | 
|---|
| 50 | res.Normalize(); | 
|---|
| 51 | return res; | 
|---|
| 52 | } | 
|---|
| 53 |  | 
|---|
| 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 |  | 
|---|
| 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{ | 
|---|
| 78 | return arg; | 
|---|
| 79 | } | 
|---|
| 80 |  | 
|---|
| 81 | /********************* Resize ********************/ | 
|---|
| 82 |  | 
|---|
| 83 | Resize_impl::Resize_impl(const Shape::impl_ptr &_arg,double _size) : | 
|---|
| 84 | ShapeOpsBase_impl(_arg), size(_size) | 
|---|
| 85 | { | 
|---|
| 86 | ASSERT(size>0,"Cannot resize a Shape to size zero or below"); | 
|---|
| 87 | } | 
|---|
| 88 |  | 
|---|
| 89 | Resize_impl::~Resize_impl(){} | 
|---|
| 90 |  | 
|---|
| 91 | bool Resize_impl::isInside(const Vector& point){ | 
|---|
| 92 | return getArg()->isInside((1/size) * point); | 
|---|
| 93 | } | 
|---|
| 94 |  | 
|---|
| 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; | 
|---|
| 105 | } | 
|---|
| 106 |  | 
|---|
| 107 | string Resize_impl::toString(){ | 
|---|
| 108 | stringstream sstr; | 
|---|
| 109 | sstr << "resize(" << getArg()->toString() << "," << size << ")"; | 
|---|
| 110 | return sstr.str(); | 
|---|
| 111 | } | 
|---|
| 112 |  | 
|---|
| 113 | std::vector<Vector> Resize_impl::getHomogeneousPointsOnSurface(const size_t N) const { | 
|---|
| 114 | std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N); | 
|---|
| 115 | for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) { | 
|---|
| 116 | *iter *= size; | 
|---|
| 117 | } | 
|---|
| 118 | return PointsOnSurface; | 
|---|
| 119 | } | 
|---|
| 120 |  | 
|---|
| 121 |  | 
|---|
| 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 |  | 
|---|
| 127 | /*************************** translate *******************/ | 
|---|
| 128 |  | 
|---|
| 129 | Translate_impl::Translate_impl(const Shape::impl_ptr &_arg, const Vector &_offset) : | 
|---|
| 130 | ShapeOpsBase_impl(_arg),offset(_offset) | 
|---|
| 131 | {} | 
|---|
| 132 |  | 
|---|
| 133 | Translate_impl::~Translate_impl(){} | 
|---|
| 134 |  | 
|---|
| 135 | bool Translate_impl::isInside(const Vector& point){ | 
|---|
| 136 | return getArg()->isInside(point-offset); | 
|---|
| 137 | } | 
|---|
| 138 |  | 
|---|
| 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; | 
|---|
| 149 | } | 
|---|
| 150 |  | 
|---|
| 151 | string Translate_impl::toString(){ | 
|---|
| 152 | stringstream sstr; | 
|---|
| 153 | sstr << "translate(" << getArg()->toString() << "," << offset << ")"; | 
|---|
| 154 | return sstr.str(); | 
|---|
| 155 | } | 
|---|
| 156 |  | 
|---|
| 157 | std::vector<Vector> Translate_impl::getHomogeneousPointsOnSurface(const size_t N) const { | 
|---|
| 158 | std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N); | 
|---|
| 159 | for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) { | 
|---|
| 160 | *iter += offset; | 
|---|
| 161 | } | 
|---|
| 162 | return PointsOnSurface; | 
|---|
| 163 | } | 
|---|
| 164 |  | 
|---|
| 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 | } | 
|---|
| 169 |  | 
|---|
| 170 | /*********************** stretch ******************/ | 
|---|
| 171 |  | 
|---|
| 172 | Stretch_impl::Stretch_impl(const Shape::impl_ptr &_arg, const Vector &_factors) : | 
|---|
| 173 | ShapeOpsBase_impl(_arg),factors(_factors) | 
|---|
| 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); | 
|---|
| 188 | return getArg()->isInside(helper); | 
|---|
| 189 | } | 
|---|
| 190 |  | 
|---|
| 191 | Vector Stretch_impl::translateIn(const Vector& point){ | 
|---|
| 192 | Vector helper=point; | 
|---|
| 193 | helper.ScaleAll(reciFactors); | 
|---|
| 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; | 
|---|
| 213 | } | 
|---|
| 214 |  | 
|---|
| 215 | string Stretch_impl::toString(){ | 
|---|
| 216 | stringstream sstr; | 
|---|
| 217 | sstr << "stretch(" << getArg()->toString() << "," << factors << ")"; | 
|---|
| 218 | return sstr.str(); | 
|---|
| 219 | } | 
|---|
| 220 |  | 
|---|
| 221 | std::vector<Vector> Stretch_impl::getHomogeneousPointsOnSurface(const size_t N) const { | 
|---|
| 222 | std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N); | 
|---|
| 223 | for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) { | 
|---|
| 224 | (*iter).ScaleAll(reciFactors); | 
|---|
| 225 | } | 
|---|
| 226 | return PointsOnSurface; | 
|---|
| 227 | } | 
|---|
| 228 |  | 
|---|
| 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 |  | 
|---|
| 236 | Transform_impl::Transform_impl(const Shape::impl_ptr &_arg, const RealSpaceMatrix &_transformation) : | 
|---|
| 237 | ShapeOpsBase_impl(_arg),transformation(_transformation) | 
|---|
| 238 | { | 
|---|
| 239 | transformationInv = transformation.invert(); | 
|---|
| 240 | } | 
|---|
| 241 |  | 
|---|
| 242 | Transform_impl::~Transform_impl(){} | 
|---|
| 243 |  | 
|---|
| 244 | bool Transform_impl::isInside(const Vector& point){ | 
|---|
| 245 | return getArg()->isInside(transformationInv * point); | 
|---|
| 246 | } | 
|---|
| 247 |  | 
|---|
| 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){ | 
|---|
| 257 | RealSpaceMatrix mat = transformation.invert().transpose(); | 
|---|
| 258 | return mat * point; | 
|---|
| 259 | } | 
|---|
| 260 |  | 
|---|
| 261 | string Transform_impl::toString(){ | 
|---|
| 262 | stringstream sstr; | 
|---|
| 263 | sstr << "transform(" << getArg()->toString() << "," << transformation << ")"; | 
|---|
| 264 | return sstr.str(); | 
|---|
| 265 | } | 
|---|
| 266 |  | 
|---|
| 267 | std::vector<Vector> Transform_impl::getHomogeneousPointsOnSurface(const size_t N) const { | 
|---|
| 268 | std::vector<Vector> PointsOnSurface = getArg()->getHomogeneousPointsOnSurface(N); | 
|---|
| 269 | for(std::vector<Vector>::iterator iter = PointsOnSurface.begin(); iter != PointsOnSurface.end(); ++iter) { | 
|---|
| 270 | *iter = transformation * (*iter); | 
|---|
| 271 | } | 
|---|
| 272 | return PointsOnSurface; | 
|---|
| 273 | } | 
|---|
| 274 |  | 
|---|
| 275 | Shape transform(const Shape &arg, const RealSpaceMatrix &transformation){ | 
|---|
| 276 | Shape::impl_ptr impl = Shape::impl_ptr(new Transform_impl(getShapeImpl(arg),transformation)); | 
|---|
| 277 | return Shape(impl); | 
|---|
| 278 | } | 
|---|