scilab/modules/elementary_functions/help/en_US/searchandsort/dsearch.xml

   1 <?xml version="1.0" encoding="utf-8"?>
   2 <!--
   3  * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
   4  * Copyright (C) Bruno Pincon
   5  * Copyright (C) 2011 - DIGITEO - Michael Baudin
   6  *
   7  * This file must be used under the terms of the CeCILL.
   8  * This source file is licensed as described in the file COPYING, which
   9  * you should have received as part of this distribution.  The terms
  10  * are also available at
  11  * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
  12  *
  13  -->
  14 <refentry xmlns="http://docbook.org/ns/docbook" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:svg="http://www.w3.org/2000/svg" xmlns:ns5="http://www.w3.org/1999/xhtml" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:db="http://docbook.org/ns/docbook"   xmlns:scilab="http://www.scilab.org"  xml:id="dsearch" xml:lang="en">
  15     <refnamediv>
  16         <refname>dsearch</refname>
  17         <refpurpose>
  18             search in ordered sets
  19         </refpurpose>
  20     </refnamediv>
  21     <refsynopsisdiv>
  22         <title>Calling Sequence</title>
  23         <synopsis>
  24             [ind, occ, info] = dsearch(X, s )
  25             [ind, occ, info] = dsearch(X, s , ch )
  26         </synopsis>
  27     </refsynopsisdiv>
  28     <refsection>
  29         <title>Arguments</title>
  30         <variablelist>
  31             <varlistentry>
  32                 <term>X</term>
  33                 <listitem>
  34                     <para>
  35                         a mx-by-nx matrix of doubles, the entries to locate.
  36                     </para>
  37                 </listitem>
  38             </varlistentry>
  39             <varlistentry>
  40                 <term>s</term>
  41                 <listitem>
  42                     <para>
  43                         a n-by-1 or 1-by-n matrix of doubles, the
  44                         intervals (if ch="c") or the set (if ch="d").
  45                         The values in s must be in strictly increasing order:
  46                         s(1) &lt; s(2) &lt; ... &lt; s(n)
  47                     </para>
  48                 </listitem>
  49             </varlistentry>
  50             <varlistentry>
  51                 <term>ch</term>
  52                 <listitem>
  53                     <para>
  54                         a 1-by-1 matrix of strings, the type of search (default ch="c").
  55                         Available values are ch="c" or ch="d".
  56                     </para>
  57                 </listitem>
  58             </varlistentry>
  59             <varlistentry>
  60                 <term>ind</term>
  61                 <listitem>
  62                     <para>
  63                         a mx-by-nx matrix of doubles.
  64                         The location of the X values in the intervals or the set defined by <literal>s</literal>.
  65                     </para>
  66                 </listitem>
  67             </varlistentry>
  68             <varlistentry>
  69                 <term>occ</term>
  70                 <listitem>
  71                     <para>
  72                         If ch="c", a (n-1)-by-1 or 1-by-(n-1) of doubles.
  73                         If ch="d", a n-by-1 or 1-by-n of doubles. The number of entries
  74                         from X in the intervals s.
  75                     </para>
  76                 </listitem>
  77             </varlistentry>
  78             <varlistentry>
  79                 <term>info</term>
  80                 <listitem>
  81                     <para>
  82                         a 1-by-1 matrix of doubles. The number of X entries which are not in <literal>s</literal>.
  83                     </para>
  84                 </listitem>
  85             </varlistentry>
  86         </variablelist>
  87     </refsection>
  88     <refsection>
  89         <title>Description</title>
  90         <para>
  91             This function locates the indices of the entries in X in the intervals or the set defined by s.
  92         </para>
  93         <para>
  94             If <literal>ch="c"</literal> (default), we consider the
  95             intervals <literal>I1 = [s(1), s(2)]</literal>
  96             and <literal>Ik = (s(k), s(k+1)]</literal> for <literal>k=2,...,n-1</literal>.
  97             Notice that the bounds of <literal>I1</literal> are closed,
  98             while the left bounds of <literal>I2</literal>, ..., <literal>In</literal> are
  99             opened.
 100             For each entry <literal>X(i)</literal>,
 101             we search the interval <literal>Ik</literal> which contains X(i), i.e. we search k such
 102             that <literal>s(k)&lt;X(i)&lt;=s(k+1)</literal>.
 103         </para>
 104         <para>
 105             More precisely,
 106         </para>
 107         <variablelist>
 108             <varlistentry>
 109                 <term>ind(i)</term>
 110                 <listitem>
 111                     <para>
 112                         is k such that <literal>Ik</literal> contains X(i) or
 113                         0 if X(i) is not in any interval <literal>Ik</literal>.
 114                     </para>
 115                 </listitem>
 116             </varlistentry>
 117             <varlistentry>
 118                 <term>occ(k)</term>
 119                 <listitem>
 120                     <para>
 121                         is the number of components of <literal>X</literal> which are in <literal>Ik</literal>
 122                     </para>
 123                 </listitem>
 124             </varlistentry>
 125             <varlistentry>
 126                 <term>info</term>
 127                 <listitem>
 128                     <para>
 129                         is the number of components of <literal>X</literal> which are not in
 130                         any interval <literal>Ik</literal>.
 131                     </para>
 132                 </listitem>
 133             </varlistentry>
 134         </variablelist>
 135         <para>
 136             If <literal>ch="d"</literal>, we consider the set {s(1),s(2),...,s(n)}.
 137             For each X(i), search k such that X(i)=s(k).
 138         </para>
 139         <para>
 140             More precisely,
 141         </para>
 142         <variablelist>
 143             <varlistentry>
 144                 <term>ind(i)</term>
 145                 <listitem>
 146                     <para>
 147                         is k if X(i)=s(k) or 0 if X(i) is not in s.
 148                     </para>
 149                 </listitem>
 150             </varlistentry>
 151             <varlistentry>
 152                 <term>occ(k)</term>
 153                 <listitem>
 154                     <para>
 155                         is the number of entries in X equal to s(k)
 156                     </para>
 157                 </listitem>
 158             </varlistentry>
 159             <varlistentry>
 160                 <term>info</term>
 161                 <listitem>
 162                     <para>
 163                         is the number of entries in X which are not in the set
 164                         {s(1),...,s(n)}.
 165                     </para>
 166                 </listitem>
 167             </varlistentry>
 168         </variablelist>
 169         <para>
 170             The ch="c" option can be used to compute the empirical histogram of a function given a dataset.
 171             The ch="d" option can be used to compute the entries in X which are present in the set s.
 172         </para>
 173     </refsection>
 174     <refsection>
 175         <title>Examples</title>
 176         <para>
 177             In the following example, we consider 3 intervals <literal>I1=[5,11]</literal>,
 178             <literal>I2=(11,15]</literal> and <literal>I3=(15,20]</literal>.
 179             We are looking for the location of the entries of <literal>X=[11 13 1 7 5 2 9]</literal>
 180             in these intervals.
 181         </para>
 182         <programlisting role="example"><![CDATA[
 183 [ind, occ, info] = dsearch([11 13 1 7 5 2 9], [5 11 15 20])
 184     ]]></programlisting>
 185         <para>
 186             The previous example produces the following output.
 187         </para>
 188         <screen><![CDATA[
 189 -->[ind, occ, info] = dsearch([11 13 1 7 5 2 9], [5 11 15 20])
 190  info  =
 191     2.
 192  occ  =
 193     4.    1.    0.
 194  ind  =
 195     1.    2.    0.    1.    1.    0.    1.
 196     ]]></screen>
 197         <para>
 198             We now explain these results.
 199         </para>
 200         <itemizedlist>
 201             <listitem>
 202                 <para>X(1)=11 is in the interval I1, hence ind(1)=1.</para>
 203             </listitem>
 204             <listitem>
 205                 <para> X(2)=13 is in the interval I2, hence ind(2)=2.</para>
 206             </listitem>
 207             <listitem>
 208                 <para> X(3)=1 is not in any interval, hence ind(3)=0.</para>
 209             </listitem>
 210             <listitem>
 211                 <para> X(4)=7 is in the interval I1, hence ind(4)=1.</para>
 212             </listitem>
 213             <listitem>
 214                 <para> X(5)=5 is in the interval I1, hence ind(5)=1.</para>
 215             </listitem>
 216             <listitem>
 217                 <para> X(6)=2 is not in any interval, hence ind(6)=0.</para>
 218             </listitem>
 219             <listitem>
 220                 <para> X(7)=9 is in the interval I1, hence ind(7)=1.</para>
 221             </listitem>
 222             <listitem>
 223                 <para> There are four X entries (5, 7, 9 and 11) in I1, hence occ(1)=4.</para>
 224             </listitem>
 225             <listitem>
 226                 <para> There is one X entry (i.e. 13) in I2, hence occ(2)=1.</para>
 227             </listitem>
 228             <listitem>
 229                 <para> There is no X entry in I3, hence occ(3)=0.</para>
 230             </listitem>
 231             <listitem>
 232                 <para> There are two X entries (i.e. 1, 2) which are not in any interval, hence info=2.</para>
 233             </listitem>
 234         </itemizedlist>
 235         <para>
 236             In the following example, we consider the set [5 11 15 20] and
 237             are searching the location of the X entries in this set.
 238         </para>
 239         <programlisting role="example"><![CDATA[
 240 [ind, occ, info] = dsearch([11 13 1 7 5 2 9], [5 11 15 20],"d" )
 241     ]]></programlisting>
 242         <para>
 243             The previous example produces the following output.
 244         </para>
 245         <screen><![CDATA[
 246 -->[ind, occ, info] = dsearch([11 13 1 7 5 2 9], [5 11 15 20],"d" )
 247  info  =
 248     5.
 249  occ  =
 250     1.    1.    0.    0.
 251  ind  =
 252     2.    0.    0.    0.    1.    0.    0.
 253     ]]></screen>
 254         <para>
 255             The following is a detailed explanation for the previous results.
 256         </para>
 257         <itemizedlist>
 258             <listitem>
 259                 <para>
 260                     X(1)=11 is in the set <literal>s</literal> at position #2, hence ind(1)=2.
 261                 </para>
 262             </listitem>
 263             <listitem>
 264                 <para>
 265                     X(2)=13 is not in the set <literal>s</literal>, hence ind(2)=0.
 266                 </para>
 267             </listitem>
 268             <listitem>
 269                 <para>
 270                     X(3)=1 is not in the set <literal>s</literal>, hence ind(3)=0.
 271                 </para>
 272             </listitem>
 273             <listitem>
 274                 <para>
 275                     X(4)=7 is not in the set <literal>s</literal>, hence ind(4)=0.
 276                 </para>
 277             </listitem>
 278             <listitem>
 279                 <para>
 280                     X(5)=5 is in the set <literal>s</literal> at position #1, hence ind(5)=1.
 281                 </para>
 282             </listitem>
 283             <listitem>
 284                 <para>
 285                     X(6)=2 is not in the set <literal>s</literal>, hence ind(6)=0.
 286                 </para>
 287             </listitem>
 288             <listitem>
 289                 <para>
 290                     X(7)=9 is not in the set <literal>s</literal>, hence ind(7)=0.
 291                 </para>
 292             </listitem>
 293             <listitem>
 294                 <para>
 295                     There is one entry X (i.e. 5) equal to 5, hence occ(1)=1.
 296                 </para>
 297             </listitem>
 298             <listitem>
 299                 <para>
 300                     There is one entry X (i.e. 11) equal to 1, hence occ(2)=1.
 301                 </para>
 302             </listitem>
 303             <listitem>
 304                 <para>
 305                     There are no entries matching <literal>s(3)</literal>, hence occ(3)=0.
 306                 </para>
 307             </listitem>
 308             <listitem>
 309                 <para>
 310                     There are no entries matching <literal>s(4)</literal>, hence occ(4)=0.
 311                 </para>
 312             </listitem>
 313             <listitem>
 314                 <para>
 315                     There are five X entries (i.e. 13, 1, 7, 2, 9) which are not in the set <literal>s</literal>, hence info=5.
 316                 </para>
 317             </listitem>
 318         </itemizedlist>
 319         <para>
 320             The values in <literal>s</literal> must be in increasing order, whatever the
 321             value of the <literal>ch</literal> option.
 322             If this is not true, an error is generated, as in the following session.
 323         </para>
 324         <screen><![CDATA[
 325 -->dsearch([11 13 1 7 5 2 9], [2 1])
 326 !--error 999
 327 dsearch   : the array s (arg 2) is not well ordered
 328 -->dsearch([11 13 1 7 5 2 9], [2 1],"d")
 329 !--error 999
 330 dsearch   : the array s (arg 2) is not well ordered
 331     ]]></screen>
 332     </refsection>
 333     <refsection>
 334         <title>Advanced Examples</title>
 335         <para>
 336             In the following example, we compare the empirical histogram of uniform
 337             random numbers in [0,1) with the uniform distribution function.
 338             To perform this comparison, we use the default search algorithm based on
 339             intervals (ch="c").
 340             We generate X as a collection of m=50 000 uniform random numbers in
 341             the range [0,1).
 342             We consider the n=10 values equally equally spaced values in the [0,1] range and
 343             consider the associated intervals.
 344             Then we count the number of entries in X which fall in the intervals:
 345             this is the empirical histogram of the uniform distribution function.
 346             The expectation for occ/m is equal to 1/(n-1).
 347         </para>
 348         <programlisting role="example"><![CDATA[
 349 m = 50000 ;
 350 n = 10;
 351 X = grand(m,1,"def");
 352 s = linspace(0,1,n)';
 353 [ind, occ] = dsearch(X, s);
 354 e = 1/(n-1)*ones(1,n-1);
 355 scf() ;
 356 plot(s(1:n-1), occ/m,"bo");
 357 plot(s(1:n-1), e,"r-");
 358 legend(["Experiment","Expectation"]);
 359 xtitle("Uniform random numbers","X","P(X)");
 360  ]]></programlisting>
 361         <scilab:image>
 362             m = 50000 ;
 363             n = 10;
 364             X = grand(m,1,"def");
 365             s = linspace(0,1,n)';
 366             [ind, occ] = dsearch(X, s);
 367             e = 1/(n-1)*ones(1,n-1);
 368             scf() ;
 369             plot(s(1:n-1), occ/m,"bo");
 370             plot(s(1:n-1), e,"r-");
 371             legend(["Experiment","Expectation"]);
 372             xtitle("Uniform random numbers","X","P(X)");
 373         </scilab:image>
 374         <para>
 375             In the following example, we compare the histogram of binomially distributed
 376             random numbers with the binomial probability distribution function B(N,p),
 377             with N=8 and p=0.5.
 378             To perform this comparison, we use the discrete search algorithm based on
 379             a set (ch="d").
 380         </para>
 381         <programlisting role="example"><![CDATA[
 382 N = 8 ;
 383 p = 0.5;
 384 m = 50000;
 385 X = grand(m,1,"bin",N,p);
 386 s = (0:N)';
 387 [ind, occ] = dsearch(X, s, "d");
 388 Pexp = occ/m;
 389 Pexa = binomial(p,N);
 390 scf() ;
 391 plot(s,Pexp,"bo");
 392 plot(s,Pexa',"r-");
 393 xtitle("Binomial distribution B(8,0.5)","X","P(X)");
 394 legend(["Experiment","Expectation"]);
 395 ]]></programlisting>
 396         <scilab:image>
 397             N = 8 ;
 398             p = 0.5;
 399             m = 50000;
 400             X = grand(m,1,"bin",N,p);
 401             s = (0:N)';
 402             [ind, occ] = dsearch(X, s, "d");
 403             Pexp = occ/m;
 404             Pexa = binomial(p,N);
 405             scf() ;
 406             plot(s,Pexp,"bo");
 407             plot(s,Pexa',"r-");
 408             xtitle("Binomial distribution B(8,0.5)","X","P(X)");
 409             legend(["Experiment","Expectation"]);
 410         </scilab:image>
 411         <para>
 412             In the following example, we use piecewise Hermite polynomials to
 413             interpolate a dataset.
 414         </para>
 415         <programlisting role="example"><![CDATA[
 416
 417       // define Hermite base functions
 418       function y=Ll(t,k,x)
 419         // Lagrange left on Ik
 420         y=(t-x(k+1))./(x(k)-x(k+1))
 421       endfunction
 422       function y=Lr(t,k,x)
 423         // Lagrange right on Ik
 424         y=(t-x(k))./(x(k+1)-x(k))
 425       endfunction
 426       function y=Hl(t,k,x)
 427         y=(1-2*(t-x(k))./(x(k)-x(k+1))).*Ll(t,k,x).^2
 428       endfunction
 429       function y=Hr(t,k,x)
 430         y=(1-2*(t-x(k+1))./(x(k+1)-x(k))).*Lr(t,k,x).^2
 431       endfunction
 432       function y=Kl(t,k,x)
 433         y=(t-x(k)).*Ll(t,k,x).^2
 434       endfunction
 435       function y=Kr(t,k,x)
 436         y=(t-x(k+1)).*Lr(t,k,x).^2
 437       endfunction
 438
 439       x = [0 ; 0.2 ; 0.35 ; 0.5 ; 0.65 ; 0.8 ;  1];
 440       y = [0 ; 0.1 ;-0.1  ; 0   ; 0.4  ;-0.1 ;  0];
 441       d = [1 ; 0   ; 0    ; 1   ; 0    ; 0   ; -1];
 442       X = linspace(0, 1, 200)';
 443       ind = dsearch(X, x);
 444
 445       // plot the curve
 446       Y = y(ind).*Hl(X,ind) + y(ind+1).*Hr(X,ind) + d(ind).*Kl(X,ind) + d(ind+1).*Kr(X,ind);
 447       scf();
 448       plot(X,Y,"k-");
 449       plot(x,y,"bo")
 450       xtitle("Hermite piecewise polynomial");
 451       legend(["Polynomial","Data"]);
 452       // NOTE : it can be verified by adding these ones :
 453       // YY = interp(X,x,y,d); plot2d(X,YY,3,"000")
 454       ]]></programlisting>
 455         <scilab:image>
 456
 457             function y=Ll(t,k,x)
 458             y=(t-x(k+1))./(x(k)-x(k+1))
 459             endfunction
 460             function y=Lr(t,k,x)
 461             y=(t-x(k))./(x(k+1)-x(k))
 462             endfunction
 463             function y=Hl(t,k,x)
 464             y=(1-2*(t-x(k))./(x(k)-x(k+1))).*Ll(t,k,x).^2
 465             endfunction
 466             function y=Hr(t,k,x)
 467             y=(1-2*(t-x(k+1))./(x(k+1)-x(k))).*Lr(t,k,x).^2
 468             endfunction
 469             function y=Kl(t,k,x)
 470             y=(t-x(k)).*Ll(t,k,x).^2
 471             endfunction
 472             function y=Kr(t,k,x)
 473             y=(t-x(k+1)).*Lr(t,k,x).^2
 474             endfunction
 475
 476             x = [0 ; 0.2 ; 0.35 ; 0.5 ; 0.65 ; 0.8 ;  1];
 477             y = [0 ; 0.1 ;-0.1  ; 0   ; 0.4  ;-0.1 ;  0];
 478             d = [1 ; 0   ; 0    ; 1   ; 0    ; 0   ; -1];
 479             X = linspace(0, 1, 200)';
 480             ind = dsearch(X, x);
 481
 482             Y = y(ind).*Hl(X,ind) + y(ind+1).*Hr(X,ind) + d(ind).*Kl(X,ind) + d(ind+1).*Kr(X,ind);
 483             scf();
 484             plot(X,Y,"k-");
 485             plot(x,y,"bo")
 486             xtitle("Hermite piecewise polynomial");
 487             legend(["Polynomial","Data"]);
 488         </scilab:image>
 489     </refsection>
 490     <refsection role="see also">
 491         <title>See Also</title>
 492         <simplelist type="inline">
 493             <member>
 494                 <link linkend="find">find</link>
 495             </member>
 496             <member>
 497                 <link linkend="tabul">tabul</link>
 498             </member>
 499         </simplelist>
 500     </refsection>
 501 </refentry>