Mercurial > octave-nkf

comparison src/dldfcn/__delaunayn__.cc @ 15076:000587f92082

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
rename src/DLD-FUNCTIONS directory to src/dldfcn * src/dldfcn: Rename from src/DLD-FUNCTIONS. * autogen.sh, src/Makefile.am, src/dldfcn/config-module.awk, src/dldfcn/config-module.sh: Change all uses of DLD-FUNCTIONS to be dldfcn. Change all uses of DLD_FUNCTIONS to be DLDFCN.
author John W. Eaton <jwe@octave.org>
date Tue, 31 Jul 2012 21:57:58 -0400
parents src/DLD-FUNCTIONS/__delaunayn__.cc@460a3c6d8bf1
children
comparison
equal deleted inserted replaced
15075:b62b0b85369c 15076:000587f92082
1 /*
2
3 Copyright (C) 2000-2012 Kai Habel
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
9 Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
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
18 along with Octave; see the file COPYING. If not, see
19 <http://www.gnu.org/licenses/>.
20
21 */
22
23 /*
24 16. July 2000 - Kai Habel: first release
25
26 25. September 2002 - Changes by Rafael Laboissiere <rafael@laboissiere.net>
27
28 * Added Qbb option to normalize the input and avoid crashes in Octave.
29 * delaunayn accepts now a second (optional) argument that must be a string
30 containing extra options to the qhull command.
31 * Fixed doc string. The dimension of the result matrix is [m, dim+1], and
32 not [n, dim-1].
33
34 6. June 2006: Changes by Alexander Barth <abarth@marine.usf.edu>
35
36 * triangulate non-simplicial facets
37 * allow options to be specified as cell array of strings
38 * change the default options (for compatibility with matlab)
39 */
40
41 #ifdef HAVE_CONFIG_H
42 #include <config.h>
43 #endif
44
45 #include <iostream>
46 #include <string>
47
48 #include "Cell.h"
49 #include "defun-dld.h"
50 #include "error.h"
51 #include "oct-obj.h"
52 #include "unwind-prot.h"
53
54 #if defined (HAVE_QHULL)
55 # include "oct-qhull.h"
56 # if defined (NEED_QHULL_VERSION)
57 char qh_version[] = "__delaunayn__.oct 2007-08-21";
58 # endif
59 #endif
60
61 static void
62 close_fcn (FILE *f)
63 {
64 gnulib::fclose (f);
65 }
66
67 DEFUN_DLD (__delaunayn__, args, ,
68 "-*- texinfo -*-\n\
69 @deftypefn {Loadable Function} {@var{T} =} __delaunayn__ (@var{pts})\n\
70 @deftypefnx {Loadable Function} {@var{T} =} __delaunayn__ (@var{pts}, @var{options})\n\
71 Undocumented internal function.\n\
72 @end deftypefn")
73
74 {
75 octave_value_list retval;
76
77 #if defined (HAVE_QHULL)
78
79 retval(0) = 0.0;
80
81 int nargin = args.length ();
82 if (nargin < 1 || nargin > 2)
83 {
84 print_usage ();
85 return retval;
86 }
87
88 Matrix p (args(0).matrix_value ());
89 const octave_idx_type dim = p.columns ();
90 const octave_idx_type n = p.rows ();
91
92 // Default options
93 std::string options;
94 if (dim <= 3)
95 options = "Qt Qbb Qc Qz";
96 else
97 options = "Qt Qbb Qc Qx";
98
99 if (nargin == 2)
100 {
101 if (args(1).is_string ())
102 options = args(1).string_value ();
103 else if (args(1).is_empty ())
104 ; // Use default options
105 else if (args(1).is_cellstr ())
106 {
107 options = "";
108 Array<std::string> tmp = args(1).cellstr_value ();
109
110 for (octave_idx_type i = 0; i < tmp.numel (); i++)
111 options += tmp(i) + " ";
112 }
113 else
114 {
115 error ("__delaunayn__: OPTIONS argument must be a string, cell array of strings, or empty");
116 return retval;
117 }
118 }
119
120 if (n > dim + 1)
121 {
122 p = p.transpose ();
123 double *pt_array = p.fortran_vec ();
124 boolT ismalloc = false;
125
126 // Qhull flags argument is not const char*
127 OCTAVE_LOCAL_BUFFER (char, flags, 9 + options.length ());
128
129 sprintf (flags, "qhull d %s", options.c_str ());
130
131 unwind_protect frame;
132
133 // Replace the outfile pointer with stdout for debugging information.
134 #if defined (OCTAVE_HAVE_WINDOWS_FILESYSTEM) && ! defined (OCTAVE_HAVE_POSIX_FILESYSTEM)
135 FILE *outfile = gnulib::fopen ("NUL", "w");
136 #else
137 FILE *outfile = gnulib::fopen ("/dev/null", "w");
138 #endif
139 FILE *errfile = stderr;
140
141 if (outfile)
142 frame.add_fcn (close_fcn, outfile);
143 else
144 {
145 error ("__delaunayn__: unable to create temporary file for output");
146 return retval;
147 }
148
149 int exitcode = qh_new_qhull (dim, n, pt_array,
150 ismalloc, flags, outfile, errfile);
151 if (! exitcode)
152 {
153 // triangulate non-simplicial facets
154 qh_triangulate ();
155
156 facetT *facet;
157 vertexT *vertex, **vertexp;
158 octave_idx_type nf = 0, i = 0;
159
160 FORALLfacets
161 {
162 if (! facet->upperdelaunay)
163 nf++;
164
165 // Double check. Non-simplicial facets will cause segfault below
166 if (! facet->simplicial)
167 {
168 error ("__delaunayn__: Qhull returned non-simplicial facets -- try delaunayn with different options");
169 exitcode = 1;
170 break;
171 }
172 }
173
174 if (! exitcode)
175 {
176 Matrix simpl (nf, dim+1);
177
178 FORALLfacets
179 {
180 if (! facet->upperdelaunay)
181 {
182 octave_idx_type j = 0;
183
184 FOREACHvertex_ (facet->vertices)
185 {
186 simpl(i, j++) = 1 + qh_pointid(vertex->point);
187 }
188 i++;
189 }
190 }
191
192 retval(0) = simpl;
193 }
194 }
195 else
196 error ("__delaunayn__: qhull failed");
197
198 // Free memory from Qhull
199 qh_freeqhull (! qh_ALL);
200
201 int curlong, totlong;
202 qh_memfreeshort (&curlong, &totlong);
203
204 if (curlong || totlong)
205 warning ("__delaunay__: did not free %d bytes of long memory (%d pieces)",
206 totlong, curlong);
207 }
208 else if (n == dim + 1)
209 {
210 // one should check if nx points span a simplex
211 // I will look at this later.
212 RowVector vec (n);
213 for (octave_idx_type i = 0; i < n; i++)
214 vec(i) = i + 1.0;
215
216 retval(0) = vec;
217 }
218
219 #else
220 error ("__delaunayn__: not available in this version of Octave");
221 #endif
222
223 return retval;
224 }
225
226 /*
227 ## No test needed for internal helper function.
228 %!assert (1)
229 */