]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Patch by Scott Miller: Implement writing in VTU format.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 13 May 2010 03:01:41 +0000 (03:01 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 13 May 2010 03:01:41 +0000 (03:01 +0000)
git-svn-id: https://svn.dealii.org/trunk@21120 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/data_out_base.h
deal.II/base/source/data_out_base.cc

index 969300459fa1cfa5368287d4ac966827f22d32ad..7f4d162b0b09a868a8a4117b37dbe0df4ac0a715 100644 (file)
@@ -1246,9 +1246,6 @@ class DataOutBase
                                          */
        unsigned int memory_consumption () const;
     };
-
-
-
                                     /**
                                      * Flags controlling the details
                                      * of output in deal.II
@@ -1395,6 +1392,12 @@ class DataOutBase
                                            */
          vtk,
 
+                                          /**
+                                           * Output in @ref
+                                           * SoftwareVTU format.
+                                           */
+         vtu,
+         
                                           /**
                                            * Output in deal.II
                                            * intermediate format.
@@ -1674,11 +1677,9 @@ class DataOutBase
                           std::ostream                            &out);
 
 /**
- * Write the given list of patches to the output stream in @ref
- * SoftwareVTK format.
- *
- * This is the file format used by @ref SoftwareVTK, as described in
- * their manual, section 14.3.
+ * Write the given list of patches to the output stream in @ref SoftwareVTK
+ * format. The data is written in the traditional VTK format as opposed to the
+ * XML-based format that write_vtu() produces.
  *
  * The vector_data_ranges argument denotes ranges of components in the
  * output that are considered a vector, rather than simply a
@@ -1693,6 +1694,26 @@ class DataOutBase
                           const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
                           const VtkFlags                          &flags,
                           std::ostream                            &out);
+                          
+
+/**
+ * Write the given list of patches to the output stream in @ref SoftwareVTK
+ * format. The data is written in the XML-based VTK format as opposed to the
+ * traditional format that write_vtk() produces.
+ *
+ * The vector_data_ranges argument denotes ranges of components in the
+ * output that are considered a vector, rather than simply a
+ * collection of scalar fields. The VTK output format has special
+ * provisions that allow these components to be output by a single
+ * name rather than having to group several scalar fields into a
+ * vector later on in the visualization program.
+ */
+    template <int dim, int spacedim>
+    static void write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
+                          const std::vector<std::string>          &data_names,
+                          const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
+                          const VtkFlags                          &flags,
+                          std::ostream                            &out);                  
 
 /**
  * Write the given list of patches to the output stream in deal.II
@@ -1816,7 +1837,7 @@ class DataOutBase
                                      * <li> <tt>povray</tt>: <tt>.pov</tt>
                                      * <li> <tt>eps</tt>: <tt>.eps</tt>
                                      * <li> <tt>gmv</tt>: <tt>.gmv</tt>
-                                     * <li> <tt>vtk</tt>: <tt>.vtk</tt>.
+                                     * <li> <tt>vtk</tt>: <tt>.vtk</tt>
                                      * <li> <tt>deal_II_intermediate</tt>: <tt>.d2</tt>.
                                      * </ul>
                                      */
@@ -2105,6 +2126,7 @@ class DataOutInterface : private DataOutBase
     using DataOutBase::tecplot;
     using DataOutBase::tecplot_binary;
     using DataOutBase::vtk;
+    using DataOutBase::vtu;
     using DataOutBase::deal_II_intermediate;
     using DataOutBase::parse_output_format;
     using DataOutBase::get_output_format_names;
@@ -2197,6 +2219,14 @@ class DataOutInterface : private DataOutBase
                                      * DataOut::write_vtk.
                                      */
     void write_vtk (std::ostream &out) const;
+    
+                                    /**
+                                     * Obtain data through get_patches()
+                                     * and write it to <tt>out</tt>
+                                     * in Vtu (VTK's XML) format. See
+                                     * DataOut::write_vtu.
+                                     */
+    void write_vtu (std::ostream &out) const;
 
                                     /**
                                      * Obtain data through get_patches()
index 8e62c1eb2ba8acb25bcd617b916b43ce96abc717..5ceaa164e97577a8f2a7d63bd80ce4a4d5e7dbb0 100644 (file)
@@ -589,6 +589,69 @@ namespace
   };
 
 
+  class VtuStream
+  {
+    public:
+                                      /**
+                                       * Constructor, storing
+                                       * persistent values for
+                                       * later use.
+                                       */
+      VtuStream (std::ostream& stream,
+                const DataOutBase::VtkFlags flags);
+
+                                      /**
+                                       * Output operator for points.
+                                       */
+      template <int dim>
+      void write_point (const unsigned int index,
+                       const Point<dim>&);
+
+                                      /**
+                                       * Write dim-dimensional cell
+                                       * with first vertex at
+                                       * number start and further
+                                       * vertices offset by the
+                                       * specified values. Values
+                                       * not needed are ignored.
+                                       *
+                                       * The order of vertices for
+                                       * these cells in different
+                                       * dimensions is
+                                       * <ol>
+                                       * <li> [0,1]
+                                       * <li> []
+                                       * <li> []
+                                       * </ol>
+                                       */
+      template <int dim>
+      void write_cell(const unsigned int index,
+                     const unsigned int start,
+                     const unsigned int x_offset,
+                     const unsigned int y_offset,
+                     const unsigned int z_offset);
+
+                                      /**
+                                       * Forwarding of output stream
+                                       */
+      template <typename T>
+      std::ostream& operator<< (const T&);
+
+    private:
+                                      /**
+                                       * The ostream to use. Since
+                                       * the life span of these
+                                       * objects is small, we use a
+                                       * very simple storage
+                                       * technique.
+                                       */
+      std::ostream& stream;
+
+                                      /**
+                                       * The flags controlling the output
+                                       */
+      const DataOutBase::VtkFlags flags;
+  };
 
 
 //----------------------------------------------------------------------//
@@ -913,8 +976,52 @@ namespace
     stream << '\n';
   }
 
+  VtuStream::VtuStream(std::ostream& out, const DataOutBase::VtkFlags f)
+                 :
+                 stream(out), flags(f)
+  {}
 
 
+  template<int dim>
+  void
+  VtuStream::write_point (const unsigned int,
+                         const Point<dim>& p)
+  {
+                                    // write out coordinates
+    stream << p;
+                                    // fill with zeroes
+    for (unsigned int i=dim; i<3; ++i)
+      stream << " 0";
+    stream << '\n';
+  }
+
+
+  template<int dim>
+  void
+  VtuStream::write_cell(
+    unsigned int,
+    unsigned int start,
+    unsigned int d1,
+    unsigned int d2,
+    unsigned int d3)
+  {   
+    stream << start << '\t'
+          << start+d1;
+    if (dim>=2)
+      {
+       stream << '\t' << start+d2+d1
+              << '\t' << start+d2;
+       if (dim>=3)
+         {
+           stream << '\t' << start+d3
+                  << '\t' << start+d3+d1
+                  << '\t' << start+d3+d2+d1
+                  << '\t' << start+d3+d2;
+         }
+      }
+    stream << '\n';
+  }
+
 
   template <typename T>
   std::ostream&
@@ -3997,6 +4104,268 @@ DataOutBase::write_vtk (const std::vector<Patch<dim,spacedim> > &patches,
 }
 
 
+template <int dim, int spacedim>
+void
+DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
+                       const std::vector<std::string>          &data_names,
+                       const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
+                       const VtkFlags                          &flags,
+                       std::ostream                            &out)
+{
+  AssertThrow (out, ExcIO());
+
+  Assert (patches.size() > 0, ExcNoPatches());
+
+  VtuStream vtu_out(out, flags);
+
+  const unsigned int n_data_sets = data_names.size();
+                                  // check against # of data sets in
+                                  // first patch. checks against all
+                                  // other patches are made in
+                                  // write_gmv_reorder_data_vectors
+  Assert ((patches[0].data.n_rows() == n_data_sets && !patches[0].points_are_available) ||
+         (patches[0].data.n_rows() == n_data_sets+spacedim && patches[0].points_are_available),
+         ExcDimensionMismatch (patches[0].points_are_available
+                               ?
+                               (n_data_sets + spacedim)
+                               :
+                               n_data_sets,
+                               patches[0].data.n_rows()));
+
+                                  ///////////////////////
+                                  // preamble
+  if (true)
+    {
+      std::time_t  time1= std::time (0);
+      std::tm     *time = std::localtime(&time1);
+      out << "<?xml version=\"1.0\" ?> \n";
+      out << "<!-- \n";
+      out << "# vtk DataFile Version 3.0"
+         << '\n'
+         << "#This file was generated by the deal.II library on "
+         << time->tm_year+1900 << "/"
+         << time->tm_mon+1 << "/"
+         << time->tm_mday << " at "
+         << time->tm_hour << ":"
+         << std::setw(2) << time->tm_min << ":"
+         << std::setw(2) << time->tm_sec
+         << "\n --> \n";
+                                      // LittleEndian vs BigEndian should not
+                                      // matter; we are writing an ascii
+                                      // file.
+      out << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">";
+      out << '\n';
+      out << "<UnstructuredGrid>";
+      out << '\n';
+    }
+
+
+                                  // first count the number of cells
+                                  // and cells for later use
+  unsigned int n_nodes;
+  unsigned int n_cells;
+  compute_sizes<dim,spacedim>(patches, n_nodes, n_cells);
+                                  // in gmv format the vertex
+                                  // coordinates and the data have an
+                                  // order that is a bit unpleasant
+                                  // (first all x coordinates, then
+                                  // all y coordinate, ...; first all
+                                  // data of variable 1, then
+                                  // variable 2, etc), so we have to
+                                  // copy the data vectors a bit around
+                                  //
+                                  // note that we copy vectors when
+                                  // looping over the patches since we
+                                  // have to write them one variable
+                                  // at a time and don't want to use
+                                  // more than one loop
+                                  //
+                                  // this copying of data vectors can
+                                  // be done while we already output
+                                  // the vertices, so do this on a
+                                  // separate task and when wanting
+                                  // to write out the data, we wait
+                                  // for that task to finish
+  Table<2,double> data_vectors (n_data_sets, n_nodes);
+
+  void (*fun_ptr) (const std::vector<Patch<dim,spacedim> > &,
+                  Table<2,double> &)
+    = &DataOutBase::template write_gmv_reorder_data_vectors<dim,spacedim>;
+  Threads::Task<> reorder_task = Threads::new_task (fun_ptr, patches, data_vectors);
+
+                                  ///////////////////////////////
+                                  // first make up a list of used
+                                  // vertices along with their
+                                  // coordinates
+                                  //
+                                  // note that we have to print
+                                  // d=1..3 dimensions
+  out << "<Piece NumberOfPoints=\"" << n_nodes <<"\" NumberOfCells=\"" << n_cells << "\" > \n";
+  out << "  <Points> \n";
+  out << "    <DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\"ascii\"> \n";
+  write_nodes(patches, vtu_out);
+  out << "    </DataArray> \n";
+  out << "  </Points> \n \n";
+                                  /////////////////////////////////
+                                  // now for the cells
+  out << "  <Cells> \n";
+  out << "    <DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\"> \n";                              
+  write_cells(patches, vtu_out);
+  out << "    </DataArray> \n";
+  
+                                  // XML VTU format uses offsets; this is
+                                  // different than the VTK format, which
+                                  // puts the number of nodes per cell in
+                                  // front of the connectivity list.
+  out << "    <DataArray type=\"Int32\" Name=\"offsets\" format=\"ascii\"> \n";
+  for(unsigned int i=0; i<n_cells; ++i)
+    out << ' ' << ( (i+1)*GeometryInfo<dim>::vertices_per_cell);
+       
+  out << "    </DataArray> \n";
+  
+                                  // next output the types of the
+                                  // cells. since all cells are
+                                  // the same, this is simple
+  out << "    <DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\"> \n";
+  
+  for (unsigned int i=0; i<n_cells; ++i)
+    out << ' ' << vtk_cell_type[dim];
+   
+  out << "    </DataArray> \n"; 
+  out << "  </Cells> \n";
+  
+
+                                  ///////////////////////////////////////
+                                  // data output.
+
+                                  // now write the data vectors to
+                                  // @p{out} first make sure that all
+                                  // data is in place
+  reorder_task.join ();
+
+                                  // then write data.  the
+                                  // 'POINT_DATA' means: node data
+                                  // (as opposed to cell data, which
+                                  // we do not support explicitly
+                                  // here). all following data sets
+                                  // are point data
+  out << "  <PointData Scalars=\"scalars\"> \n";
+
+                                  // when writing, first write out
+                                  // all vector data, then handle the
+                                  // scalar data sets that have been
+                                  // left over
+  std::vector<bool> data_set_written (n_data_sets, false);
+  for (unsigned int n_th_vector=0; n_th_vector<vector_data_ranges.size(); ++n_th_vector)
+    {
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) >=
+                  std_cxx1x::get<0>(vector_data_ranges[n_th_vector]),
+                  ExcLowerRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]),
+                                 std_cxx1x::get<0>(vector_data_ranges[n_th_vector])));
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) < n_data_sets,
+                  ExcIndexRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]),
+                                 0, n_data_sets));
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) + 1
+                  - std_cxx1x::get<0>(vector_data_ranges[n_th_vector]) <= 3,
+                  ExcMessage ("Can't declare a vector with more than 3 components "
+                              "in VTK"));
+
+                                      // mark these components as already
+                                      // written:
+      for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
+          i<=std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
+          ++i)
+       data_set_written[i] = true;
+
+                                      // write the
+                                      // header. concatenate all the
+                                      // component names with double
+                                      // underscores unless a vector
+                                      // name has been specified
+      out << "    <DataArray type=\"Float32\" Name=\"";
+      
+      if (std_cxx1x::get<2>(vector_data_ranges[n_th_vector]) != "")
+       out << std_cxx1x::get<2>(vector_data_ranges[n_th_vector]);
+      else
+       {
+         for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
+              i<std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
+              ++i)
+           out << data_names[i] << "__";
+         out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])];
+       }
+
+      out << "\" NumberOfComponents=\"3\" format=\"ascii\"> \n";
+
+                                      // now write data. pad all
+                                      // vectors to have three
+                                      // components
+      for (unsigned int n=0; n<n_nodes; ++n)
+       {
+         switch (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) -
+                 std_cxx1x::get<0>(vector_data_ranges[n_th_vector]))
+           {
+             case 0:
+                   out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n) << " 0 0"
+                       << '\n';
+                   break;
+
+             case 1:
+                   out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]),   n) << ' '
+                       << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n) << " 0"
+                       << '\n';
+                   break;
+             case 2:
+                   out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]),   n) << ' '
+                       << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n) << ' '
+                       << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+2, n)
+                       << '\n';
+                   break;
+
+             default:
+                                                    // VTK doesn't
+                                                    // support
+                                                    // anything else
+                                                    // than vectors
+                                                    // with 1, 2, or
+                                                    // 3 components
+                   Assert (false, ExcInternalError());
+           }
+       }
+    
+      out << "    </DataArray> \n";
+    }
+
+                                  // now do the left over scalar data sets
+  for (unsigned int data_set=0; data_set<n_data_sets; ++data_set)
+    if (data_set_written[data_set] == false)
+      {
+       out << "    <DataArray type=\"Float32\" Name=\""
+           << data_names[data_set]
+           << "\" format=\"ascii\"> \n";
+           
+       std::copy (data_vectors[data_set].begin(),
+                  data_vectors[data_set].end(),
+                  std::ostream_iterator<double>(out, " "));
+       out << '\n';
+       out << "    </DataArray> \n";
+      }
+
+  out << "  </PointData> \n";
+
+                                  // Finish up writing a valid XML file
+  out << " </Piece> \n";
+  out << " </UnstructuredGrid> \n";
+  out << "</VTKFile> \n";
+                                  // make sure everything now gets to
+                                  // disk
+  out.flush ();
+
+                                  // assert the stream is still ok
+  AssertThrow (out, ExcIO());
+}
+
+
 
 template <int dim, int spacedim>
 void
@@ -4213,6 +4582,13 @@ void DataOutInterface<dim,spacedim>::write_vtk (std::ostream &out) const
                          vtk_flags, out);
 }
 
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::write_vtu (std::ostream &out) const
+{
+  DataOutBase::write_vtu (get_patches(), get_dataset_names(),
+                         get_vector_data_ranges(),
+                         vtk_flags, out);
+}
 
 
 template <int dim, int spacedim>
@@ -4275,6 +4651,10 @@ DataOutInterface<dim,spacedim>::write (std::ostream &out,
       case vtk:
            write_vtk (out);
            break;
+           
+         case vtu:
+           write_vtu (out);
+           break;  
 
       case deal_II_intermediate:
            write_deal_II_intermediate (out);

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.