scilab/modules/graphics/macros/contourf.sci

   1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   2 // Copyright (C) INRIA
   3 // Copyright (C) DIGITEO - 2012 - Allan CORNET
   4 // Copyright (C) 2012 - 2016 - Scilab Enterprises
   5 //
   6 // This file is hereby licensed under the terms of the GNU GPL v2.0,
   7 // pursuant to article 5.3.4 of the CeCILL v.2.1.
   8 // This file was originally licensed under the terms of the CeCILL v2.1,
   9 // and continues to be available under such terms.
  10 // For more information, see the COPYING file which you should have received
  11 // along with this program.
  12
  13 function contourf(x, y, z, nv, style, strf, leg, rect, nax)
  14
  15     [nout, nin] = argn(0);
  16
  17     if nin == 0 then   // demo
  18         t = -%pi:0.1:%pi;
  19         m = sin(t)' * cos(t);
  20         contourf(t,t,m);
  21         return;
  22     end
  23
  24     if nin <= 0 then
  25         x=1:10;
  26     end
  27     if nin <= 1 then
  28         y=1:10;
  29     end
  30     if nin <= 2 then
  31         z=rand(size(x,"*"), size(y,"*"));
  32     end
  33     if nin <= 3 then
  34         zmin=min(z);
  35         zmax=max(z);
  36         nv = zmin + (1:10) * (zmax-zmin)/(11);
  37     end
  38     if nin <= 5 then
  39         strf="021";
  40     end
  41     if nin <= 6 then
  42         leg=" ";
  43     end
  44     if nin <= 7 then
  45         rect=[0,0,1,1];
  46     end
  47     if nin <= 8 then
  48         nax=[1,10,1,10];
  49     end
  50     if x==[] then
  51         x=1:size(z,"r");
  52     end
  53     if y==[] then
  54         y=1:size(z,"c");
  55     end
  56
  57     if type(x) <> 1 then
  58         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real vector expected.\n"), "contourf", 1));
  59     end
  60
  61     if type(y) <> 1 then
  62         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real vector expected.\n"), "contourf", 2));
  63     end
  64
  65     if type(z) <> 1 then
  66         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real matrix expected.\n"), "contourf", 3));
  67     end
  68
  69     if type(nv) <> 1 then
  70         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real matrix expected.\n"), "contourf", 4));
  71     end
  72
  73     if type(strf) <> 10 then
  74         error(msprintf(gettext("%s: Wrong type for input argument #%d: String expected.\n"), "contourf", 6));
  75     end
  76
  77     if type(leg) <> 10 then
  78         error(msprintf(gettext("%s: Wrong type for input argument #%d: String expected.\n"), "contourf", 7));
  79     end
  80
  81     if type(rect) <> 1 then
  82         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real matrix expected.\n"), "contourf", 8));
  83     end
  84
  85     if type(nax) <> 1 then
  86         error(msprintf(gettext("%s: Wrong type for input argument #%d:  Real matrix expected.\n"), "contourf", 9));
  87     end
  88
  89
  90     nvs=size(nv,"*") ;
  91     if nvs==1 then
  92         nvs=nv;
  93         zmin=min(z);
  94         zmax=max(z);
  95         nv = zmin + (1:nvs)*(zmax-zmin)/(nvs+1);
  96     end;
  97
  98     if nin <= 4 then
  99         style = -1*ones(1, nvs);
 100     end
 101
 102     if type(style) <> 1 then
 103         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real matrix expected.\n"), "contourf", 5));
 104     end
 105
 106     if nin <= 7 then
 107         rect=[min(x), min(y), max(x), max(y)];
 108     end
 109
 110     if ~isvector(x) then
 111         error(msprintf(gettext("%s: Wrong size for input argument #%d: Real vector expected.\n"), "contourf", 1));
 112     end
 113
 114     if ~isvector(y) then
 115         error(msprintf(gettext("%s: Wrong size for input argument #%d: Real vector expected.\n"), "contourf", 2));
 116     end
 117
 118     if size(strf, "*") <> 1 then
 119         error(msprintf(gettext("%s: Wrong size for input argument #%d: string expected.\n"), "contourf", 6));
 120     end
 121
 122     if size(leg, "*") <> 1 then
 123         error(msprintf(gettext("%s: Wrong size for input argument #%d: string expected.\n"), "contourf", 7));
 124     end
 125
 126     nv1 = nv;
 127     [mz,nz] = size(z);
 128     minz = min(z);
 129     maxz = max(z);
 130
 131     // Surround the matrix by a very low region to get closed contours, and
 132     // replace any NaN with low numbers as well.
 133     zz=[ %nan * ones(1,nz+2) + %nan;
 134     %nan * ones(mz,1)   + %nan, z, %nan * ones(mz,1)   + %nan;
 135     %nan * ones(1,nz+2) + %nan];
 136
 137     kk=find(isnan(zz(:)));
 138
 139     zz(kk)=minz-1e4*(maxz-minz)+zeros(kk);
 140
 141     xx = [2 * x(1) - x(2); x(:); 2 * x(mz) - x(mz - 1)];
 142     yy = [2 * y(1) - y(2); y(:); 2 * y(nz) - y(nz - 1)];
 143
 144     // Internal call to get the contours
 145     [x1,y1]=contour2di(xx,yy,zz,nv);
 146     CS=[x1;y1];
 147     // Find the indices of the curves in the c matrix, and get the
 148     // area of closed curves in order to draw patches correctly.
 149     ii = 1;
 150     ncurves = 0;
 151     I = [];
 152     Area=[];
 153
 154     while (ii < size(CS,2)),
 155         nl=CS(2,ii);
 156         ncurves = ncurves + 1;
 157         I(ncurves) = ii;
 158         xp=CS(1,ii+(1:nl));  // First patch
 159         yp=CS(2,ii+(1:nl));
 160         Area(ncurves)=sum(diff(xp).*(yp(1:nl-1)+yp(2:nl))/2);
 161         ii = ii + nl + 1;
 162     end
 163
 164     lp = size(gcf().color_map, 1);
 165
 166     if size(nv,"*") > 1 // case where nv is a vector defining the level curve values
 167         if  size(nv,"*") > lp
 168             error(msprintf(gettext("%s: Colormap too small"),"contourf"));
 169         end
 170     else
 171         if nv > lp
 172             error(msprintf(gettext("%s: Colormap too small"),"contourf"));
 173             return ;
 174         end
 175     end
 176
 177     min_nv=min(nv);
 178     max_nv=max(nv);
 179
 180     initDrawingMode = gcf().immediate_drawing;
 181     gcf().immediate_drawing = "off";
 182
 183     plot2d([min(xx);max(xx)],[min(yy);max(yy)],0,strf,leg,rect,nax);
 184
 185     // Plot patches in order of decreasing size. This makes sure that
 186     // all the lev1es get drawn, not matter if we are going up a hill or
 187     // down into a hole. When going down we shift levels though, you can
 188     // tell whether we are going up or down by checking the sign of the
 189     // area (since curves are oriented so that the high side is always
 190     // the same side). Lowest curve is largest and encloses higher data
 191     // always.
 192
 193     draw_min=1;
 194     H=[];
 195     [FA,IA]=gsort(abs(Area));
 196
 197     ax = gca();
 198     old_foreground = ax.foreground;
 199     pat = old_foreground;
 200     for jj=IA',
 201         nl=CS(2,I(jj));
 202         lev1=CS(1,I(jj));
 203         if (lev1 ~= minz | draw_min) then
 204             xp=CS(1,I(jj)+(1:nl));
 205             yp=CS(2,I(jj)+(1:nl));
 206             pat = size(find( nv <= lev1),"*");
 207             ax.foreground = pat;
 208             xfpoly(xp,yp)
 209         end
 210     end
 211
 212     if style(1)<>-1 then
 213         contour2d(xx,yy,zz,nv,style,"000",leg,rect,nax);
 214     end
 215     ax.foreground = old_foreground;
 216
 217     gcf().immediate_drawing = initDrawingMode;
 218 endfunction
 219