## Copyright (C) 2007, 2008, 2009 David Bateman

##

## This file is part of Octave.

##

## Octave 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 3 of the License, or (at

## your option) any later version.

##

## Octave 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 Octave; see the file COPYING.  If not, see

## <http://www.gnu.org/licenses/>.

## -*- texinfo -*-

## @deftypefn {Function File} {[@var{c}, @var{hg}] =} __contour__ (@dots{})

## Undocumented internal function.

## @end deftypefn

function [c, hg] = __contour__ (varargin)

  ax = varargin{1};

  zlevel = varargin{2};

  filled = "off";

  linespec.linestyle = "-";

  linespec.color = "auto";

  edgecolor = "flat";

  for i = 3 : nargin

    arg = varargin {i};

    if ((ischar (arg) || iscell (arg)))

      [linespec, valid] = __pltopt__ ("__contour__", arg, false);

      if (isempty (linespec.color))

	linespec.color = "auto";

      endif

      if (valid)

	have_line_spec = true;

	varargin(i) = [];

	break;

      endif

    endif

  endfor

  opts = {};

  i = 3;

  while (i < length (varargin))

    if (ischar (varargin {i}))

      if (strcmpi (varargin{i}, "fill"))

	filled = varargin {i + 1};

	varargin(i:i+1) = [];

      elseif (strcmpi (varargin{i}, "linecolor"))

	linespec.color = varargin {i + 1};

	edgecolor = linespec.color;

	if (ischar (edgecolor) && strcmpi (edgecolor, "auto"))

	  edgecolor = "flat";

	endif

	varargin(i:i+1) = [];

      elseif (strcmpi (varargin{i}, "edgecolor"))

	linespec.color = varargin {i + 1};

	edgecolor = linespec.color;

	if (ischar (edgecolor) && strcmpi (edgecolor, "flat"))

	  linespec.color = "auto";

	endif

	varargin(i:i+1) = [];

      else

	opts{end+1} = varargin{i};

        varargin(i) = [];

        opts{end+1} = varargin{i};

        varargin(i) = [];

      endif

    else

      i++;

    endif

  endwhile

  if (length(varargin) < 5)

    z1 = varargin{3};

    x1 = 1 : rows(z1);

    y1 = 1 : columns(z1);

  else

    x1 = varargin{3};

    y1 = varargin{4};

    z1 = varargin{5};

  endif

  if (length (varargin) == 4 || length (varargin) == 6)

    vn = varargin {end};

    vnauto = false;

  else

    vnauto = true;

    vn = 10; 

  endif

  if (isscalar (vn))

    lvl = linspace (min (z1(!isinf(z1))), max (z1(!isinf(z1))), 

		    vn + 2)(1:end-1);

  else

    lvl = vn;

  endif

  if (strcmpi (filled, "on"))

    if (isvector (x1) || isvector (y1))

      [x1, y1] = meshgrid (x1, y1);

    endif

    [nr, nc] = size (z1);

    x0 = prepad(x1, nc+1, 2 * x1(1, 1) - x1(1, 2), 2);

    x0 = postpad(x0, nc+2, 2 * x1(1, nc) - x1(1, nc - 1), 2);

    x0 = [x0(1, :); x0; x0(1, :)];

    y0 = prepad(y1, nr+1, 2 * y1(1, 1) - y1(2, 1), 1);

    y0 = postpad(y0, nr+2, 2 * y1(nr, 1) - y1(nr - 1, 1));

    y0 = [y0(:, 1), y0, y0(:, 1)];

    z0 = -Inf(nr+2, nc+2);

    z0(2:nr+1, 2:nc+1) = z1;

    [c, lev] = contourc (x0, y0, z0, lvl);

  else

    [c, lev] = contourc (x1, y1, z1, lvl);

  endif

  hg = hggroup ();

  opts = __add_datasource__ ("__countour__", hg, {"x", "y", "z"}, opts{:});

  addproperty ("xdata", hg, "data", x1);

  addproperty ("ydata", hg, "data", y1);

  addproperty ("zdata", hg, "data", z1);

  addproperty ("contourmatrix", hg, "data", c);

  addlistener (hg, "xdata", @update_data);

  addlistener (hg, "ydata", @update_data);

  addlistener (hg, "zdata", @update_data);

  addlistener (hg, "contourmatrix", @update_data);

  addproperty ("fill", hg, "radio", "on|{off}", filled);

  ## The properties zlevel and zlevelmode don't exist in matlab, but

  ## allow the use of contourgroups with the contour3, meshc and surfc 

  ## functions. 

  if (isnumeric (zlevel))

    addproperty ("zlevelmode", hg, "radio", "{none}|auto|manual", "manual")

    addproperty ("zlevel", hg, "data", zlevel);

  else

    addproperty ("zlevelmode", hg, "radio", "{none}|auto|manual", zlevel)

    if (ischar (zlevel) && strcmpi (zlevel, "manual"))

      z = varargin{3};

      z = 2 * (min (z(:)) - max (z(:)));

      addproperty ("zlevel", hg, "data", z);

    else

      addproperty ("zlevel", hg, "data", 0.);

    endif

  endif

  lvlstep = sum (abs (diff (lvl))) / (length (lvl) - 1);

  addproperty ("levellist", hg, "data", lev);

  addproperty ("levelstep", hg, "double", lvlstep);

  if (vnauto)

    addproperty ("levellistmode", hg, "radio", "{auto}|manual", "auto");

    addproperty ("levelstepmode", hg, "radio", "{auto}|manual", "auto");

  elseif (isscalar (vn))

    addproperty ("levellistmode", hg, "radio", "{auto}|manual", "auto");

    addproperty ("levelstepmode", hg, "radio", "{auto}|manual", "manual");

  else

    addproperty ("levellistmode", hg, "radio", "{auto}|manual", "manual");

    addproperty ("levelstepmode", hg, "radio", "{auto}|manual", "auto");

  endif

  addproperty ("labelspacing", hg, "double", 144);

  addproperty ("textlist", hg, "data", lev);

  addproperty ("textlistmode", hg, "radio", "{auto}|manual", "auto");

  addproperty ("textstep", hg, "double", lvlstep);

  addproperty ("textstepmode", hg, "radio", "{auto}|manual", "auto");

  addproperty ("showtext", hg, "radio", "on|{off}", "off");

  addproperty ("linecolor", hg, "color", linespec.color, "{auto}|none");

  addproperty ("linestyle", hg, "linelinestyle", linespec.linestyle);

  addproperty ("linewidth", hg, "linelinewidth", 0.5);

  ## FIXME It would be good to hide this property which is just an undocumented

  ## alias for linecolor

  addproperty ("edgecolor", hg, "color", edgecolor, "{flat}|none");

  addlistener (hg, "fill", @update_data);

  addlistener (hg, "zlevelmode", @update_zlevel);

  addlistener (hg, "zlevel", @update_zlevel);

  addlistener (hg, "levellist", @update_data);

  addlistener (hg, "levelstep", @update_data);

  addlistener (hg, "levellistmode", @update_data);

  addlistener (hg, "levelstepmode", @update_data);

  addlistener (hg, "labelspacing", @update_text);

  addlistener (hg, "textlist", @update_text);

  addlistener (hg, "textlistmode", @update_text);

  addlistener (hg, "textstep", @update_text);

  addlistener (hg, "textstepmode", @update_text);

  addlistener (hg, "showtext", @update_text);

  addlistener (hg, "linecolor", @update_line);

  addlistener (hg, "linestyle", @update_line);

  addlistener (hg, "linewidth", @update_line);

  addlistener (hg, "edgecolor", @update_edgecolor);

  add_patch_children (hg);

  if (!isempty (opts))

    set (hg, opts{:});

  endif

endfunction

function add_patch_children (hg)

  c = get (hg, "contourmatrix");

  lev = get (hg, "levellist");

  fill = get (hg, "fill");

  z = get (hg, "zlevel");

  zmode = get (hg, "zlevelmode");

  lc = get (hg, "linecolor");

  lw = get (hg, "linewidth");

  ls = get (hg, "linestyle");

  filled = get (hg, "fill");

  if (strcmpi (lc, "auto"))

    lc = "flat";

  endif

  if (strcmpi (filled, "on"))

    if (diff (lev) < 10*eps) 

      lvl_eps = 1e-6;

    else

      lvl_eps = min (diff (lev)) / 1000.0;

    endif

    ## 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);

      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 = [];

    for n = len:(-1):1

      idx = svec(n);

      ctmp = c(:, cont_idx(idx):cont_idx(idx) + cont_len(idx) - 1);

      if (all (ctmp(:,1) == ctmp(:,end)))

        ctmp(:, end) = [];

      else

	## Special case unclosed contours

      endif

      h = [h; patch(ctmp(1, :), ctmp(2, :), cont_lev(idx), "edgecolor", lc, 

		    "linestyle", ls, "linewidth", lw, "parent", hg)];

    endfor

    if (min (lev) == max (lev))

      set (gca (), "clim", [min(lev)-1, max(lev)+1]);

    else

      set (gca(), "clim", [min(lev), max(lev)]);

    endif

    set (gca (), "layer", "top");

  else

    ## Decode contourc output format.

    i1 = 1;

    h = [];

    while (i1 < length (c))

      clev = c(1,i1);

      clen = c(2,i1);

      if (all (c(:,i1+1) == c(:,i1+clen)))

	p = c(:, i1+1:i1+clen-1);

      else

	p = [c(:, i1+1:i1+clen), NaN(2, 1)];

      endif

      switch (zmode)

	case "none"

	  h = [h; patch(p(1,:), p(2,:), "facecolor", "none", 

			"edgecolor", lc, "linestyle", ls, "linewidth", lw,

			"cdata", clev, "parent", hg)]; 

	case "auto"

	  h = [h; patch(p(1,:), p(2,:), clev * ones (1, columns (p)),

			"facecolor", "none", "edgecolor", lc, 

			"linestyle", ls, "linewidth", lw, "cdata", clev, 

			"parent", hg)];

	otherwise

	  h = [h; patch(p(1,:), p(2,:), z * ones (1, columns (p)),

			"facecolor", "none", "edgecolor", lc,

			"linestyle", ls, "linewidth", lw, "cdata", clev,

			"parent", hg)];

      endswitch

      i1 += clen+1;

    endwhile

  endif

endfunction

function update_zlevel (h, d)

  z = get (h, "zlevel");

  zmode = get (h, "zlevelmode");

  kids = get (h, "children");

  switch (zmode)

    case "none"

      set (kids, "zdata", []);

    case "auto"

      for i = 1 : length (kids)

	set (kids(i), "zdata", get (kids (i), "cdata") .* 

	     ones (size (get (kids (i), "xdata"))));

      endfor

    otherwise

      for i = 1 : length (kids)

	set (kids(i), "zdata", z .* ones (size (get (kids (i), "xdata"))));

      endfor

  endswitch

endfunction

function update_edgecolor (h, d)

  ec = get (h, "edgecolor");

  lc = get (h, "linecolor");

  if (ischar (ec) && strcmpi (ec, "flat"))

    if (! strcmpi (lc, "auto"))

      set (h, "linecolor", "auto");

    endif

  elseif (! isequal (ec, lc))

    set (h, "linecolor", ec);

  endif

endfunction

function update_line (h, d)

  lc = get (h, "linecolor");

  ec = get (h, "edgecolor");

  if (strcmpi (lc, "auto"))

    lc = "flat";

  endif

  if (! isequal (ec, lc))

    set (h, "edgecolor", lc);

  endif

  set (findobj (h, "type", "patch"), "edgecolor", lc,

       "linewidth", get (h, "linewidth"), "linestyle", get (h, "linestyle"));

endfunction

function update_data (h, d)

  persistent recursive = false;

  if (!recursive)

    recursive = true;

    delete (get (h, "children"));

    if (strcmpi (get (h, "levellistmode"), "manual"))

      lvl = get (h, "levellist");

    elseif (strcmpi (get (h, "levelstepmode"), "manual"))

      z = get (h, "zdata");

      lvl = ceil ((max(z(:)) - min (z(:)) ./ get (h, "levelstep")));

    else

      lvl = 10;

    endif

    if (strcmpi (get (h, "fill"), "on"))

      X = get (h, "xdata");

      Y = get (h, "ydata");

      Z = get (h, "zdata");

      if (isvector (X) || isvector (Y))

	[X, Y] = meshgrid (X, Y);

      endif

      [nr, nc] = size (Z);

      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);

    else

      [c, lev] = contourc (get (h, "xdata"), get (h, "ydata"), 

			   get (h, "zdata"), lvl); 

    endif

    set (h, "contourmatrix", c);

    if (strcmpi (get (h, "levellistmode"), "manual"))

      ## Do nothing

    elseif (strcmpi (get (h, "levelstepmode"), "manual"))

      set (h, "levellist", lev);

    else

      set (h, "levellist", lev);

      z = get (h, "zdata");

      lvlstep = (max(z(:)) - min(z(:))) / 10;

      set (h, "levelstep", lvlstep);

    endif

    add_patch_children (h);

    update_text (h, d);

  endif

  recursive = false;

endfunction

function update_text (h, d)

  persistent recursive = false;

  if (!recursive)

    recursive = true;

    delete (findobj (h, "type", "text"));

    if (strcmpi (get (h, "textlistmode"), "manual"))

      lvl = get (h, "textlist");

    elseif (strcmpi (get (h, "textstepmode"), "manual"))

      lev = get (h, "levellist");

      if (diff (lev) < 10*eps) 

	lvl_eps = 1e-6;

      else

	lvl_eps = min (abs (diff (lev))) / 1000.0;

      endif

      stp = get (h, "textstep");

      t = [0, floor(cumsum(diff (lev)) / (abs(stp) - lvl_eps))];

      lvl = lev([true, t(1:end-1) != t(2:end)]);

      set (h, "textlist", lvl);

    else

      lvl = get (h, "levellist");

      set (h, "textlist", lvl, "textstep", get (h, "levelstep"));

    endif

    if (strcmpi (get (h, "showtext"), "on"))

      switch (get (h, "zlevelmode"))

	case "manual"

	  __clabel__ (get (h, "contourmatrix"), lvl, h, 

		      get (h, "labelspacing"), get (h, "zlevel"));

	case "auto"

	  __clabel__ (get (h, "contourmatrix"), lvl, h,

		      get (h, "labelspacing"), "auto");

	otherwise

	  __clabel__ (get (h, "contourmatrix"), lvl, h,

		      get (h, "labelspacing"), []);

      endswitch

    endif

  endif

  recursive = false;

endfunction