1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
|
# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
include <math/curfit.h>
include "names.h"
# IC_REJECT -- Reject points with large residuals from the fit.
#
# The sigma of the fit residuals is calculated. The rejection thresholds
# are set at low_reject*sigma and high_reject*sigma. Points outside the
# rejection threshold are rejected from the fit and flagged in the rejpts
# array. Finally, the remaining points are refit.
procedure ic_rejectd (cv, x, y, w, rejpts, npts, low_reject, high_reject,
niterate, grow, nreject)
pointer cv # Curve descriptor
double x[npts] # Input ordinates
double y[npts] # Input data values
double w[npts] # Weights
int rejpts[npts] # Points rejected
int npts # Number of input points
real low_reject, high_reject # Rejection threshold
int niterate # Number of rejection iterations
real grow # Rejection radius
int nreject # Number of points rejected
int i, ierr, nit, newreject
errchk ic_deviantd
begin
# Initialize rejection.
nreject = 0
call amovki (NO, rejpts, npts)
if (niterate <= 0)
return
# Find deviant points. If an error occurs reduce the number of
# iterations and start again.
iferr {
nit = 0
do i = 1, niterate {
call ic_deviantd (cv, x, y, w, rejpts, npts, low_reject,
high_reject, grow, YES, nreject, newreject)
nit = nit + 1
if (newreject == 0)
break
}
} then {
call dcvfit (cv, x, y, w, npts, WTS_USER, ierr)
nreject = 0
call amovki (NO, rejpts, npts)
do i = 1, nit
call ic_deviantd (cv, x, y, w, rejpts, npts, low_reject,
high_reject, grow, YES, nreject, newreject)
}
end
|
