Mercurial > octave
annotate src/DLD-FUNCTIONS/__voronoi__.cc @ 9786:2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
author | John W. Eaton <jwe@octave.org> |
---|---|
date | Mon, 09 Nov 2009 12:12:57 -0500 |
parents | 0631d397fbe0 |
children | 40dfc0c99116 |
rev | line source |
---|---|
6823 | 1 /* |
2 | |
8920 | 3 Copyright (C) 2000, 2007, 2008 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 20. Augiust 2000 - Kai Habel: first release | |
25 */ | |
26 | |
27 /* | |
28 2003-12-14 Rafael Laboissiere <rafael@laboissiere.net> | |
29 Added optional second argument to pass options to the underlying | |
30 qhull command | |
31 */ | |
32 | |
33 #ifdef HAVE_CONFIG_H | |
34 #include <config.h> | |
35 #endif | |
9786
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
36 |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
37 #include <cstdio> |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
38 |
9003
0631d397fbe0
replace lo_ieee_isnan by xisnan, add missing includes
Jaroslav Hajek <highegg@gmail.com>
parents:
8920
diff
changeset
|
39 #include "lo-ieee.h" |
9786
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
40 |
6823 | 41 #include "Cell.h" |
9786
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
42 #include "defun-dld.h" |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
43 #include "error.h" |
2c279308f6ab
fix includes in some src/DLD-FUNCTIONS files
John W. Eaton <jwe@octave.org>
parents:
9003
diff
changeset
|
44 #include "oct-obj.h" |
6823 | 45 |
46 #ifdef HAVE_QHULL | |
47 extern "C" { | |
48 #include <qhull/qhull_a.h> | |
49 } | |
50 | |
51 #ifdef NEED_QHULL_VERSION | |
52 char qh_version[] = "__voronoi__.oct 2007-07-24"; | |
53 #endif | |
54 #endif | |
55 | |
56 DEFUN_DLD (__voronoi__, args, , | |
57 "-*- texinfo -*-\n\ | |
58 @deftypefn {Loadable Function} {@var{tri} =} __voronoi__ (@var{point})\n\ | |
59 @deftypefnx {Loadable Function} {@var{tri} =} __voronoi__ (@var{point}, @var{options})\n\ | |
6945 | 60 Undocumented internal function.\n\ |
6823 | 61 @end deftypefn") |
62 { | |
63 octave_value_list retval; | |
6945 | 64 |
6823 | 65 #ifdef HAVE_QHULL |
6945 | 66 |
6823 | 67 retval(0) = 0.0; |
68 | |
6945 | 69 int nargin = args.length (); |
6823 | 70 if (nargin < 1 || nargin > 2) |
71 { | |
72 print_usage (); | |
73 return retval; | |
74 } | |
75 | |
6880 | 76 const char *options; |
77 | |
6823 | 78 if (nargin == 2) |
79 { | |
6945 | 80 if (! args (1).is_string ()) |
6823 | 81 { |
82 error ("__voronoi__: second argument must be a string"); | |
83 return retval; | |
84 } | |
6945 | 85 |
86 options = args (1).string_value().c_str (); | |
6823 | 87 } |
88 else | |
89 options = ""; | |
90 | |
6945 | 91 Matrix p (args(0).matrix_value ()); |
6823 | 92 |
6945 | 93 const octave_idx_type dim = p.columns (); |
94 const octave_idx_type np = p.rows (); | |
95 p = p.transpose (); | |
6823 | 96 |
6945 | 97 double *pt_array = p.fortran_vec (); |
98 | |
6823 | 99 //double pt_array[dim * np]; |
100 //for (int i = 0; i < np; i++) | |
101 // { | |
102 // for (int j = 0; j < dim; j++) | |
103 // { | |
104 // pt_array[j+i*dim] = p(i,j); | |
105 // } | |
106 // } | |
107 | |
108 boolT ismalloc = false; | |
109 | |
6945 | 110 OCTAVE_LOCAL_BUFFER (char, flags, 250); |
6880 | 111 |
6823 | 112 // hmm lot's of options for qhull here |
6945 | 113 sprintf (flags, "qhull v Fv T0 %s", options); |
6880 | 114 |
6945 | 115 // If you want some debugging information replace the 0 pointer |
116 // with stdout or some other file open for writing. | |
117 | |
118 FILE *outfile = 0; | |
6880 | 119 FILE *errfile = stderr; |
120 | |
6945 | 121 if (! qh_new_qhull (dim, np, pt_array, ismalloc, flags, outfile, errfile)) |
6823 | 122 { |
123 facetT *facet; | |
124 vertexT *vertex; | |
125 | |
126 octave_idx_type i = 0, n = 1, k = 0, m = 0, fidx = 0, j = 0, r = 0; | |
6945 | 127 OCTAVE_LOCAL_BUFFER (octave_idx_type, ni, np); |
6823 | 128 |
129 for (i = 0; i < np; i++) | |
130 ni[i] = 0; | |
6945 | 131 qh_setvoronoi_all (); |
6823 | 132 bool infinity_seen = false; |
6945 | 133 facetT *neighbor, **neighborp; |
6823 | 134 coordT *voronoi_vertex; |
6945 | 135 |
6823 | 136 FORALLfacets |
137 { | |
138 facet->seen = false; | |
139 } | |
6945 | 140 |
6823 | 141 FORALLvertices |
142 { | |
143 if (qh hull_dim == 3) | |
6945 | 144 qh_order_vertexneighbors (vertex); |
6823 | 145 infinity_seen = false; |
6945 | 146 |
147 FOREACHneighbor_ (vertex) | |
6823 | 148 { |
6945 | 149 if (! neighbor->upperdelaunay) |
6823 | 150 { |
6945 | 151 if (! neighbor->seen) |
6823 | 152 { |
153 n++; | |
6945 | 154 neighbor->seen = true; |
6823 | 155 } |
156 ni[k]++; | |
157 } | |
6945 | 158 else if (! infinity_seen) |
6823 | 159 { |
160 infinity_seen = true; | |
161 ni[k]++; | |
162 } | |
163 } | |
164 k++; | |
165 } | |
166 | |
6945 | 167 Matrix v (n, dim); |
6823 | 168 for (octave_idx_type d = 0; d < dim; d++) |
169 v(0,d) = octave_Inf; | |
170 | |
6945 | 171 boolMatrix AtInf (np, 1); |
6823 | 172 for (i = 0; i < np; i++) |
173 AtInf(i) = false; | |
8415
fa12c67683d3
Reszie third arg of __voronoi__ before retuning it
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
174 octave_value_list F (np, octave_value ()); |
6823 | 175 k = 0; |
176 i = 0; | |
6945 | 177 |
6823 | 178 FORALLfacets |
179 { | |
180 facet->seen = false; | |
181 } | |
6945 | 182 |
6823 | 183 FORALLvertices |
184 { | |
185 if (qh hull_dim == 3) | |
186 qh_order_vertexneighbors(vertex); | |
187 infinity_seen = false; | |
6945 | 188 RowVector facet_list (ni[k++]); |
6823 | 189 m = 0; |
6945 | 190 |
6823 | 191 FOREACHneighbor_(vertex) |
192 { | |
193 if (neighbor->upperdelaunay) | |
194 { | |
6945 | 195 if (! infinity_seen) |
6823 | 196 { |
197 infinity_seen = true; | |
198 facet_list(m++) = 1; | |
199 AtInf(j) = true; | |
200 } | |
201 } | |
202 else | |
203 { | |
6945 | 204 if (! neighbor->seen) |
6823 | 205 { |
206 voronoi_vertex = neighbor->center; | |
207 fidx = neighbor->id; | |
208 i++; | |
209 for (octave_idx_type d = 0; d < dim; d++) | |
210 { | |
211 v(i,d) = *(voronoi_vertex++); | |
212 } | |
213 neighbor->seen = true; | |
214 neighbor->visitid = i; | |
215 } | |
6945 | 216 facet_list(m++) = neighbor->visitid + 1; |
6823 | 217 } |
218 } | |
219 F(r++) = facet_list; | |
220 j++; | |
221 } | |
222 | |
223 Cell C(r, 1); | |
224 for (i = 0; i < r; i++) | |
225 C.elem (i) = F(i); | |
226 | |
227 retval(0) = v; | |
228 retval(1) = C; | |
8415
fa12c67683d3
Reszie third arg of __voronoi__ before retuning it
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
229 AtInf.resize (r, 1); |
6823 | 230 retval(2) = AtInf; |
231 | |
6945 | 232 // free long memory |
233 qh_freeqhull (! qh_ALL); | |
6823 | 234 |
6945 | 235 // free short memory and memory allocator |
6823 | 236 int curlong, totlong; |
237 qh_memfreeshort (&curlong, &totlong); | |
238 | |
6945 | 239 if (curlong || totlong) |
240 warning ("__voronoi__: did not free %d bytes of long memory (%d pieces)", totlong, curlong); | |
6823 | 241 } |
6945 | 242 else |
243 error ("__voronoi__: qhull failed"); | |
244 | |
6823 | 245 #else |
246 error ("__voronoi__: not available in this version of Octave"); | |
247 #endif | |
6945 | 248 |
6823 | 249 return retval; |
250 } |