songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Another test...

This commit is contained in:
Laurent Rineau 2006-06-06 13:26:02 +00:00
parent 0744238abe
commit 361b29b976
29 changed files with 2496 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
.gitattributes vendored
View File

@ -1002,6 +1002,8 @@ Mesh_2/doc_tex/Mesh_2_ref/part_of_a_cluster.gif -text svneol=unset#unset
Mesh_2/test/Mesh_2/fish-and-rectangle.poly -text
Mesh_2/test/Mesh_2/fish.edg -text
Mesh_2/test/Mesh_2/fish.poly -text
Mesh_3/applications/GNUmakefile -text
Mesh_3/applications/medit_to_cgal.C -text
Min_sphere_of_spheres_d/test_extensive/stability/maple/balls-on-boundary-2.mws -text
Min_sphere_of_spheres_d/test_extensive/stability/maple/balls-on-boundary-3.mws -text
Min_sphere_of_spheres_d/web/figs/heuristic/excess.eps -text

18
.gitignore vendored
View File

@ -73,6 +73,24 @@ Mesh_2/test/Mesh_2/test_double_map
Mesh_2/test/Mesh_2/test_filtred_container
Mesh_2/test/Mesh_2/test_meshing
Mesh_3/*.tags.xml
Mesh_3/applications/*.cgal
Mesh_3/applications/*.exe
Mesh_3/applications/*.maillage
Mesh_3/applications/*.mesh
Mesh_3/applications/*.png
Mesh_3/applications/*.surface*
Mesh_3/applications/cgal_to_medit
Mesh_3/applications/depends
Mesh_3/applications/display_distribution
Mesh_3/applications/lanteri
Mesh_3/applications/lanteri_output_tet_mesh
Mesh_3/applications/medit_to_cgal
Mesh_3/applications/my_makefile
Mesh_3/applications/off_to_ghs
Mesh_3/applications/off_to_medit
Mesh_3/applications/read_mesh
Mesh_3/applications/slivers_exuder
Mesh_3/applications/stat_mesh
Mesh_3/doxygen
Mesh_3/examples/Mesh_3/chair-after.mesh
Mesh_3/examples/Mesh_3/chair-after.png

18
Mesh_3/applications/Distribution_displayer.h Normal file
View File

@ -0,0 +1,18 @@
#ifndef DISTRIBUTION_DISPLAYER_H
#define DISTRIBUTION_DISPLAYER_H
#include <CGAL/IO/Color.h>
struct Distribution_displayer
{
virtual void fill_rectangle(double x1, double y1,
double x2, double y2,
CGAL::Color c) = 0;
virtual void segment(double x1, double y1,
double x2, double y2,
CGAL::Color c) = 0;
virtual ~Distribution_displayer() {};
};
#endif

28
Mesh_3/applications/GNUmakefile Normal file
View File

@ -0,0 +1,28 @@
include makefile
#---------------------------------------------------------------------#
#
# dependencies
#
# if you want deps, create a file my_makefile that contains
# -include depends
#
#---------------------------------------------------------------------#
-include my_makefile
dep:
rm depends; $(MAKE) depends
.PHONY: dep
depends: *.C
cat /dev/null > depends
for f in *.C; do \
echo >> depends; \
echo >> depends; \
echo "$${f%.C}$(OBJ_EXT): \\" >> depends; \
$(CGAL_CXX) $(CXXFLAGS) -M -MG $$f \
| grep '\.\..*/include/CGAL' >> depends; \
done; \
test -f depends

79
Mesh_3/applications/Gd_displayer.C Normal file
View File

@ -0,0 +1,79 @@
#include "Gd_displayer.h"
Gd_displayer::Gd_displayer(int sx, int sy)
{
im = gdImageCreate(sx, sy);
gdImageColorAllocate(im, 255, 255, 255);
set_window(0., 1., 0., 1.);
}
Gd_displayer::Gd_displayer(gdImagePtr gd)
{
im = gd;
set_window(0., 1., 0., 1.);
}
Gd_displayer::~Gd_displayer()
{
gdImageDestroy(im);
}
void Gd_displayer::fill_rectangle(double x1, double y1,
double x2, double y2,
CGAL::Color c)
{
gdImageFilledRectangle(im,
x_pixel(x1), y_pixel(y2),
x_pixel(x2), y_pixel(y1),
gd_color(c));
}
void Gd_displayer::segment(double x1, double y1,
double x2, double y2,
CGAL::Color c)
{
gdImageLine(im,
x_pixel(x1), y_pixel(y1),
x_pixel(x2), y_pixel(y2),
gd_color(c));
}
void Gd_displayer::set_window(const double x_min, const double x_max,
const double y_min, const double y_max,
int x, int y, int sx, int sy)
{
xmin = x_min;
xmax = x_max;
ymin = y_min;
ymax = y_max;
xclip = x;
yclip = y;
sxclip = sx;
syclip = sy;
gdImageSetClip(im, xclip, yclip, xclip+sxclip, yclip+syclip);
xscal = sxclip / (xmax-xmin);
yscal = syclip / (ymax-ymin);
}
int Gd_displayer::x_pixel(double x) const
{
return( xclip + static_cast<int>((x-xmin)*xscal) );
}
int Gd_displayer::y_pixel(double y) const
{
return( yclip - static_cast<int>((y-ymax)*yscal) );
}
bool Gd_displayer::save_png(const char* filename)
{
FILE *pngout = fopen(filename, "wb");
if (pngout != NULL)
{
gdImagePng(im, pngout);
fclose(pngout);
return true;
}
else
return false;
}

104
Mesh_3/applications/Gd_displayer.h Normal file
View File

@ -0,0 +1,104 @@
#include <gd.h>
#include "Distribution_displayer.h"
struct Gd_displayer : public Distribution_displayer
{
// --- CONSTRUCTORS ---
/** Create a Gd_displayer using with a new GD image.
/param sx width of the image
/param sy height of the image
*/
Gd_displayer(int sx, int sy);
/** Create a Gd_displayer using an existing GD Image. */
Gd_displayer(gdImagePtr gd);
~Gd_displayer();
/** \name VIRTUAL FUNCTIONS FROM Distribution_displayer */
//@{
void fill_rectangle(double x1, double y1,
double x2, double y2,
CGAL::Color c);
void segment(double x1, double y1,
double x2, double y2,
CGAL::Color c);
//@}
/** \name FUNCTIONS SPECIFIC TO Gd_displayer */
//@{
/** Set world coordinates of the image.*/
inline
void set_window(const double x_min,
const double x_max,
const double y_min,
const double y_max)
{
set_window(x_min, x_max, y_min, y_max,
0, 0, gdImageSX(im), gdImageSY(im));
}
/** Set world coordinates and a clip zone.
\param x x coordinate (in pixel) of the upper left point of
the clip zone
\param y y coordinate (in pixel) of the upper left point of
the clip zone
\param sx width of the clip zone
\param sy height of the clip zone
*/
void set_window(const double x_min,
const double x_max,
const double y_min,
const double y_max,
int x, int y,
int sx, int sy);
/** Returns the x coordinate, in pixel, of the world-coordinate x. */
int x_pixel(double x) const;
/** Returns the y coordinate, in pixel, of the world-coordinate y. */
int y_pixel(double y) const;
/** Save the GD image to a PNG file. */
bool save_png(const char* filename);
/** Returns the GD image. */
inline gdImagePtr image()
{
return im;
}
/** Returns the index of the color c in the image palette. */
inline
int gd_color(CGAL::Color c)
{
int i = gdImageColorExact(im,
c.red(),
c.green(),
c.blue());
if( i < 0 )
return gdImageColorAllocate(im,
c.red(),
c.green(),
c.blue());
else
return i;
}
//@}
// end specific functions
private:
/** \name real dimensions
@{ */
double xmin, xmax, ymin, ymax; //@}
/** \name clip zone
@{ */
int xclip, yclip, sxclip, syclip; //@}
/** \name scaling factors
@{ */
double xscal, yscal; //@}
gdImagePtr im; /**< the GD image */
};

30
Mesh_3/applications/Qt_widget_displayer.C Normal file
View File

@ -0,0 +1,30 @@
#include <CGAL/IO/Qt_widget.h>
#include <CGAL/Simple_cartesian.h>
#include <CGAL/IO/Qt_widget.h>
#include "Qt_widget_displayer.h"
typedef CGAL::Simple_cartesian<double> Simple_kernel;
typedef Simple_kernel::Iso_rectangle_2 Rectangle_2;
typedef Simple_kernel::Segment_2 Segment_2;
typedef Simple_kernel::Point_2 Point_2;
Qt_widget_displayer::Qt_widget_displayer(CGAL::Qt_widget* w) : widget(w) {}
void Qt_widget_displayer::fill_rectangle(double x1, double y1,
double x2, double y2,
CGAL::Color color)
{
*widget << CGAL::FillColor(color)
<< color
<< Rectangle_2(Point_2(x1, y1), Point_2(x2, y2));
}
void Qt_widget_displayer::segment(double x1, double y1,
double x2, double y2,
CGAL::Color color)
{
*widget << color << Segment_2(Point_2(x1, y1),
Point_2(x2, y2));
}

20
Mesh_3/applications/Qt_widget_displayer.h Normal file
View File

@ -0,0 +1,20 @@
#include "Distribution_displayer.h"
namespace CGAL {
class Qt_widget;
}
struct Qt_widget_displayer : public Distribution_displayer
{
Qt_widget_displayer(CGAL::Qt_widget* widget);
void fill_rectangle(double x1, double y1,
double x2, double y2,
CGAL::Color c);
void segment(double x1, double y1,
double x2, double y2,
CGAL::Color c);
private:
CGAL::Qt_widget* widget;
};

193
Mesh_3/applications/README Normal file
View File

@ -0,0 +1,193 @@
This directory is a repository of tools, related to Mesh_3, but which
are not demos nor examples. Some of them (like slivers_exuder) will
perhaps be moved to demo/Mesh_3/ in the future.
Here is the list of tools:
cgal_to_medit
medit_to_cgal
display_distribution
filter_remove_tets_from_medit
lanteri
lanteri_output_tet_mesh
off_to_ghs
read_mesh
slivers_exuder
stat_mesh
Details follows...
Files formats:
-------------
[too be written]
filter_remove_tets_from_medit:
-----------------------------
This is a filter, written in Perl.
Use it like that:
./filter_remove_tets_from_medit < file.mesh > file-notets.mesh
It takes a medit file as standard input, and returns to standard
output the same medit file, without its tetrahedra.
cgal_to_medit:
-------------
Usage:
./cgal_to_medit INPUT OUPUT
INPUT must be a file of format produced by the output operator of
CGAL::Triangulation_3, with points
CGAL::Weighted_point_with_surface_index and cells
CGAL::Mesh_3::Complex_2_in_triangulation_cell_base_3.
OUTPUT will be a medit file.
medit_to_cgal:
-------------
Usage:
./medit_to_cgal INPUT OUPUT
INPUT must be a medit file.
OUTPUT will be a file of format produced by the output operator of
CGAL::Triangulation_3, with points
CGAL::Weighted_point_with_surface_index and cells
CGAL::Mesh_3::Complex_2_in_triangulation_cell_base_3.
lanteri:
--------
This tool is designed specificaly for answering to a request of
Stephane Lanteri, INRIA. It uses several files:
- utils.h contains templated functions
'display_faces_counts',
'display_facets_by_surface_indices_statistics',
'display_vertices_by_surface_indices_statistics',
'display_cells_by_volume_indices_statistics'.
- distribution.h and distribution.C declare and define functions
'compute_distribution',
'display_distribution'.
They use Distribution_displayer.h.
- Distribution_displayer.h defines a base class
'Distribution_displayer', with two pure virtual functions:
'fill_rectangle',
'segment'
- Gd_displayer.h Qt_widget_displayer.h
Gd_displayer.C Qt_widget_displayer.C
implementes two different classes derived of
'Distribution_displayer' (two implementation for displaying).
- lanteri_process_results.C defines three functions (declared in
lanteri.C):
'process_cells',
'process_volume_edges',
'process_surface_edges'.
- lanteri_utils.h defines two templated functions:
'scan_edges_and_process',
'scan_cells_and_process'.
These two functions use functions of lanteri_process_results.C.
lanteri_output_tet_mesh:
-----------------------
This tool is designed specificaly for answering to a request of
Stephane Lanteri, INRIA.
display_distribution:
--------------------
Usage:
./display_distribution [--scale x] (--tets|--noboite|--mesh) FILE
Options:
--scale SCALE SCALE is a real number, which will be the
the vertical scaling of the histogram.
--tets FILE
--noboite FILE
--mesh FILE
Read input file FILE.
FILE must a .tets file, or a .noboite file,
or a .mesh file.
This tool displays the radius-radius ratio distribution of a mesh file,
and save it in a png file named FILE-scale-SCALE.png.
off_to_ghs:
----------
Usage:
./off_to_ghs [FILE_BASENAME]
will convert FILE_BASENAME.off to FILE_BASENAME.faces
and FILE_BASENAME.points.
off_to_medit:
----------
This tool creates a medit files with 0 tetrahedra, from an .off file.
Usage:
./off_to_medit [FILE_BASENAME]
will convert FILE_BASENAME.off to FILE_BASENAME.mesh.
read_mesh:
---------
Usage:
./read_mesh FILE
FILE must be a file of format produced by the output operator of
CGAL::Triangulation_3, with points
CGAL::Weighted_point_with_surface_index and cells
CGAL::Mesh_3::Complex_2_in_triangulation_cell_base_3.
This tools read a .mesh file and displays faces counts, such as in:
Vertices: 3247
Facets on surface: 6474
Cells: 19670
Cells in volume: 18335
This file uses utils.h, and types.h.
stat_mesh:
---------
Usage:
./stat_mesh FILE
FILE must be a file of format produced by the output operator of
CGAL::Triangulation_3, with points
CGAL::Weighted_point_with_surface_index and cells
CGAL::Mesh_3::Complex_2_in_triangulation_cell_base_3.
Same as read_mesh, but displays more informations:
Vertices: 3247
Facets on surface: 6474
Cells: 19670
Cells in volume: 18335
Statistics:
(vertices)
Vertices with index #0: 0
Vertices with index #1: 450
Vertices with index #2: 563
Vertices with index #3: 705
Vertices with index #4: 863
Vertices with index #5: 666
(facets)
Hybrid facets: 0
Facets on surface #0: 0
Facets on surface #1: 896
Facets on surface #2: 1122
Facets on surface #3: 1406
Facets on surface #4: 1722
Facets on surface #5: 1328
slivers_exuder:
--------------
This tools exudes slivers from a cgal triangulation.
Usage:
./slivers_exuder INPUT OUTPUT bound
INPUT must be a file of format produced by the output operator of
CGAL::Triangulation_3, with points
CGAL::Weighted_point_with_surface_index and cells
CGAL::Mesh_3::Complex_2_in_triangulation_cell_base_3.
OUPUT is the name of the ouput file. It will be the same format
0<=bound<=1.0 is a bound on pumping process.
Uses utils.h and weighted_types.h.

46
Mesh_3/applications/cgal_to_medit.C Normal file
View File

@ -0,0 +1,46 @@
#include "weighted_types.h"
// ios
#include <iostream>
#include <fstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include <CGAL/IO/File_medit.h>
#include "utils.h"
int main(int , char** argv)
{
Tr tr;
C2T3 c2t3(tr);
std::ifstream ifs(argv[1]);
std::ofstream ofs(argv[2]);
if( !ifs || !ofs )
{
std::cerr << "Usage:\n"
<< " " << argv[0] << " INPUT OUPUT\n" <<
"\n"
" INPUT must be " << format_cgal_description <<
"\n"
" OUTPUT will be a medit file.\n"
"\n";
return EXIT_FAILURE;
}
std::cout << " Reading " << argv[1] << std::endl;
if( CGAL::Mesh_3::input_mesh(ifs, c2t3,
true,
&std::cerr) )
// if( CGAL::input_pslg_from_medit(ifs,
// c2t3,
// true, // debug
// &std::cout) ) // debug to cout
{
std::cout << " Writing " << argv[2] << std::endl;
CGAL::output_to_medit(ofs, c2t3);
return EXIT_SUCCESS;
}
else
return EXIT_FAILURE;
}

412
Mesh_3/applications/display_distribution.C Normal file
View File

@ -0,0 +1,412 @@
#include <CGAL/Simple_cartesian.h>
#include <iostream>
#include <fstream>
#include <sstream>
#include <CGAL/Mesh_3/Slivers_exuder.h>
#ifndef CGAL_USE_QT
int main()
{
std::cerr << "This tool requires Qt.\n";
exit(1);
}
#else
#include <CGAL/IO/Qt_widget.h>
#include <qpainter.h>
#include <qapplication.h>
struct K : public CGAL::Simple_cartesian<double> {};
typedef K::Point_3 Point_3;
typedef K::Tetrahedron_3 Tetrahedron_3;
typedef K::Point_2 Point_2;
typedef K::Segment_2 Segment_2;
typedef CGAL::Iso_rectangle_2<K> Rectangle_2;
/// Global variables
typedef std::map<std::string, std::string> String_options;
typedef std::map<std::string, double> Double_options;
String_options string_options;
Double_options double_options;
void init_options()
{
string_options["tets"] = "";
string_options["noboite"] = "";
string_options["mesh"] = "";
double_options["scale"] = 1;
}
void usage(char *argv0, std::string error = "")
{
if( error != "" )
std:: cerr << "Error: " << error << std::endl;
std::cerr << "Usage:\n "
<< argv0
<< " [--scale x] (--tets|--noboite|--mesh) FILE\n"
<< "Options:\n"
<< " --scale SCALE "
<< "SCALE is a real number, which will be the\n"
<< " "
<< "the vertical scaling of the histogram.\n"
<< " --tets FILE\n"
<< " --noboite FILE\n"
<< " --mesh FILE\n"
<< " "
<< "Read input file FILE.\n"
<< " "
<< "FILE must a .tets file, or a .noboite file,\n"
<< " or a .mesh file."
<< std::endl;
exit(1);
}
void parse_argv(int argc, char** argv, int extra_args = 0)
{
if (argc >=(2 + extra_args))
{
std::string arg = argv[1+extra_args];
if( arg.substr(0, 2) == "--" )
{
Double_options::iterator opt_it =
double_options.find(arg.substr(2, arg.length()-2));
if( opt_it != double_options.end() )
{
if( argc < (3 + extra_args) )
usage(argv[0],
(arg + " must be followed by a double!").c_str());
std::stringstream s;
double val;
s << argv[extra_args + 2];
s >> val;
if( !s )
usage(argv[0], ("Bad double after " + arg + "!").c_str());
opt_it->second = val;
parse_argv(argc, argv, extra_args + 2);
}
else
{
String_options::iterator opt_it =
string_options.find(arg.substr(2, arg.length()-2));
if( opt_it != string_options.end() )
{
if( argc < (3 + extra_args) )
usage(argv[0],
(arg + " must be followed by a string!").c_str());
std::string s = argv[extra_args + 2];
opt_it->second = s;
parse_argv(argc, argv, extra_args + 2);
}
else
usage(argv[0], ("Invalid option: " + arg + "!").c_str());
}
}
}
} // end parse_argv
void output_distribution_to_png(std::vector<double>& elements,
double max,
const int number_of_classes,
std::string filename)
{
const int number_of_cells = elements.size();
std::vector<int> distribution(number_of_classes);
for(int j=0;j<number_of_cells;j++)
{ // This block is (c) Pierre Alliez 2005
int index = number_of_classes-1;
// saturate highest value to last bin
if(elements[j] < max)
{
double dratio = elements[j]/max;
index = static_cast<int>(dratio*(double)number_of_classes);
}
distribution[index]++;
}
// const int max_occurrence = *std::max_element(distribution.begin(),
// distribution.end());
CGAL::Qt_widget *widget = new CGAL::Qt_widget();
qApp->setMainWidget(widget);
widget->resize(400, 400);
widget->set_window(0, 1, 0, 1, true); // x_min, x_max,
// y_min, y_max.
widget->show();
widget->lock();
*widget << CGAL::FillColor(CGAL::Color(200, 200, 200))
<< CGAL::Color(200, 200, 200)
<< Rectangle_2(Point_2(0, 0), Point_2(1,1));
if( number_of_classes == 0 ) return;
const double width = 1.0 / number_of_classes;
const double scale = double_options["scale"];
*widget << CGAL::FillColor(CGAL::BLACK);
// *widget << Segment_2(Point_2(0., 0.), Point_2(1., 0.));
for(int k=0;k<number_of_classes;k++)
if(distribution[k]>0)
{
double height = ( (distribution[k]+0.)/number_of_cells ) * scale;
*widget << CGAL::BLACK
<< Rectangle_2(Point_2(k*width, 0),
Point_2((k+1)*width, height));
}
else
*widget << CGAL::RED << Segment_2(Point_2(k*width, 0),
Point_2((k+1)*width, 0));
widget->unlock();
if( widget->get_pixmap().save( QString(filename.c_str()),
"PNG") )
std::cerr << "Distribution saved to file " << filename
<< std::endl;
else
{
std::cerr << "Error: cannot save distribution to file "
<< filename << std::endl;
exit(1);
}
qApp->exec();
}
bool failed(const char* error)
{
std::cerr << error << std::endl;
return false;
}
bool read_tets(std::vector<double>& elements, std::istream& in)
{
// read header
int nb_vertices = 0;
int nb_tets = 0;
std::string head;
in >> nb_vertices >> head;
if( !in || head != "vertices" )
return false;
in >> nb_tets >> head;
if( !in || head != "tets" )
return false;
std::vector<Point_3> points;
points.reserve(nb_vertices);
// read points
for(int i=0;i<nb_vertices;i++)
{
float x,y,z;
in >> x >> y >> z;
if( !in )
return false;
points.push_back(Point_3(x,y,z));
}
// read tets
for(int i=0;i<nb_tets;i++)
{
int dummy, i0,i1,i2,i3;
in >> dummy >> i0 >> i1 >> i2 >> i3;
if( dummy != 4 || !in )
return false;
Tetrahedron_3 t = Tetrahedron_3(points[i0],
points[i1],
points[i2],
points[i3]);
elements.push_back(radius_ratio(t));
}
return true;
}
bool read_mesh(std::vector<double>& elements, std::istream& in)
{
// Header.
std::string head;
int version;
in >> head >> version;
if( head != "MeshVersionFormatted" ||
version != 1 ||
! in)
return failed("read_mesh: bad version");
int dimension;
in >> head >> dimension;
if( head != "Dimension" ||
dimension!= 3 ||
! in)
return failed("read_mesh: bad dimension");
// Vertices
int number_of_vertices;
in >> head >> number_of_vertices;
if( head != "Vertices" || ! in )
return failed("read_mesh: bad file (missing Vertices)");
std::vector<Point_3> points;
points.reserve(number_of_vertices);
for(int i = 0; i < number_of_vertices; ++i)
{
int dummy_i;
in >> points[i] >> dummy_i;
if( !in )
return failed("read_mesh: bad file (reading of vertices)");
}
// Facets
int number_of_facets_on_surface;
in >> head >> number_of_facets_on_surface;
if( !in || head != "Triangles" )
return failed("read_mesh: bad file (missing Triangles)");
for(int i = 0; i < 4 * number_of_facets_on_surface; ++i)
{
double dummy;
in >> dummy;
}
// Tetrahedra
int number_of_cells;
in >> head >> number_of_cells;
if( !in || head != "Tetrahedra")
return failed("read_mesh: bad file (missing Tetrahedra)");
for(int i = 0; i < number_of_cells; ++i)
{
int i0, i1, i2, i3, dummy;
in >> i0 >> i1 >> i2 >> i3 >> dummy;
if( !in )
return failed("read_mesh: bad file (reading of cells)");
const Tetrahedron_3 t = Tetrahedron_3(points[i0-1],
points[i1-1],
points[i2-1],
points[i3-1]);
elements.push_back(radius_ratio(t));
}
in >> head;
if ( !in || head != "End")
return failed("read_mesh: bad file (missing End)");
else
return true;
}
bool read_noboite(std::vector<double>& elements, std::istream& in)
{
int nb_vertices = 0;
int nb_tets = 0;
int nb_input_points;
int dummy;
in >> nb_tets >> nb_vertices >> nb_input_points
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy
>> dummy;
if( !in )
return false;
std::vector<Point_3> points;
points.reserve(nb_vertices);
elements.clear();
elements.reserve(nb_tets);
// read tets
std::vector<int> tets;
tets.reserve(4 * nb_tets);
for(int i=0;i<nb_tets;i++)
{
int i0,i1,i2,i3;
in >> i0 >> i1 >> i2 >> i3;
if( !in )
return false;
tets.push_back(i0-1);
tets.push_back(i1-1);
tets.push_back(i2-1);
tets.push_back(i3-1);
}
// read points
for(int i=0;i<nb_vertices;i++)
{
double x,y,z;
in >> x >> y >> z;
if( !in )
return false;
points.push_back(Point_3(x,y,z));
}
// compute elements
for(int i = 0; i < 4 * nb_tets; i += 4)
{
const Tetrahedron_3 t = Tetrahedron_3(points[tets[i]],
points[tets[i+1]],
points[tets[i+2]],
points[tets[i+3]]);
elements.push_back(radius_ratio(t));
}
return true;
}
int main(int argc, char** argv)
{
QApplication app(argc, argv);
init_options();
parse_argv(argc, argv, 0);
bool ghs = false;
bool tets = false;
bool mesh = false;
std::string filename = "";
std::string ghs_filename = string_options["noboite"];
std::string tets_filename = string_options["tets"];
std::string mesh_filename = string_options["mesh"];
if(ghs_filename != "")
{
ghs = true;
filename = ghs_filename;
}
if(mesh_filename != "")
{
mesh = true;
filename = mesh_filename;
}
if(tets_filename != "")
{
tets = true;
filename = tets_filename;
}
std::vector<double> elements;
std::ifstream in_file(filename.c_str());
if(tets)
tets = read_tets(elements, in_file);
if(mesh)
mesh = read_mesh(elements, in_file);
if(ghs)
ghs = read_noboite(elements, in_file);
std::stringstream png_name;
png_name << filename << "-scale-" << double_options["scale"]
<< ".png";
if(tets || mesh || ghs)
output_distribution_to_png(elements, 1., 100, png_name.str());
else
usage(argv[0], "Error: cannot read file " + filename + "!");
}
#endif // CGAL_USE_QT

43
Mesh_3/applications/distribution.C Normal file
View File

@ -0,0 +1,43 @@
#include "distribution.h"
void compute_distribution(const Qualities& qualities,
const double max_quality,
Distribution& distribution)
{
const int number_of_classes = distribution.size();
const int qualities_size = qualities.size();
for(int j = 0; j < qualities_size; ++j)
{
if(qualities[j] < max_quality)
{
++distribution[static_cast<int>((qualities[j]/max_quality)
*number_of_classes)];
}
}
}
void display_distribution(Distribution_displayer* display,
const Distribution& distribution,
const double echelle)
{
const int number_of_classes = distribution.size();
if( number_of_classes == 0 ) return;
const double width = 1.0 / number_of_classes;
display->fill_rectangle(0., 0., 1., 1., CGAL::Color(200, 200, 200));
for(int k = 0; k < number_of_classes; ++k)
if(distribution[k]>0)
{
const double height = ( distribution[k]+0. ) * echelle;
display->fill_rectangle(k * width, 0,
(k+1)* width, height,
CGAL::BLACK);
}
else
display->segment(k * width, 0.,
(k+1) * width, 0.,
CGAL::RED);
}

13
Mesh_3/applications/distribution.h Normal file
View File

@ -0,0 +1,13 @@
#include <vector>
#include "Distribution_displayer.h"
typedef std::vector<double> Qualities;
typedef std::vector<int> Distribution;
void compute_distribution(const Qualities& qualities,
const double max_quality,
Distribution& distribution);
void display_distribution(Distribution_displayer* display,
const Distribution& distribution,
const double echelle);

11
Mesh_3/applications/dump_to_ghs_points Executable file
View File

@ -0,0 +1,11 @@
#!/usr/bin/perl -w
my $n = <>;
chop $n;
print "$n\n";
while(<>)
{
chop;
print; print " 0\n";
}

14
Mesh_3/applications/filter_remove_tets_from_medit Executable file
View File

@ -0,0 +1,14 @@
#!/usr/bin/perl -W -T
$stop = 0;
while( $stop == 0 )
{
$line = <>;
print $line;
if($line =~ /Tetrahedra/)
{
$stop = 1;
print "0\n"
}
}

229
Mesh_3/applications/lanteri.C Normal file
View File

@ -0,0 +1,229 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
// regular
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Regular_triangulation_euclidean_traits_3.h>
#include <CGAL/Regular_triangulation_cell_base_3.h>
// IO.h must be included before vertex and cell bases.
#include <CGAL/Mesh_3/IO.h>
// vertex and cell bases
#include <CGAL/Surface_mesh_vertex_base_3.h>
#include <CGAL/Surface_mesh_cell_base_3.h>
#include <CGAL/Volume_mesher_cell_base_3.h>
// c2t3
#include <CGAL/Complex_2_in_triangulation_3.h>
// traits class for reading meshes
#include <CGAL/Weighted_point_with_surface_index_geom_traits.h>
// radius_ratio
#include <CGAL/Mesh_3/Slivers_exuder.h>
#include <sstream>
// traits class
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_euclidean_traits_3<K> Regular_traits;
typedef CGAL::Weighted_point_with_surface_index_geom_traits<Regular_traits> My_traits;
// vertex and cell types
typedef CGAL::Surface_mesh_vertex_base_3<My_traits> Vb;
typedef CGAL::Triangulation_cell_base_3<My_traits> Cb1;
typedef CGAL::Surface_mesh_cell_base_3<My_traits, Cb1> Cb2;
typedef CGAL::Volume_mesher_cell_base_3<My_traits, Cb2> Cb3;
template <class GT, class Cb = CGAL::Triangulation_ds_cell_base_3 <> >
class Cell_with_volume_index : public Cb
{
private:
int volume;
public:
typedef typename GT::Point_3 Point;
typedef typename Cb::Triangulation_data_structure Tds;
typedef typename Tds::Vertex_handle Vertex_handle;
typedef typename Tds::Cell_handle Cell_handle;
template < class TDS3 >
struct Rebind_TDS {
typedef typename Cb::template Rebind_TDS<TDS3>::Other Cb3;
typedef Cell_with_volume_index<GT, Cb3> Other;
};
// Constructors
Cell_with_volume_index() : Cb(), volume(-1)
{
}
Cell_with_volume_index (Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Vertex_handle v3)
: Cb (v0, v1, v2, v3), volume(-1)
{
}
Cell_with_volume_index (Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Vertex_handle v3,
Cell_handle n0,
Cell_handle n1,
Cell_handle n2,
Cell_handle n3)
: Cb (v0, v1, v2, v3, n0, n1, n2, n3), volume(-1)
{
}
// volume index
int volume_index() const
{
return volume;
}
void set_volume_index(const int i)
{
volume = i;
}
#ifdef CGAL_MESH_3_IO_H
static
std::string io_signature()
{
return CGAL::Get_io_signature<Cb>()();
}
#endif
}; // end template class Cell_with_volume_index
template <class Gt, class Cb>
inline
std::istream&
operator>>(std::istream &is, Cell_with_volume_index<Gt, Cb>& c)
{
return is >> static_cast<Cb&>(c);
}
template <class Gt, class Cb>
inline
std::ostream&
operator<<(std::ostream &os, const Cell_with_volume_index<Gt, Cb>& c)
{
return os << static_cast<const Cb&>(c);
}
typedef Cell_with_volume_index<My_traits, Cb3> Cb;
// triangulation
typedef CGAL::Triangulation_data_structure_3<Vb, Cb> Tds;
typedef CGAL::Regular_triangulation_3<My_traits, Tds> Tr;
// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2T3;
// ios
#include <iostream>
#include <fstream>
#include <string>
#include <vector>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include "utils.h"
#include "distribution.h"
bool process_cells(const std::vector<Qualities>& volume_cells_quality,
const std::string filename_prefix);
bool process_volume_edges(const std::vector<Qualities>& volume_edges_length,
const std::vector<double>& length_bounds,
const std::string filename_prefix);
bool process_surface_edges(const std::vector<Qualities>& surface_edges_length,
const std::vector<double>& length_bounds,
const std::string filename_prefix);
#include "lanteri_utils.h"
typedef Tr::Vertex_handle Vertex_handle;
int main(int , char**)
{
Tr tr;
C2T3 c2t3(tr);
double r1, r2, r3, r4, r5;
std::vector<double> size_bounds(5);
std::vector<double> radii(5);
std::cout << "Input r1, r2, r3, r4, r5:" << std::endl;
std::cin >> r1 >> r2 >> r3 >> r4 >> r5;
std::cout << "Input the corresponding 5 size bounds:" << std::endl;
std::cin >> size_bounds[0]
>> size_bounds[1]
>> size_bounds[2]
>> size_bounds[3]
>> size_bounds[4];
if(!std::cin)
return EXIT_FAILURE;
std::string filename;
std::cout << "Input filename (with extension):" << std::endl;
std::cin >> filename;
std::ifstream ifs((filename+".cgal").c_str());
if( !ifs || !std::cin)
{
return EXIT_FAILURE;
}
std::cout << " Reading " << (filename+".cgal") << std::endl;
if( ! CGAL::Mesh_3::input_mesh(ifs,
c2t3,
true, // debug
&std::cerr) ) // debug to cerr
return EXIT_FAILURE;
display_faces_counts(tr, " ", &std::cout);
std::cout << "\n Statistics:\n";
std::cout << "(vertices)\n";
display_vertices_by_surface_indices_statistics(tr, " ", &std::cout);
std:: cout << "(facets)\n";
display_facets_by_surface_indices_statistics(c2t3, " ", &std::cout);
// sets volume indices
for(Tr::Finite_cells_iterator cit = tr.finite_cells_begin();
cit != tr.finite_cells_end();
++cit)
if(cit->is_in_domain())
{
const double sq_r =
squared_distance(K::Point_3(0, 0, 0),
static_cast<K::Point_3>(tr.dual(cit)));
if( sq_r < r1*r1 )
cit->set_volume_index(1);
else if( sq_r < r2*r2 )
cit->set_volume_index(2);
else if( sq_r < r3*r3 )
cit->set_volume_index(3);
else if( sq_r < r4*r4 )
cit->set_volume_index(4);
else if( sq_r < r5*r5 )
cit->set_volume_index(5);
}
std::cout << "(cells)\n";
display_cells_by_volume_indices_statistics(tr, " ", &std::cout);
std::cout << "\n"
<< " Edges lengths:\n";
if(!scan_edges_and_process(tr, size_bounds, filename, " ", &std::cout))
return EXIT_FAILURE;
if(!scan_cells_and_process(tr, filename))
return EXIT_FAILURE;
return EXIT_SUCCESS;
}

315
Mesh_3/applications/lanteri_output_tet_mesh.C Normal file
View File

@ -0,0 +1,315 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
// regular
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Regular_triangulation_euclidean_traits_3.h>
#include <CGAL/Regular_triangulation_cell_base_3.h>
// IO.h must be included before vertex and cell bases.
#include <CGAL/Mesh_3/IO.h>
// vertex and cell bases
#include <CGAL/Surface_mesh_vertex_base_3.h>
#include <CGAL/Surface_mesh_cell_base_3.h>
#include <CGAL/Volume_mesher_cell_base_3.h>
// c2t3
#include <CGAL/Complex_2_in_triangulation_3.h>
// traits class for reading meshes
#include <CGAL/Weighted_point_with_surface_index_geom_traits.h>
// radius_ratio
#include <CGAL/Mesh_3/Slivers_exuder.h>
// Point_traits
#include <CGAL/Point_traits.h>
#include <sstream>
// traits class
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_euclidean_traits_3<K> Regular_traits;
typedef CGAL::Weighted_point_with_surface_index_geom_traits<Regular_traits> My_traits;
// vertex and cell types
typedef CGAL::Surface_mesh_vertex_base_3<My_traits> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<My_traits> Cb1;
typedef CGAL::Surface_mesh_cell_base_3<My_traits, Cb1> Cb2;
typedef CGAL::Volume_mesher_cell_base_3<My_traits, Cb2> Cb3;
template <class GT, class Cb = CGAL::Triangulation_ds_cell_base_3 <> >
class Cell_with_volume_index : public Cb
{
private:
int volume;
public:
typedef typename GT::Point_3 Point;
typedef typename Cb::Triangulation_data_structure Tds;
typedef typename Tds::Vertex_handle Vertex_handle;
typedef typename Tds::Cell_handle Cell_handle;
template < class TDS3 >
struct Rebind_TDS {
typedef typename Cb::template Rebind_TDS<TDS3>::Other Cb3;
typedef Cell_with_volume_index<GT, Cb3> Other;
};
// Constructors
Cell_with_volume_index() : Cb(), volume(0)
{
}
Cell_with_volume_index (Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Vertex_handle v3)
: Cb (v0, v1, v2, v3), volume(0)
{
}
Cell_with_volume_index (Vertex_handle v0,
Vertex_handle v1,
Vertex_handle v2,
Vertex_handle v3,
Cell_handle n0,
Cell_handle n1,
Cell_handle n2,
Cell_handle n3)
: Cb (v0, v1, v2, v3, n0, n1, n2, n3), volume(0)
{
}
// volume index
int volume_index() const
{
return volume;
}
void set_volume_index(const int i)
{
volume = i;
}
#ifdef CGAL_MESH_3_IO_H
static
std::string io_signature()
{
return CGAL::Get_io_signature<Cb>()();
}
#endif
}; // end template class Cell_with_volume_index
template <class Gt, class Cb>
inline
std::istream&
operator>>(std::istream &is, Cell_with_volume_index<Gt, Cb>& c)
{
return is >> static_cast<Cb&>(c);
}
template <class Gt, class Cb>
inline
std::ostream&
operator<<(std::ostream &os, const Cell_with_volume_index<Gt, Cb>& c)
{
return os << static_cast<const Cb&>(c);
}
typedef Cell_with_volume_index<My_traits, Cb3> Cb;
// triangulation
typedef CGAL::Triangulation_data_structure_3<Vb, Cb> Tds;
typedef CGAL::Regular_triangulation_3<My_traits, Tds> Tr;
// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2T3;
// ios
#include <iostream>
#include <fstream>
#include <string>
#include <sstream>
#include <vector>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
int main(int , char**)
{
Tr tr;
C2T3 c2t3(tr);
double r1, r2, r3, r4, r5;
std::vector<double> size_bounds(5);
std::vector<double> radii(5);
std::cout << "Input r1, r2, r3, r4, r5:" << std::endl;
std::cin >> r1 >> r2 >> r3 >> r4 >> r5;
std::cout << "Input the corresponding 5 size bounds:" << std::endl;
std::cin >> size_bounds[0]
>> size_bounds[1]
>> size_bounds[2]
>> size_bounds[3]
>> size_bounds[4];
if(!std::cin)
return EXIT_FAILURE;
std::string file_prefix;
std::cout << "Input filename prefix:" << std::endl;
std::cin >> file_prefix;
// input file
std::ifstream ifs((file_prefix + ".cgal").c_str());
// ouput tets mesh file
std::ofstream ofs_maillage((file_prefix + ".maillage").c_str());
// output five surface meshes files
std::vector<std::ofstream*> ofs_surfaces(6);
// ofs_surfaces[0] will not be used.
for(int i = 1; i <= 5; ++i)
{
std::stringstream str_stream;
str_stream << file_prefix << ".surface" << i;
ofs_surfaces[i] = new std::ofstream(str_stream.str().c_str());
}
std::cout << " Reading file " << (file_prefix + ".cgal") << std::endl;
if( CGAL::Mesh_3::input_mesh(ifs, c2t3,
true,
&std::cerr) )
{
std::cout << " Writing file " << (file_prefix + ".maillage") << std::endl
<< " and files " << (file_prefix + ".surface*") << std::endl;
typedef C2T3::Triangulation Tr;
typedef Tr::Finite_cells_iterator Finite_cells_iterator;
typedef Tr::Finite_facets_iterator Finite_facets_iterator;
typedef Tr::Finite_vertices_iterator Finite_vertices_iterator;
typedef Tr::Vertex_handle Vertex_handle;
typedef Tr::Point Point;
typedef CGAL::Point_traits<Point> P_traits;
typedef P_traits::Bare_point Bare_point;
// sets volume indices
for(Tr::Finite_cells_iterator cit = tr.finite_cells_begin();
cit != tr.finite_cells_end();
++cit)
if(cit->is_in_domain())
{
const double sq_r =
squared_distance(K::Point_3(0, 0, 0),
static_cast<K::Point_3>(tr.dual(cit)));
if( sq_r < r1*r1 )
cit->set_volume_index(1); // brain
else if( sq_r < r2*r2 )
cit->set_volume_index(2); // LCR
else if( sq_r < r3*r3 )
cit->set_volume_index(3); // head skull
else if( sq_r < r4*r4 )
cit->set_volume_index(4); // skin
else if( sq_r < r5*r5 )
cit->set_volume_index(-1);// air
}
const Tr& tr = c2t3.triangulation();
// Headers
ofs_maillage << tr.number_of_vertices() << " "
<< CGAL::number_of_cells_in_domain(tr) << " "
<< CGAL::number_of_facets_on_surface_with_index(c2t3, 5)
<< std::endl;
for(int i = 1; i <= 5; ++i)
{
*ofs_surfaces[i] << tr.number_of_vertices() << " "
<< CGAL::number_of_facets_on_surface_with_index(c2t3, i)
<< std::endl;
}
// precision
ofs_maillage << std::setprecision(20);
for(int i = 1; i <= 5; ++i)
{
*ofs_surfaces[i] << std::setprecision(20);
}
// Vertices
std::map<Vertex_handle, int> V; // vertices are counter from 1
int inum = 1;
for( Finite_vertices_iterator vit = tr.finite_vertices_begin();
vit != tr.finite_vertices_end();
++vit)
{
V[vit] = inum++;
Point p = static_cast<Point>(vit->point());
ofs_maillage << p.x() << " " << p.y() << " " << p.z()
<< std::endl;
for(int i = 1; i <= 5; ++i)
{
*ofs_surfaces[i] << p.x() << " " << p.y() << " " << p.z()
<< std::endl;
}
}
// Tetrahedra
for( Finite_cells_iterator cit = tr.finite_cells_begin();
cit != tr.finite_cells_end(); ++cit)
if( cit->is_in_domain() )
{
for (int i=0; i<4; i++)
ofs_maillage << V[cit->vertex(i)] << " ";
ofs_maillage << cit->volume_index()
<< std::endl;
}
// Facets
for( Finite_facets_iterator fit = tr.finite_facets_begin();
fit != tr.finite_facets_end(); ++fit)
{
int surface_index = 0;
if (c2t3.face_status(fit->first,fit->second)
!= C2T3::NOT_IN_COMPLEX)
{
for (int i=0; i<4; i++)
if (i != (*fit).second)
{
const Vertex_handle& vh = (*fit).first->vertex(i);
surface_index = vh->point().surface_index();
if(surface_index == 5)
{
ofs_maillage << V[vh] << " ";
}
*ofs_surfaces[surface_index] << V[vh] << " ";
}
if(surface_index == 5)
ofs_maillage << "4"
<< std::endl;
*ofs_surfaces[surface_index] << surface_index
<< std::endl;
}
}
for(int i = 1; i <= 5; ++i)
if( ofs_surfaces[i]->bad() )
return EXIT_FAILURE;
if( ofs_maillage.good() )
{
for(int i = 1; i <= 5; ++i)
delete ofs_surfaces[i];
return EXIT_SUCCESS;
}
else
return EXIT_FAILURE;
}
else
return EXIT_FAILURE;
}

117
Mesh_3/applications/lanteri_process_results.C Normal file
View File

@ -0,0 +1,117 @@
#include "distribution.h"
#include "Gd_displayer.h"
#include <vector>
#include <string>
const int global_number_of_classes = 20;
bool process_aux_1(const std::vector<Qualities>& qualities,
const std::string filename,
const int number_of_classes)
{
Gd_displayer main_display(1000, 200);
main_display.set_window(-1.1, 1.1, -1.1, 1.1,
0, 0, 200, 200);
std::vector<Gd_displayer*> displays(5);
std::vector<Distribution> distributions(5);
for(int i = 0; i < 5; ++i)
{
distributions[i].resize(number_of_classes);
compute_distribution(qualities[i+1],
1.,
distributions[i]);
displays[i] = new Gd_displayer(main_display.image());
displays[i]->set_window(-.1, 1.1, -.1, 1.1,
i * 200, 0, 200, 200);
const int max = *(std::max_element(distributions[i].begin(),
distributions[i].end()));
display_distribution(displays[i],
distributions[i],
1. / max);
}
return main_display.save_png(filename.c_str());
}
bool process_aux_2(const std::vector<Qualities>& qualities,
const std::vector<double>& length_bounds,
const std::string filename,
const int number_of_classes)
{
Gd_displayer main_display(1000, 200);
main_display.set_window(-1.1, 1.1, -1.1, 1.1,
0, 0, 200, 200);
std::vector<Gd_displayer*> displays(5);
std::vector<Distribution> distributions(5);
for(int i = 0; i < 5; ++i)
{
distributions[i].resize(number_of_classes);
double max_quality = *(std::max_element(qualities[i+1].begin(),
qualities[i+1].end()));
max_quality = std::max(max_quality, length_bounds[i]);
compute_distribution(qualities[i+1],
max_quality,
distributions[i]);
displays[i] = new Gd_displayer(main_display.image());
displays[i]->set_window(-.1, 1.1, -.1, 1.1,
i * 200, 0, 200, 200);
const int max = *(std::max_element(distributions[i].begin(),
distributions[i].end()));
display_distribution(displays[i],
distributions[i],
1. / max);
double x_position_of_length_bound = length_bounds[i] / max_quality;
displays[i]->segment(x_position_of_length_bound, 0.0,
x_position_of_length_bound, -0.05,
CGAL::BLUE);
}
return main_display.save_png(filename.c_str());
}
bool process_cells(const std::vector<Qualities>& volume_cells_quality,
const std::string filename_prefix)
{
return process_aux_1(volume_cells_quality,
filename_prefix + "_cells_radius_radius_ratio.png",
global_number_of_classes);
}
bool process_volume_edges(const std::vector<Qualities>& volume_edges_length,
const std::vector<double>& length_bounds,
const std::string filename_prefix)
{
return process_aux_2(volume_edges_length,
length_bounds,
filename_prefix + "_volume_edges_lengths.png",
global_number_of_classes);
}
bool process_surface_edges(const std::vector<Qualities>& surface_edges_length,
const std::vector<double>& length_bounds,
const std::string filename_prefix)
{
return process_aux_2(surface_edges_length,
length_bounds,
filename_prefix + "_surface_edges_lengths.png",
global_number_of_classes);
}

119
Mesh_3/applications/lanteri_utils.h Normal file
View File

@ -0,0 +1,119 @@
#include <string>
template <class Tr>
bool
scan_edges_and_process(const Tr& tr,
std::vector<double> length_bounds,
std::string filename_prefix,
std::string prefix = "",
// prefix to each line output
std::ostream* out_stream =
&std::cout
// output stream
)
{
// reminder: Qualities is std::vector<double> (voir "distribution.h")
std::vector<Qualities> surface_edges_length;
std::vector<Qualities> volume_edges_length;
std::vector<double> surface_edges_length_max;
std::vector<double> volume_edges_length_max;
for(typename Tr::Finite_edges_iterator fit = tr.finite_edges_begin();
fit!=tr.finite_edges_end();
++fit)
{
const typename Tr::Vertex_handle& va = fit->first->vertex(fit->second);
const typename Tr::Vertex_handle& vb = fit->first->vertex(fit->third);
const double length =
CGAL::sqrt(CGAL::to_double(squared_distance(va->point(),
vb->point())));
const unsigned int& index_a = va->point().surface_index();
const unsigned int& index_b = vb->point().surface_index();
if( index_a != 0 && index_a == index_b ) // surface edge
{
if( surface_edges_length.size() <= index_a )
{
surface_edges_length.resize(index_a+1);
surface_edges_length_max.resize(index_a+1);
}
if( length > surface_edges_length_max[index_a] )
surface_edges_length_max[index_a] = length;
surface_edges_length[index_a].push_back(length);
}
else // volume edge
{
const int index = fit->first->volume_index();
if(index >=0)
{
const unsigned int positive_index = index;
if( volume_edges_length.size() <= positive_index )
{
volume_edges_length.resize(positive_index+1);
volume_edges_length_max.resize(positive_index+1);
}
if( length > volume_edges_length_max[positive_index] )
volume_edges_length_max[positive_index] = length;
volume_edges_length[positive_index].push_back(length);
}
}
}
const typename Qualities::size_type surface_vector_size =
surface_edges_length_max.size();
for(unsigned int i = 0; i < surface_vector_size; ++i)
{
*out_stream << prefix
<< "Maximum length for edges on surface #" << i << ": "
<< surface_edges_length_max[i] << std::endl;
}
const typename Qualities::size_type volume_vector_size =
volume_edges_length_max.size();
for(unsigned int i = 0; i < volume_vector_size; ++i)
{
*out_stream << prefix
<< "Maximum length for edges in volume #" << i << ": "
<< volume_edges_length_max[i] << std::endl;
}
return process_surface_edges(surface_edges_length,
length_bounds,
filename_prefix) &&
process_volume_edges(volume_edges_length,
length_bounds,
filename_prefix);
}
template <class Tr>
bool
scan_cells_and_process(const Tr& tr, std::string filename_prefix)
{
std::vector<Qualities> volume_cells_quality;
for(typename Tr::Finite_cells_iterator cit = tr.finite_cells_begin();
cit != tr.finite_cells_end();
++cit)
if(cit->is_in_domain())
{
const double quality =
CGAL::to_double(radius_ratio(tr.tetrahedron(cit)));
// radius ratio is in common namespace, in Slivers_exuder.h
int index = cit->volume_index();
if(index < 0)
index = 0;
const unsigned int positive_index = index;
if( volume_cells_quality.size() <= positive_index )
volume_cells_quality.resize(positive_index+1);
volume_cells_quality[positive_index].push_back(quality);
}
return process_cells(volume_cells_quality, filename_prefix);
}

116
Mesh_3/applications/makefile Normal file
View File

@ -0,0 +1,116 @@
# Created by the script cgal_create_makefile
# This is the makefile for compiling a CGAL application.
#---------------------------------------------------------------------#
# include platform specific settings
#---------------------------------------------------------------------#
# Choose the right include file from the <cgalroot>/make directory.
# CGAL_MAKEFILE = ENTER_YOUR_INCLUDE_MAKEFILE_HERE
include $(CGAL_MAKEFILE)
#---------------------------------------------------------------------#
# compiler flags
#---------------------------------------------------------------------#
CXXFLAGS = -I../../Triangulation_3/include \
-I../../Mesh_2/include \
-I../include \
-I../../Data_structure_for_queries_3/include \
-I../../Surface_mesher/include \
$(CGAL_CXXFLAGS) \
$(LONG_NAME_PROBLEM_CXXFLAGS) \
$(OTHER_CXXFLAGS)
#---------------------------------------------------------------------#
# linker flags
#---------------------------------------------------------------------#
LIBPATH = \
$(CGAL_LIBPATH)
LDFLAGS = \
$(LONG_NAME_PROBLEM_LDFLAGS) \
$(CGAL_QT_LDFLAGS)
#---------------------------------------------------------------------#
# target entries
#---------------------------------------------------------------------#
all: \
display_distribution$(EXE_EXT) \
off_to_ghs$(EXE_EXT) \
off_to_medit$(EXE_EXT) \
lanteri$(EXE_EXT) \
lanteri_output_tet_mesh$(EXE_EXT) \
read_mesh$(EXE_EXT) \
slivers_exuder$(EXE_EXT) \
cgal_to_medit$(EXE_EXT) \
medit_to_cgal$(EXE_EXT) \
stat_mesh$(EXE_EXT)
cgal_to_medit$(EXE_EXT): cgal_to_medit$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)cgal_to_medit cgal_to_medit$(OBJ_EXT) $(LDFLAGS)
medit_to_cgal$$(EXE_EXT): medit_to_cgal$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)medit_to_cgal medit_to_cgal$(OBJ_EXT) $(LDFLAGS)
display_distribution$(EXE_EXT): display_distribution$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)display_distribution display_distribution$(OBJ_EXT) $(LDFLAGS)
off_to_ghs$(EXE_EXT): off_to_ghs$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)off_to_ghs off_to_ghs$(OBJ_EXT) $(LDFLAGS)
off_to_medit$(EXE_EXT): off_to_medit$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)off_to_medit off_to_medit$(OBJ_EXT) $(LDFLAGS)
lanteri_output_tet_mesh$(EXE_EXT): lanteri_output_tet_mesh$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)lanteri_output_tet_mesh lanteri_output_tet_mesh$(OBJ_EXT) $(LDFLAGS)
lanteri$(EXE_EXT): lanteri$(OBJ_EXT) \
distribution$(OBJ_EXT) \
Qt_widget_displayer$(OBJ_EXT) \
Gd_displayer$(OBJ_EXT) \
lanteri_process_results$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)lanteri \
lanteri$(OBJ_EXT) \
distribution$(OBJ_EXT) \
Qt_widget_displayer$(OBJ_EXT) \
Gd_displayer$(OBJ_EXT) \
lanteri_process_results$(OBJ_EXT) \
$(LDFLAGS) -lgd
read_mesh$(EXE_EXT): read_mesh$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)read_mesh read_mesh$(OBJ_EXT) $(LDFLAGS)
slivers_exuder$(OBJ_EXT): ../include/CGAL/Mesh_3/Slivers_exuder.h
slivers_exuder$(EXE_EXT): slivers_exuder$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)slivers_exuder slivers_exuder$(OBJ_EXT) $(LDFLAGS)
stat_mesh$(EXE_EXT): stat_mesh$(OBJ_EXT)
$(CGAL_CXX) $(LIBPATH) $(EXE_OPT)stat_mesh stat_mesh$(OBJ_EXT) $(LDFLAGS)
clean: \
distribution.clean \
Gd_displayer.clean \
cgal_to_medit.clean \
medit_to_cgal.clean \
lanteri.clean \
lanteri_output_tet_mesh.clean \
lanteri_process_results.clean \
Qt_widget_displayer.clean \
read_mesh.clean \
slivers_exuder.clean \
stat_mesh.clean \
display_distribution.clean \
off_to_ghs.clean \
off_to_medit.clean
#---------------------------------------------------------------------#
# suffix rules
#---------------------------------------------------------------------#
.C$(OBJ_EXT):
$(CGAL_CXX) $(CXXFLAGS) $(OBJ_OPT) $<

43
Mesh_3/applications/medit_to_cgal.C Normal file
View File

@ -0,0 +1,43 @@
#include "weighted_types.h"
// ios
#include <iostream>
#include <fstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include <CGAL/IO/File_medit.h>
#include "utils.h"
int main(int , char** argv)
{
Tr tr;
C2T3 c2t3(tr);
std::ifstream ifs(argv[1]);
std::ofstream ofs(argv[2]);
if( !ifs || !ofs )
{
std::cerr << "Usage:\n"
<< " " << argv[0] << " INPUT OUPUT\n" <<
"\n"
" INPUT must be a medit file.\n"
"\n"
" OUTPUT will be " << format_cgal_description <<
"\n";
return EXIT_FAILURE;
}
std::cout << " Reading " << argv[1] << std::endl;
if( CGAL::input_from_medit(ifs,
c2t3,
true, // debug
&std::cout) ) // debug to cout
{
std::cout << " Writing " << argv[2] << std::endl;
CGAL::Mesh_3::output_mesh(ofs, c2t3);
return EXIT_SUCCESS;
}
else
return EXIT_FAILURE;
}

71
Mesh_3/applications/off_to_ghs.C Normal file
View File

@ -0,0 +1,71 @@
#include <iostream>
#include <fstream>
#include <string>
void usage(char* name)
{
std::cerr << "Usage:\n"
<< name << " [FILE_BASENAME]\n"
<< " will convert FILE_BASENAME.off to FILE_BASENAME.faces \n"
<< " and FILE_BASENAME.points.\n";
exit(1);
}
int main(int argc, char** argv)
{
if(argc < 2)
usage(argv[0]);
std::string file_base = argv[1];
std::ifstream off((file_base + ".off").c_str());
std::ofstream points((file_base + ".points").c_str());
std::ofstream faces((file_base + ".faces").c_str());
std::string dummy_s;
int number_of_points;
int number_of_faces;
int dummy_i;
off >> dummy_s; // OFF
off >> number_of_points >> number_of_faces >> dummy_i;
points << number_of_points << "\n";
faces << number_of_faces << "\n";
while( number_of_points > 0 && ! off.eof() )
{
std::string x, y, z;
off >> x >> y >> z;
if(off)
{
points << x << " " << y << " " << z << " 0\n";
--number_of_points;
}
}
if( number_of_points != 0)
{
std::cerr << "Error: bad number of points in OFF file "
<< file_base << ".off\n";
exit(2);
}
points.close();
while( ! off.eof() )
{
int i0, i1, i2, i3;
off >> i0 >> i1 >> i2 >> i3;
if(off)
{
faces << i0 << " " << (i1 + 1) << " "
<< (i2 + 1) << " " << (i3 + 1) << " 0 0 0 0\n";
--number_of_faces;
}
}
faces.close();
if( number_of_faces != 0)
{
std::cerr << "Error: bad number of faces in OFF file "
<< file_base << ".off\n";
exit(2);
}
}

80
Mesh_3/applications/off_to_medit.C Normal file
View File

@ -0,0 +1,80 @@
#include <iostream>
#include <fstream>
#include <string>
void usage(char* name)
{
std::cerr << "Usage:\n"
<< name << " [FILE_BASENAME]\n"
<< " will convert FILE_BASENAME.off to FILE_BASENAME.mesh.\n";
exit(1);
}
int main(int argc, char** argv)
{
if(argc < 2)
usage(argv[0]);
std::string file_base = argv[1];
std::ifstream off((file_base + ".off").c_str());
std::ofstream mesh((file_base + ".mesh").c_str());
std::string dummy_s;
int number_of_points;
int number_of_faces;
int dummy_i;
mesh << "MeshVersionFormatted 1\n"
<< "Dimension 3\n"
<< "Vertices\n";
off >> dummy_s; // OFF
off >> number_of_points >> number_of_faces >> dummy_i;
mesh << number_of_points << "\n";
while( number_of_points > 0 && ! off.eof() )
{
std::string x, y, z;
off >> x >> y >> z;
if(off)
{
mesh << x << " " << y << " " << z << " 1\n";
--number_of_points;
}
}
if( number_of_points != 0)
{
std::cerr << "Error: bad number of points in OFF file "
<< file_base << ".off\n";
exit(2);
}
mesh << "Triangles" << "\n"
<< number_of_faces << "\n";
while( ! off.eof() )
{
int i0, i1, i2, i3;
off >> i0 >> i1 >> i2 >> i3;
if(off)
{
mesh << (i1 + 1) << " "
<< (i2 + 1) << " "
<< (i3 + 1) << " 0\n";
--number_of_faces;
}
}
mesh << "Tetrahedra\n0\nEnd\n";
mesh.close();
if( number_of_faces != 0)
{
std::cerr << "Error: bad number of faces in OFF file "
<< file_base << ".off\n";
exit(2);
}
}

37
Mesh_3/applications/read_mesh.C Normal file
View File

@ -0,0 +1,37 @@
#include "weighted_types.h"
// ios
#include <iostream>
#include <fstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include "utils.h"
int main(int , char** argv)
{
Tr tr;
C2T3 c2t3(tr);
std::ifstream ifs(argv[1]);
if( !ifs )
{
std::cerr << "Usage:\n"
<< " " << argv[0] << " FILE\n"
<< "\n"
<< " FILE must be " << format_cgal_description
<< "\n";
return EXIT_FAILURE;
}
std::cout << " Reading " << argv[1] << std::endl;
if( CGAL::Mesh_3::input_mesh(ifs, c2t3,
true,
&std::cerr) )
{
display_faces_counts(tr, " ", &std::cout);
return EXIT_SUCCESS;
}
else
return EXIT_FAILURE;
}

58
Mesh_3/applications/slivers_exuder.C Normal file
View File

@ -0,0 +1,58 @@
#include "weighted_types.h"
// ios
#include <iostream>
#include <fstream>
#include <sstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include <CGAL/Mesh_3/Slivers_exuder.h>
#include "utils.h"
int main(int argc, char** argv)
{
Tr tr;
C2T3 c2t3(tr);
std::ifstream ifs(argv[1]);
std::ofstream ofs(argv[2]);
if( argc != 4 || !ifs || !ofs )
{
std::cerr <<
"Usage:\n" <<
" " << argv[0] << " INPUT OUTPUT bound\n" <<
"\n" <<
" INPUT must be " << format_cgal_description <<
"\n" <<
" OUPUT is the name of the ouput file. It will be in the same format.\n" <<
"\n"
" 0<=bound<=1.0 is a bound on pumping process.\n";
return EXIT_FAILURE;
}
std::stringstream str_stream(argv[3]);
double pumping_bound;
str_stream >> pumping_bound;
std::cout << " Reading " << argv[1] << std::endl;
if(! CGAL::Mesh_3::input_mesh(ifs,
c2t3,
true, // debug
&std::cerr) ) // debug to cerr
return EXIT_FAILURE;
ifs.close();
CGAL::Mesh_3::Slivers_exuder<Tr> exuder(tr);
std::cout << " Pumping" << std::endl;
exuder.init(pumping_bound);
exuder.pump_vertices(pumping_bound);
std::cout << " Writing " << argv[2] << std::endl;
CGAL::Mesh_3::output_mesh(ofs, c2t3);
return EXIT_SUCCESS;
}

51
Mesh_3/applications/stat_mesh.C Normal file
View File

@ -0,0 +1,51 @@
#include "weighted_types.h"
// ios
#include <iostream>
#include <fstream>
#include <CGAL/IO/Complex_2_in_triangulation_3_file_writer.h>
#include <CGAL/Mesh_3/IO.h>
#include "utils.h"
int main(int , char** argv)
{
Tr tr;
C2T3 c2t3(tr);
std::ifstream ifs(argv[1]);
if( !ifs )
{
std::cerr << "Usage:\n"
<< " " << argv[0] << " FILE\n"
<< "\n"
<< " FILE must be " << format_cgal_description
<< "\n";
return EXIT_FAILURE;
}
std::cout << " Reading " << argv[1] << std::endl;
if( CGAL::Mesh_3::input_mesh(ifs, c2t3,
true,
&std::cerr) )
// if( CGAL::input_pslg_from_medit(ifs,
// c2t3,
// true, // debug
// &std::cout) ) // debug to cout
{
display_faces_counts(tr, " ", &std::cout);
std::cout << "\n Statistics:\n";
std::cout << "(vertices)\n";
display_vertices_by_surface_indices_statistics(tr, " ", &std::cout);
std:: cout << "(facets)\n";
display_facets_by_surface_indices_statistics(c2t3, " ", &std::cout);
return EXIT_SUCCESS;
}
else
return EXIT_FAILURE;
}

32
Mesh_3/applications/types.h Normal file
View File

@ -0,0 +1,32 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
// vertex and cell bases
#include <CGAL/Surface_mesher_vertex_base_3.h>
#include <CGAL/Surface_mesher_cell_base_3.h>
#include <CGAL/Volume_mesher_cell_base_3.h>
// c2t3
#include <CGAL/Complex_2_in_triangulation_3.h>
// delaunay
#include <CGAL/Delaunay_triangulation_3.h>
// traits class for reading meshes
#include <CGAL/Point_with_surface_index_geom_traits.h>
// traits class
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Point_with_surface_index_geom_traits<K> My_traits;
// vertex and cell types
typedef CGAL::Surface_mesh_vertex_base_3<My_traits> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<My_traits> Cb1;
typedef CGAL::Surface_mesh_cell_base_3<My_traits, Cb1> Cb2;
typedef CGAL::Volume_mesher_cell_base_3<My_traits, Cb2> Cb;
// triangulation
typedef CGAL::Triangulation_data_structure_3<Vb, Cb> Tds;
typedef CGAL::Delaunay_triangulation_3<My_traits, Tds> Tr;
// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2T3;

153
Mesh_3/applications/utils.h Normal file
View File

@ -0,0 +1,153 @@
#include <iostream>
#include <string>
#include <vector>
template <class Tr>
void display_faces_counts(const Tr& tr,
std::string prefix = "",
// prefix to each line output
std::ostream* out_stream = &std::cout
// output stream
)
{
*out_stream << prefix << "Vertices: " << tr.number_of_vertices() << std::endl
<< prefix << "Facets on surface: "
<< CGAL::number_of_facets_on_surface(tr) << std::endl
<< prefix << "Cells: "
<< tr.number_of_cells() << std::endl
<< prefix << "Cells in volume: "
<< CGAL::number_of_cells_in_domain(tr) << std::endl;
}
template <class C2T3>
void
display_facets_by_surface_indices_statistics(C2T3& c2t3,
std::string prefix = "",
// prefix to each line output
std::ostream* out_stream =
&std::cout
// output stream
)
{
typedef typename C2T3::Triangulation Tr;
const Tr& tr = c2t3.triangulation();
std::vector<int> number_of_facets_by_surface_index;
int number_of_hybrid_facets = 0;
for(typename Tr::Finite_facets_iterator fit = tr.finite_facets_begin();
fit != tr.finite_facets_end();
++fit)
{
const typename Tr::Cell_handle& cell = fit->first;
const int index = fit->second;
if(c2t3.face_status(cell, index) != C2T3::NOT_IN_COMPLEX)
{
const typename Tr::Vertex_handle& va = cell->vertex((index+1)&3);
const typename Tr::Vertex_handle& vb = cell->vertex((index+1)&3);
const typename Tr::Vertex_handle& vc = cell->vertex((index+1)&3);
const unsigned int va_index = va->point().surface_index();
const unsigned int vb_index = vb->point().surface_index();
const unsigned int vc_index = vc->point().surface_index();
if(va_index != vb_index || va_index != vc_index )
++number_of_hybrid_facets;
else
{
if(number_of_facets_by_surface_index.size() <= va_index)
number_of_facets_by_surface_index.resize(va_index+1);
++number_of_facets_by_surface_index[va_index];
}
}
}
*out_stream << prefix << "Hybrid facets: "
<< number_of_hybrid_facets << std::endl;
std::vector<int>::size_type vector_size =
number_of_facets_by_surface_index.size();
for(unsigned int i = 0; i < vector_size; ++i)
{
*out_stream << prefix << "Facets on surface #" << i << ": "
<< number_of_facets_by_surface_index[i] << std::endl;
}
}
template <class Tr>
void
display_vertices_by_surface_indices_statistics(const Tr& tr,
std::string prefix = "",
// prefix to each line output
std::ostream* out_stream =
&std::cout
// output stream
)
{
std::vector<int> number_of_vertices_by_surface_index;
for(typename Tr::Finite_vertices_iterator vit = tr.finite_vertices_begin();
vit != tr.finite_vertices_end();
++vit)
{
const unsigned int index = vit->point().surface_index();
if(number_of_vertices_by_surface_index.size() <= index)
number_of_vertices_by_surface_index.resize(index+1);
++number_of_vertices_by_surface_index[index];
}
std::vector<int>::size_type vector_size =
number_of_vertices_by_surface_index.size();
for(unsigned int i = 0; i < vector_size; ++i)
{
*out_stream << prefix << "Vertices with index #" << i << ": "
<< number_of_vertices_by_surface_index[i] << std::endl;
}
}
template <class Tr>
void
display_cells_by_volume_indices_statistics(const Tr& tr,
std::string prefix = "",
// prefix to each line output
std::ostream* out_stream =
&std::cout
// output stream
)
{
std::vector<int> number_of_cells_by_volume_index;
int number_of_cells_with_default_volume_index = 0;
for(typename Tr::Finite_cells_iterator cit = tr.finite_cells_begin();
cit != tr.finite_cells_end();
++cit)
if(cit->is_in_domain())
{
const int index = cit->volume_index();
if(index < 0)
++number_of_cells_with_default_volume_index;
else
{
const unsigned int positive_index = index;
if(number_of_cells_by_volume_index.size() <= positive_index)
number_of_cells_by_volume_index.resize(index+1);
++number_of_cells_by_volume_index[index];
}
}
*out_stream << prefix << "Cells with default index: "
<< number_of_cells_with_default_volume_index << std::endl;
std::vector<int>::size_type vector_size =
number_of_cells_by_volume_index.size();
for(unsigned int i = 0; i < vector_size; ++i)
{
*out_stream << prefix << "Cells with index #" << i << ": "
<< number_of_cells_by_volume_index[i] << std::endl;
}
}

44
Mesh_3/applications/weighted_types.h Normal file
View File

@ -0,0 +1,44 @@
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
// regular
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Regular_triangulation_euclidean_traits_3.h>
#include <CGAL/Regular_triangulation_cell_base_3.h>
// IO.h must be included before vertex and cell bases.
#include <CGAL/Mesh_3/IO.h>
// vertex and cell bases
#include <CGAL/Surface_mesh_vertex_base_3.h>
#include <CGAL/Surface_mesh_cell_base_3.h>
#include <CGAL/Volume_mesher_cell_base_3.h>
// c2t3
#include <CGAL/Complex_2_in_triangulation_3.h>
// traits class for reading meshes
#include <CGAL/Weighted_point_with_surface_index_geom_traits.h>
// traits class
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_euclidean_traits_3<K> Regular_traits;
typedef CGAL::Weighted_point_with_surface_index_geom_traits<Regular_traits> My_traits;
// vertex and cell types
typedef CGAL::Surface_mesh_vertex_base_3<My_traits> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<My_traits> Cb1;
typedef CGAL::Surface_mesh_cell_base_3<My_traits, Cb1> Cb2;
typedef CGAL::Volume_mesher_cell_base_3<My_traits, Cb2> Cb;
// triangulation
typedef CGAL::Triangulation_data_structure_3<Vb, Cb> Tds;
typedef CGAL::Regular_triangulation_3<My_traits, Tds> Tr;
// c2t3
typedef CGAL::Complex_2_in_triangulation_3<Tr> C2T3;
const std::string format_cgal_description =
"a file of format produced by the output operator of\n"
" CGAL::Triangulation_3, with points\n"
" CGAL::Weighted_point_with_surface_index and cells\n"
" CGAL::Volume_mesher_cell_base_3.\n";