]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
checkin implementation of povray output
authorrichter <richter@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Nov 1999 15:24:58 +0000 (15:24 +0000)
committerrichter <richter@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 23 Nov 1999 15:24:58 +0000 (15:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@1935 0785d39b-7218-0410-832d-ea1e28bc413d

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

index f459eb677fec886ca0e32b64c196e4246a4a1ac8..c50333670c9b183a54d0daa46dcbabecd4bf49b8 100644 (file)
  *
  * \subsection{POVRAY format}
  *
- * No information presently available.
+ * Output in this format creates a povray source file, include standard
+ * camera and light source definition for rendering with povray 3.1
+ * At present, this format only supports output for two-dimensional data,
+ * with values in the third direction taken from a data vector.
  *
+ * The output uses two different povray-objects:
+
+ * \begin{itemize}
+ *
+ * \item BICUBIC_PATCH
+ * A bicubic_patch is a 3-dimensional Bezier patch. It consists of 16 Points
+ * describing the surface. The 4 corner points are touched by the object,
+ * while the other 12 points pull and stretch the patch into shape.
+ * One bicubic_patch is generated on each patch. Therefor the number of 
+ * subdivisions has to be 3 to provide the patch with 16 points.
+ * A bicubic patch is not exact but generates very smooth images.
+ *
+ * \item MESH
+ * The mesh object is used to store large number of triangles.
+ * Every square of the patch data is split into one upper-left and one 
+ * lower-right triangle. If the number of subdivisions is three, 32 triangle
+ * are generated for every patch.
+ * 
+ * Using the smooth flag povray interpolates the normals on the triangles,
+ * imitating a curved surface
+ *
+ * \end{itemize}
+ *
+ * All objects get one texture definition called Tex. This texture has to be
+ * declared somewhere before the object data. This may be in an external 
+ * data file or at the beginning of the output file.
+ * Setting the external_data flag to false, an standard camera, light and
+ * texture (scaled to fit the scene) is added to the outputfile. Set to true
+ * an include file "data.inc" is included. This file is not generated by deal
+ * and has to include camera, light and the texture definition Tex.
+ *
+ * You need povray (>=3.0) to render the scene. The minimum Options for povray
+ * are:
+ * \begin{verbatim}
+ *   povray +I<inputfile> +W<horiz. size> +H<ver. size> +L<include path>
+ * \end{verbatim}
+ * If the external file "data.inc" is used, the path to this file has to be
+ * included in the povray options.
  *
  * \subsection{EPS format}
  *
@@ -401,18 +442,39 @@ class DataOutBase
                                      */
     struct PovrayFlags 
     {
-      private:
-                                        /**
-                                         * Dummy entry to suppress compiler
-                                         * warnings when copying an empty
-                                         * structure. Remove this member
-                                         * when adding the first flag to
-                                         * this structure (and remove the
-                                         * #private# as well).
-                                         */
-       int dummy;
+                                     /**
+                                      * Normal vector interpolation,
+                                      * if set to true
+                                      *
+                                      * default = false
+                                      */
+      bool smooth;
+
+                                      /**
+                                      * Use bicubic patches (b-splines)
+                                      * instead of triangles.
+                                      *
+                                      * default = false
+                                      */
+      bool bicubic_patch;
+
+                                      /**
+                                      * include external "data.inc"
+                                      * with camera, light and
+                                      * texture definition for the
+                                      * scene.
+                                      *
+                                      * default = false
+                                      */
+      bool external_data;
+
+                                      /**
+                                      * Constructor.
+                                      */
+      PovrayFlags (const bool smooth = false,
+                  const bool bicubic_patch = false,
+                  const bool external_data = false);
 
-      public:
                                         /**
                                          * Declare all flags with name
                                          * and type as offered by this
@@ -822,6 +884,13 @@ class DataOutBase
                                     /**
                                      * Exception
                                      */
+    DeclException2 (ExcUnexpectedNumberOfSubdivisions,
+                   int, int,
+                   << "The number of subdivisions is " << arg1
+                   << ", but we excepted "  << arg2<< " for bicubic patch");
+                                    /**
+                                     * Exception
+                                     */
     DeclException0 (ExcIO);
 
   private:
index 85107ef6ce3eccf8259083920f69811ee4a60be1..9ae7c6e81de80a14f462661b5f4f00f3e0914de8 100644 (file)
@@ -9,7 +9,13 @@ DataOutBase::UcdFlags::UcdFlags (const bool write_preamble) :
                write_preamble (write_preamble)
 {};
 
-
+DataOutBase::PovrayFlags::PovrayFlags (const bool smooth,
+                                      const bool bicubic_patch,
+                                      const bool external_data) :
+               smooth (smooth),
+               bicubic_patch(bicubic_patch),
+               external_data(external_data)
+{};
 
 void DataOutBase::UcdFlags::declare_parameters (ParameterHandler &prm)
 {
@@ -37,14 +43,23 @@ void DataOutBase::GnuplotFlags::parse_parameters (ParameterHandler &/*prm*/)
 
 
 
-void DataOutBase::PovrayFlags::declare_parameters (ParameterHandler &/*prm*/)
+void DataOutBase::PovrayFlags::declare_parameters (ParameterHandler &prm)
 {
+  prm.declare_entry ("Use smooth triangles", "false",
+                    Patterns::Bool());
+  prm.declare_entry ("Use bicubic patches", "false",
+                    Patterns::Bool());
+  prm.declare_entry ("Include external file", "true",
+                    Patterns::Bool ());
 };
 
 
 
-void DataOutBase::PovrayFlags::parse_parameters (ParameterHandler &/*prm*/)
+void DataOutBase::PovrayFlags::parse_parameters (ParameterHandler &prm)
 {
+  smooth        = prm.get_bool ("Use smooth triangles");
+  bicubic_patch = prm.get_bool ("Use bicubic patches");
+  external_data = prm.get_bool ("Include external file");
 };
 
 
index fb26d35bc806ab73c86b5a4fa102ed5f09302ebf..f2f12611c1c26fcbb7f3f7dfa7067852537985d3 100644 (file)
@@ -682,16 +682,306 @@ void DataOutBase::write_gnuplot (const vector<Patch<dim> > &patches,
 
 
 template <int dim>
-void DataOutBase::write_povray (const vector<Patch<dim> > &/*patches*/,
-                               const vector<string>      &/*data_names*/,
-                               const PovrayFlags         &/*flags*/,
-                               ostream                   &/*out*/
+void DataOutBase::write_povray (const vector<Patch<dim> > &patches,
+                               const vector<string>      &data_names,
+                               const PovrayFlags         &flags,
+                               ostream                   &out
 {
-  Assert (false, ExcNotImplemented());
+
+  AssertThrow (out, ExcIO());
+  Assert (dim==2, ExcNotImplemented());        // only for 2-D  
+
+  const unsigned int n_data_sets = data_names.size();
+  
+                                    // write preamble
+  if (true) 
+    {
+                                    // block this to have local
+                                    // variables destroyed after use
+      const time_t  time1= time (0);
+      const tm     *time = localtime(&time1); 
+      out << "/* This file was generated by the deal.II library." << endl
+         << "   Date =  "
+         << time->tm_year+1900 << "/"
+         << time->tm_mon+1 << "/"
+         << time->tm_mday << endl
+         << "   Time =  "
+         << time->tm_hour << ":"
+         << setw(2) << time->tm_min << ":"
+         << setw(2) << time->tm_sec << endl
+         << endl
+         << "   For a description of the POVRAY format see the POVRAY manual."
+         << endl
+         << "*/ " << endl;
+      
+                                    // include files
+      out << "#include \"colors.inc\" " << endl
+         << "#include \"textures.inc\" " << endl;
+
+
+                                    // use external include file for textures,
+                                    // camera and light
+      if (flags.external_data)
+       out << "#include \"data.inc\" " << endl;
+      else                          // all definitions in data file
+       {  
+                                    // camera
+         out << endl << endl
+             << "camera {"            << endl
+             << "  location <1,4,-7>" << endl
+             << "  look_at <0,0,0>"   << endl
+             << "  angle 30"          << endl
+             << "}"                   << endl;
+      
+                                    // light
+         out << endl 
+             << "light_source {"      << endl
+             << "  <1,4,-7>"      << endl
+             << "  color Grey"        << endl
+             << "}"                   << endl;
+         out << endl 
+             << "light_source {"      << endl
+             << "  <0,20,0>"      << endl
+             << "  color White"       << endl
+             << "}"                   << endl;
+       }
+    };
+
+  double hmin=0,hmax=0;                             // max. and min. heigth of solution 
+
+  for (typename vector<Patch<dim> >::const_iterator patch=patches.begin();
+       patch != patches.end(); ++patch)
+    {
+      const unsigned int n_subdivisions = patch->n_subdivisions;
+      
+      Assert (patch->data.m() == n_data_sets,
+             ExcUnexpectedNumberOfDatasets (patch->data.m(), n_data_sets));
+      Assert (patch->data.n() == (dim==1 ?
+                                 n_subdivisions+1 :
+                                 (dim==2 ?
+                                  (n_subdivisions+1)*(n_subdivisions+1) :
+                                  (dim==3 ?
+                                   (n_subdivisions+1)*(n_subdivisions+1)*(n_subdivisions+1) :
+                                   0))),
+             ExcInvalidDatasetSize (patch->data.n(), n_subdivisions+1));
+      
+      for (unsigned int i=0; i<n_subdivisions; ++i)
+       for (unsigned int j=0; j<n_subdivisions; ++j)
+         {
+           const int dl = i*(n_subdivisions+1)+j;
+           if ((hmin==0)||(patch->data(0,dl)<hmin)) hmin=patch->data(0,dl);
+           if ((hmax==0)||(patch->data(0,dl)>hmax)) hmax=patch->data(0,dl);
+         }
+    }
+
+  out << "#declare HMIN=" << hmin << ";" << endl
+      << "#declare HMAX=" << hmax << ";" << endl << endl;
+
+  if (!flags.external_data)
+    {
+                                    // texture with scaled niveau lines
+                                    // 10 lines in the surface
+      out << "#declare Tex=texture{" << endl
+         << "  pigment {" << endl
+         << "    gradient y" << endl
+         << "    scale y*(HMAX-HMIN)*" << 0.1 << endl
+         << "    color_map {" << endl
+         << "      [0.00 color Light_Purple] " << endl
+         << "      [0.95 color Light_Purple] " << endl
+         << "      [1.00 color White]    " << endl
+         << "} } }" << endl << endl;
+    }
+
+  if (!flags.bicubic_patch)
+    {                                  // start of mesh header
+      out << endl
+         << "mesh {" << endl;
+    }
+
+                                     // loop over all patches
+  for (typename vector<Patch<dim> >::const_iterator patch=patches.begin();
+       patch != patches.end(); ++patch)
+    {
+      const unsigned int n_subdivisions = patch->n_subdivisions;
+      
+      Assert (patch->data.m() == n_data_sets,
+             ExcUnexpectedNumberOfDatasets (patch->data.m(), n_data_sets));
+      Assert (patch->data.n() == (dim==1 ?
+                                 n_subdivisions+1 :
+                                 (dim==2 ?
+                                  (n_subdivisions+1)*(n_subdivisions+1) :
+                                  (dim==3 ?
+                                   (n_subdivisions+1)*(n_subdivisions+1)*(n_subdivisions+1) :
+                                   0))),
+             ExcInvalidDatasetSize (patch->data.n(), n_subdivisions+1));
+      
+      
+      Point<dim> ver[16];                          // value for all points in this patch
+      
+      for (unsigned int i=0; i<n_subdivisions+1; ++i)
+       {
+         for (unsigned int j=0; j<n_subdivisions+1; ++j)
+           {
+             const double x_frac = i * 1./n_subdivisions,
+               y_frac = j * 1./n_subdivisions;
+                                                  // compute coordinates for
+                                                  // this patch point, storing in ver
+             ver[i*(n_subdivisions+1)+j]= (((patch->vertices[1] * x_frac) +
+                                            (patch->vertices[0] * (1-x_frac))) * (1-y_frac) +
+                                           ((patch->vertices[2] * x_frac) +
+                                            (patch->vertices[3] * (1-x_frac))) * y_frac);
+           };
+       };
+      
+      if (!flags.bicubic_patch)
+       {                                    // setting up triangles
+         for (unsigned int i=0; i<n_subdivisions; ++i)
+           {
+             for (unsigned int j=0; j<n_subdivisions; ++j)
+               {
+                                                    // down/left vertex of triangle
+                 const int dl = i*(n_subdivisions+1)+j;
+                 if (flags.smooth)               // only if smooth triangles are used
+                   {
+                     Point<3> nrml[16];     // aproximate normal vectors in patch
+                     for (unsigned int i=0; i<n_subdivisions+1;++i)
+                       {
+                         for (unsigned int j=0; j<n_subdivisions+1;++j)
+                           {
+                             if (i==0) 
+                               nrml[i*(n_subdivisions+1)+j](0)=-patch->data(0,i*(n_subdivisions+1)+j);
+                             else 
+                               nrml[i*(n_subdivisions+1)+j](0)=-patch->data(0,(i-1)*(n_subdivisions+1)+j);
+                             if (i==n_subdivisions) 
+                               nrml[i*(n_subdivisions+1)+j](0)+=patch->data(0,i*(n_subdivisions+1)+j);
+                             else 
+                               nrml[i*(n_subdivisions+1)+j](0)+=patch->data(0,(i+1)*(n_subdivisions+1)+j);
+                             nrml[i*(n_subdivisions+1)+j](0)/=(1./n_subdivisions);
+                             if (i!=n_subdivisions&&i!=0) nrml[i*(n_subdivisions+1)+j](0)/=2;
+                             if (j==0)
+                               nrml[i*(n_subdivisions+1)+j](2)=-patch->data(0,i*(n_subdivisions+1)+j);
+                             else
+                               nrml[i*(n_subdivisions+1)+j](2)=-patch->data(0,i*(n_subdivisions+1)+j-1);
+                             if (j==n_subdivisions)
+                               nrml[i*(n_subdivisions+1)+j](2)+=patch->data(0,i*(n_subdivisions+1)+j);
+                             else
+                               nrml[i*(n_subdivisions+1)+j](2)+=patch->data(0,i*(n_subdivisions+1)+j+1);
+                             nrml[i*(n_subdivisions+1)+j](2)/=-(1./n_subdivisions);
+                             if (j!=n_subdivisions&&j!=0) nrml[i*(n_subdivisions+1)+j](2)/=2;
+                             nrml[i*(n_subdivisions+1)+j](1)=1.;
+
+                                                               // normalize Vektor
+                             double norm=sqrt(
+                                             pow(nrml[i*(n_subdivisions+1)+j](0),2.)+
+                                             pow(nrml[i*(n_subdivisions+1)+j](1),2.)+1.);
+                             
+                             for (unsigned int k=0;k<3;++k)  nrml[i*(n_subdivisions+1)+j](k)/=norm;
+                           }
+                       }
+                                                            // writing smooth_triangles
+                               
+                                                            // down/right triangle
+                     out << "smooth_triangle {" << endl << "\t<" 
+                         << ver[dl](0) << ","   
+                         << patch->data(0,dl) << ","
+                         << ver[dl](1) << ">, <"
+                         << nrml[dl] << ">," << endl
+                         << " \t<" 
+                         << ver[dl+n_subdivisions+1](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+1)  << ","
+                         << ver[dl+n_subdivisions+1](1) << ">, <"
+                         << nrml[dl+n_subdivisions+1] << ">," << endl 
+                         << "\t<" 
+                         << ver[dl+n_subdivisions+2](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+2)  << ","
+                         << ver[dl+n_subdivisions+2](1) << ">, <"
+                         << nrml[dl+n_subdivisions+2] << ">}" << endl; 
+                     
+                                                    // upper/left triangle
+                     out << "smooth_triangle {" << endl << "\t<" 
+                         << ver[dl](0) << "," 
+                         << patch->data(0,dl) << ","
+                         << ver[dl](1) << ">, <"
+                         << nrml[dl] << ">," << endl 
+                         << "\t<" 
+                         << ver[dl+n_subdivisions+2](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+2)  << ","
+                         << ver[dl+n_subdivisions+2](1) << ">, <"
+                         << nrml[dl+n_subdivisions+2] << ">," << endl 
+                         << "\t<" 
+                         << ver[dl+1](0) << "," 
+                         << patch->data(0,dl+1)  << ","
+                         << ver[dl+1](1) << ">, <"
+                         << nrml[dl+1] << ">}" << endl;
+                   }
+                 else
+                   {           
+                                           // writing standard triangles
+                                            // down/right triangle
+                     out << "triangle {" << endl << "\t<" 
+                         << ver[dl](0) << "," 
+                         << patch->data(0,dl) << ","
+                         << ver[dl](1) << ">," << endl
+                         << "\t<" 
+                         << ver[dl+n_subdivisions+1](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+1)  << ","
+                         << ver[dl+n_subdivisions+1](1) << ">," << endl 
+                         << "\t<" 
+                         << ver[dl+n_subdivisions+2](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+2)  << ","
+                         << ver[dl+n_subdivisions+2](1) << ">}" << endl; 
+                           
+                                                    // upper/left triangle
+                     out << "triangle {" << endl << "\t<" 
+                         << ver[dl](0) << "," 
+                         << patch->data(0,dl) << ","
+                         << ver[dl](1) << ">," << endl 
+                         << "\t<"
+                         << ver[dl+n_subdivisions+2](0) << "," 
+                         << patch->data(0,dl+n_subdivisions+2)  << ","
+                         << ver[dl+n_subdivisions+2](1) << ">," << endl 
+                         << "\t<" 
+                         << ver[dl+1](0) << "," 
+                         << patch->data(0,dl+1)  << ","
+                         << ver[dl+1](1) << ">}" << endl;
+                   };
+               };
+           };
+       }
+      else
+       {                                    // writing bicubic_patch
+         Assert (n_subdivisions==3, ExcUnexpectedNumberOfSubdivisions(n_subdivisions,3));
+         out << endl
+             << "bicubic_patch {" << endl
+             << "  type 0" << endl
+             << "  flatness 0" << endl
+             << "  u_steps 0" << endl
+             << "  v_steps 0" << endl;
+         for (int i=0;i<16;++i)
+           {
+             out << "\t<" << ver[i](0) << "," << patch->data(0,i) << "," << ver[i](1) << ">";
+             if (i!=15) out << ",";
+             out << endl;
+           };
+         out << "  texture {Tex}" <<  endl
+             << "}" << endl;
+       };
+    };
+  
+  if (!flags.bicubic_patch) 
+    {                                   // the end of the mesh
+      out << "  texture {Tex}" << endl
+         << "}" << endl
+         << endl;
+    }
+  
+  AssertThrow (out, ExcIO());
 };
 
 
 
+
 template <int dim>
 void DataOutBase::write_eps (const vector<Patch<dim> > &patches,
                             const vector<string>      &/*data_names*/,

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.