From 94d3e79a8617f88dc0219cfdeedfa3147833719d Mon Sep 17 00:00:00 2001 From: Adrian Kummerlaender Date: Mon, 24 Jun 2019 14:43:36 +0200 Subject: Initialize at openlb-1-3 --- src/io/stlReader.h | 290 +++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 290 insertions(+) create mode 100644 src/io/stlReader.h (limited to 'src/io/stlReader.h') diff --git a/src/io/stlReader.h b/src/io/stlReader.h new file mode 100644 index 0000000..0e44cf1 --- /dev/null +++ b/src/io/stlReader.h @@ -0,0 +1,290 @@ +/* This file is part of the OpenLB library + * + * Copyright (C) 2010-2015 Thomas Henn, Mathias J. Krause + * E-mail contact: info@openlb.net + * The most recent release of OpenLB can be downloaded at + * + * + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public + * License along with this program; if not, write to the Free + * Software Foundation, Inc., 51 Franklin Street, Fifth Floor, + * Boston, MA 02110-1301, USA. +*/ + +/** \file + * Input in STL format -- header file. + */ + +#ifndef STL_READER_H +#define STL_READER_H + +#include +#include +#include +#include +#include + +#include "communication/loadBalancer.h" +#include "geometry/cuboidGeometry3D.h" +#include "functors/analytical/indicator/indicatorF3D.h" +#include "functors/analytical/indicator/indicatorBaseF3D.h" +#include "utilities/vectorHelpers.h" +#include "octree.h" +#include "core/vector.h" + + +using namespace olb::util; + +/// All OpenLB code is contained in this namespace. +namespace olb { + +template +class Octree; + +template +struct STLpoint { + /// Constructor constructs + STLpoint() : r() {}; + /// Operator= equals + STLpoint& operator=(STLpoint const& rhs) + { + r = rhs.r; + return *this; + }; + /// CopyConstructor copies + STLpoint(STLpoint const& rhs):r(rhs.r) {}; + + /// Point coordinates in SI units + Vector r; +}; + +template +struct STLtriangle { + /** Test intersection between ray and triangle + * \param pt Raypoint + * \param dir Direction + * \param q Point of intersection (if intersection occurs) + * \param alpha Explained in http://www.uninformativ.de/bin/RaytracingSchnitttests-76a577a-CC-BY.pdf page 7-12 + * q = pt + alpha * dir + * \param rad It's complicated. Imagine you have a sphere with radius rad moving a long the ray. Then q becomes the first point of the sphere to touch the triangle. + */ + bool testRayIntersect(const Vector& pt,const Vector& dir, Vector& q, T& alpha, const T& rad = T(), bool print = false); + Vector closestPtPointTriangle(const Vector& pt) const; + + /// A triangle contains 3 Points + std::vector > point; + + /// normal of triangle + Vector normal; + + /// variables explained in http://www.uninformativ.de/bin/RaytracingSchnitttests-76a577a-CC-BY.pdf page 7-12 + /// precomputed for speedup + Vector uBeta, uGamma; + T d, kBeta, kGamma; + +public: + /// Constructor constructs + STLtriangle():point(3, STLpoint()), normal(T()), uBeta(T()), uGamma(T()), d(T()), kBeta(T()), kGamma(T()) {}; + /// CopyConstructor copies + STLtriangle(STLtriangle const& tri):point(tri.point), normal(tri.normal), uBeta(tri.uBeta), uGamma(tri.uGamma), d(tri.d), kBeta(tri.kBeta), kGamma(tri.kGamma) {}; + /// Operator= equals + STLtriangle& operator=(STLtriangle const& tri) + { + point = tri.point; + normal = tri.normal; + uBeta = tri.uBeta; + uGamma = tri.uGamma; + d = tri.d; + kBeta = tri.kBeta; + kGamma = tri.kGamma; + return *this; + }; + + ~STLtriangle() {}; + + /// Initializes triangle and precomputes member variables. + void init(); + /// Return write access to normal + inline Vector& getNormal() + { + return normal; + } + /// Returns Pt0-Pt1 + std::vector getE0(); + /// Returns Pt0-Pt2 + std::vector getE1(); +}; + +template +class STLmesh { + /// Computes distance squared betwenn p1 and p2 + T distPoints(STLpoint& p1, STLpoint& p2); + /// Filename + const std::string _fName; + /// Vector of Triangles + std::vector > _triangles; + /// Min and Max points of axis aligned bounding box coordinate in SI units + Vector _min, _max; + /// largest squared length of edge of all triangles + T _maxDist2; + /// OstreamManager + mutable OstreamManager clout; + +public: + /** + * Constructs a new STLmesh from a file + * \param Filename - Filename + * \param stlSize - Conversion factor for STL (e.g. STL in mm stlSize=10^-3) + */ + STLmesh(std::string, T stlSize = 1.); + + /** + * Constructs a new STLmesh from a file + * \param Filename - Filename + * \param stlSize - Conversion factor for STL (e.g. STL in mm stlSize=10^-3) + */ + STLmesh(const std::vector> meshPoints, T stlSize = 1.); + + /// Returns reference to a triangle + inline STLtriangle& getTri(unsigned int i) + { + return _triangles[i]; + } + /// Returns number of triangles + inline unsigned int triangleSize() const + { + return _triangles.size(); + } + /// Returns _min + inline Vector& getMin() + { + return _min; + }; + /// Returns _max + inline Vector& getMax() + { + return _max; + }; + /// Returns maxDist squared + inline float maxDist2() const + { + return _maxDist2; + } + /// Prints console output + void print(bool full = false); + /// Writes STL mesh in Si units + void write(std::string fName); + /// Compute intersection between Ray and set of triangles; returns true if intersection is found + bool testRayIntersect(const std::set& tris, const Vector& pt,const Vector& dir, Vector& q, T& alpha); +}; + + +template +class STLreader : public IndicatorF3D { +private: + /* + * Old indicate function (slower, more stable) + * Define three rays (X-, Y-, Z-direction) for each leaf and count intersections + * with STL for each ray. Odd number of intersection means inside (Majority vote). + */ + void indicate1(); + /* + * New indicate function (faster, less stable) + * Define ray in Z-direction for each Voxel in XY-layer. Indicate all nodes on the fly. + */ + void indicate2(); + /* + * Double ray approach: two times (X-, Y-, Z-direction) for each leaf. + * Could be use to deal with double layer triangles and face intersections. + */ + void indicate3(); + /// Size of the smallest voxel + T _voxelSize; + /// Factor to get Si unit (m), i.e. "0.001" means mm + T _stlSize; + /// Overlap increases Octree radius by _overlap + T _overlap; + /// Pointer to tree + Octree* _tree; + /// The filename + const std::string _fName; + /// The mesh + STLmesh _mesh; + /// Variable for output + bool _verbose; + /// The OstreamManager + mutable OstreamManager clout; + +public: + /** + * Constructs a new STLreader from a file + * \param fName The STL file name + * \param voxelSize Voxelsize in SI units + * \param stlSize Conversion factor for STL (e.g. STL in mm stlSize=10^-3) + * \param method Choose indication method + * 0: fast, less stable + * 1: slow, more stable (for untight STLs) + * \param verbose Get additional information. + */ + + STLreader(const std::string fName, T voxelSize, T stlSize=1, unsigned short int method=2, + bool verbose = false, T overlap=0., T max=0.); + /** + * Constructs a new STLreader from a file + * \param fName The STL file name + * \param voxelSize Voxelsize in SI units + * \param stlSize Conversion factor for STL (e.g. STL in mm stlSize=10^-3) + * \param method Choose indication method + * 0: fast, less stable + * 1: slow, more stable (for untight STLs) + * \param verbose Get additional information. + */ + STLreader(const std::vector> meshPoints, T voxelSize, T stlSize=1, unsigned short int method=2, + bool verbose = false, T overlap=0., T max=0.); + + ~STLreader() override; + /// Returns whether node is inside or not. + bool operator() (bool output[], const T input[]) override; + + /// Computes distance to closest triangle intersection + bool distance(T& distance,const Vector& origin, const Vector& direction, int iC=-1) override; + + /// Prints console output + void print(); + + /// Writes STL mesh in Si units + void writeSTL(std::string stlName=""); + + /// Writes Octree + void writeOctree(); + + /// Rearranges normals of triangles to point outside of geometry + void setNormalsOutside(); + + /// Returns tree + inline Octree* getTree() const + { + return _tree; + }; + + + /// Returns mesh + inline STLmesh& getMesh() + { + return _mesh; + }; +}; + +} // namespace olb + +#endif -- cgit v1.2.3