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

NAME

6       PDL::Slatec - PDL interface to the slatec numerical programming library
7

SYNOPSIS

9        use PDL::Slatec;
10
11        ($ndeg, $r, $ierr, $c) = polyfit($x, $y, $w, $maxdeg, $eps);
12

DESCRIPTION

14       This module serves the dual purpose of providing an interface to parts
15       of the slatec library and showing how to interface PDL to an external
16       library.  Using this library requires a fortran compiler; the source
17       for the routines is provided for convenience.
18
19       Currently available are routines to: manipulate matrices; calculate
20       FFT's; fit data using polynomials; and interpolate/integrate data using
21       piecewise cubic Hermite interpolation.
22
23   Piecewise cubic Hermite interpolation (PCHIP)
24       PCHIP is the slatec package of routines to perform piecewise cubic
25       Hermite interpolation of data.  It features software to produce a
26       monotone and "visually pleasing" interpolant to monotone data.
27       According to Fritsch & Carlson ("Monotone piecewise cubic
28       interpolation", SIAM Journal on Numerical Analysis 17, 2 (April 1980),
29       pp. 238-246), such an interpolant may be more reasonable than a cubic
30       spline if the data contains both "steep" and "flat" sections.
31       Interpolation of cumulative probability distribution functions is
32       another application.  These routines are cryptically named (blame
33       FORTRAN), beginning with 'ch', and accept either float or double
34       ndarrays.
35
36       Most of the routines require an integer parameter called "check"; if
37       set to 0, then no checks on the validity of the input data are made,
38       otherwise these checks are made.  The value of "check" can be set to 0
39       if a routine such as "chim" has already been successfully called.
40
41       •   If not known, estimate derivative values for the points using the
42           "chim", "chic", or "chsp" routines (the following routines require
43           both the function ("f") and derivative ("d") values at a set of
44           points ("x")).
45
46       •   Evaluate, integrate, or differentiate the resulting PCH function
47           using the routines: "chfd"; "chfe"; "chia"; "chid".
48
49       •   If desired, you can check the monotonicity of your data using
50           "chcm".
51

FUNCTIONS

53   eigsys
54       Eigenvalues and eigenvectors of a real positive definite symmetric
55       matrix.
56
57        ($eigvals,$eigvecs) = eigsys($mat)
58
59       Note: this function should be extended to calculate only eigenvalues if
60       called in scalar context!
61
62   matinv
63       Inverse of a square matrix
64
65        ($inv) = matinv($mat)
66
67   polyfit
68       Convenience wrapper routine about the "polfit" "slatec" function.
69       Separates supplied arguments and return values.
70
71       Fit discrete data in a least squares sense by polynomials in one
72       variable.  Handles broadcasting correctly--one can pass in a 2D PDL (as
73       $y) and it will pass back a 2D PDL, the rows of which are the
74       polynomial regression results (in $r corresponding to the rows of $y.
75
76        ($ndeg, $r, $ierr, $c, $coeffs, $rms) = polyfit($x, $y, $w, $maxdeg, [$eps]);
77
78        $coeffs = polyfit($x,$y,$w,$maxdeg,[$eps]);
79
80       where on input:
81
82       $x and $y are the values to fit to a polynomial.  $w are weighting
83       factors $maxdeg is the maximum degree of polynomial to use and $eps is
84       the required degree of fit.
85
86       and the output switches on list/scalar context.
87
88       In list context:
89
90       $ndeg is the degree of polynomial actually used $r is the values of the
91       fitted polynomial $ierr is a return status code, and $c is some working
92       array or other (preserved for historical purposes) $coeffs is the
93       polynomial coefficients of the best fit polynomial.  $rms is the rms
94       error of the fit.
95
96       In scalar context, only $coeffs is returned.
97
98       Historically, $eps was modified in-place to be a return value of the
99       rms error.  This usage is deprecated, and $eps is an optional parameter
100       now.  It is still modified if present.
101
102       $c is a working array accessible to Slatec - you can feed it to several
103       other Slatec routines to get nice things out.  It does not broadcast
104       correctly and should probably be fixed by someone.  If you are reading
105       this, that someone might be you.
106
107       This version of polyfit handles bad values correctly.  Bad values in $y
108       are ignored for the fit and give computed values on the fitted curve in
109       the return.  Bad values in $x or $w are ignored for the fit and result
110       in bad elements in the output.
111
112   polycoef
113       Convenience wrapper routine around the "pcoef" "slatec" function.
114       Separates supplied arguments and return values.
115
116       Convert the "polyfit"/"polfit" coefficients to Taylor series form.
117
118        $tc = polycoef($l, $c, $x);
119
120   polyvalue
121       Convenience wrapper routine around the "pvalue" "slatec" function.
122       Separates supplied arguments and return values.
123
124       For multiple input x positions, a corresponding y position is
125       calculated.
126
127       The derivatives PDL is one dimensional (of size "nder") if a single x
128       position is supplied, two dimensional if more than one x position is
129       supplied.
130
131       Use the coefficients "c" generated by "polyfit" (or "polfit") to
132       evaluate the polynomial fit of degree "l", along with the first "nder"
133       of its derivatives, at a specified point "x".
134
135        ($yfit, $yp) = polyvalue($l, $nder, $x, $c);
136
137   detslatec
138       compute the determinant of an invertible matrix
139
140         $mat = zeroes(5,5); $mat->diagonal(0,1) .= 1; # unity matrix
141         $det = detslatec $mat;
142
143       Usage:
144
145         $determinant = detslatec $matrix;
146
147         Signature: detslatec(mat(n,m); [o] det())
148
149       "detslatec" computes the determinant of an invertible matrix and barfs
150       if the matrix argument provided is non-invertible. The matrix
151       broadcasts as usual.
152
153       This routine was previously known as "det" which clashes now with det
154       which is provided by PDL::MatrixOps.
155
156   fft
157       Fast Fourier Transform
158
159         $v_in = pdl(1,0,1,0);
160         ($azero,$x,$y) = PDL::Slatec::fft($v_in);
161
162       "PDL::Slatec::fft" is a convenience wrapper for "ezfftf", and performs
163       a Fast Fourier Transform on an input vector $v_in. The return values
164       are the same as for "ezfftf".
165
166   rfft
167       reverse Fast Fourier Transform
168
169         $v_out = PDL::Slatec::rfft($azero,$x,$y);
170         print $v_in, $vout
171         [1 0 1 0] [1 0 1 0]
172
173       "PDL::Slatec::rfft" is a convenience wrapper for "ezfftb", and performs
174       a reverse Fast Fourier Transform. The input is the same as the output
175       of "PDL::Slatec::fft", and the output of "rfft" is a data vector,
176       similar to what is input into "PDL::Slatec::fft".
177
178   svdc
179         Signature: (x(n,p);[o]s(p);[o]e(p);[o]u(n,p);[o]v(p,p);[o]work(n);longlong job();longlong [o]info())
180
181       singular value decomposition of a matrix
182
183       svdc does not process bad values.  It will set the bad-value flag of
184       all output ndarrays if the flag is set for any of the input ndarrays.
185
186   poco
187         Signature: (a(n,n);rcond();[o]z(n);longlong [o]info())
188
189       Factor a real symmetric positive definite matrix and estimate the
190       condition number of the matrix.
191
192       poco does not process bad values.  It will set the bad-value flag of
193       all output ndarrays if the flag is set for any of the input ndarrays.
194
195   geco
196         Signature: (a(n,n);longlong [o]ipvt(n);[o]rcond();[o]z(n))
197
198       Factor a matrix using Gaussian elimination and estimate the condition
199       number of the matrix.
200
201       geco does not process bad values.  It will set the bad-value flag of
202       all output ndarrays if the flag is set for any of the input ndarrays.
203
204   gefa
205         Signature: (a(n,n);longlong [o]ipvt(n);longlong [o]info())
206
207       Factor a matrix using Gaussian elimination.
208
209       gefa does not process bad values.  It will set the bad-value flag of
210       all output ndarrays if the flag is set for any of the input ndarrays.
211
212   podi
213         Signature: (a(n,n);[o]det(two=2);longlong job())
214
215       Compute the determinant and inverse of a certain real symmetric
216       positive definite matrix using the factors computed by "poco".
217
218       podi does not process bad values.  It will set the bad-value flag of
219       all output ndarrays if the flag is set for any of the input ndarrays.
220
221   gedi
222         Signature: (a(n,n);longlong [o]ipvt(n);[o]det(two=2);[o]work(n);longlong job())
223
224       Compute the determinant and inverse of a matrix using the factors
225       computed by "geco" or "gefa".
226
227       gedi does not process bad values.  It will set the bad-value flag of
228       all output ndarrays if the flag is set for any of the input ndarrays.
229
230   gesl
231         Signature: (a(lda,n);longlong ipvt(n);b(n);longlong job())
232
233       Solve the real system "A*X=B" or "TRANS(A)*X=B" using the factors
234       computed by "geco" or "gefa".
235
236       gesl does not process bad values.  It will set the bad-value flag of
237       all output ndarrays if the flag is set for any of the input ndarrays.
238
239   rs
240         Signature: (a(n,n);[o]w(n);longlong matz();[o]z(n,n);[t]fvone(n);[t]fvtwo(n);longlong [o]ierr())
241
242       This subroutine calls the recommended sequence of subroutines from the
243       eigensystem subroutine package (EISPACK) to find the eigenvalues and
244       eigenvectors (if desired) of a REAL SYMMETRIC matrix.
245
246       rs does not process bad values.  It will set the bad-value flag of all
247       output ndarrays if the flag is set for any of the input ndarrays.
248
249   ezffti
250         Signature: (longlong n();[o]wsave(foo))
251
252       Subroutine ezffti initializes the work array wsave() which is used in
253       both "ezfftf" and "ezfftb".  The prime factorization of "n" together
254       with a tabulation of the trigonometric functions are computed and
255       stored in wsave().
256
257       ezffti does not process bad values.  It will set the bad-value flag of
258       all output ndarrays if the flag is set for any of the input ndarrays.
259
260   ezfftf
261         Signature: (r(n);[o]azero();[o]a(n);[o]b(n);wsave(foo))
262
263       ezfftf does not process bad values.  It will set the bad-value flag of
264       all output ndarrays if the flag is set for any of the input ndarrays.
265
266   ezfftb
267         Signature: ([o]r(n);azero();a(n);b(n);wsave(foo))
268
269       ezfftb does not process bad values.  It will set the bad-value flag of
270       all output ndarrays if the flag is set for any of the input ndarrays.
271
272   pcoef
273         Signature: (longlong l();c();[o]tc(bar);a(foo))
274
275       Convert the "polfit" coefficients to Taylor series form.  "c" and a()
276       must be of the same type.
277
278       pcoef does not process bad values.  It will set the bad-value flag of
279       all output ndarrays if the flag is set for any of the input ndarrays.
280
281   pvalue
282         Signature: (longlong l();x();[o]yfit();[o]yp(nder);a(foo))
283
284       Use the coefficients generated by "polfit" to evaluate the polynomial
285       fit of degree "l", along with the first "nder" of its derivatives, at a
286       specified point. "x" and "a" must be of the same type.
287
288       pvalue does not process bad values.  It will set the bad-value flag of
289       all output ndarrays if the flag is set for any of the input ndarrays.
290
291   chim
292         Signature: (x(n);f(n);[o]d(n);longlong [o]ierr())
293
294       Calculate the derivatives of (x,f(x)) using cubic Hermite
295       interpolation.
296
297       Calculate the derivatives at the given set of points ("$x,$f", where $x
298       is strictly increasing).  The resulting set of points - "$x,$f,$d",
299       referred to as the cubic Hermite representation - can then be used in
300       other functions, such as "chfe", "chfd", and "chia".
301
302       The boundary conditions are compatible with monotonicity, and if the
303       data are only piecewise monotonic, the interpolant will have an
304       extremum at the switch points; for more control over these issues use
305       "chic".
306
307       Error status returned by $ierr:
308
309       •   0 if successful.
310
311       •   > 0 if there were "ierr" switches in the direction of monotonicity
312           (data still valid).
313
314       •   -1 if "nelem($x) < 2".
315
316       •   -3 if $x is not strictly increasing.
317
318       chim does not process bad values.  It will set the bad-value flag of
319       all output ndarrays if the flag is set for any of the input ndarrays.
320
321   chic
322         Signature: (longlong ic(two=2);vc(two=2);mflag();x(n);f(n);[o]d(n);wk(nwk);longlong [o]ierr())
323
324       Calculate the derivatives of (x,f(x)) using cubic Hermite
325       interpolation.
326
327       Calculate the derivatives at the given points ("$x,$f", where $x is
328       strictly increasing).  Control over the boundary conditions is given by
329       the $ic and $vc ndarrays, and the value of $mflag determines the
330       treatment of points where monotoncity switches direction. A simpler,
331       more restricted, interface is available using "chim".
332
333       The first and second elements of $ic determine the boundary conditions
334       at the start and end of the data respectively.  If the value is 0, then
335       the default condition, as used by "chim", is adopted.  If greater than
336       zero, no adjustment for monotonicity is made, otherwise if less than
337       zero the derivative will be adjusted.  The allowed magnitudes for ic(0)
338       are:
339
340       •   1 if first derivative at x(0) is given in vc(0).
341
342       •   2 if second derivative at x(0) is given in vc(0).
343
344       •   3 to use the 3-point difference formula for d(0).  (Reverts to the
345           default b.c. if "n < 3")
346
347       •   4 to use the 4-point difference formula for d(0).  (Reverts to the
348           default b.c. if "n < 4")
349
350       •   5 to set d(0) so that the second derivative is continuous at x(1).
351           (Reverts to the default b.c. if "n < 4")
352
353       The values for ic(1) are the same as above, except that the first-
354       derivative value is stored in vc(1) for cases 1 and 2.  The values of
355       $vc need only be set if options 1 or 2 are chosen for $ic.
356
357       Set "$mflag = 0" if interpolant is required to be monotonic in each
358       interval, regardless of the data. This causes $d to be set to 0 at all
359       switch points. Set $mflag to be non-zero to use a formula based on the
360       3-point difference formula at switch points. If "$mflag > 0", then the
361       interpolant at swich points is forced to not deviate from the data by
362       more than "$mflag*dfloc", where "dfloc" is the maximum of the change of
363       $f on this interval and its two immediate neighbours.  If "$mflag < 0",
364       no such control is to be imposed.
365
366       The ndarray $wk is only needed for work space. However, I could not get
367       it to work as a temporary variable, so you must supply it; it is a 1D
368       ndarray with "2*n" elements.
369
370       Error status returned by $ierr:
371
372       •   0 if successful.
373
374       •   1 if "ic(0) < 0" and d(0) had to be adjusted for monotonicity.
375
376       •   2 if "ic(1) < 0" and d(n-1) had to be adjusted for monotonicity.
377
378       •   3 if both 1 and 2 are true.
379
380       •   -1 if "n < 2".
381
382       •   -3 if $x is not strictly increasing.
383
384       •   -4 if "abs(ic(0)) > 5".
385
386       •   -5 if "abs(ic(1)) > 5".
387
388       •   -6 if both -4 and -5  are true.
389
390       •   -7 if "nwk < 2*(n-1)".
391
392       chic does not process bad values.  It will set the bad-value flag of
393       all output ndarrays if the flag is set for any of the input ndarrays.
394
395   chsp
396         Signature: (longlong ic(two=2);vc(two=2);x(n);f(n);[o]d(n);wk(nwk);longlong [o]ierr())
397
398       Calculate the derivatives of (x,f(x)) using cubic spline interpolation.
399
400       Calculate the derivatives, using cubic spline interpolation, at the
401       given points ("$x,$f"), with the specified boundary conditions.
402       Control over the boundary conditions is given by the $ic and $vc
403       ndarrays.  The resulting values - "$x,$f,$d" - can be used in all the
404       functions which expect a cubic Hermite function.
405
406       The first and second elements of $ic determine the boundary conditions
407       at the start and end of the data respectively.  The allowed values for
408       ic(0) are:
409
410       •   0 to set d(0) so that the third derivative is continuous at x(1).
411
412       •   1 if first derivative at x(0) is given in "vc(0").
413
414       •   2 if second derivative at "x(0") is given in vc(0).
415
416       •   3 to use the 3-point difference formula for d(0).  (Reverts to the
417           default b.c. if "n < 3".)
418
419       •   4 to use the 4-point difference formula for d(0).  (Reverts to the
420           default b.c. if "n < 4".)
421
422       The values for ic(1) are the same as above, except that the first-
423       derivative value is stored in vc(1) for cases 1 and 2.  The values of
424       $vc need only be set if options 1 or 2 are chosen for $ic.
425
426       The ndarray $wk is only needed for work space. However, I could not get
427       it to work as a temporary variable, so you must supply it; it is a 1D
428       ndarray with "2*n" elements.
429
430       Error status returned by $ierr:
431
432       •   0 if successful.
433
434       •   -1  if "nelem($x) < 2".
435
436       •   -3  if $x is not strictly increasing.
437
438       •   -4  if "ic(0) < 0" or "ic(0) > 4".
439
440       •   -5  if "ic(1) < 0" or "ic(1) > 4".
441
442       •   -6  if both of the above are true.
443
444       •   -7  if "nwk < 2*n".
445
446       •   -8  in case of trouble solving the linear system for the interior
447           derivative values.
448
449       chsp does not process bad values.  It will set the bad-value flag of
450       all output ndarrays if the flag is set for any of the input ndarrays.
451
452   chfd
453         Signature: (x(n);f(n);d(n);longlong check();xe(ne);[o]fe(ne);[o]de(ne);longlong [o]ierr())
454
455       Interpolate function and derivative values.
456
457       Given a piecewise cubic Hermite function - such as from "chim" -
458       evaluate the function ($fe) and derivative ($de) at a set of points
459       ($xe).  If function values alone are required, use "chfe".  Set "check"
460       to 0 to skip checks on the input data.
461
462       Error status returned by $ierr:
463
464       •   0 if successful.
465
466       •   >0 if extrapolation was performed at "ierr" points (data still
467           valid).
468
469       •   -1 if "nelem($x) < 2"
470
471       •   -3 if $x is not strictly increasing.
472
473       •   -4 if "nelem($xe) < 1".
474
475       •   -5 if an error has occurred in a lower-level routine, which should
476           never happen.
477
478       chfd does not process bad values.  It will set the bad-value flag of
479       all output ndarrays if the flag is set for any of the input ndarrays.
480
481   chfe
482         Signature: (x(n);f(n);d(n);longlong check();xe(ne);[o]fe(ne);longlong [o]ierr())
483
484       Interpolate function values.
485
486       Given a piecewise cubic Hermite function - such as from "chim" -
487       evaluate the function ($fe) at a set of points ($xe).  If derivative
488       values are also required, use "chfd".  Set "check" to 0 to skip checks
489       on the input data.
490
491       Error status returned by $ierr:
492
493       •   0 if successful.
494
495       •   >0 if extrapolation was performed at "ierr" points (data still
496           valid).
497
498       •   -1 if "nelem($x) < 2"
499
500       •   -3 if $x is not strictly increasing.
501
502       •   -4 if "nelem($xe) < 1".
503
504       chfe does not process bad values.  It will set the bad-value flag of
505       all output ndarrays if the flag is set for any of the input ndarrays.
506
507   chia
508         Signature: (x(n);f(n);d(n);longlong check();la();lb();[o]ans();longlong [o]ierr())
509
510       Integrate (x,f(x)) over arbitrary limits.
511
512       Evaluate the definite integral of a piecewise cubic Hermite function
513       over an arbitrary interval, given by "[$la,$lb]". $d should contain the
514       derivative values, computed by "chim".  See "chid" if the integration
515       limits are data points.  Set "check" to 0 to skip checks on the input
516       data.
517
518       The values of $la and $lb do not have to lie within $x, although the
519       resulting integral value will be highly suspect if they are not.
520
521       Error status returned by $ierr:
522
523       •   0 if successful.
524
525       •   1 if $la lies outside $x.
526
527       •   2 if $lb lies outside $x.
528
529       •   3 if both 1 and 2 are true.
530
531       •   -1 if "nelem($x) < 2"
532
533       •   -3 if $x is not strictly increasing.
534
535       •   -4 if an error has occurred in a lower-level routine, which should
536           never happen.
537
538       chia does not process bad values.  It will set the bad-value flag of
539       all output ndarrays if the flag is set for any of the input ndarrays.
540
541   chid
542         Signature: (x(n);f(n);d(n);longlong check();longlong ia();longlong ib();[o]ans();longlong [o]ierr())
543
544       Integrate (x,f(x)) between data points.
545
546       Evaluate the definite integral of a a piecewise cubic Hermite function
547       between x($ia) and x($ib).
548
549       See "chia" for integration between arbitrary limits.
550
551       Although using a fortran routine, the values of $ia and $ib are zero
552       offset.  $d should contain the derivative values, computed by "chim".
553       Set "check" to 0 to skip checks on the input data.
554
555       Error status returned by $ierr:
556
557       •   0 if successful.
558
559       •   -1 if "nelem($x) < 2".
560
561       •   -3 if $x is not strictly increasing.
562
563       •   -4 if $ia or $ib is out of range.
564
565       chid does not process bad values.  It will set the bad-value flag of
566       all output ndarrays if the flag is set for any of the input ndarrays.
567
568   chcm
569         Signature: (x(n);f(n);d(n);longlong check();longlong [o]ismon(n);longlong [o]ierr())
570
571       Check the given piecewise cubic Hermite function for monotonicity.
572
573       The outout ndarray $ismon indicates over which intervals the function
574       is monotonic.  Set "check" to 0 to skip checks on the input data.
575
576       For the data interval "[x(i),x(i+1)]", the values of ismon(i) can be:
577
578       •   -3 if function is probably decreasing
579
580       •   -1 if function is strictly decreasing
581
582       •   0  if function is constant
583
584       •   1  if function is strictly increasing
585
586       •   2  if function is non-monotonic
587
588       •   3  if function is probably increasing
589
590       If "abs(ismon(i)) == 3", the derivative values are near the boundary of
591       the monotonicity region. A small increase produces non-monotonicity,
592       whereas a decrease produces strict monotonicity.
593
594       The above applies to "i = 0 .. nelem($x)-1". The last element of $ismon
595       indicates whether the entire function is monotonic over $x.
596
597       Error status returned by $ierr:
598
599       •   0 if successful.
600
601       •   -1 if "n < 2".
602
603       •   -3 if $x is not strictly increasing.
604
605       chcm does not process bad values.  It will set the bad-value flag of
606       all output ndarrays if the flag is set for any of the input ndarrays.
607
608   chbs
609         Signature: (x(n);f(n);d(n);longlong knotyp();longlong nknots();t(tsize);[o]bcoef(bsize);longlong [o]ndim();longlong [o]kord();longlong [o]ierr())
610
611       Piecewise Cubic Hermite function to B-Spline converter.
612
613       The resulting B-spline representation of the data (i.e. "nknots", "t",
614       "bcoeff", "ndim", and "kord") can be evaluated by "bvalu" (which is
615       currently not available).
616
617       Array sizes: "tsize = 2*n + 4", "bsize = 2*n", and "ndim = 2*n".
618
619       "knotyp" is a flag which controls the knot sequence.  The knot sequence
620       "t" is normally computed from $x by putting a double knot at each "x"
621       and setting the end knot pairs according to the value of "knotyp"
622       (where "m = ndim = 2*n"):
623
624       •   0 -   Quadruple knots at the first and last points.
625
626       •   1 -   Replicate lengths of extreme subintervals: "t( 0 ) = t( 1 ) =
627           x(0) - (x(1)-x(0))" and "t(m+3) = t(m+2) = x(n-1) +
628           (x(n-1)-x(n-2))"
629
630       •   2 -   Periodic placement of boundary knots: "t( 0 ) = t( 1 ) = x(0)
631           - (x(n-1)-x(n-2))" and "t(m+3) = t(m+2) = x(n) + (x(1)-x(0))"
632
633       •   <0 - Assume the "nknots" and "t" were set in a previous call.
634
635       "nknots" is the number of knots and may be changed by the routine.  If
636       "knotyp >= 0", "nknots" will be set to "ndim+4", otherwise it is an
637       input variable, and an error will occur if its value is not equal to
638       "ndim+4".
639
640       "t" is the array of "2*n+4" knots for the B-representation and may be
641       changed by the routine.  If "knotyp >= 0", "t" will be changed so that
642       the interior double knots are equal to the x-values and the boundary
643       knots set as indicated above, otherwise it is assumed that "t" was set
644       by a previous call (no check is made to verify that the data forms a
645       legitimate knot sequence).
646
647       Error status returned by $ierr:
648
649       •   0 if successful.
650
651       •   -4 if "knotyp > 2".
652
653       •   -5 if "knotyp < 0" and "nknots != 2*n + 4".
654
655       chbs does not process bad values.  It will set the bad-value flag of
656       all output ndarrays if the flag is set for any of the input ndarrays.
657
658   polfit
659         Signature: (x(n); y(n); w(n); longlong maxdeg(); longlong [o]ndeg(); [o]eps(); [o]r(n); longlong [o]ierr(); [o]a(foo); [o]coeffs(bar);[t]xtmp(n);[t]ytmp(n);[t]wtmp(n);[t]rtmp(n))
660
661       Fit discrete data in a least squares sense by polynomials
662                 in one variable. x(), y() and w() must be of the same type.
663              This version handles bad values appropriately
664
665       polfit processes bad values.  It will set the bad-value flag of all
666       output ndarrays if the flag is set for any of the input ndarrays.
667

AUTHOR

669       Copyright (C) 1997 Tuomas J. Lukka.  Copyright (C) 2000 Tim Jenness,
670       Doug Burke.  All rights reserved. There is no warranty. You are allowed
671       to redistribute this software / documentation under certain conditions.
672       For details, see the file COPYING in the PDL distribution. If this file
673       is separated from the PDL distribution, the copyright notice should be
674       included in the file.
675
676
677
678perl v5.38.0                      2023-07-21                         Slatec(3)
Impressum