Mercurial > octave
annotate src/DLD-FUNCTIONS/convhulln.cc @ 11553:01f703952eff
Improve docstrings for functions in DLD-FUNCTIONS directory.
Use same variable names in error() strings and in documentation.
author | Rik <octave@nomad.inbox5.com> |
---|---|
date | Sun, 16 Jan 2011 22:13:23 -0800 |
parents | fd0a3ac60b0e |
children | 12df7854fa7c |
rev | line source |
---|---|
6823 | 1 /* |
2 | |
11523 | 3 Copyright (C) 2000-2011 Kai Habel |
6823 | 4 |
5 This file is part of Octave. | |
6 | |
7 Octave is free software; you can redistribute it and/or modify it | |
8 under the terms of the GNU General Public License as published by the | |
7016 | 9 Free Software Foundation; either version 3 of the License, or (at your |
10 option) any later version. | |
6823 | 11 |
12 Octave is distributed in the hope that it will be useful, but WITHOUT | |
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
15 for more details. | |
16 | |
17 You should have received a copy of the GNU General Public License | |
7016 | 18 along with Octave; see the file COPYING. If not, see |
19 <http://www.gnu.org/licenses/>. | |
6823 | 20 |
21 */ | |
22 | |
23 /* | |
24 29. July 2000 - Kai Habel: first release | |
25 2002-04-22 Paul Kienzle | |
26 * Use warning(...) function rather than writing to cerr | |
27 2006-05-01 Tom Holroyd | |
28 * add support for consistent winding in all dimensions; output is | |
29 * guaranteed to be simplicial. | |
30 */ | |
31 | |
32 #ifdef HAVE_CONFIG_H | |
33 #include <config.h> | |
34 #endif | |
9786
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
35 |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
36 #include <sstream> |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
37 |
6823 | 38 #include "Cell.h" |
9786
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
39 #include "defun-dld.h" |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
40 #include "error.h" |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
41 #include "oct-obj.h" |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
8920
diff
changeset
|
42 #include "parse.h" |
6823 | 43 |
10474
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
44 #ifdef HAVE_QHULL |
6823 | 45 extern "C" { |
46 #include <qhull/qhull_a.h> | |
47 } | |
48 | |
10474
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
49 # ifdef NEED_QHULL_VERSION |
6823 | 50 char qh_version[] = "convhulln.oct 2007-07-24"; |
10474
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
51 # endif |
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
52 #endif /* HAVE_QHULL */ |
6823 | 53 |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
54 DEFUN_DLD (convhulln, args, nargout, |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
55 "-*- texinfo -*-\n\ |
10840 | 56 @deftypefn {Loadable Function} {@var{h} =} convhulln (@var{p})\n\ |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
57 @deftypefnx {Loadable Function} {@var{h} =} convhulln (@var{p}, @var{opt})\n\ |
8788 | 58 @deftypefnx {Loadable Function} {[@var{h}, @var{v}] =} convhulln (@dots{})\n\ |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
59 Return an index vector to the points of the enclosing convex hull.\n\ |
6823 | 60 The input matrix of size [n, dim] contains n points of dimension dim.\n\n\ |
61 If a second optional argument is given, it must be a string or cell array\n\ | |
62 of strings containing options for the underlying qhull command. (See\n\ | |
63 the Qhull documentation for the available options.) The default options\n\ | |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
64 are \"s Qci Tcv\".\n\ |
11553
01f703952eff
Improve docstrings for functions in DLD-FUNCTIONS directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
65 If the second output @var{v} is requested the volume of the convex hull is\n\ |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
66 calculated.\n\n\ |
6823 | 67 @seealso{convhull, delaunayn}\n\ |
68 @end deftypefn") | |
69 { | |
70 octave_value_list retval; | |
6959 | 71 |
6823 | 72 #ifdef HAVE_QHULL |
73 std::string options; | |
74 | |
6959 | 75 int nargin = args.length (); |
6823 | 76 if (nargin < 1 || nargin > 2) |
77 { | |
78 print_usage (); | |
79 return retval; | |
80 } | |
81 | |
82 if (nargin == 2) | |
83 { | |
6959 | 84 if (args (1).is_string ()) |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
85 options = args(1).string_value (); |
6959 | 86 else if (args(1).is_cell ()) |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
87 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
88 Cell c = args(1).cell_value (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
89 options = ""; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
90 for (octave_idx_type i = 0; i < c.numel (); i++) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
91 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
92 if (! c.elem(i).is_string ()) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
93 { |
11553
01f703952eff
Improve docstrings for functions in DLD-FUNCTIONS directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
94 error ("convhulln: OPT must be a string or cell array of strings"); |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
95 return retval; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
96 } |
6823 | 97 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
98 options = options + c.elem(i).string_value() + " "; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
99 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
100 } |
6823 | 101 else |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
102 { |
11553
01f703952eff
Improve docstrings for functions in DLD-FUNCTIONS directory.
Rik <octave@nomad.inbox5.com>
parents:
11523
diff
changeset
|
103 error ("convhulln: OPT must be a string or cell array of strings"); |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
104 return retval; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
105 } |
6823 | 106 } |
107 else | |
108 // turn on some consistency checks | |
109 options = "s Qci Tcv"; | |
110 | |
6959 | 111 Matrix p (args(0).matrix_value ()); |
6823 | 112 |
6959 | 113 const octave_idx_type dim = p.columns (); |
114 const octave_idx_type n = p.rows (); | |
115 p = p.transpose (); | |
6823 | 116 |
6959 | 117 double *pt_array = p.fortran_vec (); |
6823 | 118 |
119 boolT ismalloc = False; | |
120 | |
10257
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
121 std::ostringstream buf; |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
122 |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
123 buf << "qhull QJ " << options; |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
124 |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
125 std::string buf_string = buf.str (); |
6880 | 126 |
10257
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
127 // FIXME -- we can't just pass buf_string.c_str () to qh_new_qhull |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
128 // because the argument is not declared const. Ugh. Unless |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
129 // qh_new_qhull really needs to modify this argument, someone should |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
130 // fix QHULL. |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
131 |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
132 OCTAVE_LOCAL_BUFFER (char, flags, buf_string.length () + 1); |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
133 |
cd550069240e
assume vsnprintf from gnulib; use sstream instead of snprintf
John W. Eaton <jwe@octave.org>
parents:
10154
diff
changeset
|
134 strcpy (flags, buf_string.c_str ()); |
6823 | 135 |
7520 | 136 if (! qh_new_qhull (dim, n, pt_array, ismalloc, flags, 0, stderr)) |
6823 | 137 { |
138 // If you want some debugging information replace the NULL | |
139 // pointer with stdout | |
140 | |
6959 | 141 vertexT *vertex, **vertexp; |
6823 | 142 facetT *facet; |
143 setT *vertices; | |
144 unsigned int nf = qh num_facets; | |
145 | |
6959 | 146 Matrix idx (nf, dim); |
6823 | 147 |
148 octave_idx_type j, i = 0; | |
149 FORALLfacets | |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
150 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
151 j = 0; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
152 if (! facet->simplicial) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
153 // should never happen with QJ |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
154 error ("convhulln: non-simplicial facet"); |
6823 | 155 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
156 if (dim == 3) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
157 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
158 vertices = qh_facet3vertex (facet); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
159 FOREACHvertex_ (vertices) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
160 idx(i, j++) = 1 + qh_pointid(vertex->point); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
161 qh_settempfree (&vertices); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
162 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
163 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
164 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
165 if (facet->toporient ^ qh_ORIENTclock) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
166 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
167 FOREACHvertex_ (facet->vertices) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
168 idx(i, j++) = 1 + qh_pointid(vertex->point); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
169 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
170 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
171 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
172 FOREACHvertexreverse12_ (facet->vertices) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
173 idx(i, j++) = 1 + qh_pointid(vertex->point); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
174 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
175 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
176 if (j < dim) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
177 // likewise but less fatal |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
178 warning ("facet %d only has %d vertices", i, j); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
179 i++; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
180 } |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
181 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
182 if (nargout == 2) |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
183 // calculate volume of convex hull |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
184 // taken from qhull src/geom2.c |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
185 { |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
186 realT area; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
187 realT dist; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
188 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
189 FORALLfacets |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
190 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
191 if (! facet->normal) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
192 continue; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
193 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
194 if (facet->upperdelaunay && qh ATinfinity) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
195 continue; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
196 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
197 facet->f.area = area = qh_facetarea (facet); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
198 facet->isarea = True; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
199 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
200 if (qh DELAUNAY) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
201 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
202 if (facet->upperdelaunay == qh UPPERdelaunay) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
203 qh totarea += area; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
204 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
205 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
206 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
207 qh totarea += area; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
208 qh_distplane (qh interior_point, facet, &dist); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
209 qh totvol += -dist * area/ qh hull_dim; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
210 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
211 } |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
212 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
213 retval(1) = octave_value (qh totvol); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
214 } |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
215 |
6959 | 216 retval(0) = octave_value (idx); |
6823 | 217 } |
6959 | 218 |
219 // free long memory | |
220 qh_freeqhull (! qh_ALL); | |
221 | |
222 // free short memory and memory allocator | |
6823 | 223 int curlong, totlong; |
224 qh_memfreeshort (&curlong, &totlong); | |
225 | |
226 if (curlong || totlong) | |
6959 | 227 warning ("convhulln: did not free %d bytes of long memory (%d pieces)", |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
228 totlong, curlong); |
6823 | 229 #else |
230 error ("convhulln: not available in this version of Octave"); | |
231 #endif | |
6959 | 232 |
6823 | 233 return retval; |
234 } | |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
235 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
236 /* |
10474
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
237 %!testif HAVE_QHULL |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
238 %! cube = [0 0 0;1 0 0;1 1 0;0 1 0;0 0 1;1 0 1;1 1 1;0 1 1]; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
239 %! [h, v] = convhulln(cube,'Pp'); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
240 %! assert (v, 1.0, 1e6*eps); |
10474
b2143d97c002
Neither include qhull/qhull_a.h nor run tests unless we HAVE_QHULL.
David Grundberg <davidg@cs.umu.se>
parents:
10257
diff
changeset
|
241 %!testif HAVE_QHULL |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
242 %! tetrahedron = [1 1 1;-1 -1 1;-1 1 -1;1 -1 -1]; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
243 %! [h, v] = convhulln(tetrahedron,'Pp'); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
244 %! assert (v, 8/3, 1e6*eps); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
245 */ |