1 ## Copyright (C) 2008 Andreas Stahel <Andreas.Stahel@bfh.ch>
3 ## This program is free software; you can redistribute it and/or modify
4 ## it under the terms of the GNU General Public License as published by
5 ## the Free Software Foundation; either version 2 of the License, or
6 ## (at your option) any later version.
8 ## This program is distributed in the hope that it will be useful,
9 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
10 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 ## GNU General Public License for more details.
13 ## You should have received a copy of the GNU General Public License
14 ## along with this program; if not, write to the Free Software
15 ## Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
19 ## @deftypefn {Function File} {} tricontour (@var{tri}, @var{x}, @var{y}, @var{z}, @var{levels})
20 ## @deftypefnx {Function File} {} tricontour (@var{tri}, @var{x}, @var{y}, @var{z}, @var{levels}, @var{linespec})
21 ## Plot level curves for the values of @code{@var{z}} on a triangular mesh in 2D.
23 ## The variable @var{tri} is the triangular meshing of the points
24 ## @code{(@var{x}, @var{y})} which is returned from @code{delaunay}. The
25 ## variable @var{levels} is a vector with the values of the countour levels. If
26 ## @var{levels} is a scalar, then it corresponds to the number of
27 ## level curves to be drawn. If exactly one level curve is desired, list
28 ## the level twice in the vector @var{levels}.
30 ## If given, @var{linespec} determines the properties to use for the
33 ## The output argument @var{h} is the graphic handle to the plot.
34 ## @seealso{plot, trimesh, delaunay}
37 function h = tricontour (tri, x, y, z, levels, varargin)
44 dom=mean(dom)+0.99*(dom-mean(dom));
45 levels=linspace(dom(1),dom(2),levels);
50 lmax=levels(length(levels));
52 pData=[]; %% no preallocation
53 %% pData=zeros(12000,2); %% preallocation
57 values=[z(tri(el,1)),z(tri(el,2)),z(tri(el,3))];
60 locallevel=levels(minval<=levels+eps); # select the levels to be plotted
62 locallevel=locallevel(locallevel<=maxval+eps);
69 if (abs(dl(1))<=10*eps)
70 points(npoints,:)=[x(tri(el,1)),y(tri(el,1))];npoints++;endif
71 if (abs(dl(2))<=10*eps)
72 points(npoints,:)=[x(tri(el,2)),y(tri(el,2))];npoints++;endif
73 if (abs(dl(3))<=10*eps)
74 points(npoints,:)=[x(tri(el,3)),y(tri(el,3))];npoints++;endif
77 if ((dl(1)*dl(2)) < 0) # intersection between 1st and 2nd point
78 points(npoints,:)= ( dl(2)*[x(tri(el,1)),y(tri(el,1))]...
79 -dl(1)*[x(tri(el,2)),y(tri(el,2))])/(dl(2)-dl(1));
82 if ((dl(1)*dl(3)) < 0) # intersection between 1st and 3rd point
83 points(npoints,:)= ( dl(3)*[x(tri(el,1)),y(tri(el,1))]...
84 -dl(1)*[x(tri(el,3)),y(tri(el,3))])/(dl(3)-dl(1));
87 if ((dl(3)*dl(2)) < 0) # intersection between 2nd and 3rd point
88 points(npoints,:)= ( dl(2)*[x(tri(el,3)),y(tri(el,3))]...
89 -dl(3)*[x(tri(el,2)),y(tri(el,2))])/(dl(2)-dl(3));
93 pData=[pData;points; NaN,NaN ]; %% no preallocation
94 %% pData(pPoints+1:pPoints+npoints,1:2)=[points; NaN,NaN ]; %% preallocation
99 pData=pData(1:pPoints-1,:);
102 h= plot(pData(:,1),pData(:,2),varargin(:));
104 plot(pData(:,1),pData(:,2),varargin(:));
111 %! x = rand (100, 1)-0.5;
112 %! y = rand (100, 1)-0.5;
114 %! tri = delaunay (x, y);
115 %! tricontour (tri, x, y, z, [-0.25:0.05:0.25]);