1Primitive(3)          User Contributed Perl Documentation         Primitive(3)
2
3
4

NAME

6       PDL::Primitive - primitive operations for pdl
7

DESCRIPTION

9       This module provides some primitive and useful functions defined using
10       PDL::PP and able to use the new indexing tricks.
11
12       See PDL::Indexing for how to use indices creatively.  For explanation
13       of the signature format, see PDL::PP.
14

SYNOPSIS

16        # Pulls in PDL::Primitive, among other modules.
17        use PDL;
18
19        # Only pull in PDL::Primitive:
20        use PDL::Primitive;
21

FUNCTIONS

23   inner
24         Signature: (a(n); b(n); [o]c())
25
26       Inner product over one dimension
27
28        c = sum_i a_i * b_i
29
30       If "a() * b()" contains only bad data, c() is set bad. Otherwise c()
31       will have its bad flag cleared, as it will not contain any bad values.
32
33   outer
34         Signature: (a(n); b(m); [o]c(n,m))
35
36       outer product over one dimension
37
38       Naturally, it is possible to achieve the effects of outer product
39       simply by broadcasting over the ""*"" operator but this function is
40       provided for convenience.
41
42       outer processes bad values.  It will set the bad-value flag of all
43       output ndarrays if the flag is set for any of the input ndarrays.
44
45   x
46        Signature: (a(i,z), b(x,i),[o]c(x,z))
47
48       Matrix multiplication
49
50       PDL overloads the "x" operator (normally the repeat operator) for
51       matrix multiplication.  The number of columns (size of the 0 dimension)
52       in the left-hand argument must normally equal the number of rows (size
53       of the 1 dimension) in the right-hand argument.
54
55       Row vectors are represented as (N x 1) two-dimensional PDLs, or you may
56       be sloppy and use a one-dimensional PDL.  Column vectors are
57       represented as (1 x N) two-dimensional PDLs.
58
59       Broadcasting occurs in the usual way, but as both the 0 and 1 dimension
60       (if present) are included in the operation, you must be sure that you
61       don't try to broadcast over either of those dims.
62
63       Of note, due to how Perl v5.14.0 and above implement operator
64       overloading of the "x" operator, the use of parentheses for the left
65       operand creates a list context, that is
66
67        pdl> ( $x * $y ) x $z
68        ERROR: Argument "..." isn't numeric in repeat (x) ...
69
70       treats $z as a numeric count for the list repeat operation and does not
71       call the scalar form of the overloaded operator. To use the operator in
72       this case, use a scalar context:
73
74        pdl> scalar( $x * $y ) x $z
75
76       or by calling "matmult" directly:
77
78        pdl> ( $x * $y )->matmult( $z )
79
80       EXAMPLES
81
82       Here are some simple ways to define vectors and matrices:
83
84        pdl> $r = pdl(1,2);                # A row vector
85        pdl> $c = pdl([[3],[4]]);          # A column vector
86        pdl> $c = pdl(3,4)->(*1);          # A column vector, using NiceSlice
87        pdl> $m = pdl([[1,2],[3,4]]);      # A 2x2 matrix
88
89       Now that we have a few objects prepared, here is how to matrix-multiply
90       them:
91
92        pdl> print $r x $m                 # row x matrix = row
93        [
94         [ 7 10]
95        ]
96
97        pdl> print $m x $r                 # matrix x row = ERROR
98        PDL: Dim mismatch in matmult of [2x2] x [2x1]: 2 != 1
99
100        pdl> print $m x $c                 # matrix x column = column
101        [
102         [ 5]
103         [11]
104        ]
105
106        pdl> print $m x 2                  # Trivial case: scalar mult.
107        [
108         [2 4]
109         [6 8]
110        ]
111
112        pdl> print $r x $c                 # row x column = scalar
113        [
114         [11]
115        ]
116
117        pdl> print $c x $r                 # column x row = matrix
118        [
119         [3 6]
120         [4 8]
121        ]
122
123       INTERNALS
124
125       The mechanics of the multiplication are carried out by the "matmult"
126       method.
127
128   matmult
129         Signature: (a(t,h); b(w,t); [o]c(w,h))
130
131       Matrix multiplication
132
133       Notionally, matrix multiplication $x x $y is equivalent to the
134       broadcasting expression
135
136           $x->dummy(1)->inner($y->xchg(0,1)->dummy(2),$c);
137
138       but for large matrices that breaks CPU cache and is slow.  Instead,
139       matmult calculates its result in 32x32x32 tiles, to keep the memory
140       footprint within cache as long as possible on most modern CPUs.
141
142       For usage, see "x", a description of the overloaded 'x' operator
143
144       matmult ignores the bad-value flag of the input ndarrays.  It will set
145       the bad-value flag of all output ndarrays if the flag is set for any of
146       the input ndarrays.
147
148   innerwt
149         Signature: (a(n); b(n); c(n); [o]d())
150
151       Weighted (i.e. triple) inner product
152
153        d = sum_i a(i) b(i) c(i)
154
155       innerwt processes bad values.  It will set the bad-value flag of all
156       output ndarrays if the flag is set for any of the input ndarrays.
157
158   inner2
159         Signature: (a(n); b(n,m); c(m); [o]d())
160
161       Inner product of two vectors and a matrix
162
163        d = sum_ij a(i) b(i,j) c(j)
164
165       Note that you should probably not broadcast over "a" and "c" since that
166       would be very wasteful. Instead, you should use a temporary for "b*c".
167
168       inner2 processes bad values.  It will set the bad-value flag of all
169       output ndarrays if the flag is set for any of the input ndarrays.
170
171   inner2d
172         Signature: (a(n,m); b(n,m); [o]c())
173
174       Inner product over 2 dimensions.
175
176       Equivalent to
177
178        $c = inner($x->clump(2), $y->clump(2))
179
180       inner2d processes bad values.  It will set the bad-value flag of all
181       output ndarrays if the flag is set for any of the input ndarrays.
182
183   inner2t
184         Signature: (a(j,n); b(n,m); c(m,k); [t]tmp(n,k); [o]d(j,k)))
185
186       Efficient Triple matrix product "a*b*c"
187
188       Efficiency comes from by using the temporary "tmp". This operation only
189       scales as "N**3" whereas broadcasting using "inner2" would scale as
190       "N**4".
191
192       The reason for having this routine is that you do not need to have the
193       same broadcast-dimensions for "tmp" as for the other arguments, which
194       in case of large numbers of matrices makes this much more memory-
195       efficient.
196
197       It is hoped that things like this could be taken care of as a kind of
198       closures at some point.
199
200       inner2t processes bad values.  It will set the bad-value flag of all
201       output ndarrays if the flag is set for any of the input ndarrays.
202
203   crossp
204         Signature: (a(tri=3); b(tri); [o] c(tri))
205
206       Cross product of two 3D vectors
207
208       After
209
210        $c = crossp $x, $y
211
212       the inner product "$c*$x" and "$c*$y" will be zero, i.e. $c is
213       orthogonal to $x and $y
214
215       crossp does not process bad values.  It will set the bad-value flag of
216       all output ndarrays if the flag is set for any of the input ndarrays.
217
218   norm
219         Signature: (vec(n); [o] norm(n))
220
221       Normalises a vector to unit Euclidean length
222
223       norm processes bad values.  It will set the bad-value flag of all
224       output ndarrays if the flag is set for any of the input ndarrays.
225
226   indadd
227         Signature: (input(n); indx ind(n); [io] sum(m))
228
229       Broadcasting index add: add "input" to the "ind" element of "sum", i.e:
230
231        sum(ind) += input
232
233       Simple example:
234
235         $x = 2;
236         $ind = 3;
237         $sum = zeroes(10);
238         indadd($x,$ind, $sum);
239         print $sum
240         #Result: ( 2 added to element 3 of $sum)
241         # [0 0 0 2 0 0 0 0 0 0]
242
243       Broadcasting example:
244
245         $x = pdl( 1,2,3);
246         $ind = pdl( 1,4,6);
247         $sum = zeroes(10);
248         indadd($x,$ind, $sum);
249         print $sum."\n";
250         #Result: ( 1, 2, and 3 added to elements 1,4,6 $sum)
251         # [0 1 0 0 2 0 3 0 0 0]
252
253       The routine barfs on bad indices, and bad inputs set target outputs
254       bad.
255
256   conv1d
257         Signature: (a(m); kern(p); [o]b(m); int reflect)
258
259       1D convolution along first dimension
260
261       The m-th element of the discrete convolution of an input ndarray $a of
262       size $M, and a kernel ndarray $kern of size $P, is calculated as
263
264                                     n = ($P-1)/2
265                                     ====
266                                     \
267         ($a conv1d $kern)[m]   =     >      $a_ext[m - n] * $kern[n]
268                                     /
269                                     ====
270                                     n = -($P-1)/2
271
272       where $a_ext is either the periodic (or reflected) extension of $a so
273       it is equal to $a on " 0..$M-1 " and equal to the corresponding
274       periodic/reflected image of $a outside that range.
275
276         $con = conv1d sequence(10), pdl(-1,0,1);
277
278         $con = conv1d sequence(10), pdl(-1,0,1), {Boundary => 'reflect'};
279
280       By default, periodic boundary conditions are assumed (i.e. wrap
281       around).  Alternatively, you can request reflective boundary conditions
282       using the "Boundary" option:
283
284         {Boundary => 'reflect'} # case in 'reflect' doesn't matter
285
286       The convolution is performed along the first dimension. To apply it
287       across another dimension use the slicing routines, e.g.
288
289         $y = $x->mv(2,0)->conv1d($kernel)->mv(0,2); # along third dim
290
291       This function is useful for broadcasted filtering of 1D signals.
292
293       Compare also conv2d, convolve, fftconvolve, fftwconv, rfftwconv
294
295       WARNING: "conv1d" processes bad values in its inputs as the numeric
296       value of "$pdl->badvalue" so it is not recommended for processing pdls
297       with bad values in them unless special care is taken.
298
299       conv1d ignores the bad-value flag of the input ndarrays.  It will set
300       the bad-value flag of all output ndarrays if the flag is set for any of
301       the input ndarrays.
302
303   in
304         Signature: (a(); b(n); [o] c())
305
306       test if a is in the set of values b
307
308          $goodmsk = $labels->in($goodlabels);
309          print pdl(3,1,4,6,2)->in(pdl(2,3,3));
310         [1 0 0 0 1]
311
312       "in" is akin to the is an element of of set theory. In principle, PDL
313       broadcasting could be used to achieve its functionality by using a
314       construct like
315
316          $msk = ($labels->dummy(0) == $goodlabels)->orover;
317
318       However, "in" doesn't create a (potentially large) intermediate and is
319       generally faster.
320
321       in does not process bad values.  It will set the bad-value flag of all
322       output ndarrays if the flag is set for any of the input ndarrays.
323
324   uniq
325       return all unique elements of an ndarray
326
327       The unique elements are returned in ascending order.
328
329         PDL> p pdl(2,2,2,4,0,-1,6,6)->uniq
330         [-1 0 2 4 6]     # 0 is returned 2nd (sorted order)
331
332         PDL> p pdl(2,2,2,4,nan,-1,6,6)->uniq
333         [-1 2 4 6 nan]   # NaN value is returned at end
334
335       Note: The returned pdl is 1D; any structure of the input ndarray is
336       lost.  "NaN" values are never compare equal to any other values, even
337       themselves.  As a result, they are always unique. "uniq" returns the
338       NaN values at the end of the result ndarray.  This follows the Matlab
339       usage.
340
341       See "uniqind" if you need the indices of the unique elements rather
342       than the values.
343
344       Bad values are not considered unique by uniq and are ignored.
345
346        $x=sequence(10);
347        $x=$x->setbadif($x%3);
348        print $x->uniq;
349        [0 3 6 9]
350
351   uniqind
352       Return the indices of all unique elements of an ndarray The order is in
353       the order of the values to be consistent with uniq. "NaN" values never
354       compare equal with any other value and so are always unique.  This
355       follows the Matlab usage.
356
357         PDL> p pdl(2,2,2,4,0,-1,6,6)->uniqind
358         [5 4 1 3 6]     # the 0 at index 4 is returned 2nd, but...
359
360         PDL> p pdl(2,2,2,4,nan,-1,6,6)->uniqind
361         [5 1 3 6 4]     # ...the NaN at index 4 is returned at end
362
363       Note: The returned pdl is 1D; any structure of the input ndarray is
364       lost.
365
366       See "uniq" if you want the unique values instead of the indices.
367
368       Bad values are not considered unique by uniqind and are ignored.
369
370   uniqvec
371       Return all unique vectors out of a collection
372
373         NOTE: If any vectors in the input ndarray have NaN values
374         they are returned at the end of the non-NaN ones.  This is
375         because, by definition, NaN values never compare equal with
376         any other value.
377
378         NOTE: The current implementation does not sort the vectors
379         containing NaN values.
380
381       The unique vectors are returned in lexicographically sorted ascending
382       order. The 0th dimension of the input PDL is treated as a dimensional
383       index within each vector, and the 1st and any higher dimensions are
384       taken to run across vectors. The return value is always 2D; any
385       structure of the input PDL (beyond using the 0th dimension for vector
386       index) is lost.
387
388       See also "uniq" for a unique list of scalars; and qsortvec for sorting
389       a list of vectors lexicographcally.
390
391       If a vector contains all bad values, it is ignored as in "uniq".  If
392       some of the values are good, it is treated as a normal vector. For
393       example, [1 2 BAD] and [BAD 2 3] could be returned, but [BAD BAD BAD]
394       could not.  Vectors containing BAD values will be returned after any
395       non-NaN and non-BAD containing vectors, followed by the NaN vectors.
396
397   hclip
398         Signature: (a(); b(); [o] c())
399
400       clip (threshold) $a by $b ($b is upper bound)
401
402       hclip processes bad values.  It will set the bad-value flag of all
403       output ndarrays if the flag is set for any of the input ndarrays.
404
405   lclip
406         Signature: (a(); b(); [o] c())
407
408       clip (threshold) $a by $b ($b is lower bound)
409
410       lclip processes bad values.  It will set the bad-value flag of all
411       output ndarrays if the flag is set for any of the input ndarrays.
412
413   clip
414       Clip (threshold) an ndarray by (optional) upper or lower bounds.
415
416        $y = $x->clip(0,3);
417        $c = $x->clip(undef, $x);
418
419       clip handles bad values since it is just a wrapper around "hclip" and
420       "lclip".
421
422   clip
423         Signature: (a(); l(); h(); [o] c())
424
425       info not available
426
427       clip processes bad values.  It will set the bad-value flag of all
428       output ndarrays if the flag is set for any of the input ndarrays.
429
430   wtstat
431         Signature: (a(n); wt(n); avg(); [o]b(); int deg)
432
433       Weighted statistical moment of given degree
434
435       This calculates a weighted statistic over the vector "a".  The formula
436       is
437
438        b() = (sum_i wt_i * (a_i ** degree - avg)) / (sum_i wt_i)
439
440       Bad values are ignored in any calculation; $b will only have its bad
441       flag set if the output contains any bad data.
442
443   statsover
444         Signature: (a(n); w(n); float+ [o]avg(); float+ [o]prms(); int+ [o]median(); int+ [o]min(); int+ [o]max(); float+ [o]adev(); float+ [o]rms())
445
446       Calculate useful statistics over a dimension of an ndarray
447
448         ($mean,$prms,$median,$min,$max,$adev,$rms) = statsover($ndarray, $weights);
449
450       This utility function calculates various useful quantities of an
451       ndarray. These are:
452
453       •  the mean:
454
455            MEAN = sum (x)/ N
456
457          with "N" being the number of elements in x
458
459       •  the population RMS deviation from the mean:
460
461            PRMS = sqrt( sum( (x-mean(x))^2 )/(N-1)
462
463          The population deviation is the best-estimate of the deviation of
464          the population from which a sample is drawn.
465
466       •  the median
467
468          The median is the 50th percentile data value.  Median is found by
469          medover, so WEIGHTING IS IGNORED FOR THE MEDIAN CALCULATION.
470
471       •  the minimum
472
473       •  the maximum
474
475       •  the average absolute deviation:
476
477            AADEV = sum( abs(x-mean(x)) )/N
478
479       •  RMS deviation from the mean:
480
481            RMS = sqrt(sum( (x-mean(x))^2 )/N)
482
483          (also known as the root-mean-square deviation, or the square root of
484          the variance)
485
486       This operator is a projection operator so the calculation will take
487       place over the final dimension. Thus if the input is N-dimensional each
488       returned value will be N-1 dimensional, to calculate the statistics for
489       the entire ndarray either use clump(-1) directly on the ndarray or call
490       "stats".
491
492       Bad values are simply ignored in the calculation, effectively reducing
493       the sample size.  If all data are bad then the output data are marked
494       bad.
495
496   stats
497       Calculates useful statistics on an ndarray
498
499        ($mean,$prms,$median,$min,$max,$adev,$rms) = stats($ndarray,[$weights]);
500
501       This utility calculates all the most useful quantities in one call.  It
502       works the same way as "statsover", except that the quantities are
503       calculated considering the entire input PDL as a single sample, rather
504       than as a collection of rows. See "statsover" for definitions of the
505       returned quantities.
506
507       Bad values are handled; if all input values are bad, then all of the
508       output values are flagged bad.
509
510   histogram
511         Signature: (in(n); int+[o] hist(m); double step; double min; int msize => m)
512
513       Calculates a histogram for given stepsize and minimum.
514
515        $h = histogram($data, $step, $min, $numbins);
516        $hist = zeroes $numbins;  # Put histogram in existing ndarray.
517        histogram($data, $hist, $step, $min, $numbins);
518
519       The histogram will contain $numbins bins starting from $min, each $step
520       wide. The value in each bin is the number of values in $data that lie
521       within the bin limits.
522
523       Data below the lower limit is put in the first bin, and data above the
524       upper limit is put in the last bin.
525
526       The output is reset in a different broadcastloop so that you can take a
527       histogram of "$a(10,12)" into $b(15) and get the result you want.
528
529       For a higher-level interface, see hist.
530
531        pdl> p histogram(pdl(1,1,2),1,0,3)
532        [0 2 1]
533
534       histogram processes bad values.  It will set the bad-value flag of all
535       output ndarrays if the flag is set for any of the input ndarrays.
536
537   whistogram
538         Signature: (in(n); float+ wt(n);float+[o] hist(m); double step; double min; int msize => m)
539
540       Calculates a histogram from weighted data for given stepsize and
541       minimum.
542
543        $h = whistogram($data, $weights, $step, $min, $numbins);
544        $hist = zeroes $numbins;  # Put histogram in existing ndarray.
545        whistogram($data, $weights, $hist, $step, $min, $numbins);
546
547       The histogram will contain $numbins bins starting from $min, each $step
548       wide. The value in each bin is the sum of the values in $weights that
549       correspond to values in $data that lie within the bin limits.
550
551       Data below the lower limit is put in the first bin, and data above the
552       upper limit is put in the last bin.
553
554       The output is reset in a different broadcastloop so that you can take a
555       histogram of "$a(10,12)" into $b(15) and get the result you want.
556
557        pdl> p whistogram(pdl(1,1,2), pdl(0.1,0.1,0.5), 1, 0, 4)
558        [0 0.2 0.5 0]
559
560       whistogram processes bad values.  It will set the bad-value flag of all
561       output ndarrays if the flag is set for any of the input ndarrays.
562
563   histogram2d
564         Signature: (ina(n); inb(n); int+[o] hist(ma,mb); double stepa; double mina; int masize => ma;
565                            double stepb; double minb; int mbsize => mb;)
566
567       Calculates a 2d histogram.
568
569        $h = histogram2d($datax, $datay, $stepx, $minx,
570              $nbinx, $stepy, $miny, $nbiny);
571        $hist = zeroes $nbinx, $nbiny;  # Put histogram in existing ndarray.
572        histogram2d($datax, $datay, $hist, $stepx, $minx,
573              $nbinx, $stepy, $miny, $nbiny);
574
575       The histogram will contain $nbinx x $nbiny bins, with the lower limits
576       of the first one at "($minx, $miny)", and with bin size "($stepx,
577       $stepy)".  The value in each bin is the number of values in $datax and
578       $datay that lie within the bin limits.
579
580       Data below the lower limit is put in the first bin, and data above the
581       upper limit is put in the last bin.
582
583        pdl> p histogram2d(pdl(1,1,1,2,2),pdl(2,1,1,1,1),1,0,3,1,0,3)
584        [
585         [0 0 0]
586         [0 2 2]
587         [0 1 0]
588        ]
589
590       histogram2d processes bad values.  It will set the bad-value flag of
591       all output ndarrays if the flag is set for any of the input ndarrays.
592
593   whistogram2d
594         Signature: (ina(n); inb(n); float+ wt(n);float+[o] hist(ma,mb); double stepa; double mina; int masize => ma;
595                            double stepb; double minb; int mbsize => mb;)
596
597       Calculates a 2d histogram from weighted data.
598
599        $h = whistogram2d($datax, $datay, $weights,
600              $stepx, $minx, $nbinx, $stepy, $miny, $nbiny);
601        $hist = zeroes $nbinx, $nbiny;  # Put histogram in existing ndarray.
602        whistogram2d($datax, $datay, $weights, $hist,
603              $stepx, $minx, $nbinx, $stepy, $miny, $nbiny);
604
605       The histogram will contain $nbinx x $nbiny bins, with the lower limits
606       of the first one at "($minx, $miny)", and with bin size "($stepx,
607       $stepy)".  The value in each bin is the sum of the values in $weights
608       that correspond to values in $datax and $datay that lie within the bin
609       limits.
610
611       Data below the lower limit is put in the first bin, and data above the
612       upper limit is put in the last bin.
613
614        pdl> p whistogram2d(pdl(1,1,1,2,2),pdl(2,1,1,1,1),pdl(0.1,0.2,0.3,0.4,0.5),1,0,3,1,0,3)
615        [
616         [  0   0   0]
617         [  0 0.5 0.9]
618         [  0 0.1   0]
619        ]
620
621       whistogram2d processes bad values.  It will set the bad-value flag of
622       all output ndarrays if the flag is set for any of the input ndarrays.
623
624   fibonacci
625         Signature: (i(n); indx [o]x(n))
626
627       Constructor - a vector with Fibonacci's sequence
628
629       fibonacci does not process bad values.  It will set the bad-value flag
630       of all output ndarrays if the flag is set for any of the input
631       ndarrays.
632
633   append
634         Signature: (a(n); b(m); [o] c(mn))
635
636       append two ndarrays by concatenating along their first dimensions
637
638        $x = ones(2,4,7);
639        $y = sequence 5;
640        $c = $x->append($y);  # size of $c is now (7,4,7) (a jumbo-ndarray ;)
641
642       "append" appends two ndarrays along their first dimensions. The rest of
643       the dimensions must be compatible in the broadcasting sense. The
644       resulting size of the first dimension is the sum of the sizes of the
645       first dimensions of the two argument ndarrays - i.e. "n + m".
646
647       Similar functions include "glue" (below), which can append more than
648       two ndarrays along an arbitrary dimension, and cat, which can append
649       more than two ndarrays that all have the same sized dimensions.
650
651       append does not process bad values.  It will set the bad-value flag of
652       all output ndarrays if the flag is set for any of the input ndarrays.
653
654   glue
655         $c = $x->glue(<dim>,$y,...)
656
657       Glue two or more PDLs together along an arbitrary dimension (N-D
658       "append").
659
660       Sticks $x, $y, and all following arguments together along the specified
661       dimension.  All other dimensions must be compatible in the broadcasting
662       sense.
663
664       Glue is permissive, in the sense that every PDL is treated as having an
665       infinite number of trivial dimensions of order 1 -- so "$x->glue(3,$y)"
666       works, even if $x and $y are only one dimensional.
667
668       If one of the PDLs has no elements, it is ignored.  Likewise, if one of
669       them is actually the undefined value, it is treated as if it had no
670       elements.
671
672       If the first parameter is a defined perl scalar rather than a pdl, then
673       it is taken as a dimension along which to glue everything else, so you
674       can say "$cube = PDL::glue(3,@image_list);" if you like.
675
676       "glue" is implemented in pdl, using a combination of xchg and "append".
677       It should probably be updated (one day) to a pure PP function.
678
679       Similar functions include "append" (above), which appends only two
680       ndarrays along their first dimension, and cat, which can append more
681       than two ndarrays that all have the same sized dimensions.
682
683   cmpvec
684         Signature: (a(n); b(n); sbyte [o]c())
685
686       Compare two vectors lexicographically.
687
688       Returns -1 if a is less, 1 if greater, 0 if equal.
689
690       The output is bad if any input values up to the point of inequality are
691       bad - any after are ignored.
692
693   eqvec
694         Signature: (a(n); b(n); sbyte [o]c())
695
696       Compare two vectors, returning 1 if equal, 0 if not equal.
697
698       The output is bad if any input values are bad.
699
700   enumvec
701         Signature: (v(M,N); indx [o]k(N))
702
703       Enumerate a list of vectors with locally unique keys.
704
705       Given a sorted list of vectors $v, generate a vector $k containing
706       locally unique keys for the elements of $v (where an "element" is a
707       vector of length $M occurring in $v).
708
709       Note that the keys returned in $k are only unique over a run of a
710       single vector in $v, so that each unique vector in $v has at least one
711       0 (zero) index in $k associated with it.  If you need global keys, see
712       enumvecg().
713
714       Contributed by Bryan Jurish <moocow@cpan.org>.
715
716       enumvec does not process bad values.  It will set the bad-value flag of
717       all output ndarrays if the flag is set for any of the input ndarrays.
718
719   enumvecg
720         Signature: (v(M,N); indx [o]k(N))
721
722       Enumerate a list of vectors with globally unique keys.
723
724       Given a sorted list of vectors $v, generate a vector $k containing
725       globally unique keys for the elements of $v (where an "element" is a
726       vector of length $M occurring in $v).  Basically does the same thing
727       as:
728
729        $k = $v->vsearchvec($v->uniqvec);
730
731       ... but somewhat more efficiently.
732
733       Contributed by Bryan Jurish <moocow@cpan.org>.
734
735       enumvecg does not process bad values.  It will set the bad-value flag
736       of all output ndarrays if the flag is set for any of the input
737       ndarrays.
738
739   vsearchvec
740         Signature: (find(M); which(M,N); indx [o]found())
741
742       Routine for searching N-dimensional values - akin to vsearch() for
743       vectors.
744
745        $found   = vsearchvec($find, $which);
746        $nearest = $which->dice_axis(1,$found);
747
748       Returns for each row-vector in $find the index along dimension N of the
749       least row vector of $which greater or equal to it.  $which should be
750       sorted in increasing order.  If the value of $find is larger than any
751       member of $which, the index to the last element of $which is returned.
752
753       See also: "vsearch".  Contributed by Bryan Jurish <moocow@cpan.org>.
754
755       vsearchvec does not process bad values.  It will set the bad-value flag
756       of all output ndarrays if the flag is set for any of the input
757       ndarrays.
758
759   unionvec
760         Signature: (a(M,NA); b(M,NB); [o]c(M,NC); indx [o]nc())
761
762       Union of two vector-valued PDLs.
763
764       Input PDLs $a() and $b() MUST be sorted in lexicographic order.  On
765       return, $nc() holds the actual number of vector-values in the union.
766
767       In scalar context, slices $c() to the actual number of elements in the
768       union and returns the sliced PDL.
769
770       Contributed by Bryan Jurish <moocow@cpan.org>.
771
772       unionvec does not process bad values.  It will set the bad-value flag
773       of all output ndarrays if the flag is set for any of the input
774       ndarrays.
775
776   intersectvec
777         Signature: (a(M,NA); b(M,NB); [o]c(M,NC); indx [o]nc())
778
779       Intersection of two vector-valued PDLs.  Input PDLs $a() and $b() MUST
780       be sorted in lexicographic order.  On return, $nc() holds the actual
781       number of vector-values in the intersection.
782
783       In scalar context, slices $c() to the actual number of elements in the
784       intersection and returns the sliced PDL.
785
786       Contributed by Bryan Jurish <moocow@cpan.org>.
787
788       intersectvec does not process bad values.  It will set the bad-value
789       flag of all output ndarrays if the flag is set for any of the input
790       ndarrays.
791
792   setdiffvec
793         Signature: (a(M,NA); b(M,NB); [o]c(M,NC); indx [o]nc())
794
795       Set-difference ($a() \ $b()) of two vector-valued PDLs.
796
797       Input PDLs $a() and $b() MUST be sorted in lexicographic order.  On
798       return, $nc() holds the actual number of vector-values in the computed
799       vector set.
800
801       In scalar context, slices $c() to the actual number of elements in the
802       output vector set and returns the sliced PDL.
803
804       Contributed by Bryan Jurish <moocow@cpan.org>.
805
806       setdiffvec does not process bad values.  It will set the bad-value flag
807       of all output ndarrays if the flag is set for any of the input
808       ndarrays.
809
810   union_sorted
811         Signature: (a(NA); b(NB); [o]c(NC); indx [o]nc())
812
813       Union of two flat sorted unique-valued PDLs.  Input PDLs $a() and $b()
814       MUST be sorted in lexicographic order and contain no duplicates.  On
815       return, $nc() holds the actual number of values in the union.
816
817       In scalar context, reshapes $c() to the actual number of elements in
818       the union and returns it.
819
820       Contributed by Bryan Jurish <moocow@cpan.org>.
821
822       union_sorted does not process bad values.  It will set the bad-value
823       flag of all output ndarrays if the flag is set for any of the input
824       ndarrays.
825
826   intersect_sorted
827         Signature: (a(NA); b(NB); [o]c(NC); indx [o]nc())
828
829       Intersection of two flat sorted unique-valued PDLs.  Input PDLs $a()
830       and $b() MUST be sorted in lexicographic order and contain no
831       duplicates.  On return, $nc() holds the actual number of values in the
832       intersection.
833
834       In scalar context, reshapes $c() to the actual number of elements in
835       the intersection and returns it.
836
837       Contributed by Bryan Jurish <moocow@cpan.org>.
838
839       intersect_sorted does not process bad values.  It will set the bad-
840       value flag of all output ndarrays if the flag is set for any of the
841       input ndarrays.
842
843   setdiff_sorted
844         Signature: (a(NA); b(NB); [o]c(NC); indx [o]nc())
845
846       Set-difference ($a() \ $b()) of two flat sorted unique-valued PDLs.
847
848       Input PDLs $a() and $b() MUST be sorted in lexicographic order and
849       contain no duplicate values.  On return, $nc() holds the actual number
850       of values in the computed vector set.
851
852       In scalar context, reshapes $c() to the actual number of elements in
853       the difference set and returns it.
854
855       Contributed by Bryan Jurish <moocow@cpan.org>.
856
857       setdiff_sorted does not process bad values.  It will set the bad-value
858       flag of all output ndarrays if the flag is set for any of the input
859       ndarrays.
860
861   srand
862         Signature: (a())
863
864       Seed random-number generator with a 64-bit int. Will generate seed data
865       for a number of threads equal to the return-value of "online_cpus" in
866       PDL::Core.
867
868        srand(); # uses current time
869        srand(5); # fixed number e.g. for testing
870
871       srand does not process bad values.  It will set the bad-value flag of
872       all output ndarrays if the flag is set for any of the input ndarrays.
873
874   random
875         Signature: (a())
876
877       Constructor which returns ndarray of random numbers
878
879        $x = random([type], $nx, $ny, $nz,...);
880        $x = random $y;
881
882       etc (see zeroes).
883
884       This is the uniform distribution between 0 and 1 (assumedly excluding 1
885       itself). The arguments are the same as "zeroes" (q.v.) - i.e. one can
886       specify dimensions, types or give a template.
887
888       You can use the PDL function "srand" to seed the random generator.  If
889       it has not been called yet, it will be with the current time.
890
891       random does not process bad values.  It will set the bad-value flag of
892       all output ndarrays if the flag is set for any of the input ndarrays.
893
894   randsym
895         Signature: (a())
896
897       Constructor which returns ndarray of random numbers
898
899        $x = randsym([type], $nx, $ny, $nz,...);
900        $x = randsym $y;
901
902       etc (see zeroes).
903
904       This is the uniform distribution between 0 and 1 (excluding both 0 and
905       1, cf "random"). The arguments are the same as "zeroes" (q.v.) - i.e.
906       one can specify dimensions, types or give a template.
907
908       You can use the PDL function "srand" to seed the random generator.  If
909       it has not been called yet, it will be with the current time.
910
911       randsym does not process bad values.  It will set the bad-value flag of
912       all output ndarrays if the flag is set for any of the input ndarrays.
913
914   grandom
915       Constructor which returns ndarray of Gaussian random numbers
916
917        $x = grandom([type], $nx, $ny, $nz,...);
918        $x = grandom $y;
919
920       etc (see zeroes).
921
922       This is generated using the math library routine "ndtri".
923
924       Mean = 0, Stddev = 1
925
926       You can use the PDL function "srand" to seed the random generator.  If
927       it has not been called yet, it will be with the current time.
928
929   vsearch
930         Signature: ( vals(); xs(n); [o] indx(); [\%options] )
931
932       Efficiently search for values in a sorted ndarray, returning indices.
933
934         $idx = vsearch( $vals, $x, [\%options] );
935         vsearch( $vals, $x, $idx, [\%options ] );
936
937       vsearch performs a binary search in the ordered ndarray $x, for the
938       values from $vals ndarray, returning indices into $x.  What is a
939       "match", and the meaning of the returned indices, are determined by the
940       options.
941
942       The "mode" option indicates which method of searching to use, and may
943       be one of:
944
945       "sample"
946           invoke vsearch_sample, returning indices appropriate for sampling
947           within a distribution.
948
949       "insert_leftmost"
950           invoke vsearch_insert_leftmost, returning the left-most possible
951           insertion point which still leaves the ndarray sorted.
952
953       "insert_rightmost"
954           invoke vsearch_insert_rightmost, returning the right-most possible
955           insertion point which still leaves the ndarray sorted.
956
957       "match"
958           invoke vsearch_match, returning the index of a matching element,
959           else -(insertion point + 1)
960
961       "bin_inclusive"
962           invoke vsearch_bin_inclusive, returning an index appropriate for
963           binning on a grid where the left bin edges are inclusive of the
964           bin. See below for further explanation of the bin.
965
966       "bin_exclusive"
967           invoke vsearch_bin_exclusive, returning an index appropriate for
968           binning on a grid where the left bin edges are exclusive of the
969           bin. See below for further explanation of the bin.
970
971       The default value of "mode" is "sample".
972
973         use PDL;
974
975         my @modes = qw( sample insert_leftmost insert_rightmost match
976                         bin_inclusive bin_exclusive );
977
978         # Generate a sequence of 3 zeros, 3 ones, ..., 3 fours.
979         my $x = zeroes(3,5)->yvals->flat;
980
981         for my $mode ( @modes ) {
982           # if the value is in $x
983           my $contained = 2;
984           my $idx_contained = vsearch( $contained, $x, { mode => $mode } );
985           my $x_contained = $x->copy;
986           $x_contained->slice( $idx_contained ) .= 9;
987
988           # if the value is not in $x
989           my $not_contained = 1.5;
990           my $idx_not_contained = vsearch( $not_contained, $x, { mode => $mode } );
991           my $x_not_contained = $x->copy;
992           $x_not_contained->slice( $idx_not_contained ) .= 9;
993
994           print sprintf("%-23s%30s\n", '$x', $x);
995           print sprintf("%-23s%30s\n",   "$mode ($contained)", $x_contained);
996           print sprintf("%-23s%30s\n\n", "$mode ($not_contained)", $x_not_contained);
997         }
998
999         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1000         # sample (2)             [0 0 0 1 1 1 9 2 2 3 3 3 4 4 4]
1001         # sample (1.5)           [0 0 0 1 1 1 9 2 2 3 3 3 4 4 4]
1002         #
1003         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1004         # insert_leftmost (2)    [0 0 0 1 1 1 9 2 2 3 3 3 4 4 4]
1005         # insert_leftmost (1.5)  [0 0 0 1 1 1 9 2 2 3 3 3 4 4 4]
1006         #
1007         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1008         # insert_rightmost (2)   [0 0 0 1 1 1 2 2 2 9 3 3 4 4 4]
1009         # insert_rightmost (1.5) [0 0 0 1 1 1 9 2 2 3 3 3 4 4 4]
1010         #
1011         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1012         # match (2)              [0 0 0 1 1 1 2 9 2 3 3 3 4 4 4]
1013         # match (1.5)            [0 0 0 1 1 1 2 2 9 3 3 3 4 4 4]
1014         #
1015         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1016         # bin_inclusive (2)      [0 0 0 1 1 1 2 2 9 3 3 3 4 4 4]
1017         # bin_inclusive (1.5)    [0 0 0 1 1 9 2 2 2 3 3 3 4 4 4]
1018         #
1019         # $x                     [0 0 0 1 1 1 2 2 2 3 3 3 4 4 4]
1020         # bin_exclusive (2)      [0 0 0 1 1 9 2 2 2 3 3 3 4 4 4]
1021         # bin_exclusive (1.5)    [0 0 0 1 1 9 2 2 2 3 3 3 4 4 4]
1022
1023       Also see vsearch_sample, vsearch_insert_leftmost,
1024       vsearch_insert_rightmost, vsearch_match, vsearch_bin_inclusive, and
1025       vsearch_bin_exclusive
1026
1027   vsearch_sample
1028         Signature: (vals(); x(n); indx [o]idx())
1029
1030       Search for values in a sorted array, return index appropriate for
1031       sampling from a distribution
1032
1033         $idx = vsearch_sample($vals, $x);
1034
1035       $x must be sorted, but may be in decreasing or increasing order.
1036
1037       vsearch_sample returns an index I for each value V of $vals appropriate
1038       for sampling $vals
1039
1040       I has the following properties:
1041
1042       •   if $x is sorted in increasing order
1043
1044                     V <= x[0]  : I = 0
1045             x[0]  < V <= x[-1] : I s.t. x[I-1] < V <= x[I]
1046             x[-1] < V          : I = $x->nelem -1
1047
1048       •   if $x is sorted in decreasing order
1049
1050                      V > x[0]  : I = 0
1051             x[0]  >= V > x[-1] : I s.t. x[I] >= V > x[I+1]
1052             x[-1] >= V         : I = $x->nelem - 1
1053
1054       If all elements of $x are equal, I = $x->nelem - 1.
1055
1056       If $x contains duplicated elements, I is the index of the leftmost (by
1057       position in array) duplicate if V matches.
1058
1059       This function is useful e.g. when you have a list of probabilities for
1060       events and want to generate indices to events:
1061
1062        $x = pdl(.01,.86,.93,1); # Barnsley IFS probabilities cumulatively
1063        $y = random 20;
1064        $c = vsearch_sample($y, $x); # Now, $c will have the appropriate distr.
1065
1066       It is possible to use the cumusumover function to obtain cumulative
1067       probabilities from absolute probabilities.
1068
1069       needs major (?) work to handles bad values
1070
1071   vsearch_insert_leftmost
1072         Signature: (vals(); x(n); indx [o]idx())
1073
1074       Determine the insertion point for values in a sorted array, inserting
1075       before duplicates.
1076
1077         $idx = vsearch_insert_leftmost($vals, $x);
1078
1079       $x must be sorted, but may be in decreasing or increasing order.
1080
1081       vsearch_insert_leftmost returns an index I for each value V of $vals
1082       equal to the leftmost position (by index in array) within $x that V may
1083       be inserted and still maintain the order in $x.
1084
1085       Insertion at index I involves shifting elements I and higher of $x to
1086       the right by one and setting the now empty element at index I to V.
1087
1088       I has the following properties:
1089
1090       •   if $x is sorted in increasing order
1091
1092                     V <= x[0]  : I = 0
1093             x[0]  < V <= x[-1] : I s.t. x[I-1] < V <= x[I]
1094             x[-1] < V          : I = $x->nelem
1095
1096       •   if $x is sorted in decreasing order
1097
1098                      V >  x[0]  : I = -1
1099             x[0]  >= V >= x[-1] : I s.t. x[I] >= V > x[I+1]
1100             x[-1] >= V          : I = $x->nelem -1
1101
1102       If all elements of $x are equal,
1103
1104           i = 0
1105
1106       If $x contains duplicated elements, I is the index of the leftmost (by
1107       index in array) duplicate if V matches.
1108
1109       needs major (?) work to handles bad values
1110
1111   vsearch_insert_rightmost
1112         Signature: (vals(); x(n); indx [o]idx())
1113
1114       Determine the insertion point for values in a sorted array, inserting
1115       after duplicates.
1116
1117         $idx = vsearch_insert_rightmost($vals, $x);
1118
1119       $x must be sorted, but may be in decreasing or increasing order.
1120
1121       vsearch_insert_rightmost returns an index I for each value V of $vals
1122       equal to the rightmost position (by index in array) within $x that V
1123       may be inserted and still maintain the order in $x.
1124
1125       Insertion at index I involves shifting elements I and higher of $x to
1126       the right by one and setting the now empty element at index I to V.
1127
1128       I has the following properties:
1129
1130       •   if $x is sorted in increasing order
1131
1132                      V < x[0]  : I = 0
1133             x[0]  <= V < x[-1] : I s.t. x[I-1] <= V < x[I]
1134             x[-1] <= V         : I = $x->nelem
1135
1136       •   if $x is sorted in decreasing order
1137
1138                     V >= x[0]  : I = -1
1139             x[0]  > V >= x[-1] : I s.t. x[I] >= V > x[I+1]
1140             x[-1] > V          : I = $x->nelem -1
1141
1142       If all elements of $x are equal,
1143
1144           i = $x->nelem - 1
1145
1146       If $x contains duplicated elements, I is the index of the leftmost (by
1147       index in array) duplicate if V matches.
1148
1149       needs major (?) work to handles bad values
1150
1151   vsearch_match
1152         Signature: (vals(); x(n); indx [o]idx())
1153
1154       Match values against a sorted array.
1155
1156         $idx = vsearch_match($vals, $x);
1157
1158       $x must be sorted, but may be in decreasing or increasing order.
1159
1160       vsearch_match returns an index I for each value V of $vals.  If V
1161       matches an element in $x, I is the index of that element, otherwise it
1162       is -( insertion_point + 1 ), where insertion_point is an index in $x
1163       where V may be inserted while maintaining the order in $x.  If $x has
1164       duplicated values, I may refer to any of them.
1165
1166       needs major (?) work to handles bad values
1167
1168   vsearch_bin_inclusive
1169         Signature: (vals(); x(n); indx [o]idx())
1170
1171       Determine the index for values in a sorted array of bins, lower bound
1172       inclusive.
1173
1174         $idx = vsearch_bin_inclusive($vals, $x);
1175
1176       $x must be sorted, but may be in decreasing or increasing order.
1177
1178       $x represents the edges of contiguous bins, with the first and last
1179       elements representing the outer edges of the outer bins, and the inner
1180       elements the shared bin edges.
1181
1182       The lower bound of a bin is inclusive to the bin, its outer bound is
1183       exclusive to it.  vsearch_bin_inclusive returns an index I for each
1184       value V of $vals
1185
1186       I has the following properties:
1187
1188       •   if $x is sorted in increasing order
1189
1190                      V < x[0]  : I = -1
1191             x[0]  <= V < x[-1] : I s.t. x[I] <= V < x[I+1]
1192             x[-1] <= V         : I = $x->nelem - 1
1193
1194       •   if $x is sorted in decreasing order
1195
1196                      V >= x[0]  : I = 0
1197             x[0]  >  V >= x[-1] : I s.t. x[I+1] > V >= x[I]
1198             x[-1] >  V          : I = $x->nelem
1199
1200       If all elements of $x are equal,
1201
1202           i = $x->nelem - 1
1203
1204       If $x contains duplicated elements, I is the index of the righmost (by
1205       index in array) duplicate if V matches.
1206
1207       needs major (?) work to handles bad values
1208
1209   vsearch_bin_exclusive
1210         Signature: (vals(); x(n); indx [o]idx())
1211
1212       Determine the index for values in a sorted array of bins, lower bound
1213       exclusive.
1214
1215         $idx = vsearch_bin_exclusive($vals, $x);
1216
1217       $x must be sorted, but may be in decreasing or increasing order.
1218
1219       $x represents the edges of contiguous bins, with the first and last
1220       elements representing the outer edges of the outer bins, and the inner
1221       elements the shared bin edges.
1222
1223       The lower bound of a bin is exclusive to the bin, its upper bound is
1224       inclusive to it.  vsearch_bin_exclusive returns an index I for each
1225       value V of $vals.
1226
1227       I has the following properties:
1228
1229       •   if $x is sorted in increasing order
1230
1231                      V <= x[0]  : I = -1
1232             x[0]  <  V <= x[-1] : I s.t. x[I] < V <= x[I+1]
1233             x[-1] <  V          : I = $x->nelem - 1
1234
1235       •   if $x is sorted in decreasing order
1236
1237                      V >  x[0]  : I = 0
1238             x[0]  >= V >  x[-1] : I s.t. x[I-1] >= V > x[I]
1239             x[-1] >= V          : I = $x->nelem
1240
1241       If all elements of $x are equal,
1242
1243           i = $x->nelem - 1
1244
1245       If $x contains duplicated elements, I is the index of the righmost (by
1246       index in array) duplicate if V matches.
1247
1248       needs major (?) work to handles bad values
1249
1250   interpolate
1251         Signature: (real xi(); real x(n); y(n); [o] yi(); int [o] err())
1252
1253       routine for 1D linear interpolation
1254
1255        ( $yi, $err ) = interpolate($xi, $x, $y)
1256
1257       Given a set of points "($x,$y)", use linear interpolation to find the
1258       values $yi at a set of points $xi.
1259
1260       "interpolate" uses a binary search to find the suspects, er...,
1261       interpolation indices and therefore abscissas (ie $x) have to be
1262       strictly ordered (increasing or decreasing).  For interpolation at lots
1263       of closely spaced abscissas an approach that uses the last index found
1264       as a start for the next search can be faster (compare Numerical Recipes
1265       "hunt" routine). Feel free to implement that on top of the binary
1266       search if you like. For out of bounds values it just does a linear
1267       extrapolation and sets the corresponding element of $err to 1, which is
1268       otherwise 0.
1269
1270       See also "interpol", which uses the same routine, differing only in the
1271       handling of extrapolation - an error message is printed rather than
1272       returning an error ndarray.
1273
1274       Note that "interpolate" can use complex values for $y and $yi but $x
1275       and $xi must be real.
1276
1277       needs major (?) work to handles bad values
1278
1279   interpol
1280        Signature: (xi(); x(n); y(n); [o] yi())
1281
1282       routine for 1D linear interpolation
1283
1284        $yi = interpol($xi, $x, $y)
1285
1286       "interpol" uses the same search method as "interpolate", hence $x must
1287       be strictly ordered (either increasing or decreasing).  The difference
1288       occurs in the handling of out-of-bounds values; here an error message
1289       is printed.
1290
1291   interpND
1292       Interpolate values from an N-D ndarray, with switchable method
1293
1294         $source = 10*xvals(10,10) + yvals(10,10);
1295         $index = pdl([[2.2,3.5],[4.1,5.0]],[[6.0,7.4],[8,9]]);
1296         print $source->interpND( $index );
1297
1298       InterpND acts like indexND, collapsing $index by lookup into $source;
1299       but it does interpolation rather than direct sampling.  The
1300       interpolation method and boundary condition are switchable via an
1301       options hash.
1302
1303       By default, linear or sample interpolation is used, with constant value
1304       outside the boundaries of the source pdl.  No dataflow occurs, because
1305       in general the output is computed rather than indexed.
1306
1307       All the interpolation methods treat the pixels as value-centered, so
1308       the "sample" method will return $a->(0) for coordinate values on the
1309       set [-0.5,0.5), and all methods will return $a->(1) for a coordinate
1310       value of exactly 1.
1311
1312       Recognized options:
1313
1314       method
1315          Values can be:
1316
1317          •  0, s, sample, Sample (default for integer source types)
1318
1319             The nearest value is taken. Pixels are regarded as centered on
1320             their respective integer coordinates (no offset from the linear
1321             case).
1322
1323          •  1, l, linear, Linear (default for floating point source types)
1324
1325             The values are N-linearly interpolated from an N-dimensional cube
1326             of size 2.
1327
1328          •  3, c, cube, cubic, Cubic
1329
1330             The values are interpolated using a local cubic fit to the data.
1331             The fit is constrained to match the original data and its
1332             derivative at the data points.  The second derivative of the fit
1333             is not continuous at the data points.  Multidimensional datasets
1334             are interpolated by the successive-collapse method.
1335
1336             (Note that the constraint on the first derivative causes a small
1337             amount of ringing around sudden features such as step functions).
1338
1339          •  f, fft, fourier, Fourier
1340
1341             The source is Fourier transformed, and the interpolated values
1342             are explicitly calculated from the coefficients.  The boundary
1343             condition option is ignored -- periodic boundaries are imposed.
1344
1345             If you pass in the option "fft", and it is a list (ARRAY) ref,
1346             then it is a stash for the magnitude and phase of the source FFT.
1347             If the list has two elements then they are taken as already
1348             computed; otherwise they are calculated and put in the stash.
1349
1350       b, bound, boundary, Boundary
1351          This option is passed unmodified into indexND, which is used as the
1352          indexing engine for the interpolation.  Some current allowed values
1353          are 'extend', 'periodic', 'truncate', and 'mirror' (default is
1354          'truncate').
1355
1356       bad
1357          contains the fill value used for 'truncate' boundary.  (default 0)
1358
1359       fft
1360          An array ref whose associated list is used to stash the FFT of the
1361          source data, for the FFT method.
1362
1363   one2nd
1364       Converts a one dimensional index ndarray to a set of ND coordinates
1365
1366        @coords=one2nd($x, $indices)
1367
1368       returns an array of ndarrays containing the ND indexes corresponding to
1369       the one dimensional list indices. The indices are assumed to correspond
1370       to array $x clumped using clump(-1). This routine is used in the old
1371       vector form of "whichND", but is useful on its own occasionally.
1372
1373       Returned ndarrays have the indx datatype.  $indices can have values
1374       larger than "$x->nelem" but negative values in $indices will not give
1375       the answer you expect.
1376
1377        pdl> $x=pdl [[[1,2],[-1,1]], [[0,-3],[3,2]]]; $c=$x->clump(-1)
1378        pdl> $maxind=maximum_ind($c); p $maxind;
1379        6
1380        pdl> print one2nd($x, maximum_ind($c))
1381        0 1 1
1382        pdl> p $x->at(0,1,1)
1383        3
1384
1385   which
1386         Signature: (mask(n); indx [o] inds(n); indx [o]lastout())
1387
1388       Returns indices of non-zero values from a 1-D PDL
1389
1390        $i = which($mask);
1391
1392       returns a pdl with indices for all those elements that are nonzero in
1393       the mask. Note that the returned indices will be 1D. If you feed in a
1394       multidimensional mask, it will be flattened before the indices are
1395       calculated.  See also "whichND" for multidimensional masks.
1396
1397       If you want to index into the original mask or a similar ndarray with
1398       output from "which", remember to flatten it before calling index:
1399
1400         $data = random 5, 5;
1401         $idx = which $data > 0.5; # $idx is now 1D
1402         $bigsum = $data->flat->index($idx)->sum;  # flatten before indexing
1403
1404       Compare also "where" for similar functionality.
1405
1406       SEE ALSO:
1407
1408       "which_both" returns separately the indices of both nonzero and zero
1409       values in the mask.
1410
1411       "where_both" returns separately slices of both nonzero and zero values
1412       in the mask.
1413
1414       "where" returns associated values from a data PDL, rather than indices
1415       into the mask PDL.
1416
1417       "whichND" returns N-D indices into a multidimensional PDL.
1418
1419        pdl> $x = sequence(10); p $x
1420        [0 1 2 3 4 5 6 7 8 9]
1421        pdl> $indx = which($x>6); p $indx
1422        [7 8 9]
1423
1424       which processes bad values.  It will set the bad-value flag of all
1425       output ndarrays if the flag is set for any of the input ndarrays.
1426
1427   which_both
1428         Signature: (mask(n); indx [o] inds(n); indx [o]notinds(n); indx [o]lastout(); indx [o]lastoutn())
1429
1430       Returns indices of nonzero and zero values in a mask PDL
1431
1432        ($i, $c_i) = which_both($mask);
1433
1434       This works just as "which", but the complement of $i will be in $c_i.
1435
1436        pdl> p $x = sequence(10)
1437        [0 1 2 3 4 5 6 7 8 9]
1438        pdl> ($big, $small) = which_both($x >= 5); p "$big\n$small"
1439        [5 6 7 8 9]
1440        [0 1 2 3 4]
1441
1442       which_both processes bad values.  It will set the bad-value flag of all
1443       output ndarrays if the flag is set for any of the input ndarrays.
1444
1445   where
1446       Use a mask to select values from one or more data PDLs
1447
1448       "where" accepts one or more data ndarrays and a mask ndarray.  It
1449       returns a list of output ndarrays, corresponding to the input data
1450       ndarrays.  Each output ndarray is a 1-dimensional list of values in its
1451       corresponding data ndarray. The values are drawn from locations where
1452       the mask is nonzero.
1453
1454       The output PDLs are still connected to the original data PDLs, for the
1455       purpose of dataflow.
1456
1457       "where" combines the functionality of "which" and index into a single
1458       operation.
1459
1460       BUGS:
1461
1462       While "where" works OK for most N-dimensional cases, it does not
1463       broadcast properly over (for example) the (N+1)th dimension in data
1464       that is compared to an N-dimensional mask.  Use "whereND" for that.
1465
1466        $i = $x->where($x+5 > 0); # $i contains those elements of $x
1467                                  # where mask ($x+5 > 0) is 1
1468        $i .= -5;  # Set those elements (of $x) to -5. Together, these
1469                   # commands clamp $x to a maximum of -5.
1470
1471       It is also possible to use the same mask for several ndarrays with the
1472       same call:
1473
1474        ($i,$j,$k) = where($x,$y,$z, $x+5>0);
1475
1476       Note: $i is always 1-D, even if $x is >1-D.
1477
1478       WARNING: The first argument (the values) and the second argument (the
1479       mask) currently have to have the exact same dimensions (or horrible
1480       things happen). You *cannot* broadcast over a smaller mask, for
1481       example.
1482
1483   where_both
1484       Returns slices (non-zero in mask, zero) of an ndarray according to a
1485       mask
1486
1487        ($match_vals, $non_match_vals) = where_both($pdl, $mask);
1488
1489       This works like "which_both", but (flattened) data-flowing slices
1490       rather than index-sets are returned.
1491
1492        pdl> p $x = sequence(10) + 2
1493        [2 3 4 5 6 7 8 9 10 11]
1494        pdl> ($big, $small) = where_both($x, $x > 5); p "$big\n$small"
1495        [6 7 8 9 10 11]
1496        [2 3 4 5]
1497        pdl> p $big += 2, $small -= 1
1498        [8 9 10 11 12 13] [1 2 3 4]
1499        pdl> p $x
1500        [1 2 3 4 8 9 10 11 12 13]
1501
1502   whereND
1503       "where" with support for ND masks and broadcasting
1504
1505       "whereND" accepts one or more data ndarrays and a mask ndarray.  It
1506       returns a list of output ndarrays, corresponding to the input data
1507       ndarrays.  The values are drawn from locations where the mask is
1508       nonzero.
1509
1510       "whereND" differs from "where" in that the mask dimensionality is
1511       preserved which allows for proper broadcasting of the selection
1512       operation over higher dimensions.
1513
1514       As with "where" the output PDLs are still connected to the original
1515       data PDLs, for the purpose of dataflow.
1516
1517         $sdata = whereND $data, $mask
1518         ($s1, $s2, ..., $sn) = whereND $d1, $d2, ..., $dn, $mask
1519
1520         where
1521
1522           $data is M dimensional
1523           $mask is N < M dimensional
1524           dims($data) 1..N == dims($mask) 1..N
1525           with broadcasting over N+1 to M dimensions
1526
1527         $data   = sequence(4,3,2);   # example data array
1528         $mask4  = (random(4)>0.5);   # example 1-D mask array, has $n4 true values
1529         $mask43 = (random(4,3)>0.5); # example 2-D mask array, has $n43 true values
1530         $sdat4  = whereND $data, $mask4;   # $sdat4 is a [$n4,3,2] pdl
1531         $sdat43 = whereND $data, $mask43;  # $sdat43 is a [$n43,2] pdl
1532
1533       Just as with "where", you can use the returned value in an assignment.
1534       That means that both of these examples are valid:
1535
1536         # Used to create a new slice stored in $sdat4:
1537         $sdat4 = $data->whereND($mask4);
1538         $sdat4 .= 0;
1539         # Used in lvalue context:
1540         $data->whereND($mask4) .= 0;
1541
1542       SEE ALSO:
1543
1544       "whichND" returns N-D indices into a multidimensional PDL, from a mask.
1545
1546   whichND
1547       Return the coordinates of non-zero values in a mask.
1548
1549       WhichND returns the N-dimensional coordinates of each nonzero value in
1550       a mask PDL with any number of dimensions.  The returned values arrive
1551       as an array-of-vectors suitable for use in indexND or range.
1552
1553        $coords = whichND($mask);
1554
1555       returns a PDL containing the coordinates of the elements that are non-
1556       zero in $mask, suitable for use in "indexND" in PDL::Slices. The 0th
1557       dimension contains the full coordinate listing of each point; the 1st
1558       dimension lists all the points.  For example, if $mask has rank 4 and
1559       100 matching elements, then $coords has dimension 4x100.
1560
1561       If no such elements exist, then whichND returns a structured empty PDL:
1562       an Nx0 PDL that contains no values (but matches, broadcasting-wise,
1563       with the vectors that would be produced if such elements existed).
1564
1565       DEPRECATED BEHAVIOR IN LIST CONTEXT:
1566
1567       whichND once delivered different values in list context than in scalar
1568       context, for historical reasons.  In list context, it returned the
1569       coordinates transposed, as a collection of 1-PDLs (one per dimension)
1570       in a list.  This usage is deprecated in PDL 2.4.10, and will cause a
1571       warning to be issued every time it is encountered.  To avoid the
1572       warning, you can set the global variable "$PDL::whichND" to 's' to get
1573       scalar behavior in all contexts, or to 'l' to get list behavior in list
1574       context.
1575
1576       In later versions of PDL, the deprecated behavior will disappear.
1577       Deprecated list context whichND expressions can be replaced with:
1578
1579           @list = $x->whichND->mv(0,-1)->dog;
1580
1581       SEE ALSO:
1582
1583       "which" finds coordinates of nonzero values in a 1-D mask.
1584
1585       "where" extracts values from a data PDL that are associated with
1586       nonzero values in a mask PDL.
1587
1588       "indexND" in PDL::Slices can be fed the coordinates to return the
1589       values.
1590
1591        pdl> $s=sequence(10,10,3,4)
1592        pdl> ($x, $y, $z, $w)=whichND($s == 203); p $x, $y, $z, $w
1593        [3] [0] [2] [0]
1594        pdl> print $s->at(list(cat($x,$y,$z,$w)))
1595        203
1596
1597   setops
1598       Implements simple set operations like union and intersection
1599
1600          Usage: $set = setops($x, <OPERATOR>, $y);
1601
1602       The operator can be "OR", "XOR" or "AND". This is then applied to $x
1603       viewed as a set and $y viewed as a set. Set theory says that a set may
1604       not have two or more identical elements, but setops takes care of this
1605       for you, so "$x=pdl(1,1,2)" is OK. The functioning is as follows:
1606
1607       "OR"
1608           The resulting vector will contain the elements that are either in
1609           $x or in $y or both. This is the union in set operation terms
1610
1611       "XOR"
1612           The resulting vector will contain the elements that are either in
1613           $x or $y, but not in both. This is
1614
1615                Union($x, $y) - Intersection($x, $y)
1616
1617           in set operation terms.
1618
1619       "AND"
1620           The resulting vector will contain the intersection of $x and $y, so
1621           the elements that are in both $x and $y. Note that for convenience
1622           this operation is also aliased to "intersect".
1623
1624       It should be emphasized that these routines are used when one or both
1625       of the sets $x, $y are hard to calculate or that you get from a
1626       separate subroutine.
1627
1628       Finally IDL users might be familiar with Craig Markwardt's
1629       "cmset_op.pro" routine which has inspired this routine although it was
1630       written independently However the present routine has a few less
1631       options (but see the examples)
1632
1633       You will very often use these functions on an index vector, so that is
1634       what we will show here. We will in fact something slightly silly. First
1635       we will find all squares that are also cubes below 10000.
1636
1637       Create a sequence vector:
1638
1639         pdl> $x = sequence(10000)
1640
1641       Find all odd and even elements:
1642
1643         pdl> ($even, $odd) = which_both( ($x % 2) == 0)
1644
1645       Find all squares
1646
1647         pdl> $squares= which(ceil(sqrt($x)) == floor(sqrt($x)))
1648
1649       Find all cubes (being careful with roundoff error!)
1650
1651         pdl> $cubes= which(ceil($x**(1.0/3.0)) == floor($x**(1.0/3.0)+1e-6))
1652
1653       Then find all squares that are cubes:
1654
1655         pdl> $both = setops($squares, 'AND', $cubes)
1656
1657       And print these (assumes that "PDL::NiceSlice" is loaded!)
1658
1659         pdl> p $x($both)
1660          [0 1 64 729 4096]
1661
1662       Then find all numbers that are either cubes or squares, but not both:
1663
1664         pdl> $cube_xor_square = setops($squares, 'XOR', $cubes)
1665
1666         pdl> p $cube_xor_square->nelem()
1667          112
1668
1669       So there are a total of 112 of these!
1670
1671       Finally find all odd squares:
1672
1673         pdl> $odd_squares = setops($squares, 'AND', $odd)
1674
1675       Another common occurrence is to want to get all objects that are in $x
1676       and in the complement of $y. But it is almost always best to create the
1677       complement explicitly since the universe that both are taken from is
1678       not known. Thus use "which_both" if possible to keep track of
1679       complements.
1680
1681       If this is impossible the best approach is to make a temporary:
1682
1683       This creates an index vector the size of the universe of the sets and
1684       set all elements in $y to 0
1685
1686         pdl> $tmp = ones($n_universe); $tmp($y) .= 0;
1687
1688       This then finds the complement of $y
1689
1690         pdl> $C_b = which($tmp == 1);
1691
1692       and this does the final selection:
1693
1694         pdl> $set = setops($x, 'AND', $C_b)
1695
1696   intersect
1697       Calculate the intersection of two ndarrays
1698
1699          Usage: $set = intersect($x, $y);
1700
1701       This routine is merely a simple interface to "setops". See that for
1702       more information
1703
1704       Find all numbers less that 100 that are of the form 2*y and 3*x
1705
1706        pdl> $x=sequence(100)
1707        pdl> $factor2 = which( ($x % 2) == 0)
1708        pdl> $factor3 = which( ($x % 3) == 0)
1709        pdl> $ii=intersect($factor2, $factor3)
1710        pdl> p $x($ii)
1711        [0 6 12 18 24 30 36 42 48 54 60 66 72 78 84 90 96]
1712

AUTHOR

1714       Copyright (C) Tuomas J. Lukka 1997 (lukka@husc.harvard.edu).
1715       Contributions by Christian Soeller (c.soeller@auckland.ac.nz), Karl
1716       Glazebrook (kgb@aaoepp.aao.gov.au), Craig DeForest
1717       (deforest@boulder.swri.edu) and Jarle Brinchmann (jarle@astro.up.pt)
1718       All rights reserved. There is no warranty. You are allowed to
1719       redistribute this software / documentation under certain conditions.
1720       For details, see the file COPYING in the PDL distribution. If this file
1721       is separated from the PDL distribution, the copyright notice should be
1722       included in the file.
1723
1724       Updated for CPAN viewing compatibility by David Mertens.
1725
1726
1727
1728perl v5.36.0                      2023-01-20                      Primitive(3)
Impressum