Structural Bioinformatics Library
Template C++ / Python API for developping structural bioinformatics applications.
User Manual

Union_of_balls_surface_volume_3

Authors: F. Cazals and T. Dreyfus and S. Loriot

Introduction

This package provides a generic geometric data structure and algorithm for computing the surface area and the volume of the union of a family of 3D balls. Compared to other algorithms, it warranties:

  • time efficiency (...)
  • bounded errors.

Calculations are handled in a "as most as possible" exact (due to trigonometric functions, the volume and surface area cannot be exact), or using interval arithmetics.

Algorithms

Surfaces and Volumes

Having commented above on the decomposition of a union of balls into restrictions, it is sufficient to comment on the geometry of a restriction. As established in [47] , a partition of a restriction can be defined using two types of pyramids, as shown on Fig. Pyramides:

  • Type I: pyramids with flat base: the apex of such a pyramid is the center of the ball processed, and its base is (a portion of) the Voronoi facet shared with a neighboring ball in the Voronoi diagram.
  • Type II pyramids with curved base: the apex of such a pyramid is the center of the ball processed, and its based is a spherical polygon on the sphere processed.

As established using Gauss' divergence theorem, computing the volume of a pyramid is easily done from the surface area of its base [47] , be it flat or curvy:

  • For a type I pyramid, whose flat base is bounded by Vononoi line-segments and possibly circle arcs, the surface area of the base is computed using Gauss' divergence theorem.
  • For a type II pyramid, the surface area of the corresponding spherical cap is computed using the local version of the Gauss-Bonnet theorem [74] . This requires computing the geodesic curvature of circle arcs, and also angle tilts when two circle arcs bounding the spherical cap meet. Note that this latter operation involves an arc cosine which is computed using arithmetic interval from the Boost Interval Arithmetic Library. (Note that the result is stored as a pair of doubles representing an interval containing the true result.)

Thus, the primitive algorithm is the computation of the boundary of the union [4] , which is implemented in the package Union_of_balls_boundary_3.

Partitioning a restriction. The two types of pyramides

Implementation

Main Class

The main class is SBL::GT::T_Union_of_Balls_surface_volume_3< AlphaComplex3 , SphericalKernel3 , Union_of_balls_boundary_3 > : it takes as input a family of 3D balls, and stores for each 3D ball its contribution to the surface area and volume of the union of the whole family. The three template parameters are:

  • SphericalKernel3 : geometric kernel defining predicates and data structures for computations over a 3D sphere; a model for SphericalKernel3 is the CGAL data structure CGAL::Spherical_kernel_3, that is the default type; CGAL::Spherical_kernel_3 is instantiated such that it is compliant with the template parameter AlphaComplex3

An example of use of SBL::GT::T_Union_of_Balls_surface_volume_3 is shown in section Computing the Volume and Surface Area of Input 3D Balls .

Serialization

The class SBL::GT::T_Union_of_Balls_surface_volume_3 is partially serializable – see the section Partial Serialization. The partial serialization becomes available by including the file SBL/IO/Union_of_balls_surface_volume_3_oserialization.hpp in the source file where it is required.

An example of partial serialization of SBL::GT::T_Union_of_Balls_surface_volume_3 is shown in section Serializing the Contributions of 3D Balls to the Volume and Surface Area of their Union .

Module

This package provides also the module SBL::CSB::T_Union_of_balls_surface_volume_3_module< ModuleTraits , OutputArchive > that computes the surface area and the volume of the union of a family of 3D balls that are represented by:

  • a $\alpha$-complex data structure or,
  • a boundary data structure.

In the first case, the computation of the boundary has to be done before computing the surface area and volume. The two template parameters are:

  • ModuleTraits : traits class of the module defining the type Alpha_complex,
  • OutputArchive : type of the output archive for partial serialization; it can be the Boost output archives (default is boost::archive::xml_oarchive, for a XML archive), or the custom archives provided by the package Multiple_archives_serialization.

The module provides also options for switching between different computation modes (exact or interval arithmetics, number of floating points, etc...)

An example of usage of this module is shown in section Using Modules for Computing the Volume and Surface Area .

Examples

Computing the Volume and Surface Area of Input 3D Balls

This example shows how to compute the surface area and the volume of the union of a family of 3D balls:

  • (i) the 2D balls are loaded from a given file,
  • (ii) the surface are and volume is computed at the creation of the corresponding data structure,
  • (iii) and the volume and surface area are printed in the standard output.
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Fixed_alpha_shape_vertex_base_3.h>
#include <CGAL/Fixed_alpha_shape_cell_base_3.h>
#include <CGAL/Fixed_alpha_shape_3.h>
#include <CGAL/Exact_spherical_kernel_3.h>
#include <SBL/GT/Union_of_balls_surface_volume_3.hpp>
#include <iostream>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_vertex_base_3<K> Vbb;
typedef CGAL::Fixed_alpha_shape_vertex_base_3<K,Vbb> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<K> Rcb;
typedef CGAL::Fixed_alpha_shape_cell_base_3<K,Rcb> Cb;
typedef CGAL::Triangulation_data_structure_3<Vb,Cb> Tds;
typedef CGAL::Regular_triangulation_3<K,Tds> Triangulation_3;
typedef CGAL::Fixed_alpha_shape_3<Triangulation_3> Alpha_complex;
#ifndef NDEBUG
typedef CGAL::Exact_spherical_kernel_3 SK;
#else
typedef CGAL::Simple_cartesian<double> KE;
typedef CGAL::Algebraic_kernel_for_spheres_2_3<double> AK;
typedef CGAL::Spherical_kernel_3<KE, AK> SK;
#endif
typedef SBL::GT::T_Union_of_balls_surface_volume_3<Alpha_complex, SK> Union_of_balls_surface_volume_3;
int main(int argc, char *argv[])
{
if(argc < 2) return -1;
std::ifstream in(argv[1]);
std::vector<K::Sphere_3> spheres;
while(!in.eof())
{
K::Point_3 p;
K::FT r;
in >> p >> r;
if(!in.eof())
spheres.push_back(K::Sphere_3(p, r*r));
}
in.close();
Union_of_balls_surface_volume_3 surface_volume(spheres.begin(), spheres.end());
std::cout << "Total volume:" << surface_volume.volume() << std::endl;
std::cout << "Total area:" << surface_volume.area() << std::endl;
return 0;
}
Algorithm computing the surface area and the volume of the union of 3D balls.
Definition Union_of_balls_surface_volume_3.hpp:113

Serializing the Contributions of 3D Balls to the Volume and Surface Area of their Union

This example shows how to partially serialize the data structure SBL::CSB::T_Union_of_balls_surface_volume_3. The steps are the same as the example of section Computing the Volume and Surface Area of Input 3D Balls, except the third step, for which a Boost XML archive is created, and then the contribution of each 3D ball is partially serialized in this archive.

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Fixed_alpha_shape_vertex_base_3.h>
#include <CGAL/Fixed_alpha_shape_cell_base_3.h>
#include <CGAL/Fixed_alpha_shape_3.h>
#include <CGAL/Exact_spherical_kernel_3.h>
#include <SBL/GT/Union_of_balls_surface_volume_3.hpp>
//serialize the vertices of the alpha-complex, for the volume serialization
#include <SBL/IO/Alpha_complex_3_serialize.hpp>
#include <SBL/IO/Union_of_balls_surface_volume_3_oserialization.hpp>
#include <boost/serialization/list.hpp>
#include <boost/archive/xml_oarchive.hpp>
#include <iostream>
#include <fstream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Regular_triangulation_vertex_base_3<K> Vbb;
typedef CGAL::Fixed_alpha_shape_vertex_base_3<K,Vbb> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<K> Rcb;
typedef CGAL::Fixed_alpha_shape_cell_base_3<K,Rcb> Cb;
typedef CGAL::Triangulation_data_structure_3<Vb,Cb> Tds;
typedef CGAL::Regular_triangulation_3<K,Tds> Triangulation_3;
typedef CGAL::Fixed_alpha_shape_3<Triangulation_3> Alpha_complex;
#ifndef NDEBUG
typedef CGAL::Exact_spherical_kernel_3 SK;
#else
typedef CGAL::Simple_cartesian<double> KE;
typedef CGAL::Algebraic_kernel_for_spheres_2_3<double> AK;
typedef CGAL::Spherical_kernel_3<KE, AK> SK;
#endif
typedef SBL::GT::T_Union_of_balls_surface_volume_3<Alpha_complex, SK> Union_of_balls_surface_volume_3;
typedef SBL::IO::T_Union_of_balls_surface_volume_3_restriction_oserialization
<Union_of_balls_surface_volume_3> Restriction_oserialization;
int main(int argc, char *argv[])
{
if(argc < 2) return -1;
std::ifstream in(argv[1]);
std::vector<K::Sphere_3> spheres;
while(!in.eof())
{
K::Point_3 p;
K::FT r;
in >> p >> r;
if(!in.eof())
spheres.push_back(K::Sphere_3(p, r*r));
}
in.close();
Union_of_balls_surface_volume_3 surface_volume(spheres.begin(), spheres.end());
std::ofstream out("surface_volume.xml");
boost::archive::xml_oarchive* ar = new boost::archive::xml_oarchive(out);
std::list<Restriction_oserialization> restrictions;
for(Alpha_complex::Finite_vertices_iterator it = surface_volume.get_alpha_complex().finite_vertices_begin();
it != surface_volume.get_alpha_complex().finite_vertices_end(); it++)
restrictions.push_back(Restriction_oserialization(&surface_volume, it, true, 5));
*ar << boost::serialization::make_nvp("surface_volume_restrictions", restrictions);
delete ar;
out.close();
return 0;
}

Using Modules for Computing the Volume and Surface Area

This example presents the module SBL::CSB::T_Union_of_balls_surface_volume_3_module:

  • (i) the class Module_traits defines only the Alpha_complex type,
  • (iii) the $\alpha$-complex of the input balls is computed,
  • (iv) and the module is initialized and run, with statistics dumping at the end of the run.
//Load a file listing the spheres
#include <SBL/Models/Spheres_3_file_loader.hpp>
//CGAL alpha complex
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Regular_triangulation_3.h>
#include <CGAL/Fixed_alpha_shape_vertex_base_3.h>
#include <CGAL/Fixed_alpha_shape_cell_base_3.h>
#include <CGAL/Fixed_alpha_shape_3.h>
//Boundary
#include <SBL/GT/Union_of_balls_boundary_3_data_structure.hpp>
//Modules
#include <SBL/Modules/Alpha_complex_of_molecular_model_module.hpp>
#include <SBL/IO/Alpha_complex_of_molecular_model_molecular_view.hpp>
#include <SBL/Modules/Union_of_balls_surface_volume_3_module.hpp>
#include <SBL/Modules/Module_based_workflow.hpp>
//serialize the vertices of the alpha-complex, for the surface volume module
#include <SBL/IO/Alpha_complex_3_serialize.hpp>
//Traits and Modules
class Traits
{
public:
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef SBL::Models::T_Spheres_3_file_loader<CGAL::Weighted_point_3<K>, K::Point_3> Molecular_geometry_loader;
typedef CGAL::Regular_triangulation_vertex_base_3<K> Vbb;
typedef CGAL::Fixed_alpha_shape_vertex_base_3<K,Vbb> Vb;
typedef CGAL::Regular_triangulation_cell_base_3<K> Rcb;
typedef CGAL::Fixed_alpha_shape_cell_base_3<K,Rcb> Cb;
typedef CGAL::Triangulation_data_structure_3<Vb,Cb> Tds;
typedef CGAL::Regular_triangulation_3<K,Tds> Triangulation_3;
typedef CGAL::Fixed_alpha_shape_3<Triangulation_3> Alpha_complex;
typedef std::vector<CGAL::Weighted_point_3<K> > Particles_container;
};//end class Traits
//Joints
void set_alpha_complex(Alpha_complex_module* alpha_complex_module, Surface_volume_module* surface_volume_module)
{
surface_volume_module->get_alpha_complex() = &alpha_complex_module->get_alpha_complex();
}
int main(int argc, char *argv[]){
Traits::Molecular_geometry_loader* spheres_loader = workflow.add_loader<Traits::Molecular_geometry_loader>("File Loader");
SBL::Modules::Module_based_workflow::Vertex alpha_complex_vertex = workflow.register_module<Alpha_complex_module>("Alpha Complex");
workflow.set_start_module(alpha_complex_vertex);
SBL::Modules::Module_based_workflow::Vertex surface_volume_vertex = workflow.register_module<Surface_volume_module>("Surface / Volume");
workflow.make_module_flow(alpha_complex_vertex, surface_volume_vertex, &set_alpha_complex, "alpha-complex");
workflow.set_end_module(surface_volume_vertex);
workflow.parse_command_line(argc, argv);
workflow.load();
Alpha_complex_module& alpha_complex_module = (Alpha_complex_module&)workflow.get_module(alpha_complex_vertex);
alpha_complex_module.get_particles() = new Traits::Particles_container();
std::copy(spheres_loader->get_geometric_model().begin(), spheres_loader->get_geometric_model().end(),
std::back_inserter(*alpha_complex_module.get_particles()));
workflow.start();
return 0;
}
Representation of the boundary of the union of balls.
Definition Union_of_balls_boundary_3_data_structure.hpp:543
Loader for one or more txt files listing 3D spheres
Definition Spheres_3_file_loader.hpp:97
const Molecular_geometric_model & get_geometric_model(unsigned i) const
ith container of the loaded 3D spheres from the ith file (const).
Definition Spheres_3_file_loader.hpp:242
Representaiont of the workflow of an application enriched with default options.
Definition Module_based_workflow.hpp:626
Module building the alpha-complex of an input set of particles. Module building the alpha-complex of ...
Definition Alpha_complex_of_molecular_model_module.hpp:120
Particles_container *& get_particles(void)
Reference to the input particles container.
Definition Alpha_complex_of_molecular_model_module.hpp:328
const Alpha_complex & get_alpha_complex(void) const
Const reference to the output -complex.
Definition Alpha_complex_of_molecular_model_module.hpp:344
Workflow_graph_traits::vertex_descriptor Vertex
Representation of a verttex in the workflow that could be a module, or the Start or End vertices.
Definition Module_based_workflow.hpp:213
void parse_command_line(int argc, char *argv[])
Parses the command line options and store the values of the options: the check_options method is call...
Definition Module_based_workflow.hpp:2016
void make_module_flow(Vertex u, Vertex v, const Transformer &join, const std::string &name="")
Connect the modules represented by the vertices u (source output) and v (target input); if any output...
Definition Module_based_workflow.hpp:1907
void set_start_module(Vertex u)
Modules connected to the start vertex will be first executed; in case of multiple vertices,...
Definition Module_based_workflow.hpp:1890
Vertex register_module(const std::string &name="")
Registers a module with an input: note that modules are executed in a FIFO manner.
Definition Module_based_workflow.hpp:1754
void set_end_module(Vertex u)
The results of all modules connected to the end vertex are reported.
Definition Module_based_workflow.hpp:1897
Loader * add_loader(const std::string &name="")
Registers a loader into the automaton: note that loaders are executed before the modules.
Definition Module_based_workflow.hpp:1727
void start(void)
Runs the application by pushing into the stack the Start vertex.
Definition Module_based_workflow.hpp:2154
void load(std::size_t i)
Load the data using the loaders.
Definition Module_based_workflow.hpp:2127
Module computing the surface area and the volume of the union of an input set of 3D balls....
Definition Union_of_balls_surface_volume_3_module.hpp:175
Alpha_complex *& get_alpha_complex(void)
Reference to a pointer over the input -complex.
Definition Union_of_balls_surface_volume_3_module.hpp:479