]> https://gitweb.dealii.org/ - dealii.git/commitdiff
GridIn::read_assimp + test.
authorLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 14:24:53 +0000 (16:24 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Sat, 16 Sep 2017 14:57:03 +0000 (16:57 +0200)
include/deal.II/base/exceptions.h
include/deal.II/grid/grid_in.h
source/grid/grid_in.cc
tests/grid/grid_in_assimp_01.cc [new file with mode: 0644]
tests/grid/grid_in_assimp_01.with_assimp=on.output [new file with mode: 0644]
tests/grid/grids/torus.obj [new file with mode: 0644]

index d70c4aa50dbd7a75ab695cb3d8a471bbea399bb9..2f36c5320d2aac3101eda646e49d7f49c4709538 100644 (file)
@@ -1113,6 +1113,14 @@ namespace StandardExceptions
                     "find a valid muparser library on your system and also did "
                     "not choose the one that comes bundled with deal.II.");
 
+  /**
+   * This function requires support for the Assimp library.
+   */
+  DeclExceptionMsg (ExcNeedsAssimp,
+                    "You are attempting to use functionality that is only available "
+                    "if deal.II was configured to use Assimp, but cmake did not "
+                    "find a valid Assimp library.");
+
 #ifdef DEAL_II_WITH_CUDA
   /**
    * This exception is raised if an error happened in a CUDA kernel.
index 87f25f377df58238574e1107fa78b2708c1b8bf1..a7908864b86803996c2fe819fb0844d343b8c09c 100644 (file)
@@ -33,16 +33,13 @@ template <int dim> struct CellData;
 /**
  * This class implements an input mechanism for grid data. It allows to read a
  * grid structure into a triangulation object. At present, UCD (unstructured
- * cell data), DB Mesh, XDA, Gmsh, Tecplot, NetCDF, UNV, VTK, and Cubit are
- * supported as input format for grid data. Any numerical data other than
+ * cell data), DB Mesh, XDA, Gmsh, Tecplot, NetCDF, UNV, VTK, ASSIMP, and Cubit
+ * are supported as input format for grid data. Any numerical data other than
  * geometric (vertex locations) and topological (how vertices form cells,
  * faces, and edges) information is ignored, but the readers for the various
  * formats generally do read information that associates material ids or
- * boundary ids to cells or faces (see
- * @ref GlossMaterialId "this"
- * and
- * @ref GlossBoundaryIndicator "this"
- * glossary entry for more information).
+ * boundary ids to cells or faces (see @ref GlossMaterialId "this" and @ref
+ * GlossBoundaryIndicator "this" glossary entry for more information).
  *
  * @note Since deal.II only supports line, quadrilateral and hexahedral
  * meshes, the functions in this class can only read meshes that consist
@@ -325,7 +322,9 @@ public:
     /// Use read_tecplot()
     tecplot,
     /// Use read_vtk()
-    vtk
+    vtk,
+    /// Use read_assimp()
+    assimp,
   };
 
   /**
@@ -455,6 +454,44 @@ public:
    */
   void read_tecplot (std::istream &in);
 
+  /**
+   * Read in a file supported by Assimp, and generate a Triangulation
+   * out of it.  If you specify a @p mesh_index, only the mesh with
+   * the given index will be extracted, otherwise all meshes which are
+   * present in the file will be used to generate the Triangulation.
+   *
+   * This function can only be used to read two-dimensional meshes (possibly
+   * embedded in three dimensions), as this is the only format supported
+   * by the supported graphical formats.
+   *
+   * If @p remove_duplicates is set to true (the default), then
+   * duplicated vertices will be removed if their distance is lower
+   * than @p tol.
+   *
+   * Only the elements compatible with the given dimension and spacedimension
+   * will be extracted from the mesh, and only those elements that are
+   * compatible with deal.II are supported. If you set
+   * `ignore_unsupported_element_types`, all the other element types are simply
+   * ignored by this algorithm. If your mesh contains a mixture of triangles
+   * and quadrilaterals, for example, only the quadrilaterals will be
+   * extracted. Careful that your mesh may not make any sense if you are mixing
+   * compatible and incompatible element types. If
+   * `ignore_unsupported_element_types` is set to `false`, then an exception is
+   * thrown when an unsupporte type is encountered.
+   *
+   * @param filename The file to read from
+   * @param mesh_index Index of the mesh within the file
+   * @param remove_duplicates Remove duplicated vertices
+   * @param tol Tolerance to use when removing vertices
+   * @param ignore_unsupported_element_types Don't throw exceptions if we
+   *        encounter unsupported types during parsing
+   */
+  void read_assimp(const std::string &filename,
+                   const unsigned int mesh_index=numbers::invalid_unsigned_int,
+                   const bool remove_duplicates=true,
+                   const double tol = 1e-12,
+                   const bool ignore_unsupported_element_types = true);
+
   /**
    * Return the standard suffix for a file in this format.
    */
@@ -649,7 +686,6 @@ void
 GridIn<3>::debug_output_grid (const std::vector<CellData<3> > &cells,
                               const std::vector<Point<3> >    &vertices,
                               std::ostream                    &out);
-
 #endif // DOXYGEN
 
 DEAL_II_NAMESPACE_CLOSE
index 4ab71c303fdbce2fb3ae5ee66647d17cde59f282..5b65eed9cfe9c4cbfd6e15c5908886c25dd249b6 100644 (file)
 #include <netcdfcpp.h>
 #endif
 
+#ifdef DEAL_II_WITH_ASSIMP
+DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Waddress-of-packed-member"
+#include <assimp/Importer.hpp>      // C++ importer interface
+#include <assimp/scene.h>           // Output data structure
+#include <assimp/postprocess.h>     // Post processing flags
+#pragma GCC diagnostic pop
+DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
+
+#undef AI_CONFIG_PP_RVC_FLAGS
+#define AI_CONFIG_PP_RVC_FLAGS           \
+  aiComponent_NORMALS |            \
+  aiComponent_TANGENTS_AND_BITANGENTS |        \
+  aiComponent_COLORS    |        \
+  aiComponent_TEXCOORDS   |        \
+  aiComponent_BONEWEIGHTS |        \
+  aiComponent_ANIMATIONS |           \
+  aiComponent_TEXTURES  |          \
+  aiComponent_LIGHTS |             \
+  aiComponent_CAMERAS |            \
+  aiComponent_MATERIALS
+#endif
+
 
 DEAL_II_NAMESPACE_OPEN
 
@@ -2542,6 +2566,153 @@ void GridIn<dim, spacedim>::read_tecplot(std::istream &)
 }
 
 
+
+template <int dim, int spacedim>
+void GridIn<dim, spacedim>::read_assimp(const std::string &filename,
+                                        const unsigned int mesh_index,
+                                        const bool remove_duplicates,
+                                        const double tol,
+                                        const bool ignore_unsupported_types)
+{
+  // Only good for surface grids.
+  AssertThrow(dim<3, ExcImpossibleInDim(dim));
+
+  // Create an istance of the Importer class
+  Assimp::Importer importer;
+
+  // And have it read the given file with some  postprocessing
+  const aiScene *scene = importer.ReadFile( filename.c_str(),
+                                            aiProcess_RemoveComponent   |
+                                            aiProcess_JoinIdenticalVertices   |
+                                            aiProcess_ImproveCacheLocality  |
+                                            aiProcess_SortByPType   |
+                                            aiProcess_OptimizeGraph   |
+                                            aiProcess_OptimizeMeshes);
+
+  // If the import failed, report it
+  AssertThrow(scene, ExcMessage(importer.GetErrorString()))
+
+  AssertThrow(scene->mNumMeshes, ExcMessage("Input file contains no meshes."));
+
+  AssertThrow((mesh_index == numbers::invalid_unsigned_int) ||
+              (mesh_index < scene->mNumMeshes),
+              ExcMessage("Too few meshes in the file."));
+
+  int start_mesh = (mesh_index == numbers::invalid_unsigned_int ?
+                    0 : mesh_index);
+  int end_mesh = (mesh_index == numbers::invalid_unsigned_int ?
+                  scene->mNumMeshes : mesh_index+1);
+
+  // Deal.II objects are created empty, and then filled with imported file.
+  std::vector<Point<spacedim> > vertices;
+  std::vector<CellData<dim> > cells;
+  SubCellData subcelldata;
+
+  // A series of counters to merge cells.
+  unsigned int v_offset=0;
+  unsigned int c_offset=0;
+
+  // The index of the mesh will be used as a material index.
+  for (int m=start_mesh; m<end_mesh; ++m)
+    {
+      const aiMesh *mesh = scene->mMeshes[m];
+
+      // Check that we know what to do with this mesh, otherwise just
+      // ignore it
+      if ( (dim == 2) && mesh->mPrimitiveTypes != aiPrimitiveType_POLYGON)
+        {
+          std::stringstream s;
+          s << "Incompatible mesh "  << m
+            << "/" << scene->mNumMeshes;
+          AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+          continue;
+        }
+      else if ( (dim == 1) && mesh->mPrimitiveTypes != aiPrimitiveType_LINE)
+        {
+          std::stringstream s;
+          s << "Incompatible mesh "  << m
+            << "/" << scene->mNumMeshes;
+          AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+          continue;
+        }
+      // Vertices
+      const unsigned int n_vertices = mesh->mNumVertices;
+      const aiVector3D *mVertices = mesh->mVertices;
+
+      // Faces
+      const unsigned int n_faces = mesh->mNumFaces;
+      const aiFace *mFaces = mesh->mFaces;
+
+      vertices.resize(v_offset+n_vertices);
+      cells.resize(c_offset+n_faces);
+
+      for (unsigned int i=0; i<n_vertices; ++i)
+        for (unsigned int d=0; d<spacedim; ++d)
+          vertices[i+v_offset][d] = mVertices[i][d];
+
+      unsigned int valid_cell = c_offset;
+      for (unsigned int i=0; i<n_faces; ++i)
+        {
+          if (mFaces[i].mNumIndices == GeometryInfo<dim>::vertices_per_cell)
+            {
+              for (unsigned int f=0; f<GeometryInfo<dim>::vertices_per_cell; ++f)
+                {
+                  cells[valid_cell].vertices[f] = mFaces[i].mIndices[f]+v_offset;
+                }
+              cells[valid_cell].material_id = (types::material_id) m;
+              ++valid_cell;
+            }
+          else
+            {
+              std::stringstream s;
+              s << "Face " << i << " of mesh "
+                << m << " has " << mFaces[i].mNumIndices
+                << " vertices. We expected only "
+                << GeometryInfo<dim>::vertices_per_cell;
+              AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+            }
+        }
+      cells.resize(valid_cell);
+
+      // The vertices are added all at once. Cells are check for
+      // validity, so only valid_cells are now present in the deal.II
+      // list of cells.
+      v_offset += n_vertices;
+      c_offset = valid_cell;
+    }
+
+  // No cells were read
+  if (cells.size() == 0)
+    return;
+
+  if (remove_duplicates)
+    {
+      // The function delete_duplicated_vertices() needs to be called more
+      // than once if a vertex is duplicated more than once. So we keep
+      // calling it until the number of vertices does not change any more.
+      unsigned int n_verts = 0;
+      while (n_verts != vertices.size())
+        {
+          n_verts = vertices.size();
+          std::vector<unsigned int> considered_vertices;
+          GridTools::delete_duplicated_vertices(vertices, cells, subcelldata,
+                                                considered_vertices, tol);
+        }
+    }
+
+  GridTools::delete_unused_vertices(vertices, cells, subcelldata);
+  if (dim == spacedim)
+    GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid
+    (vertices,cells);
+
+  GridReordering<dim, spacedim>::reorder_cells(cells);
+  if (dim == 2)
+    tria->create_triangulation_compatibility(vertices, cells, subcelldata);
+  else
+    tria->create_triangulation(vertices, cells, subcelldata);
+}
+
+
 template <int dim, int spacedim>
 void GridIn<dim, spacedim>::skip_empty_lines (std::istream &in)
 {
@@ -2820,7 +2991,7 @@ void GridIn<dim, spacedim>::read (std::istream &in,
 
     case netcdf:
       Assert(false, ExcMessage("There is no read_netcdf(istream &) function. "
-                               "Use the read(_netcdf)(string &filename) "
+                               "Use the read_netcdf)(string &filename) "
                                "functions, instead."));
       return;
 
@@ -2828,6 +2999,12 @@ void GridIn<dim, spacedim>::read (std::istream &in,
       read_tecplot (in);
       return;
 
+    case assimp:
+      Assert(false, ExcMessage("There is no read_assimp(istream &) function. "
+                               "Use the read_netcdf)(string &filename, ...) "
+                               "functions, instead."));
+      return;
+
     case Default:
       break;
     }
@@ -2928,7 +3105,7 @@ GridIn<dim, spacedim>::parse_format (const std::string &format_name)
 template <int dim, int spacedim>
 std::string GridIn<dim, spacedim>::get_format_names ()
 {
-  return "dbmesh|msh|unv|vtk|ucd|abaqus|xda|netcdf|tecplot";
+  return "dbmesh|msh|unv|vtk|ucd|abaqus|xda|netcdf|tecplot|assimp";
 }
 
 namespace
diff --git a/tests/grid/grid_in_assimp_01.cc b/tests/grid/grid_in_assimp_01.cc
new file mode 100644 (file)
index 0000000..a68efa6
--- /dev/null
@@ -0,0 +1,33 @@
+//-----------------------------------------------------------
+//
+//    Copyright (C) 2017 by the deal.II authors
+//
+//    This file is part of the deal.II library.
+//
+//    The deal.II library is free software; you can use it, redistribute
+//    it, and/or modify it under the terms of the GNU Lesser General
+//    Public License as published by the Free Software Foundation; either
+//    version 2.1 of the License, or (at your option) any later version.
+//    The full text of the license can be found in the file LICENSE at
+//    the top level of the deal.II distribution.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/grid_out.h>
+
+
+int main ()
+{
+  initlog();
+  Triangulation<2,3> tria;
+  GridIn<2,3> gridin;
+  gridin.attach_triangulation(tria);
+  gridin.read_assimp(SOURCE_DIR "/grids/torus.obj", -1, true, 1e-3);
+
+  GridOut go;
+  go.write_msh(tria, deallog.get_file_stream());
+  return 0;
+}
diff --git a/tests/grid/grid_in_assimp_01.with_assimp=on.output b/tests/grid/grid_in_assimp_01.with_assimp=on.output
new file mode 100644 (file)
index 0000000..42191dc
--- /dev/null
@@ -0,0 +1,39 @@
+
+$NOD
+16
+1  -0.364915 0.735417 -1.86297
+2  0.935085 0.735417 -0.562966
+3  0.635085 0.435417 -0.562966
+4  -0.364915 0.435417 -1.56297
+5  -0.364915 1.03542 -1.56297
+6  0.635085 1.03542 -0.562966
+7  -0.364915 0.735417 -1.26297
+8  0.335085 0.735417 -0.562966
+9  -1.66492 0.735417 -0.562966
+10  -1.36492 0.435417 -0.562966
+11  -1.36492 1.03542 -0.562966
+12  -1.06491 0.735417 -0.562966
+13  -0.364915 0.735417 0.737034
+14  -0.364915 0.435417 0.437034
+15  -0.364915 1.03542 0.437034
+16  -0.364915 0.735417 0.137034
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 3 4 1 2 
+2 3 0 0 4 6 2 1 5 
+3 3 0 0 4 8 6 5 7 
+4 3 0 0 4 3 8 7 4 
+5 3 0 0 4 4 10 9 1 
+6 3 0 0 4 5 1 9 11 
+7 3 0 0 4 7 5 11 12 
+8 3 0 0 4 4 7 12 10 
+9 3 0 0 4 10 14 13 9 
+10 3 0 0 4 11 9 13 15 
+11 3 0 0 4 12 11 15 16 
+12 3 0 0 4 10 12 16 14 
+13 3 0 0 4 14 3 2 13 
+14 3 0 0 4 15 13 2 6 
+15 3 0 0 4 16 15 6 8 
+16 3 0 0 4 14 16 8 3 
+$ENDELM
diff --git a/tests/grid/grids/torus.obj b/tests/grid/grids/torus.obj
new file mode 100644 (file)
index 0000000..f3ca4cf
--- /dev/null
@@ -0,0 +1,46 @@
+# Blender v2.74 (sub 0) OBJ File: ''
+# www.blender.org
+mtllib torus.mtl
+o SuperToroid
+v 0.935085 0.735417 -0.562966
+v 0.635085 1.035417 -0.562966
+v 0.335085 0.735417 -0.562966
+v 0.635085 0.435417 -0.562966
+v -0.364915 0.735417 -1.862966
+v -0.364915 1.035417 -1.562966
+v -0.364915 0.735417 -1.262966
+v -0.364915 0.435417 -1.562966
+v -1.664915 0.735417 -0.562966
+v -1.364915 1.035417 -0.562966
+v -1.064915 0.735417 -0.562966
+v -1.364915 0.435417 -0.562966
+v -0.364915 0.735417 0.737034
+v -0.364915 1.035417 0.437034
+v -0.364915 0.735417 0.137034
+v -0.364915 0.435417 0.437034
+vn 0.577400 -0.577400 -0.577400
+vn 0.577400 0.577400 -0.577400
+vn -0.577400 0.577400 0.577400
+vn -0.577400 -0.577400 0.577400
+vn -0.577400 -0.577400 -0.577400
+vn -0.577400 0.577400 -0.577400
+vn 0.577400 0.577400 0.577400
+vn 0.577400 -0.577400 0.577400
+usemtl None
+s off
+f 5//1 1//1 4//1 8//1
+f 1//2 5//2 6//2 2//2
+f 2//3 6//3 7//3 3//3
+f 3//4 7//4 8//4 4//4
+f 9//5 5//5 8//5 12//5
+f 5//6 9//6 10//6 6//6
+f 6//7 10//7 11//7 7//7
+f 7//8 11//8 12//8 8//8
+f 13//4 9//4 12//4 16//4
+f 9//3 13//3 14//3 10//3
+f 10//2 14//2 15//2 11//2
+f 11//1 15//1 16//1 12//1
+f 1//8 13//8 16//8 4//8
+f 13//7 1//7 2//7 14//7
+f 14//6 2//6 3//6 15//6
+f 15//5 3//5 4//5 16//5

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.