Menu
▾
▴
octave-cvsupdate — CVS update notifications
This list is closed, nobody may subscribe to it.
| 2001 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
(10) |
Aug
(5) |
Sep
(3) |
Oct
(41) |
Nov
(41) |
Dec
(33) |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 2002 |
Jan
(75) |
Feb
(10) |
Mar
(170) |
Apr
(174) |
May
(66) |
Jun
(11) |
Jul
(10) |
Aug
(44) |
Sep
(73) |
Oct
(28) |
Nov
(139) |
Dec
(52) |
| 2003 |
Jan
(35) |
Feb
(93) |
Mar
(62) |
Apr
(10) |
May
(55) |
Jun
(70) |
Jul
(37) |
Aug
(16) |
Sep
(56) |
Oct
(31) |
Nov
(57) |
Dec
(83) |
| 2004 |
Jan
(85) |
Feb
(67) |
Mar
(27) |
Apr
(37) |
May
(75) |
Jun
(85) |
Jul
(160) |
Aug
(68) |
Sep
(104) |
Oct
(25) |
Nov
(39) |
Dec
(23) |
| 2005 |
Jan
(10) |
Feb
(45) |
Mar
(43) |
Apr
(19) |
May
(108) |
Jun
(31) |
Jul
(41) |
Aug
(23) |
Sep
(65) |
Oct
(58) |
Nov
(44) |
Dec
(54) |
| 2006 |
Jan
(96) |
Feb
(27) |
Mar
(69) |
Apr
(59) |
May
(67) |
Jun
(35) |
Jul
(13) |
Aug
(461) |
Sep
(160) |
Oct
(399) |
Nov
(32) |
Dec
(72) |
| 2007 |
Jan
(316) |
Feb
(305) |
Mar
(318) |
Apr
(54) |
May
(194) |
Jun
(173) |
Jul
(282) |
Aug
(91) |
Sep
(227) |
Oct
(365) |
Nov
(168) |
Dec
(18) |
| 2008 |
Jan
(71) |
Feb
(111) |
Mar
(155) |
Apr
(173) |
May
(70) |
Jun
(67) |
Jul
(55) |
Aug
(83) |
Sep
(32) |
Oct
(68) |
Nov
(80) |
Dec
(29) |
| 2009 |
Jan
(46) |
Feb
(18) |
Mar
(95) |
Apr
(76) |
May
(140) |
Jun
(98) |
Jul
(84) |
Aug
(123) |
Sep
(94) |
Oct
(131) |
Nov
(142) |
Dec
(125) |
| 2010 |
Jan
(128) |
Feb
(158) |
Mar
(172) |
Apr
(134) |
May
(94) |
Jun
(84) |
Jul
(32) |
Aug
(127) |
Sep
(167) |
Oct
(109) |
Nov
(69) |
Dec
(78) |
| 2011 |
Jan
(39) |
Feb
(58) |
Mar
(52) |
Apr
(47) |
May
(56) |
Jun
(76) |
Jul
(55) |
Aug
(54) |
Sep
(165) |
Oct
(255) |
Nov
(328) |
Dec
(263) |
| 2012 |
Jan
(82) |
Feb
(147) |
Mar
(400) |
Apr
(216) |
May
(209) |
Jun
(160) |
Jul
(86) |
Aug
(141) |
Sep
(156) |
Oct
(6) |
Nov
|
Dec
|
| 2015 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
(1) |
Aug
|
Sep
(1) |
Oct
|
Nov
(1) |
Dec
(2) |
| 2016 |
Jan
|
Feb
(2) |
Mar
(2) |
Apr
(1) |
May
(1) |
Jun
(2) |
Jul
(1) |
Aug
(1) |
Sep
|
Oct
|
Nov
(1) |
Dec
|
| 2019 |
Jan
|
Feb
|
Mar
(1) |
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
|
Oct
|
Nov
|
Dec
|
| 2021 |
Jan
|
Feb
|
Mar
|
Apr
(3) |
May
(4) |
Jun
(8) |
Jul
(2) |
Aug
(5) |
Sep
(9) |
Oct
|
Nov
|
Dec
|
|
From: <par...@us...> - 2012-05-31 13:49:49
|
Revision: 10548
http://octave.svn.sourceforge.net/octave/?rev=10548&view=rev
Author: paramaniac
Date: 2012-05-31 13:49:40 +0000 (Thu, 31 May 2012)
Log Message:
-----------
control-devel: minor fixes
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
Modified: trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 13:34:10 UTC (rev 10547)
+++ trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 13:49:40 UTC (rev 10548)
@@ -76,7 +76,7 @@
s = [];
n = [];
conf = [];
- noise = "n"
+ noise = "n";
## handle keys and values
for k = 1 : 2 : nkv
@@ -98,9 +98,6 @@
conf = logical (val);
case "noise"
noise = val;
- % none
- % e
- % v (normalized)
otherwise
warning ("%s: invalid property name '%s' ignored", method, key);
endswitch
@@ -143,19 +140,19 @@
## assemble model
[inname, outname] = get (dat, "inname", "outname");
- if (strncmpi (noise, "e", 1))
+ if (strncmpi (noise, "e", 1)) # add error inputs e, not normalized
sys = ss (a, [b, k], c, [d, eye(p)], tsam);
in_u = __labels__ (inname, "u");
in_e = __labels__ (outname, "y");
in_e = cellfun (@(x) ["e@", x], in_e, "uniformoutput", false);
inname = [in_u; in_e];
- elseif (strncmpi (noise, "v", 1))
+ elseif (strncmpi (noise, "v", 1)) # add error inputs v, normalized
sys = ss (a, [b, k*l], c, [d, l], tsam);
in_u = __labels__ (inname, "u");
in_v = __labels__ (outname, "y");
in_v = cellfun (@(x) ["v@", x], in_v, "uniformoutput", false);
inname = [in_u; in_v];
- else
+ else # no error inputs, default
sys = ss (a, b, c, d, tsam);
endif
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-31 13:34:19
|
Revision: 10547
http://octave.svn.sourceforge.net/octave/?rev=10547&view=rev
Author: paramaniac
Date: 2012-05-31 13:34:10 +0000 (Thu, 31 May 2012)
Log Message:
-----------
control-devel: return noise inputs e or v if requested
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
Modified: trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 12:31:34 UTC (rev 10546)
+++ trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 13:34:10 UTC (rev 10547)
@@ -57,7 +57,7 @@
error ("%s: keys and values must come in pairs", method);
endif
- [ns, l, m, e] = size (dat); # dataset dimensions
+ [ns, p, m, e] = size (dat); # dataset dimensions
tsam = dat.tsam;
## multi-experiment data requires equal sampling times
@@ -76,6 +76,7 @@
s = [];
n = [];
conf = [];
+ noise = "n"
## handle keys and values
for k = 1 : 2 : nkv
@@ -95,6 +96,11 @@
rcond = val;
case "confirm"
conf = logical (val);
+ case "noise"
+ noise = val;
+ % none
+ % e
+ % v (normalized)
otherwise
warning ("%s: invalid property name '%s' ignored", method, key);
endswitch
@@ -103,7 +109,7 @@
## handle s/nobr and n
nsmp = sum (ns); # total number of samples
- nobr = fix ((nsmp+1)/(2*(m+l+1)));
+ nobr = fix ((nsmp+1)/(2*(m+p+1)));
if (e > 1)
nobr = min (nobr, fix (min (ns) / 2));
endif
@@ -132,9 +138,29 @@
## perform system identification
[a, b, c, d, q, ry, s, k, x0] = slident (dat.y, dat.u, nobr, n, meth, alg, conct, ctrl, rcond, tol);
+ ## L L' = Ry, e = L v, v becomes white noise with identity covariance matrix
+ l = chol (ry, "lower");
+
## assemble model
- sys = ss (a, b, c, d, tsam);
+ [inname, outname] = get (dat, "inname", "outname");
+ if (strncmpi (noise, "e", 1))
+ sys = ss (a, [b, k], c, [d, eye(p)], tsam);
+ in_u = __labels__ (inname, "u");
+ in_e = __labels__ (outname, "y");
+ in_e = cellfun (@(x) ["e@", x], in_e, "uniformoutput", false);
+ inname = [in_u; in_e];
+ elseif (strncmpi (noise, "v", 1))
+ sys = ss (a, [b, k*l], c, [d, l], tsam);
+ in_u = __labels__ (inname, "u");
+ in_v = __labels__ (outname, "y");
+ in_v = cellfun (@(x) ["v@", x], in_v, "uniformoutput", false);
+ inname = [in_u; in_v];
+ else
+ sys = ss (a, b, c, d, tsam);
+ endif
+ sys = set (sys, "inname", inname, "outname", outname);
+
## return x0 as vector for single-experiment data
## instead of a cell containing one vector
if (numel (x0) == 1)
@@ -146,8 +172,7 @@
## state covariance matrix Q
## output covariance matrix Ry
## state-output cross-covariance matrix S
- ## L L' = Ry, e = L v, v becomes white noise with identity covariance matrix
- info = struct ("K", k, "Q", q, "Ry", ry, "S", s, "L", chol (ry, "lower"));
+ info = struct ("K", k, "Q", q, "Ry", ry, "S", s, "L", l);
endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-31 12:31:43
|
Revision: 10546
http://octave.svn.sourceforge.net/octave/?rev=10546&view=rev
Author: paramaniac
Date: 2012-05-31 12:31:34 +0000 (Thu, 31 May 2012)
Log Message:
-----------
control-devel: return scaling L
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m
trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
Modified: trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-05-31 06:04:54 UTC (rev 10545)
+++ trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-05-31 12:31:34 UTC (rev 10546)
@@ -48,9 +48,12 @@
dat = iddata (Y, U)
-[sys, x0] = moen4 (dat, 's', 10, 'n', 5) % s=10, n=5
+[sys, x0, info] = n4sid (dat, 's', 10, 'n', 5) % s=10, n=5
+L = chol (info.Ry, 'lower')
+Be = info.K * L
+
[y, t] = lsim (sys, U, [], x0);
err = norm (Y - y, 1) / norm (Y, 1)
Modified: trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 06:04:54 UTC (rev 10545)
+++ trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-31 12:31:34 UTC (rev 10546)
@@ -146,7 +146,8 @@
## state covariance matrix Q
## output covariance matrix Ry
## state-output cross-covariance matrix S
- info = struct ("K", k, "Q", q, "Ry", ry, "S", s);
+ ## L L' = Ry, e = L v, v becomes white noise with identity covariance matrix
+ info = struct ("K", k, "Q", q, "Ry", ry, "S", s, "L", chol (ry, "lower"));
endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <be...@us...> - 2012-05-31 06:05:06
|
Revision: 10545
http://octave.svn.sourceforge.net/octave/?rev=10545&view=rev
Author: benjf5
Date: 2012-05-31 06:04:54 +0000 (Thu, 31 May 2012)
Log Message:
-----------
Assembled the first solid steps towards a fastlscomplex C++ function. Don't recommend compiling it.
Added Paths:
-----------
trunk/octave-forge/extra/lssa/fastlscomplex.cpp
Removed Paths:
-------------
trunk/octave-forge/extra/lssa/fastlscomplex.c
Deleted: trunk/octave-forge/extra/lssa/fastlscomplex.c
===================================================================
--- trunk/octave-forge/extra/lssa/fastlscomplex.c 2012-05-30 18:38:57 UTC (rev 10544)
+++ trunk/octave-forge/extra/lssa/fastlscomplex.c 2012-05-31 06:04:54 UTC (rev 10545)
@@ -1,243 +0,0 @@
-// fastlscomplex, implemented for Octave
-
-
-#include <octave/oct.h>
-#ifdef __cplusplus
-extern "C"
-{
-#endif
-#include <stdio.h>
-#include <stdlib.h>
-#include <string.h>
-#include <math.h>
-#include <complex.h>
-#ifdef __cplusplus
-}
-#endif
-
-#define MAXCOLUMN 8
-
-/* We use the ISO C99 complex facility as implemented by GCC, but only in two places */
-
-typedef _Complex double Complex;
-typedef double Real;
-
-#define RE(x_) (__real__ x_)
-#define IM(x_) (__imag__ x_)
-#define CSET(x_, r_, i_) (RE(x_) = (r_), IM(x_) = (i_))
-#define PHISET(x_, p_) CSET(x_, cos(tmp=p_), sin(tmp))
-#define SCALEPHISET(x_, f_, p_) CSET(x_, f_ cos(tmp=p_), f_ sin(tmp))
-
-// Here the Data structure is included, but it's only ever a feature of the
-// #STANDALONE cases, which we're not dealing with here.
-
-typedef struct
-{ Real x, t;
-} XTElem;
-
-inline double sqr(double x)
-{ return x*x; }
-
-/* PNUM has to match the definition of EXP_IOT_SERIES! */
-#define PNUM 12
-#define SETXT(p_, op_, x_, t_) (p_)->x op_ x_; (p_++)->t op_ t_;
-#define SETT(p_, op_, x_, t_) *p_++ op_ t_;
-#define SETX(p_, op_, x_, t_) *p_++ op_ x_;
-/* h is a complex aux. variable; it is used for assignment times I everywhere */
-#define SETIX(p_, op_, x_, t_) h = x_; RE(*(p_)) op_ -IM(h); IM(*(p_)) op_ RE(h); p_++;
-
- /* Macro that sums up the power series terms into the power series
- * element record pointed to by p_.
- * By using = and += for op_, initial setting and accumulation can be selected.
- * t_ is the expression specifying the abscissa value. set_ can be either
- * SETXT to set the x and t fields of an XTElem record, or SETT/SETX to set
- * the elements of a Real array representing alternately real and imaginary
- * values.
- */
- // -10 points, comments don't match method.
-#define EXP_IOT_SERIES(p_, el_, t_, op_, setr_, seti_) \
-{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \
- tt = -t; seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(-1.0/3.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(-1.0/5.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(-1.0/7.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(-1.0/9.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(-1.0/11.0); seti_(p_, op_, x*tt, tt) \
-}
-
-/* same as the above, but without alternating signs */
-#define EXPIOT_SERIES(p_, el_, t_, op_, setr_, seti_) \
-{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \
- seti_(p_, op_, x*t, t ) \
- tt = t*t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/3.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/5.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/7.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/9.0); seti_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \
- tt *= t*(1.0/11.0); seti_(p_, op_, x*tt, tt) \
-}
-
-# define SRCARG Real *tptr, Real *xptr, int n, double *lengthptr
-# define SRCVAR int k; Real length = *lengthptr;
-# define SRCT tptr[k]
-# define SRCX xptr[k]
-# define SRCFIRST k = 0
-# define SRCAVAIL (k<n)
-
-
-
-DEFUN_DLD (fastlscomplex, args, nargout, "Computes a rapid complex least squares transform.")
-{
- int nargs = args.length();
- if ( nargs != 7 ) print_usage();
- else {
- const RowVector tvals = args(0).row_vector_value();
- const ComplexRowVector xvals = args(1).complex_row_vector_value();
- const int nval = args(2).int_value();
- const int lenval = args(3).int_value();
- const int ncoeffval = args(4).int_value();
- const int noctaveval = args(5).int_value();
- const int omegamaxval = args(6).int_value();
- if ( !error_value ) {
- ComplexRowVector ret_series ( ( ncoeffval * noctaveval ) , 0.0 );
- octave_idx_type numt = tvals.numel(), numx = xvals.numel(),
- numret = ret_series.numel();
- }
- }
-
-}
-
-void print_usage() {
- octave_stdout << "Prints the usage string for fastlscomplex, "
- << "once I get around to it.\n";
- // Really, I'll replace this soon.
-
-}
-
-void fastlscomplex(Real *tptr, Complex *xptr, int *nptr, double *lengthptr, int *ncoeffptr, int *noctaveptr, Real *omegamaxptr, Complex *rp)
-{
- int k, n = *nptr, ncoeff = *ncoeffptr, noctave = *noctaveptr;
- Real length = *lengthptr, omegamax = *omegamaxptr;
-
- struct SumVec /* Precomputation record */
- { struct SumVec *next; /* A singly linked list */
- Complex elems[PNUM]; /* the summed power series elements */
- int cnt; /* number of samples for which the power series elements were added */
- }
- *shead, *stail, *sp, *sq; /* power series element lists */
- Real dtelems[PNUM], /* power series elements of exp(-i dtau) */
- *dte, *r, /* Pointers into dtelems */
- x, /* abscissa and ordinate value, p-th power of t */
- tau, tau0, te, /* Precomputation range centers and range end */
- tau0d, /* tau_h of first summand range at level d */
- dtau = (0.5*M_PI)/omegamax,/* initial precomputation interval radius */
- dtaud, /* precomputation interval radius at d'th merging step */
- n_1 = 1.0/n, /* reciprocal of sample count */
- ooct, o, omul, /* omega/mu for octave's top omega and per band, mult. factor */
- omegaoct, omega, /* Max. frequency of octave and current frequency */
- on_1, /* n_1*(omega/mu)^p, n_1*(2*omega/mu)^p */
- mu = (0.5*M_PI)/length, /* Frequency shift: a quarter period of exp(i mu t) on length */
- tmp;
- Complex zeta, zz, /* Accumulators for spectral coefficients */
- e, emul, /* summation factor exp(-i o tau_h) */
- h, eoelems[PNUM], oeelems[PNUM],
- *eop, *oep,
- *ep, *op, /* Pointer into the current choice of eoelems and oeelems */
- *p, *q, *pe;
- int i, j; /* Coefficient and octave counter */
-
- /* Subdivision and Precomputation */
- SRCFIRST;
- tau = SRCT+dtau; te = tau+dtau;
- tau0 = tau;
- shead = stail = sp = alloca(sizeof(*sp)); sp->next = 0;
- for(te = SRCT+2*dtau; ; )
- { x = SRCX;
- EXP_IOT_SERIES(p, sp->elems, mu*(SRCT-tau), =, SETX, SETIX); sp->cnt = 1;
- for(SRCNEXT; SRCAVAIL && SRCT<te; SRCNEXT)
- { x = SRCX;
- EXP_IOT_SERIES(p, sp->elems, mu*(SRCT-tau), +=, SETX, SETIX); sp->cnt++;
- }
- if(!SRCAVAIL) break;
- tau = te+dtau; te = tau+dtau;
- sp = alloca(sizeof(*sp)); stail->next = sp; stail = sp; sp->next = 0; sp->cnt = 0;
- }
-
- ooct = omegamax/mu;
- omul = exp(-M_LN2/ncoeff);
- omegaoct = omegamax;
- tau0d = tau0;
- dtaud = dtau;
- /*** Loop over Octaves ***/
- for(j = noctave; ; ooct *= 0.5, omegaoct *= 0.5, tau0d += dtaud, dtaud *= 2)
- { /*** Results per frequency ***/
- for(i = ncoeff, o = ooct, omega = omegaoct; i--; o *= omul, omega *= omul)
- { PHISET(e, -omega*tau0d);
- PHISET(emul, -2*omega*dtaud);
- for(zeta = 0, sp = shead; sp; sp = sp->next, e *= emul)
- if(sp->cnt)
- { for(zz = 0, p = sp->elems, pe = p+PNUM, on_1 = n_1; p < pe; )
- { zz += *p++ * on_1; on_1 *= o; }
- zeta += e * zz;
- }
- *rp++ = zeta;
- }
- if(--j<=0) break; /* avoid unnecessary merging at the end */
- /* Merging of the s_h;
- * 4 different possibilities, depending on which of the merged ranges actually contain data.
- * The computation is described in the paper; The result of a merger is stored in the
- * left precomputation record (sp). Before one power series element is stored, the
- * sum and difference of the original values *p and *q are stored in eoelems and oeelems,
- * respectively. The result is then stored in *p.
- */
- EXPIOT_SERIES(r, dtelems, mu*dtaud, =, SETT, SETT);
- for(sp = shead; sp; sp = sp->next)
- { if(!(sq = sp->next) || !sq->cnt)
- { if(sp->cnt)
- for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+PNUM; p < pe; p++, dte++)
- { ep = eop; *eop++ = *p;
- for(r = dtelems, *p = *ep * *r; ; )
- { if(++r>=dte) break; --ep; h = *ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
- if(++r>=dte) break; --ep; *p -= *ep * *r;
- if(++r>=dte) break; --ep; h = -*ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
- if(++r>=dte) break; --ep; *p += *ep * *r;
- }
- }
- if(!sq) break; /* reached the last precomputation range */
- }
- else
- if(sp->cnt)
- for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+PNUM;
- p < pe; p++, q++, dte++)
- { ep = eop; *eop++ = *p+*q; *oep++ = *p-*q; op = oep;
- for(r = dtelems, *p = *ep * *r; ; )
- { if(++r>=dte) break; op -= 2; h = *op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
- if(++r>=dte) break; ep -= 2; *p -= *ep * *r;
- if(++r>=dte) break; op -= 2; h = -*op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
- if(++r>=dte) break; ep -= 2; *p += *ep * *r;
- }
- }
- else
- for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+PNUM; q<pe; q++, dte++)
- { ep = eop; *eop++ = *q;
- for(r = dtelems, *q = *ep * *r; ; )
- { if(++r>=dte) break; --ep; h = *ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
- if(++r>=dte) break; --ep; *p -= *ep * *r;
- if(++r>=dte) break; --ep; h = -*ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
- if(++r>=dte) break; --ep; *q += *ep * *r;
- }
- }
-
- sp->cnt += sq->cnt; sp->next = sq->next; /* free(sq) if malloc'ed */
- }
- }
-}
Added: trunk/octave-forge/extra/lssa/fastlscomplex.cpp
===================================================================
--- trunk/octave-forge/extra/lssa/fastlscomplex.cpp (rev 0)
+++ trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-05-31 06:04:54 UTC (rev 10545)
@@ -0,0 +1,161 @@
+/* Copyright (C) 2012 Benjamin Lewis <be...@gm...>
+ * fastlscomplex.cpp, compiles to fastlscomplex.oct
+ * This program is free software; you can redistribute it and/or modify it under
+ * the terms of the GNU General Public License as published by the Free Software
+ * Foundation; either version 3 of the License, or (at your option) any later
+ * version.
+ *
+ * This program is distributed in the hope that it will be useful, but WITHOUT
+ * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+ * details.
+ *
+ * You should have received a copy of the GNU General Public License along with
+ * this program; if not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <octave/oct.h>
+#include <complex>
+#include <string>
+#include <stdio>
+
+#define PNUM 12;
+
+// In order to reduce memory overhead, fastlscomplex will use a reference as the last input.
+
+void fastlscomplex ( RowVector tvals, ComplexRowVector xvals, int n,
+ double length, int ncoeff, int noctaves, double omegamax,
+ ComplexRowVector& result ) {
+ /* Singly-linked list which contains each precomputation record
+ * count stores the number of samples for which power series elements
+ * were added. This will be useful for accelerating computation by ignoring
+ * unused entries.
+ */
+ struct SumVec {
+ struct SumVec *next;
+ Complex<double> elements[PNUM];
+ int count; }
+ *shead, *stail, *sp, *sq; //Names have been kept from the C, may change if I want to.
+ double dtelems[PNUM], /* power series elements of exp(-i dtau) */
+ *dte, *r, /* Pointers into dtelems */
+ x, /* abscissa and ordinate value, p-th power of t */
+ tau, tau0, te, /* Precomputation range centers and range end */
+ tau0d, /* tau_h of first summand range at level d */
+ dtau = (0.5*M_PI)/omegamax,/* initial precomputation interval radius */
+ dtaud, /* precomputation interval radius at d'th merging step */
+ n_1 = 1.0/n, /* reciprocal of sample count */
+ ooct, o, omul, /* omega/mu for octave's top omega and per band, mult. factor */
+ omegaoct, omega, /* Max. frequency of octave and current frequency */
+ on_1, /* n_1*(omega/mu)^p, n_1*(2*omega/mu)^p */
+ mu = (0.5*M_PI)/length, /* Frequency shift: a quarter period of exp(i mu t) on length */
+ tmp;
+ std::complex<double> zeta, zz, // Accumulators for spectral coefficients to place in Complex<double>
+ e, emul,
+ h, eoelems[PNUM], oeelems[PNUM],
+ *eop, *oep,
+ *ep, *op,
+ *p, *q, *pe;
+
+
+ int i , j; // Counters; coefficient and octave, respectively.
+ // probable doubles: zetar, zetai, zzr, zzi, er, ei, emulr, emuli,
+ // eoelemsr[PNUM] eoelemsi[PNUM], etc. (since I want to use Complex<double>
+ // as little as possible.)
+
+ octave_idx_type k = 0;
+
+ // Subdivision and precomputation, reinvisioned in an OOWorld.
+ tau = tvals(k) + dtau;
+ te = tau + dtau;
+ tau0 = tau;
+ shead = stail = sp = alloca(sizeof(*sp));
+ sp->next = 0;
+ { te = tvals(k) + ( 2 * dtau );
+ while ( k < n ) { //THIS makes no sense, n is the wrong type. Will need to fix that.
+ EXP_IOT_SERIES(p, sp->elems, mu*(tvals(k)-tau), =, SETX, SETIX);
+ sp->count = 1;
+ //Sum terms and show that there has been at least one precomputation.
+ // I will probably replace the macro with a better explanation.
+ for(SRCNEXT; SRCAVAIL && tvals(k) < te; SRCNEXT) {
+ x = xvals(k);
+ EXP_IOT_SERIES(p,sp->elems,mu*(tvals(k)-tau), +=, SETX, SETIX);
+ sp->count++;
+ }
+ if ( k >= n ) break;
+ tau = te + dtau;
+ te = tau + dtau;
+ sp = alloca(sizeof(*sp));
+ stail->next = sp;
+ stail = sp;
+ sp->next = 0;
+ sp->count = 0;
+ } }
+ // Now isn't that just a nicer representation of much the same control structure as that ugly for-loop?
+ // Defining starting values for the loops over octaves:
+ ooct = omegamax / mu;
+ omul = exp(-M_LN2/ncoeff);
+ omegaoct = omegamax;
+ tau0d = tau0;
+ dtaud = dtau;
+ // Looping over octaves
+ for ( j = noctave ; ; ooct *= 0.5 , omegaoct *= 0.5 , tau0d += dtaud , dtaud *= 2 ) {
+ // Looping over&results per frequency
+ for ( i = ncoeff, o = ooct, omega = omegaoct; i-- ; o *= omul, omega *= omul ) {
+ e.real() = cos( - ( omega * tau0d ) ); e.imag() = sin( - ( omega * tau0d ) );
+ // sets real, imag parts of e
+ emul.real() = cos( - 2 * ( omega * dtaud ) ); emul.imag() = sin( - 2 * ( omega * dtaud ) );
+ // sets real, imag parts of emul
+ for( zeta = 0 , sp = shead; sp; sp = sp->next, e *= emul ) {
+ if ( sp->count ) {
+ for ( zz = std::complex<double>(0.0,0.0) , p = sp->elems , pe = p + PNUM , on_1 = n_1 ; p < pe ; ) {
+ zz += *p++ * on_1;
+ on_1 *= 0;
+ }
+ zeta += e * zz;
+ }
+ *rp++ = zeta;
+ }
+ if ( --j <= 0 ) break; //Avoids excess merging
+
+ EXPIOT_SERIES(r, dtelems, mu*dtaud, =, SETT, SETT);
+ for(sp = shead; sp; sp = sp->next){
+ if(!(sq = sp->next) || !sq->count ) {
+ for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+PNUM; p < pe; p++, dte++)
+ { ep = eop; *eop++ = *p;
+ for(r = dtelems, *p = *ep * *r; ; )
+ { if(++r>=dte) break; --ep; h = *ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; --ep; *p -= *ep * *r;
+ if(++r>=dte) break; --ep; h = -*ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; --ep; *p += *ep * *r;
+ }
+ }
+ if(!sq) break; /* reached the last precomputation range */
+ }
+ else
+ if(sp->cnt)
+ for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+PNUM;
+ p < pe; p++, q++, dte++)
+ { ep = eop; *eop++ = *p+*q; *oep++ = *p-*q; op = oep;
+ for(r = dtelems, *p = *ep * *r; ; )
+ { if(++r>=dte) break; op -= 2; h = *op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; ep -= 2; *p -= *ep * *r;
+ if(++r>=dte) break; op -= 2; h = -*op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; ep -= 2; *p += *ep * *r;
+ }
+ }
+ else
+ for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+PNUM; q<pe; q++, dte++)
+ { ep = eop; *eop++ = *q;
+ for(r = dtelems, *q = *ep * *r; ; )
+ { if(++r>=dte) break; --ep; h = *ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
+ if(++r>=dte) break; --ep; *p -= *ep * *r;
+ if(++r>=dte) break; --ep; h = -*ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
+ if(++r>=dte) break; --ep; *q += *ep * *r;
+ }
+ }
+
+ sp->cnt += sq->cnt; sp->next = sq->next; /* free(sq) if malloc'ed */
+ }
+ }
+ }
+}
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 18:39:03
|
Revision: 10544
http://octave.svn.sourceforge.net/octave/?rev=10544&view=rev
Author: jpicarbajal
Date: 2012-05-30 18:38:57 +0000 (Wed, 30 May 2012)
Log Message:
-----------
system-identification: nonlinear prediction error. Tisean Wrapper
Modified Paths:
--------------
trunk/octave-forge/main/system-identification/INDEX
Modified: trunk/octave-forge/main/system-identification/INDEX
===================================================================
--- trunk/octave-forge/main/system-identification/INDEX 2012-05-30 18:37:53 UTC (rev 10543)
+++ trunk/octave-forge/main/system-identification/INDEX 2012-05-30 18:38:57 UTC (rev 10544)
@@ -2,6 +2,8 @@
TISEAN wrapper
corrdim
delayvec
+ nstatz.m
pspec_mep
pspec_bin
+ recurrplt.m
Data types
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 18:38:00
|
Revision: 10543
http://octave.svn.sourceforge.net/octave/?rev=10543&view=rev
Author: jpicarbajal
Date: 2012-05-30 18:37:53 +0000 (Wed, 30 May 2012)
Log Message:
-----------
system-identification: nonlinear prediction error. Tisean Wrapper
Modified Paths:
--------------
trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
Added Paths:
-----------
trunk/octave-forge/main/system-identification/inst/tisean/nstatz.m
trunk/octave-forge/main/system-identification/inst/tisean/recurrplt.m
Added: trunk/octave-forge/main/system-identification/inst/tisean/nstatz.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/nstatz.m (rev 0)
+++ trunk/octave-forge/main/system-identification/inst/tisean/nstatz.m 2012-05-30 18:37:53 UTC (rev 10543)
@@ -0,0 +1,109 @@
+%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...>
+%%
+%% This program is free software: you can redistribute it and/or modify
+%% it under the terms of the GNU General Public License as published by
+%% the Free Software Foundation, either version 3 of the License, or
+%% any later version.
+%%
+%% This program is distributed in the hope that it will be useful,
+%% but WITHOUT ANY WARRANTY; without even the implied warranty of
+%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+%% GNU General Public License for more details.
+%%
+%% You should have received a copy of the GNU General Public License
+%% along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+%% -*- texinfo -*-
+%% @deftypefn {Function File} {@var{m} = } nstatz (@var{data}, @var{n})
+%% @deftypefn {Function File} {[@var{f} @var{p} @var{ferr}] = } nstatz (@dots{})
+%% @deftypefnx {Function File} {@dots{} = } nststz (@dots{},@var{property},@var{value})
+%% Seeks for nonstationarity in a time series by dividing it into a number of @
+%% segments and calculating the cross-forecast errors between the different @
+%% segments. The model used for the forecast is zeroth order model as proposed by @
+%% Schreiber. This function calls @code{nstat_z} from the TISEAN package.
+%%
+%% @end deftypefn
+
+function varargout = nstatz (data, n, varargin)
+
+ data = data(:);
+ [nT nc] = size (data);
+ if nc !=1
+ warning ('Tisian:InvalidInputFormat', ...
+ 'Only first column of data will be used');
+ end
+ amp = abs((max(data)-min(data)));
+
+ # --- Parse arguments --- #
+ parser = inputParser ();
+ parser.FunctionName = "nstatz";
+ parser = addParamValue (parser,'Edim', 3, @(x)x>0);
+ parser = addParamValue (parser,'Delay', 1, @(x)x>0);
+ parser = addParamValue (parser,'Points', Inf, @(x)x>0);
+ parser = addParamValue (parser,'minNeigh', 30, @(x)x>0);
+ parser = addParamValue (parser,'Rini', amp*1e-3, @(x)x>0);
+ parser = addParamValue (parser,'Rslope', 1.2, @(x)x>0);
+ parser = addParamValue (parser,'Fstep', 1, @(x)x>0);
+ parser = addParamValue (parser,'CausW', 1, @(x)x>0);
+ parser = parse(parser,varargin{:});
+
+ edim = parser.Results.Edim;
+ delay = parser.Results.Delay;
+ points = parser.Results.Points;
+ minN = parser.Results.minNeigh;
+ rini = parser.Results.Rini;
+ rslope = parser.Results.Rslope;
+ fstep = parser.Results.Fstep;
+ causw = parser.Results.CausW;
+ if ismember("CausW", parser.UsingDefaults)
+ causw = fstep;
+ end
+ clear parser
+ # ------ #
+
+ flag.num = sprintf("-#%d", n);
+
+ flag.m = sprintf("-m%d", edim);
+ flag.d = sprintf("-d%d", delay);
+ if isinf (points)
+ flag.n = "";
+ else
+ flag.n = sprintf("-d%d", points);
+ end
+ flag.k = sprintf("-k%d", minN);
+ flag.r = sprintf("-r%f", rini);
+ flag.f = sprintf("-f%f", rslope);
+ flag.s = sprintf("-s%d", fstep);
+ flag.C = sprintf("-C%d", causw);
+
+ infile = tmpnam ();
+ outfile = tmpnam ();
+
+ %% Write data to file
+ save ('-ascii',infile, 'data');
+
+ %% Prepare format of system call
+ func = file_in_loadpath ("nstat_z");
+ syscmd = sprintf("%s %s %s %s %s %s %s %s %s %s -o%s -V0 %s", func, flag.num, ...
+ flag.m, flag.d, flag.n, flag.k, flag.r, flag.f, flag.s, flag.C, ...
+ outfile, infile);
+
+ %% Function call
+ system (syscmd);
+ s = load (outfile);
+ f = s(:,1);
+ p = s(:,2);
+ ferr = s(:,3);
+
+ if nargout == 1
+ m = NA (n);
+ idx = sub2ind ([n,n], f, p);
+ m(idx) = ferr;
+ varargout{1} = m;
+ elseif nargout == 3
+ varargout = {f,p,ferr};
+ else
+ print_usage ();
+ end
+
+endfunction
Modified: trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m 2012-05-30 12:03:18 UTC (rev 10542)
+++ trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m 2012-05-30 18:37:53 UTC (rev 10543)
@@ -15,7 +15,7 @@
%% -*- texinfo -*-
%% @deftypefn {Function File} {[@var{S} @var{f} @var{AR} @var{ARerr}] = } pspec_mep (@var{data})
-%% @deftypefnx {Function File} {[@dots{}] = } pspec_mep (@dosts{},@var{property},@var{value})
+%% @deftypefnx {Function File} {[@dots{}] = } pspec_mep (@dots{},@var{property},@var{value})
%% Estimates the power spectrum of a scalar data set on the basis of the maximum @
%% entropy principle. This function calls @code{mem_spec} from the TISEAN package.
%%
Added: trunk/octave-forge/main/system-identification/inst/tisean/recurrplt.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/recurrplt.m (rev 0)
+++ trunk/octave-forge/main/system-identification/inst/tisean/recurrplt.m 2012-05-30 18:37:53 UTC (rev 10543)
@@ -0,0 +1,91 @@
+%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...>
+%%
+%% This program is free software: you can redistribute it and/or modify
+%% it under the terms of the GNU General Public License as published by
+%% the Free Software Foundation, either version 3 of the License, or
+%% any later version.
+%%
+%% This program is distributed in the hope that it will be useful,
+%% but WITHOUT ANY WARRANTY; without even the implied warranty of
+%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+%% GNU General Public License for more details.
+%%
+%% You should have received a copy of the GNU General Public License
+%% along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+%% -*- texinfo -*-
+%% @deftypefn {Function File} {[@var{M} @var{i} @var{j}] = } recurrplt (@var{data})
+%% @deftypefnx {Function File} {[@dots{}] = } recurrplt (@dots{},@var{property},@var{value})
+%% Produces a recurrence plot of the, possibly multivariate, data set. That
+%% means, for each point in the data set it looks for all points, such that the
+%% distance between these two points is smaller than a given size in a given
+%% embedding space. This function calls @code{recurr} from the TISEAN package.
+%%
+%% Returns a pairs of integers representing the indexes of the pairs of points
+%% having a distance smaller than R.
+%%
+%% Optional arguments are:
+%% @table @samp
+%% @item Edim
+%% A 1x2 vector with number of components, embedding dimension. Deafult [1,2].
+%% @item Delay
+%% Time delay. Default 1.
+%% @item R
+%% Size of the neighbourhood. Default data-interval/1000.
+%% @item Decimate
+%% Number between 0-1. Output only the percentage a of all points found
+%% should help to keep the output size reasonably smal
+%% @end table
+%%
+%% @seealso{delayvec}
+%%
+%% @end deftypefn
+
+function [i j] = recurrplt (data, varargin)
+
+ [nT nc] = size (data);
+ amp = abs((max(data(:))-min(data(:))));
+
+ # --- Parse arguments --- #
+ parser = inputParser ();
+ parser.FunctionName = "recurrplt";
+ parser = addParamValue (parser,'Edim', [1,2], @ismatrix);
+ parser = addParamValue (parser,'Delay', 1, @(x)x>0);
+ parser = addParamValue (parser,'R', amp*1e-3, @(x)x>0);
+ parser = addParamValue (parser,'Decimate', 1, @(x)x>0 && x <=1);
+ parser = parse(parser,varargin{:});
+
+ edim = parser.Results.Edim;
+ delay = parser.Results.Delay;
+ radius = parser.Results.R;
+ percent = parser.Results.Decimate;
+
+ clear parser
+ # ------ #
+
+ flag.m = sprintf("-m%d,%d", edim);
+ flag.d = sprintf("-d%d", delay);
+ flag.r = sprintf("-r%f", radius);
+ flag.P = sprintf("-%%%f", percent*100);
+ flag.c = sprintf("-c%d", nc);
+
+ infile = tmpnam ();
+ outfile = tmpnam ();
+
+ %% Write data to file
+ save ('-ascii',infile, 'data');
+
+ %% Prepare format of system call
+ func = file_in_loadpath ("recurr");
+ syscmd = sprintf("%s %s %s %s %s %s -V0 -o%s %s", ...
+ func, flag.c, flag.m, flag.d, flag.r, flag.P, outfile, infile);
+
+ %% Function call
+ system (syscmd);
+
+ %% Read results
+ d = load(outfile);
+ i = d(:,1);
+ j = d(:,2);
+
+endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 12:03:29
|
Revision: 10542
http://octave.svn.sourceforge.net/octave/?rev=10542&view=rev
Author: jpicarbajal
Date: 2012-05-30 12:03:18 +0000 (Wed, 30 May 2012)
Log Message:
-----------
system-identification: power spectrum tool. Wrapper
Modified Paths:
--------------
trunk/octave-forge/main/system-identification/INDEX
trunk/octave-forge/main/system-identification/NEWS
trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m
trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
Modified: trunk/octave-forge/main/system-identification/INDEX
===================================================================
--- trunk/octave-forge/main/system-identification/INDEX 2012-05-30 11:41:35 UTC (rev 10541)
+++ trunk/octave-forge/main/system-identification/INDEX 2012-05-30 12:03:18 UTC (rev 10542)
@@ -1,4 +1,7 @@
system-identification >> System Identification
-TISEAN
+TISEAN wrapper
+ corrdim
+ delayvec
+ pspec_mep
+ pspec_bin
Data types
-
Modified: trunk/octave-forge/main/system-identification/NEWS
===================================================================
--- trunk/octave-forge/main/system-identification/NEWS 2012-05-30 11:41:35 UTC (rev 10541)
+++ trunk/octave-forge/main/system-identification/NEWS 2012-05-30 12:03:18 UTC (rev 10542)
@@ -1,7 +1,7 @@
Summary of important user-visible changes for releases of the system-identification package
===============================================================================
-geometry-0.1.0 Release Date: XX Release Manager: Juan Pablo Carbajal
+system-identification-0.1.0 Release Date: XX Release Manager: Juan Pablo Carbajal
===============================================================================
** First official release.
Modified: trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m 2012-05-30 11:41:35 UTC (rev 10541)
+++ trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m 2012-05-30 12:03:18 UTC (rev 10542)
@@ -24,7 +24,7 @@
%%
%% @end deftypefn
-function s = pspec_bin (data, varargin)
+function [s f] = pspec_bin (data, varargin)
data = data(:);
[nT n] = size (data);
@@ -65,6 +65,7 @@
system (syscmd);
s = load (outfile);
- s = complex (s(:,1), s(:,2));
+ f = s(:,1);
+ s = s(:,2);
endfunction
Modified: trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m 2012-05-30 11:41:35 UTC (rev 10541)
+++ trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m 2012-05-30 12:03:18 UTC (rev 10542)
@@ -94,7 +94,8 @@
ARerr = str2num (strsplit (ARerrtxt, "="){2});
else
s = load (outfile);
- s = complex (s(:,1), s(:,2));
+ f = s(:,1);
+ s = s(:,2);
AR = [];
ARerr = [];
end
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 11:41:46
|
Revision: 10541
http://octave.svn.sourceforge.net/octave/?rev=10541&view=rev
Author: jpicarbajal
Date: 2012-05-30 11:41:35 +0000 (Wed, 30 May 2012)
Log Message:
-----------
system-identification: power spectrum tool. Wrapper
Added Paths:
-----------
trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m
trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
Added: trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m (rev 0)
+++ trunk/octave-forge/main/system-identification/inst/tisean/pspec_bin.m 2012-05-30 11:41:35 UTC (rev 10541)
@@ -0,0 +1,70 @@
+%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...>
+%%
+%% This program is free software: you can redistribute it and/or modify
+%% it under the terms of the GNU General Public License as published by
+%% the Free Software Foundation, either version 3 of the License, or
+%% any later version.
+%%
+%% This program is distributed in the hope that it will be useful,
+%% but WITHOUT ANY WARRANTY; without even the implied warranty of
+%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+%% GNU General Public License for more details.
+%%
+%% You should have received a copy of the GNU General Public License
+%% along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+%% -*- texinfo -*-
+%% @deftypefn {Function File} {@var{S} = } pspec_bin (@var{data})
+%% Computes a power spectrum by binning adjacent frequencies. This function @
+%% calls @code{spectrum} from the TISEAN package.
+%%
+%% @var{S} contains the power spectrum, @var{f} the frequency vector.
+%%
+%% @seealso{pspec_mep}
+%%
+%% @end deftypefn
+
+function s = pspec_bin (data, varargin)
+
+ data = data(:);
+ [nT n] = size (data);
+ if n !=1
+ warning ('Tisian:InvalidInputFormat', ...
+ 'Only first column of data will be used');
+ end
+
+ # --- Parse arguments --- #
+ parser = inputParser ();
+ parser.FunctionName = "pspec_bin";
+ parser = addParamValue (parser,'Res', 1/nT, @(x)x>0);
+ parser = addParamValue (parser,'Sampling', 1, @(x)x>0);
+ parser = parse(parser,varargin{:});
+
+ res = parser.Results.Res;
+ sampl = parser.Results.Sampling;
+
+ clear parser
+ # ------ #
+
+
+ flag.w = sprintf("-w%d", res);
+ flag.f = sprintf("-f%d", sampl);
+
+ infile = tmpnam ();
+ outfile = tmpnam ();
+
+ %% Write data to file
+ save ('-ascii',infile, 'data');
+
+ %% Prepare format of system call
+ func = file_in_loadpath ("spectrum");
+ syscmd = sprintf("%s %s %s -o%s -V0 %s", ...
+ func, flag.f, flag.w, outfile, infile);
+
+ %% Function call
+ system (syscmd);
+
+ s = load (outfile);
+ s = complex (s(:,1), s(:,2));
+
+endfunction
Added: trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m (rev 0)
+++ trunk/octave-forge/main/system-identification/inst/tisean/pspec_mep.m 2012-05-30 11:41:35 UTC (rev 10541)
@@ -0,0 +1,104 @@
+%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...>
+%%
+%% This program is free software: you can redistribute it and/or modify
+%% it under the terms of the GNU General Public License as published by
+%% the Free Software Foundation, either version 3 of the License, or
+%% any later version.
+%%
+%% This program is distributed in the hope that it will be useful,
+%% but WITHOUT ANY WARRANTY; without even the implied warranty of
+%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+%% GNU General Public License for more details.
+%%
+%% You should have received a copy of the GNU General Public License
+%% along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+%% -*- texinfo -*-
+%% @deftypefn {Function File} {[@var{S} @var{f} @var{AR} @var{ARerr}] = } pspec_mep (@var{data})
+%% @deftypefnx {Function File} {[@dots{}] = } pspec_mep (@dosts{},@var{property},@var{value})
+%% Estimates the power spectrum of a scalar data set on the basis of the maximum @
+%% entropy principle. This function calls @code{mem_spec} from the TISEAN package.
+%%
+%% @var{S} contains the power spectrum, @var{f} the frequency vector.
+%%
+%% Optional arguments are:
+%% @table @samp
+%% @item Poles
+%% Number of poles of the AR model. Deafault 128.
+%% @item Freqs
+%% Number of frequencies to use. Defualt 2000 Hz.
+%% @item Sampling
+%% Sampling frequency of the data. Default 1 Hz
+%% @item AR
+%% This is just a swtich, it doesn't need to be followed by a value. If present
+%% the coefficients of the AR model are returned in @var{AR} and the estimation
+%% error in @var{ARerr}.
+%% @end table
+%%
+%% @seealso{pspec_bin}
+%%
+%% @end deftypefn
+
+function [s f AR ARerr] = pspec_mep (data, varargin)
+
+ # --- Parse arguments --- #
+ parser = inputParser ();
+ parser.FunctionName = "pspec_mep";
+ parser = addParamValue (parser,'Poles', 128, @(x)x>0);
+ parser = addParamValue (parser,'Freqs', 2e3, @(x)x>0);
+ parser = addParamValue (parser,'Sampling', 1, @(x)x>0);
+ parser = addSwitch (parser,'AR');
+ parser = parse(parser,varargin{:});
+
+ poles = parser.Results.Poles;
+ freqs = parser.Results.Freqs;
+ sampl = parser.Results.Sampling;
+ ar = parser.Results.AR;
+
+ clear parser
+ # ------ #
+
+ data = data(:);
+ [nT n] = size (data);
+ if n !=1
+ warning ('Tisian:InvalidInputFormat', ...
+ 'Only first column of data will be used');
+ end
+
+ flag.p = sprintf("-p%d", poles);
+ flag.P = sprintf("-P%d", freqs);
+ flag.f = sprintf("-f%d", sampl);
+ flag.V = sprintf("-V%d", 2*ar);
+
+ infile = tmpnam ();
+ outfile = tmpnam ();
+
+ %% Write data to file
+ save ('-ascii',infile, 'data');
+
+ %% Prepare format of system call
+ func = file_in_loadpath ("mem_spec");
+ syscmd = sprintf("%s %s %s %s %s -o%s %s", ...
+ func, flag.p, flag.P, flag.f, flag.V, outfile, infile);
+
+ %% Function call
+ system (syscmd);
+
+ %% Read results
+ if ar
+ ARerrtxt = textread (outfile, "%s", 1){1};
+ [~, AR] = textread (outfile, "%s%f", poles, "headerlines", 1);
+ [f,s] = textread (outfile, "%f%f", freqs, "headerlines", 1+poles);
+
+ %% Parse AR error
+ ARerr = str2num (strsplit (ARerrtxt, "="){2});
+ else
+ s = load (outfile);
+ s = complex (s(:,1), s(:,2));
+ AR = [];
+ ARerr = [];
+ end
+
+
+
+endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 11:41:18
|
Revision: 10540
http://octave.svn.sourceforge.net/octave/?rev=10540&view=rev
Author: jpicarbajal
Date: 2012-05-30 11:41:08 +0000 (Wed, 30 May 2012)
Log Message:
-----------
system-identification: power spectrum tool. Wrapper
Modified Paths:
--------------
trunk/octave-forge/main/system-identification/inst/tisean/delayvec.m
Modified: trunk/octave-forge/main/system-identification/inst/tisean/delayvec.m
===================================================================
--- trunk/octave-forge/main/system-identification/inst/tisean/delayvec.m 2012-05-30 08:19:24 UTC (rev 10539)
+++ trunk/octave-forge/main/system-identification/inst/tisean/delayvec.m 2012-05-30 11:41:08 UTC (rev 10540)
@@ -53,8 +53,9 @@
save ('-ascii',infile, 'data');
%% Prepare format of the embedding vector
- syscmd = sprintf("delay -M%d -m%d %s %s -o%s -V0 %s", ...
- M, edim, flag.F, flag.D, outfile, infile);
+ func = file_in_loadpath ("delay");
+ syscmd = sprintf("%s -M%d -m%d %s %s -o%s -V0 %s", ...
+ func, M, edim, flag.F, flag.D, outfile, infile);
%% Function call
system (syscmd);
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 08:19:33
|
Revision: 10539
http://octave.svn.sourceforge.net/octave/?rev=10539&view=rev
Author: jpicarbajal
Date: 2012-05-30 08:19:24 +0000 (Wed, 30 May 2012)
Log Message:
-----------
linear-algebra: fixing text_waitbar use in nmf_pg
Modified Paths:
--------------
trunk/octave-forge/main/linear-algebra/NEWS
Modified: trunk/octave-forge/main/linear-algebra/NEWS
===================================================================
--- trunk/octave-forge/main/linear-algebra/NEWS 2012-05-30 08:15:08 UTC (rev 10538)
+++ trunk/octave-forge/main/linear-algebra/NEWS 2012-05-30 08:19:24 UTC (rev 10539)
@@ -1,3 +1,10 @@
+Summary of important user-visible changes for linear-algebra 2.2.1:
+-------------------------------------------------------------------
+
+ ** fixed bugs:
+ nmf_pg doesn't use text_waitbar by default. Miscellanoues package is
+ not required.
+
Summary of important user-visible changes for linear-algebra 2.2.0:
-------------------------------------------------------------------
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <jpi...@us...> - 2012-05-30 08:15:19
|
Revision: 10538
http://octave.svn.sourceforge.net/octave/?rev=10538&view=rev
Author: jpicarbajal
Date: 2012-05-30 08:15:08 +0000 (Wed, 30 May 2012)
Log Message:
-----------
linear-algebra: fixing text_waitbar use in nmf_pg
Modified Paths:
--------------
trunk/octave-forge/main/geometry/inst/geom2d/closed_path.m
trunk/octave-forge/main/linear-algebra/inst/nmf_pg.m
Modified: trunk/octave-forge/main/geometry/inst/geom2d/closed_path.m
===================================================================
--- trunk/octave-forge/main/geometry/inst/geom2d/closed_path.m 2012-05-30 08:00:59 UTC (rev 10537)
+++ trunk/octave-forge/main/geometry/inst/geom2d/closed_path.m 2012-05-30 08:15:08 UTC (rev 10538)
@@ -13,7 +13,6 @@
%% You should have received a copy of the GNU General Public License
%% along with this program. If not, see <http://www.gnu.org/licenses/>.
-
%% -*- texinfo -*-
%% @deftypefn {Function File} {@var{y} =} polygon (@var{x})
%% Returns a simple closed path that passes through all the points in @var{x}.
Modified: trunk/octave-forge/main/linear-algebra/inst/nmf_pg.m
===================================================================
--- trunk/octave-forge/main/linear-algebra/inst/nmf_pg.m 2012-05-30 08:00:59 UTC (rev 10537)
+++ trunk/octave-forge/main/linear-algebra/inst/nmf_pg.m 2012-05-30 08:15:08 UTC (rev 10538)
@@ -85,6 +85,12 @@
maxsubiter = parser.Results.MaxSubIter;
verbose = parser.Results.Verbose;
+ # Check if text_waitbar is loaded
+ __txtwb__ = true;
+ if !exist ('text_waitbar')
+ __txtwb__ = false;
+ end
+
clear parser
# ------ #
@@ -141,7 +147,13 @@
r,c,Wr,Wc,Hr,Hc);
fprintf ("Initial gradient norm = %f\n", initgrad);
fflush (stdout);
- text_waitbar(0,'Please wait ...');
+
+ if __txtwb__
+ text_waitbar(0,'Please wait ...');
+ else
+ printf ('Running main loop, this may take a while.\n');
+ fflush (stdout);
+ end
end
for iter = 1:maxiter
@@ -177,9 +189,10 @@
break
end
- if verbose
+ if verbose && __txtwb__
text_waitbar (iter/maxiter);
- end
+ end
+
end
if iter == maxiter
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-30 08:01:10
|
Revision: 10537
http://octave.svn.sourceforge.net/octave/?rev=10537&view=rev
Author: paramaniac
Date: 2012-05-30 08:00:59 +0000 (Wed, 30 May 2012)
Log Message:
-----------
control-devel: remove abs
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/moen4.m
Modified: trunk/octave-forge/extra/control-devel/inst/moen4.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-30 07:59:16 UTC (rev 10536)
+++ trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-30 08:00:59 UTC (rev 10537)
@@ -2825,7 +2825,7 @@
%! De = [ -0.4997 0.0451
%! -1.0011 -0.5567 ];
%!
-%!assert (abs (SYS.A), abs (Ae), 1e-4);
-%!assert (abs (SYS.B), abs (Be), 1e-4);
-%!assert (abs (SYS.C), abs (Ce), 1e-4);
-%!assert (abs (SYS.D), abs (De), 1e-4);
+%!assert (SYS.A, Ae, 1e-4);
+%!assert (SYS.B, Be, 1e-4);
+%!assert (SYS.C, Ce, 1e-4);
+%!assert (SYS.D, De, 1e-4);
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-30 07:59:29
|
Revision: 10536
http://octave.svn.sourceforge.net/octave/?rev=10536&view=rev
Author: paramaniac
Date: 2012-05-30 07:59:16 +0000 (Wed, 30 May 2012)
Log Message:
-----------
control-devel: add another test from the DaISy library
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/moen4.m
Modified: trunk/octave-forge/extra/control-devel/inst/moen4.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-30 07:28:58 UTC (rev 10535)
+++ trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-30 07:59:16 UTC (rev 10536)
@@ -744,3 +744,2088 @@
%!assert (abs (SYS.B), abs (Be), 1e-4);
%!assert (abs (SYS.C), abs (Ce), 1e-4);
%!assert (abs (SYS.D), abs (De), 1e-4);
+
+
+## [96-007] Data of a CD-player arm
+%!shared SYS, Ae, Be, Ce, De
+%! UY = [ 0.0531 -0.0313 0.0495 0.0342
+%! 0.0526 -0.0328 0.0525 0.0342
+%! 0.0504 -0.0482 0.0485 0.0567
+%! 0.0423 -0.1048 0.0285 0.1542
+%! 0.0269 -0.1515 0.0555 0.1192
+%! 0.0072 -0.1404 0.0635 0.2067
+%! -0.0091 -0.0936 0.1885 -0.0558
+%! -0.0272 -0.0770 0.1145 0.2292
+%! -0.0602 -0.0720 0.0985 0.2767
+%! -0.0962 -0.1012 0.2485 -0.0733
+%! -0.1266 -0.2732 0.2705 0.3017
+%! -0.1548 -0.4999 0.3675 0.2567
+%! -0.1831 -0.6836 0.3665 0.3492
+%! -0.2225 -0.8477 0.3145 0.7492
+%! -0.2638 -0.9533 0.4615 0.4067
+%! -0.2865 -0.9926 0.5025 0.4917
+%! -0.2877 -0.9519 0.5355 0.5942
+%! -0.2652 -0.8100 0.5835 0.2917
+%! -0.2229 -0.6546 0.4355 0.4417
+%! -0.1690 -0.4802 0.3175 0.2667
+%! -0.1022 -0.3147 0.2265 0.1067
+%! -0.0262 -0.2181 0.1655 0.1717
+%! 0.0508 -0.1214 0.1295 0.0167
+%! 0.1181 0.0027 -0.0375 0.0642
+%! 0.1639 0.1290 -0.2015 0.0767
+%! 0.1915 0.2491 -0.2445 -0.0933
+%! 0.2141 0.3501 -0.2085 -0.1833
+%! 0.2383 0.4384 -0.2055 -0.2083
+%! 0.2591 0.5325 -0.2915 -0.1533
+%! 0.2720 0.6505 -0.3725 -0.2108
+%! 0.2855 0.7959 -0.3975 -0.4083
+%! 0.3036 0.9837 -0.4585 -0.4258
+%! 0.3180 1.2078 -0.5175 -0.4108
+%! 0.3278 1.4001 -0.4415 -0.6958
+%! 0.3333 1.5244 -0.4445 -0.7183
+%! 0.3201 1.6123 -0.6065 -0.4658
+%! 0.2920 1.6983 -0.5645 -0.8133
+%! 0.2720 1.7388 -0.3875 -1.0058
+%! 0.2556 1.6854 -0.3775 -0.6208
+%! 0.2279 1.5986 -0.4115 -0.6783
+%! 0.1954 1.4945 -0.4215 -0.7808
+%! 0.1667 1.3238 -0.3475 -0.6483
+%! 0.1384 1.0917 -0.2985 -0.5033
+%! 0.1039 0.8539 -0.2995 -0.3083
+%! 0.0697 0.6475 -0.1155 -0.4483
+%! 0.0391 0.4233 -0.0315 -0.2033
+%! 0.0027 0.2235 -0.0775 0.0317
+%! -0.0330 0.0966 0.0155 -0.2458
+%! -0.0606 -0.0817 0.1125 0.0542
+%! -0.0827 -0.3000 0.2685 0.1242
+%! -0.0939 -0.4675 0.2865 0.0642
+%! -0.0950 -0.5866 0.1685 0.3717
+%! -0.0840 -0.6702 0.2735 0.1017
+%! -0.0611 -0.7676 0.2505 0.2867
+%! -0.0432 -0.7891 0.1245 0.5342
+%! -0.0289 -0.6736 0.1945 0.0842
+%! -0.0126 -0.5516 0.1885 0.2392
+%! 0.0003 -0.4055 0.1795 0.1967
+%! 0.0159 -0.2097 0.1605 -0.1208
+%! 0.0348 -0.0408 0.0925 0.0217
+%! 0.0506 0.1208 0.1395 -0.1408
+%! 0.0627 0.3078 0.0795 -0.2183
+%! 0.0723 0.4773 0.0175 -0.3158
+%! 0.0792 0.5095 0.0215 -0.3533
+%! 0.0752 0.4314 0.0085 -0.1533
+%! 0.0590 0.4053 0.0425 -0.2683
+%! 0.0378 0.4174 0.0245 -0.2658
+%! 0.0103 0.4410 0.0295 -0.1783
+%! -0.0258 0.4940 0.0905 -0.2633
+%! -0.0652 0.4837 0.1385 -0.2758
+%! -0.1035 0.3611 0.1805 -0.2258
+%! -0.1429 0.1376 0.2035 -0.0383
+%! -0.1847 -0.2251 0.3095 0.0167
+%! -0.2340 -0.7017 0.3245 0.3692
+%! -0.3036 -1.1050 0.3355 0.7517
+%! -0.3827 -1.3696 0.5995 0.4917
+%! -0.4557 -1.6423 0.7655 0.7017
+%! -0.5309 -1.9082 0.7975 1.1367
+%! -0.6018 -2.0464 0.9805 0.8317
+%! -0.6458 -2.0844 1.1445 0.7692
+%! -0.6613 -2.0949 1.2125 0.9217
+%! -0.6589 -2.0794 1.1595 0.9317
+%! -0.6399 -1.9906 1.1525 0.7942
+%! -0.6009 -1.8632 1.1895 0.6867
+%! -0.5490 -1.7523 1.0765 0.7617
+%! -0.4912 -1.6356 0.9615 0.6342
+%! -0.4272 -1.5132 0.8845 0.5817
+%! -0.3567 -1.3945 0.8295 0.5517
+%! -0.2766 -1.3215 0.7085 0.4392
+%! -0.1868 -1.3317 0.5175 0.4792
+%! -0.0971 -1.3609 0.3765 0.5317
+%! -0.0192 -1.3473 0.2125 0.6367
+%! 0.0466 -1.2796 0.0715 0.5192
+%! 0.1089 -1.2096 -0.0385 0.4117
+%! 0.1644 -1.1582 -0.1425 0.5742
+%! 0.2056 -1.0522 -0.2015 0.4542
+%! 0.2331 -0.9216 -0.3305 0.3817
+%! 0.2478 -0.8330 -0.4485 0.4142
+%! 0.2518 -0.7613 -0.4545 0.3017
+%! 0.2470 -0.7081 -0.4315 0.2892
+%! 0.2276 -0.6960 -0.4605 0.3442
+%! 0.1865 -0.7076 -0.5155 0.3692
+%! 0.1253 -0.7396 -0.4485 0.3242
+%! 0.0472 -0.8146 -0.3375 0.4292
+%! -0.0473 -0.8961 -0.2235 0.4967
+%! -0.1521 -0.9436 -0.0515 0.4692
+%! -0.2577 -0.9516 0.1265 0.5667
+%! -0.3484 -0.8765 0.3575 0.3617
+%! -0.4103 -0.7632 0.4915 0.2542
+%! -0.4523 -0.6816 0.4825 0.4667
+%! -0.4835 -0.5694 0.5725 0.2892
+%! -0.4953 -0.4093 0.6665 0.1417
+%! -0.4786 -0.2418 0.7055 0.0917
+%! -0.4345 -0.0829 0.6615 0.0192
+%! -0.3719 0.0565 0.5455 0.0617
+%! -0.2949 0.2073 0.4655 -0.0808
+%! -0.2016 0.3566 0.3255 -0.1433
+%! -0.0946 0.3950 0.1765 -0.1383
+%! 0.0180 0.2803 0.0325 -0.1083
+%! 0.1239 0.1148 -0.1245 0.0392
+%! 0.2137 -0.0357 -0.2785 0.0642
+%! 0.2851 -0.1776 -0.4545 0.1267
+%! 0.3386 -0.2838 -0.5135 0.1392
+%! 0.3724 -0.3703 -0.5395 0.1992
+%! 0.3804 -0.4262 -0.5995 0.2692
+%! 0.3606 -0.3882 -0.6205 0.1717
+%! 0.3143 -0.2657 -0.6005 0.1667
+%! 0.2456 -0.0453 -0.4925 0.0092
+%! 0.1633 0.2726 -0.3835 -0.1683
+%! 0.0757 0.6140 -0.2705 -0.3183
+%! -0.0088 0.8954 -0.1175 -0.5333
+%! -0.0827 1.0590 -0.0395 -0.5358
+%! -0.1377 1.1857 0.0545 -0.7408
+%! -0.1707 1.3715 0.0725 -0.7933
+%! -0.1911 1.5865 0.0695 -0.7758
+%! -0.2027 1.7486 0.1445 -0.9933
+%! -0.2042 1.8322 0.1265 -0.8733
+%! -0.1901 1.9362 0.2065 -1.0933
+%! -0.1550 2.0756 0.2185 -1.2008
+%! -0.1150 2.2175 0.0835 -0.9258
+%! -0.0782 2.4161 0.1225 -1.3158
+%! -0.0382 2.5968 0.0995 -1.3508
+%! 0.0009 2.6974 0.0685 -1.2508
+%! 0.0392 2.7343 0.0955 -1.4683
+%! 0.0755 2.6634 0.0165 -1.2133
+%! 0.1076 2.5477 0.0225 -1.2458
+%! 0.1431 2.4412 -0.0305 -1.2558
+%! 0.1844 2.3062 -0.1165 -1.0908
+%! 0.2260 2.1947 -0.1755 -1.1083
+%! 0.2587 2.1622 -0.2835 -0.9508
+%! 0.2798 2.1829 -0.3055 -1.0508
+%! 0.2967 2.2071 -0.3375 -1.1083
+%! 0.3114 2.2195 -0.3515 -1.0408
+%! 0.3235 2.2042 -0.3295 -1.1133
+%! 0.3315 2.0932 -0.3775 -0.9583
+%! 0.3365 1.8795 -0.3815 -0.9233
+%! 0.3446 1.6517 -0.4025 -0.8583
+%! 0.3561 1.4464 -0.4285 -0.7183
+%! 0.3619 1.2298 -0.4995 -0.5708
+%! 0.3524 1.0610 -0.6035 -0.3683
+%! 0.3355 0.9913 -0.5375 -0.5383
+%! 0.3266 0.9109 -0.4835 -0.5358
+%! 0.3218 0.7731 -0.5415 -0.2683
+%! 0.3110 0.6149 -0.5715 -0.2408
+%! 0.2969 0.4368 -0.5525 -0.1883
+%! 0.2860 0.2525 -0.4825 -0.1483
+%! 0.2740 0.0546 -0.5115 0.0767
+%! 0.2521 -0.1300 -0.5075 0.2217
+%! 0.2261 -0.2571 -0.3805 0.1317
+%! 0.2065 -0.3238 -0.3155 0.1867
+%! 0.1931 -0.2987 -0.3305 0.1942
+%! 0.1808 -0.2329 -0.3565 0.2092
+%! 0.1728 -0.2060 -0.2655 0.1117
+%! 0.1778 -0.1696 -0.1935 0.0042
+%! 0.1933 -0.0814 -0.2805 0.1292
+%! 0.2136 0.0401 -0.3195 0.0017
+%! 0.2433 0.1480 -0.3015 -0.1308
+%! 0.2830 0.1864 -0.3395 -0.1183
+%! 0.3231 0.1801 -0.5095 -0.0508
+%! 0.3548 0.1985 -0.6255 -0.0283
+%! 0.3828 0.2212 -0.5505 -0.2358
+%! 0.4071 0.2130 -0.6185 -0.0833
+%! 0.4161 0.2078 -0.7505 -0.0233
+%! 0.4070 0.1729 -0.7735 -0.1208
+%! 0.3820 0.0580 -0.7255 0.0892
+%! 0.3472 -0.1245 -0.6195 0.0492
+%! 0.3091 -0.3883 -0.6545 0.2617
+%! 0.2664 -0.6903 -0.6355 0.5492
+%! 0.2262 -1.0013 -0.4655 0.4792
+%! 0.1930 -1.3467 -0.4925 0.8692
+%! 0.1636 -1.6417 -0.5175 0.9717
+%! 0.1427 -1.8525 -0.4415 0.9067
+%! 0.1266 -2.0179 -0.4015 1.2292
+%! 0.1118 -2.1203 -0.3485 1.0842
+%! 0.1018 -2.1909 -0.3845 1.0817
+%! 0.0912 -2.2222 -0.3635 1.2492
+%! 0.0774 -2.1564 -0.2115 0.9892
+%! 0.0614 -2.0300 -0.1995 1.0467
+%! 0.0438 -1.8759 -0.1595 0.9517
+%! 0.0385 -1.7096 0.0155 0.5367
+%! 0.0457 -1.5773 0.0455 0.6742
+%! 0.0461 -1.4479 -0.0415 0.7242
+%! 0.0371 -1.2867 -0.0325 0.3842
+%! 0.0259 -1.1849 0.0445 0.4142
+%! 0.0133 -1.1751 0.1485 0.4092
+%! -0.0023 -1.1735 0.1255 0.3992
+%! -0.0234 -1.1757 0.0705 0.5492
+%! -0.0427 -1.2288 0.1825 0.3642
+%! -0.0546 -1.3379 0.1545 0.5367
+%! -0.0693 -1.4457 0.0385 0.8092
+%! -0.0862 -1.4768 0.0965 0.6067
+%! -0.0960 -1.4632 0.1635 0.6542
+%! -0.1013 -1.4416 0.1445 0.7117
+%! -0.1098 -1.3817 0.0555 0.6417
+%! -0.1307 -1.2943 0.0305 0.7067
+%! -0.1687 -1.1853 0.1585 0.6017
+%! -0.2203 -1.0350 0.2505 0.5317
+%! -0.2805 -0.8197 0.3175 0.4542
+%! -0.3427 -0.5362 0.4545 0.3067
+%! -0.3950 -0.2783 0.6245 0.1217
+%! -0.4257 -0.1265 0.6745 -0.0083
+%! -0.4347 -0.0352 0.6115 0.0642
+%! -0.4273 0.0849 0.6455 -0.0608
+%! -0.4080 0.2669 0.6375 -0.0983
+%! -0.3808 0.4593 0.5595 -0.1183
+%! -0.3378 0.6035 0.5465 -0.4158
+%! -0.2776 0.6956 0.4895 -0.3233
+%! -0.2139 0.7767 0.4145 -0.2483
+%! -0.1485 0.8759 0.3105 -0.5033
+%! -0.0798 0.9496 0.1465 -0.4033
+%! -0.0138 0.9867 0.1065 -0.4133
+%! 0.0481 1.0353 0.0825 -0.5458
+%! 0.1014 1.0584 -0.0285 -0.3958
+%! 0.1417 1.0195 -0.1245 -0.3508
+%! 0.1810 0.9643 -0.1175 -0.5058
+%! 0.2304 0.8987 -0.1325 -0.4633
+%! 0.2784 0.8223 -0.3155 -0.1833
+%! 0.3158 0.7741 -0.4165 -0.3058
+%! 0.3477 0.7241 -0.4475 -0.3608
+%! 0.3700 0.6195 -0.5515 -0.0458
+%! 0.3809 0.5139 -0.5905 -0.1633
+%! 0.3896 0.4472 -0.5885 -0.1983
+%! 0.3935 0.3697 -0.5985 0.0417
+%! 0.3879 0.3296 -0.6335 -0.0308
+%! 0.3781 0.3844 -0.6795 -0.1183
+%! 0.3661 0.4741 -0.6555 -0.0908
+%! 0.3493 0.6039 -0.5805 -0.2358
+%! 0.3240 0.7004 -0.5635 -0.2533
+%! 0.2858 0.6189 -0.5785 -0.0933
+%! 0.2432 0.4427 -0.4435 -0.2383
+%! 0.2074 0.2354 -0.3355 -0.1458
+%! 0.1672 0.0128 -0.4075 0.2967
+%! 0.1183 -0.1060 -0.3275 0.0767
+%! 0.0720 -0.1569 -0.2185 0.0792
+%! 0.0281 -0.1837 -0.1755 0.2942
+%! -0.0110 -0.1465 -0.0965 0.0092
+%! -0.0430 -0.1124 -0.1035 0.1892
+%! -0.0694 -0.0907 0.0345 0.1267
+%! -0.0822 -0.0402 0.1395 -0.0883
+%! -0.0840 0.0366 0.0545 0.1267
+%! -0.0813 0.1044 0.1125 -0.0758
+%! -0.0744 0.0936 0.1615 -0.0783
+%! -0.0736 0.0579 0.1215 0.1442
+%! -0.0809 0.0546 0.1445 -0.0458
+%! -0.0869 0.0402 0.1665 -0.0033
+%! -0.0882 0.0025 0.2585 0.0167
+%! -0.0891 -0.0893 0.2335 0.0617
+%! -0.0989 -0.2016 0.1095 0.2567
+%! -0.1157 -0.2898 0.1755 0.1517
+%! -0.1279 -0.4203 0.2505 0.2192
+%! -0.1321 -0.5523 0.2425 0.2892
+%! -0.1353 -0.6380 0.1515 0.3892
+%! -0.1451 -0.7266 0.0955 0.5492
+%! -0.1517 -0.8032 0.1775 0.3517
+%! -0.1402 -0.8444 0.1445 0.4217
+%! -0.1103 -0.8473 0.0825 0.4467
+%! -0.0677 -0.8711 0.0725 0.3342
+%! -0.0286 -0.9469 -0.0805 0.6842
+%! 0.0023 -0.9579 -0.1245 0.4242
+%! 0.0381 -0.9479 -0.1135 0.2117
+%! 0.0664 -0.9993 -0.2065 0.7217
+%! 0.0742 -1.0048 -0.1365 0.4642
+%! 0.0735 -1.0174 -0.0945 0.2817
+%! 0.0619 -1.1243 -0.1855 0.6692
+%! 0.0339 -1.2704 -0.1235 0.5067
+%! -0.0041 -1.4761 -0.0475 0.6292
+%! -0.0512 -1.6788 0.0215 0.8992
+%! -0.0973 -1.7716 0.1645 0.6467
+%! -0.1341 -1.8061 0.1835 0.7942
+%! -0.1695 -1.7935 0.1665 0.9717
+%! -0.1988 -1.7171 0.2675 0.6517
+%! -0.2107 -1.6368 0.3285 0.6417
+%! -0.2119 -1.5489 0.3375 0.7267
+%! -0.2123 -1.4181 0.3615 0.5617
+%! -0.2172 -1.2714 0.3455 0.5342
+%! -0.2319 -1.0671 0.3335 0.4817
+%! -0.2558 -0.7123 0.3665 0.2467
+%! -0.2838 -0.2669 0.4385 0.0517
+%! -0.3123 0.1531 0.5195 -0.1483
+%! -0.3390 0.5226 0.5215 -0.3008
+%! -0.3576 0.8555 0.4985 -0.4808
+%! -0.3583 1.1684 0.5345 -0.7233
+%! -0.3436 1.4304 0.5085 -0.6883
+%! -0.3233 1.6427 0.4625 -0.7483
+%! -0.2932 1.8196 0.4685 -1.0258
+%! -0.2534 1.8815 0.3835 -0.8933
+%! -0.2164 1.8872 0.2795 -0.7883
+%! -0.1816 1.9679 0.2695 -0.9983
+%! -0.1423 2.0907 0.2485 -1.0408
+%! -0.1081 2.2150 0.1715 -0.9933
+%! -0.0923 2.2823 0.0785 -0.9883
+%! -0.0926 2.2505 0.0865 -1.0533
+%! -0.0974 2.2151 0.1815 -1.1133
+%! -0.1033 2.1746 0.2105 -1.0133
+%! -0.1130 2.0369 0.1895 -0.9683
+%! -0.1261 1.8613 0.1905 -0.8858
+%! -0.1377 1.7249 0.2185 -0.8358
+%! -0.1403 1.5793 0.2285 -0.8208
+%! -0.1310 1.4183 0.2345 -0.7683
+%! -0.1180 1.2871 0.2165 -0.6133
+%! -0.1102 1.1936 0.1955 -0.5383
+%! -0.1015 1.1864 0.2075 -0.7283
+%! -0.0886 1.2592 0.1275 -0.5833
+%! -0.0749 1.3410 0.1485 -0.5658
+%! -0.0470 1.4143 0.2895 -0.9458
+%! -0.0031 1.4604 0.1465 -0.6733
+%! 0.0402 1.5447 -0.0435 -0.6083
+%! 0.0876 1.6959 -0.0645 -1.0383
+%! 0.1412 1.7869 -0.1525 -0.8083
+%! 0.1883 1.8245 -0.2255 -0.8208
+%! 0.2252 1.8679 -0.3275 -0.9933
+%! 0.2467 1.9045 -0.4395 -0.8208
+%! 0.2469 1.9533 -0.4565 -0.9058
+%! 0.2301 1.9731 -0.5135 -0.9308
+%! 0.2034 1.9264 -0.4865 -0.9233
+%! 0.1679 1.8261 -0.4125 -0.8458
+%! 0.1191 1.7087 -0.3975 -0.7083
+%! 0.0578 1.6395 -0.3765 -0.7758
+%! -0.0120 1.5785 -0.3385 -0.6708
+%! -0.0860 1.4653 -0.1925 -0.6008
+%! -0.1515 1.3092 0.0135 -0.6958
+%! -0.2008 1.0958 0.0985 -0.5583
+%! -0.2418 0.8423 0.0795 -0.3258
+%! -0.2793 0.5898 0.1115 -0.2908
+%! -0.3053 0.3063 0.2445 -0.3208
+%! -0.3232 -0.0402 0.2745 0.0467
+%! -0.3460 -0.3642 0.2835 0.3317
+%! -0.3602 -0.6486 0.4845 0.0767
+%! -0.3510 -1.0138 0.5625 0.3067
+%! -0.3365 -1.3782 0.4035 0.8092
+%! -0.3285 -1.6311 0.3775 0.7192
+%! -0.3177 -1.8760 0.4785 0.7742
+%! -0.3025 -2.0846 0.5275 1.0067
+%! -0.2861 -2.1696 0.4745 1.0042
+%! -0.2632 -2.1988 0.4385 0.9742
+%! -0.2280 -2.2098 0.4705 0.9567
+%! -0.1835 -2.2121 0.3835 1.0492
+%! -0.1343 -2.2346 0.2365 1.0167
+%! -0.0799 -2.2940 0.2005 0.9267
+%! -0.0296 -2.4018 0.1485 1.1567
+%! 0.0054 -2.5377 0.0265 1.2892
+%! 0.0379 -2.6953 -0.0205 1.1017
+%! 0.0850 -2.8762 -0.0155 1.1667
+%! 0.1365 -3.0209 -0.1045 1.4792
+%! 0.1758 -3.0638 -0.2865 1.4992
+%! 0.2069 -3.0199 -0.3975 1.3217
+%! 0.2399 -2.9711 -0.3785 1.2517
+%! 0.2687 -2.9543 -0.4125 1.3817
+%! 0.2813 -2.9447 -0.4615 1.3142
+%! 0.2737 -2.9591 -0.4575 1.2217
+%! 0.2385 -3.0112 -0.4585 1.4592
+%! 0.1778 -3.0317 -0.3535 1.3667
+%! 0.1115 -2.9655 -0.1805 1.1367
+%! 0.0469 -2.8414 -0.0535 1.2642
+%! -0.0254 -2.7009 0.0305 1.2717
+%! -0.1049 -2.5296 0.0665 1.1567
+%! -0.1749 -2.2903 0.2005 0.8842
+%! -0.2261 -2.0294 0.3485 0.7692
+%! -0.2709 -1.8409 0.3965 0.9217
+%! -0.3160 -1.7059 0.4925 0.6592
+%! -0.3553 -1.5587 0.5645 0.5542
+%! -0.3887 -1.4311 0.6065 0.6467
+%! -0.4132 -1.3691 0.6915 0.4717
+%! -0.4202 -1.3607 0.7465 0.4992
+%! -0.4079 -1.3743 0.8165 0.5242
+%! -0.3790 -1.3781 0.7885 0.4892
+%! -0.3393 -1.3442 0.6355 0.5792
+%! -0.2925 -1.2459 0.5555 0.4817
+%! -0.2378 -1.1148 0.5695 0.3467
+%! -0.1856 -0.9931 0.4825 0.4667
+%! -0.1486 -0.8185 0.3805 0.3767
+%! -0.1216 -0.5874 0.3885 0.0717
+%! -0.1015 -0.3821 0.3235 0.1667
+%! -0.0898 -0.1778 0.2855 0.0667
+%! -0.0712 0.0346 0.3355 -0.2858
+%! -0.0379 0.2085 0.3175 -0.2483
+%! -0.0027 0.3772 0.2065 -0.2333
+%! 0.0236 0.5655 0.0175 -0.3458
+%! 0.0438 0.7329 -0.0445 -0.4733
+%! 0.0635 0.8253 0.0225 -0.5433
+%! 0.0801 0.8649 0.0045 -0.4633
+%! 0.0903 0.9308 -0.0115 -0.6133
+%! 0.0901 1.0208 -0.0875 -0.5258
+%! 0.0744 1.1279 -0.1185 -0.4933
+%! 0.0558 1.2299 -0.0115 -0.8233
+%! 0.0408 1.2587 -0.0005 -0.6233
+%! 0.0191 1.2574 0.0385 -0.5433
+%! -0.0063 1.2484 0.0785 -0.7158
+%! -0.0255 1.1921 0.0205 -0.5933
+%! -0.0336 1.1108 0.0425 -0.6808
+%! -0.0357 0.9984 0.0135 -0.4883
+%! -0.0415 0.8908 0.0155 -0.3583
+%! -0.0416 0.8448 0.0965 -0.6158
+%! -0.0291 0.8026 0.0265 -0.4233
+%! -0.0135 0.7560 -0.0525 -0.2958
+%! 0.0064 0.7406 -0.0195 -0.4958
+%! 0.0347 0.7367 -0.0165 -0.4208
+%! 0.0633 0.7939 -0.0945 -0.3233
+%! 0.0912 0.9579 -0.1195 -0.5658
+%! 0.1225 1.1168 -0.1215 -0.6283
+%! 0.1480 1.1866 -0.2055 -0.4508
+%! 0.1662 1.2186 -0.1875 -0.7083
+%! 0.1851 1.2276 -0.2065 -0.7058
+%! 0.1996 1.2388 -0.2605 -0.4633
+%! 0.2131 1.3063 -0.2195 -0.7858
+%! 0.2323 1.3793 -0.3485 -0.6333
+%! 0.2558 1.4242 -0.3975 -0.6083
+%! 0.2937 1.4094 -0.3555 -0.9133
+%! 0.3392 1.2810 -0.5475 -0.4583
+%! 0.3757 1.1218 -0.6745 -0.4483
+%! 0.4103 0.9746 -0.7225 -0.6508
+%! 0.4441 0.8090 -0.8085 -0.2808
+%! 0.4710 0.6612 -0.8195 -0.3258
+%! 0.4931 0.5291 -0.8975 -0.3208
+%! 0.5052 0.4500 -0.9735 -0.1608
+%! 0.4969 0.4825 -1.0065 -0.2058
+%! 0.4640 0.5571 -1.0395 -0.1783
+%! 0.4134 0.6007 -0.9225 -0.3358
+%! 0.3525 0.5695 -0.7705 -0.2733
+%! 0.2808 0.4708 -0.6625 -0.1683
+%! 0.2038 0.3850 -0.5825 -0.2383
+%! 0.1330 0.3329 -0.4905 -0.1808
+%! 0.0800 0.2964 -0.3315 -0.2783
+%! 0.0499 0.2154 -0.2585 -0.2158
+%! 0.0356 0.0501 -0.2945 0.0017
+%! 0.0345 -0.1225 -0.2425 -0.1633
+%! 0.0424 -0.2963 -0.2555 0.0417
+%! 0.0454 -0.4329 -0.3185 0.2467
+%! 0.0459 -0.4041 -0.2455 -0.0408
+%! 0.0501 -0.2796 -0.1785 0.0867
+%! 0.0562 -0.1394 -0.1245 0.0392
+%! 0.0708 0.0026 -0.1095 -0.2233
+%! 0.0921 0.1000 -0.2235 0.0242
+%! 0.1173 0.2755 -0.1985 -0.2208
+%! 0.1574 0.5562 -0.1505 -0.5558
+%! 0.2034 0.7663 -0.3005 -0.3233
+%! 0.2345 0.8696 -0.3995 -0.4508
+%! 0.2466 0.8681 -0.4305 -0.5183
+%! 0.2426 0.8135 -0.4365 -0.4033
+%! 0.2294 0.8375 -0.3955 -0.5383
+%! 0.2117 0.9195 -0.3865 -0.4633
+%! 0.1889 1.0417 -0.3415 -0.5258
+%! 0.1615 1.2338 -0.3415 -0.6758
+%! 0.1260 1.4168 -0.4065 -0.6508
+%! 0.0810 1.5369 -0.3505 -0.7708
+%! 0.0342 1.5810 -0.1615 -0.8558
+%! -0.0091 1.5314 -0.0225 -0.8083
+%! -0.0503 1.4241 -0.0285 -0.7283
+%! -0.0876 1.2739 0.0045 -0.7208
+%! -0.1141 1.0825 0.0965 -0.7008
+%! -0.1308 0.9130 0.1345 -0.5608
+%! -0.1446 0.7945 0.1385 -0.5108
+%! -0.1608 0.7183 0.1425 -0.4833
+%! -0.1846 0.7185 0.1855 -0.4508
+%! -0.2169 0.7908 0.2355 -0.5133
+%! -0.2502 0.8577 0.2985 -0.5833
+%! -0.2747 0.8571 0.4505 -0.6633
+%! -0.2911 0.7965 0.5145 -0.5358
+%! -0.3055 0.7070 0.5125 -0.5008
+%! -0.3153 0.5460 0.5205 -0.5583
+%! -0.3242 0.2987 0.4875 -0.2383
+%! -0.3392 0.0706 0.5795 -0.1583
+%! -0.3528 -0.0710 0.6875 -0.2433
+%! -0.3649 -0.1759 0.6485 0.0217
+%! -0.3827 -0.3188 0.6335 0.1142
+%! -0.3987 -0.5033 0.7005 0.0142
+%! -0.4061 -0.6932 0.7185 0.2192
+%! -0.4107 -0.8736 0.7015 0.4317
+%! -0.4089 -1.0089 0.7395 0.3192
+%! -0.3942 -1.1045 0.7075 0.4167
+%! -0.3738 -1.2074 0.5835 0.6192
+%! -0.3490 -1.3199 0.5075 0.5067
+%! -0.3125 -1.4543 0.4785 0.5492
+%! -0.2685 -1.5904 0.4425 0.6667
+%! -0.2324 -1.6779 0.2805 0.7667
+%! -0.2177 -1.7222 0.1215 0.8817
+%! -0.2177 -1.7124 0.2035 0.6517
+%! -0.2191 -1.6440 0.2945 0.6017
+%! -0.2262 -1.5442 0.2725 0.7342
+%! -0.2446 -1.4283 0.2725 0.5967
+%! -0.2675 -1.3620 0.3145 0.5392
+%! -0.2852 -1.3847 0.3815 0.5242
+%! -0.2915 -1.4221 0.3885 0.5442
+%! -0.2859 -1.4241 0.3955 0.6142
+%! -0.2662 -1.4117 0.4345 0.5817
+%! -0.2261 -1.4121 0.3795 0.5092
+%! -0.1681 -1.4674 0.1875 0.5967
+%! -0.1067 -1.5856 0.0025 0.7967
+%! -0.0484 -1.7123 0.0105 0.7017
+%! 0.0055 -1.8149 -0.0275 0.7892
+%! 0.0441 -1.8852 -0.1755 1.0267
+%! 0.0662 -1.9856 -0.1775 0.8417
+%! 0.0769 -2.2193 -0.1235 1.0092
+%! 0.0690 -2.4788 -0.1095 1.2917
+%! 0.0433 -2.5547 -0.1075 1.1892
+%! 0.0121 -2.4502 -0.0515 1.1292
+%! -0.0163 -2.3174 0.0695 1.0167
+%! -0.0431 -2.2535 0.0355 1.0567
+%! -0.0690 -2.2237 -0.0085 1.0267
+%! -0.0830 -2.1826 0.0905 0.8217
+%! -0.0802 -2.1610 0.1245 0.9717
+%! -0.0644 -2.1314 0.0895 0.8942
+%! -0.0389 -2.0083 -0.0235 0.7767
+%! -0.0144 -1.7807 -0.0925 0.8117
+%! 0.0039 -1.5045 -0.0555 0.5167
+%! 0.0181 -1.2185 -0.0975 0.3967
+%! 0.0257 -0.9135 -0.1055 0.2917
+%! 0.0253 -0.5901 -0.0405 0.0192
+%! 0.0119 -0.2804 -0.0405 -0.0208
+%! -0.0223 0.0284 -0.0525 -0.1183
+%! -0.0722 0.3668 0.0375 -0.3658
+%! -0.1248 0.6931 0.2185 -0.5433
+%! -0.1762 0.9826 0.3255 -0.5983
+%! -0.2279 1.2387 0.3505 -0.7408
+%! -0.2801 1.4081 0.3595 -0.7783
+%! -0.3313 1.4824 0.4385 -0.7633
+%! -0.3714 1.4586 0.5705 -0.8483
+%! -0.3940 1.3294 0.5915 -0.6908
+%! -0.4058 1.2356 0.5915 -0.5733
+%! -0.4124 1.2560 0.5965 -0.6108
+%! -0.4123 1.3012 0.5405 -0.5908
+%! -0.3988 1.2608 0.5535 -0.6983
+%! -0.3752 1.0714 0.5255 -0.4733
+%! -0.3563 0.8491 0.5215 -0.2658
+%! -0.3428 0.6776 0.5805 -0.4083
+%! -0.3337 0.4690 0.4925 -0.1158
+%! -0.3337 0.2571 0.5005 0.0117
+%! -0.3316 0.0998 0.6345 -0.1208
+%! -0.3140 -0.0275 0.7025 -0.0458
+%! -0.2848 -0.1338 0.6455 0.0317
+%! -0.2575 -0.2528 0.4625 0.1967
+%! -0.2342 -0.3685 0.4365 0.0917
+%! -0.2083 -0.4770 0.4665 0.0767
+%! -0.1875 -0.6038 0.3695 0.3942
+%! -0.1767 -0.6891 0.3645 0.2767
+%! -0.1647 -0.7274 0.4225 0.1342
+%! -0.1496 -0.7568 0.3685 0.2642
+%! -0.1427 -0.8146 0.2315 0.3767
+%! -0.1438 -0.9273 0.2285 0.3517
+%! -0.1384 -1.0284 0.4105 0.2117
+%! -0.1249 -1.1491 0.4185 0.4717
+%! -0.1147 -1.3249 0.2665 0.6917
+%! -0.0989 -1.4528 0.2795 0.4867
+%! -0.0619 -1.5434 0.3525 0.5667
+%! -0.0081 -1.6087 0.2615 0.7242
+%! 0.0507 -1.6067 0.0855 0.6692
+%! 0.1051 -1.5640 -0.0405 0.7317
+%! 0.1482 -1.4225 -0.0795 0.6492
+%! 0.1804 -1.1516 -0.1505 0.4267
+%! 0.2027 -0.8082 -0.2635 0.3492
+%! 0.2135 -0.4165 -0.2375 0.1192
+%! 0.2148 -0.0383 -0.1745 -0.1108
+%! 0.2076 0.3144 -0.3045 -0.1433
+%! 0.2004 0.6244 -0.3415 -0.4108
+%! 0.2110 0.7909 -0.2105 -0.6583
+%! 0.2361 0.8213 -0.2295 -0.3833
+%! 0.2602 0.7620 -0.3495 -0.3533
+%! 0.2841 0.6144 -0.4395 -0.4783
+%! 0.3056 0.4209 -0.4835 -0.1683
+%! 0.3176 0.2084 -0.4535 -0.0783
+%! 0.3250 -0.0407 -0.4695 -0.0783
+%! 0.3287 -0.3328 -0.5005 0.1792
+%! 0.3187 -0.6266 -0.5085 0.4117
+%! 0.2961 -0.8894 -0.5195 0.4617
+%! 0.2740 -1.0959 -0.4995 0.4492
+%! 0.2542 -1.1684 -0.4795 0.6467
+%! 0.2321 -1.0800 -0.3895 0.5542
+%! 0.2110 -0.8808 -0.3495 0.3067
+%! 0.1900 -0.6415 -0.3965 0.3167
+%! 0.1682 -0.4753 -0.3005 0.0467
+%! 0.1458 -0.4414 -0.2055 0.0917
+%! 0.1163 -0.4386 -0.1425 0.2617
+%! 0.0891 -0.3733 0.0045 -0.0933
+%! 0.0760 -0.3217 0.0415 -0.0408
+%! 0.0669 -0.3010 -0.0135 0.1467
+%! 0.0560 -0.2167 -0.0305 -0.0708
+%! 0.0476 -0.1060 -0.0405 -0.0458
+%! 0.0458 -0.0230 0.0035 -0.1208
+%! 0.0540 0.0155 -0.0085 -0.1733
+%! 0.0668 -0.0147 -0.1005 0.0117
+%! 0.0792 -0.0398 -0.1425 -0.0508
+%! 0.0933 -0.0347 -0.1735 -0.0508
+%! 0.1089 -0.0374 -0.2045 -0.0283
+%! 0.1195 -0.0380 -0.2625 0.0292
+%! 0.1195 -0.0179 -0.2725 0.0267
+%! 0.1110 0.0673 -0.2335 -0.1483
+%! 0.0916 0.2159 -0.2865 -0.0383
+%! 0.0596 0.3960 -0.2405 -0.1733
+%! 0.0267 0.6076 -0.1025 -0.4333
+%! -0.0017 0.8519 -0.0395 -0.3883
+%! -0.0248 1.1211 0.0015 -0.6533
+%! -0.0414 1.3389 -0.0415 -0.6833
+%! -0.0505 1.4805 -0.0145 -0.7208
+%! -0.0418 1.5473 0.0785 -1.0008
+%! -0.0182 1.5086 -0.0465 -0.7158
+%! 0.0047 1.4786 -0.1455 -0.6858
+%! 0.0286 1.4618 -0.1135 -0.8933
+%! 0.0544 1.3645 -0.1525 -0.6733
+%! 0.0735 1.2553 -0.1945 -0.6583
+%! 0.0804 1.1186 -0.2285 -0.5433
+%! 0.0768 0.9336 -0.1345 -0.4483
+%! 0.0700 0.7648 -0.0355 -0.4333
+%! 0.0642 0.5964 -0.0615 -0.2583
+%! 0.0608 0.4208 -0.0495 -0.2258
+%! 0.0602 0.2307 -0.0155 -0.0408
+%! 0.0659 0.0514 0.0435 -0.0433
+%! 0.0820 -0.0945 0.0105 -0.0583
+%! 0.0956 -0.2256 -0.1235 0.2667
+%! 0.0948 -0.2807 -0.1055 0.1717
+%! 0.0868 -0.2945 -0.0515 0.0367
+%! 0.0754 -0.3494 -0.0825 0.1942
+%! 0.0574 -0.4033 -0.0815 0.1867
+%! 0.0357 -0.4587 -0.0115 0.1717
+%! 0.0156 -0.5346 0.0505 0.1867
+%! -0.0023 -0.6237 0.0135 0.3017
+%! -0.0165 -0.7228 0.0525 0.2717
+%! -0.0230 -0.7928 0.1605 0.2417
+%! -0.0284 -0.8219 0.1615 0.4642
+%! -0.0327 -0.8036 0.1625 0.2917
+%! -0.0210 -0.7466 0.1915 0.1567
+%! 0.0111 -0.7242 0.2005 0.3042
+%! 0.0592 -0.7573 0.1465 0.2667
+%! 0.1231 -0.7844 -0.0095 0.2467
+%! 0.1999 -0.7794 -0.1375 0.2692
+%! 0.2804 -0.7513 -0.2825 0.3117
+%! 0.3560 -0.6774 -0.4315 0.2467
+%! 0.4198 -0.5854 -0.5625 0.2242
+%! 0.4650 -0.5277 -0.6685 0.3242
+%! 0.4975 -0.4888 -0.6675 0.1642
+%! 0.5329 -0.4856 -0.6705 0.0967
+%! 0.5714 -0.5286 -0.7535 0.2292
+%! 0.5974 -0.5528 -0.8835 0.3717
+%! 0.6036 -0.5237 -0.9285 0.3542
+%! 0.6018 -0.4513 -0.8545 0.1417
+%! 0.5995 -0.3414 -0.8395 0.2067
+%! 0.5916 -0.2010 -0.8205 0.1942
+%! 0.5754 -0.0626 -0.8065 0.0892
+%! 0.5542 0.0470 -0.8155 0.0967
+%! 0.5375 0.0963 -0.7775 -0.0808
+%! 0.5306 0.0637 -0.8045 0.0367
+%! 0.5296 0.0065 -0.7795 0.0492
+%! 0.5319 -0.0821 -0.7785 -0.0208
+%! 0.5239 -0.2330 -0.9285 0.2767
+%! 0.4911 -0.3534 -0.9875 0.2492
+%! 0.4368 -0.4041 -0.9135 0.2767
+%! 0.3774 -0.4396 -0.6965 0.2317
+%! 0.3299 -0.5189 -0.5085 0.1617
+%! 0.2946 -0.5876 -0.5285 0.4617
+%! 0.2732 -0.5118 -0.4495 0.1917
+%! 0.2735 -0.3676 -0.4335 0.1092
+%! 0.2897 -0.2488 -0.4875 0.2317
+%! 0.3188 -0.1107 -0.4445 -0.0758
+%! 0.3580 0.0310 -0.4655 0.0367
+%! 0.3999 0.2099 -0.4825 -0.0933
+%! 0.4416 0.4537 -0.6305 -0.2708
+%! 0.4758 0.6744 -0.8075 -0.1533
+%! 0.5033 0.8167 -0.7235 -0.5058
+%! 0.5249 0.8812 -0.7385 -0.3783
+%! 0.5302 0.9226 -0.8115 -0.2858
+%! 0.5213 0.9572 -0.8195 -0.5783
+%! 0.4967 0.9528 -0.8725 -0.2483
+%! 0.4535 0.9563 -0.7685 -0.3808
+%! 0.4022 0.9923 -0.6685 -0.5258
+%! 0.3435 1.0462 -0.6205 -0.2433
+%! 0.2816 1.1181 -0.4365 -0.5508
+%! 0.2276 1.1811 -0.3685 -0.5858
+%! 0.1752 1.2622 -0.3765 -0.4508
+%! 0.1188 1.4145 -0.3025 -0.6808
+%! 0.0624 1.5738 -0.1975 -0.7083
+%! 0.0118 1.6516 -0.0645 -0.8358
+%! -0.0333 1.6513 -0.0625 -0.8383
+%! -0.0826 1.5769 -0.1065 -0.6158
+%! -0.1318 1.4702 0.0125 -0.7833
+%! -0.1670 1.4121 0.1115 -0.7908
+%! -0.1917 1.3594 0.0985 -0.5433
+%! -0.2125 1.2994 0.1295 -0.6108
+%! -0.2249 1.2767 0.2415 -0.6833
+%! -0.2319 1.2354 0.2795 -0.5283
+%! -0.2430 1.1914 0.2665 -0.4958
+%! -0.2598 1.1599 0.3615 -0.5608
+%! -0.2818 1.0552 0.4915 -0.4008
+%! -0.3076 0.8983 0.5765 -0.2833
+%! -0.3210 0.7477 0.6535 -0.4083
+%! -0.3101 0.6146 0.6765 -0.3383
+%! -0.2900 0.5163 0.5625 0.0067
+%! -0.2723 0.4557 0.4575 -0.0683
+%! -0.2465 0.3869 0.4035 -0.2233
+%! -0.2099 0.2499 0.3495 0.0092
+%! -0.1739 0.1037 0.3055 0.0267
+%! -0.1512 0.0119 0.1695 0.1517
+%! -0.1484 -0.0519 0.0805 0.2842
+%! -0.1510 -0.0699 0.2045 -0.0533
+%! -0.1495 -0.0575 0.2485 0.0517
+%! -0.1572 -0.0494 0.2205 0.2867
+%! -0.1715 -0.0612 0.3355 -0.0658
+%! -0.1810 -0.1391 0.3715 0.0442
+%! -0.1916 -0.2317 0.3255 0.2367
+%! -0.2045 -0.3066 0.3375 0.0842
+%! -0.2152 -0.4015 0.3625 0.2167
+%! -0.2228 -0.5104 0.4315 0.2292
+%! -0.2264 -0.6669 0.3695 0.2817
+%! -0.2322 -0.8304 0.2035 0.5217
+%! -0.2424 -0.9202 0.2115 0.4992
+%! -0.2438 -0.9549 0.3045 0.3742
+%! -0.2288 -0.9232 0.2735 0.3792
+%! -0.2076 -0.8649 0.1505 0.5142
+%! -0.1881 -0.8226 0.1345 0.4142
+%! -0.1670 -0.7409 0.1615 0.2617
+%! -0.1451 -0.6425 0.1015 0.3617
+%! -0.1230 -0.5638 0.1275 0.1917
+%! -0.0984 -0.5363 0.2195 0.1492
+%! -0.0803 -0.5953 0.1945 0.3917
+%! -0.0729 -0.6555 0.1085 0.3267
+%! -0.0632 -0.6865 0.1145 0.1717
+%! -0.0489 -0.7249 0.1485 0.3442
+%! -0.0439 -0.7028 0.0995 0.4717
+%! -0.0451 -0.6042 0.0625 0.2117
+%! -0.0406 -0.5269 0.0245 0.1367
+%! -0.0383 -0.4894 -0.0435 0.3392
+%! -0.0485 -0.4438 -0.0455 0.2067
+%! -0.0679 -0.4331 -0.0395 0.2042
+%! -0.0895 -0.4768 0.0445 0.2217
+%! -0.1077 -0.4958 0.1265 0.0817
+%! -0.1329 -0.4951 0.0275 0.3717
+%! -0.1737 -0.5221 0.0445 0.3217
+%! -0.2132 -0.5739 0.2165 0.0742
+%! -0.2403 -0.6589 0.2895 0.3267
+%! -0.2575 -0.7592 0.3165 0.3442
+%! -0.2639 -0.8295 0.2885 0.3417
+%! -0.2594 -0.8865 0.2775 0.3992
+%! -0.2440 -0.9962 0.2605 0.3367
+%! -0.2267 -1.1301 0.0985 0.6167
+%! -0.2183 -1.1910 0.0655 0.6417
+%! -0.2123 -1.1804 0.1735 0.4192
+%! -0.2063 -1.1308 0.1145 0.5367
+%! -0.2103 -1.0522 -0.0135 0.6292
+%! -0.2177 -0.9723 0.0595 0.3867
+%! -0.2085 -0.8974 0.2295 0.2317
+%! -0.1816 -0.8520 0.1905 0.4017
+%! -0.1471 -0.8541 0.0715 0.3917
+%! -0.1048 -0.8622 0.0615 0.2592
+%! -0.0576 -0.8463 0.0615 0.3767
+%! -0.0136 -0.8192 -0.0185 0.3392
+%! 0.0235 -0.7788 -0.1415 0.3292
+%! 0.0551 -0.7021 -0.1095 0.2867
+%! 0.0856 -0.6467 -0.0675 0.1667
+%! 0.1120 -0.6449 -0.2705 0.3217
+%! 0.1297 -0.6291 -0.4115 0.2392
+%! 0.1433 -0.6098 -0.3755 0.1167
+%! 0.1493 -0.5649 -0.4015 0.2967
+%! 0.1413 -0.4200 -0.4795 0.1742
+%! 0.1251 -0.2349 -0.5045 -0.0333
+%! 0.1035 -0.0723 -0.4195 -0.0283
+%! 0.0714 0.0871 -0.3485 -0.0533
+%! 0.0300 0.2866 -0.3345 -0.1883
+%! -0.0115 0.5233 -0.2235 -0.3508
+%! -0.0474 0.7916 -0.0725 -0.4258
+%! -0.0765 1.1562 -0.0105 -0.6133
+%! -0.0981 1.5229 -0.0365 -0.7958
+%! -0.1162 1.7436 -0.0585 -0.7983
+%! -0.1316 1.8822 0.0025 -0.9008
+%! -0.1349 1.9916 0.0605 -1.0358
+%! -0.1202 2.0391 0.0065 -1.0233
+%! -0.0932 2.0123 -0.0625 -0.9483
+%! -0.0609 1.9210 -0.0805 -0.9583
+%! -0.0293 1.8059 -0.1785 -0.8458
+%! -0.0052 1.7253 -0.2625 -0.7733
+%! 0.0091 1.7317 -0.2195 -0.8883
+%! 0.0063 1.7910 -0.1925 -0.7333
+%! -0.0175 1.8525 -0.1535 -0.7433
+%! -0.0460 1.9114 -0.0405 -1.0208
+%! -0.0656 1.9435 0.0365 -0.9383
+%! -0.0793 1.9967 0.1075 -0.9433
+%! -0.0890 2.0655 0.1275 -1.0708
+%! -0.0971 2.0547 0.0985 -0.9808
+%! -0.1056 1.9936 0.1425 -1.0508
+%! -0.1133 1.9289 0.1445 -1.0233
+%! -0.1202 1.8841 0.1445 -1.0283
+%! -0.1278 1.8653 0.1605 -1.1133
+%! -0.1456 1.8273 0.0905 -0.9433
+%! -0.1796 1.7730 0.0805 -0.9833
+%! -0.2191 1.6584 0.1575 -1.0358
+%! -0.2544 1.4743 0.2565 -0.9083
+%! -0.2854 1.3089 0.3565 -0.9083
+%! -0.3208 1.1494 0.3295 -0.6758
+%! -0.3642 0.9927 0.3655 -0.6683
+%! -0.4052 0.8263 0.5005 -0.7658
+%! -0.4443 0.5780 0.5355 -0.4133
+%! -0.4891 0.3103 0.6355 -0.3558
+%! -0.5373 0.1107 0.7275 -0.2808
+%! -0.5849 0.0334 0.7645 -0.1008
+%! -0.6212 0.1017 0.8415 -0.3308
+%! -0.6375 0.2141 0.8055 -0.2608
+%! -0.6349 0.3187 0.8195 -0.2683
+%! -0.6100 0.4173 0.8785 -0.5558
+%! -0.5747 0.4404 0.6735 -0.2058
+%! -0.5453 0.3679 0.5525 -0.1733
+%! -0.5116 0.2496 0.6265 -0.4183
+%! -0.4674 0.1314 0.6485 -0.1408
+%! -0.4214 0.0316 0.6365 -0.1358
+%! -0.3811 -0.0349 0.5035 -0.0783
+%! -0.3498 -0.0504 0.4445 0.0092
+%! -0.3165 -0.0772 0.5055 -0.2258
+%! -0.2736 -0.1911 0.4325 -0.0983
+%! -0.2333 -0.3845 0.3495 0.0917
+%! -0.2116 -0.6080 0.2855 0.1917
+%! -0.2140 -0.7873 0.2095 0.3642
+%! -0.2278 -0.9055 0.1795 0.2642
+%! -0.2368 -1.0029 0.1665 0.3317
+%! -0.2350 -1.0601 0.2255 0.3617
+%! -0.2232 -1.0793 0.2145 0.2592
+%! -0.2166 -1.0744 -0.0305 0.4942
+%! -0.2288 -0.9663 -0.1795 0.4292
+%! -0.2482 -0.7630 -0.0595 0.1117
+%! -0.2598 -0.5967 0.1215 0.0967
+%! -0.2618 -0.4925 0.2025 0.0492
+%! -0.2590 -0.4010 0.1605 0.0367
+%! -0.2550 -0.3328 0.1185 0.0867
+%! -0.2431 -0.2946 0.1605 -0.0808
+%! -0.2166 -0.2787 0.1645 -0.0383
+%! -0.1817 -0.3122 0.1755 0.1067
+%! -0.1408 -0.4076 0.2485 -0.0058
+%! -0.0955 -0.5121 0.1335 0.1192
+%! -0.0580 -0.5406 -0.0965 0.2292
+%! -0.0371 -0.4996 -0.1945 0.1292
+%! -0.0357 -0.4631 -0.1495 0.2142
+%! -0.0497 -0.4436 -0.0065 0.0717
+%! -0.0737 -0.4519 -0.0165 0.0992
+%! -0.1119 -0.4687 -0.0865 0.2892
+%! -0.1561 -0.4543 0.0375 -0.0208
+%! -0.1950 -0.4876 0.1145 0.0817
+%! -0.2343 -0.5827 0.1375 0.2967
+%! -0.2695 -0.6924 0.2895 0.0017
+%! -0.2983 -0.8815 0.3295 0.3067
+%! -0.3325 -1.1117 0.3055 0.5492
+%! -0.3656 -1.2956 0.3705 0.2967
+%! -0.3906 -1.4858 0.4285 0.5417
+%! -0.4198 -1.6773 0.4905 0.8242
+%! -0.4544 -1.8112 0.5555 0.6917
+%! -0.4817 -1.9092 0.5635 0.7317
+%! -0.4984 -1.9624 0.5865 0.8517
+%! -0.5052 -1.9578 0.6725 0.7392
+%! -0.5043 -1.9750 0.6675 0.8467
+%! -0.5005 -2.0244 0.6645 0.8942
+%! -0.4895 -2.0492 0.7325 0.7292
+%! -0.4713 -2.0865 0.6905 0.9517
+%! -0.4485 -2.1329 0.7045 0.9417
+%! -0.4108 -2.1318 0.7525 0.7467
+%! -0.3562 -2.1064 0.6665 0.9792
+%! -0.2897 -2.0510 0.5845 0.8717
+%! -0.2097 -1.9414 0.4475 0.7142
+%! -0.1240 -1.8255 0.2975 0.8092
+%! -0.0458 -1.7261 0.1825 0.7217
+%! 0.0176 -1.5932 0.0225 0.7017
+%! 0.0662 -1.4228 -0.0535 0.5917
+%! 0.1049 -1.3030 -0.0845 0.4792
+%! 0.1335 -1.2418 -0.1745 0.5192
+%! 0.1483 -1.1914 -0.2385 0.4542
+%! 0.1466 -1.1992 -0.2885 0.5092
+%! 0.1282 -1.2520 -0.2755 0.5542
+%! 0.0992 -1.2430 -0.2205 0.5267
+%! 0.0699 -1.1510 -0.1345 0.4517
+%! 0.0478 -1.0541 -0.0455 0.3867
+%! 0.0322 -0.9852 -0.0195 0.4267
+%! 0.0207 -0.8964 0.0005 0.2817
+%! 0.0094 -0.8166 -0.0155 0.3042
+%! -0.0108 -0.7687 0.0055 0.3492
+%! -0.0421 -0.6628 0.0855 0.1892
+%! -0.0816 -0.5166 0.1335 0.2067
+%! -0.1240 -0.4531 0.1985 0.1517
+%! -0.1587 -0.4955 0.2175 0.1317
+%! -0.1765 -0.5726 0.2375 0.2242
+%! -0.1697 -0.6386 0.2575 0.1567
+%! -0.1396 -0.6945 0.1815 0.2542
+%! -0.0986 -0.7609 0.0965 0.3417
+%! -0.0577 -0.8911 0.0185 0.4417
+%! -0.0179 -1.0590 0.0145 0.4917
+%! 0.0260 -1.1811 0.0785 0.4217
+%! 0.0683 -1.2444 0.0305 0.6217
+%! 0.0949 -1.1920 -0.0365 0.6392
+%! 0.1073 -1.0508 -0.0215 0.4067
+%! 0.1155 -0.9975 -0.0025 0.4092
+%! 0.1214 -1.0199 0.0115 0.4092
+%! 0.1192 -0.9884 0.0035 0.4167
+%! 0.1011 -0.8955 -0.0175 0.4642
+%! 0.0706 -0.7734 -0.0055 0.3317
+%! 0.0435 -0.6574 0.0385 0.1417
+%! 0.0254 -0.5681 0.0495 0.1792
+%! 0.0072 -0.4466 0.0325 0.2567
+%! -0.0080 -0.2535 0.0875 -0.0533
+%! -0.0115 -0.0873 0.0575 -0.1208
+%! -0.0112 -0.0130 -0.0395 0.0817
+%! -0.0089 0.0369 0.0145 -0.1758
+%! 0.0030 0.0533 0.0625 -0.1708
+%! 0.0194 -0.0025 0.0185 0.0467
+%! 0.0374 -0.0338 -0.0445 -0.0783
+%! 0.0607 0.0134 -0.0985 -0.0633
+%! 0.0869 0.0666 -0.1175 -0.0358
+%! 0.1116 0.1083 -0.1945 -0.0683
+%! 0.1329 0.1776 -0.2645 -0.0883
+%! 0.1515 0.2269 -0.2495 -0.1808
+%! 0.1643 0.2195 -0.3015 -0.0533
+%! 0.1677 0.2351 -0.3745 -0.0333
+%! 0.1713 0.2860 -0.3445 -0.2508
+%! 0.1813 0.2769 -0.3195 -0.1283
+%! 0.1898 0.1852 -0.3195 -0.0083
+%! 0.1937 0.0800 -0.3415 -0.0683
+%! 0.1960 0.0125 -0.3745 0.0767
+%! 0.2058 -0.0679 -0.2675 -0.0358
+%! 0.2313 -0.2067 -0.2205 -0.0433
+%! 0.2604 -0.2975 -0.3315 0.2442
+%! 0.2806 -0.2723 -0.3645 0.0692
+%! 0.2919 -0.2029 -0.3955 0.0117
+%! 0.2899 -0.0986 -0.4325 0.0517
+%! 0.2731 0.0407 -0.3855 -0.1783
+%! 0.2420 0.1252 -0.3585 -0.1008
+%! 0.1990 0.1201 -0.2645 -0.1433
+%! 0.1551 0.0612 -0.2165 -0.2383
+%! 0.1148 -0.0002 -0.2735 -0.0133
+%! 0.0804 -0.0214 -0.2015 -0.1033
+%! 0.0626 -0.0484 -0.0865 -0.2383
+%! 0.0590 -0.1404 -0.1335 0.0367
+%! 0.0604 -0.2218 -0.2005 0.0567
+%! 0.0731 -0.2085 -0.1645 -0.1383
+%! 0.0997 -0.1695 -0.1575 -0.0058
+%! 0.1295 -0.1504 -0.2145 0.0242
+%! 0.1556 -0.1147 -0.2605 -0.0783
+%! 0.1749 -0.1237 -0.2505 0.0192
+%! 0.1857 -0.2165 -0.2005 -0.0008
+%! 0.1859 -0.2956 -0.2755 0.1167
+%! 0.1726 -0.2906 -0.3025 0.1892
+%! 0.1556 -0.2457 -0.1625 -0.0583
+%! 0.1417 -0.2264 -0.1545 0.0942
+%! 0.1269 -0.2038 -0.1905 0.0842
+%! 0.1129 -0.1753 -0.1785 -0.0658
+%! 0.0992 -0.1663 -0.1525 0.1442
+%! 0.0875 -0.1169 -0.0275 -0.1183
+%! 0.0790 -0.0582 -0.0545 -0.0283
+%! 0.0663 -0.0174 -0.0605 0.0767
+%! 0.0562 0.0665 0.0935 -0.3133
+%! 0.0531 0.1512 0.0685 -0.1408
+%! 0.0461 0.2151 0.0115 -0.1208
+%! 0.0374 0.2488 0.0475 -0.3083
+%! 0.0370 0.2115 0.1245 -0.2358
+%! 0.0492 0.1377 0.1595 -0.3158
+%! 0.0650 0.0180 -0.0255 -0.0283
+%! 0.0750 -0.1500 -0.0805 0.0317
+%! 0.0891 -0.3396 0.0055 -0.0733
+%! 0.1133 -0.5416 -0.0395 0.1742
+%! 0.1401 -0.7139 -0.1485 0.2417
+%! 0.1603 -0.9009 -0.2485 0.4742
+%! 0.1780 -1.1065 -0.1815 0.4742
+%! 0.2099 -1.2319 -0.1665 0.3667
+%! 0.2553 -1.2788 -0.3375 0.6817
+%! 0.3076 -1.2666 -0.3695 0.5117
+%! 0.3685 -1.2250 -0.3785 0.4217
+%! 0.4305 -1.1970 -0.5295 0.5917
+%! 0.4856 -1.1524 -0.7095 0.4792
+%! 0.5373 -1.1079 -0.7375 0.4042
+%! 0.5897 -1.0977 -0.6755 0.3842
+%! 0.6362 -1.0851 -0.8095 0.4792
+%! 0.6699 -1.0496 -0.9775 0.4067
+%! 0.6894 -0.9976 -1.0365 0.3367
+%! 0.6901 -0.9060 -1.0275 0.4142
+%! 0.6719 -0.7677 -0.9995 0.2142
+%! 0.6371 -0.6139 -1.0315 0.2267
+%! 0.5849 -0.4463 -0.9335 0.1717
+%! 0.5219 -0.2783 -0.7995 0.0167
+%! 0.4573 -0.1421 -0.7785 0.0567
+%! 0.4020 -0.0451 -0.6615 -0.1833
+%! 0.3559 0.0389 -0.5805 -0.0758
+%! 0.3050 0.2057 -0.5415 0.0042
+%! 0.2521 0.4271 -0.4775 -0.3708
+%! 0.2065 0.5643 -0.4755 -0.3033
+%! 0.1686 0.6331 -0.3475 -0.3783
+%! 0.1358 0.6722 -0.2895 -0.3708
+%! 0.1037 0.6797 -0.3195 -0.2408
+%! 0.0778 0.6585 -0.2575 -0.4533
+%! 0.0624 0.5950 -0.2515 -0.2558
+%! 0.0517 0.5407 -0.2285 -0.1908
+%! 0.0473 0.5332 -0.1815 -0.3383
+%! 0.0465 0.6105 -0.1795 -0.2558
+%! 0.0382 0.7590 -0.1585 -0.3383
+%! 0.0154 0.8366 -0.1745 -0.3208
+%! -0.0179 0.8604 -0.0875 -0.3733
+%! -0.0470 0.8872 0.1005 -0.5033
+%! -0.0626 0.8752 0.1875 -0.3783
+%! -0.0602 0.8813 0.2065 -0.5008
+%! -0.0402 0.8863 0.1245 -0.4583
+%! -0.0124 0.8391 0.0235 -0.3433
+%! 0.0192 0.7919 -0.0165 -0.4783
+%! 0.0554 0.7361 -0.0555 -0.3733
+%! 0.0929 0.6751 -0.0475 -0.3683
+%! 0.1304 0.6184 -0.0965 -0.3933
+%! 0.1639 0.5700 -0.2405 -0.2308
+%! 0.1967 0.6093 -0.2875 -0.4483
+%! 0.2338 0.7093 -0.3125 -0.4933
+%! 0.2645 0.8005 -0.3845 -0.3158
+%! 0.2869 0.8915 -0.3885 -0.5983
+%! 0.3091 0.9385 -0.4365 -0.6158
+%! 0.3251 0.9464 -0.5195 -0.4483
+%! 0.3311 0.9586 -0.5255 -0.6083
+%! 0.3313 0.9269 -0.5245 -0.5558
+%! 0.3224 0.8837 -0.5195 -0.4508
+%! 0.3000 0.9051 -0.5285 -0.4733
+%! 0.2679 0.9557 -0.5235 -0.4883
+%! 0.2360 1.0002 -0.4485 -0.6033
+%! 0.2076 0.9830 -0.4105 -0.5083
+%! 0.1760 0.8692 -0.4285 -0.3633
+%! 0.1395 0.7508 -0.4055 -0.3758
+%! 0.1052 0.6499 -0.3275 -0.3558
+%! 0.0774 0.5699 -0.2655 -0.3033
+%! 0.0517 0.6138 -0.2495 -0.3008
+%! 0.0213 0.7287 -0.1975 -0.3383
+%! -0.0184 0.7977 -0.1255 -0.3633
+%! -0.0668 0.8589 -0.0445 -0.4283
+%! -0.1192 0.9156 0.0605 -0.5133
+%! -0.1747 0.9026 0.1815 -0.4483
+%! -0.2323 0.8126 0.3325 -0.4383
+%! -0.2839 0.6700 0.4405 -0.4483
+%! -0.3233 0.5204 0.4685 -0.3258
+%! -0.3485 0.4290 0.5085 -0.3033
+%! -0.3550 0.4219 0.5735 -0.3533
+%! -0.3426 0.4701 0.5585 -0.3108
+%! -0.3172 0.5523 0.4785 -0.3108
+%! -0.2791 0.6432 0.4565 -0.4733
+%! -0.2315 0.6805 0.4115 -0.4533
+%! -0.1894 0.6917 0.2995 -0.3708
+%! -0.1591 0.7669 0.2685 -0.5033
+%! -0.1352 0.8755 0.2685 -0.5383
+%! -0.1136 0.9348 0.2695 -0.5183
+%! -0.0901 0.9419 0.2325 -0.6033
+%! -0.0643 0.9225 0.1315 -0.5033
+%! -0.0407 0.8522 0.0945 -0.4833
+%! -0.0225 0.7295 0.0265 -0.4108
+%! -0.0137 0.6081 -0.0405 -0.3008
+%! -0.0151 0.4741 -0.0395 -0.2883
+%! -0.0274 0.3231 -0.0395 -0.1408
+%! -0.0520 0.2134 -0.0215 -0.0958
+%! -0.0864 0.1254 -0.0295 -0.0233
+%! -0.1263 0.0256 0.0055 0.0442
+%! -0.1653 -0.0778 0.0825 0.0617
+%! -0.1988 -0.1712 0.1285 0.1667
+%! -0.2238 -0.2647 0.1945 0.1292
+%! -0.2412 -0.3741 0.1865 0.2617
+%! -0.2549 -0.4006 0.1935 0.3267
+%! -0.2573 -0.3240 0.2765 0.0892
+%! -0.2433 -0.2505 0.2645 0.1417
+%! -0.2231 -0.1534 0.2415 0.1142
+%! -0.2069 -0.0410 0.2085 0.0342
+%! -0.2000 0.0238 0.1645 0.0642
+%! -0.2001 0.0544 0.1985 -0.0808
+%! -0.2037 -0.0099 0.2065 0.0317
+%! -0.2104 -0.1788 0.2635 0.1042
+%! -0.2160 -0.3564 0.2935 0.1217
+%! -0.2186 -0.5024 0.2335 0.2917
+%! -0.2186 -0.5871 0.1745 0.3217
+%! -0.2108 -0.5916 0.1735 0.2567
+%! -0.1901 -0.5845 0.2265 0.1692
+%! -0.1679 -0.6241 0.1165 0.4067
+%! -0.1573 -0.6377 0.0325 0.3817
+%! -0.1530 -0.6189 0.0765 0.1342
+%! -0.1569 -0.6339 0.0265 0.4667
+%! -0.1734 -0.6051 0.0855 0.3067
+%! -0.1951 -0.5319 0.1385 0.1567
+%! -0.2237 -0.4474 0.1275 0.3742
+%! -0.2571 -0.3082 0.2075 0.0317
+%! -0.2879 -0.2211 0.1775 0.0992
+%! -0.3188 -0.1793 0.1985 0.1917
+%! -0.3429 -0.0881 0.3475 -0.1158
+%! -0.3541 -0.0297 0.3825 -0.0333
+%! -0.3618 0.0134 0.3265 0.0117
+%! -0.3716 0.0785 0.2805 -0.0033
+%! -0.3719 0.0910 0.3965 -0.1883
+%! -0.3546 0.0535 0.4395 -0.1708
+%! -0.3345 0.0320 0.3055 0.1492
+%! -0.3160 0.0821 0.3255 -0.1533
+%! -0.2892 0.1773 0.3235 -0.2208
+%! -0.2587 0.2922 0.2255 -0.0133
+%! -0.2214 0.4446 0.2095 -0.3608
+%! -0.1671 0.5549 0.2005 -0.4208
+%! -0.1050 0.5845 0.1275 -0.2683
+%! -0.0486 0.5966 -0.0525 -0.3083
+%! 0.0032 0.5948 -0.1795 -0.4133
+%! 0.0516 0.5806 -0.2155 -0.3683
+%! 0.0863 0.5808 -0.2465 -0.2633
+%! 0.1019 0.5570 -0.2995 -0.3258
+%! 0.1007 0.4895 -0.3555 -0.2833
+%! 0.0823 0.4475 -0.3445 -0.2083
+%! 0.0486 0.4517 -0.3085 -0.1983
+%! 0.0106 0.4461 -0.2495 -0.2433
+%! -0.0204 0.3968 -0.1525 -0.2433
+%! -0.0427 0.3374 -0.1145 -0.1308
+%! -0.0570 0.3425 -0.1035 -0.2283
+%! -0.0662 0.4094 -0.1785 -0.1583
+%! -0.0749 0.4980 -0.1755 -0.1458
+%! -0.0767 0.5981 -0.0495 -0.3908
+%! -0.0689 0.6888 -0.0465 -0.2683
+%! -0.0577 0.7881 -0.0605 -0.3833
+%! -0.0469 0.8619 -0.1055 -0.4733
+%! -0.0422 0.8321 -0.1485 -0.3083
+%! -0.0387 0.7466 -0.0695 -0.5408
+%! -0.0309 0.6634 -0.0885 -0.3933
+%! -0.0268 0.5794 -0.0795 -0.2908
+%! -0.0307 0.5002 -0.0445 -0.4458
+%! -0.0508 0.4760 -0.1505 -0.1608
+%! -0.0894 0.5679 -0.0875 -0.3133
+%! -0.1283 0.7078 0.1175 -0.6508
+%! -0.1663 0.7919 0.1845 -0.3833
+%! -0.2133 0.8265 0.2345 -0.4483
+%! -0.2599 0.8374 0.3125 -0.6758
+%! -0.2969 0.8162 0.3435 -0.4808
+%! -0.3210 0.7420 0.4165 -0.5058
+%! -0.3215 0.5730 0.4795 -0.5808
+%! -0.3013 0.2979 0.4085 -0.1883
+%! -0.2680 0.0143 0.3635 -0.1408
+%! -0.2172 -0.1931 0.2695 -0.1933
+%! -0.1568 -0.3482 0.0795 0.1542
+%! -0.0990 -0.5055 0.0565 0.0917
+%! -0.0475 -0.7178 -0.0035 0.2642
+%! -0.0025 -0.9341 -0.0295 0.3692
+%! 0.0428 -1.0834 -0.0325 0.2067
+%! 0.0814 -1.2097 -0.2115 0.6567
+%! 0.1090 -1.3175 -0.1885 0.4517
+%! 0.1336 -1.4127 -0.1795 0.3767
+%! 0.1445 -1.5065 -0.3165 0.8592
+%! 0.1424 -1.5894 -0.2155 0.4442
+%! 0.1379 -1.7094 -0.2005 0.5842
+%! 0.1194 -1.7628 -0.2475 0.9192
+%! 0.0871 -1.6474 -0.1865 0.5092
+%! 0.0541 -1.4894 -0.1495 0.4967
+%! 0.0205 -1.3807 -0.0495 0.4817
+%! -0.0189 -1.3325 -0.0365 0.4492
+%! -0.0630 -1.3452 -0.0175 0.4417
+%! -0.1055 -1.4067 0.0835 0.4317
+%! -0.1473 -1.4711 0.1405 0.6042
+%! -0.1886 -1.4663 0.2055 0.4892
+%! -0.2252 -1.4405 0.2435 0.5292
+%! -0.2545 -1.4530 0.3345 0.5367
+%! -0.2726 -1.4790 0.3915 0.4667
+%! -0.2804 -1.5291 0.3805 0.6292
+%! -0.2800 -1.6064 0.4145 0.5967
+%! -0.2711 -1.6577 0.4125 0.6642
+%! -0.2543 -1.6503 0.4445 0.6142
+%! -0.2332 -1.6413 0.3945 0.6267
+%! -0.2186 -1.6598 0.2995 0.8167
+%! -0.2105 -1.6905 0.3655 0.6167
+%! -0.2016 -1.7773 0.3605 0.6767
+%! -0.1990 -1.8767 0.3035 0.8442
+%! -0.2109 -1.9359 0.2905 0.8542
+%! -0.2309 -1.9978 0.3625 0.8267
+%! -0.2464 -2.0725 0.4795 0.7117
+%! -0.2598 -2.1566 0.3915 1.0417
+%! -0.2729 -2.2060 0.4415 0.9092
+%! -0.2729 -2.2146 0.5805 0.6867
+%! -0.2651 -2.1907 0.4845 1.0592
+%! -0.2563 -2.0544 0.4705 0.7217
+%! -0.2420 -1.8793 0.4675 0.5817
+%! -0.2313 -1.7438 0.4625 0.7617
+%! -0.2284 -1.5795 0.5095 0.4092
+%! -0.2306 -1.4005 0.4165 0.4717
+%! -0.2380 -1.2226 0.4445 0.3817
+%! -0.2426 -1.0615 0.5195 0.1817
+%! -0.2387 -0.9443 0.4795 0.2317
+%! -0.2262 -0.8374 0.4195 0.1167
+%! -0.2090 -0.7585 0.3505 0.1867
+%! -0.1902 -0.6965 0.3525 0.1167
+%! -0.1702 -0.6160 0.2955 0.0917
+%! -0.1530 -0.5535 0.2285 0.1892
+%! -0.1379 -0.5331 0.2845 0.0267
+%! -0.1232 -0.5745 0.2955 0.1692
+%! -0.1141 -0.6693 0.2535 0.3117
+%! -0.1051 -0.7692 0.2395 0.2017
+%! -0.0888 -0.8714 0.2465 0.3167
+%! -0.0706 -0.9784 0.2345 0.4467
+%! -0.0576 -1.0470 0.1525 0.4467
+%! -0.0533 -1.0386 0.0805 0.4942
+%! -0.0564 -0.9563 0.0985 0.4417
+%! -0.0595 -0.8278 0.1525 0.2792
+%! -0.0590 -0.6936 0.1335 0.2642
+%! -0.0571 -0.6112 0.1255 0.2042
+%! -0.0560 -0.6221 0.1465 0.1867
+%! -0.0570 -0.6806 0.1365 0.2192
+%! -0.0627 -0.7056 0.1195 0.1917
+%! -0.0791 -0.6703 0.1115 0.3067
+%! -0.1047 -0.5488 0.2345 0.1217
+%! -0.1276 -0.3868 0.3215 -0.0508
+%! -0.1452 -0.2989 0.2885 0.0592
+%! -0.1579 -0.3142 0.3435 -0.1133
+%! -0.1698 -0.3930 0.3155 0.1142
+%! -0.1872 -0.4033 0.3395 0.1592
+%! -0.1990 -0.3194 0.4595 -0.2333
+%! -0.2078 -0.2826 0.3695 0.1642
+%! -0.2264 -0.2406 0.3875 0.1117
+%! -0.2433 -0.1958 0.4715 -0.2008
+%! -0.2512 -0.2855 0.4245 0.1467
+%! -0.2515 -0.4154 0.4715 0.0942
+%! -0.2396 -0.5006 0.4585 0.1692
+%! -0.2129 -0.5207 0.4565 0.2067
+%! -0.1705 -0.4377 0.3985 0.0617
+%! -0.1208 -0.3532 0.2195 0.2492
+%! -0.0686 -0.3256...
[truncated message content] |
|
From: <par...@us...> - 2012-05-30 07:29:10
|
Revision: 10535
http://octave.svn.sourceforge.net/octave/?rev=10535&view=rev
Author: paramaniac
Date: 2012-05-30 07:28:58 +0000 (Wed, 30 May 2012)
Log Message:
-----------
control-devel: add a test from the DaISy library
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/moen4.m
Modified: trunk/octave-forge/extra/control-devel/inst/moen4.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-29 21:54:39 UTC (rev 10534)
+++ trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-05-30 07:28:58 UTC (rev 10535)
@@ -303,3 +303,444 @@
%!assert (INFO.Ry, Rye, 1e-4);
%!assert (INFO.S, Se, 1e-4);
%!assert (X0, X0e, 1e-4);
+
+
+## [96-003] Data of a 120 MW power plant (Pont-sur-Sambre, France)
+%!shared SYS, Ae, Be, Ce, De
+%! U = [ -811 -592 421 -680 -681
+%! -812 -619 477 -685 -651
+%! -817 -565 538 -678 -677
+%! -695 -725 536 -674 -702
+%! -697 -571 531 -676 -685
+%! -697 -618 533 -681 -721
+%! -702 -579 549 -677 -699
+%! -703 -487 575 -677 -694
+%! -705 -449 561 -679 -678
+%! -705 -431 563 -680 -692
+%! -707 -502 561 -679 -686
+%! -707 -583 530 -676 -751
+%! -710 -458 540 -677 -700
+%! -713 -469 543 -679 -731
+%! -715 -506 549 -684 -635
+%! -713 -590 532 -681 -715
+%! -714 -582 528 -676 -696
+%! -713 -575 538 -679 -690
+%! -716 -382 557 -682 -701
+%! -716 -476 558 -679 -690
+%! -718 -425 565 -678 -686
+%! -719 -409 562 -679 -615
+%! -719 -508 523 -677 -737
+%! -721 -569 523 -679 -722
+%! -723 -434 542 -681 -735
+%! -723 -395 544 -676 -704
+%! -723 -428 542 -677 -729
+%! -722 -402 537 -677 -706
+%! -725 -380 534 -681 -696
+%! -726 -324 549 -676 -701
+%! -726 -211 578 -675 -614
+%! -727 -113 569 -677 -738
+%! -727 -208 554 -676 -737
+%! -727 -320 523 -684 -697
+%! -727 944 605 -680 -587
+%! -729 396 654 -681 -708
+%! -729 754 637 -682 -685
+%! -727 419 522 -677 -715
+%! -729 378 494 -668 -703
+%! -728 363 493 -669 -682
+%! -729 390 496 -665 -713
+%! -729 106 483 -664 -711
+%! -729 32 495 -661 -718
+%! -729 517 585 -661 -641
+%! -729 455 625 -659 -703
+%! -730 521 649 -687 -669
+%! -730 540 627 -689 -705
+%! -731 514 605 -694 -682
+%! -585 525 558 -685 -611
+%! -586 534 520 -680 -668
+%! -586 539 531 -681 -679
+%! -585 519 507 -682 -663
+%! -588 513 505 -667 -668
+%! -587 509 539 -680 -616
+%! -587 512 535 -668 -628
+%! -588 514 557 -667 -648
+%! -588 553 563 -676 -613
+%! -589 519 559 -684 -638
+%! -589 521 563 -682 -652
+%! -588 518 547 -678 -597
+%! -589 552 549 -688 -630
+%! -589 520 535 -685 -623
+%! -589 547 542 -678 -619
+%! -589 549 531 -684 -524
+%! -588 544 522 -1540 -580
+%! -588 564 555 -1538 -584
+%! -588 684 545 -1541 -564
+%! -590 558 546 -1541 -609
+%! -589 552 537 -1550 -601
+%! -591 532 526 -1548 -580
+%! -590 544 524 -1542 -565
+%! -591 559 535 -1538 -604
+%! -592 555 542 -1548 -629
+%! -591 577 532 -1549 -587
+%! -593 581 530 -1543 -585
+%! -592 562 540 -1548 -583
+%! -591 568 546 -1536 -587
+%! -593 550 557 -1533 -569
+%! -592 550 537 -1518 -527
+%! -593 568 551 -1533 -582
+%! -590 528 540 -1529 -492
+%! -590 542 532 -1525 -585
+%! -590 556 535 -1522 -606
+%! -591 637 535 -1516 -571
+%! -591 608 539 -1512 -582
+%! -591 545 527 -1510 -577
+%! -591 603 534 -1507 -548
+%! -592 567 521 -1507 -539
+%! -594 560 530 -1503 -583
+%! -422 549 534 -1487 -542
+%! 5 619 550 -1488 -609
+%! 5 572 541 -1487 -545
+%! 185 564 529 -1488 -528
+%! 185 571 531 -1497 -512
+%! 187 590 545 -1489 -472
+%! 186 658 544 -1485 -501
+%! 185 639 553 -1486 -512
+%! 187 604 532 -1486 -490
+%! 186 584 535 -1486 -505
+%! 187 572 541 -1488 -479
+%! 186 570 531 -1486 -477
+%! 187 579 539 -1488 -482
+%! 187 694 537 -1487 -512
+%! 187 727 546 -1494 -474
+%! 186 838 545 -1493 -470
+%! 186 703 527 -1492 -510
+%! 185 634 502 -1490 -486
+%! 185 659 502 -1490 -446
+%! 185 569 483 -1488 -472
+%! 184 552 494 -1487 -468
+%! 183 541 492 -1487 -411
+%! 184 544 519 -1486 -471
+%! 184 577 538 -1487 -420
+%! 183 536 546 -1490 -455
+%! 183 586 567 -1487 -421
+%! 183 554 566 -1485 -449
+%! 183 564 558 -1486 -455
+%! 182 558 546 -1486 -426
+%! 183 659 543 -1486 -398
+%! 182 545 531 -1487 -425
+%! 182 542 513 -1487 -297
+%! 181 549 502 -1486 -407
+%! 180 558 963 -1487 -391
+%! 182 546 1105 -1487 -388
+%! 179 540 1166 -1487 -387
+%! 181 519 1221 -1487 -390
+%! 179 510 1254 -1486 -382
+%! 176 528 1289 -1487 -379
+%! 178 510 1310 -1486 -389
+%! 178 504 1313 -1484 -357
+%! 178 507 1328 -1484 -383
+%! 178 519 1333 -1489 -426
+%! 183 509 1317 -1483 -377
+%! 177 585 1331 -1483 -352
+%! 176 765 1334 -1484 -375
+%! 196 689 1328 -1484 -381
+%! 178 570 1306 -1486 -398
+%! 176 997 1338 -1483 -360
+%! 176 839 1318 -1487 -366
+%! 176 879 1324 -90 -394
+%! 175 763 1302 -92 -436
+%! 175 739 1307 -92 -430
+%! 176 791 1336 -86 -446
+%! 175 774 1301 -91 -405
+%! 174 782 1304 -91 -414
+%! 189 722 1335 -92 -360
+%! 173 871 1341 -91 -417
+%! 173 825 1337 -86 -383
+%! 193 836 1332 -86 -449
+%! 174 832 1325 -87 -432
+%! 175 834 1326 -88 -383
+%! 176 899 1332 -86 -341
+%! 195 827 1296 -84 -371
+%! 174 897 1311 -86 -366
+%! 195 848 1320 -86 -416
+%! 192 777 1286 -83 126
+%! 171 788 1309 -84 153
+%! 171 810 1307 -84 173
+%! 174 758 1297 -85 182
+%! 188 910 1329 -1036 236
+%! 189 944 1344 -1669 245
+%! 196 859 1330 -1668 249
+%! 190 797 1320 -1669 231
+%! 191 784 1341 -1665 248
+%! 184 737 1335 -1656 240
+%! 182 733 1342 -1655 245
+%! 182 696 1350 -1653 253
+%! 185 614 1341 -1658 236
+%! 188 708 1370 -1654 238
+%! 191 729 1366 -1655 237
+%! 183 714 1380 -1651 262
+%! 186 695 1378 -1651 257
+%! 189 758 1387 -1649 276
+%! 164 572 1364 -1649 257
+%! 163 587 1343 -1650 246
+%! 161 683 1353 -1658 261
+%! 160 572 1359 -1654 261
+%! 160 607 1355 -1651 267
+%! 158 580 1349 -1655 273
+%! 161 631 1362 -1652 272
+%! 160 706 1382 -1657 298
+%! 161 601 1356 -1654 277
+%! 159 570 1354 -1653 274
+%! 158 547 1347 -1657 261
+%! 158 582 1332 -1657 257
+%! 157 570 1340 -1657 267
+%! 154 556 1343 -1657 268
+%! 157 537 1345 -1657 -425
+%! 158 555 1331 -1653 -581
+%! 158 551 1315 -1654 -643
+%! 159 590 1322 -1656 -687
+%! 160 566 1315 -1657 -737
+%! 160 553 1315 -1653 -767
+%! 161 644 1327 -1396 -731
+%! 159 640 1335 -577 -639
+%! 161 726 1334 -577 -730
+%! 175 729 1310 -573 -711
+%! 175 854 1330 -576 -690
+%! 161 729 1313 -573 -636
+%! 176 706 1314 -583 -745 ];
+%!
+%!
+%! Y = [ 117 129 -47
+%! 113 141 -42
+%! 83 150 -37
+%! 144 174 -3
+%! 156 196 19
+%! 174 192 6
+%! 171 193 3
+%! 169 224 14
+%! 155 211 14
+%! 137 175 4
+%! 130 165 13
+%! 145 172 21
+%! 128 173 30
+%! 119 194 30
+%! 127 188 11
+%! 147 176 5
+%! 162 183 19
+%! 173 178 14
+%! 131 161 4
+%! 112 139 -11
+%! 91 133 -10
+%! 75 127 -7
+%! 80 115 -10
+%! 87 132 6
+%! 51 151 6
+%! 28 167 3
+%! 22 167 -3
+%! 9 171 4
+%! 7 158 6
+%! -11 155 5
+%! -87 149 -8
+%! -153 146 -12
+%! -184 116 -20
+%! -167 118 8
+%! -563 78 -4
+%! -694 -3 -155
+%! -1000 -39 -168
+%! -1101 -28 -135
+%! -1080 55 6
+%! -1075 162 89
+%! -1035 254 117
+%! -894 329 148
+%! -747 409 218
+%! -774 416 231
+%! -781 361 189
+%! -825 288 171
+%! -897 220 166
+%! -960 175 169
+%! -994 157 196
+%! -1009 202 237
+%! -1003 267 258
+%! -981 326 267
+%! -947 361 276
+%! -921 369 280
+%! -894 369 292
+%! -876 359 302
+%! -869 356 318
+%! -857 334 316
+%! -857 336 320
+%! -868 331 323
+%! -873 316 324
+%! -870 301 326
+%! -867 294 323
+%! -878 293 322
+%! -877 300 309
+%! -891 317 233
+%! -919 308 160
+%! -919 296 111
+%! -898 304 94
+%! -884 336 101
+%! -876 344 92
+%! -877 347 77
+%! -859 354 79
+%! -852 368 84
+%! -852 362 79
+%! -846 337 67
+%! -846 322 71
+%! -836 325 82
+%! -831 330 81
+%! -848 332 84
+%! -834 334 89
+%! -830 314 87
+%! -838 322 93
+%! -859 313 73
+%! -864 317 59
+%! -852 318 48
+%! -850 313 59
+%! -858 327 65
+%! -849 332 82
+%! -822 332 77
+%! -721 377 121
+%! -517 492 193
+%! -306 596 226
+%! -117 683 244
+%! 22 716 212
+%! 87 684 157
+%! 151 642 129
+%! 198 599 105
+%! 233 549 93
+%! 244 512 83
+%! 255 487 70
+%! 247 453 49
+%! 225 437 32
+%! 175 410 4
+%! 105 386 -12
+%! 86 379 -35
+%! 101 382 -16
+%! 111 412 3
+%! 158 451 31
+%! 217 486 37
+%! 259 504 29
+%! 301 503 24
+%! 317 495 7
+%! 354 472 -3
+%! 361 441 -8
+%! 368 413 -32
+%! 362 381 -34
+%! 354 370 -42
+%! 310 360 -36
+%! 318 369 -36
+%! 313 372 -36
+%! 298 377 -28
+%! 311 341 -34
+%! 339 292 -33
+%! 362 223 -52
+%! 380 161 -73
+%! 397 113 -90
+%! 402 64 -117
+%! 398 30 -136
+%! 396 12 -148
+%! 390 -22 -189
+%! 377 -47 -211
+%! 366 -54 -227
+%! 331 -52 -240
+%! 227 -74 -289
+%! 168 -117 -355
+%! 149 -104 -357
+%! 10 -101 -344
+%! -77 -115 -374
+%! -116 -88 -330
+%! -106 -53 -223
+%! -82 -4 -98
+%! -99 5 -40
+%! -89 22 3
+%! -79 34 26
+%! -74 34 42
+%! -101 28 61
+%! -102 28 73
+%! -124 22 74
+%! -125 19 90
+%! -132 9 100
+%! -167 8 102
+%! -161 7 90
+%! -185 19 105
+%! -207 29 110
+%! -295 24 159
+%! -361 17 293
+%! -411 12 436
+%! -401 48 540
+%! -441 86 601
+%! -475 118 553
+%! -448 146 472
+%! -401 172 428
+%! -366 197 400
+%! -304 231 391
+%! -252 253 387
+%! -183 261 399
+%! -123 255 405
+%! -107 247 397
+%! -76 228 392
+%! -71 205 389
+%! -66 192 381
+%! -92 171 371
+%! -40 178 356
+%! -15 185 353
+%! -27 170 359
+%! 15 173 351
+%! 23 175 356
+%! 44 178 338
+%! 47 171 323
+%! 4 153 294
+%! 23 141 295
+%! 39 141 306
+%! 52 139 313
+%! 53 142 313
+%! 69 142 312
+%! 74 136 306
+%! 196 164 277
+%! 347 208 133
+%! 461 230 -58
+%! 526 204 -251
+%! 581 161 -385
+%! 588 119 -458
+%! 549 63 -528
+%! 497 5 -550
+%! 420 -20 -498
+%! 327 -49 -464
+%! 198 -78 -438
+%! 154 -80 -409
+%! 130 -60 -377 ];
+%!
+%! DAT = iddata (Y, U, 1228.8);
+%!
+%! [SYS, X0] = moen4 (DAT, "s", 10, "n", 8, "rcond", 0.0, "tol", -1.0);
+%!
+%! Ae = [ 0.9811 0.0574 0.3270 0.0003 0.0358 0.0403 -0.1366 0.0276
+%! 0.1043 0.7634 -0.1308 0.2252 0.0203 -0.0117 -0.2328 -0.2516
+%! -0.0612 0.1437 0.8378 -0.2400 0.0367 -0.3205 -0.0367 -0.0978
+%! -0.0213 0.0149 0.0706 0.8902 0.2415 -0.1329 0.3080 -0.0143
+%! -0.0076 0.0680 0.0557 -0.0469 0.7084 0.2852 0.2565 0.1508
+%! -0.0016 0.0603 0.0071 -0.0448 -0.0870 0.8608 0.1495 -0.1148
+%! 0.0046 0.0120 0.0040 -0.0272 0.0117 -0.0654 0.8348 -0.4557
+%! -0.0058 0.0133 -0.0112 -0.0171 -0.0353 -0.0892 0.3328 0.7650 ];
+%!
+%! Ce = [ -0.1635 0.6294 0.1157 -0.2129 0.0812 -0.0238 0.0745 0.2027
+%! 0.1775 0.1736 -0.2839 -0.0750 -0.4986 -0.1254 0.3740 0.1598
+%! 0.2012 -0.0045 -0.4447 -0.3289 0.4767 -0.3377 0.0200 -0.0978 ];
+%!
+%! Be = [ -0.0436 0.1911 -0.1345 0.0409 0.2828
+%! 0.5541 -0.4223 0.0001 -0.0151 -0.1778
+%! -0.3644 0.2868 0.0754 -0.0457 -0.1691
+%! -0.0372 0.1270 -0.0772 -0.0714 -0.0255
+%! -0.1251 -0.0021 0.1094 0.0576 0.2121
+%! -0.1372 0.1139 0.0030 0.0141 0.0493
+%! -0.0187 0.0712 -0.0042 0.0075 -0.0221
+%! -0.0834 0.0509 0.0322 -0.0089 0.0009 ];
+%!
+%! De = [ 0.1581 -0.3111 -0.0350 0.0179 -0.1403
+%! -0.0037 -0.0461 -0.1177 0.0092 -0.0242
+%! 0.0476 -0.0237 -0.0159 0.0174 0.0464 ];
+%!
+%!assert (abs (SYS.A), abs (Ae), 1e-4);
+%!assert (abs (SYS.B), abs (Be), 1e-4);
+%!assert (abs (SYS.C), abs (Ce), 1e-4);
+%!assert (abs (SYS.D), abs (De), 1e-4);
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <nir...@us...> - 2012-05-29 21:54:45
|
Revision: 10534
http://octave.svn.sourceforge.net/octave/?rev=10534&view=rev
Author: nir-krakauer
Date: 2012-05-29 21:54:39 +0000 (Tue, 29 May 2012)
Log Message:
-----------
fixed cspas_sel to calculate trend uncertainty when the input y is an array
Modified Paths:
--------------
trunk/octave-forge/main/splines/inst/csaps_sel.m
Modified: trunk/octave-forge/main/splines/inst/csaps_sel.m
===================================================================
--- trunk/octave-forge/main/splines/inst/csaps_sel.m 2012-05-29 19:26:00 UTC (rev 10533)
+++ trunk/octave-forge/main/splines/inst/csaps_sel.m 2012-05-29 21:54:39 UTC (rev 10534)
@@ -108,7 +108,7 @@
Hd = diag(H);
end
- sigma2 = MSR * (n / (n-Ht)); #estimated data error variance (wahba83)
+ sigma2 = mean(MSR(:)) * (n / (n-Ht)); #estimated data error variance (wahba83)
unc_y = sqrt(sigma2 * Hd ./ w); #uncertainty (SD) of fitted curve at each input x-value (hutchinson86)
## solve for the scaled second derivatives u and for the function values a at the knots (if p = 1, a = y)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <be...@us...> - 2012-05-29 19:26:07
|
Revision: 10533
http://octave.svn.sourceforge.net/octave/?rev=10533&view=rev
Author: benjf5
Date: 2012-05-29 19:26:00 +0000 (Tue, 29 May 2012)
Log Message:
-----------
Fixed a typo, switched a demo to a test, added fastlscomplex.c. Don't recommend running it.
Modified Paths:
--------------
trunk/octave-forge/extra/lssa/cubicwgt.m
trunk/octave-forge/extra/lssa/lscomplex.m
trunk/octave-forge/extra/lssa/lscorrcoeff.m
Added Paths:
-----------
trunk/octave-forge/extra/lssa/fastlscomplex.c
Modified: trunk/octave-forge/extra/lssa/cubicwgt.m
===================================================================
--- trunk/octave-forge/extra/lssa/cubicwgt.m 2012-05-29 15:51:05 UTC (rev 10532)
+++ trunk/octave-forge/extra/lssa/cubicwgt.m 2012-05-29 19:26:00 UTC (rev 10533)
@@ -19,15 +19,12 @@
## 1 + ( x ^ 2 * ( 2 x - 3 ) ), assuming x is in [-1,1].
## @end deftypefn
-%!demo
-%! h = 2;
-%! hcw = cubicwgt(h)
-%! m = 0.01;
-%! mcw = cubicwgt(m)
+%!shared h, m, k
+%! h = 2; m = 0.01;
%! k = [ 0 , 3 , 1.5, -1, -0.5, -0.25, 0.75 ];
-%! kcw = cubicwgt(k)
-%! kt = k';
-%! ktcw = cubicwgt(kt);
+%!assert( cubicwgt(h), 0 );
+%!assert( cubicwgt(m), 1 + m ^ 2 * ( 2 * m - 3 ));
+%!assert( cubicwgt(k), [ 1.00000 0.00000 0.00000 0.00000 0.50000 0.84375 0.15625], 1E-6);
%! ## Tests cubicwgt on two scalars and two vectors; cubicwgt will work on any array input.
Added: trunk/octave-forge/extra/lssa/fastlscomplex.c
===================================================================
--- trunk/octave-forge/extra/lssa/fastlscomplex.c (rev 0)
+++ trunk/octave-forge/extra/lssa/fastlscomplex.c 2012-05-29 19:26:00 UTC (rev 10533)
@@ -0,0 +1,243 @@
+// fastlscomplex, implemented for Octave
+
+
+#include <octave/oct.h>
+#ifdef __cplusplus
+extern "C"
+{
+#endif
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <complex.h>
+#ifdef __cplusplus
+}
+#endif
+
+#define MAXCOLUMN 8
+
+/* We use the ISO C99 complex facility as implemented by GCC, but only in two places */
+
+typedef _Complex double Complex;
+typedef double Real;
+
+#define RE(x_) (__real__ x_)
+#define IM(x_) (__imag__ x_)
+#define CSET(x_, r_, i_) (RE(x_) = (r_), IM(x_) = (i_))
+#define PHISET(x_, p_) CSET(x_, cos(tmp=p_), sin(tmp))
+#define SCALEPHISET(x_, f_, p_) CSET(x_, f_ cos(tmp=p_), f_ sin(tmp))
+
+// Here the Data structure is included, but it's only ever a feature of the
+// #STANDALONE cases, which we're not dealing with here.
+
+typedef struct
+{ Real x, t;
+} XTElem;
+
+inline double sqr(double x)
+{ return x*x; }
+
+/* PNUM has to match the definition of EXP_IOT_SERIES! */
+#define PNUM 12
+#define SETXT(p_, op_, x_, t_) (p_)->x op_ x_; (p_++)->t op_ t_;
+#define SETT(p_, op_, x_, t_) *p_++ op_ t_;
+#define SETX(p_, op_, x_, t_) *p_++ op_ x_;
+/* h is a complex aux. variable; it is used for assignment times I everywhere */
+#define SETIX(p_, op_, x_, t_) h = x_; RE(*(p_)) op_ -IM(h); IM(*(p_)) op_ RE(h); p_++;
+
+ /* Macro that sums up the power series terms into the power series
+ * element record pointed to by p_.
+ * By using = and += for op_, initial setting and accumulation can be selected.
+ * t_ is the expression specifying the abscissa value. set_ can be either
+ * SETXT to set the x and t fields of an XTElem record, or SETT/SETX to set
+ * the elements of a Real array representing alternately real and imaginary
+ * values.
+ */
+ // -10 points, comments don't match method.
+#define EXP_IOT_SERIES(p_, el_, t_, op_, setr_, seti_) \
+{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \
+ tt = -t; seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(-1.0/3.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(-1.0/5.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(-1.0/7.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(-1.0/9.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(-1.0/11.0); seti_(p_, op_, x*tt, tt) \
+}
+
+/* same as the above, but without alternating signs */
+#define EXPIOT_SERIES(p_, el_, t_, op_, setr_, seti_) \
+{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \
+ seti_(p_, op_, x*t, t ) \
+ tt = t*t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/3.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/5.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/7.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/9.0); seti_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \
+ tt *= t*(1.0/11.0); seti_(p_, op_, x*tt, tt) \
+}
+
+# define SRCARG Real *tptr, Real *xptr, int n, double *lengthptr
+# define SRCVAR int k; Real length = *lengthptr;
+# define SRCT tptr[k]
+# define SRCX xptr[k]
+# define SRCFIRST k = 0
+# define SRCAVAIL (k<n)
+
+
+
+DEFUN_DLD (fastlscomplex, args, nargout, "Computes a rapid complex least squares transform.")
+{
+ int nargs = args.length();
+ if ( nargs != 7 ) print_usage();
+ else {
+ const RowVector tvals = args(0).row_vector_value();
+ const ComplexRowVector xvals = args(1).complex_row_vector_value();
+ const int nval = args(2).int_value();
+ const int lenval = args(3).int_value();
+ const int ncoeffval = args(4).int_value();
+ const int noctaveval = args(5).int_value();
+ const int omegamaxval = args(6).int_value();
+ if ( !error_value ) {
+ ComplexRowVector ret_series ( ( ncoeffval * noctaveval ) , 0.0 );
+ octave_idx_type numt = tvals.numel(), numx = xvals.numel(),
+ numret = ret_series.numel();
+ }
+ }
+
+}
+
+void print_usage() {
+ octave_stdout << "Prints the usage string for fastlscomplex, "
+ << "once I get around to it.\n";
+ // Really, I'll replace this soon.
+
+}
+
+void fastlscomplex(Real *tptr, Complex *xptr, int *nptr, double *lengthptr, int *ncoeffptr, int *noctaveptr, Real *omegamaxptr, Complex *rp)
+{
+ int k, n = *nptr, ncoeff = *ncoeffptr, noctave = *noctaveptr;
+ Real length = *lengthptr, omegamax = *omegamaxptr;
+
+ struct SumVec /* Precomputation record */
+ { struct SumVec *next; /* A singly linked list */
+ Complex elems[PNUM]; /* the summed power series elements */
+ int cnt; /* number of samples for which the power series elements were added */
+ }
+ *shead, *stail, *sp, *sq; /* power series element lists */
+ Real dtelems[PNUM], /* power series elements of exp(-i dtau) */
+ *dte, *r, /* Pointers into dtelems */
+ x, /* abscissa and ordinate value, p-th power of t */
+ tau, tau0, te, /* Precomputation range centers and range end */
+ tau0d, /* tau_h of first summand range at level d */
+ dtau = (0.5*M_PI)/omegamax,/* initial precomputation interval radius */
+ dtaud, /* precomputation interval radius at d'th merging step */
+ n_1 = 1.0/n, /* reciprocal of sample count */
+ ooct, o, omul, /* omega/mu for octave's top omega and per band, mult. factor */
+ omegaoct, omega, /* Max. frequency of octave and current frequency */
+ on_1, /* n_1*(omega/mu)^p, n_1*(2*omega/mu)^p */
+ mu = (0.5*M_PI)/length, /* Frequency shift: a quarter period of exp(i mu t) on length */
+ tmp;
+ Complex zeta, zz, /* Accumulators for spectral coefficients */
+ e, emul, /* summation factor exp(-i o tau_h) */
+ h, eoelems[PNUM], oeelems[PNUM],
+ *eop, *oep,
+ *ep, *op, /* Pointer into the current choice of eoelems and oeelems */
+ *p, *q, *pe;
+ int i, j; /* Coefficient and octave counter */
+
+ /* Subdivision and Precomputation */
+ SRCFIRST;
+ tau = SRCT+dtau; te = tau+dtau;
+ tau0 = tau;
+ shead = stail = sp = alloca(sizeof(*sp)); sp->next = 0;
+ for(te = SRCT+2*dtau; ; )
+ { x = SRCX;
+ EXP_IOT_SERIES(p, sp->elems, mu*(SRCT-tau), =, SETX, SETIX); sp->cnt = 1;
+ for(SRCNEXT; SRCAVAIL && SRCT<te; SRCNEXT)
+ { x = SRCX;
+ EXP_IOT_SERIES(p, sp->elems, mu*(SRCT-tau), +=, SETX, SETIX); sp->cnt++;
+ }
+ if(!SRCAVAIL) break;
+ tau = te+dtau; te = tau+dtau;
+ sp = alloca(sizeof(*sp)); stail->next = sp; stail = sp; sp->next = 0; sp->cnt = 0;
+ }
+
+ ooct = omegamax/mu;
+ omul = exp(-M_LN2/ncoeff);
+ omegaoct = omegamax;
+ tau0d = tau0;
+ dtaud = dtau;
+ /*** Loop over Octaves ***/
+ for(j = noctave; ; ooct *= 0.5, omegaoct *= 0.5, tau0d += dtaud, dtaud *= 2)
+ { /*** Results per frequency ***/
+ for(i = ncoeff, o = ooct, omega = omegaoct; i--; o *= omul, omega *= omul)
+ { PHISET(e, -omega*tau0d);
+ PHISET(emul, -2*omega*dtaud);
+ for(zeta = 0, sp = shead; sp; sp = sp->next, e *= emul)
+ if(sp->cnt)
+ { for(zz = 0, p = sp->elems, pe = p+PNUM, on_1 = n_1; p < pe; )
+ { zz += *p++ * on_1; on_1 *= o; }
+ zeta += e * zz;
+ }
+ *rp++ = zeta;
+ }
+ if(--j<=0) break; /* avoid unnecessary merging at the end */
+ /* Merging of the s_h;
+ * 4 different possibilities, depending on which of the merged ranges actually contain data.
+ * The computation is described in the paper; The result of a merger is stored in the
+ * left precomputation record (sp). Before one power series element is stored, the
+ * sum and difference of the original values *p and *q are stored in eoelems and oeelems,
+ * respectively. The result is then stored in *p.
+ */
+ EXPIOT_SERIES(r, dtelems, mu*dtaud, =, SETT, SETT);
+ for(sp = shead; sp; sp = sp->next)
+ { if(!(sq = sp->next) || !sq->cnt)
+ { if(sp->cnt)
+ for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+PNUM; p < pe; p++, dte++)
+ { ep = eop; *eop++ = *p;
+ for(r = dtelems, *p = *ep * *r; ; )
+ { if(++r>=dte) break; --ep; h = *ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; --ep; *p -= *ep * *r;
+ if(++r>=dte) break; --ep; h = -*ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; --ep; *p += *ep * *r;
+ }
+ }
+ if(!sq) break; /* reached the last precomputation range */
+ }
+ else
+ if(sp->cnt)
+ for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+PNUM;
+ p < pe; p++, q++, dte++)
+ { ep = eop; *eop++ = *p+*q; *oep++ = *p-*q; op = oep;
+ for(r = dtelems, *p = *ep * *r; ; )
+ { if(++r>=dte) break; op -= 2; h = *op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; ep -= 2; *p -= *ep * *r;
+ if(++r>=dte) break; op -= 2; h = -*op * *r; RE(*p) -= IM(h); IM(*p) += RE(h);
+ if(++r>=dte) break; ep -= 2; *p += *ep * *r;
+ }
+ }
+ else
+ for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+PNUM; q<pe; q++, dte++)
+ { ep = eop; *eop++ = *q;
+ for(r = dtelems, *q = *ep * *r; ; )
+ { if(++r>=dte) break; --ep; h = *ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
+ if(++r>=dte) break; --ep; *p -= *ep * *r;
+ if(++r>=dte) break; --ep; h = -*ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h);
+ if(++r>=dte) break; --ep; *q += *ep * *r;
+ }
+ }
+
+ sp->cnt += sq->cnt; sp->next = sq->next; /* free(sq) if malloc'ed */
+ }
+ }
+}
Modified: trunk/octave-forge/extra/lssa/lscomplex.m
===================================================================
--- trunk/octave-forge/extra/lssa/lscomplex.m 2012-05-29 15:51:05 UTC (rev 10532)
+++ trunk/octave-forge/extra/lssa/lscomplex.m 2012-05-29 19:26:00 UTC (rev 10533)
@@ -20,6 +20,7 @@
## series, considering frequencies up to @var{maxfreq}, over @var{numoctaves}
## octaves and @var{numcoeff} coefficients.
##
+## @seealso{lsreal}
## @end deftypefn
Modified: trunk/octave-forge/extra/lssa/lscorrcoeff.m
===================================================================
--- trunk/octave-forge/extra/lssa/lscorrcoeff.m 2012-05-29 15:51:05 UTC (rev 10532)
+++ trunk/octave-forge/extra/lssa/lscorrcoeff.m 2012-05-29 19:26:00 UTC (rev 10533)
@@ -57,7 +57,7 @@
rx2 = x2(mask);
ry2 = y2(mask);
## I've used the same mask for all of these as it's an otherwise unimportant variable ... can this leak memory?
- lnength(rx1) ##printing this length is probably used as a warning if 0 is returned; I inculded it
+ length(rx1) ##printing this length is probably used as a warning if 0 is returned; I included it
## in particular to maintain an exact duplicate of the R function.
s = sum( wgt( ( rx1 - t ) .* so ) ) * sum( wgt( ( rx2 - t ) .* so ) );
coeff = ifelse( s != 0 , ( sum( wgt( ( rx1 - t ) .* so ) .* exp( i .* o .* rx1 ) .* ry1 )
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-29 15:51:16
|
Revision: 10532
http://octave.svn.sourceforge.net/octave/?rev=10532&view=rev
Author: paramaniac
Date: 2012-05-29 15:51:05 +0000 (Tue, 29 May 2012)
Log Message:
-----------
control: apply makefile patch by Olaf Till, add some remarks about temporary nature of the included reference lapack distribution
Modified Paths:
--------------
trunk/octave-forge/main/control/src/Makefile
Modified: trunk/octave-forge/main/control/src/Makefile
===================================================================
--- trunk/octave-forge/main/control/src/Makefile 2012-05-28 11:31:45 UTC (rev 10531)
+++ trunk/octave-forge/main/control/src/Makefile 2012-05-29 15:51:05 UTC (rev 10532)
@@ -1,7 +1,9 @@
-# LAPACK_LIBS := $(shell mkoctfile -p LAPACK_LIBS)
-# BLAS_LIBS := $(shell mkoctfile -p BLAS_LIBS)
-FLIBS := $(shell mkoctfile -p FLIBS)
+MKOCTFILE ?= mkoctfile
+# LAPACK_LIBS := $(shell $(MKOCTFILE) -p LAPACK_LIBS)
+# BLAS_LIBS := $(shell $(MKOCTFILE) -p BLAS_LIBS)
+FLIBS := $(shell $(MKOCTFILE) -p FLIBS)
+
all: control_slicot_functions.oct \
is_real_scalar.oct \
is_real_vector.oct \
@@ -20,36 +22,42 @@
mv slicot/src/*.f .
mv slicot/src_aux/*.f .
cp TG04BX.fortran TG04BX.f
- mkoctfile -c *.f
+ $(MKOCTFILE) -c *.f
ar -rc slicotlibrary.a *.o
rm -rf *.o *.f slicot
+# reference lapack is just included for debugging purposes
+# it will be removed again before an official release of the control package
lapacklibrary.a: lapack-3.4.1.tgz
tar -xzf lapack-3.4.1.tgz
mv lapack-3.4.1/BLAS/SRC/*.f .
mv lapack-3.4.1/INSTALL/*.f .
mv lapack-3.4.1/SRC/*.f .
- mkoctfile -c *.f
+ $(MKOCTFILE) -c *.f
ar -rc lapacklibrary.a *.o
rm -rf *.o *.f lapack-3.4.1
# slicot functions
control_slicot_functions.oct: control_slicot_functions.cc slicotlibrary.a lapacklibrary.a
- mkoctfile control_slicot_functions.cc common.cc slicotlibrary.a lapacklibrary.a \
+ $(MKOCTFILE) control_slicot_functions.cc common.cc slicotlibrary.a lapacklibrary.a \
${FLIBS}
+# slicot functions
+#control_slicot_functions.oct: control_slicot_functions.cc slicotlibrary.a
+# $(MKOCTFILE) control_slicot_functions.cc common.cc slicotlibrary.a \
+# ${LAPACK_LIBS} ${BLAS_LIBS} ${FLIBS}
# helpers
is_real_scalar.oct: is_real_scalar.cc
- mkoctfile is_real_scalar.cc
+ $(MKOCTFILE) is_real_scalar.cc
is_real_vector.oct: is_real_vector.cc
- mkoctfile is_real_vector.cc
+ $(MKOCTFILE) is_real_vector.cc
is_real_matrix.oct: is_real_matrix.cc
- mkoctfile is_real_matrix.cc
+ $(MKOCTFILE) is_real_matrix.cc
is_real_square_matrix.oct: is_real_square_matrix.cc
- mkoctfile is_real_square_matrix.cc
+ $(MKOCTFILE) is_real_square_matrix.cc
clean:
rm -rf *.o core octave-core *.oct *~ *.f slicot lapack-3.4.1
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <car...@us...> - 2012-05-28 11:31:56
|
Revision: 10531
http://octave.svn.sourceforge.net/octave/?rev=10531&view=rev
Author: carandraug
Date: 2012-05-28 11:31:45 +0000 (Mon, 28 May 2012)
Log Message:
-----------
GeSHi highlighting:
* comments on file
* added URL
* made keywords match case insensitive
* improved pkg regexp
* swapped keyword #1 with #4 to to workaround GeSHi bug
* removed octave deprecated functions
* removed attempt to highlight complex numbers
* fixed indentation problem
Modified Paths:
--------------
trunk/octave-forge/admin/langfiles/octave.php
Modified: trunk/octave-forge/admin/langfiles/octave.php
===================================================================
--- trunk/octave-forge/admin/langfiles/octave.php 2012-05-28 10:35:46 UTC (rev 10530)
+++ trunk/octave-forge/admin/langfiles/octave.php 2012-05-28 11:31:45 UTC (rev 10531)
@@ -6,17 +6,20 @@
* Juan Pablo Carbajal (car...@if...)
* Copyright: (c) 2012 Carnë Draug
* (c) 2012 Juan Pablo Carbajal
- * Release Version: 1.0.0
+ * Release Version: 1.0.8.11
* Date Started: 2012/05/22
*
* GNU/Octave M-file language file for GeSHi.
*
+ * This file was heavily based on octave.lang from gtksourceview. If bugs are
+ * found and/or fixed on this file, please send them to the gtksourceview
+ * project or e-mail them to this file authors. Thanks in advance
+ *
* CHANGES
* -------
- * 2012/05/22 (1.0.0)
+ * 2012/05/22 (1.0.8.11)
* - First Release
*
- *
*************************************************************************************
*
* This file is part of GeSHi.
@@ -39,20 +42,21 @@
$language_data = array (
'LANG_NAME' => 'GNU/Octave',
- 'COMMENT_SINGLE' => array(1=> '#', 2 => '%'),
- // we really can't use this since these characters need to be alone (or with
- // whitespace) on a line. We will have to use regexp
+ 'COMMENT_SINGLE' => array(1 => '#', 2 => '%'),
+ // we can't use COMMENT_MULTI since start and end of block comments need to
+ // be alone on the line (optionally, with whitespace). See COMMENT_REGEXP
'COMMENT_MULTI' => array(),
+ // we can't use QUOTEMARKS, not even HARDQUOTE, see COMMENT_REGEXP
'QUOTEMARKS' => array(),
+ 'ESCAPE_CHAR' => '',
'COMMENT_REGEXP' => array(
- // Single quote strings: we really can't use QUOTEMARKS here since new
+ // Single quote strings: we can't use QUOTEMARKS here since new
// lines will break the string. Plus, single quote strings do not even
// allow for continuation markers, only double quote strings allow it.
// Also, to do not misdetect the transpose operator ' as the start of a
// string we assert to not follow a variable name (letters, digits and
// underscores) or a closing bracket (round, square or curly) or a dot
// (to form the array transpose operator ".'" ).
- // see the source of octave.lang of gtksourceview
3 => "/(?<![0-9a-zA-Z_\)\]}\.])'.*?'/",
// Double quote strings: we also can't use QUOTEMARKS here (see single
// line quotes). However, with double quote strings both \ and ... can
@@ -64,6 +68,12 @@
// comments on the first and last line of source) and make . also match
// a newline
5 => "/^\s*?[%#]{\s*?$.*?^\s*?[%#]}\s*?$/ms",
+ // Packaging system: comes here so that pkg can also be used in the
+ // function form. The list of pkg commands is optional to the match so
+ // that at least pkg is highlighted if new commands are implemented
+ 6 => "/\bpkg(?!\s*\()\s+((un)?install|update|(un)?load|list|(global|local)_list|describe|prefix|(re)?build)?\b/",
+ // Function handles
+ 7 => "/@([a-z_][a-z1-9_]*)?/i",
),
'NUMBERS' =>
GESHI_NUMBER_INT_BASIC |
@@ -71,12 +81,15 @@
GESHI_NUMBER_HEX_PREFIX |
GESHI_NUMBER_FLT_SCI_ZERO,
'CASE_KEYWORDS' => GESHI_CAPS_NO_CHANGE,
- 'ESCAPE_CHAR' => '',
'KEYWORDS' => array(
- // Data types
+ // Reserved words
1 => array(
- 'cell', 'char', 'double', 'uint8', 'uint16', 'uint32', 'uint64',
- 'int8','int16', 'int32', 'int64', 'logical', 'single', 'struct'
+ 'break', 'case', 'catch', 'continue', 'do', 'else', 'elseif', 'end',
+ 'end_try_catch', 'end_unwind_protect', 'endfor', 'endfunction',
+ 'endif', 'endparfor', 'endswitch', 'endwhile', 'for', 'function',
+ 'if', 'otherwise', 'parfor', 'return',
+ 'switch', 'try', 'until', 'unwind_protect',
+ 'unwind_protect_cleanup', 'varargin', 'varargout', 'while'
),
// Storage type
2 => array(
@@ -84,14 +97,10 @@
),
// Internal variable
3 => array('ans'),
- // Reserved words
+ // Data types
4 => array(
- 'break', 'case', 'catch', 'continue', 'do', 'else', 'elseif', 'end',
- 'end_try_catch', 'end_unwind_protect', 'endfor', 'endfunction',
- 'endif', 'endparfor', 'endswitch', 'endwhile', 'for', 'function',
- 'if', 'otherwise', 'parfor', 'return',
- 'switch', 'try', 'until', 'unwind_protect',
- 'unwind_protect_cleanup', 'varargin', 'varargout', 'while'
+ 'cell', 'char', 'double', 'uint8', 'uint16', 'uint32', 'uint64',
+ 'int8','int16', 'int32', 'int64', 'logical', 'single', 'struct'
),
// Built in
5 => array(
@@ -187,38 +196,38 @@
'acscd', 'acsch', 'addpref', 'addtodate', 'allchild', 'amd',
'ancestor', 'anova', 'arch_fit', 'arch_rnd', 'arch_test',
'area', 'arma_rnd', 'asctime', 'asec', 'asecd', 'asech', 'asind',
- 'assert', 'atand', 'autocor', 'autocov', 'autoreg_matrix', 'autumn',
+ 'assert', 'atand', 'autoreg_matrix', 'autumn',
'axes', 'axis', 'balance', 'bar', 'barh', 'bartlett', 'bartlett_test',
- 'base2dec', 'beep', 'bessel', 'besselj', 'beta', 'betacdf', 'betai',
+ 'base2dec', 'beep', 'bessel', 'besselj', 'beta', 'betacdf',
'betainc', 'betainv', 'betaln', 'betapdf', 'betarnd', 'bicg',
'bicgstab', 'bicubic', 'bin2dec', 'bincoeff', 'binocdf', 'binoinv',
'binopdf', 'binornd', 'bitcmp', 'bitget', 'bitset', 'blackman',
'blanks', 'blkdiag', 'bone', 'box', 'brighten', 'bsxfun',
'bug_report', 'bunzip2', 'bzip2', 'calendar', 'cart2pol', 'cart2sph',
'cast', 'cauchy_cdf', 'cauchy_inv', 'cauchy_pdf', 'cauchy_rnd',
- 'caxis', 'ccolamd', 'cell2mat', 'celldisp', 'cellfun', 'cellidx',
+ 'caxis', 'ccolamd', 'cell2mat', 'celldisp', 'cellfun',
'center', 'cgs', 'chi2cdf', 'chi2inv', 'chi2pdf', 'chi2rnd',
'chisquare_test_homogeneity', 'chisquare_test_independence', 'chol',
- 'chop', 'circshift', 'cla', 'clabel', 'clf', 'clg', 'clock',
+ 'chop', 'circshift', 'cla', 'clabel', 'clf', 'clock',
'cloglog', 'close', 'closereq', 'colamd', 'colloc', 'colon',
'colorbar', 'colormap', 'colperm', 'colstyle', 'comet', 'comet3',
'comma', 'common_size', 'commutation_matrix', 'compan',
'compare_versions', 'compass', 'computer', 'cond', 'condest',
'contour', 'contour3', 'contourc', 'contourf', 'contrast', 'conv',
- 'conv2', 'convhull', 'convhulln', 'cool', 'copper', 'copyfile', 'cor',
- 'cor_test', 'corr', 'corrcoef', 'cosd', 'cot', 'cotd', 'coth', 'cov',
- 'cplxpair', 'cquad', 'cross', 'csc', 'cscd', 'csch', 'cstrcat',
- 'csvread', 'csvwrite', 'ctime', 'cumtrapz', 'curl', 'cut', 'cylinder',
+ 'conv2', 'convhull', 'convhulln', 'cool', 'copper', 'copyfile',
+ 'cor_test', 'corr', 'cosd', 'cot', 'cotd', 'coth', 'cov',
+ 'cplxpair', 'cross', 'csc', 'cscd', 'csch', 'cstrcat',
+ 'csvread', 'csvwrite', 'ctime', 'cumtrapz', 'curl', 'cylinder',
'daspect', 'daspk', 'dasrt', 'dassl', 'date', 'datenum', 'datestr',
'datetick', 'datevec', 'dblquad', 'deal', 'deblank', 'debug',
'dec2base', 'dec2bin', 'dec2hex', 'deconv', 'del2', 'delaunay',
'delaunay3', 'delaunayn', 'delete', 'demo', 'det', 'detrend',
'diffpara', 'diffuse', 'dir', 'discrete_cdf', 'discrete_inv',
- 'discrete_pdf', 'discrete_rnd', 'dispatch', 'display', 'divergence',
+ 'discrete_pdf', 'discrete_rnd', 'display', 'divergence',
'dlmread', 'dlmwrite', 'dmperm', 'doc', 'dos', 'dot', 'dsearch',
'dsearchn', 'dump_prefs', 'duplication_matrix', 'durbinlevinson',
'edit', 'eig', 'eigs', 'ellipsoid', 'empirical_cdf', 'empirical_inv',
- 'empirical_pdf', 'empirical_rnd', 'eomday', 'error_text', 'errorbar',
+ 'empirical_pdf', 'empirical_rnd', 'eomday', 'errorbar',
'etime', 'etreeplot', 'example', 'expcdf', 'expinv', 'expm', 'exppdf',
'exprnd', 'ezcontour', 'ezcontourf', 'ezmesh', 'ezmeshc', 'ezplot',
'ezplot3', 'ezpolar', 'ezsurf', 'ezsurfc', 'f_test_regression',
@@ -227,13 +236,13 @@
'fileattrib', 'fileparts', 'fileread', 'fill', 'filter', 'filter2',
'find', 'findall', 'findobj', 'findstr', 'finv', 'flag', 'flipdim',
'fliplr', 'flipud', 'fminbnd', 'fminunc', 'fpdf', 'fplot',
- 'fractdiff', 'freqz', 'freqz_plot', 'frnd', 'fsolve', 'fstat',
- 'fullfile', 'fzero', 'gamcdf', 'gaminv', 'gammai', 'gammainc',
+ 'fractdiff', 'freqz', 'freqz_plot', 'frnd', 'fsolve',
+ 'fullfile', 'fzero', 'gamcdf', 'gaminv', 'gammainc',
'gampdf', 'gamrnd', 'gca', 'gcbf', 'gcbo', 'gcd', 'gcf',
'gen_doc_cache', 'genvarname', 'geocdf', 'geoinv', 'geopdf', 'geornd',
'get_first_help_sentence', 'getappdata', 'getfield', 'getgrent',
'getpref', 'getpwent', 'getrusage', 'ginput', 'givens', 'glpk',
- 'glpkmex', 'gls', 'gmap40', 'gmres', 'gnuplot_binary', 'gplot',
+ 'gls', 'gmap40', 'gmres', 'gnuplot_binary', 'gplot',
'gradient', 'graphics_toolkit', 'gray', 'gray2ind', 'grid',
'griddata', 'griddata3', 'griddatan', 'gtext', 'guidata',
'guihandles', 'gunzip', 'gzip', 'hadamard', 'hamming', 'hankel',
@@ -244,16 +253,16 @@
'imagesc', 'imfinfo', 'imread', 'imshow', 'imwrite', 'ind2gray',
'ind2rgb', 'index', 'info', 'inpolygon', 'inputname', 'int2str',
'interp1', 'interp1q', 'interp2', 'interp3', 'interpft', 'interpn',
- 'intersect', 'intwarning', 'inv', 'invhilb', 'iqr',
- 'is_duplicate_entry', 'is_global', 'is_leap_year', 'is_valid_file_id',
+ 'intersect', 'inv', 'invhilb', 'iqr',
+ 'is_leap_year', 'is_valid_file_id',
'isa', 'isappdata', 'iscolumn', 'isdefinite', 'isdeployed', 'isdir',
'isequal', 'isequalwithequalnans', 'isfigure', 'ishermitian',
'ishghandle', 'ishold', 'isletter', 'ismac', 'ismember', 'isocolors',
'isonormals', 'isosurface', 'ispc', 'ispref', 'isprime', 'isprop',
- 'isrow', 'isscalar', 'issquare', 'isstr', 'isstrprop', 'issymmetric',
+ 'isrow', 'isscalar', 'issquare', 'isstrprop', 'issymmetric',
'isunix', 'isvector', 'jet', 'kendall', 'kolmogorov_smirnov_cdf',
'kolmogorov_smirnov_test', 'kolmogorov_smirnov_test_2', 'kron',
- 'kruskal_wallis_test', 'krylov', 'krylovb', 'kurtosis', 'laplace_cdf',
+ 'kruskal_wallis_test', 'krylov', 'kurtosis', 'laplace_cdf',
'laplace_inv', 'laplace_pdf', 'laplace_rnd', 'lcm', 'legend',
'legendre', 'license', 'lin2mu', 'line', 'linkprop', 'list_primes',
'loadaudio', 'loadobj', 'logistic_cdf', 'logistic_inv',
@@ -272,11 +281,11 @@
'ocean', 'ols', 'onenormest', 'optimget', 'optimset', 'orderfields',
'orient', 'orth', 'pack', 'paren', 'pareto', 'parseparams', 'pascal',
'patch', 'pathdef', 'pbaspect', 'pcg', 'pchip', 'pcolor', 'pcr',
- 'peaks', 'periodogram', 'perl', 'perms', 'perror', 'pie', 'pie3',
- 'pink', 'pinv', 'planerot', 'playaudio', 'plot', 'plot3',
+ 'peaks', 'periodogram', 'perl', 'perms', 'pie', 'pie3',
+ 'pink', 'pinv', 'pkg', 'planerot', 'playaudio', 'plot', 'plot3',
'plotmatrix', 'plotyy', 'poisscdf', 'poissinv', 'poisspdf',
'poissrnd', 'pol2cart', 'polar', 'poly', 'polyaffine', 'polyarea',
- 'polyder', 'polyderiv', 'polyfit', 'polygcd', 'polyint', 'polyout',
+ 'polyder', 'polyfit', 'polygcd', 'polyint', 'polyout',
'polyreduce', 'polyval', 'polyvalm', 'postpad', 'pow2', 'powerset',
'ppder', 'ppint', 'ppjumps', 'ppplot', 'ppval', 'pqpnonneg',
'prctile', 'prepad', 'primes', 'print', 'printAllBuiltins',
@@ -286,15 +295,15 @@
'qz', 'qzhess', 'rainbow', 'rand', 'randi', 'range', 'rank', 'ranks',
'rat', 'rcond', 'reallog', 'realpow', 'realsqrt', 'record',
'rectangle', 'rectint', 'recycle', 'refresh', 'refreshdata', 'regexp',
- 'regexptranslate', 'replot', 'repmat', 'residue', 'rgb2hsv',
+ 'regexptranslate', 'repmat', 'residue', 'rgb2hsv',
'rgb2ind', 'rgb2ntsc', 'ribbon', 'rindex', 'rmappdata', 'rmpref',
'roots', 'rose', 'rosser', 'rot90', 'rotdim', 'rref', 'run',
'run_count', 'run_test', 'rundemos', 'runlength', 'runtests',
- 'saveas', 'saveaudio', 'saveimage', 'saveobj', 'savepath', 'scatter',
+ 'saveas', 'saveaudio', 'saveobj', 'savepath', 'scatter',
'scatter3', 'schur', 'sec', 'secd', 'sech', 'semicolon', 'semilogx',
'semilogxerr', 'semilogy', 'semilogyerr', 'setappdata', 'setaudio',
- 'setdiff', 'setfield', 'setpref', 'setstr', 'setxor', 'shading',
- 'shell_cmd', 'shg', 'shift', 'shiftdim', 'sign_test', 'sinc', 'sind',
+ 'setdiff', 'setfield', 'setpref', 'setxor', 'shading',
+ 'shg', 'shift', 'shiftdim', 'sign_test', 'sinc', 'sind',
'sinetone', 'sinewave', 'skewness', 'slice', 'sombrero', 'sortrows',
'spaugment', 'spconvert', 'spdiags', 'spearman', 'spectral_adf',
'spectral_xdf', 'specular', 'speed', 'spencer', 'speye', 'spfun',
@@ -302,12 +311,12 @@
'sprand', 'sprandn', 'sprandsym', 'spring', 'spstats', 'spy', 'sqp',
'sqrtm', 'stairs', 'statistics', 'std', 'stdnormal_cdf',
'stdnormal_inv', 'stdnormal_pdf', 'stdnormal_rnd', 'stem', 'stem3',
- 'stft', 'str2double', 'str2num', 'strcat', 'strchr', 'strerror',
+ 'stft', 'str2double', 'str2num', 'strcat', 'strchr',
'strfind', 'strjust', 'strmatch', 'strread', 'strsplit', 'strtok',
- 'strtrim', 'strtrunc', 'structfun', 'studentize', 'sub2ind',
+ 'strtrim', 'strtrunc', 'structfun', 'sub2ind',
'subplot', 'subsindex', 'subspace', 'substr', 'substruct', 'summer',
'surf', 'surface', 'surfc', 'surfl', 'surfnorm', 'svd', 'svds',
- 'swapbytes', 'syl', 'sylvester_matrix', 'symbfact', 'symrcm',
+ 'swapbytes', 'syl', 'symbfact', 'symrcm',
'symvar', 'synthesis', 't_test', 't_test_2', 't_test_regression',
'table', 'tand', 'tar', 'tcdf', 'tempdir', 'tempname', 'test', 'text',
'textread', 'textscan', 'time', 'tinv', 'title', 'toeplitz', 'tpdf',
@@ -319,11 +328,11 @@
'unidcdf', 'unidinv', 'unidpdf', 'unidrnd', 'unifcdf', 'unifinv',
'unifpdf', 'unifrnd', 'unimplemented', 'union', 'unique', 'unix',
'unmkpp', 'unpack', 'untabify', 'untar', 'unwrap', 'unzip',
- 'urlwrite', 'usejava', 'validatestring', 'values', 'vander', 'var',
+ 'urlwrite', 'usejava', 'validatestring', 'vander', 'var',
'var_test', 'vech', 'ver', 'version', 'view', 'voronoi', 'voronoin',
'waitbar', 'waitforbuttonpress', 'warning_ids', 'wavread', 'wavwrite',
- 'wblcdf', 'wblinv', 'wblpdf', 'wblrnd', 'weekday', 'weibcdf',
- 'weibinv', 'weibpdf', 'weibrnd', 'welch_test', 'what', 'which',
+ 'wblcdf', 'wblinv', 'wblpdf', 'wblrnd', 'weekday',
+ 'welch_test', 'what', 'which',
'white', 'whitebg', 'wienrnd', 'wilcoxon_test', 'wilkinson', 'winter',
'xlabel', 'xlim', 'xor', 'ylabel', 'ylim', 'yulewalker', 'z_test',
'z_test_2', 'zip', 'zlabel', 'zlim', 'zscore', 'airy', 'arrayfun',
@@ -363,7 +372,7 @@
// Private Octave functions
8 => array(
'__all_opts__', '__contourc__', '__delaunayn__', '__dispatch__',
- '__dsearchn__', '__error_text__', '__finish__', '__fltk_uigetfile__',
+ '__dsearchn__', '__finish__', '__fltk_uigetfile__',
'__glpk__', '__gnuplot_drawnow__', '__init_fltk__',
'__init_gnuplot__', '__lin_interpn__', '__magick_read__',
'__makeinfo__', '__pchip_deriv__', '__plt_get_axis_arg__', '__qp__',
@@ -387,7 +396,10 @@
'WIFSTOPPED', 'WNOHANG', 'WSTOPSIG', 'WTERMSIG', 'WUNTRACED'
),
// Constant functions
- 10 => array ('e','eps','inf','nan','NA','pi','i','j','true','false'),
+ 10 => array (
+ 'e', 'eps', 'inf', 'Inf', 'nan', 'NaN', 'NA', 'pi', 'i', 'I', 'j',
+ 'J', 'true', 'false'
+ ),
),
'SYMBOLS' => array(
// Comparison & logical
@@ -408,121 +420,86 @@
4 => array(':'),
// Delimiters
5 => array(',', '...', ';'),
- ),
+ ),
'CASE_SENSITIVE' => array(
- GESHI_COMMENTS => false,
- 1 => false,
- 2 => false,
- 3 => false,
- 4 => false,
- 5 => false,
- 6 => false,
- 7 => false,
- 8 => false,
- 9 => false,
- 10 => false
+ GESHI_COMMENTS => true,
+ 1 => true,
+ 2 => true,
+ 3 => true,
+ 4 => true,
+ 5 => true,
+ 6 => true,
+ 7 => true,
+ 8 => true,
+ 9 => true,
+ 10 => true,
),
'URLS' => array(
1 => '',
2 => '',
3 => '',
- 4 => '',
- 5 => '',
- 6 => '',
+ 4 => 'http://octave.sourceforge.net/octave/function/{FNAME}.html',
+ 5 => 'http://octave.sourceforge.net/octave/function/{FNAME}.html',
+ 6 => 'http://octave.sourceforge.net/octave/function/{FNAME}.html',
7 => '',
8 => '',
- 9 => '',
- 10 => ''
+ 9 => 'http://octave.sourceforge.net/octave/function/{FNAME}.html',
+ 10 => 'http://octave.sourceforge.net/octave/function/{FNAME}.html',
),
'OOLANG' => true,
'OBJECT_SPLITTERS' => array(
1 => '.',
- 2 => '::'
),
- 'REGEXPS' => array(
- // Complex numbers
- 1 => array(
- GESHI_SEARCH => '([+-]?\d+(\.\d+)?(e\d+)?\s*[+-]\s*([ij]\*?)?\d+(\.\d+)?(e\d+)?(\*?[ij])?)',
- GESHI_REPLACE => '\\1',
- GESHI_MODIFIERS => 'si',
- GESHI_BEFORE => '',
- GESHI_AFTER => ''
- ),
- //Function handle
- 2 => array(
- GESHI_SEARCH => '(@([A-Za-z_][A-Za-z1-9_]*)?)',
- GESHI_REPLACE => '\\1',
- GESHI_MODIFIERS => '',
- GESHI_BEFORE => '',
- GESHI_AFTER => ''
- ),
- //Packaging system
- 3 => array(
- GESHI_SEARCH => '((\b)pkg(?!(\s)*\()(\s)+(((un)?install|(un)?load|list|(global|local)_list|describe|prefix|(re)?build)(\b))?)',
- GESHI_REPLACE => '\\1',
- GESHI_MODIFIERS => '',
- GESHI_BEFORE => '',
- GESHI_AFTER => ''
- )
- ),
+ 'REGEXPS' => array(),
'STRICT_MODE_APPLIES' => GESHI_NEVER,
'SCRIPT_DELIMITERS' => array(),
'HIGHLIGHT_STRICT_BLOCK' => array(),
'STYLES' => array(
- 'KEYWORDS' => array(
- 1 => 'color: #2E8B57; font-weight:bold;', // Data types
- 2 => 'color: #2E8B57;', // Storage type
- 3 => 'color: #0000FF; font-weight:bold;', // Internal variable
- 4 => 'color: #990000; font-weight:bold;', // Reserved words
- 5 => 'color: #008A8C; font-weight:bold;', // Built-in
- 6 => 'color: #008A8C;', // Octave functions
- 9 => 'color: #000000; font-weight:bold;', // Builtin Global Variables
- 10 => 'color: #008A8C; font-weight:bold;' // Constant functions
- ),
'COMMENTS' => array(
- 1 => 'color: #0000FF; font-style: italic;', // single quote strings
- 2 => 'color: #0000FF; font-style: italic;', // double quote strings
+ 1 => 'color: #0000FF; font-style: italic;', // comments with #
+ 2 => 'color: #0000FF; font-style: italic;', // comments with %
3 => 'color: #FF00FF; font-style: italic;', // single quote strings
4 => 'color: #FF00FF; font-style: italic;', // double quote strings
5 => 'color: #0000FF; font-style: italic;', // block comments
- 'MULTI' => 'color: #0000FF; font-style: italic;'
+ 6 => 'color: #996600; font-weight:bold;', // packaging system
+ 7 => 'color: #006600; font-weight:bold;', // function handles
+ 'MULTI' => 'color: #0000FF; font-style: italic;',
),
- 'ESCAPE_CHAR' => array(
- 0 => ''
+ 'KEYWORDS' => array(
+ 1 => 'color: #990000; font-weight:bold;', // Reserved words
+ 2 => 'color: #2E8B57;', // Storage type
+ 3 => 'color: #0000FF; font-weight:bold;', // Internal variable
+ 4 => 'color: #2E8B57; font-weight:bold;', // Data types
+ 5 => 'color: #008A8C; font-weight:bold;', // Built-in
+ 6 => 'color: #008A8C;', // Octave functions
+ 9 => 'color: #000000; font-weight:bold;', // Builtin Global Variables
+ 10 => 'color: #008A8C; font-weight:bold;', // Constant functions
),
+ 'ESCAPE_CHAR' => array(),
'BRACKETS' => array(
- 0 => 'color: #080;'
+ 0 => 'color: #080;',
),
'STRINGS' => array(
- 0 => 'color: #A020F0;'
+ // strings were specified on the COMMENT_REGEXP section
),
'NUMBERS' => array(
0 => 'color: #cc66cc;',
GESHI_NUMBER_OCT_PREFIX => 'color: #208080;',
GESHI_NUMBER_HEX_PREFIX => 'color: #208080;',
- GESHI_NUMBER_FLT_SCI_ZERO => 'color:#800080;'
+ GESHI_NUMBER_FLT_SCI_ZERO => 'color:#800080;',
),
- 'METHODS' => array(
- 1 => '',
- 2 => ''
- ),
+ 'METHODS' => array(),
'SYMBOLS' => array(
- 0 => 'color: #FF9696; font-weight:bold;', // Comparison & logical
- 1 => 'color: #CC0000; font-weight:bold;', // Aritmethical
- 2 => 'color: #993333; font-weight:bold;', // Elementwise arithmetical
- 3 => 'color: #FF0000; font-weight:bold;', // Arithmetical & assignation
- 4 => 'color: #33F;', // Indexer
- 5 => 'color: #33F;' // Delimiters
+ 0 => 'color: #FF9696; font-weight:bold;', // Comparison & logical
+ 1 => 'color: #CC0000; font-weight:bold;', // Aritmethical
+ 2 => 'color: #993333; font-weight:bold;', // Elementwise arithmetical
+ 3 => 'color: #FF0000; font-weight:bold;', // Arithmetical & assignation
+ 4 => 'color: #33F;', // Indexer
+ 5 => 'color: #33F;', // Delimiters
),
- 'REGEXPS' => array(
- 1 => 'color: #9966FF; font-weight:bold;', // Complex number
- 2 => 'color: #006600; font-weight:bold;', //Function handle
- 3 => 'color: #996600; font-weight:bold;', // Package manager
- ),
- 'SCRIPT' => array(
- 0 => ''
- )
- )
+ 'REGEXPS' => array(),
+ 'SCRIPT' => array(),
+ ),
);
?>
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <car...@us...> - 2012-05-28 10:35:53
|
Revision: 10530
http://octave.svn.sourceforge.net/octave/?rev=10530&view=rev
Author: carandraug
Date: 2012-05-28 10:35:46 +0000 (Mon, 28 May 2012)
Log Message:
-----------
bessel_zero_Jnu: add binding for the bessel_zero_Jnu function. Patch by Rafael Laboissiere <ra...@la...>
Modified Paths:
--------------
trunk/octave-forge/main/gsl/INDEX
trunk/octave-forge/main/gsl/PKG_ADD
trunk/octave-forge/main/gsl/src/buildgsl_sf.sh
Added Paths:
-----------
trunk/octave-forge/main/gsl/src/double_int_to_double.cc.template
Modified: trunk/octave-forge/main/gsl/INDEX
===================================================================
--- trunk/octave-forge/main/gsl/INDEX 2012-05-28 07:55:30 UTC (rev 10529)
+++ trunk/octave-forge/main/gsl/INDEX 2012-05-28 10:35:46 UTC (rev 10530)
@@ -36,6 +36,7 @@
bessel_yl
bessel_zero_J0
bessel_zero_J1
+ bessel_zero_Jnu
beta_gsl
clausen
conicalP_0
Modified: trunk/octave-forge/main/gsl/PKG_ADD
===================================================================
--- trunk/octave-forge/main/gsl/PKG_ADD 2012-05-28 07:55:30 UTC (rev 10529)
+++ trunk/octave-forge/main/gsl/PKG_ADD 2012-05-28 10:35:46 UTC (rev 10530)
@@ -30,6 +30,7 @@
autoload ("bessel_Ynu", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
autoload ("bessel_zero_J0", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
autoload ("bessel_zero_J1", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
+autoload ("bessel_zero_Jnu", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
autoload ("beta_gsl", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
autoload ("Chi", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
autoload ("Ci", fullfile (fileparts (mfilename ("fullpath")), "gsl_sf.oct"));
Modified: trunk/octave-forge/main/gsl/src/buildgsl_sf.sh
===================================================================
--- trunk/octave-forge/main/gsl/src/buildgsl_sf.sh 2012-05-28 07:55:30 UTC (rev 10529)
+++ trunk/octave-forge/main/gsl/src/buildgsl_sf.sh 2012-05-28 10:35:46 UTC (rev 10530)
@@ -1183,6 +1183,17 @@
./replace_template.sh int_int_double_to_double.cc.template >> gsl_sf.cc
+# (double, int) to double
+
+export octave_name=bessel_zero_Jnu
+export funcname=gsl_sf_bessel_zero_Jnu
+cat << \EOF > docstring.txt
+These routines compute the location of the n-th positive zero of the
+Bessel function J_x().
+EOF
+./replace_template.sh double_int_to_double.cc.template >> gsl_sf.cc
+
+
export octave_name=hyperg_U
export funcname=gsl_sf_hyperg_U
cat << \EOF > docstring.txt
Added: trunk/octave-forge/main/gsl/src/double_int_to_double.cc.template
===================================================================
--- trunk/octave-forge/main/gsl/src/double_int_to_double.cc.template (rev 0)
+++ trunk/octave-forge/main/gsl/src/double_int_to_double.cc.template 2012-05-28 10:35:46 UTC (rev 10530)
@@ -0,0 +1,112 @@
+DEFUN_DLD(GSL_OCTAVE_NAME, args, nargout, "\
+ -*- texinfo -*-\n\
+@deftypefn {Loadable Function} {@var{y} =} GSL_OCTAVE_NAME (@var{x}, @var{n})\n\
+@deftypefnx {Loadable Function} {[@var{y}, @var{err}] =} GSL_OCTAVE_NAME (@dots{})\n\
+\n\
+GSL_FUNC_DOCSTRING
+\n\
+@var{err} contains an estimate of the absolute error in the value @var{y}.\n\
+\n\
+This function is from the GNU Scientific Library,\n\
+see @url{http://www.gnu.org/software/gsl/} for documentation.\n\
+@end deftypefn\n\
+")
+{
+ int i;
+ dim_vector dv;
+
+ gsl_set_error_handler (octave_gsl_errorhandler);
+
+ if(args.length() != 2) {
+ print_usage ();
+ return octave_value();
+ }
+ if(!args(0).is_real_type() || !args(1).is_real_type()) {
+ error("The arguments must be real.");
+ print_usage ();
+ return octave_value();
+ }
+
+ // Nice combinatorial explosion here
+ NDArray x = args(0).array_value();
+ NDArray n = args(1).array_value();
+ if(n.length() == x.length()) {
+ dv = x.dims();
+ NDArray y(dv);
+ int len = x.length();
+ if(nargout < 2) {
+ for(i = 0; i < len; i++) {
+ y.xelem(i) = GSL_FUNC_NAME (x.xelem(i),
+ static_cast<int>(n.xelem(i)));
+ }
+ return octave_value(y);
+ } else {
+ NDArray err(dv);
+ gsl_sf_result result;
+ octave_value_list retval;
+ for(i = 0; i < len; i++) {
+ GSL_FUNC_NAME_e (x.xelem(i),
+ static_cast<int>(n.xelem(i)), &result);
+ y.xelem(i) = result.val;
+ err.xelem(i) = result.err;
+ }
+ retval(1) = octave_value(err);
+ retval(0) = octave_value(y);
+ return retval;
+ }
+ } else if(n.length() == 1) {
+ dv = x.dims();
+ NDArray y(dv);
+ int len = x.length();
+ int nint = static_cast<int>(n.xelem(0));
+ if(nargout < 2) {
+ for(i = 0; i < len; i++) {
+ y.xelem(i) = GSL_FUNC_NAME (x.xelem(i), nint);
+ }
+ return octave_value(y);
+ } else {
+ NDArray err(dv);
+ gsl_sf_result result;
+ octave_value_list retval;
+ for(i = 0; i < len; i++) {
+ GSL_FUNC_NAME_e (x.xelem(i), nint, &result);
+ y.xelem(i) = result.val;
+ err.xelem(i) = result.err;
+ }
+ retval(1) = octave_value(err);
+ retval(0) = octave_value(y);
+ return retval;
+ }
+ } else if(x.length() == 1) {
+ dv = n.dims();
+ NDArray y(dv);
+ int len = n.length();
+ double xdouble = x.xelem(0);
+ if(nargout < 2) {
+ for(i = 0; i < len; i++) {
+ y.xelem(i) = GSL_FUNC_NAME (xdouble,
+ static_cast<int>(n.xelem(i)));
+ }
+ return octave_value(y);
+ } else {
+ NDArray err(dv);
+ gsl_sf_result result;
+ octave_value_list retval;
+ for(i = 0; i < len; i++) {
+ GSL_FUNC_NAME_e (xdouble,
+ static_cast<int>(n.xelem(i)), &result);
+ y.xelem(i) = result.val;
+ err.xelem(i) = result.err;
+ }
+ retval(1) = octave_value(err);
+ retval(0) = octave_value(y);
+ return retval;
+ }
+ } else {
+ error("First and second argument must either have the same size, or one of them must be scalar.");
+ print_usage ();
+ }
+
+ return octave_value();
+
+}
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <cd...@us...> - 2012-05-28 07:55:36
|
Revision: 10529
http://octave.svn.sourceforge.net/octave/?rev=10529&view=rev
Author: cdf
Date: 2012-05-28 07:55:30 +0000 (Mon, 28 May 2012)
Log Message:
-----------
makefile change to support linking apple veclib
Modified Paths:
--------------
trunk/octave-forge/main/odepkg/src/Makefile
Modified: trunk/octave-forge/main/odepkg/src/Makefile
===================================================================
--- trunk/octave-forge/main/odepkg/src/Makefile 2012-05-27 20:31:31 UTC (rev 10528)
+++ trunk/octave-forge/main/odepkg/src/Makefile 2012-05-28 07:55:30 UTC (rev 10529)
@@ -9,7 +9,9 @@
# MKOCTFILE = $(MKOCTFILE) has already been defined in Makeconf
# MKMEXFILE = $(MKOCTFILE) --mex
-MKOCTFILE = mkoctfile
+MKOCTFILE ?= mkoctfile
+
+FFLAGS := $(shell $(MKOCTFILE) -p FFLAGS)
ifeq (gfortran,$(findstring gfortran,$(F77)))
MKF77FILE = FFLAGS="-fno-automatic $(FFLAGS)" $(MKOCTFILE)
endif
@@ -20,8 +22,8 @@
MKF77FILE = FFLAGS="$(FFLAGS)" $(MKOCTFILE)
endif
-LAPACK_LIBS := $(shell mkoctfile -p LAPACK_LIBS)
-FLIBS := $(shell mkoctfile -p FLIBS)
+LAPACK_LIBS := $(shell $(MKOCTFILE) -p BLAS_LIBS) $(shell $(MKOCTFILE) -p LAPACK_LIBS)
+FLIBS := $(shell $(MKOCTFILE) -p FLIBS)
EXTERNALDIRS = hairer cash daskr
EXTERNALPACKS = $(patsubst %, %.tgz, $(EXTERNALDIRS))
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-27 20:31:37
|
Revision: 10528
http://octave.svn.sourceforge.net/octave/?rev=10528&view=rev
Author: paramaniac
Date: 2012-05-27 20:31:31 +0000 (Sun, 27 May 2012)
Log Message:
-----------
control-devel: handle key/value pairs for are
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m
trunk/octave-forge/extra/control-devel/devel/DestillationMEarx.m
trunk/octave-forge/extra/control-devel/devel/GlassFurnaceARX.m
trunk/octave-forge/extra/control-devel/devel/HeatingSystem.m
trunk/octave-forge/extra/control-devel/devel/LakeErieARX.m
trunk/octave-forge/extra/control-devel/devel/pHarx.m
trunk/octave-forge/extra/control-devel/inst/arx.m
Modified: trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -51,7 +51,7 @@
dat = iddata (Y, U)
% [sys, x0] = ident (dat, 15, 8) % s=15, n=8
-[sys, x0] = arx (dat, 8, 8)
+[sys, x0] = arx (dat, 'na', 8, 'nb', 8)
[y, t] = lsim (sys, U, [], x0);
Modified: trunk/octave-forge/extra/control-devel/devel/DestillationMEarx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/DestillationMEarx.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/DestillationMEarx.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -72,8 +72,8 @@
dat = iddata (Y, U)
[sys, x0] = moen4 (dat, 's', 5, 'n', 4) % s=5, n=4
-sys2 = arx (dat, 4, 4);
-[sys2, x02] = arx (dat, 4, 4);
+sys2 = arx (dat, 'na', 4, 'nb', 4);
+[sys2, x02] = arx (dat, 'na', 4, 'nb', 4);
x0=x0{1};
x02=x02{1};
Modified: trunk/octave-forge/extra/control-devel/devel/GlassFurnaceARX.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/GlassFurnaceARX.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/GlassFurnaceARX.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -49,7 +49,7 @@
dat = iddata (Y, U)
%[sys, x0] = ident (dat, 10, 5) % s=10, n=5
-sys = arx (dat, 5, 5)
+sys = arx (dat, 5)
%[y, t] = lsim (sys, U, [], x0);
[y, t] = lsim (sys(:, 1:3), U);
Modified: trunk/octave-forge/extra/control-devel/devel/HeatingSystem.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/HeatingSystem.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/HeatingSystem.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -77,7 +77,7 @@
% s=15, n=7
[sys1, x0] = moen4 (dat, 's', 15, 'n', 7)
%sys2 = arx (dat, 7, 7) % normally na = nb
-[sys2, x02] = arx (dat, 7, 7);
+[sys2, x02] = arx (dat, 7);
[y1, t1] = lsim (sys1, U, [], x0);
%[y2, t] = lsim (sys2(:, 1), U);
Modified: trunk/octave-forge/extra/control-devel/devel/LakeErieARX.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/LakeErieARX.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/LakeErieARX.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -78,7 +78,7 @@
[sys, x0] = moen4 (dat, 's', 5, 'n', 4) % s=5, n=4
% sys2 = arx (dat, 4, 4)
-[sys2, x02] = arx (dat, 4, 4)
+[sys2, x02] = arx (dat, 4)
x0=x0{1};
x02=x02{1};
Modified: trunk/octave-forge/extra/control-devel/devel/pHarx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/pHarx.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/devel/pHarx.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -52,7 +52,7 @@
dat = iddata (Y, U)
% [sys, x0] = ident (dat, 15, 6) % s=15, n=6
-sys = arx (dat, 6, 6) % normally na = nb
+sys = arx (dat, 6) % normally na = nb
% [y, t] = lsim (sys, U, [], x0);
[y, t] = lsim (sys(:, 1:2), U);
Modified: trunk/octave-forge/extra/control-devel/inst/arx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/arx.m 2012-05-27 20:03:58 UTC (rev 10527)
+++ trunk/octave-forge/extra/control-devel/inst/arx.m 2012-05-27 20:31:31 UTC (rev 10528)
@@ -24,11 +24,11 @@
## Created: April 2012
## Version: 0.1
-function [sys, varargout] = arx (dat, na, nb)
+function [sys, varargout] = arx (dat, varargin)
## TODO: delays
- if (nargin != 3)
+ if (nargin < 2)
print_usage ();
endif
@@ -36,8 +36,24 @@
error ("arx: first argument must be an iddata dataset");
endif
+% if (nargin > 2) # arx (dat, ...)
+ if (is_real_scalar (varargin{1})) # arx (dat, n, ...)
+ varargin = horzcat (varargin(2:end), {"na"}, varargin(1), {"nb"}, varargin(1));
+ endif
+ if (isstruct (varargin{1})) # arx (dat, opt, ...), arx (dat, n, opt, ...)
+ varargin = horzcat (__opt2cell__ (varargin{1}), varargin(2:end));
+ endif
+% endif
+
+ nkv = numel (varargin); # number of keys and values
+
+ if (rem (nkv, 2))
+ error ("arx: keys and values must come in pairs");
+ endif
+
+
## p: outputs, m: inputs, ex: experiments
- [~, p, m, ex] = size (dat);
+ [~, p, m, ex] = size (dat); # dataset dimensions
## extract data
Y = dat.y;
@@ -50,8 +66,31 @@
else
tsam = tsam{1};
endif
+
+
+ ## default arguments
+ na = [];
+ nb = []; % ???
+ nk = [];
-
+ ## handle keys and values
+ for k = 1 : 2 : nkv
+ key = lower (varargin{k});
+ val = varargin{k+1};
+ switch (key)
+ ## TODO: proper argument checking
+ case "na"
+ na = val;
+ case "nb"
+ nb = val;
+ case "nk"
+ error ("nk");
+ otherwise
+ warning ("arx: invalid property name '%s' ignored", key);
+ endswitch
+ endfor
+
+
if (is_real_scalar (na, nb))
na = repmat (na, p, 1); # na(p-by-1)
nb = repmat (nb, p, m); # nb(p-by-m)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-27 20:04:04
|
Revision: 10527
http://octave.svn.sourceforge.net/octave/?rev=10527&view=rev
Author: paramaniac
Date: 2012-05-27 20:03:58 +0000 (Sun, 27 May 2012)
Log Message:
-----------
control-devel: add missing semicolons
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/devel/armax.m
trunk/octave-forge/extra/control-devel/devel/rarx.m
trunk/octave-forge/extra/control-devel/inst/arx.m
Modified: trunk/octave-forge/extra/control-devel/devel/armax.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/armax.m 2012-05-27 19:41:24 UTC (rev 10526)
+++ trunk/octave-forge/extra/control-devel/devel/armax.m 2012-05-27 20:03:58 UTC (rev 10527)
@@ -113,7 +113,7 @@
## this makes only sense for state-space models, therefore convert TF to SS
if (nargout > 1)
sys = prescale (ss (sys(:,1:m)));
- x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0)
+ x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0);
## return x0 as vector for single-experiment data
## instead of a cell containing one vector
if (numel (x0) == 1)
Modified: trunk/octave-forge/extra/control-devel/devel/rarx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/rarx.m 2012-05-27 19:41:24 UTC (rev 10526)
+++ trunk/octave-forge/extra/control-devel/devel/rarx.m 2012-05-27 20:03:58 UTC (rev 10527)
@@ -113,7 +113,7 @@
## this makes only sense for state-space models, therefore convert TF to SS
if (nargout > 1)
sys = prescale (ss (sys(:,1:m)));
- x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0)
+ x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0);
## return x0 as vector for single-experiment data
## instead of a cell containing one vector
if (numel (x0) == 1)
Modified: trunk/octave-forge/extra/control-devel/inst/arx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/arx.m 2012-05-27 19:41:24 UTC (rev 10526)
+++ trunk/octave-forge/extra/control-devel/inst/arx.m 2012-05-27 20:03:58 UTC (rev 10527)
@@ -113,7 +113,7 @@
## this makes only sense for state-space models, therefore convert TF to SS
if (nargout > 1)
sys = prescale (ss (sys(:,1:m)));
- x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0)
+ x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0);
## return x0 as vector for single-experiment data
## instead of a cell containing one vector
if (numel (x0) == 1)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-27 19:41:30
|
Revision: 10526
http://octave.svn.sourceforge.net/octave/?rev=10526&view=rev
Author: paramaniac
Date: 2012-05-27 19:41:24 +0000 (Sun, 27 May 2012)
Log Message:
-----------
control-devel: add armax draft code
Added Paths:
-----------
trunk/octave-forge/extra/control-devel/devel/armax.m
Added: trunk/octave-forge/extra/control-devel/devel/armax.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/armax.m (rev 0)
+++ trunk/octave-forge/extra/control-devel/devel/armax.m 2012-05-27 19:41:24 UTC (rev 10526)
@@ -0,0 +1,215 @@
+## Copyright (C) 2012 Lukas F. Reichlin
+##
+## This file is part of LTI Syncope.
+##
+## LTI Syncope is free software: you can redistribute it and/or modify
+## it under the terms of the GNU General Public License as published by
+## the Free Software Foundation, either version 3 of the License, or
+## (at your option) any later version.
+##
+## LTI Syncope is distributed in the hope that it will be useful,
+## but WITHOUT ANY WARRANTY; without even the implied warranty of
+## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+## GNU General Public License for more details.
+##
+## You should have received a copy of the GNU General Public License
+## along with LTI Syncope. If not, see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{sys} =} arx (@var{dat}, @var{na}, @var{nb})
+## ARX
+## @end deftypefn
+
+## Author: Lukas Reichlin <luk...@gm...>
+## Created: April 2012
+## Version: 0.1
+
+function [sys, varargout] = armax (dat, na, nb, nc)
+
+ ## TODO: delays
+
+ if (nargin != 3)
+ print_usage ();
+ endif
+
+ if (! isa (dat, "iddata"))
+ error ("arx: first argument must be an iddata dataset");
+ endif
+
+ ## p: outputs, m: inputs, ex: experiments
+ [~, p, m, ex] = size (dat);
+
+ ## extract data
+ Y = dat.y;
+ U = dat.u;
+ tsam = dat.tsam;
+
+ ## multi-experiment data requires equal sampling times
+ if (ex > 1 && ! isequal (tsam{:}))
+ error ("arx: require equally sampled experiments");
+ else
+ tsam = tsam{1};
+ endif
+
+
+ if (is_real_scalar (na, nb))
+ na = repmat (na, p, 1); # na(p-by-1)
+ nb = repmat (nb, p, m); # nb(p-by-m)
+ elseif (! (is_real_vector (na) && is_real_matrix (nb) \
+ && rows (na) == p && rows (nb) == p && columns (nb) == m))
+ error ("arx: require na(%dx1) instead of (%dx%d) and nb(%dx%d) instead of (%dx%d)", \
+ p, rows (na), columns (na), p, m, rows (nb), columns (nb));
+ endif
+
+ max_nb = max (nb, [], 2); # one maximum for each row/output, max_nb(p-by-1)
+ n = max (na, max_nb); # n(p-by-1)
+
+ ## create empty cells for numerator and denominator polynomials
+ num = cell (p, m+p);
+ den = cell (p, m+p);
+
+ ## MIMO (p-by-m) models are identified as p MISO (1-by-m) models
+ ## For multi-experiment data, minimize the trace of the error
+ for i = 1 : p # for every output
+ Phi = cell (ex, 1); # one regression matrix per experiment
+ for e = 1 : ex # for every experiment
+ ## avoid warning: toeplitz: column wins anti-diagonal conflict
+ ## therefore set first row element equal to y(1)
+ PhiY = toeplitz (Y{e}(1:end-1, i), [Y{e}(1, i); zeros(na(i)-1, 1)]);
+ ## create MISO Phi for every experiment
+ PhiU = arrayfun (@(x) toeplitz (U{e}(1:end-1, x), [U{e}(1, x); zeros(nb(i,x)-1, 1)]), 1:m, "uniformoutput", false);
+ Phi{e} = (horzcat (-PhiY, PhiU{:}))(n(i):end, :);
+ endfor
+
+ ## compute parameter vector Theta
+ Theta = __theta__ (Phi, Y, i, n);
+
+ ## extract polynomial matrices A and B from Theta
+ ## A is a scalar polynomial for output i, i=1:p
+ ## B is polynomial row vector (1-by-m) for output i
+ A = [1; Theta(1:na(i))]; # a0 = 1, a1 = Theta(1), an = Theta(n)
+ ThetaB = Theta(na(i)+1:end); # all polynomials from B are in one column vector
+ B = mat2cell (ThetaB, nb(i,:)); # now separate the polynomials, one for each input
+ B = reshape (B, 1, []); # make B a row cell (1-by-m)
+ B = cellfun (@(x) [0; x], B, "uniformoutput", false); # b0 = 0 (leading zero required by filt)
+
+ ## add error inputs
+ Be = repmat ({0}, 1, p); # there are as many error inputs as system outputs (p)
+ Be(i) = 1; # inputs m+1:m+p are zero, except m+i which is one
+ num(i, :) = [B, Be]; # numerator polynomials for output i, individual for each input
+ den(i, :) = repmat ({A}, 1, m+p); # in a row (output i), all inputs have the same denominator polynomial
+ endfor
+
+ ## A(q) y(t) = B(q) u(t) + e(t)
+ ## there is only one A per row
+ ## B(z) and A(z) are a Matrix Fraction Description (MFD)
+ ## y = A^-1(q) B(q) u(t) + A^-1(q) e(t)
+ ## since A(q) is a diagonal polynomial matrix, its inverse is trivial:
+ ## the corresponding transfer function has common row denominators.
+
+ sys = filt (num, den, tsam); # filt creates a transfer function in z^-1
+
+ ## compute initial state vector x0 if requested
+ ## this makes only sense for state-space models, therefore convert TF to SS
+ if (nargout > 1)
+ sys = prescale (ss (sys(:,1:m)));
+ x0 = slib01cd (Y, U, sys.a, sys.b, sys.c, sys.d, 0.0)
+ ## return x0 as vector for single-experiment data
+ ## instead of a cell containing one vector
+ if (numel (x0) == 1)
+ x0 = x0{1};
+ endif
+ varargout{1} = x0;
+ endif
+
+endfunction
+
+
+%function theta = __theta__ (phi, y, i, n)
+function Theta = __theta__ (Phi, Y, i, n)
+
+
+ if (numel (Phi) == 1) # single-experiment dataset
+ % recursive pseudolinear regression, naive formula
+ gamma = ?
+ Theta = ?
+ R = ?
+ for t = 1 : rows (Phi{1})
+ phi = Phi{1}(t,:); # note that my phi is Ljung's phi.'
+ y = Y{1}(t+n(i), :);
+ epsilon = y - phi*Theta;
+ R += gamma * (phi.'*phi - R)
+ Theta += gamma * R \ phi.' * epsilon;
+ endfor
+%{
+ % recursive least-squares with efficient matrix inversion
+ [pr, pc] = size (Phi{1});
+ lambda = 1; % default 1
+ Theta = zeros (pc, 1);
+ P = 10 * eye (pc);
+ for t = 1 : pr
+ phi = Phi{1}(t,:); # note that my phi is Ljung's phi.'
+ y = Y{1}(t+n(i), :);
+ den = lambda + phi*P*phi.';
+ L = P * phi.' / den;
+ P = (P - (P * phi.' * phi * P) / den) / lambda;
+ Theta += L * (y - phi*Theta);
+ endfor
+%}
+%{
+ ## use "square-root algorithm"
+ A = horzcat (phi{1}, y{1}(n(i)+1:end, i)); # [Phi, Y]
+ R0 = triu (qr (A, 0)); # 0 for economy-size R (without zero rows)
+ R1 = R0(1:end-1, 1:end-1); # R1 is triangular - can we exploit this in R1\R2?
+ R2 = R0(1:end-1, end);
+ theta = __ls_svd__ (R1, R2); # R1 \ R2
+
+ ## Theta = Phi \ Y(n+1:end, :); # naive formula
+ ## theta = __ls_svd__ (phi{1}, y{1}(n(i)+1:end, i));
+%}
+ else # multi-experiment dataset
+ ## TODO: find more sophisticated formula than
+ ## Theta = (Phi1' Phi + Phi2' Phi2 + ...) \ (Phi1' Y1 + Phi2' Y2 + ...)
+
+ ## covariance matrix C = (Phi1' Phi + Phi2' Phi2 + ...)
+ tmp = cellfun (@(Phi) Phi.' * Phi, phi, "uniformoutput", false);
+ rc = cellfun (@rcond, tmp); # C auch noch testen? QR oder SVD?
+ C = plus (tmp{:});
+
+ ## PhiTY = (Phi1' Y1 + Phi2' Y2 + ...)
+ tmp = cellfun (@(Phi, Y) Phi.' * Y(n(i)+1:end, i), phi, y, "uniformoutput", false);
+ PhiTY = plus (tmp{:});
+
+ ## pseudoinverse Theta = C \ Phi'Y
+ theta = __ls_svd__ (C, PhiTY);
+ endif
+
+endfunction
+
+
+function x = __ls_svd__ (A, b)
+
+ ## solve the problem Ax=b
+ ## x = A\b would also work,
+ ## but this way we have better control and warnings
+
+ ## solve linear least squares problem by pseudoinverse
+ ## the pseudoinverse is computed by singular value decomposition
+ ## M = U S V* ---> M+ = V S+ U*
+ ## Th = Ph \ Y = Ph+ Y
+ ## Th = V S+ U* Y, S+ = 1 ./ diag (S)
+
+ [U, S, V] = svd (A, 0); # 0 for "economy size" decomposition
+ S = diag (S); # extract main diagonal
+ r = sum (S > eps*S(1));
+ if (r < length (S))
+ warning ("arx: rank-deficient coefficient matrix");
+ warning ("sampling time too small");
+ warning ("persistence of excitation");
+ endif
+ V = V(:, 1:r);
+ S = S(1:r);
+ U = U(:, 1:r);
+ x = V * (S .\ (U' * b)); # U' is the conjugate transpose
+
+endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <cd...@us...> - 2012-05-27 18:53:51
|
Revision: 10525
http://octave.svn.sourceforge.net/octave/?rev=10525&view=rev
Author: cdf
Date: 2012-05-27 18:53:45 +0000 (Sun, 27 May 2012)
Log Message:
-----------
vectorized nrbcrvderiveval
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbcrvderiveval.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbcrvderiveval.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbcrvderiveval.m 2012-05-27 18:22:21 UTC (rev 10524)
+++ trunk/octave-forge/extra/nurbs/inst/nrbcrvderiveval.m 2012-05-27 18:53:45 UTC (rev 10525)
@@ -1,6 +1,3 @@
-function ck = nrbcrvderiveval (crv, u, d)
-
-%
% NRBCRVDERIVEVAL: Evaluate n-th order derivatives of a NURBS curve.
%
% usage: skl = nrbcrvderiveval (crv, u, d)
@@ -34,35 +31,48 @@
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+function ck = nrbcrvderiveval (crv, u, d)
+ ck = arrayfun (@(x) nrbcrvderiveval__ (crv, x, d), u, 'UniformOutput', false);
+ ck = cat (3, ck{:});
+end
+
+function ck = nrbcrvderiveval__ (crv, u, d)
+
+ persistent nc;
+ if isempty (nc)
+ nc = [0 0 0 0 0;
+ 1 0 0 0 0;
+ 2 1 0 0 0;
+ 3 3 1 0 0;
+ 4 6 4 1 0];
+ end
+
+ ck = zeros (3, d+1);
+ wders = curvederiveval (crv.number-1, crv.order-1, crv.knots, squeeze (crv.coefs(4, :)), u, d);
+
+ for idim = 1:3
- ck = zeros (3, d+1, numel(u));
-
- for iu = 1:numel(u);
- wders = squeeze (curvederiveval (crv.number-1, crv.order-1, ...
- crv.knots, squeeze (crv.coefs(4, :)), u(iu), d));
+ Aders = curvederiveval (crv.number-1, crv.order-1, crv.knots, squeeze (crv.coefs(idim, :)), u, d);
+
+ ck(idim, 1) = Aders(1) / wders(1);
+ for k = 1:d
+ ck(idim, k+1) = (Aders(k+1) - sum (nc(k+1, 1:k) .* wders(2:k+1).' .* squeeze (ck(idim, k:-1:1)))) / wders(1);
+ end
+
+ end
- for idim = 1:3
- Aders = squeeze (curvederiveval (crv.number-1, crv.order-1, ...
- crv.knots, squeeze (crv.coefs(idim, :)), u(iu), d));
- for k=0:d
- v = Aders(k+1);
- for i=1:k
- v = v - nchoosek(k,i)*wders(i+1)*ck(idim, k-i+1, iu);
- end
- ck(idim, k+1, iu) = v/wders(1);
- end
- end
- end
end
+
%!test
%! knots = [0 0 0 1 1 1];
%! coefs(:,1) = [0; 0; 0; 1];
%! coefs(:,2) = [1; 0; 1; 1];
%! coefs(:,3) = [1; 1; 1; 2];
%! crv = nrbmak (coefs, knots);
-%! u = linspace (0, 1, 10);
-%! ck = nrbcrvderiveval (crv, u, 2);
+%! u = linspace (0, 1, 100);
+%! ck = nrbcrvderiveval (crv, u, 2);
%! w = @(x) 1 + x.^2;
%! dw = @(x) 2*x;
%! F1 = @(x) (2*x - x.^2)./w(x);
@@ -76,4 +86,4 @@
%! d2F3 = @(x) -2./w(x) - 2*x.*(2-2*x)./w(x).^2 - (8*x-6*x.^2)./w(x).^2 + 8*x.^2.*(2*x-x.^2)./w(x).^3;
%! assert ([F1(u); F2(u); F3(u)], squeeze(ck(:, 1, :)), 1e2*eps);
%! assert ([dF1(u); dF2(u); dF3(u)], squeeze(ck(:, 2, :)), 1e2*eps);
-%! assert ([d2F1(u); d2F2(u); d2F3(u)], squeeze(ck(:, 3, :)), 1e2*eps);
+%! assert ([d2F1(u); d2F2(u); d2F3(u)], squeeze(ck(:, 3, :)), 1e2*eps);
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-05-27 18:22:28
|
Revision: 10524
http://octave.svn.sourceforge.net/octave/?rev=10524&view=rev
Author: paramaniac
Date: 2012-05-27 18:22:21 +0000 (Sun, 27 May 2012)
Log Message:
-----------
control-devel: recursive least squares seems to work
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m
trunk/octave-forge/extra/control-devel/devel/rarx.m
trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
Added Paths:
-----------
trunk/octave-forge/extra/control-devel/devel/HeatingSystemRLS.m
Modified: trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m 2012-05-26 01:07:07 UTC (rev 10523)
+++ trunk/octave-forge/extra/control-devel/devel/CDplayerARX.m 2012-05-27 18:22:21 UTC (rev 10524)
@@ -51,8 +51,11 @@
dat = iddata (Y, U)
% [sys, x0] = ident (dat, 15, 8) % s=15, n=8
-sys = arx (dat, 4, 4)
+[sys, x0] = arx (dat, 8, 8)
+[y, t] = lsim (sys, U, [], x0);
+
+%{
%[y, t] = lsim (sys, U, [], x0);
%[y, t] = lsim (sys(:,1:2), U);
@@ -65,6 +68,7 @@
y = [y1, y2];
t = 0:length(U)-1;
+%}
err = norm (Y - y, 1) / norm (Y, 1)
Added: trunk/octave-forge/extra/control-devel/devel/HeatingSystemRLS.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/HeatingSystemRLS.m (rev 0)
+++ trunk/octave-forge/extra/control-devel/devel/HeatingSystemRLS.m 2012-05-27 18:22:21 UTC (rev 10524)
@@ -0,0 +1,102 @@
+%{
+1. Contributed by:
+
+ Roy Smith
+ Dept. of Electrical & Computer Engineering
+ University of California,
+ Santa Barbara, CA 93106
+ U.S.A.
+ ro...@ec...
+
+2. Process/Description:
+
+ The experiment is a simple SISO heating system.
+ The input drives a 300 Watt Halogen lamp, suspended
+ several inches above a thin steel plate. The output
+ is a thermocouple measurement taken from the back of
+ the plate.
+
+3. Sampling interval:
+
+ 2.0 seconds
+
+4. Number of samples
+
+ 801
+
+5. Inputs:
+
+ u: input drive voltage
+ ...
+6. Outputs:
+
+ y: temperature (deg. C)
+ ...
+7. References:
+
+ The use of this experiment and data for robust
+ control model validation is described in:
+
+ "Sampled Data Model Validation: an Algorithm and
+ Experimental Application," Geir Dullerud & Roy Smith,
+ International Journal of Robust and Nonlinear Control,
+ Vol. 6, No. 9/10, pp. 1065-1078, 1996.
+
+8. Known properties/peculiarities
+
+ The data (and nominal model) is the above paper have the
+ output expressed in 10's deg. C. This has been rescaled
+ to the original units of deg. C. in the DaISy data set.
+ There is also a -1 volt offset in u in the data shown plotted
+ in the original paper. This has been removed in the
+ DaISy dataset.
+
+ The data shows evidence of discrepancies. One of the
+ issues studied in the above paper is the size of these
+ discrepancies - measured in this case in terms of the norm
+ of the smallest perturbation required to account for the
+ difference between the nominal model and the data.
+
+ The steady state input (prior to the start of the experiment)
+ is u = 6.0 Volts.
+
+%}
+
+clear all, close all, clc
+
+load heating_system.dat
+U=heating_system(:,2);
+Y=heating_system(:,3);
+
+
+dat = iddata (Y, U, 2.0, 'inname', 'input drive voltage', \
+ 'inunit', 'Volt', \
+ 'outname', 'temperature', \
+ 'outunit', '°C')
+
+% s=15, n=7
+[sys1, x0] = moen4 (dat, 's', 15, 'n', 7)
+%sys2 = arx (dat, 7, 7) % normally na = nb
+[sys2, x02] = arx (dat, 7, 7);
+[sys3, x03] = rarx (dat, 7, 7);
+
+[y1, t1] = lsim (sys1, U, [], x0);
+%[y2, t] = lsim (sys2(:, 1), U);
+[y2, t] = lsim (sys2, U, [], x02);
+[y3, t] = lsim (sys3, U, [], x02);
+
+
+err1 = norm (Y - y1, 1) / norm (Y, 1)
+err2 = norm (Y - y2, 1) / norm (Y, 1)
+err2 = norm (Y - y3, 1) / norm (Y, 1)
+
+figure (1)
+plot (t, Y, t, y1, t, y2, t, y3)
+title ('DaISy: Heating System [99-001]')
+xlim ([t(1), t(end)])
+xlabel ('Time [s]')
+ylabel ('Temperature [Degree Celsius]')
+legend ('measured', 'simulated subspace', 'simulated ARX', 'simulated RLS', 'location', 'southeast')
+
+
+
Modified: trunk/octave-forge/extra/control-devel/devel/rarx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/rarx.m 2012-05-26 01:07:07 UTC (rev 10523)
+++ trunk/octave-forge/extra/control-devel/devel/rarx.m 2012-05-27 18:22:21 UTC (rev 10524)
@@ -131,16 +131,16 @@
if (numel (Phi) == 1) # single-experiment dataset
% recursive least-squares with efficient matrix inversion
- lambda = default 1
- Theta = ?
- P = ?
- for t =
- phi = Phi{1}(t,:);
- y = Y{1}(t?offset n(i), :);
- ## note phi != phi.'
+ [pr, pc] = size (Phi{1});
+ lambda = 1; % default 1
+ Theta = zeros (pc, 1);
+ P = 10 * eye (pc);
+ for t = 1 : pr
+ phi = Phi{1}(t,:); # note that my phi is Ljung's phi.'
+ y = Y{1}(t+n(i), :);
den = lambda + phi*P*phi.';
- L = P * phi / den;
- P = (P - (P * phi.' * phi * P) / den) / lambda
+ L = P * phi.' / den;
+ P = (P - (P * phi.' * phi * P) / den) / lambda;
Theta += L * (y - phi*Theta);
endfor
%{
Modified: trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-26 01:07:07 UTC (rev 10523)
+++ trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-05-27 18:22:21 UTC (rev 10524)
@@ -58,7 +58,7 @@
endif
[ns, l, m, e] = size (dat); # dataset dimensions
- tsam = dat.tsam
+ tsam = dat.tsam;
## multi-experiment data requires equal sampling times
if (e > 1 && ! isequal (tsam{:}))
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
48 messages has been excluded from this view by a project administrator.
