scilab/modules/graphics/macros/histplot.sci

   1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   2 // Copyright (C) INRIA
   3 // Copyright (C) Bruno Pincon
   4 // Copyright (C) 2010 - Samuel Gougeon
   5 // Copyright (C) 2012 - Scilab Enterprises - Adeline CARNIS
   6 // Copyright (C) 2013 - A. Khorshidi (new option)
   7 // Copyright (C) 2013 - Scilab Enterprises - Paul Bignier: added output
   8 // Copyright (C) 2016 - Samuel GOUGEON
   9 //
  10 // Copyright (C) 2012 - 2016 - Scilab Enterprises
  11 //
  12 // This file is hereby licensed under the terms of the GNU GPL v2.0,
  13 // pursuant to article 5.3.4 of the CeCILL v.2.1.
  14 // This file was originally licensed under the terms of the CeCILL v2.1,
  15 // and continues to be available under such terms.
  16 // For more information, see the COPYING file which you should have received
  17 // along with this program.
  18
  19 function [y, ind] = histplot(n,data,style,strf,leg,rect,nax,logflag,frameflag,axesflag,normalization,polygon)
  20     // histplot(n,data,<opt_arg_seq>)
  21     // draws histogram of entries in  data put into n classes
  22     //
  23     // histplot(xi,data,<opt_arg_seq>)
  24     // generates the histogram of entries in data put into classes
  25     // [xi(1),xi(2)], (xi(k) xi(k+1)], k=2,..,n.
  26     // xi's are assumed st. increasing (this point is verified now).
  27     //
  28     // optional args:
  29     //     1/ the same than for a plot2d:
  30     //        style,strf,leg,rect,nax,logflag,frameflag,axesflag
  31     //     2/ normalization flag (default value %t). When true the
  32     //        histogram is normalized so that to approach a density:
  33     //          xmax
  34     //         /
  35     //         |  h(x) dx = 1   (true if xmin <= min(data) and max(data) <= xmax)
  36     //         /
  37     //         xmin
  38     //
  39     // Example : enter histplot()
  40     //
  41     // modifs to use dsearch (Bruno Pincon 10/12/2001)
  42     // others modifs from Bruno (feb 2005):
  43     //    - may be a cleaner and shorter way to deal with optional arg ?
  44     //    - now the histgram is drawn using only one polyline
  45     //      (so properties (color, thickness,...) are easier
  46     //       to change with new graphics).
  47     //    - removed computation of nax and rect if they are not
  48     //      passed (let plot2d doing it)
  49     //    - modify a little the demo
  50     //    - add some checking on n|x and data
  51     //
  52
  53     [lhs, rhs] = argn()
  54
  55     y = [];
  56     ind = [];
  57     if rhs == 0 then   // demo
  58         histplot([-4.5:0.25:4.5],rand(1,20000,"n"),style=2,axesflag=1,..
  59         frameflag=1,rect=[-4.5 0 4.5 0.47]);
  60         deff("[y]=f(x)","y=exp(-x.*x/2)/sqrt(2*%pi);");
  61         x=-4.5:0.125:4.5;
  62         x=x';
  63         plot2d(x,f(x),26,"000");
  64         titre= gettext("histplot() : (normalized) histogram plot");
  65         xtitle(titre,"C (Classes)","N(C) / (Nmax length(C))");  // Not clear
  66         legend(gettext("Gaussian random sample histogram"), ..
  67         gettext("Exact gaussian density"));
  68         return
  69     end
  70
  71     if rhs < 2
  72         error(msprintf(gettext("%s: Wrong number of input argument(s): At least %d expected.\n"),"histplot",2));
  73     end
  74
  75     if type(n) ~= 1 |  ~isreal(n)
  76         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real expected.\n"),"histplot",1));
  77     elseif type(data) ~= 1 | ~isreal(data)
  78         error(msprintf(gettext("%s: Wrong type for input argument #%d: Real expected.\n"),"histplot",2));
  79     end
  80
  81     // This is the only specific optional argument for histplot
  82     if exists("normalization","l")==0, normalization=%t,end
  83     if exists("polygon","l")==0, polygon=%f,end
  84
  85     // Now parse optional arguments to be sent to plot2d
  86     opt_arg_seq = "";
  87     opt_arg_list = ["style","strf","leg","rect","nax","logflag","frameflag","axesflag"]
  88     for opt_arg = opt_arg_list
  89         if exists(opt_arg,"local") then
  90             opt_arg_seq = opt_arg_seq +","+ opt_arg + "=" + opt_arg
  91         end
  92     end
  93
  94     if normalization then
  95         tmp = "densityNorm,normWith: all"
  96     else
  97         tmp = "counts"
  98     end
  99     [y, tmp, tmp, ind] = histc(data, n, tmp);
 100
 101     if length(n) == 1 then  // The number of classes is provided
 102         x = linspace(min(data), max(data), n+1); // Class boundary
 103     else // The classes are provided
 104         x = matrix(n,1,-1);   // Force row form
 105         n = length(x)-1
 106     end
 107
 108     if polygon then
 109         xmid=(x(1:$-1)+x(2:$))/2;...
 110         xp=[x(1)-(x(2)-x(1))/2 xmid x($)+(x($)-x($-1))/2];...
 111         yp=[0 y 0];
 112     end // new lines
 113
 114     // now form the polyline
 115     //    X = [x1 x1 x2 x2 x2 x3 x3 x3  x4 ...   xn xn+1 xn+1]'
 116     //    Y = [0  y1 y1 0  y2 y2 0  y3  y3 ... 0 yn yn   0 ]'
 117     X = [x(1);x(1);matrix([1;1;1]*x(2:n),-1,1);x(n+1);x(n+1)]
 118     // BUG#1885
 119     // We start the histplot line to %eps rather than 0
 120     // So when switching to logarithmic mode we do not fall
 121     // in log(0) special behaviour.
 122     Y = [matrix([%eps;1;1]*y,-1,1);%eps]
 123
 124     if opt_arg_seq == [] then
 125         plot2d(X,Y)
 126         if polygon then plot(xp,yp,"b-o"), end // new line
 127     else
 128         execstr("plot2d(X,Y"+opt_arg_seq+")")
 129         if polygon then plot(xp,yp,"r-o"), end // new line
 130     end
 131     e = gca();
 132     if polygon then
 133         e = e.children(2).children
 134     else
 135         e = e.children(1).children
 136     end
 137     e.fill_mode = "on";
 138     e.data(:,3) = -0.1;  // unmask the x-axis
 139     c = e.foreground;
 140     f = gcf();
 141     c = f.color_map(max(1,c),:)
 142     e.background = addcolor(1-(1-c)/20); // default filling color = edges one but fainter
 143 endfunction