From: christian.wuelker Date: Mon, 15 Apr 2013 11:21:24 +0000 (+0000) Subject: Added a perspective view to GridOut::write_svg X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=933eb59c20059ab7cf13b9de15832bbe31aec598;p=dealii-svn.git Added a perspective view to GridOut::write_svg git-svn-id: https://svn.dealii.org/trunk@29274 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/grid/grid_out.h b/deal.II/include/deal.II/grid/grid_out.h index e9c75bd23f..0382975532 100644 --- a/deal.II/include/deal.II/grid/grid_out.h +++ b/deal.II/include/deal.II/grid/grid_out.h @@ -618,40 +618,46 @@ namespace GridOutFlags unsigned int boundary_line_thickness; /// Margin around the plotted area - unsigned int margin_in_percent; + bool margin; /** * Background style. */ - enum Background { + enum Background{ /// Use transparent value of SVG transparent, /// Use white background white, - /// Use a gradient from white (top) to steelblue (bottom), and add date and time plus a "deal.II" logo. - dealii }; + /// Use a gradient from white (top) to steelblue (bottom), and add date and time plus a deal.II logo. Automatically draws a margin. + dealii}; Background background; + // View angles for the perspective view of the grid; Default is 0, 0 (top view). + int azimuth_angle, polar_angle; + /** * Cell coloring. */ - enum Coloring { + enum Coloring{ /// No cell coloring none, - /// Convert the material id into the cell color + /// Convert the material id into the cell color (default) material_id, /// Convert the level number into the cell color level_number, /// Convert the subdomain id into the cell color - subdomain_id }; + subdomain_id}; Coloring coloring; + // Interpret the level number of the cells as altitude over the x-y-plane (may be useful in the perpspective view). + bool convert_level_number_to_height; + /** * Cell labeling (fixed order). * - * The following booleans determine which property of the cell + * The following booleans determine which properties of the cell * shall be displayed as text in the middle of each cell. */ bool label_level_number; /// default: true @@ -668,17 +674,20 @@ namespace GridOutFlags /** * Constructor. */ - Svg (const unsigned int line_thickness = 3, - const unsigned int boundary_line_thickness = 7, - const unsigned int margin_in_percent = 7, - const Background background = dealii, - const Coloring coloring = material_id, - const bool label_level_number = true, - const bool label_cell_index = true, - const bool label_material_id = false, - const bool label_subdomain_id = false, - const bool draw_colorbar = true, - const bool draw_legend = true); + Svg(const unsigned int line_thickness = 3, + const unsigned int boundary_line_thickness = 7, + bool margin = true, + const Background background = dealii, + const int azimuth_angle = 0, + const int polar_angle = 0, + const Coloring coloring = material_id, + const bool convert_level_number_to_height = false, + const bool label_level_number = true, + const bool label_cell_index = true, + const bool label_material_id = false, + const bool label_subdomain_id = false, + const bool draw_colorbar = true, + const bool draw_legend = true); /** * Declare parameters in @@ -796,7 +805,8 @@ namespace GridOutFlags * * @ingroup grid * @ingroup output - * @author Wolfgang Bangerth, Guido Kanschat, Luca Heltai, 1999, 2003, 2006; SVG by Christian Wülker, 2013; postscript format based on an implementation by Stefan Nauber, 1999 + * @author Wolfgang Bangerth, Guido Kanschat, Luca Heltai, Stefan Nauber, Christian Wülker + * @date 1999 - 2013 */ class GridOut { @@ -1028,10 +1038,13 @@ public: /** * Write the triangulation in the SVG format. * - * SVG (Scalable Vector Graphics) is an XML-based vector image - * format recommended by the World Wide Web Consortium (W3C). This - * function conforms to the latest specification SVG 1.1, released - * on August 16, 2011. + * SVG (Scalable Vector Graphics) is + * an XML-based vector image format + * developed and maintained by the + * World Wide Web Consortium (W3C). + * This function conforms to the + * latest specification SVG 1.1, + * released on August 16, 2011. * * The cells of the triangulation are written as polygons with * additional lines at the boundary of the triangulation. A coloring @@ -1040,8 +1053,9 @@ public: * colorbar can be drawn to encode the chosen coloring. Moreover, a * cell label can be added, showing level index, etc. * - * @note Only implemented for two-dimensional grids in two space - * dimensions. + * @note Yet only implemented for + * two-dimensional grids in two + * space dimensions. * */ template @@ -1478,12 +1492,36 @@ private: const unsigned int starting_index, std::ostream &out) const; + + /** + * This function projects a three-dimensional + * point (Point<3> point) onto a two-dimensional + * image plane, specified by the position of + * the camera viewing system + * (Point<3> camera_position), + * camera direction + * (Point<3> camera_position), + * camera horizontal + * (Point<3> camera_horizontal, necessary for + * the correct alignment of the later images), + * and the focus of the camera + * (float camera_focus). + * + * For SVG output of grids. + */ + Point<2> svg_project_point(Point<3> point, + Point<3> camera_position, + Point<3> camera_direction, + Point<3> camera_horizontal, + float camera_focus) const; + /** * Return the number of faces in the triangulation which have a * boundary indicator not equal to zero. Only these faces are * explicitly printed in the write_* functions; all faces * with indicator numbers::internal_face_boundary_id are interior * ones and an indicator with value zero for faces at the boundary + * are considered default. * * This function always returns an empty list in one dimension. diff --git a/deal.II/source/grid/grid_out.cc b/deal.II/source/grid/grid_out.cc index 717e0d8607..b2db17f04b 100644 --- a/deal.II/source/grid/grid_out.cc +++ b/deal.II/source/grid/grid_out.cc @@ -342,9 +342,12 @@ namespace GridOutFlags Svg::Svg(const unsigned int line_thickness, const unsigned int boundary_line_thickness, - const unsigned int margin_in_percent, + bool margin, const Background background, + const int azimuth_angle, + const int polar_angle, const Coloring coloring, + const bool convert_level_number_to_height, const bool label_level_number, const bool label_cell_index, const bool label_material_id, @@ -353,9 +356,12 @@ namespace GridOutFlags const bool draw_legend) : line_thickness(line_thickness), boundary_line_thickness(boundary_line_thickness), - margin_in_percent(margin_in_percent), + margin(margin), background(background), + azimuth_angle(azimuth_angle), + polar_angle(polar_angle), coloring(coloring), + convert_level_number_to_height(convert_level_number_to_height), label_level_number(label_level_number), label_cell_index(label_cell_index), label_material_id(label_material_id), @@ -625,6 +631,7 @@ GridOut::memory_consumption () const sizeof(eps_flags_2) + sizeof(eps_flags_3) + sizeof(xfig_flags) + + sizeof(svg_flags) + sizeof(mathgl_flags)); } @@ -1341,6 +1348,7 @@ void GridOut::write_xfig ( } + template void GridOut::write_svg(const Triangulation &tria, std::ostream &out) const { @@ -1350,23 +1358,32 @@ void GridOut::write_svg(const Triangulation &tria, std::ostream &o void GridOut::write_svg(const Triangulation<2,2> &tria, std::ostream &out) const { + + unsigned int n_materials = 0; + unsigned int n_levels = 0; + unsigned int n_subdomains = 0; + //unsigned int n_level_subdomains = 0; // TODO [CW]: CellAccessor does not seem to possess this attribute ... - const int dim = 2; + unsigned int n = 0; + + unsigned int min_level, max_level; + + // Svg files require an underlying drawing grid. + // We set the number of height elements of this + // grid to 4000 and adapt the number of width + // elements with respect to the bounding box of + // the given triangulation. + const unsigned int height = 4000; // TODO [CW]: consider other options, especially if the polar angle approaches 90° + unsigned int width; + + unsigned int margin_in_percent = 0; + if(svg_flags.margin || svg_flags.background == GridOutFlags::Svg::Background::dealii) margin_in_percent = 8.5; -//------------------ determine the bounding box of the given triangulation, -//------------------ the materials being used, the levels being used, -//------------------ the subdomains being used (, and the level subdomains being used) - unsigned int n_materials = 0; - unsigned int n_levels = 0; - unsigned int n_subdomains = 0; - //unsigned int n_level_subdomains = 0; // TODO [CW]: CellAccessor does not - // seem to possess this attribute ... - - const unsigned int height = 4000; // TODO [CW]: consider other options ... -// the following + // initial font size for cell labels + unsigned int cell_label_font_size; - unsigned int cell_label_font_size = static_cast(0.5+ (height/100.) * 2.25); // initial font size for cell labels - unsigned int font_size = static_cast(0.5+ (height/100.) * 2.); // font size for date, time, legend, and colorbar + // font size for date, time, legend, and colorbar + unsigned int font_size = static_cast(.5 + (height/100.) * 2.); // get date and time time_t time_stamp; @@ -1374,82 +1391,108 @@ void GridOut::write_svg(const Triangulation<2,2> &tria, std::ostream &out) const time_stamp = time(0); now = localtime(&time_stamp); + // vectors and variables for the perspective view + Point<3> camera_position; + Point<3> camera_direction; + Point<3> camera_horizontal; + float camera_focus; + + Point<3> point; + Point<2> projection_decomposition; + + float x_max_perspective, x_min_perspective; + float y_max_perspective, y_min_perspective; + + float x_dimension_perspective, y_dimension_perspective; + + Triangulation<2,2>::active_cell_iterator cell = tria.begin_active(), endc = tria.end(); + // auxiliary variables for the bounding box and the range of cell levels double x_min = cell->vertex(0)[0]; double x_max = cell->vertex(0)[0]; double y_min = cell->vertex(0)[1]; double y_max = cell->vertex(0)[1]; - unsigned int min_level = cell->level(); - unsigned int max_level = cell->level(); + double x_dimension, y_dimension; - // array for the materials being used + min_level = cell->level(); + max_level = cell->level(); + + // auxiliary array for the materials being used (material ids 255 max.) unsigned int materials[256]; - for (unsigned int material_index = 0; material_index < 256; material_index++) materials[material_index] = 0; + for(unsigned int material_index = 0; material_index < 256; material_index++) + materials[material_index] = 0; - // array for the levels being used + // auxiliary array for the levels being used (level number 255 max.) unsigned int levels[256]; - for (unsigned int level_index = 0; level_index < 256; level_index++) levels[level_index] = 0; + for(unsigned int level_index = 0; level_index < 256; level_index++) + levels[level_index] = 0; - // array for the subdomains being used + // auxiliary array for the subdomains being used (subdomain id 255 max.) unsigned int subdomains[256]; - for (unsigned int subdomain_index = 0; subdomain_index < 256; subdomain_index++) subdomains[subdomain_index] = 0; + for(unsigned int subdomain_index = 0; subdomain_index < 256; subdomain_index++) + subdomains[subdomain_index] = 0; + // TODO [CW]: not yet implemented ... /* - // array for the level subdomains being used + // auxiliary array for the level subdomains being used int level_subdomains[256]; for(int level_subdomain_index = 0; level_subdomain_index < 256; level_subdomain_index++) level_subdomains[level_subdomain_index] = 0; */ - // determine the bounding box of the triangulation and check the cells for material and subdomain - for (; cell != endc; ++cell) + // We use an active cell iterator to determine the + // bounding box of the given triangulation and check + // the cells for material id, level number, subdomain id + // (, and level subdomain id). + for(; cell != endc; ++cell) + { + for(unsigned int vertex_index = 0; vertex_index < 4; vertex_index++) { - for (unsigned int vertex_index = 0; vertex_index < GeometryInfo::vertices_per_cell; vertex_index++) - { - if (cell->vertex(vertex_index)[0] < x_min) x_min = cell->vertex(vertex_index)[0]; - if (cell->vertex(vertex_index)[0] > x_max) x_max = cell->vertex(vertex_index)[0]; - if (cell->vertex(vertex_index)[1] < y_min) y_min = cell->vertex(vertex_index)[1]; - if (cell->vertex(vertex_index)[1] > y_max) y_max = cell->vertex(vertex_index)[1]; - } - - if ((unsigned int)cell->level() < min_level) min_level = cell->level(); - if ((unsigned int)cell->level() > max_level) max_level = cell->level(); - - if ((unsigned int)cell->material_id()) materials[(unsigned int)cell->material_id()] = 1; - else if ((unsigned int)cell->material_id() == 0) materials[0] = 1; - - if ((int)cell->level()) levels[(unsigned int)cell->level()] = 1; - else if ((unsigned int)cell->level() == 0) levels[0] = 1; + if(cell->vertex(vertex_index)[0] < x_min) x_min = cell->vertex(vertex_index)[0]; + if(cell->vertex(vertex_index)[0] > x_max) x_max = cell->vertex(vertex_index)[0]; + + if(cell->vertex(vertex_index)[1] < y_min) y_min = cell->vertex(vertex_index)[1]; + if(cell->vertex(vertex_index)[1] > y_max) y_max = cell->vertex(vertex_index)[1]; + } - if ((unsigned int)cell->subdomain_id()) subdomains[(unsigned int)cell->subdomain_id()] = 1; - else if ((unsigned int)cell->subdomain_id() == 0) subdomains[0] = 1; + if((unsigned int)cell->level() < min_level) min_level = cell->level(); + if((unsigned int)cell->level() > max_level) max_level = cell->level(); - // if((unsigned int)(cell->level_subdomain_id())) level_subdomains[(unsigned int)cell->level_subdomain_id()] = 1; - // else if((unsigned int)(cell->level_subdomain_id()) == 0) level_subdomains[0] = 1; - } + if((unsigned int)cell->material_id()) materials[(unsigned int)cell->material_id()] = 1; + else if((unsigned int)cell->material_id() == 0) materials[0] = 1; + + if((int)cell->level()) levels[(unsigned int)cell->level()] = 1; + else if((unsigned int)cell->level() == 0) levels[0] = 1; + + if((unsigned int)cell->subdomain_id()) subdomains[(unsigned int)cell->subdomain_id()] = 1; + else if((unsigned int)cell->subdomain_id() == 0) subdomains[0] = 1; + + // if((unsigned int)(cell->level_subdomain_id())) level_subdomains[(unsigned int)cell->level_subdomain_id()] = 1; + // else if((unsigned int)(cell->level_subdomain_id()) == 0) level_subdomains[0] = 1; + } + + x_dimension = x_max - x_min; + y_dimension = y_max - y_min; - const double x_dimension = x_max - x_min; - const double y_dimension = y_max - y_min; - // count the materials being used - for (unsigned int material_index = 0; material_index < 256; material_index++) - { - if (materials[material_index]) n_materials++; - } - + for(unsigned int material_index = 0; material_index < 256; material_index++) + { + if(materials[material_index]) n_materials++; + } + // count the levels being used - for (unsigned int level_index = 0; level_index < 256; level_index++) - { - if (levels[level_index]) n_levels++; - } - + for(unsigned int level_index = 0; level_index < 256; level_index++) + { + if(levels[level_index]) n_levels++; + } + // count the subdomains being used - for (unsigned int subdomain_index = 0; subdomain_index < 256; subdomain_index++) - { - if (subdomains[subdomain_index]) n_subdomains++; - } - + for(unsigned int subdomain_index = 0; subdomain_index < 256; subdomain_index++) + { + if(subdomains[subdomain_index]) n_subdomains++; + } + /* // count the level subdomains being used for(int level_subdomain_index = 0; level_subdomain_index < 256; level_subdomain_index++) @@ -1458,522 +1501,690 @@ void GridOut::write_svg(const Triangulation<2,2> &tria, std::ostream &out) const } */ -//------------------ create the svg file with internal style sheet - const unsigned int width = static_cast(0.5+ height * (x_dimension/y_dimension)); - unsigned int additional_width = 0; + switch(svg_flags.coloring) + { + case GridOutFlags::Svg::Coloring::material_id: n = n_materials; + break; + case GridOutFlags::Svg::Coloring::level_number: n = n_levels; + break; + case GridOutFlags::Svg::Coloring::subdomain_id: n = n_subdomains; + break; + // case GridOutFlags::Svg::Coloring::level_subdomain_id: n = n_level_subdomains; + // break; + default: + break; + } - if (svg_flags.draw_legend && - (svg_flags.label_level_number || - svg_flags.label_cell_index || - svg_flags.label_material_id || - svg_flags.label_subdomain_id)) // || svg_flags.label_level_subdomain_id )) - { - additional_width = static_cast(0.5+ height * .4); // additional width for legend - } - else if (svg_flags.draw_colorbar && svg_flags.coloring) - { - additional_width = static_cast(0.5+ height * .175); // additional width for colorbar // TODO [CW]: not enough / automatic adjustment ... - } + // set the camera position to top view, targeting at the origin + camera_position[0] = 0; + camera_position[1] = 0; + camera_position[2] = 2. * std::max(x_dimension, y_dimension); -//TODO: [CW] Only output date and time optionally -//TODO: [CW] Use std::endl instead of '\n' - - out << "" - << '\n' << '\n'; + camera_direction[0] = 0; + camera_direction[1] = 0; + camera_direction[2] = -1; + + camera_horizontal[0] = 1; + camera_horizontal[1] = 0; + camera_horizontal[2] = 0; + + camera_focus = .5 * std::max(x_dimension, y_dimension); - // basic svg header and inclusion of the corresponding style sheet - out << "" - << '\n' << '\n'; + float* camera_position_temp = new float[3]; + float* camera_direction_temp = new float[3]; + float* camera_horizontal_temp = new float[3]; + + const float angle_factor = 3.14159265 / 180.; + + // (I) rotate the camera to the chosen polar angle + camera_position_temp[1] = cos(angle_factor * svg_flags.polar_angle) * camera_position[1] - sin(angle_factor * svg_flags.polar_angle) * camera_position[2]; + camera_position_temp[2] = sin(angle_factor * svg_flags.polar_angle) * camera_position[1] + cos(angle_factor * svg_flags.polar_angle) * camera_position[2]; + + camera_direction_temp[1] = cos(angle_factor * svg_flags.polar_angle) * camera_direction[1] - sin(angle_factor * svg_flags.polar_angle) * camera_direction[2]; + camera_direction_temp[2] = sin(angle_factor * svg_flags.polar_angle) * camera_direction[1] + cos(angle_factor * svg_flags.polar_angle) * camera_direction[2]; + + camera_horizontal_temp[1] = cos(angle_factor * svg_flags.polar_angle) * camera_horizontal[1] - sin(angle_factor * svg_flags.polar_angle) * camera_horizontal[2]; + camera_horizontal_temp[2] = sin(angle_factor * svg_flags.polar_angle) * camera_horizontal[1] + cos(angle_factor * svg_flags.polar_angle) * camera_horizontal[2]; + + camera_position[1] = camera_position_temp[1]; + camera_position[2] = camera_position_temp[2]; + + camera_direction[1] = camera_direction_temp[1]; + camera_direction[2] = camera_direction_temp[2]; + + camera_horizontal[1] = camera_horizontal_temp[1]; + camera_horizontal[2] = camera_horizontal_temp[2]; + + // (II) rotate the camera to the chosen azimuth angle + camera_position_temp[0] = cos(angle_factor * svg_flags.azimuth_angle) * camera_position[0] - sin(angle_factor * svg_flags.azimuth_angle) * camera_position[1]; + camera_position_temp[1] = sin(angle_factor * svg_flags.azimuth_angle) * camera_position[0] + cos(angle_factor * svg_flags.azimuth_angle) * camera_position[1]; + + camera_direction_temp[0] = cos(angle_factor * svg_flags.azimuth_angle) * camera_direction[0] - sin(angle_factor * svg_flags.azimuth_angle) * camera_direction[1]; + camera_direction_temp[1] = sin(angle_factor * svg_flags.azimuth_angle) * camera_direction[0] + cos(angle_factor * svg_flags.azimuth_angle) * camera_direction[1]; + + camera_horizontal_temp[0] = cos(angle_factor * svg_flags.azimuth_angle) * camera_horizontal[0] - sin(angle_factor * svg_flags.azimuth_angle) * camera_horizontal[1]; + camera_horizontal_temp[1] = sin(angle_factor * svg_flags.azimuth_angle) * camera_horizontal[0] + cos(angle_factor * svg_flags.azimuth_angle) * camera_horizontal[1]; + + camera_position[0] = camera_position_temp[0]; + camera_position[1] = camera_position_temp[1]; + + camera_direction[0] = camera_direction_temp[0]; + camera_direction[1] = camera_direction_temp[1]; + + camera_horizontal[0] = camera_horizontal_temp[0]; + camera_horizontal[1] = camera_horizontal_temp[1]; + + delete camera_position_temp; + delete camera_direction_temp; + delete camera_horizontal_temp; + + // translate the camera to the given triangulation + camera_position[0] = x_min + .5 * x_dimension; + camera_position[1] = y_min + .5 * y_dimension; - if (svg_flags.background == 2) - { - out << " " << '\n' - << " " << '\n' - << " " << '\n' - << " " - << '\n'; + camera_position[0] += 2. * std::max(x_dimension, y_dimension) * sin(angle_factor * svg_flags.polar_angle) * sin(angle_factor * svg_flags.azimuth_angle); + camera_position[1] -= 2. * std::max(x_dimension, y_dimension) * sin(angle_factor * svg_flags.polar_angle) * cos(angle_factor * svg_flags.azimuth_angle); + + + // determine the bounding box of the given triangulation on the projection plane of the camera viewing system + cell = tria.begin_active(); + endc = tria.end(); + + point[0] = cell->vertex(0)[0]; + point[1] = cell->vertex(0)[1]; + point[2] = 0; + + float min_level_min_vertex_distance = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); + } + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + x_max_perspective = projection_decomposition[0]; + x_min_perspective = projection_decomposition[0]; + + y_max_perspective = projection_decomposition[1]; + y_min_perspective = projection_decomposition[1]; + + for(; cell != endc; ++cell) + { + point[0] = cell->vertex(0)[0]; + point[1] = cell->vertex(0)[1]; + point[2] = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); } - + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + if(x_max_perspective < projection_decomposition[0]) x_max_perspective = projection_decomposition[0]; + if(x_min_perspective > projection_decomposition[0]) x_min_perspective = projection_decomposition[0]; + + if(y_max_perspective < projection_decomposition[1]) y_max_perspective = projection_decomposition[1]; + if(y_min_perspective > projection_decomposition[1]) y_min_perspective = projection_decomposition[1]; + + point[0] = cell->vertex(1)[0]; + point[1] = cell->vertex(1)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + if(x_max_perspective < projection_decomposition[0]) x_max_perspective = projection_decomposition[0]; + if(x_min_perspective > projection_decomposition[0]) x_min_perspective = projection_decomposition[0]; + + if(y_max_perspective < projection_decomposition[1]) y_max_perspective = projection_decomposition[1]; + if(y_min_perspective > projection_decomposition[1]) y_min_perspective = projection_decomposition[1]; + + point[0] = cell->vertex(2)[0]; + point[1] = cell->vertex(2)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + if(x_max_perspective < projection_decomposition[0]) x_max_perspective = projection_decomposition[0]; + if(x_min_perspective > projection_decomposition[0]) x_min_perspective = projection_decomposition[0]; + + if(y_max_perspective < projection_decomposition[1]) y_max_perspective = projection_decomposition[1]; + if(y_min_perspective > projection_decomposition[1]) y_min_perspective = projection_decomposition[1]; + + point[0] = cell->vertex(3)[0]; + point[1] = cell->vertex(3)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + if(x_max_perspective < projection_decomposition[0]) x_max_perspective = projection_decomposition[0]; + if(x_min_perspective > projection_decomposition[0]) x_min_perspective = projection_decomposition[0]; + + if(y_max_perspective < projection_decomposition[1]) y_max_perspective = projection_decomposition[1]; + if(y_min_perspective > projection_decomposition[1]) y_min_perspective = projection_decomposition[1]; + + if((unsigned int)cell->level() == min_level) min_level_min_vertex_distance = cell->minimum_vertex_distance(); + } + + x_dimension_perspective = x_max_perspective - x_min_perspective; + y_dimension_perspective = y_max_perspective - y_min_perspective; + + cell_label_font_size = static_cast(.5 + (height/100.) * 2.25) * 9. * (min_level_min_vertex_distance / std::min(x_dimension, y_dimension)); + + +// create the svg file with an internal style sheet + width = static_cast(.5 + height * (x_dimension_perspective / y_dimension_perspective)); + unsigned int additional_width = 0; + + if(svg_flags.draw_legend && (svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id)) // || svg_flags.label_level_subdomain_id )) + { + additional_width = static_cast(.5 + height * .4); // additional width for legend + } + else if(svg_flags.draw_colorbar && svg_flags.coloring) + { + additional_width = static_cast(.5 + height * .175); // additional width for colorbar + } + + //out << "" << '\n'; + + // basic svg header + out << "" + << '\n' << '\n'; + + + if (svg_flags.background == GridOutFlags::Svg::Background::dealii) + { + out << " " << '\n' + << " " << '\n' + << " " << '\n' + << " " << '\n'; + } + out << '\n'; - - out << "" - << '\n' << '\n'; + out << "]]>" << '\n' << '\n'; - out << " " - << '\n'; - - if (svg_flags.background == 2) - { - unsigned int x_offset = 0; - - if (svg_flags.margin_in_percent) - x_offset = static_cast(0.5+ (height/100.) * (svg_flags.margin_in_percent/2.)); - else - x_offset = static_cast(0.5+ height * .025); - - out << " (0.5+ height * .045) << '\"' - << " style=\"font-weight:100; fill:lightsteelblue; text-anchor:start; font-family:Courier; font-size:" - << static_cast(0.5+ height * .035) << "\">" - << "deal.II" << "" - << '\n'; - - out << " (0.5+ height * .045 * 3.75) << "\" y=\"" << static_cast(0.5+ height * .045) << '\"' - << " style=\"fill:lightsteelblue; text-anchor:start; font-size:" << static_cast(0.5+ font_size * .75) << "\">" - << now->tm_mday << '/' << now->tm_mon + 1 << '/' << now->tm_year + 1900 - << " - " << now->tm_hour << ':'; - - if (now->tm_min < 10) out << '0'; - - out << now->tm_min - << "" - << '\n' << '\n'; - } + // background rectangle + out << " " << '\n'; -//------------------ draw the cells - out << " " - << '\n'; + if(svg_flags.background == GridOutFlags::Svg::Background::dealii) + { + unsigned int x_offset = 0; + + if(svg_flags.margin) x_offset = static_cast(.5 + (height/100.) * (margin_in_percent/2.)); + else x_offset = static_cast(.5 + height * .025); + + out << " (.5 + height * .0525) << '\"' + << " style=\"font-weight:100; fill:lightsteelblue; text-anchor:start; font-family:Courier; font-size:" << static_cast(.5 + height * .045) << "\">" + << "deal.II" << "" << '\n'; + + // out << " (.5 + height * .045 * 4.75) << "\" y=\"" << static_cast(.5 + height * .0525) << '\"' + // << " style=\"fill:lightsteelblue; text-anchor:start; font-size:" << font_size << "\">" + // << now->tm_mday << '/' << now->tm_mon + 1 << '/' << now->tm_year + 1900 + // << " - " << now->tm_hour << ':'; + // + // if(now->tm_min < 10) out << '0'; + // + // out << now->tm_min + // << ""<< '\n' << '\n'; + } - cell = tria.begin_active(); - endc = tria.end(); +// draw the cells, starting out from the minimal level (in order to guaranty a correct perspective view) + out << " " << '\n'; - for (; cell != endc; ++cell) + for(unsigned int level_index = min_level; level_index <= max_level; level_index++) + { + Triangulation<2,2>::cell_iterator cell = tria.begin(level_index), endc = tria.end(level_index); + + for(; cell != endc; ++cell) { - // draw the current cell - out << " active()) continue; + + // draw the current cell + out << " active() && svg_flags.convert_level_number_to_height) out << 's'; + switch (svg_flags.coloring) - { - case 1: - out << (unsigned int)cell->material_id(); + { + case GridOutFlags::Svg::Coloring::material_id: out << (unsigned int)cell->material_id(); break; - case 2: - out << (unsigned int)cell->level(); + case GridOutFlags::Svg::Coloring::level_number: out << (unsigned int)cell->level(); break; - case 3: - out << (unsigned int)cell->subdomain_id(); + case GridOutFlags::Svg::Coloring::subdomain_id: out << (unsigned int)cell->subdomain_id(); break; - // case 4: out << (unsigned int)cell->level_subdomain_id(); break; - default: + // case GridOutFlags::Svg::Coloring::level_subdomain_id: out << (unsigned int)cell->level_subdomain_id(); + // break; + default: break; - } + } out << '\"'; } + + out << " d=\"M "; - out << " d=\"" - << "M " << static_cast(0.5+ ((cell->vertex(0)[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << ' ' << static_cast(0.5+ ((cell->vertex(0)[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << " L " << static_cast(0.5+ ((cell->vertex(1)[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << ' ' << static_cast(0.5+ ((cell->vertex(1)[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << " L " << static_cast(0.5+ ((cell->vertex(3)[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << ' ' << static_cast(0.5+ ((cell->vertex(3)[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << " L " << static_cast(0.5+ ((cell->vertex(2)[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << ' ' << static_cast(0.5+ ((cell->vertex(2)[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << " L " << static_cast(0.5+ ((cell->vertex(0)[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << ' ' << static_cast(0.5+ ((cell->vertex(0)[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\"/>" - << '\n'; + point[0] = cell->vertex(0)[0]; + point[1] = cell->vertex(0)[1]; + point[2] = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); + } + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) << ' ' + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + out << " L "; + + point[0] = cell->vertex(1)[0]; + point[1] = cell->vertex(1)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) << ' ' + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + out << " L "; + + point[0] = cell->vertex(3)[0]; + point[1] = cell->vertex(3)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); - // label the current cell - if (svg_flags.label_level_number || - svg_flags.label_cell_index || - svg_flags.label_material_id || - svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) + out << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) << ' ' + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + out << " L "; + + point[0] = cell->vertex(2)[0]; + point[1] = cell->vertex(2)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) << ' ' + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + out << " L "; + + point[0] = cell->vertex(0)[0]; + point[1] = cell->vertex(0)[1]; + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) << ' ' + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + out << "\"/>" << '\n'; + + // label the current cell + if(svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) { - out << " level() - << "\" x=\"" << static_cast(0.5+ ((cell->center()[0] - x_min) / x_dimension) * (width - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\" y=\"" << static_cast(0.5+ ((cell->center()[1] - y_min) / y_dimension) * (height - (height/100.) * 2. * svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) + point[0] = cell->center()[0]; + point[1] = cell->center()[1]; + point[2] = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); + } + + float distance_to_camera = sqrt(pow(point[0] - camera_position[0], 2.) + pow(point[1] - camera_position[1], 2.) + pow(point[2] - camera_position[2], 2.)); + float distance_factor = distance_to_camera / (2. * std::max(x_dimension, y_dimension)); + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << " (.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) + << "\" y=\"" << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)) + << "\" style=\"font-size:" << static_cast(.5 + cell_label_font_size * pow(.5, (float)cell->level() - 4. + 3.5 * distance_factor)) << "\">"; - if (svg_flags.label_level_number) - { - out << cell->level(); - } + if(svg_flags.label_level_number) + { + out << cell->level(); + } - if (svg_flags.label_cell_index) - { - if (svg_flags.label_level_number) out << ','; - out << cell->index(); - } + if(svg_flags.label_cell_index) + { + if(svg_flags.label_level_number) out << ','; + out << cell->index(); + } - if (svg_flags.label_material_id) - { - if (svg_flags.label_level_number - || svg_flags.label_cell_index) out << ','; - out << (int)cell->material_id(); - } + if(svg_flags.label_material_id) + { + if(svg_flags.label_level_number || svg_flags.label_cell_index) out << ','; + out << (int)cell->material_id(); + } - if (svg_flags.label_subdomain_id) - { - if (svg_flags.label_level_number - || svg_flags.label_cell_index - || svg_flags.label_material_id) out << ','; - out << cell->subdomain_id(); - } - /* - if(svg_flags.label_level_subdomain_id) - { - if(svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id) out << ','; - out << cell->level_subdomain_id(); - } - */ - out << "" - << '\n'; + if(svg_flags.label_subdomain_id) + { + if(svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id) out << ','; + out << cell->subdomain_id(); + } + + // if(svg_flags.label_level_subdomain_id) + // { + // if(svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id) out << ','; + // out << cell->level_subdomain_id(); + // } + + out << "" << '\n'; } -//------------------ draw the boundary - if (svg_flags.boundary_line_thickness) + // if the current cell lies at the boundary of the triangulation, draw the additional boundary line + if(svg_flags.boundary_line_thickness) { - for (unsigned int faceIndex = 0; faceIndex < GeometryInfo::faces_per_cell; faceIndex++) + for(unsigned int faceIndex = 0; faceIndex < 4; faceIndex++) + { + if(cell->at_boundary(faceIndex)) { - if (cell->at_boundary(faceIndex)) - { - out << " (0.5+ ((cell->face(faceIndex)->vertex(0)[0] - x_min)/x_dimension) * (width - (height/100.) *2.* svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\" y1=\"" << static_cast(0.5+ ((cell->face(faceIndex)->vertex(0)[1] - y_min)/y_dimension) * (height - (height/100.) *2.* svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\" x2=\"" << static_cast(0.5+ ((cell->face(faceIndex)->vertex(1)[0] - x_min)/x_dimension) * (width - (height/100.) *2.* svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\" y2=\"" << static_cast(0.5+ ((cell->face(faceIndex)->vertex(1)[1] - y_min)/y_dimension) * (height - (height/100.) *2.* svg_flags.margin_in_percent) + ((height/100.) * svg_flags.margin_in_percent)) - << "\"/>" - << '\n'; - } + + point[0] = cell->face(faceIndex)->vertex(0)[0]; + point[1] = cell->face(faceIndex)->vertex(0)[1]; + point[2] = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); + } + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << " (.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) + << "\" y1=\"" + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)); + + point[0] = cell->face(faceIndex)->vertex(1)[0]; + point[1] = cell->face(faceIndex)->vertex(1)[1]; + point[2] = 0; + + if(svg_flags.convert_level_number_to_height) + { + point[2] = .3 * ((float)cell->level() / (float)n_levels) * std::max(x_dimension, y_dimension); + } + + projection_decomposition = GridOut::svg_project_point(point, camera_position, camera_direction, camera_horizontal, camera_focus); + + out << "\" x2=\"" + << static_cast(.5 + ((projection_decomposition[0] - x_min_perspective) / x_dimension_perspective) * (width - (width/100.) * 2. * margin_in_percent) + ((width/100.) * margin_in_percent)) + << "\" y2=\"" + << static_cast(.5 + height - (height/100.) * margin_in_percent - ((projection_decomposition[1] - y_min_perspective) / y_dimension_perspective) * (height - (height/100.) * 2. * margin_in_percent)) + << "\"/>" << '\n'; } + } } + } } -//------------------ draw the legend - if (svg_flags.draw_legend && - (svg_flags.label_level_number || - svg_flags.label_cell_index || - svg_flags.label_material_id || - svg_flags.label_subdomain_id)) // || svg_flags.label_level_subdomain_id )) - { - out << '\n'; - out << " " - << '\n'; +// draw the legend + if(svg_flags.draw_legend) out << '\n' << " " << '\n'; - unsigned int additional_width = 0; - if (!svg_flags.margin_in_percent) additional_width = static_cast(0.5+ (height/100.) * 2.5); + unsigned int line_offset = 0; - out << " (0.5+ (height/100.) * svg_flags.margin_in_percent) - << "\" width=\"" << static_cast(0.5+ (height/100.) * (40. - svg_flags.margin_in_percent)) - << "\" height=\"" << static_cast(0.5+ height*.2) << "\"/>" - << '\n'; + additional_width = 0; + if(!svg_flags.margin) additional_width = static_cast(.5 + (height/100.) * 2.5); - unsigned int line_offset = 0; + // explanation of the cell labeling + if(svg_flags.draw_legend && (svg_flags.label_level_number || svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id)) // || svg_flags.label_level_subdomain_id )) + { + out << " (.5 + (height/100.) * margin_in_percent) + << "\" width=\"" << static_cast(.5 + (height/100.) * (40. - margin_in_percent)) << "\" height=\"" << static_cast(.5 + height * .165) << "\"/>" << '\n'; - out << " (0.5+ (height/100.) * 1.25) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size) - << "\" style=\"text-anchor:start; font-weight:bold; font-size:" << font_size - << "\">" - << "cell label" - << "" - << '\n'; + out << " (.5 + (height/100.) * 1.25) + << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size) + << "\" style=\"text-anchor:start; font-weight:bold; font-size:" << font_size + << "\">" << "cell label" + << "" << '\n'; - if (svg_flags.label_level_number) - { - out << " (0.5+ (height/100.) * 2.) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size) - << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size - << "\">" - << "level_number"; + if(svg_flags.label_level_number) + { + out << " (.5 + (height/100.) * 2.) + << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size) + << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size + << "\">" << "level_number"; - if (svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) - out << ','; + if(svg_flags.label_cell_index || svg_flags.label_material_id || svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) + out << ','; - out << "" - << '\n'; - } + out << "" << '\n'; + } - if (svg_flags.label_cell_index) - { - out << " (0.5+ (height/100.) * 2.) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size ) - << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size - << "\">" - << "cell_index"; + if(svg_flags.label_cell_index) + { + out << " (.5 + (height/100.) * 2.) + << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size ) + << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size + << "\">" + << "cell_index"; - if (svg_flags.label_material_id || svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) - out << ','; + if(svg_flags.label_material_id || svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) + out << ','; - out << "" - << '\n'; - } + out << "" << '\n'; + } - if (svg_flags.label_material_id) - { - out << " (0.5+ (height/100.) * 2.) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size ) - << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size - << "\">" - << "material_id"; + if(svg_flags.label_material_id) + { + out << " (.5 + (height/100.) * 2.) + << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size ) + << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size + << "\">" + << "material_id"; - if (svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) - out << ','; + if(svg_flags.label_subdomain_id) // || svg_flags.label_level_subdomain_id) + out << ','; - out << "" - << '\n'; - } + out << "" << '\n'; + } - if (svg_flags.label_subdomain_id) - { - out << " (0.5+ (height/100.) * 2.) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size ) - << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size - << "\">" - << "subdomain_id"; + if(svg_flags.label_subdomain_id) + { + out << " (.5 + (height/100.) * 2.) + << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size ) + << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size + << "\">" + << "subdomain_id"; - // if(svg_flags.label_level_subdomain_id) - // out << ','; + // if(svg_flags.label_level_subdomain_id) + // out << ','; - out << "" - << '\n'; - } - /* - if(svg_flags.label_level_subdomain_id) - { - out << " (0.5+ (height/100.) * 2.) - << "\" y=\"" << static_cast(0.5+ (height/100.) * svg_flags.margin_in_percent + (++line_offset) * 1.5 * font_size ) - << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size - << "\">" - << "level_subdomain_id" - << "" << '\n'; - } - */ + out << "" << '\n'; } + + // if(svg_flags.label_level_subdomain_id) + // { + // out << " (.5 + (height/100.) * 2.) + // << "\" y=\"" << static_cast(.5 + (height/100.) * margin_in_percent + (++line_offset) * 1.5 * font_size ) + // << "\" style=\"text-anchor:start; font-style:oblique; font-size:" << font_size + // << "\">" + // << "level_subdomain_id" + // << "" << '\n'; + // } + } -//------------------ draw the colorbar - if (svg_flags.draw_colorbar && svg_flags.coloring) - { - out << '\n'; + // show azimuth angle and polar angle as text below the explanation of the cell labeling + if(svg_flags.draw_legend) + { + out << " (.5 + (height/100.) * margin_in_percent + 10.75 * font_size) + << "\" style=\"text-anchor:start; font-size:" << font_size << "\">" + << "azimuth: " << svg_flags.azimuth_angle << "°, polar: " << svg_flags.polar_angle << "°" << '\n'; + } - out << " " - << '\n'; - unsigned int additional_width = 0; - if (!svg_flags.margin_in_percent) additional_width = static_cast(0.5+ (height/100.) * 2.5); +// draw the colorbar + if(svg_flags.draw_colorbar && svg_flags.coloring) + { + out << '\n' << " " << '\n'; + + out << " (.5 + (height/100.) * (margin_in_percent + 29.) - (font_size/1.25)) + << "\" style=\"text-anchor:start; font-weight:bold; font-size:" << font_size << "\">"; + + switch (svg_flags.coloring) + { + case 1: out << "material_id"; + break; + case 2: out << "level_number"; + break; + case 3: out << "subdomain_id"; + break; + // case 4: out << "level_subdomain_id"; + // break; + default: + break; + } - unsigned int n = 0; + out << "" << '\n'; - out << " (0.5+ (height/100.) * (svg_flags.margin_in_percent+27.5) - (font_size/1.25)) - << "\" style=\"text-anchor:start; font-weight:bold; font-size:" << font_size << "\">"; + unsigned int element_height = static_cast(((height/100.) * (71. - 2.*margin_in_percent)) / n); + unsigned int element_width = static_cast(.5 + (height/100.) * 2.5); - switch (svg_flags.coloring) - { - case 1: - out << "material_id", n = n_materials; + int labeling_index = 0; + + for(unsigned int index = 0; index < n; index++) + { + switch(svg_flags.coloring) + { + case GridOutFlags::Svg::Coloring::material_id: while(!materials[labeling_index]) labeling_index++; break; - case 2: - out << "level_number", n = n_levels; + case GridOutFlags::Svg::Coloring::level_number: while(!levels[labeling_index]) labeling_index++; break; - case 3: - out << "subdomain_id", n = n_subdomains; + case GridOutFlags::Svg::Coloring::subdomain_id: while(!subdomains[labeling_index]) labeling_index++; break; - // case 4: out << "level_subdomain_id", n = n_level_subdomains; break; - default: + // case GridOutFlags::Svg::Coloring::level_subdomain_id: while(!level_subdomains[labeling_index]) labeling_index++; + // break; + default: break; - } + } - out << "" - << '\n'; + out << " (.5 + (height/100.) * (margin_in_percent + 29)) + (n-index-1) * element_height + << "\" width=\"" << element_width + << "\" height=\"" << element_height + << "\"/>" << '\n'; - unsigned int element_height = static_cast(((height/100.) * (72.5 - 2.*svg_flags.margin_in_percent)) / n); - unsigned int element_width = static_cast(0.5+ (height/100.) * 2.5); + out << " (.5 + (height/100.) * (margin_in_percent + 29)) + (n-index-1 + .5) * element_height + static_cast(.5 + font_size * .35) << "\"" + << " style=\"text-anchor:start; font-size:" << static_cast(.5 + font_size); - int labeling_index = 0; + if(index == 0 || index == n-1) out << "; font-weight:bold"; - for (unsigned int index = 0; index < n; index++) - { - switch (svg_flags.coloring) - { - case 1: - while (!materials[labeling_index]) labeling_index++; - break; - case 2: - while (!levels[labeling_index]) labeling_index++; - break; - case 3: - while (!subdomains[labeling_index]) labeling_index++; - break; - // case 4: while(!level_subdomains[labeling_index]) labeling_index++; break; - default: - break; - } + out << "\">" << labeling_index; - out << " (0.5+ (height/100.) * (svg_flags.margin_in_percent + 27.5)) + (n-index-1) * element_height - << "\" width=\"" << element_width - << "\" height=\"" << element_height - << "\"/>" - << '\n'; - - unsigned int additional_width = 0; - if (!svg_flags.margin_in_percent) additional_width = static_cast(0.5+ (height/100.) * 2.5); - - out << " (0.5+ (height/100.) * (svg_flags.margin_in_percent + 27.5)) + (n-index-1 + .5) * element_height + static_cast(0.5+ font_size * .35) << "\"" - << " style=\"text-anchor:start; font-size:" << static_cast(0.5+ font_size); + if(index == n-1) out << " max"; + if(index == 0) out << " min"; - if (index == 0 || index == n-1) out << "; font-weight:bold"; + out << "" << '\n'; - out << "\">" << labeling_index; - - if (index == n-1) out << " max"; - if (index == 0) out << " min"; - - out << "" - << '\n'; - - labeling_index++; - } + labeling_index++; } + } -//------------------ finalize the svg file - out << '\n'; - out << ""; +// finalize the svg file + out << '\n' << ""; out.flush(); } @@ -2423,7 +2634,6 @@ void GridOut::write_ucd_lines (const Triangulation<2,3> &, - template void GridOut::write_ucd_faces (const Triangulation &tria, const unsigned int starting_index, @@ -2512,6 +2722,37 @@ void GridOut::write_ucd_lines (const Triangulation &tria, } +Point<2> GridOut::svg_project_point(Point<3> point, Point<3> camera_position, Point<3> camera_direction, Point<3> camera_horizontal, float camera_focus) const +{ + // ... + Point<3> camera_vertical; + camera_vertical[0] = camera_horizontal[1] * camera_direction[2] - camera_horizontal[2] * camera_direction[1]; + camera_vertical[1] = camera_horizontal[2] * camera_direction[0] - camera_horizontal[0] * camera_direction[2]; + camera_vertical[2] = camera_horizontal[0] * camera_direction[1] - camera_horizontal[1] * camera_direction[0]; + + float phi; + phi = camera_focus; + phi /= (point[0] - camera_position[0]) * camera_direction[0] + (point[1] - camera_position[1]) * camera_direction[1] + (point[2] - camera_position[2]) * camera_direction[2]; + + Point<3> projection; + projection[0] = camera_position[0] + phi * (point[0] - camera_position[0]); + projection[1] = camera_position[1] + phi * (point[1] - camera_position[1]); + projection[2] = camera_position[2] + phi * (point[2] - camera_position[2]); + + Point<2> projection_decomposition; + projection_decomposition[0] = (projection[0] - camera_position[0] - camera_focus * camera_direction[0]) * camera_horizontal[0]; + projection_decomposition[0] += (projection[1] - camera_position[1] - camera_focus * camera_direction[1]) * camera_horizontal[1]; + projection_decomposition[0] += (projection[2] - camera_position[2] - camera_focus * camera_direction[2]) * camera_horizontal[2]; + + projection_decomposition[1] = (projection[0] - camera_position[0] - camera_focus * camera_direction[0]) * camera_vertical[0]; + projection_decomposition[1] += (projection[1] - camera_position[1] - camera_focus * camera_direction[1]) * camera_vertical[1]; + projection_decomposition[1] += (projection[2] - camera_position[2] - camera_focus * camera_direction[2]) * camera_vertical[2]; + + return projection_decomposition; +} + + + namespace internal { namespace