Poster of Linux kernelThe best gift for a Linux geek
PostscriptIO

PostscriptIO

Section: C Library Functions (3) Updated: Thu Apr 7 2011
Local index Up
 

NAME

PostscriptIO -  

SYNOPSIS


#include <postscript_io.h>

Inherits MeshOutput< MeshBase >.  

Public Member Functions


PostscriptIO (const MeshBase &mesh)

virtual ~PostscriptIO ()

virtual void write (const std::string &)

void plot_quadratic_elem (const Elem *elem)

void plot_linear_elem (const Elem *elem)

virtual void write_equation_systems (const std::string &, const EquationSystems &)

virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 

Public Attributes


Real shade_value

Real line_width
 

Protected Member Functions


const MeshBase & mesh () const
 

Private Member Functions


void _compute_edge_bezier_coeffs (const Elem *elem)
 

Private Attributes


std::vector< Point > _bezier_coeffs

Point _offset

Real _scale

Point _current_point

std::ostringstream _cell_string

std::ofstream _out
 

Static Private Attributes


static const float _bezier_transform [3][3]
 

Detailed Description

This class implements writing 2D meshes in Postscript. It borrows several ideas from, and is a more simple-minded version of, the DataOutBase::write_eps() function from Deal II. Only output is supported here, and only the Mesh (none of the data) is written. The main use I imagined for this class is creating nice Mesh images for publications, since I didn't find/don't know of a free visualization program which would do this.

Author:

John W. Peterson, 2008

Definition at line 49 of file postscript_io.h.  

Constructor & Destructor Documentation

 

PostscriptIO::PostscriptIO (const MeshBase &mesh)Constructor.

Definition at line 38 of file postscript_io.C.

References _bezier_coeffs.

  : MeshOutput<MeshBase> (mesh),
    shade_value(0.0),
    line_width(0.5),
    //_M(3,3),
    _offset(0., 0.),
    _scale(1.0),
    _current_point(0., 0.)
{
  // This code is still undergoing some development.
  libmesh_experimental();

  // Entries of transformation matrix from physical to Bezier coords.
  // _M(0,0) = -1./6.;    _M(0,1) = 1./6.;    _M(0,2) = 1.; 
  // _M(1,0) = -1./6.;    _M(1,1) = 0.5  ;    _M(1,2) = 1./6.; 
  // _M(2,0) = 0.    ;    _M(2,1) = 1.   ;    _M(2,2) = 0.;

  // Make sure there is enough room to store Bezier coefficients.
  _bezier_coeffs.resize(3);
}
 

PostscriptIO::~PostscriptIO () [virtual]Destructor.

Definition at line 61 of file postscript_io.C.

{
}
 

Member Function Documentation

 

void PostscriptIO::_compute_edge_bezier_coeffs (const Elem *elem) [private]Given a quadratic edge Elem which lies in the x-y plane, computes the Bezier coefficients. These may be passed to the Postscript routine 'curveto'.

Definition at line 263 of file postscript_io.C.

References _bezier_coeffs, _bezier_transform, _offset, _scale, libMeshEnums::EDGE3, Elem::point(), and Elem::type().

Referenced by plot_quadratic_elem().

{
  // I only know how to do this for an Edge3!
  libmesh_assert (elem->type() == EDGE3);

  // Get x-coordinates into an array, transform them,
  // and repeat for y.
  float
    phys_coords[3] = {0., 0., 0.},
    bez_coords[3]  = {0., 0., 0.};
  
  for (unsigned int i=0; i<2; ++i)
    {
      // Initialize vectors.  Physical coordinates are initialized
      // by their postscript-scaled values. 
      for (unsigned int j=0; j<3; ++j)
        {
          phys_coords[j] = (elem->point(j)(i) - _offset(i)) * _scale;
          bez_coords[j] = 0.; // zero out result vector
        }

      // Multiply matrix times vector
      for (unsigned int j=0; j<3; ++j)
        for (unsigned int k=0; k<3; ++k)
          bez_coords[j] += _bezier_transform[j][k]*phys_coords[k];

      // Store result in _bezier_coeffs
      for (unsigned int j=0; j<3; ++j)
        _bezier_coeffs[j](i) = phys_coords[j];
    }
}
 

const MeshBase & MeshOutput< MeshBase >::mesh () const [protected, inherited]Returns the object as a read-only reference.

Referenced by write(), FroIO::write(), MEDITIO::write_ascii(), EnsightIO::write_geometry_ascii(), EnsightIO::write_scalar_ascii(), GnuPlotIO::write_solution(), DivaIO::write_stream(), and EnsightIO::write_vector_ascii().  

void PostscriptIO::plot_linear_elem (const Elem *elem)Draws an element with straight lines

Definition at line 189 of file postscript_io.C.

References _cell_string, _current_point, _offset, _out, _scale, Elem::n_vertices(), Elem::point(), and shade_value.

Referenced by write().

{
  // Clear the string contents.  Yes, this really is how you do that...
  _cell_string.str('');
          
  // The general strategy is:
  // 1.) Use m  := {moveto} to go to vertex 0.
  // 2.) Use l  := {lineto} commands to draw lines to vertex 1, 2, ... N-1.
  // 3a.) Use lx := {lineto closepath stroke} command at  vertex N to draw the last line.
  // 3b.)     lf := {lineto closepath fill} command to shade the cell just drawn 
  // All of our 2D elements' vertices are numbered in counterclockwise order,
  // so we can just draw them in the same order.

  // 1.)
  _current_point = (elem->point(0) - _offset) * _scale;
  _cell_string << _current_point(0) << ' ' << _current_point(1) << ' '; // write x y 
  _cell_string << 'm ';

  // 2.)
  const unsigned int nv=elem->n_vertices();
  for (unsigned int v=1; v<nv-1; ++v)
    {
      _current_point = (elem->point(v) - _offset) * _scale;
      _cell_string << _current_point(0) << ' ' << _current_point(1) << ' '; // write x y 
      _cell_string << 'l ';
    }
          
  // 3.)
  _current_point = (elem->point(nv-1) - _offset) * _scale;
  _cell_string << _current_point(0) << ' ' << _current_point(1) << ' '; // write x y 

  // We draw the shaded (interior) parts first, if applicable.
  if (shade_value > 0.0)
    _out << shade_value << ' sg ' << _cell_string.str() << 'lf;
          
  // Draw the black lines (I guess we will always do this)
  _out << '0 sg ' << _cell_string.str() << 'lx;
}
 

void PostscriptIO::plot_quadratic_elem (const Elem *elem)Draws an element with Bezier curves

Definition at line 231 of file postscript_io.C.

References _bezier_coeffs, _compute_edge_bezier_coeffs(), _current_point, _offset, _out, _scale, Elem::build_side(), libMeshEnums::EDGE3, AutoPtr< Tp >::get(), and Elem::n_sides().

{
  for (unsigned int ns=0; ns<elem->n_sides(); ++ns)
    {
      // Build the quadratic side
      AutoPtr<Elem> side = elem->build_side(ns);

      // Be sure it's quadratic (Edge2).  Eventually we could
      // handle cubic elements as well...
      libmesh_assert( side->type() == EDGE3 );

      _out << '0 sg ';

      // Move to the first point on this side.
      _current_point = (side->point(0) - _offset) * _scale;
      _out << _current_point(0) << ' ' << _current_point(1) << ' '; // write x y 
      _out << 'm ';
      
      // Compute _bezier_coeffs for this edge.  This fills up
      // the _bezier_coeffs vector.
      this->_compute_edge_bezier_coeffs(side.get());
      
      // Print curveto path to file
      for (unsigned int i=0; i<_bezier_coeffs.size(); ++i)
        _out << _bezier_coeffs[i](0) << ' ' << _bezier_coeffs[i](1) << ' ';
      _out << ' cs;
    }
}
 

void PostscriptIO::write (const std::string &fname) [virtual]This method implements writing a mesh to a specified file.

Implements MeshOutput< MeshBase >.

Definition at line 67 of file postscript_io.C.

References _offset, _out, _scale, MeshBase::active_elements_begin(), MeshBase::active_elements_end(), MeshTools::bounding_box(), line_width, MeshOutput< MeshBase >::mesh(), MeshBase::mesh_dimension(), plot_linear_elem(), and libMesh::processor_id().

{
  if (libMesh::processor_id() == 0)
    {
      // Get a constant reference to the mesh.
      const MeshBase& mesh = MeshOutput<MeshBase>::mesh();

      // Only works in 2D
      libmesh_assert(mesh.mesh_dimension() == 2);

      // Create output file stream.
      // _out is now a private member of the class.
      _out.open(fname.c_str());
      
      // Make sure it opened correctly
      if (!_out.good())
        libmesh_file_error(fname.c_str());

      // The mesh bounding box gives us info about what the
      // Postscript bounding box should be.
      MeshTools::BoundingBox bbox = MeshTools::bounding_box(mesh);

      // Add a little extra padding to the 'true' bounding box so
      // that we can still see the boundary
      const Real percent_padding = 0.01;
      const Real dx=bbox.second(0)-bbox.first(0); libmesh_assert(dx > 0.0);
      const Real dy=bbox.second(1)-bbox.first(1); libmesh_assert(dy > 0.0);
      
      const Real x_min = bbox.first(0)  - percent_padding*dx;
      const Real y_min = bbox.first(1)  - percent_padding*dy;
      const Real x_max = bbox.second(0) + percent_padding*dx;
      const Real y_max = bbox.second(1) + percent_padding*dy;

      // Width of the output as given in postscript units.
      // This usually is given by the strange unit 1/72 inch. 
      // A width of 300 represents a size of roughly 10 cm.
      const Real width = 300;
      _scale = width / (x_max-x_min);
      _offset(0) = x_min;
      _offset(1) = y_min;
      
      // Header writing stuff stolen from Deal.II
      std::time_t  time1= std::time (0);
      std::tm     *time = std::localtime(&time1);
      _out << '%!PS-Adobe-2.0 EPSF-1.2' << '
        //<< '%!PS-Adobe-1.0' << ' // Lars' PS version
          << '%%Filename: ' << fname << '
          << '%%Title: LibMesh Output' << '
          << '%%Creator: LibMesh: A C++ finite element library' << '
          << '%%Creation Date: '
          << time->tm_year+1900 << '/'
          << time->tm_mon+1 << '/'
          << time->tm_mday << ' - '
          << time->tm_hour << ':'
          << std::setw(2) << time->tm_min << ':'
          << std::setw(2) << time->tm_sec << '
          << '%%BoundingBox: '
        // lower left corner
          << '0 0 '
        // upper right corner
          << static_cast<unsigned int>( rint((x_max-x_min) * _scale ))
          << ' '
          << static_cast<unsigned int>( rint((y_max-y_min) * _scale ))
          << ';

      // define some abbreviations to keep
      // the output small:
      // m=move turtle to
      // l=define a line
      // s=set rgb color
      // sg=set gray value
      // lx=close the line and plot the line
      // lf=close the line and fill the interior
      _out << '/m {moveto} bind def'      << '
          << '/l {lineto} bind def'      << '
          << '/s {setrgbcolor} bind def' << '
          << '/sg {setgray} bind def'    << '
          << '/cs {curveto stroke} bind def' << '
          << '/lx {lineto closepath stroke} bind def' << '
          << '/lf {lineto closepath fill} bind def'   << ';

      _out << '%%EndProlog' << ';
      //          << ';
      
      // Set line width in the postscript file.  
      _out << line_width << ' setlinewidth' << ';

      // Set line cap and join options
      _out << '1 setlinecap' << ';
      _out << '1 setlinejoin' << ';
  
      // allow only five digits for output (instead of the default
      // six); this should suffice even for fine grids, but reduces
      // the file size significantly
      _out << std::setprecision (5);

      // Loop over the active elements, draw lines for the edges.  We
      // draw even quadratic elements with straight sides, i.e. a straight
      // line sits between each pair of vertices.  Also we draw every edge
      // for an element regardless of the fact that it may overlap with
      // another.  This would probably be a useful optimization...
      MeshBase::const_element_iterator       el     = mesh.active_elements_begin();
      const MeshBase::const_element_iterator end_el = mesh.active_elements_end(); 
      for ( ; el != end_el; ++el)
        {
          //const Elem* elem = *el;

          this->plot_linear_elem(*el);
          //this->plot_quadratic_elem(*el); // Experimental
        }
      
      // Issue the showpage command, and we're done.
      _out << 'showpage' << std::endl;
      
    } // end if (libMesh::processor_id() == 0)
}
 

virtual void MeshOutput< MeshBase >::write_equation_systems (const std::string &, const EquationSystems &) [virtual, inherited]This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.

Reimplemented in VTKIO.  

virtual void MeshOutput< MeshBase >::write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) [inline, virtual, inherited]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.

Reimplemented in ExodusII_IO, GmshIO, GMVIO, GnuPlotIO, MEDITIO, and TecplotIO.

Definition at line 91 of file mesh_output.h.

  { libmesh_error(); }
 

Member Data Documentation

 

std::vector<Point> PostscriptIO::_bezier_coeffs [private]Vector containing 3 points corresponding to Bezier coefficients, as computed by _compute_edge_bezier_coeffs.

Definition at line 115 of file postscript_io.h.

Referenced by _compute_edge_bezier_coeffs(), plot_quadratic_elem(), and PostscriptIO().  

const float PostscriptIO::_bezier_transform [static, private]Initial value:


{
  {-1./6., 1./6.,    1.}, 
  {-1./6.,   0.5, 1./6.}, 
  {    0.,    1.,    0.} 
}
Coefficients of the transformation from physical-space edge coordinates to Bezier basis coefficients. Transforms x and y separately.

Definition at line 109 of file postscript_io.h.

Referenced by _compute_edge_bezier_coeffs().  

std::ostringstream PostscriptIO::_cell_string [private]Drawing style-independent data for a single cell. This can be used as a temporary buffer for storing data which may be sent to the output stream multiple times.

Definition at line 137 of file postscript_io.h.

Referenced by plot_linear_elem().  

Point PostscriptIO::_current_point [private]A point object used for temporary calculations

Definition at line 130 of file postscript_io.h.

Referenced by plot_linear_elem(), and plot_quadratic_elem().  

Point PostscriptIO::_offset [private]Amount to add to every (x,y) point to place it in Postscript coordinates.

Definition at line 120 of file postscript_io.h.

Referenced by _compute_edge_bezier_coeffs(), plot_linear_elem(), plot_quadratic_elem(), and write().  

std::ofstream PostscriptIO::_out [private]Output file stream which will be opened when the file name is known

Definition at line 142 of file postscript_io.h.

Referenced by plot_linear_elem(), plot_quadratic_elem(), and write().  

Real PostscriptIO::_scale [private]Amount by which to stretch each point to place it in Postscript coordinates.

Definition at line 125 of file postscript_io.h.

Referenced by _compute_edge_bezier_coeffs(), plot_linear_elem(), plot_quadratic_elem(), and write().  

Real PostscriptIO::line_widthControl the thickness of the lines used. 0.5 is a reasonable default for printed images, but you may need to decrease this value (or choose it adaptively) when there are very slim cells present in the mesh.

Definition at line 83 of file postscript_io.h.

Referenced by write().  

Real PostscriptIO::shade_valueControls greyscale shading of cells. By default this value is 0.0 (which actually corresponds to black) and this indicates 'no shading' i.e. only mesh lines will be drawn. Any other value in (0,1] will cause the cells to be grey-shaded to some degree, with higher values being lighter. A value of 0.75 gives decent results.

Definition at line 75 of file postscript_io.h.

Referenced by plot_linear_elem().

 

Author

Generated automatically by Doxygen for libMesh from the source code.


 

Index

NAME
SYNOPSIS
Public Member Functions
Public Attributes
Protected Member Functions
Private Member Functions
Private Attributes
Static Private Attributes
Detailed Description
Constructor & Destructor Documentation
PostscriptIO::PostscriptIO (const MeshBase &mesh)Constructor.
PostscriptIO::~PostscriptIO () [virtual]Destructor.
Member Function Documentation
void PostscriptIO::_compute_edge_bezier_coeffs (const Elem *elem) [private]Given a quadratic edge Elem which lies in the x-y plane, computes the Bezier coefficients. These may be passed to the Postscript routine 'curveto'.
const MeshBase & MeshOutput< MeshBase >::mesh () const [protected, inherited]Returns the object as a read-only reference.
void PostscriptIO::plot_linear_elem (const Elem *elem)Draws an element with straight lines
void PostscriptIO::plot_quadratic_elem (const Elem *elem)Draws an element with Bezier curves
void PostscriptIO::write (const std::string &fname) [virtual]This method implements writing a mesh to a specified file.
virtual void MeshOutput< MeshBase >::write_equation_systems (const std::string &, const EquationSystems &) [virtual, inherited]This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.
virtual void MeshOutput< MeshBase >::write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) [inline, virtual, inherited]This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.
Member Data Documentation
std::vector<Point> PostscriptIO::_bezier_coeffs [private]Vector containing 3 points corresponding to Bezier coefficients, as computed by _compute_edge_bezier_coeffs.
const float PostscriptIO::_bezier_transform [static, private]Initial value:
std::ostringstream PostscriptIO::_cell_string [private]Drawing style-independent data for a single cell. This can be used as a temporary buffer for storing data which may be sent to the output stream multiple times.
Point PostscriptIO::_current_point [private]A point object used for temporary calculations
Point PostscriptIO::_offset [private]Amount to add to every (x,y) point to place it in Postscript coordinates.
std::ofstream PostscriptIO::_out [private]Output file stream which will be opened when the file name is known
Real PostscriptIO::_scale [private]Amount by which to stretch each point to place it in Postscript coordinates.
Real PostscriptIO::line_widthControl the thickness of the lines used. 0.5 is a reasonable default for printed images, but you may need to decrease this value (or choose it adaptively) when there are very slim cells present in the mesh.
Real PostscriptIO::shade_valueControls greyscale shading of cells. By default this value is 0.0 (which actually corresponds to black) and this indicates 'no shading' i.e. only mesh lines will be drawn. Any other value in (0,1] will cause the cells to be grey-shaded to some degree, with higher values being lighter. A value of 0.75 gives decent results.
Author

This document was created by man2html, using the manual pages.
Time: 21:52:54 GMT, April 16, 2011