From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001
From: Joseph Hunkeler <jhunkeler@gmail.com>
Date: Wed, 8 Jul 2015 20:46:52 -0400
Subject: Initial commit

---
 pkg/xtools/inlfit/infit.gx | 99 ++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 99 insertions(+)
 create mode 100644 pkg/xtools/inlfit/infit.gx

(limited to 'pkg/xtools/inlfit/infit.gx')

diff --git a/pkg/xtools/inlfit/infit.gx b/pkg/xtools/inlfit/infit.gx
new file mode 100644
index 00000000..069bf584
--- /dev/null
+++ b/pkg/xtools/inlfit/infit.gx
@@ -0,0 +1,99 @@
+include	<math/nlfit.h>
+include	<pkg/inlfit.h>
+
+# IN_FIT -- Fit a function using non-linear least squares. The function
+# can have an arbitrary number of independent variables. This is the main
+# entry point for the non-interactive part of the INLFIT package.
+
+procedure in_fit$t (in, nl, x, y, wts, npts, nvars, wtflag, stat)
+
+pointer	in				# INLFIT pointer
+pointer	nl				# NLFIT pointer
+PIXEL	x[ARB]				# Ordinates (npts * nvars)
+PIXEL	y[npts]				# Data to be fit
+PIXEL	wts[npts]			# Weights
+int	npts				# Number of points
+int	nvars				# Number of variables
+int	wtflag				# Type of weighting
+int	stat				# Error code (output)
+
+int	i, ndeleted
+pointer	sp, wts1, str
+int	in_geti()
+PIXEL	in_get$t
+
+begin
+
+#	# Debug.
+#	call eprintf ("in_fit: in=%d, nl=%d, npts=%d, nvars=%d\n")
+#	    call pargi (in)
+#	    call pargi (nl)
+#	    call pargi (npts)
+#	    call pargi (nvars)
+
+	# Allocate string, and rejection weight space. The latter are
+	# are used to mark rejected points with a zero weight before
+	# calling NLFIT.
+
+	call smark (sp)
+	call salloc (str, SZ_LINE, TY_CHAR)
+	call salloc (wts1, npts, TY_PIXEL)
+	call amov$t (wts, Mem$t[wts1], npts)
+
+	# Initialize rejected point list, and the buffer containing
+	# the minimum and maximum variable values.
+	call in_bfinit$t (in, npts, nvars)
+
+	# Set independent variable limits.
+	call in_limit$t (in, x, npts, nvars)
+
+	# Reinitialize.
+	call in_nlinit$t (in, nl)
+
+	# Check number of data points. If no points are present
+	# set the error flag to the appropiate value, and return.
+	if (npts == 0) {
+	    stat = NO_DEG_FREEDOM
+	    call in_puti (in, INLFITERROR, NO_DEG_FREEDOM)
+	    call sfree (sp)
+	    return
+	}
+
+	# Check the number of deleted points.
+	ndeleted = 0
+	do i = 1, npts {
+	    if (wts[i] <= PIXEL(0.0))
+		ndeleted = ndeleted + 1
+	}
+	if ((npts - ndeleted) < in_geti (in, INLNFPARAMS)) {
+	    stat = NO_DEG_FREEDOM
+	    call in_puti (in, INLFITERROR, NO_DEG_FREEDOM)
+	    call sfree (sp)
+	    return
+	}
+
+	# Call NLFIT.
+	call nlfit$t (nl, x, y, wts, npts, nvars, wtflag, stat)
+
+	# Update fit status into the INLFIT structure.
+	call in_puti (in, INLFITERROR, stat)
+
+	# Do pixel rejection and refit, if at least one of the rejection
+	# limits is positive. Otherwise clear number of rejected points.
+
+	if (in_get$t (in, INLLOW)  > PIXEL (0.0) ||
+	    in_get$t (in, INLHIGH) > PIXEL (0.0)) {
+	    call in_reject$t (in, nl, x, y, Mem$t[wts1], npts, nvars, wtflag)
+	    if (in_geti (in, INLNREJPTS) > 0) {
+		do i = 1, npts {
+		    if (Mem$t[wts1+i-1] > PIXEL(0.0))
+			wts[i] = Mem$t[wts1+i-1]
+		}
+	    }
+	    stat = in_geti (in, INLFITERROR)
+	} else
+	    call in_puti (in, INLNREJPTS, 0)
+
+	# Free memory.
+	call sfree (sp)
+end
-- 
cgit