From e762573f800d840fa1b4397c79d2ddd23e0d7cd7 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Wed, 8 Oct 2008 18:52:54 +0000
Subject: [PATCH] simple VTK reader

---
 Common/OpenFile.cpp   |   2 +-
 Geo/GModel.h          |   2 +-
 Geo/GModelIO_Mesh.cpp | 113 +++++++++++++++++++++++++++++++++++++++++-
 Geo/MElement.h        |   2 +-
 4 files changed, 114 insertions(+), 5 deletions(-)

diff --git a/Common/OpenFile.cpp b/Common/OpenFile.cpp
index 751b1e154e..c995201643 100644
--- a/Common/OpenFile.cpp
+++ b/Common/OpenFile.cpp
@@ -307,7 +307,7 @@ int MergeFile(const char *name, int warn_if_missing)
     status = m->readUNV(name);
   }
   else if(!strcmp(ext, ".vtk") || !strcmp(ext, ".VTK")){
-    status = m->readVTK(name);
+    status = m->readVTK(name, CTX.big_endian);
   }
   else if(!strcmp(ext, ".wrl") || !strcmp(ext, ".WRL") || 
           !strcmp(ext, ".vrml") || !strcmp(ext, ".VRML") ||
diff --git a/Geo/GModel.h b/Geo/GModel.h
index 60dec101f1..394ccd7807 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -320,7 +320,7 @@ class GModel
 	       bool saveAll=false, double scalingFactor=1.0);
 
   // VTK format
-  int readVTK(const std::string &name);
+  int readVTK(const std::string &name, bool bigEndian=false);
   int writeVTK(const std::string &name, bool binary=false,
                bool saveAll=false, double scalingFactor=1.0,
 	       bool bigEndian=false);
diff --git a/Geo/GModelIO_Mesh.cpp b/Geo/GModelIO_Mesh.cpp
index 537d953e42..89f153d954 100644
--- a/Geo/GModelIO_Mesh.cpp
+++ b/Geo/GModelIO_Mesh.cpp
@@ -2008,8 +2008,117 @@ int GModel::writeVTK(const std::string &name, bool binary, bool saveAll,
   return 1;
 }
 
-int GModel::readVTK(const std::string &name)
+int GModel::readVTK(const std::string &name, bool bigEndian)
 {
-  Msg::Error("VTK reader is not implemented yet");
+  FILE *fp = fopen(name.c_str(), "rb");
+  if(!fp){
+    Msg::Error("Unable to open file '%s'", name.c_str());
+    return 0;
+  }
+
+  char buffer[256], buffer2[256];
+
+  fgets(buffer, sizeof(buffer), fp); // version line
+  fgets(buffer, sizeof(buffer), fp); // title
+
+  fscanf(fp, "%s", buffer); // ASCII or BINARY
+  bool binary = false;
+  if(!strcmp(buffer, "BINARY")) binary = true;
+
+  if(fscanf(fp, "%s %s", &buffer, &buffer2) != 2) return 0;
+  if(strcmp(buffer, "DATASET") || strcmp(buffer2, "UNSTRUCTURED_GRID")){
+    Msg::Error("VTK reader can only read unstructured datasets");
+    return 0; 
+  }
+  
+  // read mesh vertices
+  int numVertices;
+  if(fscanf(fp, "%s %d %s\n", &buffer, &numVertices, buffer2) != 3) return 0;
+  if(strcmp(buffer, "POINTS") || strcmp(buffer2, "double")){
+    Msg::Error("VTK reader only accepts point data in double precision");
+    return 0;
+  }
+  if(!numVertices){
+    Msg::Warning("No points in dataset");
+    return 0;
+  }
+  Msg::Info("Reading %d points", numVertices);
+  std::vector<MVertex*> vertices(numVertices);
+  for(int i = 0 ; i < numVertices; i++){
+    double xyz[3];
+    if(binary){
+      if(fread(xyz, sizeof(double), 3, fp) != 3) return 0;
+      if(!bigEndian) SwapBytes((char*)xyz, sizeof(double), 3);
+    }
+    else{
+      if(fscanf(fp, "%lf %lf %lf", &xyz[0], &xyz[1], &xyz[2]) != 3) return 0;
+    }
+    vertices[i] = new MVertex(xyz[0], xyz[1], xyz[2]);
+  }
+  
+  // read mesh elements
+  int numElements, totalNumInt;
+  if(fscanf(fp, "%s %d %d\n", &buffer, &numElements, &totalNumInt) != 3) return 0;
+  if(strcmp(buffer, "CELLS") || !numElements){
+    Msg::Warning("No cells in dataset");
+    return 0;
+  }
+  Msg::Info("Reading %d cells", numElements);
+  std::vector<std::vector<MVertex*> > cells(numElements);
+  for(unsigned int i = 0; i < cells.size(); i++){
+    int numVerts, n[100];
+    if(binary){
+      if(fread(&numVerts, sizeof(int), 1, fp) != 1) return 0;
+      if(!bigEndian) SwapBytes((char*)&numVerts, sizeof(int), 1);
+      if(fread(n, sizeof(int), numVerts, fp) != numVerts) return 0;
+      if(!bigEndian) SwapBytes((char*)n, sizeof(int), numVerts);
+    }
+    else{
+      if(fscanf(fp, "%d", &numVerts) != 1) return 0;
+      for(int j = 0; j < numVerts; j++){
+	if(fscanf(fp, "%d", &n[j]) != 1) return 0;
+      }
+    }
+    for(int j = 0; j < numVerts; j++){
+      if(n[j] >= 0 && n[j] < vertices.size())
+	cells[i].push_back(vertices[n[j]]);
+      else
+	Msg::Error("Bad vertex index");
+    }
+  }
+  if(fscanf(fp, "%s %d\n", &buffer, &numElements) != 2) return 0;
+  if(strcmp(buffer, "CELL_TYPES") || numElements != cells.size()){
+    Msg::Error("No or invalid number of cells types");
+    return 0;
+  }
+  std::map<int, std::vector<MElement*> > elements[7];
+  for(unsigned int i = 0; i < cells.size(); i++){
+    int type;
+    if(binary){
+      if(fread(&type, sizeof(int), 1, fp) != 1) return 0;
+      if(!bigEndian) SwapBytes((char*)&type, sizeof(int), 1);
+    }
+    else{
+      if(fscanf(fp, "%d", &type) != 1) return 0;
+    }
+    switch(type){
+    case 3: elements[0][1].push_back(new MLine(cells[i])); break;
+    case 5: elements[1][1].push_back(new MTriangle(cells[i])); break;
+    case 8: elements[2][1].push_back(new MQuadrangle(cells[i])); break;
+    case 10: elements[3][1].push_back(new MTetrahedron(cells[i])); break;
+    case 12: elements[4][1].push_back(new MHexahedron(cells[i])); break;
+    case 13: elements[5][1].push_back(new MPrism(cells[i])); break;
+    case 14: elements[6][1].push_back(new MPyramid(cells[i])); break;
+    default: 
+      Msg::Error("Unknown type of cell %d", type);
+      break;
+    }
+  }  
+
+  for(int i = 0; i < (int)(sizeof(elements)/sizeof(elements[0])); i++) 
+    _storeElementsInEntities(elements[i]);
+  _associateEntityWithMeshVertices();
+  _storeVerticesInEntities(vertices);
+
   return 0;
 }
diff --git a/Geo/MElement.h b/Geo/MElement.h
index f942e958e2..76fcc86c5b 100644
--- a/Geo/MElement.h
+++ b/Geo/MElement.h
@@ -344,7 +344,7 @@ class MLine : public MElement {
   virtual void getFaceRep(int num, double *x, double *y, double *z, SVector3 *n){}
   virtual int getTypeForMSH() const { return MSH_LIN_2; }
   virtual int getTypeForUNV() const { return 21; } // linear beam
-  //virtual int getTypeForVTK() const { return 3; }
+  virtual int getTypeForVTK() const { return 3; }
   virtual const char *getStringForPOS() const { return "SL"; }
   virtual const char *getStringForBDF() const { return "CBAR"; }
   virtual void revert() 
-- 
GitLab