# HG changeset patch # User jwe # Date 1192826613 0 # Node ID e54cc03d53f6010627309e173cd2d37c5e94cf9c # Parent 32da63705bf9d968a3880ffe477bcf345b946e8c [project @ 2007-10-19 20:43:32 by jwe] diff -r 32da63705bf9 -r e54cc03d53f6 scripts/ChangeLog --- a/scripts/ChangeLog Fri Oct 19 19:05:31 2007 +0000 +++ b/scripts/ChangeLog Fri Oct 19 20:43:33 2007 +0000 @@ -1,3 +1,8 @@ +2007-10-19 Kai Habel + + * plot/contourf.m: New function. + * plot/Makefile.in (SOURCES): Add it to the list. + 2007-10-19 John W. Eaton * plot/subplot.m: Doc fix. diff -r 32da63705bf9 -r e54cc03d53f6 scripts/plot/Makefile.in --- a/scripts/plot/Makefile.in Fri Oct 19 19:05:31 2007 +0000 +++ b/scripts/plot/Makefile.in Fri Oct 19 20:43:33 2007 +0000 @@ -70,6 +70,7 @@ closereq.m \ contour.m \ contourc.m \ + contourf.m \ drawnow.m \ errorbar.m \ figure.m \ diff -r 32da63705bf9 -r e54cc03d53f6 scripts/plot/contourf.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/scripts/plot/contourf.m Fri Oct 19 20:43:33 2007 +0000 @@ -0,0 +1,259 @@ +## Copyright (C) 2007 Kai Habel +## Copyright (C) 2003 Shai Ayal +## +## This program is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 2, or (at your option) +## any later version. +## +## OctPlot is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with OctPlot; see the file COPYING. If not, write to the Free +## Software Foundation, 59 Temple Place - Suite 330, Boston, MA +## 02111-1307, USA. + +## -*- texinfo -*- +## @deftypefn {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}, @var{lvl}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}, @var{n}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{x}, @var{y}, @var{z}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}, @var{n}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}, @var{lvl}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{z}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@var{ax}, @dots{}) +## @deftypefnx {Function File} {[@var{c}, @var{h}] =} contourf (@dots{}, @var{"property"}, @var{val}) +## Compute and plot filled contours of the matrix @var{z}. +## Parameters @var{x}, @var{y} and @var{n} or @var{lvl} are optional. +## +## The return value @var{c} is a 2xn matrix containing the contour lines +## as described in the help to the contourc function. +## +## The return value @var{h} is handle-vector to the patch objects creating +## the filled contours. +## +## If @var{x} and @var{y} are ommited they are taken as the row/column +## index of @var{z}. @var{n} is a scalar denoting the number of lines +## to compute. Alternatively @var{lvl} is a vector containing the +## contour levels. If only one value (e.g. lvl0) is wanted, set +## @var{lvl} to [lvl0, lvl0]. If both @var{n} or @var{lvl} are omitted +## a default value of 10 contour level is assumed. +## +## If provided, the filled contours are added to the axes object +## @var{ax} instead of the current axis. +## +## The following example plots filled contours of the @code{peaks} +## function. +## @example +## [x, y, z] = peaks (50); +## contourf (x, y, z, -7:9) +## @end example +## @seealso{contour, contourc, patch} +## @end deftypefn + +## Author: Kai Habel +## Author: shaia + +function varargout = contourf (varargin) + + [X, Y, Z, lvl, ax, patch_props] = parse_args (varargin); + [nr, nc] = size (Z); + [minx, maxx] = deal (min (X(:)), max (X(:))); + [miny, maxy] = deal (min (Y(:)), max (Y(:))); + + if (diff (lvl) < 10*eps) + lvl_eps = 1e-6; + else + lvl_eps = min (diff (lvl)) / 1000.0; + endif + + X0 = prepad(X, nc+1, 2 * X(1, 1) - X(1, 2), 2); + X0 = postpad(X0, nc+2, 2 * X(1, nc) - X(1, nc - 1), 2); + X0 = [X0(1, :); X0; X0(1, :)]; + Y0 = prepad(Y, nr+1, 2 * Y(1, 1) - Y(2, 1), 1); + Y0 = postpad(Y0, nr+2, 2 * Y(nr, 1) - Y(nr - 1, 1)); + Y0 = [Y0(:, 1), Y0, Y0(:, 1)]; + + Z0 = -Inf(nr+2, nc+2); + Z0(2:nr+1, 2:nc+1) = Z; + [c, lev] = contourc (X0, Y0, Z0, lvl); + cmap = colormap (); + + levx = linspace (min (lev), max (lev), size (cmap, 1)); + + newplot (); + + ## Decode contourc output format. + i1 = 1; + ncont = 0; + while (i1 < columns (c)) + ncont++; + cont_lev(ncont) = c(1, i1); + cont_len(ncont) = c(2, i1); + cont_idx(ncont) = i1+1; + + ii = i1+1:i1+cont_len(ncont); + cur_cont = c(:, ii); + c(:, ii); + startidx = ii(1); + stopidx = ii(end); + cont_area(ncont) = polyarea (c(1, ii), c(2, ii)); + i1 += c(2, i1) + 1; + endwhile + + ## Handle for each level the case where we have (a) hole(s) in a patch. + ## Those are to be filled with the color of level below or with the + ## background colour. + for k = 1:numel (lev) + lvl_idx = find (abs (cont_lev - lev(k)) < lvl_eps); + len = numel (lvl_idx); + if (len > 1) + ## mark = logical(zeros(size(lvl_idx))); + mark = false (size (lvl_idx)); + a = 1; + while (a < len) + # take 1st patch + b = a + 1; + pa_idx = lvl_idx(a); + # get pointer to contour start, and contour length + curr_ct_idx = cont_idx(pa_idx); + curr_ct_len = cont_len(pa_idx); + # get contour + curr_ct = c(:, curr_ct_idx:curr_ct_idx+curr_ct_len-1); + b_vec = (a+1):len; + next_ct_pt_vec = c(:, cont_idx(lvl_idx(b_vec))); + in = inpolygon (next_ct_pt_vec(1,:), next_ct_pt_vec(2,:), + curr_ct(1, :), curr_ct(2, :)); + mark(b_vec(in)) = !mark(b_vec(in)); + a++; + endwhile + if (numel (mark) > 0) + ## All marked contours describe a hole in a larger contour of + ## the same level and must be filled with colour of level below. + ma_idx = lvl_idx(mark); + if (k > 1) + ## Find color of level below. + tmp = find(abs(cont_lev - lev(k - 1)) < lvl_eps); + lvl_bel_idx = tmp(1); + ## Set color of patches found. + cont_lev(ma_idx) = cont_lev(lvl_bel_idx); + else + ## Set lowest level contour to NaN. + cont_lev(ma_idx) = NaN; + endif + endif + endif + endfor + + ## The algorithm can create patches with the size of the plotting + ## area, we would like to draw only the patch with the highest level. + del_idx = []; + max_idx = find (cont_area == max (cont_area)); + if (numel (max_idx) > 1) + # delete double entries + del_idx = max_idx(1:end-1); + cont_area(del_idx) = cont_lev(del_idx) = []; + cont_len(del_idx) = cont_idx(del_idx) = []; + endif + + ## Now we have everything together and can start plotting the patches + ## beginning with largest area. + [tmp, svec] = sort (cont_area); + len = ncont - numel (del_idx); + h = zeros (1, len); + for n = len:-1:1 + idx = svec(n); + ii = cont_idx(idx):cont_idx(idx) + cont_len(idx) - 2; + h(n) = patch (c(1, ii), c(2, ii), cont_lev(idx), patch_props{:}); + endfor + + if (min (lev) == max (lev)) + set (gca (), "clim", [min(lev)-1, max(lev)+1]); + else + set (gca(), "clim", [min(lev), max(lev)]); + endif + + if (nargout > 0) + varargout{2} = h; + varargout{1} = c; + endif + +endfunction + +function [X, Y, Z, lvl, ax, patch_props] = parse_args (arg) + + patch_props = {}; + nolvl = false; + + if (isinteger (arg{1}) && ishandle (arg{1}) + && strncmpi (get (arg{1}, "type"), "axis", 4)) + ax = arg{1}; + arg{1} = []; + else + ax = gca (); + endif + + for n = 1:numel (arg) + if (ischar (arg{n})) + patch_props = arg(n:end); + arg(n:end) = []; + break; + endif + endfor + + if (mod (numel (patch_props), 2) != 0) + error ("patch: property value is missing"); + endif + + if (numel (arg) < 3) + Z = arg{1}; + [X, Y] = meshgrid (1:columns (Z), 1:rows (Z)); + endif + + if (numel (arg) == 1) + nolvl = true; + arg(1) = []; + elseif (numel (arg) == 2) + lvl = arg{2}; + arg(1:2) = []; + elseif (numel (arg) == 3) + arg{1:3}; + [X, Y, Z] = deal (arg{1:3}); + arg(1:3) = []; + nolvl = true; + elseif (numel (arg) == 4) + [X, Y, Z, lvl] = deal (arg{1:4}); + arg(1:4) = []; + endif + + if (any (size (X) != size (Y))) + error ("patch: X and Y must be of same size") + endif + + if (isvector (X) || isvector (Y)) + [X, Y] = meshgrid (X, Y); + endif + + Z_no_inf = Z(!isinf (Z)); + [minz, maxz] = deal (min (Z_no_inf(:)), max (Z_no_inf(:))); + Z(isnan (Z)) = -Inf; + + if (nolvl) + lvl = linspace (minz, maxz, 12); + endif + + if (isscalar (lvl)) + lvl = linspace (minz, maxz, lvl + 2)(1:end-1); + else + idx1 = find(lvl < minz); + idx2 = find(lvl > maxz); + lvl(idx1(1:end-1)) = []; + lvl(idx2) = []; + if (isempty (lvl)) + lvl = [minz, minz]; + endif + endif + +endfunction diff -r 32da63705bf9 -r e54cc03d53f6 src/ChangeLog --- a/src/ChangeLog Fri Oct 19 19:05:31 2007 +0000 +++ b/src/ChangeLog Fri Oct 19 20:43:33 2007 +0000 @@ -1,3 +1,14 @@ +2007-10-19 Kai Habel + + * DLD-FUNCTIONS/__contourc__.cc (add_point): Rename from + cl_add_point. Change all uses. + (end_contour): Rename from cl_end_contour. Change all uses. + (start_contour): Rename from cl_start_contour. Change all uses. + (drawcn): Rename from cl_drawcn. New algorithm for locating contours. + (mark_facets): New function. + (cntr): Rename from cl_cntr. Change all uses. New algorithm for + locating contours. + 2007-10-19 John W. Eaton * pt-idx.cc (tree_index_expression::lvalue): Correctly compute diff -r 32da63705bf9 -r e54cc03d53f6 src/DLD-FUNCTIONS/__contourc__.cc --- a/src/DLD-FUNCTIONS/__contourc__.cc Fri Oct 19 19:05:31 2007 +0000 +++ b/src/DLD-FUNCTIONS/__contourc__.cc Fri Oct 19 20:43:33 2007 +0000 @@ -1,5 +1,6 @@ /* Contour lines for function evaluated on a grid. +Copyright (C) 2007 Kai Habel Copyright (C) 2004, 2007 Shai Ayal Adapted to an oct file from the stand alone contourl by Victro Munoz @@ -35,6 +36,8 @@ #include #endif +#include + #include "quit.h" #include "defun-dld.h" @@ -45,15 +48,13 @@ static Matrix contourc; static int elem; -// this is the quanta in which we increase this_contour +// This is the quanta in which we increase this_contour. #define CONTOUR_QUANT 50 -// cl_add_point(x,y); -// -// Add a coordinate point (x,y) to this_contour +// Add a coordinate point (x,y) to this_contour. static void -cl_add_point (double x, double y) +add_point (double x, double y) { if (elem % CONTOUR_QUANT == 0) this_contour = this_contour.append (Matrix (2, CONTOUR_QUANT, 0)); @@ -63,237 +64,238 @@ elem++; } -// cl_end_contour(); -// -// Adds contents of current contour to contourc. +// Add contents of current contour to contourc. // this_contour.cols () - 1; static void -cl_end_contour (void) +end_contour (void) { if (elem > 2) { this_contour (1, 0) = elem - 1; contourc = contourc.append (this_contour.extract_n (0, 0, 2, elem)); } + this_contour = Matrix (); elem = 0; } -// cl_start_contour(flev,x,y); -// -// Start a new contour, and adds contents of current one to contourc +// Start a new contour, and add contents of current one to contourc. static void -cl_start_contour (double flev, double x, double y) +start_contour (double lvl, double x, double y) { - cl_end_contour (); + end_contour (); this_contour.resize (2, 0); - cl_add_point (flev, flev); - cl_add_point (x, y); + add_point (lvl, 0); + add_point (x, y); } static void -cl_drawcn (RowVector & X, RowVector & Y, Matrix & Z, double flev, - int krow, int kcol, double lastx, double lasty, int startedge, - Matrix & ipts) +drawcn (const RowVector& X, const RowVector& Y, const Matrix& Z, + double lvl, int r, int c, double ct_x, double ct_y, + uint start_edge, bool first, charMatrix& mark) { - - int kx = 0, lx = Z.cols () - 1, ky = 0, ly = Z.rows () - 1; - - double f[4]; - double px[4], py[4], locx[4], locy[4]; - int iedge[4]; - int num, first, inext, kcolnext, krownext; + double px[4], py[4], pz[4], tmp; + uint stop_edge, next_edge, pt[2]; + int next_r, next_c; - px[0] = X (krow + 1); - px[1] = X (krow); - px[2] = X (krow); - px[3] = X (krow + 1); - py[0] = Y (kcol); - py[1] = Y (kcol); - py[2] = Y (kcol + 1); - py[3] = Y (kcol + 1); + //get x, y, and z - lvl for current facet + px[0] = px[3] = X(c); + px[1] = px[2] = X(c+1); - f[0] = Z (krow + 1, kcol) - flev; - f[1] = Z (krow, kcol) - flev; - f[2] = Z (krow, kcol + 1) - flev; - f[3] = Z (krow + 1, kcol + 1) - flev; + py[0] = py[1] = Y(r); + py[2] = py[3] = Y(r+1); - for (int i = 0, j = 1; i < 4; i++, j = (j + 1) % 4) - { - iedge[i] = (f[i] * f[j] > 0.0) ? -1 : ((f[i] * f[j] < 0.0) ? 1 : 0); - } - - // Mark this square as done - ipts(krow,kcol) = 1; - - // Check if no contour has been crossed i.e. iedge[i] = -1 - if (iedge[0] == -1 && iedge[1] == -1 && iedge[2] == -1 && iedge[3] == -1) - return; + pz[3] = Z(r+1, c) - lvl; + pz[2] = Z(r+1, c + 1) - lvl; + pz[1] = Z(r, c+1) - lvl; + pz[0] = Z(r, c) - lvl; - // Check if this is a completely flat square - in which case ignore it - if (f[0] == 0.0 && f[1] == 0.0 && f[2] == 0.0 && f[3] == 0.0) - return; + // Facet edge and point naming assignment. + // + // 0-----1 .-0-. + // | | | | + // | | 3 1 + // | | | | + // 3-----2 .-2-. - // Calculate intersection points - num = 0; - if (startedge < 0) + // Get mark value of current facet. + char id = static_cast (mark(r, c)); + + // Check startedge s. + if (start_edge == 255) { - first = 1; - } - else - { - locx[num] = lastx; - locy[num] = lasty; - num++; - first = 0; + // Find start edge. + for (uint k = 0; k < 4; k++) + if (static_cast (pow(2, k)) & id) + start_edge = k; } - for (int k = 0, i = (startedge < 0 ? 0 : startedge); k < 4; - k++, i = (i + 1) % 4) - { - if (i == startedge) - continue; + if (start_edge == 255) + return; + + // Decrease mark value of current facet for start edge. + mark(r, c) -= static_cast (pow(2, start_edge)); - // If the contour is an edge check it hasn't already been done - if (f[i] == 0.0 && f[(i + 1) % 4] == 0.0) - { - kcolnext = kcol; - krownext = krow; - if (i == 0) - kcolnext--; - if (i == 1) - krownext--; - if (i == 2) - kcolnext++; - if (i == 3) - krownext++; - if (kcolnext < kx || kcolnext >= lx || krownext < ky - || krownext >= ly || ipts(krownext,kcolnext) == 1) - continue; - } - if (iedge[i] == 1 || f[i] == 0.0) + // Next point (clockwise). + pt[0] = start_edge; + pt[1] = (pt[0] + 1) % 4; + + // Calculate contour segment start if first of contour. + if (first) + { + tmp = fabs (pz[pt[1]]) / fabs (pz[pt[0]]); + + if (xisnan (tmp)) + ct_x = ct_y = 0.5; + else { - int j = (i + 1) % 4; - if (f[i] != 0.0) - { - locx[num] = - (px[i] * fabs (f[j]) + px[j] * fabs (f[i])) / fabs (f[j] - - f[i]); - locy[num] = - (py[i] * fabs (f[j]) + py[j] * fabs (f[i])) / fabs (f[j] - - f[i]); - } - else - { - locx[num] = px[i]; - locy[num] = py[i]; - } - // If this is the start of the contour then move to the point - if (first == 1) - { - cl_start_contour (flev, locx[num], locy[num]); - first = 0; - } - else - { - // Link to the next point on the contour - cl_add_point (locx[num], locy[num]); - // Need to follow contour into next grid box - // Easy case where contour does not pass through corner - if (f[i] != 0.0) - { - kcolnext = kcol; - krownext = krow; - inext = (i + 2) % 4; - if (i == 0) - kcolnext--; - if (i == 1) - krownext--; - if (i == 2) - kcolnext++; - if (i == 3) - krownext++; - if (kcolnext >= kx && kcolnext < lx - && krownext >= ky && krownext < ly - && ipts(krownext,kcolnext) == 0) - { - cl_drawcn (X, Y, Z, flev, krownext, kcolnext, - locx[num], locy[num], inext, ipts); - } - } - // Hard case where contour passes through corner. This - // is still not perfect - it may lose the contour which - // won't upset the contour itself (we can find it again - // later) but might upset the labelling (which is only - // relevant for the PLPlot implementation, since we - // don't worry about labels---for now!) - else - { - kcolnext = kcol; - krownext = krow; - inext = (i + 2) % 4; - if (i == 0) - { - kcolnext--; - krownext++; - } - if (i == 1) - { - krownext--; - kcolnext--; - } - if (i == 2) - { - kcolnext++; - krownext--; - } - if (i == 3) - { - krownext++; - kcolnext++; - } - if (kcolnext >= kx && kcolnext < lx - && krownext >= ky && krownext < ly - && ipts(krownext,kcolnext) == 0) - { - cl_drawcn (X, Y, Z, flev, krownext, kcolnext, - locx[num], locy[num], inext, ipts); - } + ct_x = px[pt[0]] + (px[pt[1]] - px[pt[0]])/(1 + tmp); + ct_y = py[pt[0]] + (py[pt[1]] - py[pt[0]])/(1 + tmp); + } + + start_contour (lvl, ct_x, ct_y); + } + + // Find stop edge FIXME: control flow --> while. + for (uint k = 1; k <= 4; k++) + { + if (start_edge == 0 || start_edge == 2) + stop_edge = (start_edge + k) % 4; + else + stop_edge = (start_edge - k) % 4; + + if (static_cast (pow(2, stop_edge)) & id) + break; + } + + pt[0] = stop_edge; + pt[1] = (pt[0] + 1) % 4; + tmp = fabs (pz[pt[1]]) / fabs (pz[pt[0]]); - } - if (first == 1) - { - // Move back to first point - cl_start_contour (flev, locx[num], locy[num]); - first = 0; - } - else - { - first = 1; - } - num++; - } - } + if (xisnan (tmp)) + ct_x = ct_y = 0.5; + else + { + ct_x = px[pt[0]] + (px[pt[1]] - px[pt[0]])/(1 + tmp); + ct_y = py[pt[0]] + (py[pt[1]] - py[pt[0]])/(1 + tmp); + } + + // Add point to contour. + add_point (ct_x, ct_y); + + // Decrease id value of current facet for start edge. + mark(r, c) -= static_cast(pow(2,stop_edge)); + + // Find next facet. + next_c = c; + next_r = r; + + if (stop_edge == 0) + next_r--; + else if (stop_edge == 1) + next_c++; + else if (stop_edge == 2) + next_r++; + else if (stop_edge == 3) + next_c--; + + // Check if next facet is not done yet. + // Go to next facet. + if (next_r >= 0 && next_c >= 0 && next_r < mark.rows () + && next_c < mark.cols () && mark(next_r, next_c) > 0) + { + next_edge = (stop_edge + 2) % 4; + drawcn (X, Y, Z, lvl, next_r, next_c, ct_x, ct_y, next_edge, false, mark); } } static void -cl_cntr (RowVector & X, RowVector & Y, Matrix & Z, double flev) +mark_facets (const Matrix& Z, charMatrix& mark, double lvl) { - Matrix ipts (Z.rows (), Z.cols (), 0); + uint nr = mark.rows (); + uint nc = mark.cols (); + + double f[4]; + + for (uint c = 0; c < nc; c++) + for (uint r = 0; r < nr; r++) + { + f[0] = Z(r, c) - lvl; + f[1] = Z(r, c+1) - lvl; + f[3] = Z(r+1, c) - lvl; + f[2] = Z(r+1, c+1) - lvl; + + for (uint i = 0; i < 4; i++) + if (fabs(f[i]) < DBL_EPSILON) + f[i] = DBL_EPSILON; + + if (f[1] * f[2] < 0) + mark(r, c) += 2; + + if (f[0] * f[3] < 0) + mark(r, c) += 8; + } - for (int krow = 0; krow < Z.rows () - 1; krow++) + for (uint r = 0; r < nr; r++) + for (uint c = 0; c < nc; c++) + { + f[0] = Z(r, c) - lvl; + f[1] = Z(r, c+1) - lvl; + f[3] = Z(r+1, c) - lvl; + f[2] = Z(r+1, c+1) - lvl; + + for (uint i = 0; i < 4; i++) + if (fabs(f[i]) < DBL_EPSILON) + f[i] = DBL_EPSILON; + + if (f[0] * f[1] < 0) + mark(r, c) += 1; + + if (f[2] * f[3] < 0) + mark(r, c) += 4; + } +} + +static void +cntr (const RowVector& X, const RowVector& Y, const Matrix& Z, double lvl) +{ + uint nr = Z.rows (); + uint nc = Z.cols (); + + charMatrix mark (nr - 1, nc - 1, 0); + + mark_facets (Z, mark, lvl); + + // Find contours that start at a domain edge. + + for (uint c = 0; c < nc - 1; c++) { - for (int kcol = 0; kcol < Z.cols () - 1; kcol++) - { - if (ipts(krow,kcol) == 0) - { - cl_drawcn (X, Y, Z, flev, krow, kcol, 0.0, 0.0, -2, ipts); - } - } + // Top. + if (mark(0, c) & 1) + drawcn (X, Y, Z, lvl, 0, c, 0.0, 0.0, 0, true, mark); + + // Bottom. + if (mark(nr - 2, c) & 4) + drawcn (X, Y, Z, lvl, nr - 2, c, 0.0, 0.0, 2, true, mark); } + + for (uint r = 0; r < nr - 1; r++) + { + // Left. + if (mark(r, 0) & 8) + drawcn (X, Y, Z, lvl, r, 0, 0.0, 0.0, 3, true, mark); + + // Right. + if (mark(r, nc - 2) & 2) + drawcn (X, Y, Z, lvl, r, nc - 2, 0.0, 0.0, 1, true, mark); + } + + for (uint r = 0; r < nr - 1; r++) + for (uint c = 0; c < nc - 1; c++) + if (mark (r, c) > 0) + drawcn (X, Y, Z, lvl, r, c, 0.0, 0.0, 255, true, mark); } DEFUN_DLD (__contourc__, args, , @@ -308,7 +310,7 @@ { RowVector X = args (0).row_vector_value (); RowVector Y = args (1).row_vector_value (); - Matrix Z = args (2).matrix_value ().transpose (); + Matrix Z = args (2).matrix_value (); RowVector L = args (3).row_vector_value (); if (! error_state) @@ -316,9 +318,9 @@ contourc.resize (2, 0); for (int i = 0; i < L.length (); i++) - cl_cntr (X, Y, Z, L (i)); + cntr (X, Y, Z, L (i)); - cl_end_contour (); + end_contour (); retval = contourc; }