654b850f0e17ef16d41f317d16842946f9ad1f01
[scilab.git] / scilab / modules / elementary_functions / macros / linspace.sci
1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2 // Copyright (C) INRIA
3 // Copyright (C) DIGITEO - 2011 - Allan CORNET
4 // Copyright (C) CNES - 2011 - Guillaume AZEMA
5 // Copyright (C) 2012 - 2016 - Scilab Enterprises
6 // Copyright (C) 2016, 2018 - Samuel GOUGEON
7 //
8 // This file is hereby licensed under the terms of the GNU GPL v2.0,
9 // pursuant to article 5.3.4 of the CeCILL v.2.1.
10 // This file was originally licensed under the terms of the CeCILL v2.1,
11 // and continues to be available under such terms.
12 // For more information, see the COPYING file which you should have received
13 // along with this program.
14
15 function y = linspace(d1, d2, n)
16     // Linearly spaced vector.
17     // linspace(x1, x2) generates a row vector of 100 linearly
18     // equally spaced points between x1 and x2.
19     // linspace(x1, x2, n) generates n points between x1 and x2.
20
21     // CHECKING ARGUMENTS
22     // ------------------
23     rhs = argn(2);
24     if rhs < 2 then
25         msg = gettext("%s: Wrong number of input argument(s): %d to %d expected.\n")
26         error(msprintf(msg, "linspace", 2, 3));
27     end
28     if size(d1,2)<>1 then
29         msg = gettext("%s: Argument #%d: Column expected.\n")
30         error(msprintf(msg, "linspace", 1));
31     end
32     if ~and(size(d1) == size(d2)) then
33         msg = gettext("%s: Arguments #%d and #%d: Same sizes expected.\n")
34         error(msprintf(msg, "linspace", 1, 2));
35     end
36     msg = gettext("%s: Argument #%d: Number(s) expected.\n")
37     if ~or(type(d1)==[1 5 8]) then
38         error(msprintf(msg, "linspace", 1));
39     end
40     if ~or(type(d2)==[1 5 8]) then
41         error(msprintf(msg, "linspace", 2));
42     end
43     if type(d1)==8 & type(d2)==8 & inttype(d1)~=inttype(d2) then
44         msg = gettext("%s: Arguments #%d and #%d: Same integer types expected.\n")
45         error(msprintf(msg, "linspace", 2));
46     end
47     msg = gettext("%s: Argument #%d: %%nan and %%inf values are forbidden.\n")
48     if type(d1)~=8 & (or(isinf(d1)) | or(isnan(d1))) then
49         error(msprintf(msg, "linspace", 1));
50     end
51     if type(d2)~=8 & (or(isinf(d2)) | or(isnan(d2))) then
52         error(msprintf(msg, "linspace", 2));
53     end
54
55     if rhs == 2 then
56         n = 100;
57     else
58         if and(type(n)<>[1 8]) | size(n,"*")<>1 | int(n)<>n then
59             msg = gettext("%s: Argument #%d: An integer value expected.\n")
60             error(msprintf(msg, "linspace",3));
61         end
62         n = double(n); // Convert for the operations to come
63     end
64
65     // PROCESSING
66     // ----------
67     if n==1
68         y = d2
69     elseif n==2
70         y = [d1 d2]
71     elseif n>2
72         if or(inttype(d1)==[8 18])
73             y = linspace_integers_64(d1,d2,n)
74         else
75             if type(d1)==8
76                 span = double(d2) - double(d1)
77             else
78                 span = d2 - d1
79             end
80             y = (span * (0:n-1)) / (n-1) + d1 * ones(1,n);
81         end
82         // Forces the last value to be exactly the given d2:
83         // http://bugzilla.scilab.org/10966
84         y(:,$) = d2;
85     else
86         y = []
87     end
88 endfunction
89 // -----------------
90 function y = linspace_integers_64(d1,d2,n)
91         s = d2>=d1
92         span = zeros(d1)
93         if typeof(d1)=="int64"  then
94             span = uint64(span)
95         end
96         span(s) = d2(s) - d1(s)
97         span(~s) = d1(~s) - d2(~s)
98         step = span/(n-1)
99         y = iconvert(zeros(size(d1,1),n), inttype(d1))
100         if or(s)  // d2 > d1
101             y(s,:) = d1(s)*ones(1,n) + step(s)*(0:n-1)
102         end
103         if or(~s)
104             y(~s,:) = d1(~s)*ones(1,n) - step(~s)*(0:n-1)
105         end
106         y(:,$) = d2
107         // Computing the balancing corrections
108             // We computes actual intervals
109         if typeof(d1)=="uint64" then
110             i = int64(zeros(size(d1,1),n-1))
111             i(s, :) = y(s, 2:$)    -  y(s, 1:$-1)
112             i(~s,:) = y(~s, 1:$-1) -  y(~s, 2:$)
113         else // int64
114             i = diff(y,1,2);
115         end
116             // diff wrt the 1st interval's width
117         e = i - i(:,1) * ones(1,n-1)
118
119             // Computing corrections
120         S = mean(double(e), 'c') * ones(1,n-1)
121         cumE = cumsum(e, 'c')
122         cumS = cumsum(S, 'c')
123         delta = cumS - cumE
124         // Applying the corrections
125         y(:,2:$) = y(:,2:$) + delta
126 endfunction