Mercurial > octave
diff libinterp/dldfcn/convhulln.cc @ 20898:8da80da1ac37
maint: Use ovl() more places in the code.
author | Rik <rik@octave.org> |
---|---|
date | Mon, 14 Dec 2015 17:53:27 -0800 |
parents | 1142cf6abc0d |
children | 48b2ad5ee801 |
line wrap: on
line diff
--- a/libinterp/dldfcn/convhulln.cc Mon Dec 14 15:34:39 2015 -0800 +++ b/libinterp/dldfcn/convhulln.cc Mon Dec 14 17:53:27 2015 -0800 @@ -104,14 +104,15 @@ @seealso{convhull, delaunayn, voronoin}\n\ @end deftypefn") { - octave_value_list retval; - #if defined (HAVE_QHULL) int nargin = args.length (); + if (nargin < 1 || nargin > 2) print_usage (); + octave_value_list retval; + Matrix points (args(0).matrix_value ()); const octave_idx_type dim = points.columns (); const octave_idx_type num_points = points.rows (); @@ -159,11 +160,11 @@ #endif FILE *errfile = stderr; - if (outfile) - frame.add_fcn (close_fcn, outfile); - else + if (! outfile) error ("convhulln: unable to create temporary file for output"); + frame.add_fcn (close_fcn, outfile); + // qh_new_qhull command and points arguments are not const... std::string cmd = "qhull" + options; @@ -174,108 +175,106 @@ int exitcode = qh_new_qhull (dim, num_points, points.fortran_vec (), ismalloc, cmd_str, outfile, errfile); - if (! exitcode) - { - bool nonsimp_seen = false; + if (exitcode) + error ("convhulln: qhull failed"); - octave_idx_type nf = qh num_facets; + bool nonsimp_seen = false; + + octave_idx_type nf = qh num_facets; + + Matrix idx (nf, dim + 1); - Matrix idx (nf, dim + 1); + facetT *facet; - facetT *facet; + octave_idx_type i = 0; - octave_idx_type i = 0; + FORALLfacets + { + octave_idx_type j = 0; - FORALLfacets + if (! (nonsimp_seen || facet->simplicial || qh hull_dim == 2)) { - octave_idx_type j = 0; + nonsimp_seen = true; - if (! (nonsimp_seen || facet->simplicial || qh hull_dim == 2)) - { - nonsimp_seen = true; + if (cmd.find ("QJ") != std::string::npos) + // Should never happen with QJ. + error ("convhulln: qhull failed: option 'QJ' returned non-simplicial facet"); + } + + if (dim == 3) + { + setT *vertices = qh_facet3vertex (facet); - if (cmd.find ("QJ") != std::string::npos) - // Should never happen with QJ. - error ("convhulln: qhull failed: option 'QJ' returned non-simplicial facet"); - } + vertexT *vertex, **vertexp; + + FOREACHvertex_ (vertices) + idx(i, j++) = 1 + qh_pointid(vertex->point); - if (dim == 3) + qh_settempfree (&vertices); + } + else + { + if (facet->toporient ^ qh_ORIENTclock) { - setT *vertices = qh_facet3vertex (facet); - vertexT *vertex, **vertexp; - FOREACHvertex_ (vertices) + FOREACHvertex_ (facet->vertices) idx(i, j++) = 1 + qh_pointid(vertex->point); - - qh_settempfree (&vertices); } else { - if (facet->toporient ^ qh_ORIENTclock) - { - vertexT *vertex, **vertexp; + vertexT *vertex, **vertexp; + + FOREACHvertexreverse12_ (facet->vertices) + idx(i, j++) = 1 + qh_pointid(vertex->point); + } + } + if (j < dim) + warning ("convhulln: facet %d only has %d vertices", i, j); + + i++; + } + + // Remove extra dimension if all facets were simplicial. + + if (! nonsimp_seen) + idx.resize (nf, dim, 0.0); + + if (nargout == 2) + { + // Calculate volume of convex hull, taken from qhull src/geom2.c. - FOREACHvertex_ (facet->vertices) - idx(i, j++) = 1 + qh_pointid(vertex->point); - } - else - { - vertexT *vertex, **vertexp; + realT area; + realT dist; + + FORALLfacets + { + if (! facet->normal) + continue; + + if (facet->upperdelaunay && qh ATinfinity) + continue; - FOREACHvertexreverse12_ (facet->vertices) - idx(i, j++) = 1 + qh_pointid(vertex->point); - } + facet->f.area = area = qh_facetarea (facet); + facet->isarea = True; + + if (qh DELAUNAY) + { + if (facet->upperdelaunay == qh UPPERdelaunay) + qh totarea += area; } - if (j < dim) - warning ("convhulln: facet %d only has %d vertices", i, j); - - i++; + else + { + qh totarea += area; + qh_distplane (qh interior_point, facet, &dist); + qh totvol += -dist * area/ qh hull_dim; + } } - // Remove extra dimension if all facets were simplicial. - - if (! nonsimp_seen) - idx.resize (nf, dim, 0.0); - - if (nargout == 2) - { - // Calculate volume of convex hull, taken from qhull src/geom2.c. - - realT area; - realT dist; - - FORALLfacets - { - if (! facet->normal) - continue; - - if (facet->upperdelaunay && qh ATinfinity) - continue; + retval(1) = octave_value (qh totvol); + } - facet->f.area = area = qh_facetarea (facet); - facet->isarea = True; - - if (qh DELAUNAY) - { - if (facet->upperdelaunay == qh UPPERdelaunay) - qh totarea += area; - } - else - { - qh totarea += area; - qh_distplane (qh interior_point, facet, &dist); - qh totvol += -dist * area/ qh hull_dim; - } - } - - retval(1) = octave_value (qh totvol); - } - - retval(0) = idx; - } - else - error ("convhulln: qhull failed"); + retval(0) = idx; // Free memory from Qhull qh_freeqhull (! qh_ALL); @@ -287,11 +286,11 @@ warning ("convhulln: did not free %d bytes of long memory (%d pieces)", totlong, curlong); + return retval; + #else error ("convhulln: not available in this version of Octave"); #endif - - return retval; } /*