From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001 From: Joseph Hunkeler Date: Wed, 8 Jul 2015 20:46:52 -0400 Subject: Initial commit --- sys/vops/lz/aravr.x | 44 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 44 insertions(+) create mode 100644 sys/vops/lz/aravr.x (limited to 'sys/vops/lz/aravr.x') diff --git a/sys/vops/lz/aravr.x b/sys/vops/lz/aravr.x new file mode 100644 index 00000000..c3f0fb8f --- /dev/null +++ b/sys/vops/lz/aravr.x @@ -0,0 +1,44 @@ +# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc. + +include + +# ARAV -- Compute the mean and standard deviation of a sample array by +# iteratively rejecting points further than KSIG from the mean. If the +# value of KSIG is given as 0.0, a cutoff value will be automatically +# calculated from the standard deviation and number of points in the sample. +# The number of pixels remaining in the sample upon termination is returned +# as the function value. + +int procedure aravr (a, npix, mean, sigma, ksig) + +real a[ARB] # input data array +real mean, sigma, ksig, deviation, lcut, hcut, lgpx +int npix, ngpix, old_ngpix, awvgr() + +begin + lcut = -MAX_REAL # no rejection to start + hcut = MAX_REAL + ngpix = MAX_INT + + # Iteratively compute mean, sigma and reject outliers until no + # more pixels are rejected, or until there are no more pixels. + + repeat { + old_ngpix = ngpix + ngpix = awvgr (a, npix, mean, sigma, lcut, hcut) + if (ngpix <= 1 || sigma <= EPSILONR) + break + + if (ksig == 0.0) { # Chauvenet's relation + lgpx = log10 (real(ngpix)) + deviation = (lgpx * (-0.1042 * lgpx + 1.1695) + .8895) * sigma + } else + deviation = sigma * abs(ksig) + + lcut = mean - deviation # compute window + hcut = mean + deviation + + } until (ngpix >= old_ngpix) + + return (ngpix) +end -- cgit