|
7042
|
1 ## Copyright (C) 2007 Kai Habel |
|
|
2 ## Copyright (C) 2003 Shai Ayal |
|
|
3 ## |
|
7164
|
4 ## This file is part of Octave. |
|
|
5 ## |
|
|
6 ## Octave is free software; you can redistribute it and/or modify it |
|
7042
|
7 ## under the terms of the GNU General Public License as published by |
|
7164
|
8 ## the Free Software Foundation; either version 3 of the License, or (at |
|
|
9 ## your option) any later version. |
|
7042
|
10 ## |
|
7164
|
11 ## Octave is distributed in the hope that it will be useful, but |
|
7042
|
12 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
|
|
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
|
|
14 ## General Public License for more details. |
|
|
15 ## |
|
|
16 ## You should have received a copy of the GNU General Public License |
|
7164
|
17 ## along with Octave; see the file COPYING. If not, see |
|
|
18 ## <http://www.gnu.org/licenses/>. |
|
7042
|
19 |
|
|
20 ## -*- texinfo -*- |
|
|
21 ## @deftypefn {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}, @var{lvl}) |
|
|
22 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}, @var{n}) |
|
|
23 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}) |
|
|
24 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}, @var{n}) |
|
|
25 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}, @var{lvl}) |
|
|
26 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}) |
|
|
27 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{ax}, @dots{}) |
|
|
28 ## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@dots{}, @var{"property"}, @var{val}) |
|
|
29 ## Compute and plot filled contours of the matrix @var{z}. |
|
|
30 ## Parameters @var{x}, @var{y} and @var{n} or @var{lvl} are optional. |
|
|
31 ## |
|
|
32 ## The return value @var{c} is a 2xn matrix containing the contour lines |
|
|
33 ## as described in the help to the contourc function. |
|
|
34 ## |
|
|
35 ## The return value @var{h} is handle-vector to the patch objects creating |
|
|
36 ## the filled contours. |
|
|
37 ## |
|
|
38 ## If @var{x} and @var{y} are ommited they are taken as the row/column |
|
|
39 ## index of @var{z}. @var{n} is a scalar denoting the number of lines |
|
|
40 ## to compute. Alternatively @var{lvl} is a vector containing the |
|
|
41 ## contour levels. If only one value (e.g. lvl0) is wanted, set |
|
|
42 ## @var{lvl} to [lvl0, lvl0]. If both @var{n} or @var{lvl} are omitted |
|
|
43 ## a default value of 10 contour level is assumed. |
|
|
44 ## |
|
|
45 ## If provided, the filled contours are added to the axes object |
|
|
46 ## @var{ax} instead of the current axis. |
|
|
47 ## |
|
|
48 ## The following example plots filled contours of the @code{peaks} |
|
|
49 ## function. |
|
|
50 ## @example |
|
|
51 ## [x, y, z] = peaks (50); |
|
|
52 ## contourf (x, y, z, -7:9) |
|
|
53 ## @end example |
|
|
54 ## @seealso{contour, contourc, patch} |
|
|
55 ## @end deftypefn |
|
|
56 |
|
|
57 ## Author: Kai Habel <kai.habel@gmx.de> |
|
|
58 ## Author: shaia |
|
|
59 |
|
|
60 function varargout = contourf (varargin) |
|
|
61 |
|
|
62 [X, Y, Z, lvl, ax, patch_props] = parse_args (varargin); |
|
|
63 [nr, nc] = size (Z); |
|
|
64 [minx, maxx] = deal (min (X(:)), max (X(:))); |
|
|
65 [miny, maxy] = deal (min (Y(:)), max (Y(:))); |
|
|
66 |
|
|
67 if (diff (lvl) < 10*eps) |
|
|
68 lvl_eps = 1e-6; |
|
|
69 else |
|
|
70 lvl_eps = min (diff (lvl)) / 1000.0; |
|
|
71 endif |
|
|
72 |
|
|
73 X0 = prepad(X, nc+1, 2 * X(1, 1) - X(1, 2), 2); |
|
|
74 X0 = postpad(X0, nc+2, 2 * X(1, nc) - X(1, nc - 1), 2); |
|
|
75 X0 = [X0(1, :); X0; X0(1, :)]; |
|
|
76 Y0 = prepad(Y, nr+1, 2 * Y(1, 1) - Y(2, 1), 1); |
|
|
77 Y0 = postpad(Y0, nr+2, 2 * Y(nr, 1) - Y(nr - 1, 1)); |
|
|
78 Y0 = [Y0(:, 1), Y0, Y0(:, 1)]; |
|
|
79 |
|
|
80 Z0 = -Inf(nr+2, nc+2); |
|
|
81 Z0(2:nr+1, 2:nc+1) = Z; |
|
|
82 [c, lev] = contourc (X0, Y0, Z0, lvl); |
|
|
83 cmap = colormap (); |
|
|
84 |
|
|
85 levx = linspace (min (lev), max (lev), size (cmap, 1)); |
|
|
86 |
|
|
87 newplot (); |
|
|
88 |
|
|
89 ## Decode contourc output format. |
|
|
90 i1 = 1; |
|
|
91 ncont = 0; |
|
|
92 while (i1 < columns (c)) |
|
|
93 ncont++; |
|
|
94 cont_lev(ncont) = c(1, i1); |
|
|
95 cont_len(ncont) = c(2, i1); |
|
|
96 cont_idx(ncont) = i1+1; |
|
|
97 |
|
|
98 ii = i1+1:i1+cont_len(ncont); |
|
|
99 cur_cont = c(:, ii); |
|
|
100 c(:, ii); |
|
|
101 startidx = ii(1); |
|
|
102 stopidx = ii(end); |
|
|
103 cont_area(ncont) = polyarea (c(1, ii), c(2, ii)); |
|
|
104 i1 += c(2, i1) + 1; |
|
|
105 endwhile |
|
|
106 |
|
|
107 ## Handle for each level the case where we have (a) hole(s) in a patch. |
|
|
108 ## Those are to be filled with the color of level below or with the |
|
|
109 ## background colour. |
|
|
110 for k = 1:numel (lev) |
|
|
111 lvl_idx = find (abs (cont_lev - lev(k)) < lvl_eps); |
|
|
112 len = numel (lvl_idx); |
|
|
113 if (len > 1) |
|
|
114 ## mark = logical(zeros(size(lvl_idx))); |
|
|
115 mark = false (size (lvl_idx)); |
|
|
116 a = 1; |
|
|
117 while (a < len) |
|
|
118 # take 1st patch |
|
|
119 b = a + 1; |
|
|
120 pa_idx = lvl_idx(a); |
|
|
121 # get pointer to contour start, and contour length |
|
|
122 curr_ct_idx = cont_idx(pa_idx); |
|
|
123 curr_ct_len = cont_len(pa_idx); |
|
|
124 # get contour |
|
|
125 curr_ct = c(:, curr_ct_idx:curr_ct_idx+curr_ct_len-1); |
|
|
126 b_vec = (a+1):len; |
|
|
127 next_ct_pt_vec = c(:, cont_idx(lvl_idx(b_vec))); |
|
|
128 in = inpolygon (next_ct_pt_vec(1,:), next_ct_pt_vec(2,:), |
|
|
129 curr_ct(1, :), curr_ct(2, :)); |
|
|
130 mark(b_vec(in)) = !mark(b_vec(in)); |
|
|
131 a++; |
|
|
132 endwhile |
|
|
133 if (numel (mark) > 0) |
|
|
134 ## All marked contours describe a hole in a larger contour of |
|
|
135 ## the same level and must be filled with colour of level below. |
|
|
136 ma_idx = lvl_idx(mark); |
|
|
137 if (k > 1) |
|
|
138 ## Find color of level below. |
|
|
139 tmp = find(abs(cont_lev - lev(k - 1)) < lvl_eps); |
|
|
140 lvl_bel_idx = tmp(1); |
|
|
141 ## Set color of patches found. |
|
|
142 cont_lev(ma_idx) = cont_lev(lvl_bel_idx); |
|
|
143 else |
|
|
144 ## Set lowest level contour to NaN. |
|
|
145 cont_lev(ma_idx) = NaN; |
|
|
146 endif |
|
|
147 endif |
|
|
148 endif |
|
|
149 endfor |
|
|
150 |
|
|
151 ## The algorithm can create patches with the size of the plotting |
|
|
152 ## area, we would like to draw only the patch with the highest level. |
|
|
153 del_idx = []; |
|
|
154 max_idx = find (cont_area == max (cont_area)); |
|
|
155 if (numel (max_idx) > 1) |
|
|
156 # delete double entries |
|
|
157 del_idx = max_idx(1:end-1); |
|
|
158 cont_area(del_idx) = cont_lev(del_idx) = []; |
|
|
159 cont_len(del_idx) = cont_idx(del_idx) = []; |
|
|
160 endif |
|
|
161 |
|
|
162 ## Now we have everything together and can start plotting the patches |
|
|
163 ## beginning with largest area. |
|
|
164 [tmp, svec] = sort (cont_area); |
|
|
165 len = ncont - numel (del_idx); |
|
|
166 h = zeros (1, len); |
|
|
167 for n = len:-1:1 |
|
|
168 idx = svec(n); |
|
|
169 ii = cont_idx(idx):cont_idx(idx) + cont_len(idx) - 2; |
|
|
170 h(n) = patch (c(1, ii), c(2, ii), cont_lev(idx), patch_props{:}); |
|
|
171 endfor |
|
|
172 |
|
|
173 if (min (lev) == max (lev)) |
|
|
174 set (gca (), "clim", [min(lev)-1, max(lev)+1]); |
|
|
175 else |
|
|
176 set (gca(), "clim", [min(lev), max(lev)]); |
|
|
177 endif |
|
|
178 |
|
|
179 if (nargout > 0) |
|
|
180 varargout{2} = h; |
|
|
181 varargout{1} = c; |
|
|
182 endif |
|
|
183 |
|
|
184 endfunction |
|
|
185 |
|
|
186 function [X, Y, Z, lvl, ax, patch_props] = parse_args (arg) |
|
|
187 |
|
|
188 patch_props = {}; |
|
|
189 nolvl = false; |
|
|
190 |
|
|
191 if (isinteger (arg{1}) && ishandle (arg{1}) |
|
|
192 && strncmpi (get (arg{1}, "type"), "axis", 4)) |
|
|
193 ax = arg{1}; |
|
|
194 arg{1} = []; |
|
|
195 else |
|
|
196 ax = gca (); |
|
|
197 endif |
|
|
198 |
|
|
199 for n = 1:numel (arg) |
|
|
200 if (ischar (arg{n})) |
|
|
201 patch_props = arg(n:end); |
|
|
202 arg(n:end) = []; |
|
|
203 break; |
|
|
204 endif |
|
|
205 endfor |
|
|
206 |
|
|
207 if (mod (numel (patch_props), 2) != 0) |
|
|
208 error ("patch: property value is missing"); |
|
|
209 endif |
|
|
210 |
|
|
211 if (numel (arg) < 3) |
|
|
212 Z = arg{1}; |
|
|
213 [X, Y] = meshgrid (1:columns (Z), 1:rows (Z)); |
|
|
214 endif |
|
|
215 |
|
|
216 if (numel (arg) == 1) |
|
|
217 nolvl = true; |
|
|
218 arg(1) = []; |
|
|
219 elseif (numel (arg) == 2) |
|
|
220 lvl = arg{2}; |
|
|
221 arg(1:2) = []; |
|
|
222 elseif (numel (arg) == 3) |
|
|
223 arg{1:3}; |
|
|
224 [X, Y, Z] = deal (arg{1:3}); |
|
|
225 arg(1:3) = []; |
|
|
226 nolvl = true; |
|
|
227 elseif (numel (arg) == 4) |
|
|
228 [X, Y, Z, lvl] = deal (arg{1:4}); |
|
|
229 arg(1:4) = []; |
|
|
230 endif |
|
|
231 |
|
|
232 if (any (size (X) != size (Y))) |
|
|
233 error ("patch: X and Y must be of same size") |
|
|
234 endif |
|
|
235 |
|
|
236 if (isvector (X) || isvector (Y)) |
|
|
237 [X, Y] = meshgrid (X, Y); |
|
|
238 endif |
|
|
239 |
|
|
240 Z_no_inf = Z(!isinf (Z)); |
|
|
241 [minz, maxz] = deal (min (Z_no_inf(:)), max (Z_no_inf(:))); |
|
|
242 Z(isnan (Z)) = -Inf; |
|
|
243 |
|
|
244 if (nolvl) |
|
|
245 lvl = linspace (minz, maxz, 12); |
|
|
246 endif |
|
|
247 |
|
|
248 if (isscalar (lvl)) |
|
|
249 lvl = linspace (minz, maxz, lvl + 2)(1:end-1); |
|
|
250 else |
|
|
251 idx1 = find(lvl < minz); |
|
|
252 idx2 = find(lvl > maxz); |
|
|
253 lvl(idx1(1:end-1)) = []; |
|
|
254 lvl(idx2) = []; |
|
|
255 if (isempty (lvl)) |
|
|
256 lvl = [minz, minz]; |
|
|
257 endif |
|
|
258 endif |
|
|
259 |
|
|
260 endfunction |
