From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001 From: Joseph Hunkeler Date: Wed, 8 Jul 2015 20:46:52 -0400 Subject: Initial commit --- pkg/xtools/icfit/icreject.gx | 57 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 57 insertions(+) create mode 100644 pkg/xtools/icfit/icreject.gx (limited to 'pkg/xtools/icfit/icreject.gx') diff --git a/pkg/xtools/icfit/icreject.gx b/pkg/xtools/icfit/icreject.gx new file mode 100644 index 00000000..79965384 --- /dev/null +++ b/pkg/xtools/icfit/icreject.gx @@ -0,0 +1,57 @@ +# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc. + +include +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_reject$t (cv, x, y, w, rejpts, npts, low_reject, high_reject, + niterate, grow, nreject) + +pointer cv # Curve descriptor +PIXEL x[npts] # Input ordinates +PIXEL y[npts] # Input data values +PIXEL 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_deviant$t + +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_deviant$t (cv, x, y, w, rejpts, npts, low_reject, + high_reject, grow, YES, nreject, newreject) + nit = nit + 1 + if (newreject == 0) + break + } + } then { + call $tcvfit (cv, x, y, w, npts, WTS_USER, ierr) + nreject = 0 + call amovki (NO, rejpts, npts) + do i = 1, nit + call ic_deviant$t (cv, x, y, w, rejpts, npts, low_reject, + high_reject, grow, YES, nreject, newreject) + } +end -- cgit