remove libfit

This commit is contained in:
Matthew Flatt 2011-12-21 16:05:53 -07:00
parent 5a7f1b15ed
commit d020c75202
22 changed files with 6 additions and 5350 deletions

View File

@ -661,8 +661,6 @@ mz-extras :+= (- (package: "swindle")
;; -------------------- plot
plt-extras :+= (package: "plot")
(src: "fit")
(lib: "libfit*")
;; -------------------- mzcom
plt-extras :+= (- (package: "mzcom" #:src? #t)

View File

@ -18,8 +18,6 @@
plot-foreground plot-background
plot3d-angle plot3d-altitude))
"deprecated/renderers.rkt"
;; Curve fitting
"deprecated/fit.rkt"
;; Miscellaneous
"deprecated/math.rkt")
@ -31,9 +29,6 @@
contour shade
surface)
(only-doc-out (all-defined-out))
;; Curve fitting
(rename-out [fit-int fit])
(struct-out fit-result)
;; Miscellaneous
make-vec derivative gradient)

View File

@ -1,56 +0,0 @@
(module fit-low-level racket/base
(require mzlib/foreign mzlib/runtime-path
(for-syntax racket/base))
(unsafe!)
(define-runtime-path libfit-path '(so "libfit"))
(define libfit (ffi-lib libfit-path))
(define do-fit-int
(get-ffi-obj "do_fit" libfit
(_fun (func : (_fun _int _pointer -> _double))
(val-num : _int = (length x-values))
(x-values : (_list i _double*))
(y-values : (_list i _double*))
(z-values : (_list i _double*))
(errors : (_list i _double*))
(param-num : _int = (length params))
(params : (_list i _double*))
-> (_list o _double* param-num))))
(define (do-fit callback x-vals y-vals z-vals errors params)
(do-fit-int (lambda (argc argv)
(let ([args (cblock->list argv _double argc)])
(apply callback args)))
x-vals y-vals z-vals errors params))
(define get-asym-error
(get-ffi-obj "get_asym_error" libfit
(_fun (len : _?) ; len is only used for list conversion
-> (_list o _double* len))))
(define get-asym-error-percent
(get-ffi-obj "get_asym_error_percent" libfit
(_fun (len : _?) ; len is only used for list conversion
-> (_list o _double* len))))
(define get-rms
(get-ffi-obj "get_rms" libfit
(_fun -> _double*)))
(define get-varience
(get-ffi-obj "get_varience" libfit
(_fun -> _double*)))
(define (fit-internal f-of-x-y x-vals y-vals z-vals err-vals params)
(let* ([len (length params)]
[fit-result (do-fit f-of-x-y x-vals y-vals z-vals err-vals params)]
[asym-error (get-asym-error len)]
[asym-error-percent (get-asym-error-percent len)]
[rms (get-rms)]
[varience (get-varience)])
(list fit-result asym-error asym-error-percent rms varience)))
(provide fit-internal))

View File

@ -1,57 +0,0 @@
(module fit mzscheme
(require unstable/lazy-require
"math.rkt")
;; Require lazily so the rest of 'plot' still works without libfit:
(lazy-require ["fit-low-level.rkt" (fit-internal)])
; a structure contain a the results of a curve-fit
(define-struct fit-result (
rms
variance
names
final-params
std-error
std-error-percent
function
) (make-inspector))
; fit-int : (number* -> number) (list-of (symbol number)) (list-of (vector number [number] number number)) -> fit-result
(define (fit-int function guesses data)
(let* ((independent-vars (- (procedure-arity function) (length guesses)))
(f-of-x-y (cond
[(= 1 independent-vars)
(lambda (x y . rest)
(apply function x rest))]
[(= 2 independent-vars)
function]
[else
(error "Function provided is either not of one or two independent variables or the number of
guesses given is incorrect")]))
(x-vals (map vector-x data))
(y-vals (if (= 1 independent-vars)
x-vals
(map vector-y data)))
(z-vals (if (= 1 independent-vars)
(map vector-y data)
(map vector-z data)))
(err-vals (if (= 1 independent-vars)
(map vector-z data)
(map (lambda (vec) (vector-ref vec 4)) data)))
(result (fit-internal f-of-x-y x-vals y-vals z-vals err-vals (map cadr guesses))))
(if (null? result)
null
(begin
;(display result)
(make-fit-result
(list-ref result 3)
(list-ref result 4)
(map car guesses)
(car result)
(cadr result)
(caddr result)
(lambda args (apply function(append args (car result)))))))))
(provide fit-int
(struct fit-result (rms variance names final-params
std-error std-error-percent function))))

View File

@ -127,116 +127,6 @@ Returns @racket[#t] if @racket[v] is one of the following symbols,
@; ----------------------------------------
@section[#:tag "curve-fit"]{Curve Fitting}
@define[fit-warning]{
@para{
@bold{Do not use the @(racket fit) function. It is going to be removed in Racket 5.2.1.}
It relies on old C code that nobody understands or is willing to maintain, and that is also slightly crashy.
}}
@fit-warning
Quite independent of plotting, and for reasons lost in the sands of time,
the @racketmodname[plot] library provides a non-linear, least-squares
fit algorithm to fit parameterized functions to given data.
The code that implements the algorithm is public
domain, and is used by the @tt{gnuplot} package.
To fit a particular function to a curve:
@itemize[
@item{Set up the independent and dependent variable data. The first
item in each vector is the independent variable, the second is the
result. The last item is the weight of the error; we can leave it
as @racket[1] since all the items weigh the same.
@racketblock[
(define data '(#(0 3 1)
#(1 5 1)
#(2 7 1)
#(3 9 1)
#(4 11 1)))
]
}
@item{Set up the function to be fitted using fit. This particular
function looks like a line. The independent variables must come
before the parameters.
@racketblock[
(define fit-fun
(lambda (x m b) (+ b (* m x))))
]
}
@item{If possible, come up with some guesses for the values of the
parameters. The guesses can be left as one, but each parameter must
be named.}
@item{Do the fit.
@racketblock[
(define fitted
(fit fit-fun
'((m 1) (b 1))
data))
]
}
@item{View the resulting parameters; for example,
@racketblock[
(fit-result-final-params fitted)
]
will produce @racketresultfont{(2.0 3.0)}.
}
@item{For some visual feedback of the fit result, plot the function
with the new parameters. For convenience, the structure that is
returned by the fit command has already the function.
@racketblock[
(plot (mix (points data)
(line (fit-result-function fitted)))
#:y-max 15)
]}]
A more realistic example can be found in
@filepath{compat/tests/fit-demo-2.rkt} in the @filepath{plot} collection.
@defproc[(fit [f (real? ... . -> . real?)]
[guess-list (list/c (list symbol? real?))]
[data (or/c (list-of (vector/c real? real? real?))
(list-of (vector/c real? real? real? real?)))])
fit-result?]{
@fit-warning
Attempts to fit a @defterm{fittable function} to the data that is
given. The @racket[guess-list] should be a set of arguments and
values. The more accurate your initial guesses are, the more likely
the fit is to succeed; if there are no good values for the guesses,
leave them as @racket[1].}
@defstruct[fit-result ([rms real?]
[variance real?]
[names (listof symbol?)]
[final-params (listof real?)]
[std-error (listof real?)]
[std-error-percent (listof real?)]
[function (real? ... . -> . real?)])]{
The @racket[params] field contains an associative list of the
parameters specified in @racket[fit] and their values. Note that the
values may not be correct if the fit failed to converge. For a visual
test, use the @racket[function] field to get the function with the
parameters in place and plot it along with the original data.}
@; ----------------------------------------
@section{Miscellaneous Functions}
@defproc[(derivative [f (real? . -> . real?)] [h real? .000001])

View File

@ -18,6 +18,8 @@ The update from PLoT version 5.1.3 to 5.2 introduces a few incompatibilities:
The argument change in @(racket plot3d) is similar.
This should not affect most code because PLoT encourages regarding these data types as black boxes.}
@item{The @(racket plot-extend) module no longer exists.}
@item{The @racket[fit] function and @racket[fit-result] functions have been removed.}
]
This section of the PLoT manual will help you port code written for PLoT 5.1.3 and earlier to the most recent PLoT.

File diff suppressed because it is too large Load Diff

File diff suppressed because it is too large Load Diff

BIN
collects/plot/tests/sqr.png Normal file

Binary file not shown.

After

Width:  |  Height:  |  Size: 12 KiB

Binary file not shown.

After

Width:  |  Height:  |  Size: 69 KiB

View File

@ -39,28 +39,20 @@ all:
3m:
cd racket; $(MAKE) 3m
$(MAKE) @MAKE_GRACKET@-3m
$(MAKE) @MAKE_FIT@-3m
gracket-3m:
cd gracket; $(MAKE) 3m
fit-3m:
cd fit; $(MAKE) 3m
no-3m:
$(NOOP)
cgc:
cd racket; $(MAKE) cgc
$(MAKE) @MAKE_GRACKET@-cgc
$(MAKE) @MAKE_FIT@-cgc
gracket-cgc:
cd gracket; $(MAKE) cgc
fit-cgc:
cd fit; $(MAKE) cgc
no-cgc:
$(NOOP)
@ -73,8 +65,7 @@ both:
SETUP_ARGS = -X "$(DESTDIR)$(collectsdir)" -N "raco setup" -l- setup $(PLT_SETUP_OPTIONS) $(PLT_ISO) @INSTALL_SETUP_FLAGS@
# Pass compile and link flags to `make install' for use by any
# collection-setup actions (currently in "plot") that compile
# and link C code:
# collection-setup actions that compile and link C code:
CFLAGS = @CFLAGS@ @COMPFLAGS@ @PREFLAGS@
LDFLAGS = @LDFLAGS@
WITH_ENV_VARS = env CFLAGS="$(CFLAGS)" LDFLAGS="$(LDFLAGS)"
@ -92,7 +83,6 @@ install-common-first:
mkdir -p $(ALLDIRINFO)
install-common-middle:
$(MAKE) install-@MAKE_FIT@
$(MAKE) @MAKE_COPYTREE@-run
$(MAKE) install-@MAKE_GRACKET@-post-collects
$(MAKE) lib-finish
@ -100,9 +90,6 @@ install-common-middle:
install-common-last:
$(MAKE) fix-paths
install-fit:
cd fit; $(MAKE) install
install-no:
$(NOOP)

22
src/configure vendored
View File

@ -721,7 +721,6 @@ ICP
MRLIBINSTALL
LIBFINISH
MAKE_GRACKET
MAKE_FIT
MAKE_COPYTREE
MAKE_FINISH
WXPRECOMP
@ -1333,7 +1332,6 @@ Optional Features:
--enable-foreign support foreign calls (enabled by default)
--enable-places support places (3m only; usually enabled by default)
--enable-futures support futures (usually enabled by default)
--enable-plot support plot libraries (enabled by default)
--enable-float support single-precision floats (enabled by default)
--enable-floatinstead use single-precision by default
--enable-racket=<path> use <path> as Racket executable to build Racket
@ -1986,13 +1984,6 @@ if test "${enable_futures+set}" = set; then
enableval=$enable_futures;
fi
# Check whether --enable-plot was given.
if test "${enable_plot+set}" = set; then
enableval=$enable_plot;
else
enable_plot=yes
fi
# Check whether --enable-float was given.
if test "${enable_float+set}" = set; then
enableval=$enable_float;
@ -2314,7 +2305,6 @@ show_explicitly_enabled "${enable_xonx}" "Unix style"
show_explicitly_enabled "${enable_shared}" "Shared libraries"
show_explicitly_disabled "${enable_gracket}" GRacket
show_explicitly_disabled "${enable_plot}" Plot fit library
if test "$LIBTOOLPROG" != "" ; then
echo "=== Libtool program: $LIBTOOLPROG"
@ -8947,7 +8937,6 @@ LIBS="$LIBS $EXTRALIBS"
mk_needed_dir()
@ -9015,14 +9004,6 @@ fi
makefiles="$makefiles foreign/Makefile"
ac_configure_args="$ac_configure_args$SUB_CONFIGURE_EXTRAS"
if test -d "${srcdir}/fit" && test "${enable_plot}" = "yes" ; then
makefiles="$makefiles
fit/Makefile"
MAKE_FIT=fit
else
MAKE_FIT=no
fi
if test "${enable_gracket}" = "yes" ; then
makefiles="$makefiles
gracket/Makefile
@ -9794,7 +9775,6 @@ ICP!$ICP$ac_delim
MRLIBINSTALL!$MRLIBINSTALL$ac_delim
LIBFINISH!$LIBFINISH$ac_delim
MAKE_GRACKET!$MAKE_GRACKET$ac_delim
MAKE_FIT!$MAKE_FIT$ac_delim
MAKE_COPYTREE!$MAKE_COPYTREE$ac_delim
MAKE_FINISH!$MAKE_FINISH$ac_delim
WXPRECOMP!$WXPRECOMP$ac_delim
@ -9832,7 +9812,7 @@ LIBOBJS!$LIBOBJS$ac_delim
LTLIBOBJS!$LTLIBOBJS$ac_delim
_ACEOF
if test `sed -n "s/.*$ac_delim\$/X/p" conf$$subs.sed | grep -c X` = 46; then
if test `sed -n "s/.*$ac_delim\$/X/p" conf$$subs.sed | grep -c X` = 45; then
break
elif $ac_last_try; then
{ { echo "$as_me:$LINENO: error: could not make $CONFIG_STATUS" >&5

View File

@ -1,37 +0,0 @@
srcdir = @srcdir@
prefix = @prefix@
exec_prefix = @exec_prefix@
bindir = @bindir@
libdir = @libdir@
libpltdir = @libpltdir@
collectsdir = @collectsdir@
builddir = @builddir@
ICP=@ICP@
CC = @CC@
# See ../Makefile about RUN_RACKET_<X>, which
# typically redirects to RUN_THIS_RACKET_<X>:
RUN_THIS_RACKET_CGC = ../racket/racket@CGC@
WITH_ENV = env CC="@PLAIN_CC@" CFLAGS="@CFLAGS@ @COMPFLAGS@ @PREFLAGS@" LDFLAGS="@LDFLAGS@"
FIT_SRCS = $(srcdir)/fit.c $(srcdir)/matrix.c
XCOLLECTS = # -X ../racket/gc2/xform-collects
fit-lib: libfit@SO_SUFFIX@
3m:
$(MAKE) fit-lib
cgc:
$(MAKE) fit-lib
libfit@SO_SUFFIX@:
$(WITH_ENV) @RUN_RACKET_CGC@ -c $(srcdir)/build.rkt "libfit" $(FIT_SRCS)
install:
cd ..; $(ICP) fit/libfit@SO_SUFFIX@ "$(DESTDIR)$(libpltdir)/libfit@SO_SUFFIX@"

View File

@ -1 +0,0 @@
..\..\racket build.rkt libfit fit.c matrix.c

View File

@ -1,32 +0,0 @@
(module build racket/base
(require racket/path
racket/file
dynext/file
dynext/link
dynext/compile)
(define-values (libname c-files)
(let ([l (vector->list (current-command-line-arguments))])
(values (car l)
(cdr l))))
(define sys-subpath (system-library-subpath #f))
(define so-name (append-extension-suffix libname))
(parameterize (;; we compile a simple .so, not an extension
[current-standard-link-libraries '()])
(when (or (not (file-exists? so-name))
(let ([so-time (file-or-directory-modify-seconds so-name)])
(for/or ([f c-files])
((file-or-directory-modify-seconds f) . > . so-time))))
(let ([o-files
(for/list ([c-file c-files])
(let ([o-file (append-object-suffix (path-replace-suffix (file-name-from-path c-file) #""))])
;; first #f means not quiet (here and in link-extension)
(compile-extension #f c-file o-file null)
o-file))])
(let* ([flags (if (string=? "i386-cygwin" (path->string sys-subpath))
;; DLL needs every dependence explicit:
'("-lc" "-lm" "-lcygwin" "-lkernel32")
null)])
(link-extension #f (append o-files flags) so-name))))))

View File

@ -1,6 +0,0 @@
#if (defined(__WIN32__) || defined(WIN32) || defined(_WIN32))
# define MZ_DLLEXPORT __declspec(dllexport)
#else
# define MZ_DLLEXPORT
#endif

View File

@ -1,752 +0,0 @@
/* NOTICE: Change of Copyright Status
*
* The author of this module, Carsten Grammes, has expressed in
* personal email that he has no more interest in this code, and
* doesn't claim any copyright. He has agreed to put this module
* into the public domain.
*
* Lars Hecking 15-02-1999
*/
/*
* Nonlinear least squares fit according to the
* Marquardt-Levenberg-algorithm
*
* added as Patch to Gnuplot (v3.2 and higher)
* by Carsten Grammes
* Experimental Physics, University of Saarbruecken, Germany
*
* Internet address: cagr@rz.uni-sb.de
*
* Copyright of this module: 1993, 1998 Carsten Grammes
*
* Permission to use, copy, and distribute this software and its
* documentation for any purpose with or without fee is hereby granted,
* provided that the above copyright notice appear in all copies and
* that both that copyright notice and this permission notice appear
* in supporting documentation.
*
* This software is provided "as is" without express or implied warranty.
*
* 930726: Recoding of the Unix-like raw console I/O routines by:
* Michele Marziani (marziani@ferrara.infn.it)
* drd: start unitialised variables at 1 rather than NEARLY_ZERO
* (fit is more likely to converge if started from 1 than 1e-30 ?)
*
* HBB (Broeker@physik.rwth-aachen.de) : fit didn't calculate the errors
* in the 'physically correct' (:-) way, if a third data column containing
* the errors (or 'uncertainties') of the input data was given. I think
* I've fixed that, but I'm not sure I really understood the M-L-algo well
* enough to get it right. I deduced my change from the final steps of the
* equivalent algorithm for the linear case, which is much easier to
* understand. (I also made some minor, mostly cosmetic changes)
*
* HBB (again): added error checking for negative covar[i][i] values and
* for too many parameters being specified.
*
* drd: allow 3d fitting. Data value is now called fit_z internally,
* ie a 2d fit is z vs x, and a 3d fit is z vs x and y.
*
* Lars Hecking : review update command, for VMS in particular, where
* it is not necessary to rename the old file.
*
* HBB, 971023: lifted fixed limit on number of datapoints, and number
* of parameters.
*/
#define FIT_MAIN
#define NULL 0
//#include <scheme.h>
#include "matrix.h"
#include "fit.h"
#include <math.h>
#include <stdlib.h>
#include <string.h>
/* #define STANDARD stderr */
enum marq_res {
OK, ERROR, BETTER, WORSE
};
typedef enum marq_res marq_res_t;
#ifdef INFINITY
# undef INFINITY
#endif
#define INFINITY 1e30
#define NEARLY_ZERO 1e-30
/* Relative change for derivatives */
#define DELTA 0.001
#define MAX_DATA 2048
#define MAX_PARAMS 32
#define MAX_LAMBDA 1e20
#define MIN_LAMBDA 1e-20
#define LAMBDA_UP_FACTOR 10
#define LAMBDA_DOWN_FACTOR 10
#define TBOOLEAN int
# define PLUSMINUS "+/-"
/* HBB 971023: new, allow for dynamic adjustment of these: */
static UNUSED int max_data;
static UNUSED int max_params;
static double epsilon = 1e-5; /* convergence limit */
static int maxiter = 0; /* HBB 970304: maxiter patch */
static UNUSED char *FIXED = "# FIXED";
static UNUSED char *GNUFITLOG = "FIT_LOG";
static UNUSED char *FITLIMIT = "FIT_LIMIT";
static UNUSED char *FITSTARTLAMBDA = "FIT_START_LAMBDA";
static UNUSED char *FITLAMBDAFACTOR = "FIT_LAMBDA_FACTOR";
static UNUSED char *FITMAXITER = "FIT_MAXITER"; /* HBB 970304: maxiter patch */
static UNUSED char *FITSCRIPT = "FIT_SCRIPT";
static UNUSED char *DEFAULT_CMD = "replot"; /* if no fitscript spec. */
static int num_data, num_params;
static UNUSED int columns;
static double *fit_x;
static double *fit_y;
static double *fit_z ;
static double *err_data;
static double *a;
/* static fixstr * par_name; */
static double startup_lambda = 0;
static double lambda_down_factor = LAMBDA_DOWN_FACTOR;
static double lambda_up_factor = LAMBDA_UP_FACTOR;
static void * current_fun;
/*****************************************************************
internal vars to store results of fit
*****************************************************************/
double rms = 0;
double varience = 0;
double *asym_error;
double *asym_error_percent;
MZ_DLLEXPORT
double get_rms()
{return rms;}
MZ_DLLEXPORT
double get_varience()
{return varience;}
MZ_DLLEXPORT
double * get_asym_error()
{return asym_error;}
MZ_DLLEXPORT
double * get_asym_error_percent()
{return asym_error_percent;}
/*****************************************************************
internal Prototypes
*****************************************************************/
/*static void printmatrix __PROTO((double **C, int m, int n)); */
static UNUSED void print_matrix_and_vectors (double **C, double *d, double *r, int m, int n);
static marq_res_t marquardt (double a[], double **alpha, double *chisq,
double *lambda);
static TBOOLEAN analyze (double a[], double **alpha, double beta[],
double *chisq);
static void calculate (double *zfunc, double **dzda, double a[]);
static void call_scheme (double *par, double *data);
static TBOOLEAN regress (double a[]);
//static void show_fit (int i, double chisq, double last_chisq, double *a,
// double lambda, FILE * device);
/*****************************************************************
New utility routine: print a matrix (for debugging the alg.)
*****************************************************************/
static UNUSED void printmatrix(C, m, n)
double **C;
int m, n;
{
int i, j;
for (i = 0; i < m; i++) {
for (j = 0; j < n - 1; j++);
/* Dblf2("%.8g |", C[i][j]); */
/* Dblf2("%.8g\n", C[i][j]); */
}
/* Dblf("\n"); */
}
/**************************************************************************
Yet another debugging aid: print matrix, with diff. and residue vector
**************************************************************************/
static UNUSED void print_matrix_and_vectors(C, d, r, m, n)
double **C;
double *d, *r;
int m, n;
{
int i, j;
for (i = 0; i < m; i++) {
for (j = 0; j < n; j++);
/* Dblf2("%8g ", C[i][j]); */
/* Dblf3("| %8g | %8g\n", d[i], r[i]); */
}
/* Dblf("\n"); */
}
/*****************************************************************
Marquardt's nonlinear least squares fit
*****************************************************************/
static marq_res_t marquardt(a, C, chisq, lambda)
double a[];
double **C;
double *chisq;
double *lambda;
{
int i, j;
static double *da = 0, /* delta-step of the parameter */
*temp_a = 0, /* temptative new params set */
*d = 0, *tmp_d = 0, **tmp_C = 0, *residues = 0;
double tmp_chisq;
/* Initialization when lambda == -1 */
if (*lambda == -1) { /* Get first chi-square check */
TBOOLEAN analyze_ret;
temp_a = vec(num_params);
d = vec(num_data + num_params);
tmp_d = vec(num_data + num_params);
da = vec(num_params);
residues = vec(num_data + num_params);
tmp_C = matr(num_data + num_params, num_params);
analyze_ret = analyze(a, C, d, chisq);
/* Calculate a useful startup value for lambda, as given by Schwarz */
/* FIXME: this is doesn't turn out to be much better, really... */
if (startup_lambda != 0)
*lambda = startup_lambda;
else {
*lambda = 0;
for (i = 0; i < num_data; i++)
for (j = 0; j < num_params; j++)
*lambda += C[i][j] * C[i][j];
*lambda = sqrt(*lambda / num_data / num_params);
}
/* Fill in the lower square part of C (the diagonal is filled in on
each iteration, see below) */
for (i = 0; i < num_params; i++)
for (j = 0; j < i; j++)
C[num_data + i][j] = 0, C[num_data + j][i] = 0;
/* printmatrix(C, num_data+num_params, num_params); */
return analyze_ret ? OK : ERROR;
}
/* once converged, free dynamic allocated vars */
if (*lambda == -2) {
return OK;
}
/* Givens calculates in-place, so make working copies of C and d */
for (j = 0; j < num_data + num_params; j++)
memcpy(tmp_C[j], C[j], num_params * sizeof(double));
memcpy(tmp_d, d, num_data * sizeof(double));
/* fill in additional parts of tmp_C, tmp_d */
for (i = 0; i < num_params; i++) {
/* fill in low diag. of tmp_C ... */
tmp_C[num_data + i][i] = *lambda;
/* ... and low part of tmp_d */
tmp_d[num_data + i] = 0;
}
/* printmatrix(tmp_C, num_data+num_params, num_params); */
/* FIXME: residues[] isn't used at all. Why? Should it be used? */
Givens(tmp_C, tmp_d, da, residues, num_params + num_data, num_params, 1);
/*print_matrix_and_vectors (tmp_C, tmp_d, residues,
num_params+num_data, num_params); */
/* check if trial did ameliorate sum of squares */
for (j = 0; j < num_params; j++)
temp_a[j] = a[j] + da[j];
if (!analyze(temp_a, tmp_C, tmp_d, &tmp_chisq)) {
/* FIXME: will never be reached: always returns TRUE */
return ERROR;
}
if (tmp_chisq < *chisq) { /* Success, accept new solution */
if (*lambda > MIN_LAMBDA) {
/* (void) putc('/', stderr); */
*lambda /= lambda_down_factor;
}
*chisq = tmp_chisq;
for (j = 0; j < num_data; j++) {
memcpy(C[j], tmp_C[j], num_params * sizeof(double));
d[j] = tmp_d[j];
}
for (j = 0; j < num_params; j++)
a[j] = temp_a[j];
return BETTER;
} else { /* failure, increase lambda and return */
/* (void) putc('*', stderr); */
*lambda *= lambda_up_factor;
return WORSE;
}
}
/* FIXME: in the new code, this function doesn't really do enough to be
* useful. Maybe it ought to be deleted, i.e. integrated with
* calculate() ?
*/
/*****************************************************************
compute chi-square and numeric derivations
*****************************************************************/
static TBOOLEAN analyze(a, C, d, chisq)
double a[];
double **C;
double d[];
double *chisq;
{
/*
* used by marquardt to evaluate the linearized fitting matrix C
* and vector d, fills in only the top part of C and d
* I don't use a temporary array zfunc[] any more. Just use
* d[] instead.
*/
int i, j;
*chisq = 0;
calculate(d, C, a);
for (i = 0; i < num_data; i++) {
/* note: order reversed, as used by Schwarz */
d[i] = (d[i] - fit_z[i]) / err_data[i];
*chisq += d[i] * d[i];
for (j = 0; j < num_params; j++)
C[i][j] /= err_data[i];
}
/* FIXME: why return a value that is always TRUE ? */
return 1;
}
/* To use the more exact, but slower two-side formula, activate the
following line: */
#define TWO_SIDE_DIFFERENTIATION
/*****************************************************************
compute function values and partial derivatives of chi-square
*****************************************************************/
static void calculate(zfunc, dzda, a)
double *zfunc;
double **dzda;
double a[];
{
int k, p;
double tmp_a;
double *tmp_high, *tmp_pars;
#ifdef TWO_SIDE_DIFFERENTIATION
double *tmp_low;
#endif
tmp_high = vec(num_data); /* numeric derivations */
#ifdef TWO_SIDE_DIFFERENTIATION
tmp_low = vec(num_data);
#endif
tmp_pars = vec(num_params);
/* first function values */
call_scheme(a, zfunc);
/* then derivatives */
for (p = 0; p < num_params; p++)
tmp_pars[p] = a[p];
for (p = 0; p < num_params; p++) {
tmp_a = fabs(a[p]) < NEARLY_ZERO ? NEARLY_ZERO : a[p];
tmp_pars[p] = tmp_a * (1 + DELTA);
call_scheme(tmp_pars, tmp_high);
#ifdef TWO_SIDE_DIFFERENTIATION
tmp_pars[p] = tmp_a * (1 - DELTA);
call_scheme(tmp_pars, tmp_low);
#endif
for (k = 0; k < num_data; k++)
#ifdef TWO_SIDE_DIFFERENTIATION
dzda[k][p] = (tmp_high[k] - tmp_low[k]) / (2 * tmp_a * DELTA);
#else
dzda[k][p] = (tmp_high[k] - zfunc[k]) / (tmp_a * DELTA);
#endif
tmp_pars[p] = a[p];
}
}
/*****************************************************************
evaluate the scheme function
*****************************************************************/
static void call_scheme(par, data)
double *par;
double *data;
{
int rators = 2 + num_params;
double * rands =
(double *) malloc(rators * sizeof(double));
int i;
/* set up the constant params */
for(i = 0 ; i< num_params; i++) {
rands[i+2] = par[i];
}
/* now calculate the function at the existing points */
for (i = 0; i < num_data; i++) {
rands[0] = fit_x[i];
rands[1] = fit_y[i];
data[i] = ((double (*) (int, double *) )current_fun) // ouch!
(rators, rands);
}
free(rands);
}
/* /\***************************************************************** */
/* evaluate the scheme function */
/* *****************************************************************\/ */
/* static void call_scheme(par, data) */
/* double *par; */
/* double *data; */
/* { */
/* int rators = 2 + num_params; */
/* Scheme_Object ** rands = */
/* scheme_malloc(rators * sizeof(Scheme_Object)); */
/* int i; */
/* /\* set up the constant params *\/ */
/* for(i = 0 ; i< num_params; i++) { */
/* rands[i+2] = scheme_make_double(par[i]); */
/* } */
/* /\* now calculate the function at the existing points *\/ */
/* for (i = 0; i < num_data; i++) { */
/* rands[0] = scheme_make_double(fit_x[i]); */
/* rands[1] = scheme_make_double(fit_y[i]); */
/* data[i] = scheme_real_to_double(scheme_apply(current_fun, rators, rands)); */
/* } */
/* } */
/*****************************************************************
Frame routine for the marquardt-fit
*****************************************************************/
static TBOOLEAN regress(a)
double a[];
{
double **covar, *dpar, **C, chisq, last_chisq, lambda;
int iter, i, j;
marq_res_t res;
chisq = last_chisq = INFINITY;
C = matr(num_data + num_params, num_params);
lambda = -1; /* use sign as flag */
iter = 0; /* iteration counter */
/* Initialize internal variables and 1st chi-square check */
if ((res = marquardt(a, C, &chisq, &lambda)) == ERROR)
return 0; /* an error occurded */
res = BETTER;
/* show_fit(iter, chisq, chisq, a, lambda, STANDARD); */
/* MAIN FIT LOOP: do the regression iteration */
do {
if (res == BETTER) {
iter++;
last_chisq = chisq;
}
if ((res = marquardt(a, C, &chisq, &lambda)) == BETTER)
{};
/* show_fit(iter, chisq, last_chisq, a, lambda, STANDARD); */
} while ((res != ERROR)
&& (lambda < MAX_LAMBDA)
&& ((maxiter == 0) || (iter <= maxiter))
&& (res == WORSE
|| ((chisq > NEARLY_ZERO)
? ((last_chisq - chisq) / chisq)
: (last_chisq - chisq)) > epsilon
)
);
/* fit done */
/* save all the info that was otherwise printed out */
rms = sqrt(chisq / (num_data - num_params));
varience = chisq / (num_data - num_params);
asym_error = malloc (num_params * sizeof (double));
asym_error_percent = malloc (num_params * sizeof (double)) ;
/* don't know what the following code does... */
/* compute covar[][] directly from C */
Givens(C, 0, 0, 0, num_data, num_params, 0);
covar = C + num_data;
Invert_RtR(C, covar, num_params);
dpar = vec(num_params);
for (i = 0; i < num_params; i++) {
/* FIXME: can this still happen ? */
if (covar[i][i] <= 0.0) /* HBB: prevent floating point exception later on */
return 0; /* Eex("Calculation error: non-positive diagonal element in covar. matrix"); */
dpar[i] = sqrt(covar[i][i]);
}
/* transform covariances into correlations */
for (i = 0; i < num_params; i++) {
/* only lower triangle needs to be handled */
for (j = 0; j <= i; j++)
covar[i][j] /= dpar[i] * dpar[j];
}
/* scale parameter errors based on chisq */
chisq = sqrt(chisq / (num_data - num_params));
for (i = 0; i < num_params; i++)
dpar[i] *= chisq;
for(i = 0; i< num_params; i++)
{
double temp =
(fabs(a[i]) < NEARLY_ZERO) ? 0.0 : fabs(100.0 * dpar[i] / a[i]);
asym_error[i] = dpar[i];
asym_error_percent[i] = temp;
}
return 1;
/******** CRAP LEFT OVER FROM GNUPLOT ***********/
/* HBB 970304: the maxiter patch: */
/*
if ((maxiter > 0) && (iter > maxiter)) {
Dblf2("\nMaximum iteration count (%d) reached. Fit stopped.\n", maxiter);
} else {
Dblf2("\nAfter %d iterations the fit converged.\n", iter);
}
Dblf2("final sum of squares of residuals : %g\n", chisq);
if (chisq > NEARLY_ZERO) {
Dblf2("rel. change during last iteration : %g\n\n", (chisq - last_chisq) / chisq);
} else {
Dblf2("abs. change during last iteration : %g\n\n", (chisq - last_chisq));
}
if (res == ERROR)
// Eex("FIT: error occurred during fit");
*/
/* compute errors in the parameters */
/* if (num_data == num_params) { */
/* int i; */
/* Dblf("\nExactly as many data points as there are parameters.\n"); */
/* Dblf("In this degenerate case, all errors are zero by definition.\n\n"); */
/* Dblf("Final set of parameters \n"); */
/* Dblf("======================= \n\n"); */
/* for (i = 0; i < num_params; i++) */
/* Dblf3("%-15.15s = %-15g\n", par_name[i], a[i]); */
/* } else if (chisq < NEARLY_ZERO) { */
/* int i; */
/* Dblf("\nHmmmm.... Sum of squared residuals is zero. Can't compute errors.\n\n"); */
/* Dblf("Final set of parameters \n"); */
/* Dblf("======================= \n\n"); */
/* for (i = 0; i < num_params; i++) */
/* Dblf3("%-15.15s = %-15g\n", par_name[i], a[i]); */
/* } else { */
/* Dblf2("degrees of freedom (ndf) : %d\n", num_data - num_params); */
/* Dblf2("rms of residuals (stdfit) = sqrt(WSSR/ndf) : %g\n", sqrt(chisq / (num_data - num_params))); */
/* Dblf2("variance of residuals (reduced chisquare) = WSSR/ndf : %g\n\n", chisq / (num_data - num_params)); */
/* /\* get covariance-, Korrelations- and Kurvature-Matrix *\/ */
/* /\* and errors in the parameters *\/ */
/* /\* compute covar[][] directly from C *\/ */
/* Givens(C, 0, 0, 0, num_data, num_params, 0); */
/* /\*printmatrix(C, num_params, num_params); *\/ */
/* /\* Use lower square of C for covar *\/ */
/* covar = C + num_data; */
/* Invert_RtR(C, covar, num_params); */
/* /\*printmatrix(covar, num_params, num_params); *\/ */
/* /\* calculate unscaled parameter errors in dpar[]: *\/ */
/* dpar = vec(num_params); */
/* for (i = 0; i < num_params; i++) { */
/* /\* FIXME: can this still happen ? *\/ */
/* if (covar[i][i] <= 0.0) /\* HBB: prevent floating point exception later on *\/ */
/* Eex("Calculation error: non-positive diagonal element in covar. matrix"); */
/* dpar[i] = sqrt(covar[i][i]); */
/* } */
/* /\* transform covariances into correlations *\/ */
/* for (i = 0; i < num_params; i++) { */
/* /\* only lower triangle needs to be handled *\/ */
/* for (j = 0; j <= i; j++) */
/* covar[i][j] /= dpar[i] * dpar[j]; */
/* } */
/* /\* scale parameter errors based on chisq *\/ */
/* chisq = sqrt(chisq / (num_data - num_params)); */
/* for (i = 0; i < num_params; i++) */
/* dpar[i] *= chisq; */
/* Dblf("Final set of parameters Asymptotic Standard Error\n"); */
/* Dblf("======================= ==========================\n\n"); */
/* for (i = 0; i < num_params; i++) { */
/* double temp = */
/* (fabs(a[i]) < NEARLY_ZERO) ? 0.0 : fabs(100.0 * dpar[i] / a[i]); */
/* Dblf6("%-15.15s = %-15g %-3.3s %-12.4g (%.4g%%)\n", */
/* par_name[i], a[i], PLUSMINUS, dpar[i], temp); */
/* } */
/* Dblf("\n\ncorrelation matrix of the fit parameters:\n\n"); */
/* Dblf(" "); */
/* for (j = 0; j < num_params; j++) */
/* Dblf2("%-6.6s ", par_name[j]); */
/* Dblf("\n"); */
/* for (i = 0; i < num_params; i++) { */
/* Dblf2("%-15.15s", par_name[i]); */
/* for (j = 0; j <= i; j++) { */
/* /\* Only print lower triangle of symmetric matrix *\/ */
/* Dblf2("%6.3f ", covar[i][j]); */
/* } */
/* Dblf("\n"); */
/* } */
/* free(dpar); */
/* } */
return 1;
}
/*****************************************************************
display actual state of the fit
*****************************************************************/
/* static void show_fit(i, chisq, last_chisq, a, lambda, device) */
/* int i; */
/* double chisq; */
/* double last_chisq; */
/* double *a; */
/* double lambda; */
/* FILE *device; */
//{
/*
int k;
fprintf(device, "\n\n\
Iteration %d\n\
WSSR : %-15g delta(WSSR)/WSSR : %g\n\
delta(WSSR) : %-15g limit for stopping : %g\n\
lambda : %g\n\n%s parameter values\n\n",
i, chisq, chisq > NEARLY_ZERO ? (chisq - last_chisq) / chisq : 0.0,
chisq - last_chisq, epsilon, lambda,
(i > 0 ? "resultant" : "initial set of free"));
for (k = 0; k < num_params; k++)
fprintf(device, "%-15.15s = %g\n", par_name[k], a[k]);
*/
//}
/*****************************************************************
Interface to scheme
*****************************************************************/
MZ_DLLEXPORT
double * do_fit(void * function,
int n_values,
double * x_values,
double * y_values,
double * z_values,
double * errors,
int n_parameters,
double * parameters) {
/* reset lambda and other parameters if desired */
int i;
current_fun = function;
num_data = n_values;
fit_x = x_values;
fit_y = y_values;
fit_z = z_values; /* value is stored in z */
err_data = errors;
a = parameters;
num_params = n_parameters;
/* redim_vec(&a, num_params); */
/* par_name = (fixstr *) gp_realloc(par_name, (num_params + 1) * sizeof(fixstr), "fit param"); */
/* avoid parameters being equal to zero */
for (i = 0; i < num_params; i++) {
if (a[i] == 0) {
a[i] = NEARLY_ZERO;
}
}
if(regress(a)) {
gc_cleanup();
return a;
}
else { /* something went wrong */
gc_cleanup();
return NULL;
}
}

View File

@ -1,62 +0,0 @@
/* $Id: fit.h,v 1.5 2005/03/15 23:19:40 eli Exp $ */
/* GNUPLOT - fit.h */
/* NOTICE: Change of Copyright Status
*
* The author of this module, Carsten Grammes, has expressed in
* personal email that he has no more interest in this code, and
* doesn't claim any copyright. He has agreed to put this module
* into the public domain.
*
* Lars Hecking 15-02-1999
*/
/*
* Header file: public functions in fit.c
*
*
* Copyright of this module: Carsten Grammes, 1993
* Experimental Physics, University of Saarbruecken, Germany
*
* Internet address: cagr@rz.uni-sb.de
*
* Permission to use, copy, and distribute this software and its
* documentation for any purpose with or without fee is hereby granted,
* provided that the above copyright notice appear in all copies and
* that both that copyright notice and this permission notice appear
* in supporting documentation.
*
* This software is provided "as is" without express or implied warranty.
*/
#include "dllexport.h"
#ifdef __GNUC__
# define UNUSED __attribute__((unused))
#else
# define UNUSED
#endif
MZ_DLLEXPORT
double * do_fit(void * function,
int n_values,
double * x_values,
double * y_values,
double * z_values,
double * errors,
int n_parameters,
double * parameters);
MZ_DLLEXPORT
double get_rms();
MZ_DLLEXPORT
double get_varience();
MZ_DLLEXPORT
double * get_asym_error();
MZ_DLLEXPORT
double * get_asym_error_percent();

View File

@ -1,315 +0,0 @@
/* NOTICE: Change of Copyright Status
*
* The author of this module, Carsten Grammes, has expressed in
* personal email that he has no more interest in this code, and
* doesn't claim any copyright. He has agreed to put this module
* into the public domain.
*
* Lars Hecking 15-02-1999
*/
/*
* Matrix algebra, part of
*
* Nonlinear least squares fit according to the
* Marquardt-Levenberg-algorithm
*
* added as Patch to Gnuplot (v3.2 and higher)
* by Carsten Grammes
* Experimental Physics, University of Saarbruecken, Germany
*
* Internet address: cagr@rz.uni-sb.de
*
* Copyright of this module: Carsten Grammes, 1993
*
* Permission to use, copy, and distribute this software and its
* documentation for any purpose with or without fee is hereby granted,
* provided that the above copyright notice appear in all copies and
* that both that copyright notice and this permission notice appear
* in supporting documentation.
*
* This software is provided "as is" without express or implied warranty.
*/
#define NULL 0
#define null 0
#include "fit.h"
#include "matrix.h"
#include <math.h>
#include <stdlib.h>
// create a simple gc malloc...
typedef struct Node {
struct Node * next;
void * ptr;
} Node;
Node * head = null;
void * my_gc_malloc(int size) {
void * ptr = malloc(size);
Node * n = (Node *)malloc(sizeof (Node));
n->ptr = ptr;
n->next = head;
head = n;
return ptr;
}
void gc_cleanup(){
while(head) {
Node * current = head;
head = current->next;
free(current->ptr);
free(current);
}
}
/*****************************************************************/
#define Swap(a,b) {double temp = (a); (a) = (b); (b) = temp;}
#define WINZIG 1e-30
/*****************************************************************
internal prototypes
*****************************************************************/
static int fsign (double x);
/*****************************************************************
first straightforward vector and matrix allocation functions
*****************************************************************/
MZ_DLLEXPORT
double *vec (n)
int n;
{
/* allocates a double vector with n elements */
double *dp;
if( n < 1 )
return (double *) NULL;
dp = (double *) my_gc_malloc (n * sizeof(double));
return dp;
}
MZ_DLLEXPORT
double **matr (rows, cols)
int rows;
int cols;
{
/* allocates a double matrix */
register int i;
register double **m;
if ( rows < 1 || cols < 1 )
return NULL;
m = (double **) my_gc_malloc (rows * sizeof(double *));
m[0] = (double *) my_gc_malloc (rows * cols * sizeof(double));
for ( i = 1; i<rows ; i++ )
m[i] = m[i-1] + cols;
return m;
}
void free_matr (m)
double **m;
{
free (m[0]);
free (m);
}
MZ_DLLEXPORT
double *redim_vec (v, n)
double **v;
int n;
{
if ( n < 1 )
*v = NULL;
else
*v = (double *) my_gc_malloc( n * sizeof(double));
return *v;
}
MZ_DLLEXPORT
void redim_ivec (v, n)
int **v;
int n;
{
if ( n < 1 ) {
*v = NULL;
return;
}
*v = (int *) my_gc_malloc ( n * sizeof(int));
}
/* HBB: TODO: is there a better value for 'epsilon'? how to specify
* 'inline'? is 'fsign' really not available elsewhere? use
* row-oriented version (p. 309) instead?
*/
static int fsign(x)
double x;
{
return( x>0 ? 1 : (x < 0) ? -1 : 0) ;
}
/*****************************************************************
Solve least squares Problem C*x+d = r, |r| = min!, by Given rotations
(QR-decomposition). Direct implementation of the algorithm
presented in H.R.Schwarz: Numerische Mathematik, 'equation'
number (7.33)
If 'd == NULL', d is not accesed: the routine just computes the QR
decomposition of C and exits.
If 'want_r == 0', r is not rotated back (\hat{r} is returned
instead).
*****************************************************************/
MZ_DLLEXPORT
void Givens (C, d, x, r, N, n, want_r)
double **C;
double *d;
double *x;
double *r;
int N;
int n;
int want_r;
{
int i, j, k;
double w, gamma, sigma, rho, temp;
double epsilon = 1e-5; /* FIXME (?)*/
/*
* First, construct QR decomposition of C, by 'rotating away'
* all elements of C below the diagonal. The rotations are
* stored in place as Givens coefficients rho.
* Vector d is also rotated in this same turn, if it exists
*/
for (j = 0; j<n; j++)
for (i = j+1; i<N; i++)
if (C[i][j]) {
if (fabs(C[j][j])<epsilon*fabs(C[i][j])) { /* find the rotation parameters */
w = -C[i][j];
gamma = 0;
sigma = 1;
rho = 1;
} else {
w = fsign(C[j][j])*sqrt(C[j][j]*C[j][j] + C[i][j]*C[i][j]);
if (w == 0) {
// Eex3 ( "w = 0 in Givens(); Cjj = %g, Cij = %g", C[j][j], C[i][j]);
}
gamma = C[j][j]/w;
sigma = -C[i][j]/w;
rho = (fabs(sigma)<gamma) ? sigma : fsign(sigma)/gamma;
}
C[j][j] = w;
C[i][j] = rho; /* store rho in place, for later use */
for (k = j+1; k<n; k++) { /* rotation on index pair (i,j) */
temp = gamma*C[j][k] - sigma*C[i][k];
C[i][k] = sigma*C[j][k] + gamma*C[i][k];
C[j][k] = temp;
}
if (d) { /* if no d vector given, don't use it */
temp = gamma*d[j] - sigma*d[i]; /* rotate d */
d[i] = sigma*d[j] + gamma*d[i];
d[j] = temp;
}
}
if (!d) /* stop here if no d was specified */
return;
for (i = n-1; i >= 0; i--) { /* solve R*x+d = 0, by backsubstitution */
double s = d[i];
r[i] = 0; /* ... and also set r[i] = 0 for i<n */
for (k = i+1; k<n; k++)
s += C[i][k]*x[k];
if (C[i][i] == 0) {
//Eex ( "Singular matrix in Givens()");
}
x[i] = - s / C[i][i];
}
for (i = n; i < N; i++)
r[i] = d[i]; /* set the other r[i] to d[i] */
if (!want_r) /* if r isn't needed, stop here */
return;
/* rotate back the r vector */
for (j = n-1; j >= 0; j--)
for (i = N-1; i >= 0; i--) {
if ((rho = C[i][j]) == 1) { /* reconstruct gamma, sigma from stored rho */
gamma = 0;
sigma = 1;
} else if (fabs(rho)<1) {
sigma = rho;
gamma = sqrt(1-sigma*sigma);
} else {
gamma = 1/fabs(rho);
sigma = fsign(rho)*sqrt(1-gamma*gamma);
}
temp = gamma*r[j] + sigma*r[i]; /* rotate back indices (i,j) */
r[i] = -sigma*r[j] + gamma*r[i];
r[j] = temp;
}
}
/* Given a triangular Matrix R, compute (R^T * R)^(-1), by forward
* then back substitution
*
* R, I are n x n Matrices, I is for the result. Both must already be
* allocated.
*
* Will only calculate the lower triangle of I, as it is symmetric
*/
MZ_DLLEXPORT
void Invert_RtR ( R, I, n)
double **R;
double **I;
int n;
{
int i, j, k;
/* fill in the I matrix, and check R for regularity : */
for (i = 0; i<n; i++) {
for (j = 0; j<i; j++) /* upper triangle isn't needed */
I[i][j] = 0;
I[i][i] = 1;
if (! R[i][i])
{
// Eex ("Singular matrix in Invert_RtR");
}
}
/* Forward substitution: Solve R^T * B = I, store B in place of I */
for (k = 0; k<n; k++)
for (i = k; i<n; i++) { /* upper half needn't be computed */
double s = I[i][k];
for (j = k; j<i; j++) /* for j<k, I[j][k] always stays zero! */
s -= R[j][i] * I[j][k];
I[i][k] = s / R[i][i];
}
/* Backward substitution: Solve R * A = B, store A in place of B */
for (k = 0; k<n; k++)
for (i = n-1; i >= k; i--) { /* don't compute upper triangle of A */
double s = I[i][k];
for (j = i+1; j<n; j++)
s -= R[i][j] * I[j][k];
I[i][k] = s / R[i][i];
}
}

View File

@ -1,75 +0,0 @@
/* $Id: matrix.h,v 1.5 2005/03/15 23:23:56 eli Exp $ */
/* GNUPLOT - matrix.h */
/* NOTICE: Change of Copyright Status
*
* The author of this module, Carsten Grammes, has expressed in
* personal email that he has no more interest in this code, and
* doesn't claim any copyright. He has agreed to put this module
* into the public domain.
*
* Lars Hecking 15-02-1999
*/
/*
* Header file: public functions in matrix.c
*
*
* Copyright of this module: Carsten Grammes, 1993
* Experimental Physics, University of Saarbruecken, Germany
*
* Internet address: cagr@rz.uni-sb.de
*
* Permission to use, copy, and distribute this software and its
* documentation for any purpose with or without fee is hereby granted,
* provided that the above copyright notice appear in all copies and
* that both that copyright notice and this permission notice appear
* in supporting documentation.
*
* This software is provided "as is" without express or implied warranty.
*/
#ifndef MATRIX_H
#define MATRIX_H
#include "dllexport.h"
#ifdef EXT
#undef EXT
#endif
#ifdef MATRIX_MAIN
#define EXT
#else
#define EXT extern
#endif
/******* public functions ******/
MZ_DLLEXPORT
EXT double *vec (int n);
MZ_DLLEXPORT
EXT int *ivec (int n);
MZ_DLLEXPORT
EXT double **matr (int r, int c);
MZ_DLLEXPORT
EXT double *redim_vec (double **v, int n);
MZ_DLLEXPORT
EXT void redim_ivec (int **v, int n);
EXT void solve (double **a, int n, double **b, int m);
MZ_DLLEXPORT
EXT void Givens (double **C, double *d, double *x, double *r, int N, int n, int want_r);
MZ_DLLEXPORT
EXT void Invert_RtR (double **R, double **I, int n);
#endif
// a kludgy version of a malloc
void * my_gc_malloc(int size);
void gc_cleanup();

View File

@ -84,7 +84,6 @@
["libpangocairo-1.0-0.dll" 94625]
["libpangowin32-1.0-0.dll" 102210]
["libpangoft2-1.0-0.dll" 679322]
["libfit.dll" 73728]
,@(if (getenv "PLT_WIN_GTK")
'(["libatk-1.0-0.dll" 153763]
["libgtk-win32-2.0-0.dll" 4740156]
@ -110,9 +109,8 @@
["libgthread-2.0-0.dll" 126615]
["libpangocairo-1.0-0.dll" 185168]
["libpangowin32-1.0-0.dll" 192656]
["libpangoft2-1.0-0.dll" 1188615]
["libfit.dll" 69120]]]
;; Databse libraries
["libpangoft2-1.0-0.dll" 1188615]]]
;; Database libraries
[db
[win32/i386
["sqlite3.dll" 570947]]

View File

@ -31,7 +31,6 @@ AC_ARG_ENABLE(jit, [ --enable-jit support JIT compiler (enabled
AC_ARG_ENABLE(foreign, [ --enable-foreign support foreign calls (enabled by default)], , enable_foreign=yes)
AC_ARG_ENABLE(places, [ --enable-places support places (3m only; usually enabled by default)])
AC_ARG_ENABLE(futures, [ --enable-futures support futures (usually enabled by default)])
AC_ARG_ENABLE(plot, [ --enable-plot support plot libraries (enabled by default)], , enable_plot=yes )
AC_ARG_ENABLE(float, [ --enable-float support single-precision floats (enabled by default)], , enable_float=yes)
AC_ARG_ENABLE(floatinstead, [ --enable-floatinstead use single-precision by default])
@ -239,7 +238,6 @@ show_explicitly_enabled "${enable_xonx}" "Unix style"
show_explicitly_enabled "${enable_shared}" "Shared libraries"
show_explicitly_disabled "${enable_gracket}" GRacket
show_explicitly_disabled "${enable_plot}" Plot fit library
if test "$LIBTOOLPROG" != "" ; then
echo "=== Libtool program: $LIBTOOLPROG"
@ -1262,7 +1260,6 @@ AC_SUBST(MRLIBINSTALL)
AC_SUBST(LIBFINISH)
AC_SUBST(MAKE_GRACKET)
AC_SUBST(MAKE_FIT)
AC_SUBST(MAKE_COPYTREE)
AC_SUBST(MAKE_FINISH)
@ -1340,14 +1337,6 @@ fi
makefiles="$makefiles foreign/Makefile"
ac_configure_args="$ac_configure_args$SUB_CONFIGURE_EXTRAS"
if test -d "${srcdir}/fit" && test "${enable_plot}" = "yes" ; then
makefiles="$makefiles
fit/Makefile"
MAKE_FIT=fit
else
MAKE_FIT=no
fi
if test "${enable_gracket}" = "yes" ; then
makefiles="$makefiles
gracket/Makefile