]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Finish reasonable support for the EPS format (except documentation).
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Jul 1999 11:49:42 +0000 (11:49 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Jul 1999 11:49:42 +0000 (11:49 +0000)
git-svn-id: https://svn.dealii.org/trunk@1546 0785d39b-7218-0410-832d-ea1e28bc413d

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

index b1912649948798c3dd0ee52d60e4c1ed0c264dd7..89d5eb55a6698e8cb36ba7e2bd0c675aded7cff9 100644 (file)
  * Given the lines as described above, a cut through this data in Gnuplot
  * can then be achieved like this:
  * \begin{verbatim}
- * set data style lines
- * splot [:][:][0:] "T" using 1:2:($3==.5 ? $4 : -1)
+ *   set data style lines
+ *   splot [:][:][0:] "T" using 1:2:($3==.5 ? $4 : -1)
  * \end{verbatim}
  * This command plots data in x- and y-direction unbounded, but in z-direction
  * only those data points which are above the x-y-plane (we assume here a
  *
  * To be filled in.
  * precision=5; viewpoint=gnuplot default; no border
- * shade or not; grid or not; grid shaded; data vector
+ * shade or not; grid or not; grid shaded; data vector; memory; color=one per cell
  *
  * \subsection{GMV format}
  *
  * details of the viewpoint, light source, etc.
  *
  *
- * @author Wolfgang Bangerth 1999
+ * @author Wolfgang Bangerth 1999; EPS output based on an earlier implementation by Stefan Nauber for the old DataOut class
  */
 class DataOutBase 
 {
@@ -332,6 +332,7 @@ class DataOutBase
        int dummy;
     };
 
+    
                                     /**
                                      * Flags describing the details of
                                      * output in encapsulated postscript
@@ -345,12 +346,20 @@ class DataOutBase
                                          * for generating the height
                                          * information. By default, the
                                          * first data vector is taken,
-                                         * i.e. #height_value==0#, if
+                                         * i.e. #height_vector==0#, if
                                          * there is any data vector. If there
                                          * is no data vector, no height
                                          * information is generated.
                                          */
        unsigned int height_vector;
+
+                                        /**
+                                         * Number of the vector which is
+                                         * to be taken to colorize cells.
+                                         * The same applies as for
+                                         * #height_vector#.
+                                         */
+       unsigned int color_vector;
        
                                         /**
                                          * Enum denoting the possibilities
@@ -439,17 +448,105 @@ class DataOutBase
                                          */
        bool   draw_mesh;
 
+                                        /**
+                                         * Flag whether to fill the regions
+                                         * between the lines bounding the cells
+                                         * or not. If not, no hidden line removal
+                                         * is performed, which in this crude
+                                         * implementation is done through
+                                         * writing the cells in a back-to-front
+                                         * order, thereby hiding the cells in
+                                         * the background by cells in the
+                                         * foreground.
+                                         *
+                                         * If this flag is #false# and #draw_mesh#
+                                         * is #false# as well, nothing will be
+                                         * printed.
+                                         *
+                                         * If this flag is #true#, then the cells
+                                         * will be drawn either colored by one
+                                         * of the data sets (if #shade_cells# is
+                                         * #true#), or pure white (if
+                                         * #shade_cells# is false or if there are
+                                         * no data sets).
+                                         *
+                                         * Default is #true#.
+                                         */
+       bool   draw_cells;
+
+                                        /**
+                                         * Flag to determine whether the cells
+                                         * shall be colorized by one the data
+                                         * set denoted by #color_vector#, or
+                                         * simply be painted in white. This
+                                         * flag only makes sense if
+                                         * #draw_cells==true#. Colorization is
+                                         * done through the #color_function#.
+                                         *
+                                         * Default is #true#.
+                                         */
+       bool   shade_cells;
+
+                                        /**
+                                         * Structure keeping the three color
+                                         * values in the RGB system.
+                                         */
+       struct RgbValues { float red, green, blue; };
+
+                                        /**
+                                         * Definition of a function pointer
+                                         * type taking a value and returning
+                                         * a triple of color values in RGB
+                                         * values.
+                                         *
+                                         * Besides the actual value by which
+                                         * the color is to be computed, min
+                                         * and max values of the data to
+                                         * be colorized are given as well.
+                                         */
+       typedef RgbValues (*ColorFunction) (const double value,
+                                           const double min_value,
+                                           const double max_value);
+
+                                        /**
+                                         * This is a pointer to the function
+                                         * which is used to colorize the cells.
+                                         * By default, it points to the
+                                         * static function #default_color_function#
+                                         * which is a member of this class.
+                                         */
+       ColorFunction color_function;
+
+
+                                        /**
+                                         * Default colorization function. This
+                                         * one does what one usually wants:
+                                         * it shifts colors from black (lowest
+                                         * value) through blue, green and red
+                                         * to white (highest value).
+                                         *
+                                         * This function was originally written
+                                         * by Stefan Nauber.
+                                         */
+       static RgbValues default_color_function (const double value,
+                                                const double min_value,
+                                                const double max_value);
+       
                                         /**
                                          * Constructor.
                                          */
-       EpsFlags (const unsigned int height_vector = 0,
-                 const SizeType     size_type    = width,
-                 const unsigned int size         = 300,
-                 const double       line_width   = 0.5,
-                 const double       azimut_angle = 60,
-                 const double       turn_angle   = 30,
-                 const double       z_scaling    = 1.0,
-                 const bool         draw_mesh    = true);
+       EpsFlags (const unsigned int  height_vector = 0,
+                 const unsigned int  color_vector  = 0,
+                 const SizeType      size_type     = width,
+                 const unsigned int  size          = 300,
+                 const double        line_width    = 0.5,
+                 const double        azimut_angle  = 60,
+                 const double        turn_angle    = 30,
+                 const double        z_scaling     = 1.0,
+                 const bool          draw_mesh     = true,
+                 const bool          draw_cells    = true,
+                 const bool          shade_cells   = true,
+                 const ColorFunction color_function= &default_color_function);
     };
 
                                     /**
@@ -555,10 +652,10 @@ class DataOutBase
                                     /**
                                      * Exception
                                      */
-    DeclException2 (ExcInvalidHeightVectorNumber,
+    DeclException2 (ExcInvalidVectorNumber,
                    int, int,
                    << "The number " << arg1 << " of the vector to be used for "
-                   << "height information is invalid, since there are only "
+                   << "this information is invalid, since there are only "
                    << arg2 << " data sets.");
                                     /**
                                      * Exception
@@ -584,12 +681,13 @@ class DataOutBase
        Point<2> vertices[4];
        
                                         /**
-                                         * Color values.
+                                         * Data value from which the actual
+                                         * colors will be computed by
+                                         * the colorization function stated
+                                         * in the #EpsFlags# class.
                                          */
-       float red;
-       float green;
-       float blue;
-
+       float color_value;
+       
                                         /**
                                          * Depth into the picture, which
                                          * is defined as the distance from
index 37e373744534614b5625daddfa8ffdd3d1c88d3d..8d55b521d83ecae888b84aed3b7975f55d60cf70 100644 (file)
@@ -23,26 +23,121 @@ DataOutBase::UcdFlags::UcdFlags (const bool write_preamble) :
 {};
 
 
-DataOutBase::EpsFlags::EpsFlags (const unsigned int height_vector,
-                                const SizeType     size_type,
-                                const unsigned int size,
-                                const double       line_width,
-                                const double       azimut_angle,
-                                const double       turn_angle,
-                                const double       z_scaling,
-                                const bool         draw_mesh) :
+DataOutBase::EpsFlags::EpsFlags (const unsigned int  height_vector,
+                                const unsigned int  color_vector,
+                                const SizeType      size_type,
+                                const unsigned int  size,
+                                const double        line_width,
+                                const double        azimut_angle,
+                                const double        turn_angle,
+                                const double        z_scaling,
+                                const bool          draw_mesh,
+                                const bool          draw_cells,
+                                const bool          shade_cells,
+                                const ColorFunction color_function) :
                height_vector(height_vector),
+               color_vector(color_vector),
                size_type(size_type),
                size(size),
                line_width(line_width),
                azimut_angle(azimut_angle),
                turn_angle(turn_angle),
                z_scaling(z_scaling),
-               draw_mesh(draw_mesh)
+               draw_mesh(draw_mesh),
+               draw_cells(draw_cells),
+               shade_cells(shade_cells),
+               color_function(color_function)
 {};
 
 
 
+DataOutBase::EpsFlags::RgbValues
+DataOutBase::EpsFlags::default_color_function (const double x,
+                                              const double xmin,
+                                              const double xmax)
+{
+  RgbValues rgb_values;
+  
+// A difficult color scale:
+//     xmin          = black  (1)
+// 3/4*xmin+1/4*xmax = blue   (2)
+// 1/2*xmin+1/2*xmax = green  (3)
+// 1/4*xmin+3/4*xmax = red    (4)
+//              xmax = white  (5)
+// Makes the following color functions:
+//
+// red      green    blue
+//       __
+//      /      /\  /  /\    /
+// ____/    __/  \/  /  \__/
+
+//     { 0                                (1) - (3)
+// r = { ( 4*x-2*xmin+2*xmax)/(xmax-xmin) (3) - (4)
+//     { 1                                (4) - (5)
+//
+//     { 0                                (1) - (2)
+// g = { ( 4*x-3*xmin-  xmax)/(xmax-xmin) (2) - (3)
+//     { (-4*x+  xmin+3*xmax)/(xmax-xmin) (3) - (4)
+//     { ( 4*x-  xmin-3*xmax)/(xmax-xmin) (4) - (5)
+//
+//     { ( 4*x-4*xmin       )/(xmax-xmin) (1) - (2)
+// b = { (-4*x+2*xmin+2*xmax)/(xmax-xmin) (2) - (3)
+//     { 0                                (3) - (4)
+//     { ( 4*x-  xmin-3*xmax)/(xmax-xmin) (4) - (5)
+
+  double sum   =   xmax+  xmin;
+  double sum13 =   xmin+3*xmax;
+  double sum22 = 2*xmin+2*xmax;
+  double sum31 = 3*xmin+  xmax;
+  double dif = xmax-xmin;
+  double rezdif = 1.0/dif;
+
+  int where;
+
+  if (x<(sum31)/4)
+    where = 0;
+  else if (x<(sum22)/4)
+    where = 1;
+  else if (x<(sum13)/4)
+    where = 2;
+  else
+    where = 3;
+
+  if (dif!=0)
+    {
+      switch (where)
+       {
+         case 0:
+               rgb_values.red   = 0;
+               rgb_values.green = 0;
+               rgb_values.blue  = (x-xmin)*4.*rezdif;
+               break;
+         case 1:
+               rgb_values.red   = 0;
+               rgb_values.green = (4*x-3*xmin-xmax)*rezdif;
+               rgb_values.blue  = (sum22-4.*x)*rezdif;
+               break;
+         case 2:
+               rgb_values.red   = (4*x-2*sum)*rezdif;
+               rgb_values.green = (xmin+3*xmax-4*x)*rezdif;
+               rgb_values.blue  = 0;
+               break;
+         case 3:
+               rgb_values.red   = 1;
+               rgb_values.green = (4*x-xmin-3*xmax)*rezdif;
+               rgb_values.blue  = (4.*x-sum13)*rezdif;
+         default:
+               break;
+       };
+    }
+  else // White 
+    rgb_values.red = rgb_values.green = rgb_values.blue = 1;
+
+  return rgb_values;
+};
+
+
+
 bool DataOutBase::EpsCell2d::operator < (const EpsCell2d &e) const
 {
                                   // note the "wrong" order in
@@ -749,6 +844,12 @@ void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
                                         // for sorting
        multiset<EpsCell2d> cells;
 
+                                        // two variables in which we
+                                        // will store the minimum and
+                                        // maximum values of the field
+                                        // to be used for colorization
+       double min_color_value, max_color_value;
+       
                                         // compute the cells for output and
                                         // enter them into the set above
                                         // note that since dim==2, we
@@ -791,8 +892,8 @@ void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
 
                  Assert ((flags.height_vector < patch->data.m()) ||
                          patch->data.m() == 0,
-                         ExcInvalidHeightVectorNumber (flags.height_vector,
-                                                       patch->data.m()));
+                         ExcInvalidVectorNumber (flags.height_vector,
+                                                 patch->data.m()));
                  const double heights[4]
                    = { patch->data.m() != 0 ?
                        patch->data(flags.height_vector,i*(n_subdivisions+1) + j)       * flags.z_scaling : 0,
@@ -855,7 +956,7 @@ void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
 
 //       ( cz -sz 0 )( 1 0    0 )(x)   ( cz*x-sz*(cx*y-sx*z)+0*(sx*y+cx*z) )
 // Dxz = ( sz  cz 0 )( 0 cx -sx )(y) = ( sz*x+cz*(cx*y-sx*z)+0*(sx*y+cx*z) )
-//     (  0   0 1 )( 0 sx  cx )(z)   (  0*x+   *(cx*y-sx*z)+1*(sx*y+cx*z) )
+//      (  0   0 1 )( 0 sx  cx )(z)   (  0*x+  *(cx*y-sx*z)+1*(sx*y+cx*z) )
                    };
 
                                                   // compute coordinates of
@@ -871,6 +972,45 @@ void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
                                   -sx*cz*center_point(1)
                                   +cx*center_height;
 
+                 if (flags.draw_cells && flags.shade_cells)
+                   {
+                     Assert ((flags.color_vector < patch->data.m()) ||
+                             patch->data.m() == 0,
+                             ExcInvalidVectorNumber (flags.color_vector,
+                                                     patch->data.m()));
+                     const double color_values[4]
+                       = { patch->data.m() != 0 ?
+                           patch->data(flags.color_vector,i*(n_subdivisions+1) + j)       : 1,
+                       
+                           patch->data.m() != 0 ?
+                           patch->data(flags.color_vector,(i+1)*(n_subdivisions+1) + j)   : 1,
+                           
+                           patch->data.m() != 0 ?
+                           patch->data(flags.color_vector,(i+1)*(n_subdivisions+1) + j+1) : 1,
+                           
+                           patch->data.m() != 0 ?
+                           patch->data(flags.color_vector,i*(n_subdivisions+1) + j+1)     : 1};
+
+                                                      // set color value to average of the value
+                                                      // at the vertices
+                     eps_cell.color_value = (color_values[0] +
+                                             color_values[1] +
+                                             color_values[2] +
+                                             color_values[3]) / 4;
+
+                                                      // update bounds of color
+                                                      // field
+                     if (patch == patches.begin())
+                       min_color_value = max_color_value = eps_cell.color_value;
+                     else
+                       {
+                         min_color_value = (min_color_value < eps_cell.color_value ?
+                                            min_color_value : eps_cell.color_value);
+                         max_color_value = (max_color_value > eps_cell.color_value ?
+                                            max_color_value : eps_cell.color_value);
+                       };
+                   };
+                 
                                                   // finally add this cell
                  cells.insert (eps_cell);
                };
@@ -969,16 +1109,31 @@ void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
        for (multiset<EpsCell2d>::const_iterator cell=cells.begin();
             cell!=cells.end(); ++cell)
          {
-//         out << cell->red << ' ' << cell->green << ' ' << cell->blue << " s "
-           out << 1 << ' ' << 1 << ' ' << 1 << " s "
-               << (cell->vertices[0]-offset) * scale << " m "
-               << (cell->vertices[1]-offset) * scale << " l "
-               << (cell->vertices[2]-offset) * scale << " l "
-               << (cell->vertices[3]-offset) * scale << " lf"
-               << endl;
-
+           if (flags.draw_cells) 
+             {
+               if (flags.shade_cells)
+                 {
+                   const EpsFlags::RgbValues rgb_values
+                     = (*flags.color_function) (cell->color_value,
+                                                min_color_value,
+                                                max_color_value);
+                   
+                   out << rgb_values.red   << ' '
+                       << rgb_values.green << ' '
+                       << rgb_values.blue  << " s ";
+                 }
+               else
+                 out << "1 1 1 s ";
+
+               out << (cell->vertices[0]-offset) * scale << " m "
+                   << (cell->vertices[1]-offset) * scale << " l "
+                   << (cell->vertices[2]-offset) * scale << " l "
+                   << (cell->vertices[3]-offset) * scale << " lf"
+                   << endl;
+             };
+           
            if (flags.draw_mesh)
-             out << "0 sg " 
+             out << "0 sg "      // draw lines in black
                  << (cell->vertices[0]-offset) * scale << " m "
                  << (cell->vertices[1]-offset) * scale << " l "
                  << (cell->vertices[2]-offset) * scale << " l "

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.