scilab/modules/elementary_functions/macros/unique.sci

   1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   2 // Copyright (C) INRIA
   3 // Copyright (C) 2012 - 2016 - Scilab Enterprises
   4 // Copyright (C) 2018, 2019 - Samuel GOUGEON
   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 [x, ki, ko, nb] = unique(x, varargin)
  14     // extract unique components of a vector
  15     // varargin : orient=1|2|"r"|"c", "uniqueNan", "keepOrder"
  16     //
  17     // History:
  18     // * 2019 - S. Gougeon :
  19     //   - add uniqueNan option: http://bugzilla.scilab.org/15522
  20     //   - add keepOrder option: http://bugzilla.scilab.org/15795
  21     //   - add nb output option: http://bugzilla.scilab.org/8418
  22
  23     keepOrder = %f
  24     uniqueNan = %f
  25     orient = "*"
  26     newInf = [] // init Inf substitute in case of "uniqueNan" and or(x==%inf)
  27     ki = []
  28     ko = []
  29     nb = []
  30
  31     // CHECKING INPUT ARGUMENTS
  32     // ------------------------
  33     in = varargin
  34     i = 2;  // index of the current input argument
  35     if size(in)>0 then
  36         a = in(1)
  37         if typeof(a)=="string"
  38             a = convstr(a)
  39         end
  40         select a
  41         case 1
  42             orient = "r"
  43         case 2
  44             orient = "c"
  45         case "uniquenan"
  46             uniqueNan = %t
  47         case "keeporder"
  48             keepOrder = %t
  49         case "r"
  50             orient = "r"
  51         case "c"
  52             orient = "c"
  53         else
  54             msg = _("%s: Argument #%d: Must be in the set {%s}.\n")
  55             error(msprintf(msg, "unique", i, "1,2,""r"",""c"",""keepOrder"",""uniqueNan"""))
  56         end
  57         in(1) = null()
  58         i = 3
  59     end
  60     while size(in)>0 & i<5 then
  61         a = in(1)
  62         if typeof(a)=="string"
  63             a = convstr(a)
  64         end
  65         select a
  66         case "uniquenan"
  67             uniqueNan = %t
  68         case "keeporder"
  69             keepOrder = %t
  70         else
  71             msg = _("%s: Argument #%d: Must be in the set {%s}.\n")
  72             error(msprintf(msg, "unique", i, """keepOrder"",""uniqueNan"""))
  73         end
  74         in(1) = null()
  75         i = i+1
  76     end
  77     uniqueNan = uniqueNan & or(type(x)==[1 5])
  78
  79     sz = size(x);
  80     if size(x, orient)==1 then
  81         ki = 1
  82         return
  83     end
  84     if uniqueNan
  85         [x, newInf] = uniqueProcessNan(x, [], "removeNan")
  86     end
  87     getK = argn(1)>1 | keepOrder
  88
  89
  90     // [] trivial case
  91     // ---------------
  92     if isempty(x) then
  93         return  // ki, nb are already []. x is [] or sparse([])
  94     end
  95
  96     // PROCESSING complex numbers
  97     // --------------------------
  98     if or(type(x)==[1 5]) then
  99         if ~isreal(x)
 100             if isreal(x,0)
 101                 x = real(x);
 102             else
 103                 if orient=="*"
 104                     x = [real(x(:)) imag(x(:))]
 105                     if ~getK
 106                         x = unique(x,"r")
 107                     else
 108                         [x, ki, ko, nb] = unique(x,"r")
 109                     end
 110                     x = complex(x(:,1),x(:,2));
 111                     if sz(1)==1 // => put results in row
 112                         x = x.'
 113                         if getK
 114                             ki = ki'
 115                             nb = nb'
 116                         end
 117                     end
 118                 elseif orient=="r" | orient==1
 119                     x = [real(x) imag(x)]
 120                     if ~getK
 121                         x = unique(x,"r")
 122                     else
 123                         [x, ki, ko, nb] = unique(x,"r")
 124                     end
 125                     x = complex(x(:,1:sz(2)), x(:,sz(2)+1:$));
 126                 elseif orient=="c" | orient==2
 127                     x = [real(x) ; imag(x)]
 128                     if ~getK
 129                         x = unique(x,"c")
 130                     else
 131                         [x, ki, ko, nb] = unique(x,"c")
 132                     end
 133                     x = complex(x(1:sz(1),:), x(sz(1)+1:$,:));
 134                 end
 135                 if uniqueNan
 136                     x = uniqueProcessNan(x, newInf, "restoreNan")
 137                 end
 138                 if keepOrder
 139                     [ki, kk] = gsort(ki,"g","i")
 140                     select orient
 141                     case "*"
 142                         x = x(kk)
 143                     case "r"
 144                         x = x(kk,:)
 145                     case "c"
 146                         x = x(:,kk)
 147                     end
 148                     nb = nb(kk)
 149                 end
 150                 return
 151             end
 152         end
 153     end
 154
 155     // PROCESSING text and other numerical types
 156     // -----------------------------------------
 157     if orient=="*" then
 158         if getK then
 159             [x,ki] = gsort(x,"g","i");
 160             keq = x(2:$) == x(1:$-1);
 161             if argn(1)>2
 162                 nb = [0 find(~keq) size(x,"*")]
 163                 nb = nb(2:$) - nb(1:$-1)
 164             end
 165             keq = find(keq);
 166             if keq<>[] then keq = keq+1;end
 167             x(keq) = [];
 168             ki(keq) = [];
 169             if size(x,1)>1 | ndims(x)>2
 170                 x = x(:)
 171                 ki = ki(:)
 172                 nb = nb(:)
 173             end
 174         else
 175             x = gsort(x,"g","d");
 176             x = x($:-1:1);
 177             x( find(x(2:$) == x(1:$-1)) ) = [];
 178         end
 179     elseif  orient==1|orient=="r" then
 180         if getK then
 181             [x,ki] = gsort(x,"lr","i");
 182             keq = and(x(2:$,:) == x(1:$-1,:),"c")
 183             if argn(1)>2
 184                 nb = [0 find(~keq) size(x,1)]
 185                 nb = nb(2:$) - nb(1:$-1)
 186                 nb = nb(:)
 187             end
 188             keq = find(keq)
 189             if keq<>[] then keq = keq+1;end
 190             x(keq,:) = [];
 191             ki(keq,:) = [];
 192         else
 193             x = gsort(x,"lr","i");
 194             x( find(and(x(2:$,:) == x(1:$-1,:),"c")),:) = [];
 195         end
 196     elseif  orient==2|orient=="c" then
 197         if getK then
 198             [x,ki] = gsort(x,"lc","i");
 199             keq = and(x(:,2:$) == x(:,1:$-1),"r")
 200             if argn(1)>2
 201                 nb = [0 find(~keq) size(x,2)]
 202                 nb = nb(2:$) - nb(1:$-1)
 203             end
 204             keq = find(keq)
 205             if keq<>[] then keq = keq+1;end
 206             x(:,keq) = [];
 207             ki(:,keq) = [];
 208         else
 209             x = gsort(x,"lc","i");
 210             x(:, find(and(x(:,2:$) == x(:,1:$-1),"r")) ) = [];
 211         end
 212     end
 213     if uniqueNan
 214         x = uniqueProcessNan(x, newInf, "restoreNan")
 215     end
 216     if keepOrder
 217         [ki, kk] = gsort(ki,"g","i")
 218         select orient
 219         case "*"
 220             x = x(kk)
 221         case "r"
 222             x = x(kk,:)
 223         case "c"
 224             x = x(:,kk)
 225         end
 226         if argn(1)>2
 227             nb = nb(kk)
 228         end
 229     end
 230 endfunction
 231 // -------------------------------------------------------------------
 232
 233 // To consider Nan mutually equal, we replace all of them with a "regular" substitute.
 234 // Since Nan are sorted as > Inf, we must use anyway Inf as the Nan substitute.
 235 // If the original array have already some Inf, we must priorly replace them with
 236 // a decimal greater than the finite maximum of the array values.
 237 // After processing, we restore Inf => Nan (, and maxNum => Inf).
 238
 239 function [x, newInf] = uniqueProcessNan(x, newInf, way)
 240
 241     if way=="removeNan" & or(isnan(x)) then
 242         // Replacing Nan
 243         // -------------
 244         if isreal(x)
 245             if or(x==%inf)
 246                 b = x==%inf
 247                 m = max([1 max(x(~b))])
 248                 newInf = m*1.5
 249                 x(b) = newInf
 250             end
 251             x(x<>x) = %inf
 252         else
 253             r = real(x)
 254             i = imag(x)
 255             if or([r i]==%inf)
 256                 br = r==%inf
 257                 m = max(r(~br),1)
 258                 bi = i==%inf
 259                 m = max(i(~bi),m,1)
 260                 newInf = m*1.5
 261                 r(br) = newInf
 262                 i(bi) = newInf
 263             end
 264             r(r<>r) = %inf
 265             i(i<>i) = %inf
 266             x = complex(r,i);
 267         end
 268
 269     // Restoring  NaN
 270     // --------------
 271     elseif way=="restoreNan"
 272         if isreal(x)
 273             x(x==%inf) = %nan
 274             if newInf~=[]
 275                 x(x==newInf) = %inf
 276             end
 277         else
 278             r = real(x)
 279             r(r==%inf) = %nan
 280             i = imag(x)
 281             i(i==%inf) = %nan
 282             if newInf~=[]
 283                 r(r==newInf) = %inf
 284                 i(i==newInf) = %inf
 285             end
 286             x = complex(r, i)
 287         end
 288     end
 289 endfunction