Mercurial > octave
annotate src/DLD-FUNCTIONS/convhulln.cc @ 10154:40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Wed, 20 Jan 2010 17:33:41 -0500 |
parents | 2c279308f6ab |
children | cd550069240e |
rev | line source |
---|---|
6823 | 1 /* |
2 | |
8920 | 3 Copyright (C) 2000, 2007, 2008, 2009 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 |
44 #ifdef HAVE_QHULL | |
45 #if defined(HAVE__SNPRINTF) && !defined(HAVE_SNPRINTF) | |
46 #define snprintf _snprintf | |
47 #endif | |
48 | |
49 extern "C" { | |
50 #include <qhull/qhull_a.h> | |
51 } | |
52 | |
53 #ifdef NEED_QHULL_VERSION | |
54 char qh_version[] = "convhulln.oct 2007-07-24"; | |
55 #endif | |
56 #endif | |
57 | |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
58 DEFUN_DLD (convhulln, args, nargout, |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
59 "-*- texinfo -*-\n\ |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
60 @deftypefn {Loadable Function} {@var{h} =} convhulln (@var{p})\n\ |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
61 @deftypefnx {Loadable Function} {@var{h} =} convhulln (@var{p}, @var{opt})\n\ |
8788 | 62 @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
|
63 Return an index vector to the points of the enclosing convex hull.\n\ |
6823 | 64 The input matrix of size [n, dim] contains n points of dimension dim.\n\n\ |
65 If a second optional argument is given, it must be a string or cell array\n\ | |
66 of strings containing options for the underlying qhull command. (See\n\ | |
67 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
|
68 are \"s Qci Tcv\".\n\ |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
69 If the second output @var{V} is requested the volume of the convex hull is\n\ |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
70 calculated.\n\n\ |
6823 | 71 @seealso{convhull, delaunayn}\n\ |
72 @end deftypefn") | |
73 { | |
74 octave_value_list retval; | |
6959 | 75 |
6823 | 76 #ifdef HAVE_QHULL |
77 std::string options; | |
78 | |
6959 | 79 int nargin = args.length (); |
6823 | 80 if (nargin < 1 || nargin > 2) |
81 { | |
82 print_usage (); | |
83 return retval; | |
84 } | |
85 | |
86 if (nargin == 2) | |
87 { | |
6959 | 88 if (args (1).is_string ()) |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
89 options = args(1).string_value (); |
6959 | 90 else if (args(1).is_cell ()) |
10154
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 Cell c = args(1).cell_value (); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
93 options = ""; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
94 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
|
95 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
96 if (! c.elem(i).is_string ()) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
97 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
98 error ("convhulln: second argument must be a string or cell array of strings"); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
99 return retval; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
100 } |
6823 | 101 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
102 options = options + c.elem(i).string_value() + " "; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
103 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
104 } |
6823 | 105 else |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
106 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
107 error ("convhulln: second argument must be a string or cell array of strings"); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
108 return retval; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
109 } |
6823 | 110 } |
111 else | |
112 // turn on some consistency checks | |
113 options = "s Qci Tcv"; | |
114 | |
6959 | 115 Matrix p (args(0).matrix_value ()); |
6823 | 116 |
6959 | 117 const octave_idx_type dim = p.columns (); |
118 const octave_idx_type n = p.rows (); | |
119 p = p.transpose (); | |
6823 | 120 |
6959 | 121 double *pt_array = p.fortran_vec (); |
6823 | 122 |
123 boolT ismalloc = False; | |
124 | |
6959 | 125 OCTAVE_LOCAL_BUFFER (char, flags, 250); |
6880 | 126 |
6823 | 127 // hmm, lots of options for qhull here |
128 // QJ guarantees that the output will be triangles | |
6959 | 129 snprintf (flags, 250, "qhull QJ %s", options.c_str ()); |
6823 | 130 |
7520 | 131 if (! qh_new_qhull (dim, n, pt_array, ismalloc, flags, 0, stderr)) |
6823 | 132 { |
133 // If you want some debugging information replace the NULL | |
134 // pointer with stdout | |
135 | |
6959 | 136 vertexT *vertex, **vertexp; |
6823 | 137 facetT *facet; |
138 setT *vertices; | |
139 unsigned int nf = qh num_facets; | |
140 | |
6959 | 141 Matrix idx (nf, dim); |
6823 | 142 |
143 octave_idx_type j, i = 0; | |
144 FORALLfacets | |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
145 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
146 j = 0; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
147 if (! facet->simplicial) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
148 // should never happen with QJ |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
149 error ("convhulln: non-simplicial facet"); |
6823 | 150 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
151 if (dim == 3) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
152 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
153 vertices = qh_facet3vertex (facet); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
154 FOREACHvertex_ (vertices) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
155 idx(i, j++) = 1 + qh_pointid(vertex->point); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
156 qh_settempfree (&vertices); |
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 else |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
159 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
160 if (facet->toporient ^ qh_ORIENTclock) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
161 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
162 FOREACHvertex_ (facet->vertices) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
163 idx(i, j++) = 1 + qh_pointid(vertex->point); |
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 else |
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 FOREACHvertexreverse12_ (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 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
171 if (j < dim) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
172 // likewise but less fatal |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
173 warning ("facet %d only has %d vertices", i, j); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
174 i++; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
175 } |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
176 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
177 if (nargout == 2) |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
178 // calculate volume of convex hull |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
179 // taken from qhull src/geom2.c |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
180 { |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
181 realT area; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
182 realT dist; |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
183 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
184 FORALLfacets |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
185 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
186 if (! facet->normal) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
187 continue; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
188 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
189 if (facet->upperdelaunay && qh ATinfinity) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
190 continue; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
191 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
192 facet->f.area = area = qh_facetarea (facet); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
193 facet->isarea = True; |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
194 |
10154
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
195 if (qh DELAUNAY) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
196 { |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
197 if (facet->upperdelaunay == qh UPPERdelaunay) |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
198 qh totarea += area; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
199 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
200 else |
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 qh totarea += area; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
203 qh_distplane (qh interior_point, facet, &dist); |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
204 qh totvol += -dist * area/ qh hull_dim; |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
205 } |
40dfc0c99116
DLD-FUNCTIONS/*.cc: untabify
John W. Eaton <jwe@octave.org>
parents:
9786
diff
changeset
|
206 } |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
207 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
208 retval(1) = octave_value (qh totvol); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
209 } |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
210 |
6959 | 211 retval(0) = octave_value (idx); |
6823 | 212 } |
6959 | 213 |
214 // free long memory | |
215 qh_freeqhull (! qh_ALL); | |
216 | |
217 // free short memory and memory allocator | |
6823 | 218 int curlong, totlong; |
219 qh_memfreeshort (&curlong, &totlong); | |
220 | |
221 if (curlong || totlong) | |
6959 | 222 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
|
223 totlong, curlong); |
6823 | 224 #else |
225 error ("convhulln: not available in this version of Octave"); | |
226 #endif | |
6959 | 227 |
6823 | 228 return retval; |
229 } | |
8786
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
230 |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
231 /* |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
232 %!test |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
233 %! 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
|
234 %! [h, v] = convhulln(cube,'Pp'); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
235 %! assert (v, 1.0, 1e6*eps); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
236 %!test |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
237 %! 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
|
238 %! [h, v] = convhulln(tetrahedron,'Pp'); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
239 %! assert (v, 8/3, 1e6*eps); |
dbd428efbf56
Add calculation convex hull volume
Kai Habel <kai.habel@gmx.de>
parents:
7520
diff
changeset
|
240 */ |