Mercurial > hg > octave-lyh
annotate doc/interpreter/geometry.txi @ 8434:7ceb99b0abbf
Removed typo from geometry.txi
author | Francesco Potortì <pot@gnu.org> |
---|---|
date | Fri, 26 Dec 2008 22:56:19 +0100 |
parents | fa78cb8d8a5c |
children | 8ae26422a6ce |
rev | line source |
---|---|
7018 | 1 @c Copyright (C) 2007 John W. Eaton and David Bateman |
2 @c | |
3 @c This file is part of Octave. | |
4 @c | |
5 @c Octave is free software; you can redistribute it and/or modify it | |
6 @c under the terms of the GNU General Public License as published by the | |
7 @c Free Software Foundation; either version 3 of the License, or (at | |
8 @c your option) any later version. | |
9 @c | |
10 @c Octave is distributed in the hope that it will be useful, but WITHOUT | |
11 @c ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
12 @c FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
13 @c for more details. | |
14 @c | |
15 @c You should have received a copy of the GNU General Public License | |
16 @c along with Octave; see the file COPYING. If not, see | |
17 @c <http://www.gnu.org/licenses/>. | |
6558 | 18 |
19 @node Geometry | |
20 @chapter Geometry | |
21 | |
8347
fa78cb8d8a5c
corrections for typos
Brian Gough<bjg@network-theory.co.uk>
parents:
7984
diff
changeset
|
22 Much of geometry code in Octave is based on the QHull library@footnote{Barber, |
6832 | 23 C.B., Dobkin, D.P., and Huhdanpaa, H.T., "The Quickhull algorithm for |
24 convex hulls," ACM Trans. on Mathematical Software, 22(4):469-483, Dec | |
25 1996, @url{http://www.qhull.org}}. Some of the documentation for Qhull, | |
26 particularly for the options that can be passed to @code{delaunay}, | |
27 @code{voronoi} and @code{convhull}, etc, is relevant to Octave users. | |
28 | |
6823 | 29 @menu |
30 * Delaunay Triangulation:: | |
31 * Voronoi Diagrams:: | |
32 * Convex Hull:: | |
33 * Interpolation on Scattered Data:: | |
34 @end menu | |
35 | |
36 @node Delaunay Triangulation | |
37 @section Delaunay Triangulation | |
38 | |
6832 | 39 The Delaunay triangulation is constructed from a set of |
40 circum-circles. These circum-circles are chosen so that there are at | |
41 least three of the points in the set to triangulation on the | |
42 circumference of the circum-circle. None of the points in the set of | |
43 points falls within any of the circum-circles. | |
44 | |
45 In general there are only three points on the circumference of any | |
8434
7ceb99b0abbf
Removed typo from geometry.txi
Francesco Potortì <pot@gnu.org>
parents:
8347
diff
changeset
|
46 circum-circle. However, in some cases, and in particular for the |
6832 | 47 case of a regular grid, 4 or more points can be on a single |
48 circum-circle. In this case the Delaunay triangulation is not unique. | |
49 | |
6823 | 50 @DOCSTRING(delaunay) |
51 | |
52 The 3- and N-dimensional extension of the Delaunay triangulation are | |
53 given by @code{delaunay3} and @code{delaunayn} respectively. | |
54 @code{delaunay3} returns a set of tetrahedra that satisfy the | |
55 Delaunay circum-circle criteria. Similarly, @code{delaunayn} returns the | |
56 N-dimensional simplex satisfying the Delaunay circum-circle criteria. | |
7007 | 57 The N-dimensional extension of a triangulation is called a tessellation. |
6823 | 58 |
59 @DOCSTRING(delaunay3) | |
60 | |
61 @DOCSTRING(delaunayn) | |
62 | |
6832 | 63 An example of a Delaunay triangulation of a set of points is |
64 | |
65 @example | |
66 @group | |
67 rand ("state", 2); | |
68 x = rand (10, 1); | |
69 y = rand (10, 1); | |
70 T = delaunay (x, y); | |
71 X = [ x(T(:,1)); x(T(:,2)); x(T(:,3)); x(T(:,1)) ]; | |
72 Y = [ y(T(:,1)); y(T(:,2)); y(T(:,3)); y(T(:,1)) ]; | |
73 axis ([0, 1, 0, 1]); | |
74 plot(X, Y, "b", x, y, "r*"); | |
75 @end group | |
76 @end example | |
77 | |
6855 | 78 @ifset HAVE_QHULL |
6832 | 79 @ifnotinfo |
80 @noindent | |
81 The result of which can be seen in @ref{fig:delaunay}. | |
82 | |
83 @float Figure,fig:delaunay | |
84 @image{delaunay,8cm} | |
85 @caption{Delaunay triangulation of a random set of points} | |
86 @end float | |
87 @end ifnotinfo | |
6855 | 88 @end ifset |
6832 | 89 |
6823 | 90 @menu |
6832 | 91 * Plotting the Triangulation:: |
6823 | 92 * Identifying points in Triangulation:: |
93 @end menu | |
94 | |
6832 | 95 @node Plotting the Triangulation |
96 @subsection Plotting the Triangulation | |
97 | |
98 Octave has the functions @code{triplot} and @code{trimesh} to plot the | |
99 Delaunay triangulation of a 2-dimensional set of points. | |
100 | |
101 @DOCSTRING(triplot) | |
102 | |
103 @DOCSTRING(trimesh) | |
104 | |
105 The difference between @code{triplot} and @code{trimesh} is that the | |
106 former only plots the 2-dimensional triangulation itself, whereas the | |
107 second plots the value of some function @code{f (@var{x}, @var{y})}. | |
108 An example of the use of the @code{triplot} function is | |
109 | |
110 @example | |
111 @group | |
112 rand ("state", 2) | |
113 x = rand (20, 1); | |
114 y = rand (20, 1); | |
115 tri = delaunay (x, y); | |
116 triplot (tri, x, y); | |
117 @end group | |
118 @end example | |
119 | |
120 that plot the Delaunay triangulation of a set of random points in | |
121 2-dimensions. | |
122 @ifnotinfo | |
123 The output of the above can be seen in @ref{fig:triplot}. | |
124 | |
125 @float Figure,fig:triplot | |
126 @image{triplot,8cm} | |
127 @caption{Delaunay triangulation of a random set of points} | |
128 @end float | |
129 @end ifnotinfo | |
130 | |
6823 | 131 @node Identifying points in Triangulation |
132 @subsection Identifying points in Triangulation | |
133 | |
134 It is often necessary to identify whether a particular point in the | |
7007 | 135 N-dimensional space is within the Delaunay tessellation of a set of |
6823 | 136 points in this N-dimensional space, and if so which N-Simplex contains |
7007 | 137 the point and which point in the tessellation is closest to the desired |
6823 | 138 point. The functions @code{tsearch} and @code{dsearch} perform this |
139 function in a triangulation, and @code{tsearchn} and @code{dsearchn} in | |
7007 | 140 an N-dimensional tessellation. |
6823 | 141 |
142 To identify whether a particular point represented by a vector @var{p} | |
143 falls within one of the simplices of an N-Simplex, we can write the | |
144 Cartesian coordinates of the point in a parametric form with respect to | |
145 the N-Simplex. This parametric form is called the Barycentric | |
146 Coordinates of the point. If the points defining the N-Simplex are given | |
147 by @code{@var{N} + 1} vectors @var{t}(@var{i},:), then the Barycentric | |
8347
fa78cb8d8a5c
corrections for typos
Brian Gough<bjg@network-theory.co.uk>
parents:
7984
diff
changeset
|
148 coordinates defining the point @var{p} are given by |
6823 | 149 |
150 @example | |
151 @var{p} = sum (@var{beta}(1:@var{N}+1) * @var{t}(1:@var{N}+1),:) | |
152 @end example | |
153 | |
154 @noindent | |
155 where there are @code{@var{N} + 1} values @code{@var{beta}(@var{i})} | |
156 that together as a vector represent the Barycentric coordinates of the | |
157 point @var{p}. To ensure a unique solution for the values of | |
158 @code{@var{beta}(@var{i})} an additional criteria of | |
159 | |
160 @example | |
161 sum (@var{beta}(1:@var{N}+1)) == 1 | |
162 @end example | |
163 | |
164 @noindent | |
165 is imposed, and we can therefore write the above as | |
166 | |
167 @example | |
168 @var{p} - @var{t}(end, :) = @var{beta}(1:end-1) * (@var{t}(1:end-1, :) | |
169 - ones(@var{N}, 1) * @var{t}(end, :) | |
170 @end example | |
171 | |
172 @noindent | |
173 Solving for @var{beta} we can then write | |
174 | |
175 @example | |
176 @var{beta}(1:end-1) = (@var{p} - @var{t}(end, :)) / (@var{t}(1:end-1, :) | |
177 - ones(@var{N}, 1) * @var{t}(end, :)) | |
178 @var{beta}(end) = sum(@var{beta}(1:end-1)) | |
179 @end example | |
180 | |
181 @noindent | |
182 which gives the formula for the conversion of the Cartesian coordinates | |
183 of the point @var{p} to the Barycentric coordinates @var{beta}. An | |
184 important property of the Barycentric coordinates is that for all points | |
185 in the N-Simplex | |
186 | |
187 @example | |
188 0 <= @var{beta}(@var{i}) <= 1 | |
189 @end example | |
190 | |
191 @noindent | |
192 Therefore, the test in @code{tsearch} and @code{tsearchn} essentially | |
193 only needs to express each point in terms of the Barycentric coordinates | |
194 of each of the simplices of the N-Simplex and test the values of | |
195 @var{beta}. This is exactly the implementation used in | |
196 @code{tsearchn}. @code{tsearch} is optimized for 2-dimensions and the | |
197 Barycentric coordinates are not explicitly formed. | |
198 | |
199 @DOCSTRING(tsearch) | |
200 | |
201 @DOCSTRING(tsearchn) | |
202 | |
203 An example of the use of @code{tsearch} can be seen with the simple | |
204 triangulation | |
205 | |
206 @example | |
207 @group | |
208 @var{x} = [-1; -1; 1; 1]; | |
209 @var{y} = [-1; 1; -1; 1]; | |
210 @var{tri} = [1, 2, 3; 2, 3, 1]; | |
211 @end group | |
212 @end example | |
213 | |
214 @noindent | |
215 consisting of two triangles defined by @var{tri}. We can then identify | |
216 which triangle a point falls in like | |
217 | |
218 @example | |
219 @group | |
220 tsearch (@var{x}, @var{y}, @var{tri}, -0.5, -0.5) | |
221 @result{} 1 | |
222 tsearch (@var{x}, @var{y}, @var{tri}, 0.5, 0.5) | |
223 @result{} 2 | |
224 @end group | |
225 @end example | |
226 | |
227 @noindent | |
228 and we can confirm that a point doesn't lie within one of the triangles like | |
229 | |
230 @example | |
231 @group | |
232 tsearch (@var{x}, @var{y}, @var{tri}, 2, 2) | |
233 @result{} NaN | |
234 @end group | |
235 @end example | |
236 | |
237 The @code{dsearch} and @code{dsearchn} find the closest point in a | |
238 tessellation to the desired point. The desired point does not | |
239 necessarily have to be in the tessellation, and even if it the returned | |
6832 | 240 point of the tessellation does not have to be one of the vertexes of the |
6823 | 241 N-simplex within which the desired point is found. |
242 | |
243 @DOCSTRING(dsearch) | |
244 | |
245 @DOCSTRING(dsearchn) | |
246 | |
247 An example of the use of @code{dsearch}, using the above values of | |
248 @var{x}, @var{y} and @var{tri} is | |
249 | |
250 @example | |
251 @group | |
252 dsearch (@var{x}, @var{y}, @var{tri}, -2, -2) | |
253 @result{} 1 | |
254 @end group | |
255 @end example | |
256 | |
257 If you wish the points that are outside the tessellation to be flagged, | |
258 then @code{dsearchn} can be used as | |
259 | |
260 @example | |
261 @group | |
262 dsearchn ([@var{x}, @var{y}], @var{tri}, [-2, -2], NaN) | |
263 @result{} NaN | |
264 dsearchn ([@var{x}, @var{y}], @var{tri}, [-0.5, -0.5], NaN) | |
265 @result{} 1 | |
266 @end group | |
267 @end example | |
268 | |
269 @noindent | |
270 where the point outside the tessellation are then flagged with @code{NaN}. | |
271 | |
272 @node Voronoi Diagrams | |
273 @section Voronoi Diagrams | |
274 | |
275 A Voronoi diagram or Voronoi tessellation of a set of points @var{s} in | |
276 an N-dimensional space, is the tessellation of the N-dimensional space | |
277 such that all points in @code{@var{v}(@var{p})}, a partitions of the | |
278 tessellation where @var{p} is a member of @var{s}, are closer to @var{p} | |
279 than any other point in @var{s}. The Voronoi diagram is related to the | |
6832 | 280 Delaunay triangulation of a set of points, in that the vertexes of the |
8347
fa78cb8d8a5c
corrections for typos
Brian Gough<bjg@network-theory.co.uk>
parents:
7984
diff
changeset
|
281 Voronoi tessellation are the centers of the circum-circles of the |
6832 | 282 simplicies of the Delaunay tessellation. |
6823 | 283 |
284 @DOCSTRING(voronoi) | |
285 | |
286 @DOCSTRING(voronoin) | |
287 | |
6832 | 288 An example of the use of @code{voronoi} is |
289 | |
290 @example | |
291 @group | |
292 rand("state",9); | |
293 x = rand(10,1); | |
294 y = rand(10,1); | |
295 tri = delaunay (x, y); | |
296 [vx, vy] = voronoi (x, y, tri); | |
297 triplot (tri, x, y, "b"); | |
298 hold on; | |
299 plot (vx, vy, "r"); | |
300 @end group | |
301 @end example | |
302 | |
6855 | 303 @ifset HAVE_QHULL |
6832 | 304 @ifnotinfo |
305 @noindent | |
306 The result of which can be seen in @ref{fig:voronoi}. Note that the | |
307 circum-circle of one of the triangles has been added to this figure, to | |
308 make the relationship between the Delaunay tessellation and the Voronoi | |
309 diagram clearer. | |
310 | |
311 @float Figure,fig:voronoi | |
312 @image{voronoi,8cm} | |
313 @caption{Delaunay triangulation and Voronoi diagram of a random set of points} | |
314 @end float | |
315 @end ifnotinfo | |
6855 | 316 @end ifset |
6832 | 317 |
6847 | 318 Additional information about the size of the facets of a Voronoi |
319 diagram, and which points of a set of points is in a polygon can be had | |
320 with the @code{polyarea} and @code{inpolygon} functions respectively. | |
6832 | 321 |
322 @DOCSTRING(polyarea) | |
323 | |
324 An example of the use of @code{polyarea} might be | |
325 | |
326 @example | |
327 @group | |
328 rand ("state", 2); | |
329 x = rand (10, 1); | |
330 y = rand (10, 1); | |
331 [c, f] = voronoin ([x, y]); | |
332 af = zeros (size(f)); | |
333 for i = 1 : length (f) | |
334 af(i) = polyarea (c (f @{i, :@}, 1), c (f @{i, :@}, 2)); | |
335 endfor | |
336 @end group | |
337 @end example | |
338 | |
7984
bbaa5d7d0143
Some documentation updates
David Bateman <dbateman@free.fr>
parents:
7018
diff
changeset
|
339 Facets of the Voronoi diagram with a vertex at infinity have infinity |
bbaa5d7d0143
Some documentation updates
David Bateman <dbateman@free.fr>
parents:
7018
diff
changeset
|
340 area. A simplified version of @code{polyarea} for rectangles is |
bbaa5d7d0143
Some documentation updates
David Bateman <dbateman@free.fr>
parents:
7018
diff
changeset
|
341 available with @code{rectint} |
bbaa5d7d0143
Some documentation updates
David Bateman <dbateman@free.fr>
parents:
7018
diff
changeset
|
342 |
bbaa5d7d0143
Some documentation updates
David Bateman <dbateman@free.fr>
parents:
7018
diff
changeset
|
343 @DOCSTRING(rectint) |
6832 | 344 |
6847 | 345 @DOCSTRING(inpolygon) |
346 | |
347 An example of the use of @code{inpolygon} might be | |
348 | |
349 @example | |
350 @group | |
351 randn ("state", 2); | |
352 x = randn (100, 1); | |
353 y = randn (100, 1); | |
354 vx = cos (pi * [-1 : 0.1: 1]); | |
355 vy = sin (pi * [-1 : 0.1 : 1]); | |
356 in = inpolygon (x, y, vx, vy); | |
357 plot(vx, vy, x(in), y(in), "r+", x(!in), y(!in), "bo"); | |
358 axis ([-2, 2, -2, 2]); | |
359 @end group | |
360 @end example | |
361 | |
362 @ifnotinfo | |
363 @noindent | |
364 The result of which can be seen in @ref{fig:inpolygon}. | |
365 | |
366 @float Figure,fig:inpolygon | |
367 @image{inpolygon,8cm} | |
368 @caption{Demonstration of the @code{inpolygon} function to determine the | |
369 points inside a polygon} | |
370 @end float | |
371 @end ifnotinfo | |
372 | |
6823 | 373 @node Convex Hull |
374 @section Convex Hull | |
375 | |
7001 | 376 The convex hull of a set of points is the minimum convex envelope |
6832 | 377 containing all of the points. Octave has the functions @code{convhull} |
7007 | 378 and @code{convhulln} to calculate the convex hull of 2-dimensional and |
6832 | 379 N-dimensional sets of points. |
380 | |
6823 | 381 @DOCSTRING(convhull) |
382 | |
383 @DOCSTRING(convhulln) | |
384 | |
6832 | 385 An example of the use of @code{convhull} is |
6823 | 386 |
6832 | 387 @example |
388 @group | |
389 x = -3:0.05:3; | |
390 y = abs (sin (x)); | |
391 k = convhull (x, y); | |
392 plot (x(k), y(k), "r-", x, y, "b+"); | |
393 axis ([-3.05, 3.05, -0.05, 1.05]); | |
394 @end group | |
395 @end example | |
6823 | 396 |
6855 | 397 @ifset HAVE_QHULL |
6832 | 398 @ifnotinfo |
399 @noindent | |
400 The output of the above can be seen in @ref{fig:convhull}. | |
6823 | 401 |
6832 | 402 @float Figure,fig:convhull |
403 @image{convhull,8cm} | |
404 @caption{The convex hull of a simple set of points} | |
405 @end float | |
406 @end ifnotinfo | |
6855 | 407 @end ifset |
6823 | 408 |
409 @node Interpolation on Scattered Data | |
410 @section Interpolation on Scattered Data | |
411 | |
6832 | 412 An important use of the Delaunay tessellation is that it can be used to |
413 interpolate from scattered data to an arbitrary set of points. To do | |
414 this the N-simplex of the known set of points is calculated with | |
415 @code{delaunay}, @code{delaunay3} or @code{delaunayn}. Then the | |
416 simplicies in to which the desired points are found are | |
417 identified. Finally the vertices of the simplicies are used to | |
418 interpolate to the desired points. The functions that perform this | |
419 interpolation are @code{griddata}, @code{griddata3} and | |
420 @code{griddatan}. | |
421 | |
6823 | 422 @DOCSTRING(griddata) |
423 | |
424 @DOCSTRING(griddata3) | |
425 | |
426 @DOCSTRING(griddatan) | |
6832 | 427 |
428 An example of the use of the @code{griddata} function is | |
429 | |
430 @example | |
431 @group | |
432 rand("state",1); | |
433 x=2*rand(1000,1)-1; | |
434 y=2*rand(size(x))-1; | |
435 z=sin(2*(x.^2+y.^2)); | |
436 [xx,yy]=meshgrid(linspace(-1,1,32)); | |
437 griddata(x,y,z,xx,yy); | |
438 @end group | |
439 @end example | |
440 | |
6855 | 441 @ifset HAVE_QHULL |
6832 | 442 @noindent |
443 that interpolates from a random scattering of points, to a uniform | |
444 grid. | |
445 @ifnotinfo | |
446 The output of the above can be seen in @ref{fig:griddata}. | |
447 | |
448 @float Figure,fig:griddata | |
449 @image{griddata,8cm} | |
450 @caption{Interpolation from a scattered data to a regular grid} | |
451 @end float | |
452 @end ifnotinfo | |
6855 | 453 @end ifset |