* mln/io/vtk/save.hh: New.
---
milena/ChangeLog | 6 +
milena/mln/io/vtk/save.hh | 623 +++++++++++++++++++++++++++++++++++++++++++++
2 files changed, 629 insertions(+), 0 deletions(-)
create mode 100644 milena/mln/io/vtk/save.hh
diff --git a/milena/ChangeLog b/milena/ChangeLog
index 3d8b6bc..2cb03b7 100644
--- a/milena/ChangeLog
+++ b/milena/ChangeLog
@@ -1,3 +1,9 @@
+2010-06-24 Roland Levillain <roland(a)lrde.epita.fr>
+
+ Start a VTK output for complex-based images.
+
+ * mln/io/vtk/save.hh: New.
+
2010-07-15 Roland Levillain <roland(a)lrde.epita.fr>
mesh-complex-pinv-curv-segm and mesh-complex-pinv-curv-skel apps.
diff --git a/milena/mln/io/vtk/save.hh b/milena/mln/io/vtk/save.hh
new file mode 100644
index 0000000..10c3f0b
--- /dev/null
+++ b/milena/mln/io/vtk/save.hh
@@ -0,0 +1,623 @@
+// Copyright (C) 2008, 2009, 2010 EPITA Research and Development
+// Laboratory (LRDE)
+//
+// This file is part of Olena.
+//
+// Olena is free software: you can redistribute it and/or modify it under
+// the terms of the GNU General Public License as published by the Free
+// Software Foundation, version 2 of the License.
+//
+// Olena is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with Olena. If not, see <http://www.gnu.org/licenses/>.
+//
+// As a special exception, you may use this file as part of a free
+// software project without restriction. Specifically, if other files
+// instantiate templates or use macros or inline functions from this
+// file, or you compile this file and link it with other files to produce
+// an executable, this file does not by itself cause the resulting
+// executable to be covered by the GNU General Public License. This
+// exception does not however invalidate any other reasons why the
+// executable file might be covered by the GNU General Public License.
+
+#ifndef MLN_IO_VTK_SAVE_HH
+# define MLN_IO_VTK_SAVE_HH
+
+/// \file
+/// Input saving function for VTK files.
+///
+/// \see
http://www.vtk.org/VTK/img/file-formats.pdf
+/// \see
http://dunne.uni-hd.de/VisuSimple/documents/vtkfileformat.html
+
+# include <cstdlib>
+
+# include <iostream>
+# include <fstream>
+# include <sstream>
+
+# include <string>
+
+# include <mln/core/alias/complex_image.hh>
+# include <mln/core/image/complex_neighborhoods.hh>
+# include <mln/core/image/complex_neighborhood_piter.hh>
+
+
+namespace mln
+{
+
+ namespace io
+ {
+
+ namespace vtk
+ {
+
+ /** \brief Save a (binary) VTK image into a complex image.
+
+ \param[in] ima The image to save.
+ \param[in] filename The name of the file where to save the image.
+
+ The image is said binary since data represent only the
+ existence of faces. */
+ void save(const bin_2complex_image3df& ima,
+ const std::string& filename);
+
+ /** \brief Save an 8-bit grey-level VTK image into a complex image.
+
+ \param[in] ima The image to save.
+ \param[in] filename The name of the file where to save the image.
+
+ Only data is attached to 2-faces is saved; the VTK file
+ cannot store data attached to faces of other dimensions. */
+ void save(const int_u8_2complex_image3df& ima,
+ const std::string& filename);
+
+ /** \brief Save a floating-point value grey-level VTK image into
+ a complex image.
+
+ \param[in] ima The image to save.
+ \param[in] filename The name of the file where to save the image.
+
+ Only data is attached to 2-faces is saved; the VTK file
+ cannot store data attached to faces of other dimensions. */
+ void save(const float_2complex_image3df& ima,
+ const std::string& filename);
+
+ /** \brief Save a 3x8-bit RGB (color) VTK image into a complex image.
+
+ \param[in] ima The image to save.
+ \param[in] filename The name of the file where to save the image.
+
+ Only data is attached to 2-faces is saved; the VTK file
+ cannot store data attached to faces of other dimensions. */
+ void save(const rgb8_2complex_image3df& ima,
+ const std::string& filename);
+
+
+ namespace internal
+ {
+
+ template <typename I, typename E>
+ struct vtk_saver : public Object<E>
+ {
+ /// Type of the image.
+ typedef I image;
+
+ /// Type of the values.
+ typedef mln_value(I) value;
+
+ /// Dimension of the built complex.
+ static const unsigned D = 2;
+
+ /// \brief Constructor, with static checks.
+ vtk_saver();
+
+ /// Save an image \a ima into \a filename.
+ void operator()(const I& ima, const std::string& filename) const;
+
+ protected:
+ /// Helper factoring the task of writing scalar data
+ /// associated to faces.
+ void write_scalar_data(std::ostream& ostr, const image& ima,
+ const std::string& data_type) const;
+ };
+
+
+ struct bin_vtk_saver
+ : public vtk_saver< bin_2complex_image3df, bin_vtk_saver >
+ {
+ /// \brief Save face data.
+ void write_face_data(std::ostream& ostr, const image& ima) const;
+ };
+
+ struct int_u8_vtk_saver
+ : public vtk_saver< int_u8_2complex_image3df, int_u8_vtk_saver >
+ {
+ /// \brief Save face data.
+ void write_face_data(std::ostream& ostr, const image& ima) const;
+ };
+
+
+ struct float_vtk_saver
+ : public vtk_saver< float_2complex_image3df, float_vtk_saver >
+ {
+ /// \brief Save face data.
+ void write_face_data(std::ostream& ostr, const image& ima) const;
+ };
+
+
+ struct rgb8_vtk_saver
+ : public vtk_saver< rgb8_2complex_image3df, rgb8_vtk_saver >
+ {
+ /// \brief Save face data.
+ void write_face_data(std::ostream& ostr, const image& ima) const;
+ };
+
+ } // end of namespace mln::io::vtk::internal
+
+
+
+# ifndef MLN_INCLUDE_ONLY
+
+ /*----------.
+ | Facades. |
+ `----------*/
+
+ void
+ save(const bin_2complex_image3df& ima, const std::string& filename)
+ {
+ trace::entering("mln::io::vtk::save");
+ internal::bin_vtk_saver()(ima, filename);
+ trace::exiting("mln::io::vtk::save");
+ }
+
+ void
+ save(const int_u8_2complex_image3df& ima, const std::string& filename)
+ {
+ trace::entering("mln::io::vtk::save");
+ internal::int_u8_vtk_saver()(ima, filename);
+ trace::exiting("mln::io::vtk::save");
+ }
+
+ void
+ save(const float_2complex_image3df& ima, const std::string& filename)
+ {
+ trace::entering("mln::io::vtk::save");
+ internal::float_vtk_saver()(ima, filename);
+ trace::exiting("mln::io::vtk::save");
+ }
+
+ void
+ save(const rgb8_2complex_image3df& ima, const std::string& filename)
+ {
+ trace::entering("mln::io::vtk::save");
+ internal::rgb8_vtk_saver()(ima, filename);
+ trace::exiting("mln::io::vtk::save");
+ }
+
+
+ /*-------------------------.
+ | Actual implementations. |
+ `-------------------------*/
+
+ // -------- //
+ // Canvas. //
+ // -------- //
+
+ namespace internal
+ {
+
+ template <typename I, typename E>
+ vtk_saver<I, E>::vtk_saver()
+ {
+ // Concept checking.
+ void (E::*m1)(std::ostream&, const I&) const =
+ &E::write_face_data;
+ m1 = 0;
+ }
+
+
+ template <typename I, typename E>
+ void
+ vtk_saver<I, E>::operator()(const I& ima,
+ const std::string& filename) const
+ {
+ const std::string me = "mln::io::vtk::save";
+
+ std::ofstream ostr(filename.c_str());
+ if (!ostr)
+ {
+ std::cerr << me << ": `" << filename <<
"' invalid file."
+ << std::endl;
+ /* FIXME: Too violent. We should allow the use of
+ exceptions, at least to have Milena's code behave
+ correctly in interpreted environments (std::exit() or
+ std::abort() causes the termination of a Python
+ interpreter, for instance!). */
+ std::exit(1);
+ }
+
+ /*---------.
+ | Header. |
+ `---------*/
+
+ /* ``The legacy VTK file formats consist of five basic
+ parts.'' */
+
+ /* ``1. The first part is the file version and
+ identifier. This part contains the single line:
+
+ # vtk DataFile Version x.x.
+
+ This line must be exactly as shown with the
+ exception of the version number x.x, which will vary
+ with different releases of VTK. (Note: the current
+ version number is 3.0. Version 1.0 and 2.0 files are
+ compatible with version 3.0 files.)'' */
+ ostr << "# vtk DataFile Version 2.0" << std::endl;
+
+ /* ``2. The second part is the header. The header consists
+ of a character string terminated by end-of-line
+ character `\n'. The header is 256 characters
+ maximum. The header can be used to describe the data
+ and include any other pertinent information.'' */
+ ostr << "Generated by Milena 1.0
http://olena.lrde.epita.fr"
+ << std::endl;
+
+ /* ``3. The next part is the file format. The file format
+ describes the type of file, either ASCII or
+ binary. On this line the single word ASCII or BINARY
+ must appear.'' */
+ ostr << "ASCII" << std::endl;
+
+ /*-------.
+ | Data. |
+ `-------*/
+
+ /* ``4. The fourth part is the dataset structure. The
+ geometry part describes the geometry and topology of
+ the dataset. This part begins with a line containing
+ the keyword DATASET followed by a keyword describing
+ the type of dataset. Then, depending upon the type
+ of dataset, other keyword/data combinations define
+ the actual data.''
+
+ [...]
+
+ Dataset Format. The Visualization Toolkit supports
+ five different dataset formats: structured points,
+ structured grid, rectilinear grid, unstructured
+ grid, and polygonal data.'' */
+
+ ostr << "DATASET POLYDATA" << std::endl << std::endl;
+
+ // --------- //
+ // Complex. //
+ // --------- //
+
+ typedef mln_geom(I) G;
+
+ /* ``* Polygonal Data
+ The polygonal dataset consists of arbitrary
+ combinations of surface graphics primitives
+ vertices (and polyvertices), lines (and
+ polylines), polygons (of various types), and
+ triangle strips. Polygonal data is defined by
+ the POINTS, VERTICES, LINES, POLYGONS, or
+ TRIANGLE_STRIPS sections. The POINTS definition
+ is the same as we saw for structured grid
+ datasets.'' */
+
+ // ---------------------------- //
+ // Geometry (point locations). //
+ // ---------------------------- //
+
+ ostr << "POINTS "
+ << ima.domain().cplx().template nfaces_of_static_dim<0>()
+ << " float" << std::endl;
+ // Iterate on 0-faces (vertices).
+ p_n_faces_fwd_piter<D, G> v(ima.domain(), 0);
+ for_all(v)
+ {
+ mln_invariant(v.to_site().size() == 1);
+ ostr << v.to_site().front()[0] << ' '
+ << v.to_site().front()[1] << ' '
+ << v.to_site().front()[2] << std::endl;
+ }
+ ostr << std::endl;
+
+ /* ``The VERTICES, LINES, POLYGONS, or
+ TRIANGLE_STRIPS keywords define the polygonal
+ dataset topology. Each of these keywords
+ requires two parameters: the number of cells `n'
+ and the size of the cell list `size'. The cell
+ list size is the total number of integer values
+ required to represent the list (i.e., sum of
+ `numPoints' and connectivity indices over each
+ cell). None of the keywords VERTICES, LINES,
+ POLYGONS, or TRIANGLE_STRIPS is required.'' */
+
+ // ---------- //
+ // Vertices. //
+ // ---------- //
+
+ /* We do not use
+
+ ima.domain().cplx().template nfaces_of_static_dim<N>()
+
+ to get the number of N-faces, since the image may be
+ masked, and exhibit less N-faces than its underlying
+ complex. Iterating on the N-faces is safer. */
+ /* FIXME: Is there anything faster? See what the interface
+ of the (morphed) image can provide. */
+ unsigned nvertices = 0;
+ for_all(v)
+ ++nvertices;
+
+ if (nvertices > 0)
+ {
+ ostr << "VERTICES " << nvertices << ' '
+ /* Each vertex requires two numbers: the cardinal of
+ its ends (which is always 1) and the indices of the
+ point among the POINTS section. Hence the total
+ number of values in the VERTEX section is
+ nvertices * 2. */
+ << nvertices * 2 << std::endl;
+
+ for_all(v)
+ ostr << "1 " << v.unproxy_().face().face_id() <<
std::endl;
+ ostr << std::endl;
+ }
+
+ // ------- //
+ // Edges. //
+ // ------- //
+
+ // Same comment as above about the count of N-faces.
+ unsigned nedges = 0;
+ p_n_faces_fwd_piter<D, G> e(ima.domain(), 1);
+ for_all (e)
+ ++nedges;
+
+ if (nedges > 0)
+ {
+ ostr << "LINES " << nedges << ' '
+ /* Each edge requires three numbers: the cardinal of
+ its ends (which is always 2) and the indices of
+ these ends among the POINTS section. Hence the
+ total number of values in the LINES section is
+ nedges * 3. */
+ << nedges * 3 << std::endl;
+
+ // Vertices adjacent to edges.
+ typedef complex_lower_neighborhood<D, G> adj_vertices_nbh_t;
+ adj_vertices_nbh_t adj_vertices_nbh;
+ mln_niter(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, e);
+ // Iterate on 1-faces (edges).
+ for_all (e)
+ {
+ ostr << "2";
+ // Iterate on vertices (0-faces).
+ for_all (adj_v)
+ {
+ // FIXME: Accessing the face id is too complicated.
+ ostr << " " << adj_v.unproxy_().face().face_id();
+ }
+ ostr << std::endl;
+ }
+ ostr << std::endl;
+ }
+
+ // ---------- //
+ // Polygons. //
+ // ---------- //
+
+ // Same comment as above about the count of N-faces.
+ unsigned npolygons = 0;
+ p_n_faces_fwd_piter<D, G> p(ima.domain(), 2);
+ for_all (p)
+ ++npolygons;
+
+ if (npolygons > 0)
+ {
+ // A neighborhood where neighbors are the set of 0-faces
+ // transitively adjacent to the reference point.
+ typedef complex_m_face_neighborhood<D, G> nbh_t;
+ nbh_t nbh;
+ mln_fwd_niter(nbh_t) u(nbh, p);
+ /* FIXME: We should be able to pass this value (m)
+ either at the construction of the neighborhood or at
+ the construction of the iterator. */
+ u.iter().set_m(0);
+
+ /* Compute the number of values (`size') to be passed as
+ second parameter of the POLYGONS keyword. */
+ unsigned polygons_size = 0;
+ // Iterate on polygons (2-face).
+ for_all(p)
+ {
+ unsigned nvertices = 0;
+ /* FIXME: There may be a faster way to do this (e.g.,
+ the neighbordhood may provide a method returning
+ the number of P's neighbors. */
+ // Iterate on transitively adjacent vertices (0-face).
+ for_all(u)
+ ++nvertices;
+ // The number of values describing this polygon P is
+ // the cardinal of its set of vertices (1 value) plus
+ // the NVERTICES indices of these vertices.
+ polygons_size += 1 + nvertices;
+ }
+ ostr << "POLYGONS " << npolygons << ' '
<< polygons_size
+ << std::endl;
+
+ /* Output polygons (one per line), with their number of
+ vertices and the indices of these vertices. */
+ // Iterate on polygons (2-face).
+ for_all(p)
+ {
+ unsigned nvertices = 0;
+ std::ostringstream vertices;
+ // Iterate on transitively adjacent vertices (0-face).
+ for_all(u)
+ {
+ // FIXME: Likewise, this is a bit too long.
+ vertices << ' ' << u.unproxy_().face().face_id();
+ ++nvertices;
+ }
+ ostr << nvertices << vertices.str() << std::endl;
+ }
+ ostr << std::endl;
+ }
+
+ /* ``5. The final part describes the dataset attributes.
+ This part begins with the keywords POINT_DATA or
+ CELL_DATA,followed by an integer number specifying
+ the number of points or cells, respectively. (It
+ doesn't matter whether POINT_DATA or CELL_DATA comes
+ first.) Other keyword/data combinations then define
+ the actual dataset attribute values (i.e., scalars,
+ vectors, tensors, normals, texture coordinates, or
+ field data).'' */
+
+ unsigned nfaces = nvertices + nedges + npolygons;
+ if (nfaces > 0)
+ {
+ // We don't use POINT_DATA (to associate values to
+ // POINTs), since CELL_DATA is used to associate values
+ // to VERTICES, EDGES and POLYGONS.
+ ostr << "CELL_DATA " << nfaces << std::endl;
+ exact(this)->write_face_data(ostr, ima);
+ }
+
+ ostr.close();
+ }
+
+ // ---------------- //
+ // Specific parts. //
+ // ---------------- //
+
+ /* ``Dataset Attribute Format. The Visualization Toolkit
+ supports the following dataset attributes: scalars
+ (one to four components), vectors, normals,
+ texture coordinates (1D, 2D, and 3D), 3 x 3
+ tensors, and field data. In addition, a lookup
+ table using the RGBA color specification,
+ associated with the scalar data, can be defined as
+ well. Dataset attributes are supported for both
+ points and cells.
+ Each type of attribute data has a `dataName'
+ associated with it. This is a character string
+ (without embedded whitespace) used to identify a
+ particular data. The `dataName' is used by the VTK
+ readers to extract data. As a result, more than
+ one attribute data of the same type can be
+ included in a file. For example, two different
+ scalar fields defined on the dataset points,
+ pressure and temperature, can be contained in the
+ same file. (If the appropriate dataName is not
+ specified in the VTK reader, then the first data
+ of that type is extracted from the file.)
+
+ * Scalars
+ Scalar definition includes specification of a
+ lookup table. The definition of a lookup table
+ is optional.
+
+ [...]
+
+ SCALARS dataName dataType numComp
+ LOOKUP_TABLE tableName
+ s0
+ s1
+ ...
+ sn-1''
+
+ Note: values accepted by Paraview 3.8 for `dataType' are:
+ "bit", "char", "unsigned_char", "short",
"unsigned_short",
+ "vtkidtype", "int", "unsigned_int", "long",
+ "unsigned_long", "vtktypeint64", "vtktypeuint64",
"float",
+ "double", "string", "utf8_string" and
"variant". */
+
+ /** \{ */
+ template <typename I, typename E>
+ void
+ vtk_saver<I, E>::write_scalar_data(std::ostream& ostr,
+ const image& ima,
+ const std::string& data_type) const
+ {
+ ostr << "SCALARS values " << data_type << std::endl
+ << "LOOKUP_TABLE default" << std::endl;
+ // Iterate on all faces, dimension increasing.
+ mln_fwd_piter(image) p(ima.domain());
+ for_all(p)
+ ostr << ima(p) << std::endl;
+ }
+
+
+ void
+ bin_vtk_saver::write_face_data(std::ostream& ostr,
+ const image& ima) const
+ {
+ write_scalar_data(ostr, ima, "bit");
+ }
+
+ void
+ int_u8_vtk_saver::write_face_data(std::ostream& ostr,
+ const image& ima) const
+ {
+ write_scalar_data(ostr, ima, "unsigned_char");
+ }
+
+ void
+ float_vtk_saver::write_face_data(std::ostream& ostr,
+ const image& ima) const
+ {
+ write_scalar_data(ostr, ima, "float");
+ }
+
+ /* ``The definition of color scalars (i.e., unsigned
+ char values directly mapped to color) varies
+ depending upon the number of values (`nValues')
+ per scalar. If the file format is ASCII, the
+ color scalars are defined using nValues float
+ values between (0,1).
+
+ COLOR_SCALARS dataName nValues
+ c00 c01 ... c0(nValues-1)
+ c10 c11 ... c1(nValues-1)
+ ...
+ c(n-1)0 c(n-1)1 ... c(n-1)(nValues-1)'' */
+
+ void
+ rgb8_vtk_saver::write_face_data(std::ostream& ostr,
+ const image& ima) const
+ {
+ ostr << "COLOR_SCALARS values 4" << std::endl;
+ // Iterate on all faces, dimension increasing.
+ mln_fwd_piter_(image) p(ima.domain());
+ for_all(p)
+ // RGBA values (with alpha channel always set to 1.0).
+ ostr << float(ima(p).red()) / mln_max(value::red_t) << ' '
+ << float(ima(p).green()) / mln_max(value::green_t) << ' '
+ << float(ima(p).blue()) / mln_max(value::blue_t) << ' '
+ << 1.f
+ << std::endl;
+ }
+ /** \} */
+
+ } // end of namespace mln::io::vtk::internal
+
+# endif // ! MLN_INCLUDE_ONLY
+
+ } // end of namespace mln::io::vtk
+
+ } // end of namespace mln::io
+
+} // end of namespace mln
+
+
+#endif // ! MLN_IO_VTK_SAVE_HH
--
1.5.6.5