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: <raf...@us...> - 2012-07-17 13:00:21
|
Revision: 10748
http://octave.svn.sourceforge.net/octave/?rev=10748&view=rev
Author: rafavzqz
Date: 2012-07-17 13:00:15 +0000 (Tue, 17 Jul 2012)
Log Message:
-----------
Fixed bug in nrbplot. New version of nrbkntplot
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m
trunk/octave-forge/extra/nurbs/inst/nrbplot.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m 2012-07-17 11:30:00 UTC (rev 10747)
+++ trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m 2012-07-17 13:00:15 UTC (rev 10748)
@@ -20,7 +20,7 @@
%
% nrbctrlplot
%
-% Copyright (C) 2011 Rafael Vazquez
+% Copyright (C) 2011, 2012 Rafael Vazquez
%
% 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
@@ -80,46 +80,15 @@
elseif (size (nurbs.knots,2) == 3) % plot a NURBS volume
- nsub = 30;
- nrbplot (nurbs, [nsub nsub nsub], 'light', light, 'colormap', cmap);
- hold on
+ nsub = 100;
- % And plot the knots
- knt1 = unique (nurbs.knots{1});
- knt2 = unique (nurbs.knots{2});
- knt3 = unique (nurbs.knots{3});
- kv_face1 = nrbeval (nurbs, {knt1(1), knt2, linspace(knt3(1),knt3(end),nsub)});
- kw_face1 = nrbeval (nurbs, {knt1(1), linspace(knt2(1),knt2(end),nsub), knt3});
- kv_face2 = nrbeval (nurbs, {knt1(end), knt2, linspace(knt3(1),knt3(end),nsub)});
- kw_face2 = nrbeval (nurbs, {knt1(end), linspace(knt2(1),knt2(end),nsub), knt3});
- ku_face3 = nrbeval (nurbs, {knt1, knt2(1), linspace(knt3(1),knt3(end),nsub)});
- kw_face3 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2(1), knt3});
- ku_face4 = nrbeval (nurbs, {knt1, knt2(end), linspace(knt3(1),knt3(end),nsub)});
- kw_face4 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2(end), knt3});
- ku_face5 = nrbeval (nurbs, {knt1, linspace(knt2(1),knt2(end),nsub), knt3(1)});
- kv_face5 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2, knt3(1)});
- ku_face6 = nrbeval (nurbs, {knt1, linspace(knt2(1),knt2(end),nsub), knt3(end)});
- kv_face6 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2, knt3(end)});
-
- for ii = 1:numel(knt1)
- plot3 (squeeze (ku_face3(1,ii,:,:)), squeeze (ku_face3(2,ii,:,:)), squeeze (ku_face3(3,ii,:,:)));
- plot3 (squeeze (ku_face4(1,ii,:,:)), squeeze (ku_face4(2,ii,:,:)), squeeze (ku_face4(3,ii,:,:)));
- plot3 (squeeze (ku_face5(1,ii,:,:)), squeeze (ku_face5(2,ii,:,:)), squeeze (ku_face5(3,ii,:,:)));
- plot3 (squeeze (ku_face6(1,ii,:,:)), squeeze (ku_face6(2,ii,:,:)), squeeze (ku_face6(3,ii,:,:)));
+ % Plot the boundaries
+ bnd = nrbextract (nurbs);
+ nrbkntplot (bnd(1));
+ hold on
+ for iface = 2:6
+ nrbkntplot (bnd(iface));
end
- for ii = 1:numel(knt2)
- plot3 (squeeze (kv_face1(1,:,ii,:)), squeeze (kv_face1(2,:,ii,:)), squeeze (kv_face1(3,:,ii,:)));
- plot3 (squeeze (kv_face2(1,:,ii,:)), squeeze (kv_face2(2,:,ii,:)), squeeze (kv_face2(3,:,ii,:)));
- plot3 (squeeze (kv_face5(1,:,ii,:)), squeeze (kv_face5(2,:,ii,:)), squeeze (kv_face5(3,:,ii,:)));
- plot3 (squeeze (kv_face6(1,:,ii,:)), squeeze (kv_face6(2,:,ii,:)), squeeze (kv_face6(3,:,ii,:)));
- end
- for ii = 1:numel(knt3)
- plot3 (squeeze (kw_face1(1,:,:,ii)), squeeze(kw_face1(2,:,:,ii)), squeeze (kw_face1(3,:,:,ii)));
- plot3 (squeeze (kw_face2(1,:,:,ii)), squeeze(kw_face2(2,:,:,ii)), squeeze (kw_face2(3,:,:,ii)));
- plot3 (squeeze (kw_face3(1,:,:,ii)), squeeze(kw_face3(2,:,:,ii)), squeeze (kw_face3(3,:,:,ii)));
- plot3 (squeeze (kw_face4(1,:,:,ii)), squeeze(kw_face4(2,:,:,ii)), squeeze (kw_face4(3,:,:,ii)));
- end
- end
else % plot a NURBS curve
nsub = 1000;
Modified: trunk/octave-forge/extra/nurbs/inst/nrbplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 11:30:00 UTC (rev 10747)
+++ trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 13:00:15 UTC (rev 10748)
@@ -32,6 +32,7 @@
%
% Copyright (C) 2000 Mark Spink
% Copyright (C) 2010 Carlo de Falco, Rafael Vazquez
+% Copyright (C) 2012 Rafael Vazquez
%
% 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
@@ -109,13 +110,13 @@
elseif (size (nurbs.knots,2) == 3) % plot the boundaries of a NURBS volume
bnd = nrbextract (nurbs);
hold_flag = ishold;
- nrbplot (bnd(1), subd(2:3));
+ nrbplot (bnd(1), subd(2:3), varargin{:});
hold on
- nrbplot (bnd(2), subd(2:3));
- nrbplot (bnd(3), subd([1 3]));
- nrbplot (bnd(4), subd([1 3]));
- nrbplot (bnd(5), subd(1:2));
- nrbplot (bnd(6), subd(1:2));
+ nrbplot (bnd(2), subd(2:3), varargin{:});
+ nrbplot (bnd(3), subd([1 3]), varargin{:});
+ nrbplot (bnd(4), subd([1 3]), varargin{:});
+ nrbplot (bnd(5), subd(1:2), varargin{:});
+ nrbplot (bnd(6), subd(1:2), varargin{:});
if (~hold_flag)
hold off
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <raf...@us...> - 2012-07-17 11:30:11
|
Revision: 10747
http://octave.svn.sourceforge.net/octave/?rev=10747&view=rev
Author: rafavzqz
Date: 2012-07-17 11:30:00 +0000 (Tue, 17 Jul 2012)
Log Message:
-----------
Plot of trivariates is much easier
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbplot.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 11:08:00 UTC (rev 10746)
+++ trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 11:30:00 UTC (rev 10747)
@@ -107,30 +107,15 @@
shading faceted;
end
elseif (size (nurbs.knots,2) == 3) % plot the boundaries of a NURBS volume
- knt = nurbs.knots;
+ bnd = nrbextract (nurbs);
hold_flag = ishold;
- px = nrbeval (nurbs, {knt{1}([1 end]) linspace(knt{2}(1),knt{2}(end),subd(2)) linspace(knt{3}(1),knt{3}(end),subd(3))});
- py = nrbeval (nurbs, {linspace(knt{1}(1),knt{1}(end),subd(1)) knt{2}([1 end]) linspace(knt{3}(1),knt{3}(end),subd(3))});
- pz = nrbeval (nurbs, {linspace(knt{1}(1),knt{1}(end),subd(1)) linspace(knt{2}(1),knt{2}(end),subd(2)) knt{3}([1 end])});
- if (strcmp (light, 'on'))
- surfl (squeeze (pz(1,:,:,1)), squeeze (pz(2,:,:,1)), squeeze (pz(3,:,:,1)));
- hold on
- surfl (squeeze (pz(1,:,:,2)), squeeze (pz(2,:,:,2)), squeeze (pz(3,:,:,2)));
- surfl (squeeze (py(1,:,1,:)), squeeze (py(2,:,1,:)), squeeze (py(3,:,1,:)));
- surfl (squeeze (py(1,:,2,:)), squeeze (py(2,:,2,:)), squeeze (py(3,:,2,:)));
- surfl (squeeze (px(1,1,:,:)), squeeze (px(2,1,:,:)), squeeze (px(3,1,:,:)));
- surfl (squeeze (px(1,2,:,:)), squeeze (px(2,2,:,:)), squeeze (px(3,2,:,:)));
- shading interp;
- else
- surf (squeeze (pz(1,:,:,1)), squeeze (pz(2,:,:,1)), squeeze (pz(3,:,:,1)));
- hold on
- surf (squeeze (pz(1,:,:,2)), squeeze (pz(2,:,:,2)), squeeze (pz(3,:,:,2)));
- surf (squeeze (py(1,:,1,:)), squeeze (py(2,:,1,:)), squeeze (py(3,:,1,:)));
- surf (squeeze (py(1,:,2,:)), squeeze (py(2,:,2,:)), squeeze (py(3,:,2,:)));
- surf (squeeze (px(1,1,:,:)), squeeze (px(2,1,:,:)), squeeze (px(3,1,:,:)));
- surf (squeeze (px(1,2,:,:)), squeeze (px(2,2,:,:)), squeeze (px(3,2,:,:)));
- shading faceted;
- end
+ nrbplot (bnd(1), subd(2:3));
+ hold on
+ nrbplot (bnd(2), subd(2:3));
+ nrbplot (bnd(3), subd([1 3]));
+ nrbplot (bnd(4), subd([1 3]));
+ nrbplot (bnd(5), subd(1:2));
+ nrbplot (bnd(6), subd(1:2));
if (~hold_flag)
hold off
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <raf...@us...> - 2012-07-17 11:08:07
|
Revision: 10746
http://octave.svn.sourceforge.net/octave/?rev=10746&view=rev
Author: rafavzqz
Date: 2012-07-17 11:08:00 +0000 (Tue, 17 Jul 2012)
Log Message:
-----------
Plot for non-unit interval
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m
trunk/octave-forge/extra/nurbs/inst/nrbplot.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m 2012-07-17 08:53:26 UTC (rev 10745)
+++ trunk/octave-forge/extra/nurbs/inst/nrbkntplot.m 2012-07-17 11:08:00 UTC (rev 10746)
@@ -56,8 +56,8 @@
% And plot the knots
knt1 = unique (nurbs.knots{1});
knt2 = unique (nurbs.knots{2});
- p1 = nrbeval (nurbs, {knt1, linspace(0.0,1.0,nsub)});
- p2 = nrbeval (nurbs, {linspace(0.0,1.0,nsub), knt2});
+ p1 = nrbeval (nurbs, {knt1, linspace(knt2(1),knt2(end),nsub)});
+ p2 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2});
if (any (nurbs.coefs(3,:)))
% surface in a 3D space
@@ -88,18 +88,18 @@
knt1 = unique (nurbs.knots{1});
knt2 = unique (nurbs.knots{2});
knt3 = unique (nurbs.knots{3});
- kv_face1 = nrbeval (nurbs, {0, knt2, linspace(0.0,1.0,nsub)});
- kw_face1 = nrbeval (nurbs, {0, linspace(0.0,1.0,nsub), knt3});
- kv_face2 = nrbeval (nurbs, {1, knt2, linspace(0.0,1.0,nsub)});
- kw_face2 = nrbeval (nurbs, {1, linspace(0.0,1.0,nsub), knt3});
- ku_face3 = nrbeval (nurbs, {knt1, 0, linspace(0.0,1.0,nsub)});
- kw_face3 = nrbeval (nurbs, {linspace(0.0,1.0,nsub), 0, knt3});
- ku_face4 = nrbeval (nurbs, {knt1, 1, linspace(0.0,1.0,nsub)});
- kw_face4 = nrbeval (nurbs, {linspace(0.0,1.0,nsub), 1, knt3});
- ku_face5 = nrbeval (nurbs, {knt1, linspace(0.0,1.0,nsub), 0});
- kv_face5 = nrbeval (nurbs, {linspace(0.0,1.0,nsub), knt2, 0});
- ku_face6 = nrbeval (nurbs, {knt1, linspace(0.0,1.0,nsub), 1});
- kv_face6 = nrbeval (nurbs, {linspace(0.0,1.0,nsub), knt2, 1});
+ kv_face1 = nrbeval (nurbs, {knt1(1), knt2, linspace(knt3(1),knt3(end),nsub)});
+ kw_face1 = nrbeval (nurbs, {knt1(1), linspace(knt2(1),knt2(end),nsub), knt3});
+ kv_face2 = nrbeval (nurbs, {knt1(end), knt2, linspace(knt3(1),knt3(end),nsub)});
+ kw_face2 = nrbeval (nurbs, {knt1(end), linspace(knt2(1),knt2(end),nsub), knt3});
+ ku_face3 = nrbeval (nurbs, {knt1, knt2(1), linspace(knt3(1),knt3(end),nsub)});
+ kw_face3 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2(1), knt3});
+ ku_face4 = nrbeval (nurbs, {knt1, knt2(end), linspace(knt3(1),knt3(end),nsub)});
+ kw_face4 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2(end), knt3});
+ ku_face5 = nrbeval (nurbs, {knt1, linspace(knt2(1),knt2(end),nsub), knt3(1)});
+ kv_face5 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2, knt3(1)});
+ ku_face6 = nrbeval (nurbs, {knt1, linspace(knt2(1),knt2(end),nsub), knt3(end)});
+ kv_face6 = nrbeval (nurbs, {linspace(knt1(1),knt1(end),nsub), knt2, knt3(end)});
for ii = 1:numel(knt1)
plot3 (squeeze (ku_face3(1,ii,:,:)), squeeze (ku_face3(2,ii,:,:)), squeeze (ku_face3(3,ii,:,:)));
Modified: trunk/octave-forge/extra/nurbs/inst/nrbplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 08:53:26 UTC (rev 10745)
+++ trunk/octave-forge/extra/nurbs/inst/nrbplot.m 2012-07-17 11:08:00 UTC (rev 10746)
@@ -95,7 +95,9 @@
% plot the curve or surface
if (iscell (nurbs.knots))
if (size (nurbs.knots,2) == 2) % plot a NURBS surface
- p = nrbeval (nurbs, {linspace(0.0,1.0,subd(1)) linspace(0.0,1.0,subd(2))});
+ knt = nurbs.knots;
+ p = nrbeval (nurbs, {linspace(knt{1}(1),knt{1}(end),subd(1)) ...
+ linspace(knt{2}(1),knt{2}(end),subd(2))});
if (strcmp (light,'on'))
% light surface
surfl (squeeze(p(1,:,:)), squeeze(p(2,:,:)), squeeze(p(3,:,:)));
@@ -105,10 +107,11 @@
shading faceted;
end
elseif (size (nurbs.knots,2) == 3) % plot the boundaries of a NURBS volume
+ knt = nurbs.knots;
hold_flag = ishold;
- px = nrbeval (nurbs, {[0 1] linspace(0.0,1.0,subd(2)) linspace(0.0,1.0,subd(3))});
- py = nrbeval (nurbs, {linspace(0.0,1.0,subd(1)) [0 1] linspace(0.0,1.0,subd(3))});
- pz = nrbeval (nurbs, {linspace(0.0,1.0,subd(1)) linspace(0.0,1.0,subd(2)) [0 1]});
+ px = nrbeval (nurbs, {knt{1}([1 end]) linspace(knt{2}(1),knt{2}(end),subd(2)) linspace(knt{3}(1),knt{3}(end),subd(3))});
+ py = nrbeval (nurbs, {linspace(knt{1}(1),knt{1}(end),subd(1)) knt{2}([1 end]) linspace(knt{3}(1),knt{3}(end),subd(3))});
+ pz = nrbeval (nurbs, {linspace(knt{1}(1),knt{1}(end),subd(1)) linspace(knt{2}(1),knt{2}(end),subd(2)) knt{3}([1 end])});
if (strcmp (light, 'on'))
surfl (squeeze (pz(1,:,:,1)), squeeze (pz(2,:,:,1)), squeeze (pz(3,:,:,1)));
hold on
@@ -138,7 +141,7 @@
end
else
% plot a NURBS curve
- p = nrbeval (nurbs, linspace (0.0, 1.0, subd));
+ p = nrbeval (nurbs, linspace (nurbs.knots(1), nurbs.knots(end), subd));
if (any (nurbs.coefs(3,:)))
% 3D curve
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <raf...@us...> - 2012-07-17 08:53:37
|
Revision: 10745
http://octave.svn.sourceforge.net/octave/?rev=10745&view=rev
Author: rafavzqz
Date: 2012-07-17 08:53:26 +0000 (Tue, 17 Jul 2012)
Log Message:
-----------
Plot the points only once
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m 2012-07-17 00:18:06 UTC (rev 10744)
+++ trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m 2012-07-17 08:53:26 UTC (rev 10745)
@@ -8,7 +8,7 @@
%
% INPUT:
%
-% nurbs: NURBS curve or surface, see nrbmak.
+% nurbs: NURBS curve, surface or volume, see nrbmak.
%
% Example:
%
@@ -22,7 +22,7 @@
%
% nrbkntplot
%
-% Copyright (C) 2011 Rafael Vazquez
+% Copyright (C) 2011, 2012 Rafael Vazquez
%
% 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
@@ -54,19 +54,21 @@
nsub = 100;
nrbplot (nurbs, [nsub nsub nsub], 'light', light, 'colormap', cmap);
hold on
-% Plot the the control net
+% Plot the control points
+ coefs = bsxfun (@rdivide, nurbs.coefs(1:3,:,:,:), nurbs.coefs(4,:,:,:));
+ coefs = reshape (coefs, 3, []);
+ plot3 (coefs(1,:), coefs(2,:), coefs(3,:), 'r.','MarkerSize',20);
+% Plot the control net
for ii = 1:size (nurbs.coefs, 2)
for jj = 1:size (nurbs.coefs, 3)
coefs = reshape (nurbs.coefs(1:3,ii,jj,:), 3, []);
weights = reshape (nurbs.coefs(4,ii,jj,:), 1, []);
plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
- plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
end
for kk = 1:size (nurbs.coefs, 4)
coefs = reshape (nurbs.coefs(1:3,ii,:,kk), 3, []);
weights = reshape (nurbs.coefs(4,ii,:,kk), 1, []);
plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
- plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
end
end
for jj = 1:size (nurbs.coefs, 3)
@@ -74,7 +76,6 @@
coefs = reshape (nurbs.coefs(1:3,:,jj,kk), 3, []);
weights = reshape (nurbs.coefs(4,:,jj,kk), 1, []);
plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
- plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
end
end
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <be...@us...> - 2012-07-17 00:18:12
|
Revision: 10744
http://octave.svn.sourceforge.net/octave/?rev=10744&view=rev
Author: benjf5
Date: 2012-07-17 00:18:06 +0000 (Tue, 17 Jul 2012)
Log Message:
-----------
fastlscomplex has been mostly rewritten, and has problems with rounding errors (but generates some correct results!).
Modified Paths:
--------------
trunk/octave-forge/extra/lssa/SampleScriptWithVostokData
trunk/octave-forge/extra/lssa/fastlscomplex.cc
Modified: trunk/octave-forge/extra/lssa/SampleScriptWithVostokData
===================================================================
--- trunk/octave-forge/extra/lssa/SampleScriptWithVostokData 2012-07-15 20:13:22 UTC (rev 10743)
+++ trunk/octave-forge/extra/lssa/SampleScriptWithVostokData 2012-07-17 00:18:06 UTC (rev 10744)
@@ -35,17 +35,190 @@
ls_complex_co2_gas_age = lscomplex(co2(:,3),co2(:,4),1,100,20);
ls_real_co2_ice_age = lsreal(co2(:,2),co2(:,4),1,100,20);
ls_real_co2_gas_age = lsreal(co2(:,3),co2(:,4),1,100,20);
+ls_complex_ch4_ice_age = lscomplex(ch4(:,2),ch4(:,4),1,100,20);
+ls_complex_ch4_gas_age = lscomplex(ch4(:,3),ch4(:,4),1,100,20);
+ls_real_ch4_ice_age = lsreal(ch4(:,2),ch4(:,4),1,100,20);
+ls_real_ch4_gas_age = lsreal(ch4(:,3),ch4(:,4),1,100,20);
+ls_complex_o18_ice_age = lscomplex(o18(:,2),o18(:,4),1,100,20);
+ls_complex_o18_gas_age = lscomplex(o18(:,3),o18(:,4),1,100,20);
+ls_real_o18_ice_age = lsreal(o18(:,2),o18(:,4),1,100,20);
+ls_real_o18_gas_age = lsreal(o18(:,3),o18(:,4),1,100,20);
+ls_complex_deut = lscomplex(deut(:,2),deut(:,3),1,100,20);
+ls_real_deut = lsreal(deut(:,2),deut(:,3),1,100,20);
+ls_complex_dust = lscomplex(dust(:,2),dust(:,3),1,100,20);
+ls_real_dust = lsreal(dust(:,2),dust(:,3),1,100,20);
+x_data_axis_vector = [ -430000, 0 ];
+## Useful because all of the data extends over 430 000 years up to the
+## present.
+
+## Setting up the CO2 plots:
figure(co2_fig);
subplot(4,2,1);
+axis(x_data_axis_vector);
plot(-(co2(:,2)),co2(:,4));
title("Gas levels over ice age");
subplot(4,2,2);
+axis(x_data_axis_vector);
plot(-(co2(:,3),co2(:,4));
title("Gas levels over gas age");
subplot(4,2,3);
+plot(real(ls_complex_co2_ice_age));
+hold on;
+plot(imag(ls_complex_co2_ice_age),'r');
+title("Complex L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,4);
+plot(real(ls_complex_co2_gas_age));
+hold on;
+plot(imag(ls_complex_co2_gas_age),'r');
+title("Complex L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+subplot(4,2,5);
+plot(real(ls_real_co2_ice_age));
+hold on;
+plot(imag(ls_real_co2_ice_age),'r');
+title("Real L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,6);
+plot(real(ls_real_co2_gas_age));
+hold on;
+plot(imag(ls_real_co2_gas_age));
+title("Real L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+## At this point, we have transforms of both datasets, real and complex,
+## and just need to figure out what cool thing to do with the remaining slot.
+## Setting up the CH4 plots
+figure(ch4_fig);
+subplot(4,2,1);
+axis(x_data_axis_vector);
+plot(-(ch4(:,2)),ch4(:,4));
+title("Gas levels over ice age");
+subplot(4,2,2);
+axis(x_data_axis_vector);
+plot(-(ch4(:,3),ch4(:,4));
+title("Gas levels over gas age");
+subplot(4,2,3);
+plot(real(ls_complex_ch4_ice_age));
+hold on;
+plot(imag(ls_complex_ch4_ice_age),'r');
+title("Complex L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,4);
+plot(real(ls_complex_ch4_gas_age));
+hold on;
+plot(imag(ls_complex_ch4_gas_age),'r');
+title("Complex L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+subplot(4,2,5);
+plot(real(ls_real_ch4_ice_age));
+hold on;
+plot(imag(ls_real_ch4_ice_age),'r');
+title("Real L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,6);
+plot(real(ls_real_ch4_gas_age));
+hold on;
+plot(imag(ls_real_ch4_gas_age));
+title("Real L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+## Setting up the O18 plots:
+figure(o18_fig);
+subplot(4,2,1);
+axis(x_data_axis_vector);
+plot(-(o18(:,2)),o18(:,4));
+title("Gas levels over ice age");
+subplot(4,2,2);
+axis(x_data_axis_vector);
+plot(-(o18(:,3),o18(:,4));
+title("Gas levels over gas age");
+subplot(4,2,3);
+plot(real(ls_complex_o18_ice_age));
+hold on;
+plot(imag(ls_complex_o18_ice_age),'r');
+title("Complex L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,4);
+plot(real(ls_complex_o18_gas_age));
+hold on;
+plot(imag(ls_complex_o18_gas_age),'r');
+title("Complex L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+subplot(4,2,5);
+plot(real(ls_real_o18_ice_age));
+hold on;
+plot(imag(ls_real_o18_ice_age),'r');
+title("Real L-S transform of Gas/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,2,6);
+plot(real(ls_real_o18_gas_age));
+hold on;
+plot(imag(ls_real_o18_gas_age));
+title("Real L-S transform of Gas/gas age data");
+legend("Real part","Imaginary part");
+
+## Setting up Dust plots:
+figure(dust_fig);
+subplot(4,1,1);
+axis(x_data_axis_vector);
+plot(-(dust(:,2)),dust(:,3));
+title("Dust levels over ice age");
+subplot(4,1,2);
+plot(real(ls_complex_dust_ice_age));
+hold on;
+plot(imag(ls_complex_dust_ice_age),'r');
+title("Complex L-S transform of Dust/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,1,3);
+plot(real(ls_real_dust_ice_age));
+hold on;
+plot(imag(ls_real_dust_ice_age),'r');
+title("Real L-S transform of Dust/ice age data");
+legend("Real part","Imaginary part");
+
+## Setting up Deuterium plots:
+figure(deut_fig);
+subplot(4,1,1);
+axis(x_data_axis_vector);
+plot(-(deut(:,2)),deut(:,3));
+title("Deuterium levels over ice age");
+subplot(4,1,2);
+plot(real(ls_complex_deut_ice_age));
+hold on;
+plot(imag(ls_complex_deut_ice_age),'r');
+title("Complex L-S transform of Deuterium/ice age data");
+legend("Real part","Imaginary part");
+subplot(4,1,3);
+plot(real(ls_real_deut_ice_age));
+hold on;
+plot(imag(ls_real_deut_ice_age),'r');
+title("Real L-S transform of Deuterium/ice age data");
+legend("Real part","Imaginary part");
+
+co2_ch4_comparison_figure = figure("visible","off","name","CO2/CH4
+comparison");
+subplot(4,1,1);
+axes(x_data_axis_vector);
+plot(-(co2(:,2)),co2(:,4));
+hold on;
+plot(-(ch4(:,2)),ch4(:,4),'g');
+title("CO2 and CH4 data");
+legend("CO2","CH4");
+
+subplot(4,1,2);
+plot(abs(ls_complex_co2_ice_age));
+hold on;
+plot(abs(ls_complex_ch4_gas_age),'g');
+title("Abs. values of CO2 and CH4 L-S complex transforms");
+legend("CO2,CH4");
+
+
+
+
+
+
## to implement:
## - displays of all the data and flaws in trying to model with just
## using L-S data
Modified: trunk/octave-forge/extra/lssa/fastlscomplex.cc
===================================================================
--- trunk/octave-forge/extra/lssa/fastlscomplex.cc 2012-07-15 20:13:22 UTC (rev 10743)
+++ trunk/octave-forge/extra/lssa/fastlscomplex.cc 2012-07-17 00:18:06 UTC (rev 10744)
@@ -11,6 +11,8 @@
#include <iostream>
#include <exception>
+
+
ComplexRowVector flscomplex( RowVector tvec , ComplexRowVector xvec ,
double maxfreq , int octaves , int coefficients);
@@ -67,14 +69,13 @@
ComplexRowVector results = ComplexRowVector (coefficients * octaves );
- double tau, delta_tau, tau_0, tau_h, n_inv, mu,
+ double tau, delta_tau, tau_0, tau_h, n_inv, mu, te,
omega_oct, omega_multiplier, octavemax, omega_working,
- loop_tau_0, loop_delta_tau;
+ loop_tau_0, loop_delta_tau, on_1, n_1, o;
double length = ( tvec((tvec.numel()-1)) - tvec( octave_idx_type (0)));
int octave_iter, coeff_iter;
- std::complex<double> zeta, z_accumulator, exp_term, exp_multiplier, alpha;
+ std::complex<double> zeta, zz, z_accumulator, exp_term, exp_multiplier, alpha, h, *tpra, *temp_ptr_alpha, temp_alpha[12], *tprb, *temp_ptr_beta, temp_beta[12], temp_array[12], *p, x;
octave_idx_type n = tvec.numel();
- std::complex<double> temp_array[12];
for ( int array_iter = 0 ; array_iter < 12 ; array_iter++ ) {
temp_array[array_iter] = std::complex<double> ( 0 , 0 );
}
@@ -83,7 +84,7 @@
for ( int i = 1 ; i < 12 ; i++ ) {
factorial_array[i] = factorial_array[i-1] * i;
}
- n_inv = 1.0 / n;
+ n_1 = n_inv = 1.0 / n;
mu = (0.5 * M_PI)/length; // Per the article; this is in place to improve numerical accuracy if desired.
/* Viz. the paper, in which Dtau = c / omega_max, and c is stated as pi/2 for floating point processors,
* In the case of this computation, I'll go by the recommendation.
@@ -91,64 +92,118 @@
delta_tau = M_PI / ( 2 * maxfreq );
tau_0 = tvec(0) + delta_tau;
tau_h = tau_0;
- size_t precomp_subset_count = (size_t) ceil( ( tvec(tvec.numel()-1) - tvec(0) ) / ( 2 * delta_tau ) );
- // I've used size_t because it will work for my purposes without threatening undefined behaviour.
- const std::complex<double> im = std::complex<double> ( 0 , 1 );
+ te = tau_h + delta_tau;
octave_idx_type k ( 0 ); // Iterator for accessing xvec, tvec.
Precomputation_Record * precomp_records_head, *record_current, *record_tail, *record_ref, *record_next;
record_current = precomp_records_head = new Precomputation_Record;
- for ( size_t p_r_iter = 1 ; p_r_iter < precomp_subset_count ; p_r_iter++ ) {
+ for ( te = tvec(k) + (2 * delta_tau) ; ; ) {
+ x = std::complex<double>(xvec(k));
+ {
+ double t = mu*(tvec(k)-tau_h), tt;
+ p = record_current->power_series;
+ // p = 0
+ *p++ = std::complex<double>(x);
+ // p = 1
+ tt = -t;
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ // p = 2
+ tt *= t*(1.0/2.0);
+ *p++ = x*tt;
+ // p = 3
+ tt *= t*(-1.0/3.0);
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ // p = 4
+ tt *= t*(1.0/4.0);
+ *p++ = x*tt;
+ // p = 5
+ tt *= t*(-1.0/5.0);
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ // p = 6
+ tt *= t*(1.0/6.0);
+ *p++ = x*tt;
+ // p = 7
+ tt *= t*(-1.0/7.0);
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ // p = 8
+ tt *= t*(1.0/8.0);
+ *p++ = x*tt;
+ // p = 9
+ tt *= t*(-1.0/9.0);
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ // p = 10
+ tt *= t*(1.0/10.0);
+ *p++ = x*tt;
+ // p = 11
+ tt *= t*(-1.0/11.0);
+ h = x * tt;
+ *p++ = std::complex<double>(-h.imag(),h.real());
+ }
+ record_current->stored_data = true;
+ for(k++; ( k < n ) && tvec(k) < te ; k++ ) {
+ x = std::complex<double>(xvec(k));
+ {
+ double t = mu*(tvec(k)-tau_h), tt;
+ p = record_current->power_series;
+ // p = 0
+ *p++ += std::complex<double>(x);
+ // p = 1
+ tt = -t;
+ h = x * tt;
+ *p++ += std::complex<double>(- h.imag(), h.real());
+ // p = 2
+ tt *= t*(1.0/2.0);
+ *p++ += x*tt;
+ // p = 3
+ tt *= t*(-1.0/3.0);
+ h = x * tt;
+ *p++ += std::complex<double>(-h.imag(),h.real());
+ // p = 4
+ tt *= t*(1.0/4.0);
+ *p++ += x*tt;
+ // p = 5
+ tt *= t*(-1.0/5.0);
+ h = x * tt;
+ *p++ += std::complex<double>(-h.imag(),h.real());
+ // p = 6
+ tt *= t*(1.0/6.0);
+ *p++ += x*tt;
+ // p = 7
+ tt *= t*(-1.0/7.0);
+ h = x * tt;
+ *p++ += std::complex<double>(-h.imag(),h.real());
+ // p = 8
+ tt *= t*(1.0/8.0);
+ *p++ += x*tt;
+ // p = 9
+ tt *= t*(-1.0/9.0);
+ h = x * tt;
+ *p++ += std::complex<double>(-h.imag(),h.real());
+ // p = 10
+ tt *= t*(1.0/10.0);
+ *p++ += x*tt;
+ // p = 11
+ tt *= t*(-1.0/11.0);
+ h = x * tt;
+ *p++ += std::complex<double>(-h.imag(),h.real());
+ }
+ record_current->stored_data = true;
+ }
+ if( k >= n ) break;
+ tau_h = te + delta_tau;
+ te = tau_h + delta_tau;
record_current->next = new Precomputation_Record;
record_current = record_current->next;
}
record_tail = record_current;
record_current = precomp_records_head;
record_tail->next = 0;
- /* A test needs to be included for if there was a failure, but since
- * precomp_subset_count is of type size_t, it should be okay. */
- for( ; record_current != 0 ; record_current = record_current->next ) {
- for ( int j = 0 ; j < 12 ; j++ ) {
- record_current->power_series[j] = std::complex<double> ( 0 , 0 );
- } // To avoid any trouble down the line, although it is an annoyance.
- // Error is traced this far. Time to see if it's in this loop.
- while ( (k < n) && (abs(tvec(k)-tau_h) <= delta_tau) ) {
- double p;
- // alpha = std::complex<double> ( 0 , 0 );
- for ( int j = 0 ; j < 12 ; j++ ) {
- alpha.real() = xvec(k).real();
- alpha.imag() = xvec(k).imag();
- // alpha = std::complex<double> ( xvec(k).real() , xvec(k).imag() );
- // alpha *= ( -1 * im * mu * ( tvec(k) - tau_h ) ) * p;
- if ( !( j % 2 ) ) {
- if ( ! ( j % 4 ) ) {
- alpha.real() = xvec(k).real() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- alpha.imag() = xvec(k).imag() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- } else {
- alpha.real() = -1 * xvec(k).real() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- alpha.imag() = -1 * xvec(k).imag() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- }
- } else {
- if ( ! ( j % 3 ) ) {
- alpha.real() = -1 * xvec(k).imag() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- alpha.imag() = -1 * xvec(k).real() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- } else {
- alpha.real() = xvec(k).imag() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- alpha.imag() = xvec(k).real() * pow(mu,j) * pow(tvec(k)-tau_h,j) / factorial_array[j];
- }
- }
- record_current->power_series[j].real() += alpha.real();
- record_current->power_series[j].imag() += alpha.imag();
- }
- // Computes each next step of the power series for the given power series element.
- // j was reused since it was a handy inner-loop variable, even though I used it twice here.
- record_current->stored_data = true;
- k++;
- }
- tau_h += ( 2 * delta_tau );
- }
- // At this point all precomputation records have been exhausted; short-circuit is abused to avoid overflow errors.
/* Summation of coefficients for each frequency. As we have ncoeffs * noctaves elements,
* it makes sense to work from the top down, as we have omega_max by default (maxfreq)
@@ -162,40 +217,37 @@
octave_idx_type iter ( 0 );
- double real_part = 0, imag_part = 0, real_part_accumulator = 0, imag_part_accumulator = 0;
-
// Loops need to first travel over octaves, then coefficients;
for ( octave_iter = octaves ; ; omega_oct *= 0.5 , octavemax *= 0.5 , loop_tau_0 += loop_delta_tau , loop_delta_tau *= 2 ) {
- omega_working = omega_oct;
- exp_term = std::complex<double> ( cos( - omega_working * loop_tau_0 ) ,
- sin ( - omega_working * loop_tau_0 ) );
- exp_multiplier = std::complex<double> ( cos ( - 2 * omega_working * loop_delta_tau ) ,
- sin ( - 2 * omega_working * loop_delta_tau ) );
- for ( coeff_iter = 0 ; coeff_iter < coefficients ; coeff_iter++, omega_working *= omega_multiplier){
- real_part_accumulator = 0;
- imag_part_accumulator = 0;
- real_part = 0;
- imag_part = 0;
- for ( record_current = precomp_records_head ; record_current ;
+ o = omega_oct;
+ omega_working = octavemax;
+ for ( coeff_iter = 0 ; coeff_iter < coefficients ; coeff_iter++, o *= omega_multiplier, omega_working *= omega_multiplier){
+ exp_term = std::complex<double> ( cos( - omega_working * loop_tau_0 ) ,
+ sin ( - omega_working * loop_tau_0 ) );
+ exp_multiplier = std::complex<double> ( cos ( - 2 * omega_working * loop_delta_tau ) ,
+ sin ( - 2 * omega_working * loop_delta_tau ) );
+ for ( zeta = 0, record_current = precomp_records_head ; record_current ;
record_current = record_current->next, exp_term *= exp_multiplier ) {
- for ( int array_iter = 0 ; array_iter < 12 ; array_iter++ ) {
- z_accumulator = ( pow(omega_working,array_iter) * record_current->power_series[array_iter] );
- real_part_accumulator += z_accumulator.real();
- imag_part_accumulator += z_accumulator.imag();
+ if ( record_current->stored_data ) {
+ int p;
+ for ( zz = 0 , p = 0, on_1 = n_1 ; p < 12 ; p++ ) {
+ zz += record_current->power_series[p] * on_1 ;
+ on_1 *= o;
+ }
+ zeta += exp_term * zz;
}
- real_part = real_part + ( exp_term.real() * real_part_accumulator - ( exp_term.imag() * imag_part_accumulator ) );
- imag_part = imag_part + ( exp_term.imag() * real_part_accumulator + exp_term.real() * imag_part_accumulator );
}
- results(iter) = std::complex<double> ( n_inv * real_part , n_inv * imag_part );
+ results(iter) = std::complex<double> (zeta);
iter++;
}
if ( !(--octave_iter) ) break;
/* If we've already reached the lowest value, stop.
* Otherwise, merge with the next computation range.
*/
- double exp_power_series_elements[12];
+ double *exp_pse_ptr, *exp_ptr, exp_power_series_elements[12];
exp_power_series_elements[0] = 1;
+ exp_pse_ptr = exp_ptr = exp_power_series_elements;
for ( int r_iter = 1 ; r_iter < 12 ; r_iter++ ) {
exp_power_series_elements[r_iter] = exp_power_series_elements[r_iter-1]
* ( mu * loop_delta_tau) * ( 1.0 / ( (double) r_iter ) );
@@ -204,61 +256,92 @@
for ( record_current = precomp_records_head ; record_current ;
record_current = record_current->next ) {
if ( ! ( record_ref = record_current->next ) || ! record_ref->stored_data ) {
+ // In this case, there is no next record, but this record has data.
if ( record_current->stored_data ) {
- std::complex<double> temp[12];
- for( int array_init = 0 ; array_init < 12 ; array_init++ ) { temp[array_init] = std::complex<double>(0,0); }
- for( int p = 0 ; p < 12 ; p ++ ) {
- double step_floor_r = floor( ( (double) p ) / 2.0 );
- double step_floor_i = floor( ( (double) ( p - 1 ) ) / 2.0 );
- for( int q = 0 ; q < step_floor_r ; q++ ) {
- temp[p] += exp_power_series_elements[2*q] * pow((double)-1,q) * record_current->power_series[p - ( 2 * q )];
+ int p = 0;
+ for( exp_pse_ptr = exp_power_series_elements + 1 , temp_ptr_alpha = temp_alpha ; p < 12 ; p++ , exp_pse_ptr++ ) {
+ tpra = temp_ptr_alpha;
+ *(temp_ptr_alpha++) = std::complex<double>(record_current->power_series[p]);
+ for( exp_ptr = exp_power_series_elements, record_current->power_series[p] = *temp_ptr_alpha * *exp_ptr; ; ) {
+ /* This next block is from Mathias' code, and it does a few
+ * ... unsavoury things. First off, it uses conditionals with
+ * break in order to avoid potentially accessing null regions
+ * of memory, and then it does ... painful things with a few
+ * numbers. However, remembering that most of these will not
+ * actually be accessed for the first iterations, it's easier.
+ */
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ h = *tpra * *exp_ptr;
+ record_current->power_series[p].real() -= h.imag();
+ record_current->power_series[p].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ record_current->power_series[p] -= *tpra * *exp_ptr;
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ h = -*tpra * *exp_ptr;
+ record_current->power_series[p].real() -= h.imag();
+ record_current->power_series[p].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ record_current->power_series[p] += *tpra * *exp_ptr;
}
- for( int q = 0 ; q <= step_floor_i ; q++ ) {
- temp[p] += im * exp_power_series_elements[2 * q + 1] * pow((double)-1,q) * record_current->power_series[p - ( 2 * q ) - 1];
- }
}
- for ( int array_iter = 0 ; array_iter < 12 ; array_iter++ ) {
- record_current->power_series[array_iter].real() = temp[array_iter].real();
- record_current->power_series[array_iter].imag() = temp[array_iter].imag();
- }
if ( ! record_ref ) break; // Last record was reached
}
else {
record_next = record_ref;
if ( record_current->stored_data ) {
- std::complex<double> temp[12];
- for( int array_init = 0 ; array_init < 12 ; array_init++ ) { temp[array_init] = std::complex<double>(0,0); }
- for( int p = 0 ; p < 12 ; p ++ ) {
- double step_floor_r = floor( ( (double) p ) / 2.0 );
- double step_floor_i = floor( ( (double) ( p - 1 ) ) / 2.0 );
- for( int q = 0 ; q < step_floor_r ; q++ ) {
- temp[p] += exp_power_series_elements[2*q] * pow((double)-1,q) * ( record_current->power_series[p - ( 2 * q )] - record_next->power_series[p - (2*q)] );
+ int p = 0, q = 0;
+ for( exp_pse_ptr = exp_power_series_elements + 1, temp_ptr_alpha = temp_alpha, temp_ptr_beta = temp_beta; p < 12 ; p++, q++, exp_pse_ptr++ ) {
+ tpra = temp_ptr_alpha;
+ *temp_ptr_alpha++ = record_current->power_series[p] + record_next->power_series[q];
+ *temp_ptr_beta++ = record_current->power_series[p] - record_next->power_series[1];
+ tprb = temp_ptr_beta;
+ for( exp_ptr = exp_power_series_elements, record_current->power_series[p] = *tpra * *exp_ptr; ; ) {
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ tprb -= 2;
+ h = *tprb * *exp_ptr;
+ record_current->power_series[p].real() -= h.imag();
+ record_current->power_series[p].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ tpra -= 2;
+ record_current->power_series[p] -= *tpra * *exp_ptr;
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ tprb -= 2;
+ h = - *tprb * *exp_ptr;
+ record_current->power_series[p].real() -= h.imag();
+ record_current->power_series[p].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ tpra -= 2;
+ record_current->power_series[p] += *tpra * *exp_ptr;
}
- for( int q = 0 ; q <= step_floor_i ; q++ ) {
- temp[p] += im * exp_power_series_elements[2 * q + 1] * pow((double)-1,q) * ( record_current->power_series[p - ( 2 * q ) - 1] - record_next->power_series[p - ( 2 * q ) - 1 ] );
- }
}
- for ( int array_iter = 0 ; array_iter < 12 ; array_iter++ ) {
- record_current->power_series[array_iter].real() = temp[array_iter].real();
- record_current->power_series[array_iter].imag() = temp[array_iter].imag();
- }
} else {
- std::complex<double> temp[12];
- for( int array_init = 0 ; array_init < 12 ; array_init++ ) { temp[array_init] = std::complex<double>(0,0); }
- for( int p = 0 ; p < 12 ; p ++ ) {
- double step_floor_r = floor( ( (double) p ) / 2.0 );
- double step_floor_i = floor( ( (double) ( p - 1 ) ) / 2.0 );
- for( int q = 0 ; q < step_floor_r ; q++ ) {
- temp[p] += exp_power_series_elements[2*q] * pow((double)-1,q) * record_next->power_series[p - ( 2 * q )];
+ int q = 0;
+ for( exp_pse_ptr = exp_power_series_elements + 1, temp_ptr_alpha = temp_alpha, temp_ptr_beta = temp_beta; q < 12; q++, exp_pse_ptr++ ) {
+ tpra = temp_ptr_alpha;
+ *temp_ptr_alpha++ = std::complex<double>(record_next->power_series[q]);
+ for ( exp_ptr = exp_power_series_elements, record_next->power_series[q] = *tpra * *exp_ptr; ; ) {
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ h = *tpra * *exp_ptr;
+ record_next->power_series[q].real() -= h.imag();
+ record_next->power_series[q].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ record_next->power_series[q] -= *tpra * *exp_ptr;
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ h = -*tpra * *exp_ptr;
+ record_next->power_series[q].real() -= h.imag();
+ record_next->power_series[q].imag() += h.real();
+ if ( ++exp_ptr >= exp_pse_ptr ) break;
+ --tpra;
+ record_next->power_series[q] += *tpra * *exp_ptr;
}
- for( int q = 0 ; q <= step_floor_i ; q++ ) {
- temp[p] += im * exp_power_series_elements[2 * q + 1] * pow((double)-1,(q+1)) * record_next->power_series[p - ( 2 * q ) - 1];
- }
}
- for ( int array_iter = 0 ; array_iter < 12 ; array_iter++ ) {
- record_current->power_series[array_iter].real() = temp[array_iter].real();
- record_current->power_series[array_iter].imag() = temp[array_iter].imag();
- }
}
record_current->stored_data = true;
record_ref = record_next;
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <prn...@us...> - 2012-07-15 20:13:28
|
Revision: 10743
http://octave.svn.sourceforge.net/octave/?rev=10743&view=rev
Author: prnienhuis
Date: 2012-07-15 20:13:22 +0000 (Sun, 15 Jul 2012)
Log Message:
-----------
Updated texinfo header (multiline messages)
Modified Paths:
--------------
trunk/octave-forge/extra/java/inst/errordlg.m
Modified: trunk/octave-forge/extra/java/inst/errordlg.m
===================================================================
--- trunk/octave-forge/extra/java/inst/errordlg.m 2012-07-12 05:03:35 UTC (rev 10742)
+++ trunk/octave-forge/extra/java/inst/errordlg.m 2012-07-15 20:13:22 UTC (rev 10743)
@@ -16,9 +16,9 @@
## -*- texinfo -*-
## @deftypefn {Function file} {@var{P} =} errordlg (@var{MESSAGE} [,@var{TITLE}])
##
-## Displays the @var{MESSAGE} using an error dialog box.
-## @var{TITLE} (character string or cell string array for multi-line text) can
-## be used optionally to decorate the dialog caption.
+## Displays the @var{MESSAGE} (character string or cell string array for
+## multi-line text) using an error dialog box.
+## @var{TITLE} can be used optionally to decorate the dialog caption.
## The return value is always 1.
##
## @end deftypefn
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <be...@us...> - 2012-07-12 05:03:42
|
Revision: 10742
http://octave.svn.sourceforge.net/octave/?rev=10742&view=rev
Author: benjf5
Date: 2012-07-12 05:03:35 +0000 (Thu, 12 Jul 2012)
Log Message:
-----------
Fixed o18.csv (had the wrong data before) and added the beginning of a demonstration script.
Modified Paths:
--------------
trunk/octave-forge/extra/lssa/data/o18.csv
trunk/octave-forge/extra/lssa/lombcoeff.m
trunk/octave-forge/extra/lssa/lscomplex.m
trunk/octave-forge/extra/lssa/lsreal.m
Added Paths:
-----------
trunk/octave-forge/extra/lssa/SampleScriptWithVostokData
Added: trunk/octave-forge/extra/lssa/SampleScriptWithVostokData
===================================================================
--- trunk/octave-forge/extra/lssa/SampleScriptWithVostokData (rev 0)
+++ trunk/octave-forge/extra/lssa/SampleScriptWithVostokData 2012-07-12 05:03:35 UTC (rev 10742)
@@ -0,0 +1,53 @@
+## Copyright (C) 2012 Benjamin Lewis
+## Licensed under the GNU GPL v2
+
+## This is just a sample script to introduce the purpose and usage of
+## the Lomb-Scargle Least Squares method with experimental data, here
+## using the Vostok ice core data collected and measured by J.R. Petit
+## et. al. and published in Nature; also available from the NOAA's
+## Paleoclimatology pages here:
+## <http://www.ncdc.noaa.gov/paleo/icecore/antarctica/vostok/vostok_data.html>.
+
+co2 = csvread("./data/co2.csv")(2:end,2:end);
+ch4 = csvread("./data/ch4.csv")(2:end,2:end);
+o18 = csvread("./data/o18.csv")(2:end,2:end);
+deut = csvread("./data/deut.csv")(2:end,2:end);
+dust = csvread("./data/dust.csv")(2:end,2:end);
+## The limited ranges are to deal with artifacts from the extraction of
+## the R data, notably that it leaves an extra column on the front of 0s
+## and the first row is text that Octave refuses to process.
+
+## Columns in co2 are Depth, Ice Age, Gas Age, and CO2 Concentration.
+## Columns in ch4 are Depth, Ice Age, Gas Age, and CH4 Concentration.
+## Columns in o18 are Depth, Ice Age, Gas Age, and Atmospheric O18.
+## Columns in dust are Depth, Ice Age, Dust Concentration.
+## Columns in deut are Depth, Ice Age, D concentration, and DeltaTS.
+
+co2_fig = figure("visible","off","name","CO2");
+ch4_fig = figure("visible","off","name","CH4");
+o18_fig = figure("visible","off","name","O18");
+deut_fig = figure("visible","off","name","Deuterium");
+dust_fig = figure("visible","off","name","Dust");
+## Generates figures and attaches handles to them for easy access.
+
+## Now we need some data to display; I'll run a few functions.
+ls_complex_co2_ice_age = lscomplex(co2(:,2),co2(:,4),1,100,20);
+ls_complex_co2_gas_age = lscomplex(co2(:,3),co2(:,4),1,100,20);
+ls_real_co2_ice_age = lsreal(co2(:,2),co2(:,4),1,100,20);
+ls_real_co2_gas_age = lsreal(co2(:,3),co2(:,4),1,100,20);
+
+figure(co2_fig);
+subplot(4,2,1);
+plot(-(co2(:,2)),co2(:,4));
+title("Gas levels over ice age");
+subplot(4,2,2);
+plot(-(co2(:,3),co2(:,4));
+title("Gas levels over gas age");
+subplot(4,2,3);
+
+
+## to implement:
+## - displays of all the data and flaws in trying to model with just
+ ## using L-S data
+## - correlations of every data set with every other data set
+## - Comparing ls* results to periodogram results
\ No newline at end of file
Modified: trunk/octave-forge/extra/lssa/data/o18.csv
===================================================================
--- trunk/octave-forge/extra/lssa/data/o18.csv 2012-07-10 14:20:31 UTC (rev 10741)
+++ trunk/octave-forge/extra/lssa/data/o18.csv 2012-07-12 05:03:35 UTC (rev 10742)
@@ -1,523 +1,319 @@
-"","depth.corrected","Ice.age..GT4.","Dust.Conc..ppm."
-"1",124,4509,0.025
-"2",125,4552,0.018
-"3",128,4690,0.015
-"4",140,5259,0.033
-"5",149,5674,0.028
-"6",166,6486,0.013
-"7",174,6874,0.02
-"8",175,6924,0.018
-"9",182,7267,0.038
-"10",191,7697,0.023
-"11",192,7745,0.065
-"12",199,8091,0.038
-"13",220,9107,0.04
-"14",224,9298,0.05
-"15",244,10265,0.018
-"16",249,10515,0.03
-"17",256,10861,0.033
-"18",260,11053,0.045
-"19",274,11749,0.05
-"20",278,11973,0.025
-"21",288,12569,0.075
-"22",296,13055,0.075
-"23",299,13237,0.04
-"24",309,13828,0.033
-"25",319,14404,0.085
-"26",324,14713,0.163
-"27",325,14775,0.18
-"28",329,15032,0.17
-"29",349,16426,0.638
-"30",350,16502,0.773
-"31",352,16653,0.32
-"32",363,17544,1.098
-"33",374,18446,0.22
-"34",379,18870,1.21
-"35",399,20528,0.973
-"36",402,20777,0.71
-"37",403,20859,1.735
-"38",417,22023,0.733
-"39",424,22629,1.213
-"40",433,23412,1.745
-"41",447,24619,1.573
-"42",449,24781,0.668
-"43",451,24941,0.45
-"44",460,25697,0.968
-"45",474,26849,1.48
-"46",475,26933,1.353
-"47",492,28420,0.51
-"48",499,29006,0.415
-"49",500,29085,0.365
-"50",502,29249,0.783
-"51",524,31130,0.708
-"52",549,33209,0.495
-"53",574,35185,0.09
-"54",575,35261,0.113
-"55",598,37066,0.49
-"56",599,37148,0.603
-"57",601,37313,0.293
-"58",624,39185,0.203
-"59",625,39263,0.22
-"60",649,41191,0.188
-"61",650,41268,0.21
-"62",674,43073,0.238
-"63",675,43150,0.185
-"64",698,45055,0.19
-"65",699,45139,0.468
-"66",700,45225,0.375
-"67",724,47295,0.28
-"68",749,49350,0.128
-"69",750,49430,0.208
-"70",774,51306,0.073
-"71",775,51384,0.15
-"72",777,51549,0.295
-"73",799,53355,0.19
-"74",800,53436,0.165
-"75",802,53606,0.155
-"76",823,55377,0.278
-"77",824,55457,0.205
-"78",848,57289,0.21
-"79",849,57362,0.145
-"80",873,59262,0.833
-"81",874,59345,0.468
-"82",898,61406,0.475
-"83",899,61496,0.378
-"84",923,63688,0.98
-"85",924,63776,0.708
-"86",948,66045,1.1
-"87",949,66144,0.748
-"88",955,66708,0.425
-"89",973,68424,1.285
-"90",974,68520,0.86
-"91",984,69413,0.328
-"92",993,70243,0.305
-"93",998,70683,0.16
-"94",999,70766,0.063
-"95",1005,71264,0.085
-"96",1015,72070,0.098
-"97",1023,72753,0.155
-"98",1024,72840,0.053
-"99",1034,73706,0.19
-"100",1041,74322,0.205
-"101",1048,74903,0.063
-"102",1049,74985,0.103
-"103",1073,76854,0.073
-"104",1074,76933,0.075
-"105",1098,78937,0.178
-"106",1099,79020,0.233
-"107",1123,80869,0.28
-"108",1124,80943,0.193
-"109",1144,82352,0.078
-"110",1148,82647,0.188
-"111",1149,82721,0.223
-"112",1173,84378,0.068
-"113",1174,84448,0.058
-"114",1198,86016,0.065
-"115",1199,86083,0.078
-"116",1223,87766,0.09
-"117",1224,87837,0.09
-"118",1248,89597,0.19
-"119",1249,89672,0.208
-"120",1273,91480,0.113
-"121",1274,91547,0.12
-"122",1298,93254,0.27
-"123",1299,93325,0.238
-"124",1323,95000,0.07
-"125",1324,95070,0.103
-"126",1348,96722,0.048
-"127",1349,96791,0.043
-"128",1350,96861,0.053
-"129",1373,98505,0.075
-"130",1374,98572,0.068
-"131",1375,98640,0.005
-"132",1398,100232,0.17
-"133",1399,100301,0.14
-"134",1423,101868,0.133
-"135",1424,101930,0.11
-"136",1448,103481,0.248
-"137",1449,103553,0.18
-"138",1450,103625,0.315
-"139",1451,103697,0.258
-"140",1473,105271,0.058
-"141",1474,105342,0.058
-"142",1475,105413,0.058
-"143",1498,107048,0.128
-"144",1499,107121,0.098
-"145",1500,107194,0.085
-"146",1523,109121,0.093
-"147",1524,109206,0.05
-"148",1525,109292,0.065
-"149",1548,111159,0.07
-"150",1549,111239,0.043
-"151",1550,111319,0.01
-"152",1573,113195,0.05
-"153",1574,113272,0.078
-"154",1575,113348,0.043
-"155",1582,113891,0.068
-"156",1591,114562,0.058
-"157",1598,115077,0.035
-"158",1599,115150,0.08
-"159",1606,115617,0.028
-"160",1616,116296,0.048
-"161",1623,116762,0.033
-"162",1624,116827,0.028
-"163",1625,116891,0.073
-"164",1632,117326,0.068
-"165",1635,117514,0.015
-"166",1643,118012,0.035
-"167",1648,118314,0.038
-"168",1649,118376,0.058
-"169",1655,118735,0.11
-"170",1665,119345,0.045
-"171",1673,119831,0.04
-"172",1674,119894,0.038
-"173",1675,119955,0.073
-"174",1682,120379,0.053
-"175",1692,120959,0.085
-"176",1698,121314,0.02
-"177",1699,121373,0.013
-"178",1709,121951,0.018
-"179",1716,122345,0.213
-"180",1723,122742,0.055
-"181",1724,122801,0.053
-"182",1732,123268,0.06
-"183",1739,123658,0.04
-"184",1749,124209,0.12
-"185",1756,124599,0.048
-"186",1761,124877,0.043
-"187",1768,125271,0.028
-"188",1783,126103,0.078
-"189",1798,126904,0.048
-"190",1799,126955,0.06
-"191",1800,127006,0.26
-"192",1816,127810,0.043
-"193",1823,128161,0.028
-"194",1824,128210,0.038
-"195",1825,128259,0.045
-"196",1830,128501,0.018
-"197",1848,129428,0.135
-"198",1849,129486,0.115
-"199",1850,129545,0.02
-"200",1866,130410,0.035
-"201",1872,130759,0.015
-"202",1873,130819,0.113
-"203",1874,130880,0.038
-"204",1875,130942,0.14
-"205",1879,131190,0.045
-"206",1883,131439,0.033
-"207",1886,131641,0.023
-"208",1890,131908,0.02
-"209",1895,132264,0.03
-"210",1898,132477,0.048
-"211",1899,132548,0.045
-"212",1900,132619,0.078
-"213",1901,132691,0.055
-"214",1902,132763,0.038
-"215",1903,132836,0.055
-"216",1906,133065,0.025
-"217",1907,133142,0.06
-"218",1912,133528,0.06
-"219",1919,134091,0.105
-"220",1921,134266,0.095
-"221",1923,134447,0.5
-"222",1924,134538,0.12
-"223",1932,135308,0.298
-"224",1933,135408,0.095
-"225",1934,135507,0.105
-"226",1940,136103,0.16
-"227",1943,136411,0.123
-"228",1945,136614,0.245
-"229",1948,136923,0.883
-"230",1949,137026,0.86
-"231",1953,137446,0.158
-"232",1954,137549,0.183
-"233",1963,138532,0.288
-"234",1964,138644,0.235
-"235",1966,138868,0.553
-"236",1970,139308,0.338
-"237",1973,139643,0.848
-"238",1974,139756,0.498
-"239",1979,140319,0.405
-"240",1981,140542,0.553
-"241",1984,140876,0.688
-"242",1989,141434,0.245
-"243",1990,141547,0.438
-"244",1991,141661,0.303
-"245",1996,142218,0.565
-"246",1998,142440,1.01
-"247",1999,142551,0.753
-"248",2004,143117,0.625
-"249",2009,143675,0.74
-"250",2014,144237,0.46
-"251",2023,145226,0.973
-"252",2024,145337,0.673
-"253",2025,145449,0.738
-"254",2032,146242,0.528
-"255",2043,147459,0.423
-"256",2049,148141,0.77
-"257",2052,148471,0.158
-"258",2063,149679,0.265
-"259",2069,150333,0.55
-"260",2074,150868,0.418
-"261",2080,151502,0.325
-"262",2081,151611,0.318
-"263",2089,152513,0.313
-"264",2098,153523,0.871
-"265",2103,154094,0.305
-"266",2111,155029,1.177
-"267",2119,155983,1.228
-"268",2128,157080,0.671
-"269",2135,157909,0.715
-"270",2143,158815,0.426
-"271",2151,159651,0.13
-"272",2160,160626,0.237
-"273",2168,161545,0.728
-"274",2176,162440,0.238
-"275",2184,163318,0.22
-"276",2191,164136,0.36
-"277",2199,165148,0.618
-"278",2202,165522,0.092
-"279",2208,166284,0.271
-"280",2217,167414,0.115
-"281",2225,168391,0.154
-"282",2232,169201,0.073
-"283",2241,170221,0.151
-"284",2248,171040,0.284
-"285",2256,172004,0.307
-"286",2264,173027,0.025
-"287",2272,173988,0.105
-"288",2280,174920,0.142
-"289",2287,175746,0.232
-"290",2294,176630,0.391
-"291",2308,178417,0.341
-"292",2325,180546,0.202
-"293",2325,180546,0.067
-"294",2333,181502,0.138
-"295",2342,182694,0.236
-"296",2348,183525,0.561
-"297",2356,184641,0.715
-"298",2364,185788,0.761
-"299",2373,187067,0.274
-"300",2387,189030,0.147
-"301",2398,190469,0.05
-"302",2400,190705,0.055
-"303",2404,191156,0.016
-"304",2410,191895,0.114
-"305",2414,192382,0.067
-"306",2420,193152,0.142
-"307",2425,193816,0.141
-"308",2430,194444,0.052
-"309",2434,194926,0.104
-"310",2440,195676,0.113
-"311",2444,196202,0.221
-"312",2451,197086,0.067
-"313",2454,197462,0.111
-"314",2459,198075,0.114
-"315",2462,198426,0.074
-"316",2471,199452,0.047
-"317",2473,199682,0.052
-"318",2475,199904,0.063
-"319",2477,200124,0.048
-"320",2479,200345,0.021
-"321",2491,201635,0.03
-"322",2502,202826,0.029
-"323",2505,203155,0.024
-"324",2505.41,203200,0.024
-"325",2508.41,203527,0.023
-"326",2509,203592,0.021
-"327",2512.41,203980,0.02
-"328",2515,204278,0.029
-"329",2515.41,204325,0.028
-"330",2521,204973,0.022
-"331",2524.41,205386,0.032
-"332",2525,205458,0.024
-"333",2527.41,205760,0.032
-"334",2531,206232,0.025
-"335",2531.41,206287,0.028
-"336",2534.41,206681,0.031
-"337",2537.41,207072,0.042
-"338",2540.41,207463,0.041
-"339",2543.41,207857,0.023
-"340",2545,208061,0.029
-"341",2546.41,208244,0.023
-"342",2549.41,208623,0.031
-"343",2555.41,209364,0.028
-"344",2559.41,209864,0.021
-"345",2564.41,210466,0.023
-"346",2568.41,210941,0.032
-"347",2569.41,211057,0.015
-"348",2574.41,211651,0.022
-"349",2576.41,211892,0.04
-"350",2579.41,212237,0.029
-"351",2583.41,212700,0.03
-"352",2585.41,212928,0.022
-"353",2597.41,214297,0.031
-"354",2600.41,214639,0.018
-"355",2603.41,214978,0.06
-"356",2607.41,215432,0.022
-"357",2612.41,216003,0.036
-"358",2612.41,216003,0.054
-"359",2615.41,216369,0.029
-"360",2619.41,216839,0.032
-"361",2629.41,217989,0.022
-"362",2634.41,218602,0.035
-"363",2637.41,219000,0.082
-"364",2640.41,219393,0.041
-"365",2642.41,219660,0.053
-"366",2646.41,220226,0.093
-"367",2651.41,220967,0.077
-"368",2653.41,221272,0.069
-"369",2660.41,222426,0.119
-"370",2664.41,223135,0.138
-"371",2667.41,223681,0.194
-"372",2670.41,224240,0.158
-"373",2674.41,224959,0.029
-"374",2677.41,225469,0.085
-"375",2680.41,225971,0.069
-"376",2684.41,226653,0.036
-"377",2691.41,227958,0.076
-"378",2694.41,228488,0.07
-"379",2697.41,229052,0.045
-"380",2701.41,229802,0.073
-"381",2703.41,230161,0.031
-"382",2706.41,230691,0.076
-"383",2715.41,232206,0.036
-"384",2715.41,232206,0.043
-"385",2718.41,232682,0.039
-"386",2719.41,232838,0.007
-"387",2721.41,233146,0.033
-"388",2721.41,233146,0.055
-"389",2723.41,233453,0.056
-"390",2733.41,234896,0.096
-"391",2736.41,235329,0.038
-"392",2748.41,236970,0.058
-"393",2751.41,237338,0.033
-"394",2754.41,237690,0.083
-"395",2756.41,237912,0.079
-"396",2762.41,238578,0.012
-"397",2765.41,238932,0.015
-"398",2768.41,239298,0.014
-"399",2771.41,239698,0.011
-"400",2774.41,240141,0.029
-"401",2777.41,240614,0.014
-"402",2780.41,241118,0.06
-"403",2780.41,241118,0.067
-"404",2783.41,241674,0.029
-"405",2786.41,242267,0.083
-"406",2789.41,242850,0.039
-"407",2792.41,243406,0.023
-"408",2795.41,243999,0.06
-"409",2797.41,244424,0.131
-"410",2797.41,244424,0.052
-"411",2798.41,244639,0.156
-"412",2801.41,245283,0.104
-"413",2807.41,246574,0.378
-"414",2816.41,248545,0.191
-"415",2819.41,249219,0.116
-"416",2825.41,250501,0.078
-"417",2828.41,251104,0.169
-"418",2834.41,252422,0.208
-"419",2837.41,253075,0.12
-"420",2840.41,253763,0.115
-"421",2846.41,255056,0.219
-"422",2849.41,255711,0.053
-"423",2852.41,256387,0.223
-"424",2858.41,257827,0.439
-"425",2864.41,259340,1.012
-"426",2868.41,260380,0.64
-"427",2871.41,261179,0.942
-"428",2874.41,261973,0.167
-"429",2878.41,263043,0.554
-"430",2883.41,264383,0.791
-"431",2888.41,265707,0.558
-"432",2897.41,268050,0.093
-"433",2903.41,269532,0.168
-"434",2907.41,270507,0.077
-"435",2913.41,271920,0.017
-"436",2918.41,273037,0.085
-"437",2922.41,273918,0.123
-"438",2928.41,275216,0.086
-"439",2937.41,277120,0.045
-"440",2943.41,278456,0.08
-"441",2948.41,279632,0.146
-"442",2954.41,281175,0.184
-"443",2958.41,282253,0.14
-"444",2968.41,284925,0.204
-"445",2974.41,286639,0.041
-"446",2978.41,287726,0.059
-"447",2981.41,288530,0.059
-"448",2987.41,290214,0.066
-"449",2993.41,291830,0.054
-"450",2998.41,293076,0.042
-"451",3001.41,293893,0.063
-"452",3005.41,295035,0.045
-"453",3008.41,295862,0.09
-"454",3017.41,298327,0.257
-"455",3021.41,299383,0.097
-"456",3025.41,300397,0.078
-"457",3028.41,301156,0.053
-"458",3031.41,301917,0.032
-"459",3034.41,302691,0.013
-"460",3035.41,302937,0.029
-"461",3038.41,303657,0.022
-"462",3044.41,305288,0.021
-"463",3047.41,306200,0.111
-"464",3050.41,307131,0.102
-"465",3053.41,308031,0.026
-"466",3056.41,308907,0.02
-"467",3058.41,309483,0.009
-"468",3061.41,310330,0.036
-"469",3065.41,311411,0.022
-"470",3068.41,312208,0.019
-"471",3077.41,314543,0.021
-"472",3080.41,315297,0.045
-"473",3085.41,316524,0.02
-"474",3086.41,316765,0.027
-"475",3088.41,317255,0.044
-"476",3091.41,317999,0.037
-"477",3098.41,319708,0.019
-"478",3107.41,321875,0.032
-"479",3110.41,322513,0.014
-"480",3118.41,324218,0.008
-"481",3124.41,325629,0.012
-"482",3128.41,326761,0.01
-"483",3131.41,327730,0.046
-"484",3137.41,329848,0.117
-"485",3140.41,331050,0.142
-"486",3147.41,334306,0.586
-"487",3148.41,334805,1.035
-"488",3152.41,336766,0.395
-"489",3155.41,338233,0.366
-"490",3158.41,339706,0.439
-"491",3161.41,341187,0.164
-"492",3164.41,342594,0.034
-"493",3172.41,346523,0.277
-"494",3178.41,349620,0.401
-"495",3188.41,355187,0.38
-"496",3194.41,358466,0.24
-"497",3197.41,360071,0.155
-"498",3205.41,364333,0.076
-"499",3208.41,365962,0.268
-"500",3211.41,367507,0.058
-"501",3214.41,369034,0.091
-"502",3218.41,371188,0.224
-"503",3221.41,372991,0.208
-"504",3224.41,374711,0.071
-"505",3227.41,376381,0.085
-"506",3233.41,379524,0.019
-"507",3236.41,381227,0.125
-"508",3238.41,382318,0.033
-"509",3248.41,388065,0.122
-"510",3254.41,391928,0.114
-"511",3257.41,393908,0.157
-"512",3260.41,395792,0.023
-"513",3263.41,397482,0.023
-"514",3266.41,399134,0.019
-"515",3268.41,400228,0.021
-"516",3274.41,403501,0.033
-"517",3278.41,405617,0.016
-"518",3287.41,410200,0.018
-"519",3288.41,410683,0.012
-"520",3300.41,417095,0.013
-"521",3306.41,420527,0.015
-"522",3308.41,421761,0.012
+"","Depth.corrected","Ice.age..GT4.","gas.age..GT4.","Atm.O18"
+"1",114.8,4095,724,-0.06
+"2",116.3,4168,813,0.03
+"3",120.3,4352,987,-0.06
+"4",125.2,4561,1184,-0.02
+"5",130.1,4791,1418,0.08
+"6",135.3,5039,1675,0.04
+"7",139.7,5245,1932,-0.05
+"8",139.8,5250,1935,0.14
+"9",143.7,5429,2108,-0.09
+"10",149.6,5702,2370,-0.05
+"11",152.8,5854,2541,-0.09
+"12",156.9,6048,2745,-0.08
+"13",161.8,6279,3041,-0.11
+"14",165.1,6441,3214,-0.13
+"15",166.8,6526,3306,-0.05
+"16",169.7,6665,3453,-0.15
+"17",174.2,6884,3679,0.06
+"18",174.6,6905,3695,-0.14
+"19",177.1,7028,3814,0.04
+"20",184.3,7379,4262,0.15
+"21",184.9,7408,4299,-0.19
+"22",205.6,8403,5207,-0.18
+"23",213.1,8768,5529,-0.13
+"24",213.8,8801,5559,-0.2
+"25",224,9298,6009,-0.24
+"26",230,9589,6278,-0.05
+"27",230.4,9608,6293,-0.11
+"28",250,10564,7313,0.1
+"29",252.9,10711,7475,-0.21
+"30",270.8,11579,8285,0.04
+"31",271.3,11606,8318,0.06
+"32",280.2,12099,8808,0.04
+"33",280.9,12139,8856,-0.11
+"34",281.3,12161,8879,-0.38
+"35",290.4,12719,9443,0.24
+"36",291.3,12773,9503,0.33
+"37",294.5,12964,9655,0.25
+"38",298.5,13208,9898,0.29
+"39",299.5,13267,9952,0.24
+"40",311.1,13944,10544,0.58
+"41",311.6,13971,10570,0.32
+"42",322.6,14627,11066,0.88
+"43",322.8,14639,11073,0.93
+"44",330.5,15133,11288,0.9
+"45",336.7,15549,11499,0.72
+"46",345.8,16185,11909,0.87
+"47",355.9,16962,12609,1.11
+"48",363.9,17618,13299,1.02
+"49",372,18283,13850,1.21
+"50",372.1,18291,13853,0.97
+"51",383,19199,14302,1.06
+"52",396.1,20286,14927,1.23
+"53",400.2,20628,15123,0.91
+"54",410.2,21459,15615,0.92
+"55",422.3,22483,16185,1.16
+"56",422.4,22491,16189,1.33
+"57",431.4,23270,16719,0.97
+"58",446.5,24579,17937,0.98
+"59",456.5,25402,18797,0.89
+"60",467.6,26310,19661,0.74
+"61",476.7,27079,20435,0.98
+"62",486.7,27969,21288,0.82
+"63",495.8,28752,21999,0.67
+"64",496.8,28832,22076,0.76
+"65",502.4,29282,22590,0.9
+"66",511.9,30114,23561,0.49
+"67",515.5,30411,23928,1.02
+"68",528,31473,25057,0.44
+"69",537.1,32227,25812,0.44
+"70",547.1,33047,26477,0.41
+"71",562,34281,27904,0.37
+"72",571.9,35024,28784,0.33
+"73",576.4,35366,29056,0.47
+"74",576.5,35374,29062,0.52
+"75",597.8,37049,30942,0.35
+"76",615.8,38518,32740,0.64
+"77",627.8,39483,34133,0.58
+"78",631.5,39774,34495,0.58
+"79",636.8,40205,34893,0.57
+"80",637.8,40291,34947,0.29
+"81",645.7,40930,35449,0.49
+"82",659.7,41977,36237,0.46
+"83",671.7,42892,37241,0.49
+"84",677.7,43359,37752,0.57
+"85",686.6,44086,38599,0.52
+"86",696.6,44939,39582,0.43
+"87",709.6,46070,41011,0.57
+"88",721.5,47092,41928,0.65
+"89",734.5,48176,42711,0.46
+"90",747.5,49233,43509,0.5
+"91",759.4,50169,44336,0.43
+"92",771.4,51117,45279,0.43
+"93",784.4,52146,46655,0.29
+"94",796.3,53133,47826,0.18
+"95",809.3,54224,49156,0.31
+"96",821.3,55243,50188,0.56
+"97",834.2,56254,51152,0.37
+"98",846.2,57154,51881,0.33
+"99",859.2,58148,52763,0.22
+"100",871.1,59101,53867,0.3
+"101",884.1,60185,55442,0.26
+"102",892.2,60887,56257,0.6
+"103",897.1,61326,56713,0.31
+"104",909,62401,57529,0.77
+"105",921,63512,58373,0.71
+"106",934,64698,59346,0.8
+"107",945.9,65847,60229,0.72
+"108",959.9,67169,61370,0.77
+"109",970.9,68222,62301,0.67
+"110",983.8,69393,63515,0.83
+"111",1020.5,72532,66581,0.77
+"112",1032.6,73589,67850,0.78
+"113",1045.6,74703,69299,0.71
+"114",1056.7,75602,70607,0.65
+"115",1069.8,76605,71543,0.56
+"116",1082.8,77674,72407,0.58
+"117",1095.9,78766,73689,0.42
+"118",1108,79714,74974,0.35
+"119",1121,80719,75937,0.31
+"120",1133.1,81599,76634,0.37
+"121",1146.2,82515,77505,0.14
+"122",1158.3,83352,78596,0.07
+"123",1171.3,84259,79753,0.02
+"124",1183.4,85062,80841,-0.05
+"125",1196.5,85915,81826,-0.03
+"126",1209.5,86812,82827,-0.25
+"127",1221.6,87666,83681,-0.21
+"128",1233.7,88548,84668,0.21
+"129",1234.7,88623,84733,0.13
+"130",1245.6,89419,85433,0.37
+"131",1246.8,89512,85498,0.05
+"132",1258.8,90423,86174,0.31
+"133",1284,92259,87755,0.26
+"134",1297,93180,88597,0.37
+"135",1309.1,94033,89355,0.23
+"136",1322.2,94946,90430,0.33
+"137",1334.2,95775,91475,0.43
+"138",1347.3,96673,92352,0.47
+"139",1358.4,97469,93181,0.56
+"140",1384.5,99310,95132,0.48
+"141",1408.7,100949,96632,0.41
+"142",1432.8,102476,98365,0.02
+"143",1444.9,103266,99254,-0.11
+"144",1459,104266,100462,-0.07
+"145",1473,105271,101636,-0.11
+"146",1485.1,106121,102314,-0.01
+"147",1496.2,106917,102869,-0.14
+"148",1507.2,107784,103548,-0.14
+"149",1512.3,108224,104017,0.06
+"150",1520.3,108897,104686,0.14
+"151",1529.3,109650,105444,0.33
+"152",1539.4,110459,106102,0.46
+"153",1551.5,111439,106686,0.44
+"154",1566.2,112644,107470,0.59
+"155",1592,114639,109824,0.73
+"156",1636.3,117593,113370,0.73
+"157",1656.1,118802,115226,0.39
+"158",1683.3,120457,117267,0.28
+"159",1693.1,121024,117891,-0.07
+"160",1716.7,122384,119319,-0.21
+"161",1732.7,123308,120424,-0.19
+"162",1757,124654,121904,-0.25
+"163",1784.3,126177,123521,-0.34
+"164",1823.5,128186,125622,-0.22
+"165",1831.9,128592,126162,-0.21
+"166",1845.2,129284,127030,0.17
+"167",1858,129971,127752,0.27
+"168",1868.5,130554,128325,0.65
+"169",1882.3,131393,128995,1
+"170",1883.4,131465,129069,0.9
+"171",1907.6,133188,130379,1.13
+"172",1934.5,135557,131959,1.23
+"173",1954.6,137611,133340,1.21
+"174",1955.8,137733,133418,1.12
+"175",1965.5,138812,133941,1.16
+"176",1982.3,140689,135202,1.03
+"177",1990.8,141637,135994,0.86
+"178",2016.4,144496,138408,0.91
+"179",2058.8,149220,143362,0.59
+"180",2092.1,152863,147205,0.6
+"181",2104,154207,148765,0.74
+"182",2112.1,155163,149816,0.82
+"183",2132.2,157583,151981,0.98
+"184",2142,158702,153069,0.89
+"185",2161.1,160751,154894,0.88
+"186",2170,161776,156102,0.8
+"187",2191,164136,159270,0.87
+"188",2201.1,165412,160281,0.79
+"189",2222.4,168089,162810,0.66
+"190",2236.1,169649,163928,0.38
+"191",2254.1,171772,166309,0.09
+"192",2281.1,175045,169989,0.11
+"193",2293,176501,171357,-0.13
+"194",2301,177517,172412,-0.1
+"195",2312.2,178992,174181,0.04
+"196",2332,181382,176144,-0.01
+"197",2341.1,182570,177430,0.24
+"198",2353.1,184243,179606,0.47
+"199",2372.1,186941,181627,0.88
+"200",2391.1,189582,184076,0.87
+"201",2404.5,191213,185778,0.72
+"202",2415.4,192560,187408,0.53
+"203",2425.6,193893,189454,0.47
+"204",2435.6,195126,190962,0.38
+"205",2444.1,196215,191858,0.12
+"206",2455.2,197612,193211,-0.11
+"207",2465.2,198790,194611,-0.03
+"208",2472.2,199590,195320,-0.11
+"209",2492.2,201763,198189,-0.14
+"210",2502.4,202870,199497,-0.06
+"211",2509,203592,200370,-0.15
+"212",2512.4,203979,200759,0.04
+"213",2514.8,204256,201104,-0.02
+"214",2520.4,204903,201707,0.2
+"215",2520.8,204948,201751,0.03
+"216",2530.7,206189,202917,0.39
+"217",2531.7,206327,203049,0.15
+"218",2536.8,206995,203585,0.39
+"219",2546.4,208243,204656,0.66
+"220",2549.7,208659,204934,0.52
+"221",2554.8,209293,205426,0.51
+"222",2572.7,211450,207759,0.47
+"223",2586.7,213077,209691,0.35
+"224",2595.7,214103,210894,0.26
+"225",2607.7,215464,212419,0.06
+"226",2616.7,216523,213539,-0.06
+"227",2625.7,217550,214611,-0.12
+"228",2634.7,218639,215622,-0.14
+"229",2644.7,219976,217045,0
+"230",2649.7,220719,217606,0.13
+"231",2658.7,222130,218640,0.05
+"232",2667.2,223642,219735,0.08
+"233",2677.2,225433,221031,0.14
+"234",2686.2,226975,222104,0.44
+"235",2692.2,228100,223181,0.51
+"236",2699.2,229392,224894,0.77
+"237",2706.2,230654,225994,0.84
+"238",2712.2,231682,226779,0.68
+"239",2718.2,232649,227928,0.65
+"240",2734.2,235013,231061,0.37
+"241",2736.2,235299,231478,0.27
+"242",2739.2,235722,232087,0.24
+"243",2743.2,236277,232849,0.09
+"244",2745.2,236545,233199,-0.01
+"245",2763.9,238752,235990,0.03
+"246",2772.9,239913,237741,0.48
+"247",2778.9,240859,238490,0.76
+"248",2787.9,242559,239486,0.75
+"249",2796.9,244313,240495,0.62
+"250",2802.9,245599,241301,0.63
+"251",2811.9,247548,243128,0.84
+"252",2820.9,249542,244726,0.75
+"253",2828.9,251206,246382,0.63
+"254",2835.9,252757,247944,0.51
+"255",2844.9,254756,250337,0.45
+"256",2853.9,256735,252048,0.37
+"257",2859.9,258202,253691,0.47
+"258",2869.9,260776,255933,0.59
+"259",2877.9,262907,257646,0.69
+"260",2885.9,265028,259567,0.66
+"261",2891.9,266614,261142,0.66
+"262",2900.9,268897,263563,0.73
+"263",2909.9,271097,266106,0.56
+"264",2918.9,273148,268568,0.4
+"265",2924.9,274459,270247,0.2
+"266",2933.9,276388,272620,0
+"267",2942.9,278337,274784,-0.15
+"268",2948.9,279752,276308,-0.09
+"269",2957.9,282123,278250,0
+"270",2965.9,284181,279851,-0.03
+"271",2971.9,285948,281656,0.09
+"272",2980.9,288387,283752,0.53
+"273",2989.9,290913,286629,0.6
+"274",2995.9,292475,288134,0.59
+"275",3004.9,294900,291089,0.45
+"276",3013.9,297367,293273,0.25
+"277",3019.9,298989,295052,-0.02
+"278",3028.9,301282,297568,0.02
+"279",3037.9,303533,300215,0.19
+"280",3043.9,305141,301918,0.26
+"281",3052.9,307880,304243,0.57
+"282",3061.9,310471,306589,0.59
+"283",3067.9,312072,308551,0.56
+"284",3076.9,314415,311344,0.38
+"285",3085.9,316647,313862,0.18
+"286",3091.9,318123,315533,-0.03
+"287",3100.9,320329,317837,-0.21
+"288",3109.9,322409,319956,-0.16
+"289",3115.9,323676,321706,0.03
+"290",3124.9,325754,323802,0.77
+"291",3130.9,327565,325218,1.02
+"292",3142.9,332112,328540,1.12
+"293",3147.9,334555,330355,1.19
+"294",3154.9,337989,332841,1.21
+"295",3160.9,340942,336045,1.05
+"296",3163.9,342356,337494,0.92
+"297",3172.9,346772,342238,0.7
+"298",3178.9,349884,345278,0.67
+"299",3184.9,353220,348373,0.68
+"300",3190.9,356559,351525,0.67
+"301",3196.9,359819,355044,0.63
+"302",3199.9,361384,356764,0.41
+"303",3208.9,366230,361967,0.07
+"304",3211.9,367764,363534,0.11
+"305",3220.9,372699,368878,0.23
+"306",3226.9,376115,371710,0.33
+"307",3232.9,379228,375539,0.15
+"308",3238.9,382590,378831,0.39
+"309",3244.9,385814,382627,0.35
+"310",3250.9,389631,385590,0.65
+"311",3259.9,395516,391478,0.52
+"312",3268.9,400502,397554,0.55
+"313",3274.9,403775,401044,0.38
+"314",3277.9,405354,402666,0.24
+"315",3284.9,408947,406607,0.02
+"316",3290.9,411959,409860,-0.08
+"317",3296.9,415188,412938,0.02
+"318",3299.9,416824,414603,-0.06
Modified: trunk/octave-forge/extra/lssa/lombcoeff.m
===================================================================
--- trunk/octave-forge/extra/lssa/lombcoeff.m 2012-07-10 14:20:31 UTC (rev 10741)
+++ trunk/octave-forge/extra/lssa/lombcoeff.m 2012-07-12 05:03:35 UTC (rev 10742)
@@ -22,19 +22,19 @@
## @seealso{lombnormcoeff}
## @end deftypefn
-%!shared t, x, o, maxfreq
+%!test
+%! shared t, x, o, maxfreq
%! maxfreq = 4 / ( 2 * pi );
%! t = linspace(0,8); x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)
%! - 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t)
%! + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
-%! assert( lombcoeff(t,x,o(1)) = );
-%! assert( lombcoeff(t,x,o(2)) = );
-%! assert( lombcoeff(t,x,o(3)) = );
+%!assert( lombcoeff(t,x,o(1)),10788.9848389923,5e-10 );
+%!assert( lombcoeff(t,x,o(2)),12352.6413413457,5e-10 );
+%!assert( lombcoeff(t,x,o(3)),13673.4098969780,5e-10 );
function coeff = lombcoeff(T, X, o)
theta = atan2(sum(sin(2 .* o .* T )), sum(cos(2.*o.*T)))/ (2 * o );
- coeff = ( sum(X .* cos(o .* T - tau))**2)/(sum(cos(o.*T-tau).**2)) +
- ( sum(X .* sin(o .* T - tau))**2)/(sum(sin(o.*T-tau).**2));
+ coeff = ( sum(X .* cos(o .* T - theta))**2)/(sum(cos(o.*T-theta).**2)) + ( sum(X .* sin(o .* T - theta))**2)/(sum(sin(o.*T-theta).**2));
end function
Modified: trunk/octave-forge/extra/lssa/lscomplex.m
===================================================================
--- trunk/octave-forge/extra/lssa/lscomplex.m 2012-07-10 14:20:31 UTC (rev 10741)
+++ trunk/octave-forge/extra/lssa/lscomplex.m 2012-07-12 05:03:35 UTC (rev 10742)
@@ -26,7 +26,7 @@
%!test
%! shared t, x, o, maxfreq
%! maxfreq = 4 / ( 2 * pi ); t = [0:0.008:8]; x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)- 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t) + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
-%! assert( lscomplex(t,x,maxfreq,2,2), [-0.40075 - 2.36687i , 1.22666 - 2.24390i , 1.93643 - 1.51554i , 2.12505 - 0.95410i ] );
+%! assert( lscomplex(t,x,maxfreq,2,2), [-0.400754376933531 - 2.366871097665244i, 1.226663545950135 - 2.243899314661490i, 1.936433327880238 - 1.515538553198501i, 2.125045509991203 - 0.954100898917708i ], 6e-14 );
Modified: trunk/octave-forge/extra/lssa/lsreal.m
===================================================================
--- trunk/octave-forge/extra/lssa/lsreal.m 2012-07-10 14:20:31 UTC (rev 10741)
+++ trunk/octave-forge/extra/lssa/lsreal.m 2012-07-12 05:03:35 UTC (rev 10742)
@@ -25,9 +25,10 @@
%!shared t, x, o, maxfreq
%! maxfreq = 4 / ( 2 * pi );
-%! t = linspace(0,8); x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)
-%! - 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t)
-%! + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
+%! t = linspace(0,8);
+%! x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t) - 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t) + cos((1/4)*maxfreq.*t));
+%! assert(lsreal(t,x,maxfreq,2,2),[-1.68275915310663 + 4.70126183846743i, 1.93821553170889 + 4.95660209883437i, 4.38145452686697 + 2.14403733658600i, 5.27425332281147 - 0.73933440226597i],6e-14)
+%! #In the assert here, I've got an error bound large enough to catch individual system errors which would present no real issue.
function transform = lsreal( t, x, omegamax, ncoeff, noctave)
## the R function runs the following command:
@@ -38,7 +39,7 @@
## where rp = complex(noctave*ncoeff) so ... I can just store that as noctave*ncoeff and have no
## problems, I guess.
## Possibly throw an error if ncoeff <= 0.
- k = n = min(length(x),length(t)); ## THIS IS VECTOR-ONLY. I'd need to add another bit of code to
+ k = n = length(t); ## THIS IS VECTOR-ONLY. I'd need to add another bit of code to
## make it array-safe, and that's not knowing right now what else will be necessary.
transform = zeros(1,(noctave * ncoeff)); ## In the C code, this is rendered as a Complex, but Octave doesn't care.
od = 2 ^ ( - 1 / ncoeff ); ## this will cause a crash if ncoeff=0; prefer error & quit?
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <lm...@us...> - 2012-07-10 14:20:42
|
Revision: 10741
http://octave.svn.sourceforge.net/octave/?rev=10741&view=rev
Author: lmarkov
Date: 2012-07-10 14:20:31 +0000 (Tue, 10 Jul 2012)
Log Message:
-----------
A fuzzy-logic-toolkit/NEWS
M fuzzy-logic-toolkit/DESCRIPTION
M fuzzy-logic-toolkit/ChangeLog
M fuzzy-logic-toolkit/INDEX
Modified Paths:
--------------
trunk/octave-forge/main/fuzzy-logic-toolkit/ChangeLog
trunk/octave-forge/main/fuzzy-logic-toolkit/DESCRIPTION
trunk/octave-forge/main/fuzzy-logic-toolkit/INDEX
Added Paths:
-----------
trunk/octave-forge/main/fuzzy-logic-toolkit/NEWS
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/ChangeLog
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/ChangeLog 2012-07-10 14:20:04 UTC (rev 10740)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/ChangeLog 2012-07-10 14:20:31 UTC (rev 10741)
@@ -1,7 +1,34 @@
-2011-11-12 Version 0.3.0
+2012-07-10 L. Markowsky <lm...@us...>
+ * Version 0.4.0 released.
* ChangeLog: Updated file.
* DESCRIPTION: Updated file.
+ * INDEX: Updated file.
+ * NEWS: New file.
+ * inst/fcm.m: New file. Addition of the Fuzzy C-Means clustering algorithm
+ to the toolkit.
+ * inst/fcm_demo_1.m: New file. Addition of demo script.
+ * inst/fcm_demo_2.m: New file. Addition of demo script.
+ * inst/gustafson_kessel.m: New file. Addition of the Gustafson-Kessel
+ clustering algorithm to the toolkit.
+ * inst/gustafson_kessel_demo_1.m: New file. Addition of demo script.
+ * inst/gustafson_kessel_demo_2.m: New file. Addition of demo script.
+ * inst/partition_coeff.m: New file. Addition of a measure of cluster validity.
+ * inst/partition_entropy.m: New file. Addition of a measure of cluster validity.
+ * inst/xie_beni_index.m: New file. Addition of a measure of cluster validity.
+ * inst/private/fcm_compute_convergence_criterion.m: New file.
+ * inst/private/fcm_compute_objective_fcn.m: New file.
+ * inst/private/fcm_init_prototype.m: New file.
+ * inst/private/fcm_update_cluster_centers.m: New file.
+ * inst/private/fcm_update_membership_fcn.m: New file.
+ * inst/private/square_dist_matrix.m: New file.
+ * New demos tested under: Fedora 16/Octave 3.4.3
+
+2011-11-12 L. Markowsky <lm...@us...>
+
+ * Version 0.3.0 released.
+ * ChangeLog: Updated file.
+ * DESCRIPTION: Updated file.
* inst/*.m: Many trivial changes to comments and spacing in parameter lists.
* inst/addrule.m: Edited comment to describe use with hedges.
* inst/algebraic_product.m: New file.
@@ -56,8 +83,9 @@
* Demos tested under: Fedora 15/Octave 3.4.2
* Demos tested under: Windows 7/Octave 3.2.4
-2011-09-01 Version 0.2.4
+2011-09-01 L. Markowsky <lm...@us...>
+ * Version 0.2.4 released.
* ChangeLog: Updated file.
* DESCRIPTION: Updated file.
* INDEX: Updated file.
@@ -102,8 +130,9 @@
* Demos tested under: Fedora 15/Octave 3.2.4
* Demos tested under: Windows 7/Octave 3.2.4
-2011-07-19 Version 0.2.3
+2011-07-19 L. Markowsky <lm...@us...>
+ * Version 0.2.3 released.
* ChangeLog: Updated file.
* DESCRIPTION: Updated file.
* INDEX: Updated file.
@@ -123,8 +152,9 @@
* Demos tested under: Fedora 15/Octave 3.4.0
* Demos tested under: Fedora 15/Octave 3.2.4
-2011-06-21 Version 0.2.2
+2011-06-21 L. Markowsky <lm...@us...>
+ * Version 0.2.2 released.
* ChangeLog: New file.
* DESCRIPTION: Updated file.
* inst/addmf.m: Modified to workaround a bug in Octave 3.4.0.
@@ -140,23 +170,26 @@
* Demos tested under: Fedora 15/Octave 3.4.0
* Demos tested under: Fedora 15/Octave 3.2.4
-2011-06-08 Version 0.2.1
+2011-06-08 L. Markowsky <lm...@us...>
+ * Version 0.2.1 released.
* Initial release on Octave-Forge.
* Merged membership function demos into related function files.
* Created documentation for Octave-Forge website.
* DESCRIPTION: Updated file.
* Demos tested under: Fedora 13/Octave 3.2.4
-2011-05-25 Version 0.2
+2011-05-25 L. Markowsky <lm...@us...>
+ * Version 0.2 released.
* Moved tests/demos/* and tests/fis/* to inst/*.
* Changed indentation and spacing to conform to Octave style.
* Converted comments to texinfo.
* DESCRIPTION: Update file.
* Demos tested under: Fedora 13/Octave 3.2.4
-2011-04-19 Version 0.1
+2011-04-19 L. Markowsky <lm...@us...>
+ * Version 0.1 released.
* Initial release on SourceForge.
* Demos tested under: Fedora 13/Octave 3.2.4
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/DESCRIPTION
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/DESCRIPTION 2012-07-10 14:20:04 UTC (rev 10740)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/DESCRIPTION 2012-07-10 14:20:31 UTC (rev 10741)
@@ -1,6 +1,6 @@
Name: fuzzy-logic-toolkit
-Version: 0.3.0
-Date: 2011-11-13
+Version: 0.4.0
+Date: 2012-07-10
Author: L. Markowsky <lm...@us...>
Maintainer: L. Markowsky <lm...@us...>
Title: Octave Fuzzy Logic Toolkit
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/INDEX
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/INDEX 2012-07-10 14:20:04 UTC (rev 10740)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/INDEX 2012-07-10 14:20:31 UTC (rev 10741)
@@ -61,3 +61,13 @@
mamdani_tip_demo
sugeno_tip_demo
+Fuzzy Clustering Functions and Demos
+ fcm
+ fcm_demo_1
+ fcm_demo_2
+ gustafson_kessel
+ gustafson_kessel_demo_1
+ gustafson_kessel_demo_2
+ partition_coeff
+ partition_entropy
+ xie_beni_index
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/NEWS
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/NEWS (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/NEWS 2012-07-10 14:20:31 UTC (rev 10741)
@@ -0,0 +1,15 @@
+Summary of important user-visible changes for fuzzy-logic-toolkit 0.4.0:
+------------------------------------------------------------------------
+
+ ** The following functions are new:
+ fcm
+ gustafson_kessel
+ partition_coeff
+ partition_entropy
+ xie_beni_index
+
+ ** The following demo scripts are new:
+ fcm_demo_1
+ fcm_demo_2
+ gustafson_kessel_demo_1
+ gustafson_kessel_demo_2
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <lm...@us...> - 2012-07-10 14:20:13
|
Revision: 10740
http://octave.svn.sourceforge.net/octave/?rev=10740&view=rev
Author: lmarkov
Date: 2012-07-10 14:20:04 +0000 (Tue, 10 Jul 2012)
Log Message:
-----------
M inst/gustafson_kessel_demo_1.m
M inst/gustafson_kessel_demo_2.m
M inst/xie_beni_index.m
M inst/gustafson_kessel.m
M inst/fcm_demo_1.m
M inst/fcm_demo_2.m
M inst/fcm.m
Modified Paths:
--------------
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/xie_beni_index.m
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -24,7 +24,9 @@
## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} fcm (@var{input_data}, @var{num_clusters}, @var{options})
## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} fcm (@var{input_data}, @var{num_clusters}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
##
-## Return the soft partition of a set of unlabeled data points.
+## Using the Fuzzy C-Means algorithm, calculate and return the soft partition
+## of a set of unlabeled data points.
+##
## Also, if @var{display_intermediate_results} is true, display intermediate
## results after each iteration.
##
@@ -112,7 +114,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm
## Directory: fuzzy-logic-toolkit/inst/
## Filename: fcm.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
function [cluster_centers, soft_partition, obj_fcn_history] = ...
fcm (input_data, num_clusters, options = [2.0, 100, 1e-5, 1])
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -18,8 +18,8 @@
## -*- texinfo -*-
## @deftypefn {Script File} {} fcm_demo_1
-## Use the Fuzzy C-Means algorithm to classify unlabeled data points and
-## evaluate the quality of the resulting clusters.
+## Use the Fuzzy C-Means algorithm to classify a small set of unlabeled
+## data points and evaluate the quality of the resulting clusters.
##
## This demo:
## @itemize @minus
@@ -45,7 +45,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
## Directory: fuzzy-logic-toolkit/inst/
## Filename: fcm_demo_1.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
##------------------------------------------------------------------------------
## Note: The input_data is taken from Chapter 13, Example 17 in
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -18,8 +18,8 @@
## -*- texinfo -*-
## @deftypefn {Script File} {} fcm_demo_2
-## Use the Fuzzy C-Means algorithm to classify unlabeled data points and
-## evaluate the quality of the resulting clusters.
+## Use the Fuzzy C-Means algorithm to classify three-dimensional unlabeled
+## data points and evaluate the quality of the resulting clusters.
##
## The demo:
## @itemize @minus
@@ -45,7 +45,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
## Directory: fuzzy-logic-toolkit/inst/
## Filename: fcm_demo_2.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
##------------------------------------------------------------------------------
## Note: The input_data was selected to form three areas of different shapes.
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -26,7 +26,9 @@
## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, @var{options})
## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
##
-## Return the soft partition of a set of unlabeled data points.
+## Using the Gustafson-Kessel algorithm, calculate and return the soft partition
+## of a set of unlabeled data points.
+##
## Also, if @var{display_intermediate_results} is true, display intermediate
## results after each iteration.
##
@@ -117,7 +119,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering gustafson_kessel
## Directory: fuzzy-logic-toolkit/inst/
## Filename: gustafson_kessel.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
function [cluster_centers, soft_partition, obj_fcn_history] = ...
gustafson_kessel (input_data, num_clusters, cluster_volume = [], options = [2.0, 100, 1e-5, 1])
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -18,8 +18,8 @@
## -*- texinfo -*-
## @deftypefn {Script File} {} gustafson_kessel_demo_1
-## Use the Gustafson-Kessel algorithm to classify unlabeled data points and
-## evaluate the quality of the resulting clusters.
+## Use the Gustafson-Kessel algorithm to classify a small set of unlabeled
+## data points and evaluate the quality of the resulting clusters.
##
## The demo:
## @itemize @minus
@@ -45,7 +45,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering gustafson_kessel demo
## Directory: fuzzy-logic-toolkit/inst/
## Filename: gustafson_kessel_demo_1.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
##------------------------------------------------------------------------------
## Note: The input_data is taken from Chapter 13, Example 17 in
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -18,8 +18,8 @@
## -*- texinfo -*-
## @deftypefn {Script File} {} gustafson_kessel_demo_2
-## Use the Gustafson-Kessel algorithm to classify unlabeled data points and
-## evaluate the quality of the resulting clusters.
+## Use the Gustafson-Kessel algorithm to classify three-dimensional unlabeled
+## data points and evaluate the quality of the resulting clusters.
##
## The demo:
## @itemize @minus
@@ -45,7 +45,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
## Directory: fuzzy-logic-toolkit/inst/
## Filename: gustafson_kessel_demo_2.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
##------------------------------------------------------------------------------
## Note: The input_data was selected to form three areas of different shapes.
Modified: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/xie_beni_index.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/xie_beni_index.m 2012-07-09 13:17:11 UTC (rev 10739)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/xie_beni_index.m 2012-07-10 14:20:04 UTC (rev 10740)
@@ -19,7 +19,7 @@
## -*- texinfo -*-
## @deftypefn {Function File} {@var{vxb} =} xie_beni_index (@var{input_data}, @var{cluster_centers}, @var{soft_partition})
##
-## Return the partition coefficient for a given soft partition.
+## Return the Xie-Beni validity index for a given soft partition.
##
## The arguments to xie_beni_index are:
## @itemize @w
@@ -49,7 +49,7 @@
## Keywords: fuzzy-logic-toolkit fuzzy xie beni index cluster validity
## Directory: fuzzy-logic-toolkit/inst/
## Filename: xie_beni_index.m
-## Last-Modified: 8 July 2012
+## Last-Modified: 10 July 2012
function vxb = xie_beni_index (input_data, cluster_centers, soft_partition)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <lm...@us...> - 2012-07-09 13:17:26
|
Revision: 10739
http://octave.svn.sourceforge.net/octave/?rev=10739&view=rev
Author: lmarkov
Date: 2012-07-09 13:17:11 +0000 (Mon, 09 Jul 2012)
Log Message:
-----------
A inst/gustafson_kessel_demo_1.m
A inst/gustafson_kessel_demo_2.m
A inst/partition_coeff.m
A inst/fcm.m
A inst/xie_beni_index.m
A inst/gustafson_kessel.m
A inst/private/square_distance_matrix.m
A inst/private/fcm_init_prototype.m
A inst/private/fcm_update_cluster_centers.m
A inst/private/fcm_compute_convergence_criterion.m
A inst/private/fcm_update_membership_fcn.m
A inst/private/fcm_compute_objective_fcn.m
A inst/fcm_demo_1.m
A inst/partition_entropy.m
A inst/fcm_demo_2.m
Added Paths:
-----------
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_coeff.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_entropy.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_convergence_criterion.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_objective_fcn.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_init_prototype.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_cluster_centers.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_membership_fcn.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/square_distance_matrix.m
trunk/octave-forge/main/fuzzy-logic-toolkit/inst/xie_beni_index.m
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,204 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{cluster_centers} =} fcm (@var{input_data}, @var{num_clusters})
+## @deftypefnx {Function File} {@var{cluster_centers} =} fcm (@var{input_data}, @var{num_clusters}, @var{options})
+## @deftypefnx {Function File} {@var{cluster_centers} =} fcm (@var{input_data}, @var{num_clusters}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} fcm (@var{input_data}, @var{num_clusters})
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} fcm (@var{input_data}, @var{num_clusters}, @var{options})
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} fcm (@var{input_data}, @var{num_clusters}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
+##
+## Return the soft partition of a set of unlabeled data points.
+## Also, if @var{display_intermediate_results} is true, display intermediate
+## results after each iteration.
+##
+## The required arguments to fcm are:
+## @itemize @w
+## @item
+## @var{input_data} - a matrix of input data points; each row corresponds to one point
+## @item
+## @var{num_clusters} - the number of clusters to form
+## @end itemize
+##
+## The optional arguments to fcm are:
+## @itemize @w
+## @item
+## @var{m} - the parameter (exponent) in the objective function; default = 2.0
+## @item
+## @var{max_iterations} - the maximum number of iterations before stopping; default = 100
+## @item
+## @var{epsilon} - the stopping criteria; default = 1e-5
+## @item
+## @var{display_intermediate_results} - if 1, display results after each iteration, and if 0, do not; default = 1
+## @end itemize
+##
+## The default values are used if any of the optional arguments are missing or
+## evaluate to NaN.
+##
+## The return values are:
+## @itemize @w
+## @item
+## @var{cluster_centers} - a matrix of the cluster centers; each row corresponds to one point
+## @item
+## @var{soft_partition} - a constrained soft partition matrix
+## @item
+## @var{obj_fcn_history} - the values of the objective function after each iteration
+## @end itemize
+##
+## Three important matrices used in the calculation are X (the input points
+## to be clustered), V (the cluster centers), and Mu (the membership of each
+## data point in each cluster). Each row of X and V denotes a single point,
+## and Mu(i, j) denotes the membership degree of input point X(j, :) in the
+## cluster having center V(i, :).
+##
+## X is identical to the required argument @var{input_data}; V is identical
+## to the output @var{cluster_centers}; and Mu is identical to the output
+## @var{soft_partition}.
+##
+## If n denotes the number of input points and k denotes the number of
+## clusters to be formed, then X, V, and Mu have the dimensions:
+##
+## @example
+## @group
+## 1 2 ... #features
+## 1 [ ]
+## X = input_data = 2 [ ]
+## ... [ ]
+## n [ ]
+## @end group
+## @end example
+##
+## @example
+## @group
+## 1 2 ... #features
+## 1 [ ]
+## V = cluster_centers = 2 [ ]
+## ... [ ]
+## k [ ]
+## @end group
+## @end example
+##
+## @example
+## @group
+## 1 2 ... n
+## 1 [ ]
+## Mu = soft_partition = 2 [ ]
+## ... [ ]
+## k [ ]
+## @end group
+## @end example
+##
+## @seealso{fcm_demo_1, fcm_demo_2, gustafson_kessel, gustafson_kessel_demo_1, gustafson_kessel_demo_2, partition_coeff, partition_entropy, xie_beni_index}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: fcm.m
+## Last-Modified: 8 July 2012
+
+function [cluster_centers, soft_partition, obj_fcn_history] = ...
+ fcm (input_data, num_clusters, options = [2.0, 100, 1e-5, 1])
+
+ ## If fcm was called with an incorrect number of arguments, or the
+ ## arguments do not have the correct type, print an error message and halt.
+
+ if ((nargin != 2) && (nargin != 3))
+ puts ("Type 'help fcm' for more information.\n");
+ error ("fcm requires 2 or 3 arguments\n");
+ elseif (!is_real_matrix (input_data))
+ puts ("Type 'help fcm' for more information.\n");
+ error ("fcm's first argument must be matrix of real numbers\n");
+ elseif (!(is_int (num_clusters) && (num_clusters > 1)))
+ puts ("Type 'help fcm' for more information.\n");
+ error ("fcm's second argument must be an integer greater than 1\n");
+ elseif (!(isreal (options) && isvector (options)))
+ puts ("Type 'help fcm' for more information.\n");
+ error ("fcm's third (optional) argument must be a vector of real numbers\n");
+ endif
+
+ ## Assign options to the more readable variable names: m, max_iterations,
+ ## epsilon, and display_intermediate_results. If options are missing or
+ ## NaN (not a number), use the default values.
+
+ default_options = [2.0, 100, 1e-5, 1];
+
+ for i = 1 : 4
+ if ((length (options) < i) || isna (options(i)) || isnan (options(i)))
+ options(i) = default_options(i);
+ endif
+ endfor
+
+ m = options(1);
+ max_iterations = options(2);
+ epsilon = options(3);
+ display_intermediate_results = options(4);
+
+ ## Call a private function to compute the output.
+
+ [cluster_centers, soft_partition, obj_fcn_history] = ...
+ fcm_private (input_data, num_clusters, m, max_iterations, epsilon,
+ display_intermediate_results);
+endfunction
+
+##------------------------------------------------------------------------------
+## Note: This function (fcm_private) is an implementation of Figure 13.4 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 380 (International Edition)
+## and Algorithm 4.1 in Fuzzy and Neural Control, by Robert Babuska,
+## November 2009, p. 63.
+##------------------------------------------------------------------------------
+
+function [V, Mu, obj_fcn_history] = ...
+ fcm_private (X, k, m, max_iterations, epsilon, display_intermediate_results)
+
+ ## Initialize the prototype and the calculation.
+ V = fcm_init_prototype (X, k);
+ obj_fcn_history = zeros (max_iterations);
+ convergence_criterion = epsilon + 1;
+ iteration = 0;
+
+ ## Calculate a few numbers here to reduce redundant computation.
+ k = rows (V);
+ n = rows (X);
+ sqr_dist = square_distance_matrix (X, V);
+
+ ## Loop until the objective function is within tolerance or the maximum
+ ## number of iterations has been reached.
+ while (convergence_criterion > epsilon && ++iteration <= max_iterations)
+ V_previous = V;
+ Mu = fcm_update_membership_fcn (V, X, m, k, n, sqr_dist);
+ Mu_m = Mu .^ m;
+ V = fcm_update_cluster_centers (Mu_m, X, k);
+ sqr_dist = square_distance_matrix (X, V);
+ obj_fcn_history(iteration) = fcm_compute_objective_fcn (Mu_m, sqr_dist);
+ if (display_intermediate_results)
+ printf ("Iteration count = %d, Objective fcn = %8.6f\n", ...
+ iteration, obj_fcn_history(iteration));
+ endif
+ convergence_criterion = fcm_compute_convergence_criterion (V, V_previous);
+ endwhile
+
+ ## Remove extraneous entries from the tail of the objective function history.
+ if (convergence_criterion <= epsilon)
+ obj_fcn_history = obj_fcn_history(1 : iteration);
+ endif
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_1.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,90 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Script File} {} fcm_demo_1
+## Use the Fuzzy C-Means algorithm to classify unlabeled data points and
+## evaluate the quality of the resulting clusters.
+##
+## This demo:
+## @itemize @minus
+## @item
+## classifies a small set of unlabeled data points using the Fuzzy C-Means
+## algorithm into two fuzzy clusters
+## @item
+## calculates and prints (on standard output) three validity measures:
+## the partition coefficient, the partition entropy, and the Xie-Beni
+## validity index
+## @item
+## plots the input points together with the cluster centers
+## @end itemize
+##
+## For a description of the data structures used in the script, see the
+## documentation for fcm.
+##
+## @seealso{fcm, fcm_demo_2, partition_coeff, partition_entropy, xie_beni_index, gustafson_kessel, gustafson_kessel_demo_1, gustafson_kessel_demo_2}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: fcm_demo_1.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: The input_data is taken from Chapter 13, Example 17 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 381 (International Edition).
+##------------------------------------------------------------------------------
+
+## Use fcm to classify the data in matrix x.
+input_data = [2 12; 4 9; 7 13; 11 5; 12 7; 14 4]
+number_of_clusters = 2
+[cluster_centers, soft_partition, obj_fcn_history] = ...
+ fcm (input_data, number_of_clusters)
+
+## Calculate and print the three validity measures.
+printf ("Partition Coefficient: %f\n", partition_coeff (soft_partition));
+printf ("Partition Entropy (with a = 2): %f\n", partition_entropy (soft_partition, 2));
+printf ("Xie-Beni Index: %f\n\n", xie_beni_index (input_data, cluster_centers, soft_partition));
+
+## Plot the data points as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'FCM Demo 1');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 2), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 2), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 2');
+grid
+hold
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/fcm_demo_2.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,121 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Script File} {} fcm_demo_2
+## Use the Fuzzy C-Means algorithm to classify unlabeled data points and
+## evaluate the quality of the resulting clusters.
+##
+## The demo:
+## @itemize @minus
+## @item
+## classifies three-dimensional unlabeled data points using the Fuzzy
+## C-Means algorithm into three fuzzy clusters
+## @item
+## calculates and prints (on standard output) three validity measures:
+## the partition coefficient, the partition entropy, and the Xie-Beni
+## validity index
+## @item
+## plots the input points together with the cluster centers
+## @end itemize
+##
+## For a description of the data structures used in the script, see the
+## documentation for fcm.
+##
+## @seealso{fcm, fcm_demo_1, partition_coeff, partition_entropy, xie_beni_index, gustafson_kessel, gustafson_kessel_demo_1, gustafson_kessel_demo_2}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: fcm_demo_2.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: The input_data was selected to form three areas of different shapes.
+##------------------------------------------------------------------------------
+
+## Use fcm to classify the data in matrix x.
+input_data = [1 11 5; 1 12 6; 1 13 5; 2 11 7; 2 12 6; 2 13 7; 3 11 6; 3 12 5;
+ 3 13 7; 1 1 10; 1 3 9; 2 2 11; 3 1 9; 3 3 10; 3 5 11; 4 4 9;
+ 4 6 8; 5 5 8; 5 7 9; 6 6 10; 9 10 12; 9 12 13; 9 13 14; 10 9 13;
+ 10 13 12; 11 10 14; 11 12 13; 12 6 12; 12 7 15; 12 9 15; 14 6 14;
+ 14 8 13]
+number_of_clusters = 3
+[cluster_centers, soft_partition, obj_fcn_history] = ...
+ fcm (input_data, number_of_clusters, [NaN NaN NaN 0])
+
+## Calculate and print the three validity measures.
+printf ("Partition Coefficient: %f\n", partition_coeff (soft_partition));
+printf ("Partition Entropy (with a = 2): %f\n", partition_entropy (soft_partition, 2));
+printf ("Xie-Beni Index: %f\n\n", xie_beni_index (input_data, cluster_centers, soft_partition));
+
+## Plot the data points in two dimensions (using features 1 and 2)
+## as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'FCM Demo 2');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 2), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers in two dimensions (using features 1 and 2)
+## as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 2), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 2');
+grid
+hold
+
+## Plot the data points in two dimensions (using features 1 and 3)
+## as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'FCM Demo 2');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 3), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers in two dimensions (using features 1 and 3)
+## as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 3), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 3');
+grid
+hold
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,272 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{cluster_centers} =} gustafson_kessel (@var{input_data}, @var{num_clusters})
+## @deftypefnx {Function File} {@var{cluster_centers} =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume})
+## @deftypefnx {Function File} {@var{cluster_centers} =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, @var{options})
+## @deftypefnx {Function File} {@var{cluster_centers} =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters})
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume})
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, @var{options})
+## @deftypefnx {Function File} {[@var{cluster_centers}, @var{soft_partition}, @var{obj_fcn_history}] =} gustafson_kessel (@var{input_data}, @var{num_clusters}, @var{cluster_volume}, [@var{m}, @var{max_iterations}, @var{epsilon}, @var{display_intermediate_results}])
+##
+## Return the soft partition of a set of unlabeled data points.
+## Also, if @var{display_intermediate_results} is true, display intermediate
+## results after each iteration.
+##
+## The required arguments to gustafson_kessel are:
+## @itemize @w
+## @item
+## @var{input_data} - a matrix of input data points; each row corresponds to one point
+## @item
+## @var{num_clusters} - the number of clusters to form
+## @end itemize
+##
+## The third (optional) argument to gustafson_kessel is a vector of cluster volumes.
+## If omitted, a vector of 1's will be used as the default.
+##
+## The fourth (optional) argument to gustafson_kessel is a vector consisting of:
+## @itemize @w
+## @item
+## @var{m} - the parameter (exponent) in the objective function; default = 2.0
+## @item
+## @var{max_iterations} - the maximum number of iterations before stopping; default = 100
+## @item
+## @var{epsilon} - the stopping criteria; default = 1e-5
+## @item
+## @var{display_intermediate_results} - if 1, display results after each iteration, and if 0, do not; default = 1
+## @end itemize
+##
+## The default values are used if any of the four elements of the vector are missing or
+## evaluate to NaN.
+##
+## The return values are:
+## @itemize @w
+## @item
+## @var{cluster_centers} - a matrix of the cluster centers; each row corresponds to one point
+## @item
+## @var{soft_partition} - a constrained soft partition matrix
+## @item
+## @var{obj_fcn_history} - the values of the objective function after each iteration
+## @end itemize
+##
+## Three important matrices used in the calculation are X (the input points
+## to be clustered), V (the cluster centers), and Mu (the membership of each
+## data point in each cluster). Each row of X and V denotes a single point,
+## and Mu(i, j) denotes the membership degree of input point X(j, :) in the
+## cluster having center V(i, :).
+##
+## X is identical to the required argument @var{input_data}; V is identical
+## to the output @var{cluster_centers}; and Mu is identical to the output
+## @var{soft_partition}.
+##
+## If n denotes the number of input points and k denotes the number of
+## clusters to be formed, then X, V, and Mu have the dimensions:
+##
+## @example
+## @group
+## 1 2 ... #features
+## 1 [ ]
+## X = input_data = 2 [ ]
+## ... [ ]
+## n [ ]
+## @end group
+## @end example
+##
+## @example
+## @group
+## 1 2 ... #features
+## 1 [ ]
+## V = cluster_centers = 2 [ ]
+## ... [ ]
+## k [ ]
+## @end group
+## @end example
+##
+## @example
+## @group
+## 1 2 ... n
+## 1 [ ]
+## Mu = soft_partition = 2 [ ]
+## ... [ ]
+## k [ ]
+## @end group
+## @end example
+##
+## @seealso{fcm, fcm_demo_1, fcm_demo_2, gustafson_kessel_demo_1, gustafson_kessel_demo_2, partition_coeff, partition_entropy, xie_beni_index}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering gustafson_kessel
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: gustafson_kessel.m
+## Last-Modified: 8 July 2012
+
+function [cluster_centers, soft_partition, obj_fcn_history] = ...
+ gustafson_kessel (input_data, num_clusters, cluster_volume = [], options = [2.0, 100, 1e-5, 1])
+
+ ## If gustafson_kessel was called with an incorrect number of arguments, or the
+ ## arguments do not have the correct type, print an error message and halt.
+
+ if ((nargin < 2) || (nargin > 4))
+ puts ("Type 'help gustafson_kessel' for more information.\n");
+ error ("gustafson_kessel requires 2, 3, or 4 arguments\n");
+ elseif (!is_real_matrix (input_data))
+ puts ("Type 'help gustafson_kessel' for more information.\n");
+ error ("gustafson_kessel's first argument must be matrix of real numbers\n");
+ elseif (!(is_int (num_clusters) && (num_clusters > 1)))
+ puts ("Type 'help gustafson_kessel' for more information.\n");
+ error ("gustafson_kessel's second argument must be an integer greater than 1\n");
+ elseif (!(isequal (cluster_volume, []) || ...
+ (isreal (cluster_volume) && isvector (cluster_volume))))
+ puts ("Type 'help gustafson_kessel' for more information.\n");
+ error ("gustafson_kessel's third (optional) argument must be a vector of real numbers\n");
+ elseif (!(isreal (options) && isvector (options)))
+ puts ("Type 'help gustafson_kessel' for more information.\n");
+ error ("gustafson_kessel's fourth (optional) argument must be a vector of real numbers\n");
+ endif
+
+ ## If the cluster volume matrix was not entered, create a default value
+ ## (a vector of 1's).
+
+ if (isequal (cluster_volume, []))
+ cluster_volume = ones (1, num_clusters);
+ endif
+
+ ## Assign options to the more readable variable names: m, max_iterations,
+ ## epsilon, and display_intermediate_results. If options are missing or
+ ## NaN (not a number), use the default values.
+
+ default_options = [2.0, 100, 1e-5, 1];
+
+ for i = 1 : 4
+ if ((length (options) < i) || isna (options(i)) || isnan (options(i)))
+ options(i) = default_options(i);
+ endif
+ endfor
+
+ m = options(1);
+ max_iterations = options(2);
+ epsilon = options(3);
+ display_intermediate_results = options(4);
+
+ ## Call a private function to compute the output.
+
+ [cluster_centers, soft_partition, obj_fcn_history] = ...
+ gustafson_kessel_private (input_data, num_clusters, cluster_volume, m, ...
+ max_iterations, epsilon, display_intermediate_results);
+endfunction
+
+##------------------------------------------------------------------------------
+## Function: gustafson_kessel_private
+## Purpose: Classify unlabeled data points using the Gustafson-Kessel
+## algorithm.
+## Note: This function (gustafson_kessel_private) is an implementation of
+## Algorithm 4.2 in Fuzzy and Neural Control, by Robert Babuska,
+## November 2009, p. 69.
+##------------------------------------------------------------------------------
+
+function [V, Mu, obj_fcn_history] = ...
+ gustafson_kessel_private (X, k, cluster_volume, m, max_iterations, epsilon, display_intermediate_results)
+
+ ## Initialize the prototype and the calculation.
+ V = fcm_init_prototype (X, k);
+ obj_fcn_history = zeros (max_iterations);
+ convergence_criterion = epsilon + 1;
+ iteration = 0;
+
+ ## Calculate a few numbers here to reduce redundant computation.
+ k = rows (V);
+ n = rows (X);
+ sqr_dist = square_distance_matrix (X, V);
+
+ ## Loop until the objective function is within tolerance or the maximum
+ ## number of iterations has been reached.
+ while (convergence_criterion > epsilon && ++iteration <= max_iterations)
+ V_previous = V;
+ Mu = fcm_update_membership_fcn (V, X, m, k, n, sqr_dist);
+ Mu_m = Mu .^ m;
+ V = fcm_update_cluster_centers (Mu_m, X, k);
+ sqr_dist = gk_square_distance_matrix (X, V, Mu_m, cluster_volume);
+ obj_fcn_history(iteration) = fcm_compute_objective_fcn (Mu_m, sqr_dist);
+ if (display_intermediate_results)
+ printf ("Iteration count = %d, Objective fcn = %8.6f\n", ...
+ iteration, obj_fcn_history(iteration));
+ endif
+ convergence_criterion = fcm_compute_convergence_criterion (V, V_previous);
+ endwhile
+
+ ## Remove extraneous entries from the tail of the objective function history.
+ if (convergence_criterion <= epsilon)
+ obj_fcn_history = obj_fcn_history(1 : iteration);
+ endif
+
+endfunction
+
+##------------------------------------------------------------------------------
+## Function: gk_square_distance_matrix
+##------------------------------------------------------------------------------
+
+function sqr_dist = gk_square_distance_matrix (X, V, Mu_m, cluster_volume)
+
+ k = rows (V);
+ n = rows (X);
+ num_features = columns (X);
+ sqr_dist = zeros (k, n);
+
+ for i = 1 : k
+ Vi = V(i, :);
+ covariance_matrix = compute_covariance_matrix (X, V, Mu_m, i);
+
+ for j = 1 : n
+ Vi_to_Xj = X(j, :) - Vi;
+ A = cluster_volume(i) * ...
+ det (covariance_matrix) ^ (1.0 / num_features) * ...
+ inv (covariance_matrix);
+ sqr_dist(i, j) = sum (Vi_to_Xj .* (A * Vi_to_Xj')');
+
+ endfor
+
+ endfor
+
+endfunction
+
+##------------------------------------------------------------------------------
+## Function: compute_covariance_matrix
+##------------------------------------------------------------------------------
+
+function covariance_matrix = compute_covariance_matrix (X, V, Mu_m, i)
+
+ num_features = columns (V);
+ n = rows (X);
+ num = zeros (num_features);
+ denom = 0.0;
+ Vi = V(i, :);
+
+ for j = 1 : n
+ Vi_to_Xj = X(j, :) - Vi;
+ num += Mu_m(i, j) * Vi_to_Xj' * Vi_to_Xj;
+ denom += Mu_m(i, j);
+ endfor
+
+ covariance_matrix = num / denom;
+
+endfunction
+
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_1.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,90 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Script File} {} gustafson_kessel_demo_1
+## Use the Gustafson-Kessel algorithm to classify unlabeled data points and
+## evaluate the quality of the resulting clusters.
+##
+## The demo:
+## @itemize @minus
+## @item
+## classifies a small set of unlabeled data points using the Gustafson-Kessel
+## algorithm into two fuzzy clusters
+## @item
+## calculates and prints (on standard output) three validity measures:
+## the partition coefficient, the partition entropy, and the Xie-Beni
+## validity index
+## @item
+## plots the input points together with the cluster centers
+## @end itemize
+##
+## For a description of the data structures used in the script, see the
+## documentation for gustafson_kessel.
+##
+## @seealso{gustafson_kessel, gustafson_kessel_demo_2, partition_coeff, partition_entropy, xie_beni_index, fcm, fcm_demo_1, fcm_demo_2}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering gustafson_kessel demo
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: gustafson_kessel_demo_1.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: The input_data is taken from Chapter 13, Example 17 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 381 (International Edition).
+##------------------------------------------------------------------------------
+
+## Use gustafson_kessel to classify the data in matrix x.
+input_data = [2 12; 4 9; 7 13; 11 5; 12 7; 14 4]
+number_of_clusters = 2
+[cluster_centers, soft_partition, obj_fcn_history] = ...
+ gustafson_kessel (input_data, number_of_clusters)
+
+## Calculate and print the three validity measures.
+printf ("Partition Coefficient: %f\n", partition_coeff (soft_partition));
+printf ("Partition Entropy (with a = 2): %f\n", partition_entropy (soft_partition, 2));
+printf ("Xie-Beni Index: %f\n\n", xie_beni_index (input_data, cluster_centers, soft_partition));
+
+## Plot the data points as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'Gustafson-Kessel Demo 1');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 2), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 2), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 2');
+grid
+hold
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/gustafson_kessel_demo_2.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,120 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Script File} {} gustafson_kessel_demo_2
+## Use the Gustafson-Kessel algorithm to classify unlabeled data points and
+## evaluate the quality of the resulting clusters.
+##
+## The demo:
+## @itemize @minus
+## @item
+## classifies three-dimensional unlabeled data points using the
+## Gustafson-Kessel algorithm into three fuzzy clusters
+## @item
+## calculates and prints (on standard output) three validity measures:
+## the partition coefficient, the partition entropy, and the Xie-Beni
+## validity index
+## @item
+## plots the input points together with the cluster centers
+## @end itemize
+##
+## For a description of the data structures used in the script, see the
+## documentation for gustafson_kessel.
+##
+## @seealso{gustafson_kessel, gustafson_kessel_demo_1, partition_coeff, partition_entropy, xie_beni_index, fcm, fcm_demo_1, fcm_demo_2}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm demo
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: gustafson_kessel_demo_2.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: The input_data was selected to form three areas of different shapes.
+##------------------------------------------------------------------------------
+
+## Use gustafson_kessel to classify the data in matrix x.
+input_data = [1 11 5; 1 12 6; 1 13 5; 2 11 7; 2 12 6; 2 13 7; 3 11 6; 3 12 5;
+ 3 13 7; 1 1 10; 1 3 9; 2 2 11; 3 1 9; 3 3 10; 3 5 11; 4 4 9;
+ 4 6 8; 5 5 8; 5 7 9; 6 6 10; 9 10 12; 9 12 13; 9 13 14; 10 9 13;
+ 10 13 12; 11 10 14; 11 12 13; 12 6 12; 12 7 15; 12 9 15; 14 6 14;
+ 14 8 13]
+number_of_clusters = 3
+[cluster_centers, soft_partition, obj_fcn_history] = ...
+ gustafson_kessel (input_data, number_of_clusters, [1 1 1], [NaN NaN NaN 0])
+
+## Calculate and print the three validity measures.
+printf ("Partition Coefficient: %f\n", partition_coeff (soft_partition));
+printf ("Partition Entropy (with a = 2): %f\n", partition_entropy (soft_partition, 2));
+printf ("Xie-Beni Index: %f\n\n", xie_beni_index (input_data, cluster_centers, soft_partition));
+
+## Plot the data points in two dimensions (using features 1 and 2)
+## as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'Gustafson-Kessel Demo 2');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 2), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers in two dimensions (using features 1 and 2)
+## as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 2), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 2');
+grid
+
+## Plot the data points in two dimensions (using features 1 and 3)
+## as small blue x's.
+figure ('NumberTitle', 'off', 'Name', 'Gustafson-Kessel Demo 2');
+for i = 1 : rows (input_data)
+ plot (input_data(i, 1), input_data(i, 3), 'LineWidth', 2, 'marker', 'x', ...
+ 'color', 'b');
+ hold on;
+endfor
+
+## Plot the cluster centers in two dimensions (using features 1 and 3)
+## as larger red *'s.
+for i = 1 : number_of_clusters
+ plot (cluster_centers(i, 1), cluster_centers(i, 3), 'LineWidth', 4, ...
+ 'marker', '*', 'color', 'r');
+ hold on;
+endfor
+
+## Make the figure look a little better:
+## -- scale and label the axes
+## -- show gridlines
+xlim ([0 15]);
+ylim ([0 15]);
+xlabel ('Feature 1');
+ylabel ('Feature 3');
+grid
+hold
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_coeff.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_coeff.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_coeff.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,76 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{vpc} =} partition_coeff (@var{soft_partition})
+##
+## Return the partition coefficient for a given soft partition.
+##
+## The argument to partition_coeff is:
+## @itemize @w
+## @item
+## @var{soft_partition} - the membership degree of each input data point in each cluster
+## @end itemize
+##
+## The return value is:
+## @itemize @w
+## @item
+## @var{vpc} - the partition coefficient for the given soft partition
+## @end itemize
+##
+## For more information about the @var{soft_partition} matrix, please see the
+## documentation for function fcm.
+##
+## @seealso{fcm, fcm_demo_1, fcm_demo_2, gustafson_kessel, gustafson_kessel_demo_1, gustafson_kessel_demo_2, partition_entropy, xie_beni_index}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition coefficient cluster validity
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: partition_coeff.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: This function is an implementation of Equation 13.9 (corrected -- the
+## equation in the book omits the exponent 2) in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 384 (International Edition).
+##------------------------------------------------------------------------------
+
+function vpc = partition_coeff (soft_partition)
+
+ ## If partition_coeff was called with an incorrect number of arguments, or the
+ ## argument does not have the correct type, print an error message and halt.
+
+ if (nargin != 1)
+ puts ("Type 'help partition_coeff' for more information.\n");
+ error ("partition_coeff requires 1 argument\n");
+ elseif (!(is_real_matrix (soft_partition) &&
+ (min (min (soft_partition)) >= 0) &&
+ (max (max (soft_partition)) <= 1)))
+ puts ("Type 'help partition_coeff' for more information.\n");
+ error ("partition_coeff's argument must be a matrix of real numbers mu, with 0 <= mu <= 1\n");
+ endif
+
+ ## Compute and return the partition coefficient.
+
+ soft_part_sqr = soft_partition .* soft_partition;
+ vpc = (sum (sum (soft_part_sqr))) / columns (soft_partition);
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_entropy.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_entropy.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/partition_entropy.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,82 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{vpe} =} partition_entropy (@var{soft_partition}, @var{a})
+##
+## Return the partition entropy for a given soft partition.
+##
+## The arguments to partition_entropy are:
+## @itemize @w
+## @item
+## @var{soft_partition} - the membership degree of each input data point in each cluster
+## @item
+## @var{a} - the log base to use in the calculation; must be a real number a > 1
+## @end itemize
+##
+## The return value is:
+## @itemize @w
+## @item
+## @var{vpe} - the partition entropy for the given soft partition
+## @end itemize
+##
+## For more information about the @var{soft_partition} matrix, please see the
+## documentation for function fcm.
+##
+## @seealso{fcm, fcm_demo_1, fcm_demo_2, gustafson_kessel, gustafson_kessel_demo_1, gustafson_kessel_demo_2, partition_coeff, xie_beni_index}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition entropy cluster validity
+## Directory: fuzzy-logic-toolkit/inst/
+## Filename: partition_entropy.m
+## Last-Modified: 8 July 2012
+
+##------------------------------------------------------------------------------
+## Note: This function is an implementation of Equation 13.10 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 384 (International Edition).
+##------------------------------------------------------------------------------
+
+function vpe = partition_entropy (soft_partition, a)
+
+ ## If partition_entropy was called with an incorrect number of arguments, or the
+ ## argument does not have the correct type, print an error message and halt.
+
+ if (nargin != 2)
+ puts ("Type 'help partition_entropy' for more information.\n");
+ error ("partition_entropy requires 2 arguments\n");
+ elseif (!(is_real_matrix (soft_partition) &&
+ (min (min (soft_partition)) >= 0) &&
+ (max (max (soft_partition)) <= 1)))
+ puts ("Type 'help partition_entropy' for more information.\n");
+ error ("partition_entropy's first argument must be a matrix of real numbers mu, with 0 <= mu <= 1\n");
+ elseif (!(is_real (a) && a > 1))
+ puts ("Type 'help partition_entropy' for more information.\n");
+ error ("partition_entropy's second argument argument must be a real number a > 1\n");
+ endif
+
+ ## Compute and return the partition entropy.
+
+ n = columns (soft_partition);
+ Mu = soft_partition;
+ log_a_Mu = log (Mu) / log (a);
+ vpe = -(sum (sum (Mu .* log_a_Mu))) / n;
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_convergence_criterion.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_convergence_criterion.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_convergence_criterion.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,40 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{convergence_criterion} =} fcm_compute_convergence_criterion (@var{V}, @var{V_previous})
+##
+## Compute the sum of the changes in position (using the Euclidean
+## distance) of the cluster centers.
+##
+## @seealso{fcm, fcm_init_prototype, fcm_update_membership_fcn, fcm_update_cluster_centers, fcm_compute_objective_fcn}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm private
+## Directory: fuzzy-logic-toolkit/inst/private/
+## Filename: fcm_compute_convergence_criterion.m
+## Last-Modified: 7 July 2012
+
+function convergence_criterion = fcm_compute_convergence_criterion (V, V_previous)
+
+ V_delta = V - V_previous;
+ convergence_criterion = sum (sqrt (sum (V_delta .* V_delta)'));
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_objective_fcn.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_objective_fcn.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_compute_objective_fcn.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,44 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{obj_fcn} =} fcm_compute_objective_fcn (@var{Mu_m}, @var{sqr_dist})
+##
+## Compute the objective function for the current iteration.
+##
+## @seealso{fcm, fcm_init_prototype, fcm_update_membership_fcn, fcm_update_cluster_centers, fcm_compute_convergence_criterion}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm private
+## Directory: fuzzy-logic-toolkit/inst/private/
+## Filename: fcm_compute_objective_fcn.m
+## Last-Modified: 7 July 2012
+
+##------------------------------------------------------------------------------
+## Note: This function is an implementation of Equation 13.3 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 379 (International Edition).
+##------------------------------------------------------------------------------
+
+function obj_fcn = fcm_compute_objective_fcn (Mu_m, sqr_dist)
+
+ obj_fcn = sum (sum (Mu_m .* sqr_dist));
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_init_prototype.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_init_prototype.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_init_prototype.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,47 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{V} =} fcm_init_prototype (@var{X}, @var{k})
+##
+## Initialize k cluster centers to random locations in the ranges
+## given by the min/max values of each feature of the dataset.
+##
+## @seealso{fcm, fcm_update_membership_fcn, fcm_update_cluster_centers, fcm_compute_objective_fcn, fcm_compute_convergence_criterion}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm private
+## Directory: fuzzy-logic-toolkit/inst/private/
+## Filename: fcm_init_prototype.m
+## Last-Modified: 7 July 2012
+
+function V = fcm_init_prototype (X, k)
+
+ num_features = columns (X);
+ min_feature_value = min (X);
+ max_feature_value = max (X);
+ V = rand (k, num_features);
+
+ for i = 1 : num_features
+ V(:, i) = (max_feature_value(i) - min_feature_value(i)) * V(:, i) + ...
+ min_feature_value(i);
+ endfor
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_cluster_centers.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_cluster_centers.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_cluster_centers.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,54 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{V} =} fcm_update_cluster_centers (@var{Mu_m}, @var{X}, @var{k})
+##
+## Update the cluster centers to correspond to the given membership
+## function values.
+##
+## @seealso{fcm, fcm_init_prototype, fcm_update_membership_fcn, fcm_compute_objective_fcn, fcm_compute_convergence_criterion}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm private
+## Directory: fuzzy-logic-toolkit/inst/private/
+## Filename: fcm_update_cluster_centers.m
+## Last-Modified: 7 July 2012
+
+##------------------------------------------------------------------------------
+## Note: This function is an implementation of Equation 13.5 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 380 (International Edition).
+##------------------------------------------------------------------------------
+
+function V = fcm_update_cluster_centers (Mu_m, X, k)
+
+ V = Mu_m * X;
+ sum_Mu_m = sum (Mu_m');
+
+ if (prod (sum_Mu_m) == 0)
+ error ("division by 0 in function fcm_update_cluster_centers\n");
+ endif
+
+ for i = 1 : k
+ V(i, :) /= sum_Mu_m(i);
+ endfor
+
+endfunction
Added: trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_membership_fcn.m
===================================================================
--- trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_membership_fcn.m (rev 0)
+++ trunk/octave-forge/main/fuzzy-logic-toolkit/inst/private/fcm_update_membership_fcn.m 2012-07-09 13:17:11 UTC (rev 10739)
@@ -0,0 +1,75 @@
+## Copyright (C) 2012 L. Markowsky <lm...@us...>
+##
+## This file is part of the fuzzy-logic-toolkit.
+##
+## The fuzzy-logic-toolkit 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.
+##
+## The fuzzy-logic-toolkit 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 the fuzzy-logic-toolkit; see the file COPYING. If not,
+## see <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{Mu} =} fcm_update_membership_fcn (@var{V}, @var{X}, @var{m}, @var{k}, @var{n}, @var{sqr_dist})
+##
+## Compute Mu for each (cluster center, input point) pair.
+##
+## @seealso{fcm, fcm_init_prototype, fcm_update_cluster_centers, fcm_compute_objective_fcn, fcm_compute_convergence_criterion}
+##
+## @end deftypefn
+
+## Author: L. Markowsky
+## Keywords: fuzzy-logic-toolkit fuzzy partition clustering fcm private
+## Directory: fuzzy-logic-toolkit/inst/private/
+## Filename: fcm_update_membership_fcn.m
+## Last-Modified: 7 July 2012
+
+##------------------------------------------------------------------------------
+## Note: This function is an implementation of Equation 13.4 in
+## Fuzzy Logic: Intelligence, Control and Information, by J. Yen and
+## R. Langari, Prentice Hall, 1999, page 380 (International Edition)
+## and Step 3 of Algorithm 4.1 in Fuzzy and Neural Control,
+## by Robert Babuska, November 2009, p. 63.
+##------------------------------------------------------------------------------
+
+function Mu = fcm_update_membership_fcn (V, X, m, k, n, sqr_dist)
+
+ Mu = zeros (k, n);
+
+ if (min (min (sqr_dist)) > 0)
+ exponent = 1.0 / (m - 1);
+ for i = 1 : k
+ for j = 1 : n
+ summation = 0.0;
+ for l = 1 : k
+ summation += (sqr_dist(i, j) / sqr_dist(l, j))^exponent;
+ endfor
+ if (summation != 0)
+ Mu(i, ...
[truncated message content] |
|
From: <raf...@us...> - 2012-07-09 13:06:54
|
Revision: 10738
http://octave.svn.sourceforge.net/octave/?rev=10738&view=rev
Author: rafavzqz
Date: 2012-07-09 13:06:48 +0000 (Mon, 09 Jul 2012)
Log Message:
-----------
Minor change in the help
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbtform.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbtform.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbtform.m 2012-07-09 13:06:25 UTC (rev 10737)
+++ trunk/octave-forge/extra/nurbs/inst/nrbtform.m 2012-07-09 13:06:48 UTC (rev 10738)
@@ -22,8 +22,8 @@
% The NURBS is transform as defined a transformation matrix of size (4,4),
% such as a rotation, translation or change in scale. The transformation
% matrix can define a single transformation or multiple series of
-% transformations. The matrix can be simple constructed by the functions
-% vecscale, vectrans, vecrotx, vecroty, and vecrotz.
+% transformations. The matrix can be simply constructed by the functions
+% vecscale, vectrans and vecrot, and also vecrotx, vecroty, and vecrotz.
%
% Examples:
%
@@ -34,7 +34,7 @@
%
% See also:
%
-% vecscale, vectrans, vecrotx, vecroty, vecrotz
+% vecscale, vectrans, vecrot, vecrotx, vecroty, vecrotz
%
% Copyright (C) 2000 Mark Spink
% Copyright (C) 2010 Rafael Vazquez
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <raf...@us...> - 2012-07-09 13:06:36
|
Revision: 10737
http://octave.svn.sourceforge.net/octave/?rev=10737&view=rev
Author: rafavzqz
Date: 2012-07-09 13:06:25 +0000 (Mon, 09 Jul 2012)
Log Message:
-----------
Plot the 3D control mesh
Modified Paths:
--------------
trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m
Modified: trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m
===================================================================
--- trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m 2012-07-08 18:08:33 UTC (rev 10736)
+++ trunk/octave-forge/extra/nurbs/inst/nrbctrlplot.m 2012-07-09 13:06:25 UTC (rev 10737)
@@ -51,7 +51,33 @@
if (iscell (nurbs.knots))
if (size (nurbs.knots,2) == 3)
- error ('nrbctrlplot: not implemented for NURBS volumes');
+ nsub = 100;
+ nrbplot (nurbs, [nsub nsub nsub], 'light', light, 'colormap', cmap);
+ hold on
+% Plot the the control net
+ for ii = 1:size (nurbs.coefs, 2)
+ for jj = 1:size (nurbs.coefs, 3)
+ coefs = reshape (nurbs.coefs(1:3,ii,jj,:), 3, []);
+ weights = reshape (nurbs.coefs(4,ii,jj,:), 1, []);
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
+ end
+ for kk = 1:size (nurbs.coefs, 4)
+ coefs = reshape (nurbs.coefs(1:3,ii,:,kk), 3, []);
+ weights = reshape (nurbs.coefs(4,ii,:,kk), 1, []);
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
+ end
+ end
+ for jj = 1:size (nurbs.coefs, 3)
+ for kk = 1:size (nurbs.coefs, 4)
+ coefs = reshape (nurbs.coefs(1:3,:,jj,kk), 3, []);
+ weights = reshape (nurbs.coefs(4,:,jj,kk), 1, []);
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'k--')
+ plot3 (coefs(1,:)./weights, coefs(2,:)./weights, coefs(3,:)./weights,'r.','MarkerSize',20)
+ end
+ end
+
elseif (size (nurbs.knots,2) == 2) % plot a NURBS surface
nsub = 100;
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <prn...@us...> - 2012-07-08 18:08:39
|
Revision: 10736
http://octave.svn.sourceforge.net/octave/?rev=10736&view=rev
Author: prnienhuis
Date: 2012-07-08 18:08:33 +0000 (Sun, 08 Jul 2012)
Log Message:
-----------
pre_install.m checks JAVA_HOME before trying to compile Java pkg. Works in *nix & Windows; Mac OSX will follow later
Added Paths:
-----------
trunk/octave-forge/extra/java/pre_install.m
Added: trunk/octave-forge/extra/java/pre_install.m
===================================================================
--- trunk/octave-forge/extra/java/pre_install.m (rev 0)
+++ trunk/octave-forge/extra/java/pre_install.m 2012-07-08 18:08:33 UTC (rev 10736)
@@ -0,0 +1,81 @@
+## Copyright (C) 2012 Philip Nienhuis <prn...@us...>
+##
+## 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 Octave; see the file COPYING. If not, see
+## <http://www.gnu.org/licenses/>.
+
+## Check for JAVA_HOME setting and JDK presence before attempting to
+## install the Java package
+
+## Author: Philip Nienhuis
+## Created: 2012-06-24
+
+function [ ret ] = pre_install ()
+
+ jdk_ok = false;
+ ## Get JAVA_HOME contents
+ jh = getenv ("JAVA_HOME");
+
+ ## Check if has been set at all
+ if (isempty (jh))
+ printf ("\nError while trying to install Java package:\n");
+ printf ("environment variable 'JAVA_HOME' has not been set.\n");
+ printf (" use 'setenv (\"JAVA_HOME\", \"/full/path/to/javaJDK\")'\n");
+ else
+
+ ## Check if JAVA_HOME points to a jvm (that is, given the variety of
+ ## Java installations in the wild, merely check JAVA_HOME/jre/lib)
+ jhd = dir ([jh filesep "jre" filesep "lib"]);
+
+ if (! isempty (jhd))
+ ## Search for a subdir (hopefully <arch>/) containing a
+ ## subdir client/ (*nix) or a file jvm.cfg (Windows)
+ ijhd = find (cell2mat ({jhd.isdir}));
+ ii = 3; # Ignore current and parent dirs
+ while ii < numel (ijhd)
+ jhsd = dir ([jh filesep "jre" filesep "lib" filesep jhd(ijhd(ii)).name]);
+ ## Check if client is a subdir (hopefully of <arch>/)
+ id = strmatch ("client", {jhsd.name});
+ if ((! isempty (id)) && jhsd(id).isdir)
+ cl_dir = [jh filesep "jre" filesep "lib" filesep jhd(ijhd(ii)).name filesep "client"];
+ jhcsd = dir (cl_dir);
+ ## Check for a libjvm* file inside. Should work if it's a link too.
+ jdk_ok = ! isempty (strmatch ("libjvm", {jhcsd.name}));
+ if (! jdk_ok);
+ printf (" No libjvm library found in %s\n", cl_dir);
+ endif
+ endif
+ ## Below line especially for Windows installations
+ jdk_ok = jdk_ok || (! isempty (strmatch ("JVM.CFG", upper ({jhsd.name}))));
+ if (jdk_ok); ii += numel (ijhd); endif
+ ++ii;
+ endwhile
+ endif
+ if (! jdk_ok);
+ printf ("\nError while trying to install Java package:\n");
+ printf ("JAVA_HOME environment variable does not properly point to a JDK\n");
+ endif
+ endif
+
+ if (! jdk_ok);
+ printf (" Hint:\n");
+ printf (" JAVA_HOME should usually be set such that either:\n");
+ printf (" (on *nix:)\n");
+ printf (" <JAVA_HOME>/jre/lib/<arch>/client/ contains libjvm.so (file or symlink)\n");
+ printf (" (on Windows:) \n");
+ printf (" <JAVA_HOME>/jre/lib/<arch>/ contains a file jvm.cfg\n");
+ printf (" (<arch> depends on your system hardware, can be i386, x86_64, alpha, arm, ...)\n\n");
+ error ("Aborting pkg install");
+ endif
+
+endfunction
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <prn...@us...> - 2012-07-08 18:04:56
|
Revision: 10735
http://octave.svn.sourceforge.net/octave/?rev=10735&view=rev
Author: prnienhuis
Date: 2012-07-08 18:04:50 +0000 (Sun, 08 Jul 2012)
Log Message:
-----------
Eeference added
Modified Paths:
--------------
trunk/octave-forge/extra/java/inst/cell2mlstr.m
Modified: trunk/octave-forge/extra/java/inst/cell2mlstr.m
===================================================================
--- trunk/octave-forge/extra/java/inst/cell2mlstr.m 2012-07-08 18:04:09 UTC (rev 10734)
+++ trunk/octave-forge/extra/java/inst/cell2mlstr.m 2012-07-08 18:04:50 UTC (rev 10735)
@@ -18,6 +18,8 @@
## separated by EOL
## Author: Philip <Philip@DESKPRN>
+## Based on a suggestion by D. Bateman,
+## https://savannah.gnu.org/bugs/?func=detailitem&item_id=31468#comment4
## Created: 2012-06-29
function [ ret ] = cell2mlstr (cstr)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <prn...@us...> - 2012-07-08 18:04:15
|
Revision: 10734
http://octave.svn.sourceforge.net/octave/?rev=10734&view=rev
Author: prnienhuis
Date: 2012-07-08 18:04:09 +0000 (Sun, 08 Jul 2012)
Log Message:
-----------
Fixed wrong func name in error messages
Modified Paths:
--------------
trunk/octave-forge/extra/java/inst/listdlg.m
trunk/octave-forge/extra/java/inst/msgbox.m
Modified: trunk/octave-forge/extra/java/inst/listdlg.m
===================================================================
--- trunk/octave-forge/extra/java/inst/listdlg.m 2012-07-06 20:06:01 UTC (rev 10733)
+++ trunk/octave-forge/extra/java/inst/listdlg.m 2012-07-08 18:04:09 UTC (rev 10734)
@@ -61,7 +61,7 @@
## @end deftypefn
## @seealso{errordlg, helpdlg, inputdlg, questdlg, warndlg}
-function varargout = listdlg(varargin)
+function varargout = listdlg (varargin)
if nargin < 2
print_usage ();
Modified: trunk/octave-forge/extra/java/inst/msgbox.m
===================================================================
--- trunk/octave-forge/extra/java/inst/msgbox.m 2012-07-06 20:06:01 UTC (rev 10733)
+++ trunk/octave-forge/extra/java/inst/msgbox.m 2012-07-08 18:04:09 UTC (rev 10734)
@@ -62,7 +62,7 @@
endswitch
if (! ischar (title))
- error ("questdlg: character string expected for title");
+ error ("msgbox: character string expected for title");
endif
ret = java_invoke ('org.octave.JDialogBox', dlg, message, title );
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-06 20:06:07
|
Revision: 10733
http://octave.svn.sourceforge.net/octave/?rev=10733&view=rev
Author: abarth93
Date: 2012-07-06 20:06:01 +0000 (Fri, 06 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/full.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m 2012-07-05 14:11:58 UTC (rev 10732)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m 2012-07-06 20:06:01 UTC (rev 10733)
@@ -1,6 +1,10 @@
-% Create a BaseArray of size sz.
-% sz has at least two elements.
-
+% Create a BaseArray.
+% BA = BaseArray(SZ)
+% Create a BaseArray of size SZ. BaseArray is an abstract class.
+% Derived classes should implement the methods subsref and subsasgn.
+% BaseArray implements several reduction methods such as sum, prod and mean.
+% SZ should have least two elements.
+%
function retval = BaseArray(sz)
self.sz = sz;
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/full.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/full.m 2012-07-05 14:11:58 UTC (rev 10732)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/full.m 2012-07-06 20:06:01 UTC (rev 10733)
@@ -1,4 +1,4 @@
-% make full (dense) array
+% Make full (dense) array.
function F = full(self);
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m 2012-07-05 14:11:58 UTC (rev 10732)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m 2012-07-06 20:06:01 UTC (rev 10733)
@@ -1,6 +1,6 @@
% Reduce array using callback fundtions.
% [S,N] = reduce(SELF,FUNRED,FUNELEM,DIM)
-% reduce array using the function FUNRED applied to all elements
+% Reduce array using the function FUNRED applied to all elements
% after the function FUNELEM was applied along dimension DIM.
function [s,n] = reduce(self,funred,funelem,dim)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <be...@us...> - 2012-07-05 14:12:08
|
Revision: 10732
http://octave.svn.sourceforge.net/octave/?rev=10732&view=rev
Author: benjf5
Date: 2012-07-05 14:11:58 +0000 (Thu, 05 Jul 2012)
Log Message:
-----------
Working on tests; they are all works in progress.
Modified Paths:
--------------
trunk/octave-forge/extra/lssa/lombcoeff.m
trunk/octave-forge/extra/lssa/lscomplex.m
trunk/octave-forge/extra/lssa/lsreal.m
Modified: trunk/octave-forge/extra/lssa/lombcoeff.m
===================================================================
--- trunk/octave-forge/extra/lssa/lombcoeff.m 2012-07-05 13:51:57 UTC (rev 10731)
+++ trunk/octave-forge/extra/lssa/lombcoeff.m 2012-07-05 14:11:58 UTC (rev 10732)
@@ -22,7 +22,16 @@
## @seealso{lombnormcoeff}
## @end deftypefn
+%!shared t, x, o, maxfreq
+%! maxfreq = 4 / ( 2 * pi );
+%! t = linspace(0,8); x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)
+%! - 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t)
+%! + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
+%! assert( lombcoeff(t,x,o(1)) = );
+%! assert( lombcoeff(t,x,o(2)) = );
+%! assert( lombcoeff(t,x,o(3)) = );
+
function coeff = lombcoeff(T, X, o)
theta = atan2(sum(sin(2 .* o .* T )), sum(cos(2.*o.*T)))/ (2 * o );
coeff = ( sum(X .* cos(o .* T - tau))**2)/(sum(cos(o.*T-tau).**2)) +
Modified: trunk/octave-forge/extra/lssa/lscomplex.m
===================================================================
--- trunk/octave-forge/extra/lssa/lscomplex.m 2012-07-05 13:51:57 UTC (rev 10731)
+++ trunk/octave-forge/extra/lssa/lscomplex.m 2012-07-05 14:11:58 UTC (rev 10732)
@@ -23,7 +23,13 @@
## @seealso{lsreal}
## @end deftypefn
+%!test
+%! shared t, x, o, maxfreq
+%! maxfreq = 4 / ( 2 * pi ); t = [0:0.008:8]; x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)- 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t) + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
+%! assert( lscomplex(t,x,maxfreq,2,2), [-0.40075 - 2.36687i , 1.22666 - 2.24390i , 1.93643 - 1.51554i , 2.12505 - 0.95410i ] );
+
+
function transform = lscomplex( t , x , omegamax , ncoeff , noctave )
n = length(t); ## VECTOR ONLY, and since t and x have the same number of entries, there's no problem.
transform = zeros(1,ncoeff*noctave);
Modified: trunk/octave-forge/extra/lssa/lsreal.m
===================================================================
--- trunk/octave-forge/extra/lssa/lsreal.m 2012-07-05 13:51:57 UTC (rev 10731)
+++ trunk/octave-forge/extra/lssa/lsreal.m 2012-07-05 14:11:58 UTC (rev 10732)
@@ -23,6 +23,12 @@
##
## @end deftypefn
+%!shared t, x, o, maxfreq
+%! maxfreq = 4 / ( 2 * pi );
+%! t = linspace(0,8); x = ( 2.*sin(maxfreq.*t) + 3.*sin((3/4)*maxfreq.*t)
+%! - 0.5 .* sin((1/4)*maxfreq.*t) - 0.2 .* cos(maxfreq .* t)
+%! + cos((1/4)*maxfreq.*t)); o = [ maxfreq , 3 / 4 * maxfreq , 1 / 4 * maxfreq ];
+
function transform = lsreal( t, x, omegamax, ncoeff, noctave)
## the R function runs the following command:
## nureal( double X, double Y, int min(X,Y), int ncoeff, int noctave, double omegamax, complex rp)
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-05 13:52:06
|
Revision: 10731
http://octave.svn.sourceforge.net/octave/?rev=10731&view=rev
Author: abarth93
Date: 2012-07-05 13:51:57 +0000 (Thu, 05 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/extra/ncArray/INDEX
trunk/octave-forge/extra/ncArray/inst/@ncArray/coord.m
trunk/octave-forge/extra/ncArray/inst/@ncArray/subsasgn.m
trunk/octave-forge/extra/ncArray/inst/@ncArray/subsref.m
trunk/octave-forge/extra/ncArray/inst/ncCatArray.m
Modified: trunk/octave-forge/extra/ncArray/INDEX
===================================================================
--- trunk/octave-forge/extra/ncArray/INDEX 2012-07-05 13:39:24 UTC (rev 10730)
+++ trunk/octave-forge/extra/ncArray/INDEX 2012-07-05 13:51:57 UTC (rev 10731)
@@ -4,5 +4,26 @@
ncCatArray
nccoord
cached_decompress
+ @ncArray/coord
+ @ncArray/ncArray
+ @ncArray/subsasgn
+ @ncArray/subsref
+ncBaseArray
+ @BaseArray/BaseArray
+ @BaseArray/end
+ @BaseArray/full
+ @BaseArray/isnumeric
+ @BaseArray/max
+ @BaseArray/mean
+ @BaseArray/min
+ @BaseArray/moment
+ @BaseArray/numel
+ @BaseArray/prod
+ @BaseArray/reduce
+ @BaseArray/size
+ @BaseArray/std
+ @BaseArray/sum
+ @BaseArray/sumsq
+ @BaseArray/var
Test script
test_ncarray
\ No newline at end of file
Modified: trunk/octave-forge/extra/ncArray/inst/@ncArray/coord.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@ncArray/coord.m 2012-07-05 13:39:24 UTC (rev 10730)
+++ trunk/octave-forge/extra/ncArray/inst/@ncArray/coord.m 2012-07-05 13:51:57 UTC (rev 10731)
@@ -1,8 +1,10 @@
+% Get coordinate.
+% c = coord(A)
+% Get a array of structures with the coordinate of ncArray A.
+
function c = coord(self)
c = self.coord;
-%'her'
-
% Copyright (C) 2012 Alexander Barth <bar...@gm...>
%
% This program is free software; you can redistribute it and/or modify
Modified: trunk/octave-forge/extra/ncArray/inst/@ncArray/subsasgn.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@ncArray/subsasgn.m 2012-07-05 13:39:24 UTC (rev 10730)
+++ trunk/octave-forge/extra/ncArray/inst/@ncArray/subsasgn.m 2012-07-05 13:51:57 UTC (rev 10731)
@@ -1,3 +1,8 @@
+% Subscripted assignment.
+% subsasgn (A, idx, rhs)
+% Perform the subscripted assignment operation according to the subscript specified by idx.
+% A slice of the NetCDF variable can be saved by using A(index1,index2,...) = rhs;
+
function self = subsasgn(self,idx,x)
self = subsasgn(self.var,idx,x);
Modified: trunk/octave-forge/extra/ncArray/inst/@ncArray/subsref.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@ncArray/subsref.m 2012-07-05 13:39:24 UTC (rev 10730)
+++ trunk/octave-forge/extra/ncArray/inst/@ncArray/subsref.m 2012-07-05 13:51:57 UTC (rev 10731)
@@ -1,3 +1,13 @@
+% Subscripted element selection operation.
+% out = subsref (A, idx)
+% Perform the subscripted element selection operation according to the subscript specified by idx.
+% A slice of the NetCDF variable can be load by using A(index1,index2,...) and attributes
+% can be loaded by A.attribute_name or A.('attribute_name')
+% If index selection is followed by struct selection 'coord', then the coordinates corresponding to the
+% slice are loaded:
+% For example, the coordinate of element A(4,3) are:
+% [lon,lat] = A(4,3).coord;
+
function varargout = subsref(self,idx)
if strcmp(idx(1).type,'()')
Modified: trunk/octave-forge/extra/ncArray/inst/ncCatArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/ncCatArray.m 2012-07-05 13:39:24 UTC (rev 10730)
+++ trunk/octave-forge/extra/ncArray/inst/ncCatArray.m 2012-07-05 13:51:57 UTC (rev 10731)
@@ -7,7 +7,7 @@
% create a concatenated array from variables (varname) in a list of
% netcdf files along dimension dim.Individual elements can be accessed by
% subscribs, e.g. C(2,3) and the corrsponding subset of the appropriate file is loaded
-
+%
% This list of netcdf files can be specified as a cell array (filenames),
% shell wildcard pattern (e.g. file_*.nc) or a function handle
% filenamefun. In this later case, this i-th filename is
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-05 13:39:34
|
Revision: 10730
http://octave.svn.sourceforge.net/octave/?rev=10730&view=rev
Author: abarth93
Date: 2012-07-05 13:39:24 +0000 (Thu, 05 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/end.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/max.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/mean.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/min.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/prod.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/size.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/std.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/sum.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/sumsq.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/var.m
trunk/octave-forge/extra/ncArray/inst/nccoord.m
Added Paths:
-----------
trunk/octave-forge/extra/ncArray/inst/@BaseArray/isnumeric.m
trunk/octave-forge/extra/ncArray/inst/@BaseArray/moment.m
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/BaseArray.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,4 +1,4 @@
-% create a BaseArray of size size sz
+% Create a BaseArray of size sz.
% sz has at least two elements.
function retval = BaseArray(sz)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/end.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/end.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/end.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,3 +1,5 @@
+% Return last index along a dimension.
+
function e = end(self,k,n)
e = size(self,k);
Added: trunk/octave-forge/extra/ncArray/inst/@BaseArray/isnumeric.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/isnumeric.m (rev 0)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/isnumeric.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -0,0 +1,22 @@
+% Test if array is numeric.
+% isn = isnumeric(self)
+
+function isn = isnumeric(self)
+isn = true;
+
+
+% Copyright (C) 2012 Alexander Barth <bar...@gm...>
+%
+% 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 2 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/>.
+
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/max.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/max.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/max.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,6 @@
-% s = max (x, [], dim)
-% compute the maximum along dimension dim
-% See also
-% max
+% Compute the maximum.
+% S = max (X, [], DIM)
+% Compute the maximum along dimension DIM.
function s = max(self,B,varargin)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/mean.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/mean.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/mean.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,7 @@
+% Compute the mean.
% s = mean (X, DIM)
-% compute the mean along dimension dim
-% See also
-% mean
+% Compute the mean along dimension DIM.
+
function s = mean(self,varargin)
funred = @plus;
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/min.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/min.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/min.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,6 @@
-% s = min (x, [], dim)
-% compute the minimum along dimension dim
-% See also
-% min
+% Compute the minimum.
+% S = min (X, [], DIM)
+% Compute the minimum along dimension DIM.
function s = min(self,B,varargin)
Added: trunk/octave-forge/extra/ncArray/inst/@BaseArray/moment.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/moment.m (rev 0)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/moment.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -0,0 +1,36 @@
+% Compute the central moment.
+% M = moment (X, ORDER, DIM)
+% compute the central moment of the given order along dimension DIM.
+
+function s = moment(self,order,varargin)
+
+m = mean(self,varargin{:});
+
+funred = @plus;
+funelem = @(x) (x-m).^order;
+
+[s,n] = reduce(self,funred,funelem,varargin{:});
+
+if isempty(s)
+ s = NaN;
+else
+ s = s/n;
+end
+
+
+
+% Copyright (C) 2012 Alexander Barth <bar...@gm...>
+%
+% 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 2 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/>.
+
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,3 +1,6 @@
+% Number of elements.
+% n = numel(A)
+
function n = numel(self)
n = prod(self.sz);
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/prod.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/prod.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/prod.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,3 +1,6 @@
+% Compute the product.
+% P = prod (X, DIM)
+% Compute the product of all elements along dimension DIM.
function s = prod(self,varargin)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/reduce.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,6 +1,7 @@
-% [s,n] = reduce(self,funred,funelem,dim)
-% reduce array using the function funred applied to all elements
-% after the function funelem was applied along dimension dim
+% Reduce array using callback fundtions.
+% [S,N] = reduce(SELF,FUNRED,FUNELEM,DIM)
+% reduce array using the function FUNRED applied to all elements
+% after the function FUNELEM was applied along dimension DIM.
function [s,n] = reduce(self,funred,funelem,dim)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/size.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/size.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/size.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,3 +1,6 @@
+% Size of array.
+% sz = size(self,dim)
+
function sz = size(self,dim)
sz = self.sz;
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/std.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/std.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/std.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,6 @@
-% s = std (X, OPT, DIM)
-% compute the standard deviation
-% See also
-% std
+% Compute the standard deviation.
+% S = std (X, OPT, DIM)
+
function s = std(self,varargin)
s = sqrt(var(self,varargin{:}));
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/sum.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/sum.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/sum.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,6 @@
-% s = sum (x, dim)
-% compute the sum along dimension dim
-% See also
-% sum
+% Compute the sum.
+% S = sum (X, DIM)
+% Compute the sum along dimension DIM.
function s = sum(self,varargin)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/sumsq.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/sumsq.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/sumsq.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,6 @@
-% s = sumsq (x, dim)
-% compute the sum squared along dimension dim
-% See also
-% sum
+% Compute the sum squared.
+% S = sumsq (X, DIM)
+% Compute the sum squared along dimension DIM.
function s = sumsq(self,varargin)
Modified: trunk/octave-forge/extra/ncArray/inst/@BaseArray/var.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/var.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/var.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,7 +1,7 @@
-% s = var (x, opt, dim)
-% compute the variance along dimension dim
-% See also
-% var
+% Compute the variance.
+% V = var (X, OPT, DIM)
+% Compute the variance along dimension DIM.
+% If OPT is equal to 1, then the variance is bias-corrected.
function s = var(self,opt,varargin)
Modified: trunk/octave-forge/extra/ncArray/inst/nccoord.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/nccoord.m 2012-07-05 13:13:20 UTC (rev 10729)
+++ trunk/octave-forge/extra/ncArray/inst/nccoord.m 2012-07-05 13:39:24 UTC (rev 10730)
@@ -1,4 +1,4 @@
-% Coordinate of a NetCDF variable.
+% Coordinates of a NetCDF variable.
%
% coord = nccoord(filename,varname)
% get coordinates of the variable varname in the
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-05 13:13:29
|
Revision: 10729
http://octave.svn.sourceforge.net/octave/?rev=10729&view=rev
Author: abarth93
Date: 2012-07-05 13:13:20 +0000 (Thu, 05 Jul 2012)
Log Message:
-----------
Added Paths:
-----------
trunk/octave-forge/extra/ncArray/NEWS
Added: trunk/octave-forge/extra/ncArray/NEWS
===================================================================
--- trunk/octave-forge/extra/ncArray/NEWS (rev 0)
+++ trunk/octave-forge/extra/ncArray/NEWS 2012-07-05 13:13:20 UTC (rev 10729)
@@ -0,0 +1,4 @@
+Summary of important user-visible changes for ncArray 1.0.0:
+-------------------------------------------------------------------
+
+ ** Initial release
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-05 13:12:27
|
Revision: 10728
http://octave.svn.sourceforge.net/octave/?rev=10728&view=rev
Author: abarth93
Date: 2012-07-05 13:12:20 +0000 (Thu, 05 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/extra/ncArray/INDEX
trunk/octave-forge/extra/ncArray/inst/@ncArray/ncArray.m
trunk/octave-forge/extra/ncArray/inst/cached_decompress.m
trunk/octave-forge/extra/ncArray/inst/ncCatArray.m
trunk/octave-forge/extra/ncArray/inst/nccoord.m
trunk/octave-forge/extra/ncArray/inst/test_ncarray.m
Modified: trunk/octave-forge/extra/ncArray/INDEX
===================================================================
--- trunk/octave-forge/extra/ncArray/INDEX 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/INDEX 2012-07-05 13:12:20 UTC (rev 10728)
@@ -1,5 +1,8 @@
ncarray >> ncArray
ncArray
+ ncArray
ncCatArray
- ncArray
- nccoord
\ No newline at end of file
+ nccoord
+ cached_decompress
+Test script
+ test_ncarray
\ No newline at end of file
Modified: trunk/octave-forge/extra/ncArray/inst/@ncArray/ncArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@ncArray/ncArray.m 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/inst/@ncArray/ncArray.m 2012-07-05 13:12:20 UTC (rev 10728)
@@ -1,7 +1,9 @@
-% V = ncBaseArray(filename,varname)
-% V = ncBaseArray(filename,varname,'property',value,...)
-% create a ncBaseArray that can be accessed as a normal matlab array.
+% Create an array representing a NetCDF variable.
%
+% V = ncArray(filename,varname)
+% V = ncArray(filename,varname,'property',value,...)
+% create a ncArray that can be accessed as a normal array.
+%
% For read access filename can be compressed if it has the extensions
% ".gz" or ".bz2". It use the function cache_decompress to cache to
% decompressed files.
@@ -32,7 +34,8 @@
%
% V.('_someStrangeAttribute') = 123;
%
-% see also cache_decompress
+% see also cache_decompress, ncCatArray
+% Web: http://modb.oce.ulg.ac.be/mediawiki/index.php/ncArray
% hidded constructor signature:
% data = ncArray(filename,varname)
Modified: trunk/octave-forge/extra/ncArray/inst/cached_decompress.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/cached_decompress.m 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/inst/cached_decompress.m 2012-07-05 13:12:20 UTC (rev 10728)
@@ -1,4 +1,5 @@
% Decompress a file using a cache.
+%
% [fname,success]=cached_decompress(filename)
%
% Input:
@@ -12,11 +13,11 @@
function [fname]=cached_decompress(url)
-global CASHED_GUNZIP_DIR
-global CASHED_GUNZIP_LOG_FID
+global CACHED_DECOMPRESS_DIR
+global CACHED_DECOMPRESS_LOG_FID
-cache_dir = CASHED_GUNZIP_DIR;
+cache_dir = CACHED_DECOMPRESS_DIR;
if isempty(cache_dir)
% cache_dir = fullfile(getenv('HOME'),'tmp','Cache');
cache_dir = fullfile(getenv('HOME'),'tmp','Cache');
@@ -27,7 +28,7 @@
if exist(cache_dir,'dir') ~= 7
error(['cache directory for compressed files does not exist. '...
'Please create the directory %s or change le value of the '...
- 'global variable CASHED_GUNZIP_DIR'],cache_dir);
+ 'global variable CACHED_DECOMPRESS_DIR'],cache_dir);
end
else
fname = url;
@@ -36,7 +37,7 @@
% where to print logs? default to screen
-fid=CASHED_GUNZIP_LOG_FID;
+fid=CACHED_DECOMPRESS_LOG_FID;
if (isempty(fid))
fid=1;
Modified: trunk/octave-forge/extra/ncArray/inst/ncCatArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/ncCatArray.m 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/inst/ncCatArray.m 2012-07-05 13:12:20 UTC (rev 10728)
@@ -26,6 +26,9 @@
% shell wildcard pattern, while in Matlab rdir is used. The function rdir
% is available from Matlab exchange under BSD license
% (http://www.mathworks.com/matlabcentral/fileexchange/19550).
+%
+% see also cache_decompress, ncArray
+% Web: http://modb.oce.ulg.ac.be/mediawiki/index.php/ncArray
% Author: Alexander Barth (bar...@gm...)
%
Modified: trunk/octave-forge/extra/ncArray/inst/nccoord.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/nccoord.m 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/inst/nccoord.m 2012-07-05 13:12:20 UTC (rev 10728)
@@ -1,4 +1,5 @@
-% Coordinate of a NetCDF variable
+% Coordinate of a NetCDF variable.
+%
% coord = nccoord(filename,varname)
% get coordinates of the variable varname in the
% netcdf file called filename. The netcdf is assumed to
Modified: trunk/octave-forge/extra/ncArray/inst/test_ncarray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/test_ncarray.m 2012-07-04 15:19:39 UTC (rev 10727)
+++ trunk/octave-forge/extra/ncArray/inst/test_ncarray.m 2012-07-05 13:12:20 UTC (rev 10728)
@@ -1,4 +1,4 @@
-% test ncBaseArray, ncCatArray and ncArray
+% Test ncBaseArray, ncCatArray and ncArray.
function test_ncarray()
varname = 'SST';
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <par...@us...> - 2012-07-04 15:19:50
|
Revision: 10727
http://octave.svn.sourceforge.net/octave/?rev=10727&view=rev
Author: paramaniac
Date: 2012-07-04 15:19:39 +0000 (Wed, 04 Jul 2012)
Log Message:
-----------
control-devel: work on docstrings
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/inst/arx.m
trunk/octave-forge/extra/control-devel/inst/moen4.m
trunk/octave-forge/extra/control-devel/inst/moesp.m
trunk/octave-forge/extra/control-devel/inst/n4sid.m
Modified: trunk/octave-forge/extra/control-devel/inst/arx.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/arx.m 2012-07-04 10:23:38 UTC (rev 10726)
+++ trunk/octave-forge/extra/control-devel/inst/arx.m 2012-07-04 15:19:39 UTC (rev 10727)
@@ -20,7 +20,7 @@
## @deftypefnx {Function File} {[@var{sys}, @var{x0}] =} arx (@var{dat}, @var{n}, @var{opt}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}] =} arx (@var{dat}, @var{opt}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}] =} arx (@var{dat}, @var{'na'}, @var{na}, @var{'nb'}, @var{nb})
-## ARX
+## Estimate ARX model using QR factorization.
##
## @strong{Inputs}
## @table @var
@@ -50,16 +50,14 @@
## @end table
##
##
-##
## @strong{Option Keys and Values}
## @table @var
## @item 'na'
-## The desired order of the resulting state-space system @var{sys}.
-## @var{s} > @var{n} > 0.
+## Order of the polynomial A(q) and number of poles.
##
## @item 'nb'
-## The desired order of the resulting state-space system @var{sys}.
-## @var{s} > @var{n} > 0.
+## Order of the polynomial B(q)+1 and number of zeros+1.
+## @var{nb} <= @var{na}.
## @end table
##
##
Modified: trunk/octave-forge/extra/control-devel/inst/moen4.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-07-04 10:23:38 UTC (rev 10726)
+++ trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-07-04 15:19:39 UTC (rev 10727)
@@ -20,9 +20,9 @@
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moen4 (@var{dat}, @var{n}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moen4 (@var{dat}, @var{opt}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moen4 (@var{dat}, @var{n}, @var{opt}, @dots{})
-## Combined method: MOESP algorithm for finding the
-## matrices A and C, and N4SID algorithm for
-## finding the matrices B and D.
+## Estimate state-space model using combined subspace method:
+## MOESP algorithm for finding the matrices A and C,
+## and N4SID algorithm for finding the matrices B and D.
##
## @strong{Inputs}
## @table @var
@@ -93,19 +93,40 @@
## Hankel matrices.
## @end table
##
-##
## @item 'tol'
-## The desired order of the resulting reduced order system @var{Gr}.
+## Absolute tolerance used for determining an estimate of
+## the system order. If @var{tol} >= 0, the estimate is
+## indicated by the index of the last singular value greater
+## than or equal to @var{tol}. (Singular values less than @var{tol}
+## are considered as zero.) When @var{tol} = 0, an internally
+## computed default value, @var{tol} = @var{s}*@var{eps}*SV(1), is used,
+## where SV(1) is the maximal singular value, and @var{eps} is
+## the relative machine precision.
+## When @var{tol} < 0, the estimate is indicated by the
+## index of the singular value that has the largest
+## logarithmic gap to its successor.
##
-##
## @item 'rcond'
-## The desired order of the resulting reduced order system @var{Gr}.
+## The tolerance to be used for estimating the rank of
+## matrices. If the user sets @var{rcond} > 0, the given value
+## of @var{rcond} is used as a lower bound for the reciprocal
+## condition number; an m-by-n matrix whose estimated
+## condition number is less than 1/@var{rcond} is considered to
+## be of full rank. If the user sets @var{rcond} <= 0, then an
+## implicitly computed, default tolerance, defined by
+## @var{rcond} = m*n*@var{eps}, is used instead, where @var{eps} is the
+## relative machine precision.
##
-##
## @item 'confirm'
-## The desired order of the resulting reduced order system @var{Gr}.
+## Specifies whether or not the user's confirmation of the
+## system order estimate is desired, as follows:
+## @table @var
+## @item true
+## User's confirmation.
+## @item false
+## No confirmation. Default value.
+## @end table
##
-##
## @item 'noise'
## The desired type of noise input channels.
## @table @var
@@ -170,28 +191,28 @@
## @end iftex
## @ifnottex
## @example
-## x[k+1] = A x[k] + B u[k]
+## ^ ^ ^
+## x[k+1] = A x[k] + B u[k] + K(y[k] - y[k])
+## ^ ^
## y[k] = C x[k] + D u[k]
## @end example
## @end ifnottex
-## @end table
##
-##
-## @item 'method'
-## Specifies the computational approach to be used.
-## Valid values corresponding to this key are:
-## @table @var
-## @item 'descriptor'
-## Use the inverse free descriptor system approach.
-## @item 'standard'
-## Use the inversion based standard approach.
-## @item 'auto'
-## Switch automatically to the inverse free
-## descriptor approach in case of badly conditioned
-## feedthrough matrices in V or W. Default method.
+## @iftex
+## @tex
+## $$ \\widehat{x}_{k+1} = (A-KC) \\widehat{x}_k + (B-KD) u_k + K y_k $$
+## $$ \\widehat{y}_k = C \\widehat{x}_k + D u_k + 0 y_k $$
+## @end tex
+## @end iftex
+## @ifnottex
+## @example
+## ^ ^
+## x[k+1] = (A-KC) x[k] + (B-KD) u[k] + K y[k]
+## ^ ^
+## y[k] = C x[k] + D u[k] + 0 y[k]
+## @end example
+## @end ifnottex
## @end table
-##
-##
## @end table
##
##
Modified: trunk/octave-forge/extra/control-devel/inst/moesp.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moesp.m 2012-07-04 10:23:38 UTC (rev 10726)
+++ trunk/octave-forge/extra/control-devel/inst/moesp.m 2012-07-04 15:19:39 UTC (rev 10727)
@@ -20,6 +20,7 @@
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moesp (@var{dat}, @var{n}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moesp (@var{dat}, @var{opt}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} moesp (@var{dat}, @var{n}, @var{opt}, @dots{})
+## Estimate state-space model using MOESP algorithm.
## MOESP: Multivariable Output Error State sPace.
## @end deftypefn
Modified: trunk/octave-forge/extra/control-devel/inst/n4sid.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/n4sid.m 2012-07-04 10:23:38 UTC (rev 10726)
+++ trunk/octave-forge/extra/control-devel/inst/n4sid.m 2012-07-04 15:19:39 UTC (rev 10727)
@@ -20,6 +20,7 @@
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} n4sid (@var{dat}, @var{n}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} n4sid (@var{dat}, @var{opt}, @dots{})
## @deftypefnx {Function File} {[@var{sys}, @var{x0}, @var{info}] =} n4sid (@var{dat}, @var{n}, @var{opt}, @dots{})
+## Estimate state-space model using N4SID algorithm.
## N4SID: Numerical algorithm for Subspace State Space System IDentification.
## @end deftypefn
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <sch...@us...> - 2012-07-04 10:23:48
|
Revision: 10726
http://octave.svn.sourceforge.net/octave/?rev=10726&view=rev
Author: schloegl
Date: 2012-07-04 10:23:38 +0000 (Wed, 04 Jul 2012)
Log Message:
-----------
fix normalization when A0 is not the identity matrix
Modified Paths:
--------------
trunk/octave-forge/extra/tsa/inst/mvfilter.m
Modified: trunk/octave-forge/extra/tsa/inst/mvfilter.m
===================================================================
--- trunk/octave-forge/extra/tsa/inst/mvfilter.m 2012-07-03 15:39:22 UTC (rev 10725)
+++ trunk/octave-forge/extra/tsa/inst/mvfilter.m 2012-07-04 10:23:38 UTC (rev 10726)
@@ -1,4 +1,3 @@
-
function [x,z]=mvfilter(B,A,x,z)
% Multi-variate filter function
%
@@ -30,7 +29,7 @@
% see also: MVAR, FILTER
% $Id$
-% Copyright (C) 1996-2003 by Alois Schloegl <a.s...@ie...>
+% Copyright (C) 1996-2003,2010,2012 by Alois Schloegl <alo...@is...>
%
% 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
@@ -84,18 +83,13 @@
end;
%%%%% normalization to A{1}=I;
-if p<=q,
- for k=1:p,
- %A{k}=A{k}/A{1};
- A(:,k*M+(1:M)) = A(:,k*M+(1:M)) / A(:,1:M);
+ for k=1:p,
+ A(:,k*M+(1:M)) = A(:,1:M) \ A(:,k*M+(1:M));
end;
- A(:,1:M) = eye(M);
-else
for k=0:q,
- %B{k}=B{k}/A{1};
- B(:,k*M+(1:M)) = B(:,k*M+(1:M)) / A(:,1:M);
+ B(:,k*M+(1:M)) = A(:,1:M) \ B(:,k*M+(1:M));
end;
-end;
+ A(:,1:M) = eye(M);
for k = 1:N,
acc = B(:,1:M) * x(:,k) + z(:,1); % / A{1};
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-03 15:39:32
|
Revision: 10725
http://octave.svn.sourceforge.net/octave/?rev=10725&view=rev
Author: abarth93
Date: 2012-07-03 15:39:22 +0000 (Tue, 03 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/extra/ncArray/inst/@ncBaseArray/ncBaseArray.m
Added Paths:
-----------
trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m
Added: trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m (rev 0)
+++ trunk/octave-forge/extra/ncArray/inst/@BaseArray/numel.m 2012-07-03 15:39:22 UTC (rev 10725)
@@ -0,0 +1,3 @@
+function n = numel(self)
+
+n = prod(self.sz);
Modified: trunk/octave-forge/extra/ncArray/inst/@ncBaseArray/ncBaseArray.m
===================================================================
--- trunk/octave-forge/extra/ncArray/inst/@ncBaseArray/ncBaseArray.m 2012-07-03 15:23:26 UTC (rev 10724)
+++ trunk/octave-forge/extra/ncArray/inst/@ncBaseArray/ncBaseArray.m 2012-07-03 15:39:22 UTC (rev 10725)
@@ -44,7 +44,7 @@
prop = varargin;
for i=1:2:length(prop)
- if strcmp(prop{i},'tooBigToLoad ')
+ if strcmp(prop{i},'tooBigToLoad')
self.tooBigToLoad = prop{i+1};
end
end
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
|
From: <aba...@us...> - 2012-07-03 15:23:37
|
Revision: 10724
http://octave.svn.sourceforge.net/octave/?rev=10724&view=rev
Author: abarth93
Date: 2012-07-03 15:23:26 +0000 (Tue, 03 Jul 2012)
Log Message:
-----------
Modified Paths:
--------------
trunk/octave-forge/main/octcdf/INDEX
trunk/octave-forge/main/octcdf/inst/ncinfo.m
trunk/octave-forge/main/octcdf/inst/ncread.m
trunk/octave-forge/main/octcdf/inst/ncwrite.m
Modified: trunk/octave-forge/main/octcdf/INDEX
===================================================================
--- trunk/octave-forge/main/octcdf/INDEX 2012-07-03 14:57:21 UTC (rev 10723)
+++ trunk/octave-forge/main/octcdf/INDEX 2012-07-03 15:23:26 UTC (rev 10724)
@@ -22,6 +22,9 @@
ncvar
netcdf
ncisrecord
+ ncread
+ ncwrite
+ ncinfo
Examples
example_netcdf
example_opendap
Modified: trunk/octave-forge/main/octcdf/inst/ncinfo.m
===================================================================
--- trunk/octave-forge/main/octcdf/inst/ncinfo.m 2012-07-03 14:57:21 UTC (rev 10723)
+++ trunk/octave-forge/main/octcdf/inst/ncinfo.m 2012-07-03 15:23:26 UTC (rev 10724)
@@ -1,3 +1,5 @@
+% Information about complete netCDF file or variable.
+%
% finfo = ncinfo(filename)
% vinfo = ncinfo(filename,varname)
% return information about complete netCDF file (filename) or about
Modified: trunk/octave-forge/main/octcdf/inst/ncread.m
===================================================================
--- trunk/octave-forge/main/octcdf/inst/ncread.m 2012-07-03 14:57:21 UTC (rev 10723)
+++ trunk/octave-forge/main/octcdf/inst/ncread.m 2012-07-03 15:23:26 UTC (rev 10724)
@@ -1,3 +1,5 @@
+% Read a NetCDF variable.
+%
% x = ncread(filename,varname)
% x = ncread(filename,varname,start,count,stride)
% read the variable varname from file filename.
Modified: trunk/octave-forge/main/octcdf/inst/ncwrite.m
===================================================================
--- trunk/octave-forge/main/octcdf/inst/ncwrite.m 2012-07-03 14:57:21 UTC (rev 10723)
+++ trunk/octave-forge/main/octcdf/inst/ncwrite.m 2012-07-03 15:23:26 UTC (rev 10724)
@@ -1,3 +1,5 @@
+% Write a NetCDF variable.
+%
% ncwrite(filename,varname,x)
% ncwrite(filename,varname,x,start,stride)
% write the variable varname to file filename.
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.
