]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Check in a first draft of the ability to write compressed VTU files. This
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 25 Jan 2011 04:37:54 +0000 (04:37 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 25 Jan 2011 04:37:54 +0000 (04:37 +0000)
isn't working right now, so it is disabled, but this way other people can
"throw an eye on it".

git-svn-id: https://svn.dealii.org/trunk@23259 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/source/base/data_out_base.cc

index 27e234d0bd59a57a98aa40aa3391aced8b801588..9d3dcac2589ee3f59112551ce7ed9a181ec9cecd 100644 (file)
 
 #include <sstream>
 
+#ifdef HAVE_LIBZ
+#  include <zlib.h>
+#endif
+
 DEAL_II_NAMESPACE_OPEN
 
 
@@ -4270,6 +4274,216 @@ DataOutBase::write_vtk (const std::vector<Patch<dim,spacedim> > &patches,
 }
 
 
+namespace
+{
+                                  // the functions in this namespace are
+                                  // taken from the libb64 project, see
+                                  // http://sourceforge.net/projects/libb64
+                                  //
+                                  // libb64 has been placed in the public
+                                  // domain
+  namespace base64
+  {
+    typedef enum
+    {
+         step_A, step_B, step_C
+    } base64_encodestep;
+
+    typedef struct
+    {
+       base64_encodestep step;
+       char result;
+    } base64_encodestate;
+
+    void base64_init_encodestate(base64_encodestate* state_in)
+    {
+      state_in->step = step_A;
+      state_in->result = 0;
+    }
+
+    inline
+    char base64_encode_value(char value_in)
+    {
+      static const char* encoding
+       = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/";
+      if (value_in > 63) return '=';
+      return encoding[(int)value_in];
+    }
+
+    int base64_encode_block(const char* plaintext_in,
+                           int length_in,
+                           char* code_out,
+                           base64_encodestate *state_in)
+    {
+      const char* plainchar = plaintext_in;
+      const char* const plaintextend = plaintext_in + length_in;
+      char* codechar = code_out;
+      char result;
+      char fragment;
+       
+      result = state_in->result;
+       
+      switch (state_in->step)
+       {
+         while (1)
+           {
+             case step_A:
+                   if (plainchar == plaintextend)
+                     {
+                       state_in->result = result;
+                       state_in->step = step_A;
+                       return codechar - code_out;
+                     }
+                   fragment = *plainchar++;
+                   result = (fragment & 0x0fc) >> 2;
+                   *codechar++ = base64_encode_value(result);
+                   result = (fragment & 0x003) << 4;
+             case step_B:
+                   if (plainchar == plaintextend)
+                     {
+                       state_in->result = result;
+                       state_in->step = step_B;
+                       return codechar - code_out;
+                     }
+                   fragment = *plainchar++;
+                   result |= (fragment & 0x0f0) >> 4;
+                   *codechar++ = base64_encode_value(result);
+                   result = (fragment & 0x00f) << 2;
+             case step_C:
+                   if (plainchar == plaintextend)
+                     {
+                       state_in->result = result;
+                       state_in->step = step_C;
+                       return codechar - code_out;
+                     }
+                   fragment = *plainchar++;
+                   result |= (fragment & 0x0c0) >> 6;
+                   *codechar++ = base64_encode_value(result);
+                   result  = (fragment & 0x03f) >> 0;
+                   *codechar++ = base64_encode_value(result);
+           }
+       }
+                                      /* control should not reach here */
+      return codechar - code_out;
+    }
+
+    int base64_encode_blockend(char* code_out, base64_encodestate* state_in)
+    {
+      char* codechar = code_out;
+       
+      switch (state_in->step)
+       {
+         case step_B:
+               *codechar++ = base64_encode_value(state_in->result);
+               *codechar++ = '=';
+               *codechar++ = '=';
+               break;
+         case step_C:
+               *codechar++ = base64_encode_value(state_in->result);
+               *codechar++ = '=';
+               break;
+         case step_A:
+               break;
+       }
+      *codechar++ = '\0';
+       
+      return codechar - code_out;
+    }
+  }
+  
+    
+                                  /**
+                                   * Do a base64 encoding of the given data.
+                                   *
+                                   * The function allocates memory as
+                                   * necessary and returns a pointer to
+                                   * it. The calling function must release
+                                   * this memory again.
+                                   */
+  char *
+  encode_block (const char *data,
+               const int   data_size)
+  {
+    base64::base64_encodestate state;
+    base64::base64_init_encodestate(&state);
+
+    char *encoded_data = new char[sizeof(uint32_t)+2*data_size+1];
+
+    uint32_t int_header = data_size;
+
+    const int encoded_length_header
+      = base64::base64_encode_block ((char*)&int_header, sizeof(int_header),
+                                    encoded_data, &state);
+
+    const int encoded_length_data
+      = base64::base64_encode_block (data, data_size,
+                                    encoded_data + encoded_length_header, &state);
+    base64::base64_encode_blockend (encoded_data + encoded_length_header + encoded_length_data,
+                                   &state);
+    
+    return encoded_data;
+  }
+  
+
+  
+#ifdef HAVE_LIBZ
+                                  /**
+                                   * Do a zlib compression followed by a
+                                   * base64 encoding of the given data.
+                                   *
+                                   * The function allocates memory as
+                                   * necessary and returns a pointer to
+                                   * it. The calling function must release
+                                   * this memory again.
+                                   */
+  template <typename T>
+  std::pair<char *, char *>
+  compress_and_encode_block (const std::vector<T> &data)
+  {
+    if (data.size() != 0)
+      {
+                                        // allocate a buffer for compressing
+                                        // data and do so
+       uLongf compressed_data_length
+         = compressBound (data.size() * sizeof(T));
+       char *compressed_data = new char[compressed_data_length];
+       int err = compress2 ((Bytef *) compressed_data,
+                            &compressed_data_length,
+                            (const Bytef *) &data[0],
+                            data.size() * sizeof(T),
+                            Z_BEST_COMPRESSION);
+       Assert (err == Z_OK, ExcInternalError());
+       
+                                        // now encode the compression header
+       const uint32_t compression_header[5]
+         = { 1,                                   /* number of blocks */
+             (uint32_t)(data.size() * sizeof(T)), /* size of block */
+             0, /* size of last block */
+             (uint32_t)compressed_data_length,
+             0};  /* list of compressed sizes of blocks */
+
+       char *encoded_header = encode_block ((char*)&compression_header[0],
+                                            5 * sizeof(compression_header[0]));
+       
+       
+                                        // now do the encoding in base64
+       char *encoded_data = encode_block (compressed_data,
+                                          compressed_data_length);
+       
+                                        // release the buffer for the
+                                        // compressed data and return the
+                                        // encoded data
+       delete[] compressed_data;
+       
+       return std::make_pair (encoded_header, encoded_data);
+      }
+    else
+      return std::pair<char *, char *>(0,0);
+  }
+#endif
+}
+
+
 template <int dim, int spacedim>
 void
 DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
@@ -4334,10 +4548,17 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
          << 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 << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\"";
+#ifdef HAVE_LIBZ && 0
+      out << " compressor=\"vtkZLibDataCompressor\"";
+#endif
+#ifdef DEAL_II_WORDS_BIGENDIAN
+      out << " byte_order=\"BigEndian\"";
+#else
+      out << " byte_order=\"LittleEndian\"";
+#endif
+      out << ">";
       out << '\n';
       out << "<UnstructuredGrid>";
       out << '\n';
@@ -4375,7 +4596,8 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
   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);
+  Threads::Task<> reorder_task = Threads::new_task (fun_ptr, patches,
+                                                   data_vectors);
 
                                   ///////////////////////////////
                                   // first make up a list of used
@@ -4385,16 +4607,75 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
                                   // note that according to the standard, we
                                   // have to print d=1..3 dimensions, even if
                                   // we are in reality in 2d, for example
-  out << "<Piece NumberOfPoints=\"" << n_nodes <<"\" NumberOfCells=\"" << n_cells << "\" >\n";
+  out << "<Piece NumberOfPoints=\"" << n_nodes
+      <<"\" NumberOfCells=\"" << n_cells << "\" >\n";
   out << "  <Points>\n";
+#ifndef HAVE_LIBZ && 0
   out << "    <DataArray type=\"Float64\" NumberOfComponents=\"3\" format=\"ascii\">\n";
   write_nodes(patches, vtu_out);
   out << "    </DataArray>\n";
+#else
+  out << "    <DataArray type=\"Float64\" NumberOfComponents=\"3\" format=\"binary\">\n";
+  {
+                                    // collect vertices into one big array,
+                                    // then compress and encode it. we don't
+                                    // know up front how many vertices we
+                                    // need, so only allocate a typical guess
+    std::vector<double> vertices;
+    vertices.reserve (patches.size() * 3 * GeometryInfo<spacedim>::vertices_per_cell);
+    for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+        patch!=patches.end(); ++patch)
+      {
+       const unsigned int n_subdivisions = patch->n_subdivisions;
+       const unsigned int n = n_subdivisions+1;
+                                        // Length of loops in all
+                                        // dimensions. If a dimension
+                                        // is not used, a loop of
+                                        // length one will do the job.
+       const unsigned int n1 = (dim>0) ? n : 1;
+       const unsigned int n2 = (dim>1) ? n : 1;
+       const unsigned int n3 = (dim>2) ? n : 1;
+
+       for (unsigned int i3=0; i3<n3; ++i3)
+         for (unsigned int i2=0; i2<n2; ++i2)
+           for (unsigned int i1=0; i1<n1; ++i1)
+             {
+               Point<spacedim> node;
+               compute_node(node, &*patch,
+                            i1,
+                            i2,
+                            i3,
+                            n_subdivisions);
+               for (unsigned int d=0; d<spacedim; ++d)
+                 vertices.push_back (node[d]);
+               for (unsigned int d=spacedim; d<3; ++d)
+                 vertices.push_back (0);
+             }
+      }
+
+/*    
+    const char *
+      encoded_vertices
+      = encode_block ((char*)&vertices[0],
+                     vertices.size() * sizeof(vertices[0]));
+    out << encoded_vertices << std::endl;
+    delete[] encoded_vertices;
+*/
+    const std::pair<char *, char *>
+      data
+      = compress_and_encode_block (vertices);
+    out << data.first << data.second;
+    
+    delete[] data.first;
+    delete[] data.second;
+  }
+  out << "    </DataArray>\n";
+#endif
   out << "  </Points>\n\n";
                                   /////////////////////////////////
                                   // now for the cells
   out << "  <Cells>\n";
-  out << "    <DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\">\n";                               
+  out << "    <DataArray type=\"Int32\" Name=\"connectivity\" format=\"ascii\">\n";
   write_cells(patches, vtu_out);
   out << "    </DataArray>\n";
   

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.