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: <prn...@us...> - 2012-06-06 21:08:58
|
Revision: 10573 http://octave.svn.sourceforge.net/octave/?rev=10573&view=rev Author: prnienhuis Date: 2012-06-06 21:08:51 +0000 (Wed, 06 Jun 2012) Log Message: ----------- Improved & simplified interface check and made it less verbose Modified Paths: -------------- trunk/octave-forge/main/io/inst/xlsopen.m Modified: trunk/octave-forge/main/io/inst/xlsopen.m =================================================================== --- trunk/octave-forge/main/io/inst/xlsopen.m 2012-06-06 21:07:54 UTC (rev 10572) +++ trunk/octave-forge/main/io/inst/xlsopen.m 2012-06-06 21:08:51 UTC (rev 10573) @@ -105,293 +105,303 @@ ## rediscovery of interfaces between xlsopen calls, e.g. javaclasspath changes) ## 2011-09-08 Minor code cleanup ## 2012-01-26 Fixed "seealso" help string +## 2012-06-06 Improved interface detection logic. No more messages if same interface is +## requested & used consecutively ## -## Latest subfunction update: 2012-03-21 +## Latest subfunction update: 2012-06-06 function [ xls ] = xlsopen (filename, xwrite=0, reqinterface=[]) - persistent xlsinterfaces; persistent chkintf; - # xlsinterfaces.<intf> = [] (not yet checked), 0 (found to be unsupported) or 1 (OK) - if (isempty (chkintf)); - chkintf = 1; - xlsinterfaces = struct ('COM', [], 'POI', [], 'JXL', [], 'OXS', [], 'UNO', []); - endif + persistent xlsinterfaces; persistent chkintf; persistent lastintf; + # xlsinterfaces.<intf> = [] (not yet checked), 0 (found to be unsupported) or 1 (OK) + if (isempty (chkintf)); + chkintf = 1; + xlsinterfaces = struct ('COM', [], 'POI', [], 'JXL', [], 'OXS', [], 'UNO', []); + endif + if (isempty (lastintf)) + lastintf = "---"; + endif - xlssupport = 0; + xlssupport = 0; - if (nargout < 1) - usage ("XLS = xlsopen (Xlfile [, Rw] [, reqintf]). But no return argument specified!"); - endif + if (nargout < 1) + usage ("XLS = xlsopen (Xlfile [, Rw] [, reqintf]). But no return argument specified!"); + endif - if (~(islogical (xwrite) || isnumeric (xwrite))) - usage ("Numerical or logical value expected for arg # 2") - endif + if (~(islogical (xwrite) || isnumeric (xwrite))) + usage ("Numerical or logical value expected for arg # 2") + endif - if (~isempty (reqinterface)) - if ~(ischar (reqinterface) || iscell (reqinterface)), usage ("Arg # 3 not recognized"); endif - # Turn arg3 into cell array if needed - if (~iscell (reqinterface)), reqinterface = {reqinterface}; endif - xlsinterfaces.COM = 0; xlsinterfaces.POI = 0; xlsinterfaces.JXL = 0; - xlsinterfaces.OXS = 0; xlsinterfaces.UNO = 0; - for ii=1:numel (reqinterface) - reqintf = toupper (reqinterface {ii}); - # Try to invoke requested interface(s) for this call. Check if it - # is supported anyway by emptying the corresponding var. - if (strcmp (reqintf, 'COM')) - xlsinterfaces.COM = []; - elseif (strcmp (reqintf, 'POI')) - xlsinterfaces.POI = []; - elseif (strcmp (reqintf, 'JXL')) - xlsinterfaces.JXL = []; - elseif (strcmp (reqintf, 'OXS')) - xlsinterfaces.OXS = []; - elseif (strcmp (reqintf, 'UNO')) - xlsinterfaces.UNO = []; - else - usage (sprintf ("Unknown .xls interface \"%s\" requested. Only COM, POI, JXL, OXS or UNO supported\n", reqinterface{})); - endif - endfor - printf ("Checking requested interface(s):\n"); - xlsinterfaces = getxlsinterfaces (xlsinterfaces); - # Well, is/are the requested interface(s) supported on the system? - # FIXME check for multiple interfaces - xlsintf_cnt = 0; - for ii=1:numel (reqinterface) - if (~xlsinterfaces.(toupper (reqinterface{ii}))) - # No it aint - printf ("%s is not supported.\n", upper (reqinterface{ii})); - else - ++xlsintf_cnt; - endif - endfor - # Reset interface check indicator if no requested support found - if (~xlsintf_cnt) - chkintf = []; - xls = []; - return - endif - endif - - # Var xwrite is really used to avoid creating files when wanting to read, or - # not finding not-yet-existing files when wanting to write. - - # Check if Excel file exists. Adapt file open mode for readwrite argument - if (xwrite), fmode = 'r+b'; else fmode = 'rb'; endif - fid = fopen (filename, fmode); - if (fid < 0) # File doesn't exist... - if (~xwrite) # ...which obviously is fatal for reading... - error ( sprintf ("File %s not found\n", filename)); - else # ...but for writing, we need more info: - fid = fopen (filename, 'rb'); # Check if it exists at all... - if (fid < 0) # File didn't exist yet. Simply create it - printf ("Creating file %s\n", filename); - xwrite = 3; - else # File exists, but is not writable => Error - fclose (fid); # Do not forget to close the handle neatly - error (sprintf ("Write mode requested but file %s is not writable\n", filename)) - endif - endif - else - # Close file anyway to avoid COM or Java errors - fclose (fid); - endif - -# Check for the various Excel interfaces. No problem if they've already -# been checked, getxlsinterfaces (far below) just returns immediately then. - xlsinterfaces = getxlsinterfaces (xlsinterfaces); + if (~isempty (reqinterface)) + if ~(ischar (reqinterface) || iscell (reqinterface)), usage ("Arg # 3 not recognized"); endif + # Turn arg3 into cell array if needed + if (~iscell (reqinterface)), reqinterface = {reqinterface}; endif + ## Check if previously used interface matches a requested interface + if (isempty (regexpi (reqinterface, lastintf, 'once'){1})) + ## New interface requested + xlsinterfaces.COM = 0; xlsinterfaces.POI = 0; xlsinterfaces.JXL = 0; + xlsinterfaces.OXS = 0; xlsinterfaces.UNO = 0; + for ii=1:numel (reqinterface) + reqintf = toupper (reqinterface {ii}); + # Try to invoke requested interface(s) for this call. Check if it + # is supported anyway by emptying the corresponding var. + if (strcmpi (reqintf, 'COM')) + xlsinterfaces.COM = []; + elseif (strcmpi (reqintf, 'POI')) + xlsinterfaces.POI = []; + elseif (strcmpi (reqintf, 'JXL')) + xlsinterfaces.JXL = []; + elseif (strcmpi (reqintf, 'OXS')) + xlsinterfaces.OXS = []; + elseif (strcmpi (reqintf, 'UNO')) + xlsinterfaces.UNO = []; + else + usage (sprintf ("Unknown .xls interface \"%s\" requested. Only COM, POI, JXL, OXS or UNO supported\n", reqinterface{})); + endif + endfor + printf ("Checking requested interface(s):\n"); + xlsinterfaces = getxlsinterfaces (xlsinterfaces); + # Well, is/are the requested interface(s) supported on the system? + # FIXME check for multiple interfaces + xlsintf_cnt = 0; + for ii=1:numel (reqinterface) + if (~xlsinterfaces.(toupper (reqinterface{ii}))) + # No it aint + printf ("%s is not supported.\n", upper (reqinterface{ii})); + else + ++xlsintf_cnt; + endif + endfor + # Reset interface check indicator if no requested support found + if (~xlsintf_cnt) + chkintf = []; + xls = []; + return + endif + endif + endif -# Supported interfaces determined; Excel file type check moved to seperate interfaces. - chk1 = strcmpi (filename(end-3:end), '.xls'); # Regular (binary) BIFF - chk2 = strcmpi (filename(end-4:end-1), '.xls'); # Zipped XML / OOXML - - # Initialize file ptr struct - xls = struct ("xtype", 'NONE', "app", [], "filename", [], "workbook", [], "changed", 0, "limits", []); + # Var xwrite is really used to avoid creating files when wanting to read, or + # not finding not-yet-existing files when wanting to write. - # Keep track of which interface is selected - xlssupport = 0; + # Check if Excel file exists. Adapt file open mode for readwrite argument + if (xwrite), fmode = 'r+b'; else fmode = 'rb'; endif + fid = fopen (filename, fmode); + if (fid < 0) # File doesn't exist... + if (~xwrite) # ...which obviously is fatal for reading... + error ( sprintf ("File %s not found\n", filename)); + else # ...but for writing, we need more info: + fid = fopen (filename, 'rb'); # Check if it exists at all... + if (fid < 0) # File didn't exist yet. Simply create it + printf ("Creating file %s\n", filename); + xwrite = 3; + else # File exists, but is not writable => Error + fclose (fid); # Do not forget to close the handle neatly + error (sprintf ("Write mode requested but file %s is not writable\n", filename)) + endif + endif + else + # Close file anyway to avoid COM or Java errors + fclose (fid); + endif + + # Check for the various Excel interfaces. No problem if they've already + # been checked, getxlsinterfaces (far below) just returns immediately then. + xlsinterfaces = getxlsinterfaces (xlsinterfaces); - # Interface preference order is defined below: currently COM -> POI -> JXL -> OXS -> UNO - if (xlsinterfaces.COM && ~xlssupport) - # Excel functioning has been tested above & file exists, so we just invoke it - app = actxserver ("Excel.Application"); - try # Because Excel itself can still crash on file formats etc. - app.Application.DisplayAlerts = 0; - if (xwrite < 2) - # Open workbook - wb = app.Workbooks.Open (canonicalize_file_name (filename)); - elseif (xwrite > 2) - # Create a new workbook - wb = app.Workbooks.Add (); - ### Uncommenting the below statement can be useful in multi-user environments. - ### Be sure to uncomment correspondig stanza in xlsclose to avoid zombie Excels - # wb.SaveAs (canonicalize_file_name (filename)) - endif - xls.app = app; - xls.xtype = 'COM'; - xls.workbook = wb; - xls.filename = filename; - xlssupport += 1; - catch - warning ( sprintf ("ActiveX error trying to open or create file %s\n", filename)); - app.Application.DisplayAlerts = 1; - app.Quit (); - delete (app); - end_try_catch - endif - - if (xlsinterfaces.POI && ~xlssupport) - if ~(chk1 || chk2) - error ("Unsupported file format for Apache POI.") - endif - # Get handle to workbook - try - if (xwrite > 2) - if (chk1) - wb = java_new ('org.apache.poi.hssf.usermodel.HSSFWorkbook'); - elseif (chk2) - wb = java_new ('org.apache.poi.xssf.usermodel.XSSFWorkbook'); - endif - xls.app = 'new_POI'; - else - xlsin = java_new ('java.io.FileInputStream', filename); - wb = java_invoke ('org.apache.poi.ss.usermodel.WorkbookFactory', 'create', xlsin); - xls.app = xlsin; - endif - xls.xtype = 'POI'; - xls.workbook = wb; - xls.filename = filename; - xlssupport += 2; - catch - clear xlsin; - if (xlsinterfaces.JXL) - printf ('Couldn''t open file %s using POI; trying Excel''95 format with JXL...\n', filename); - endif - end_try_catch - endif - - if (xlsinterfaces.JXL && ~xlssupport) - if (~chk1) - error ("JXL can only read reliably from .xls files") - endif - try - xlsin = java_new ('java.io.File', filename); - if (xwrite > 2) - # Get handle to new xls-file - wb = java_invoke ('jxl.Workbook', 'createWorkbook', xlsin); - else - # Open existing file - wb = java_invoke ('jxl.Workbook', 'getWorkbook', xlsin); - endif - xls.xtype = 'JXL'; - xls.app = xlsin; - xls.workbook = wb; - xls.filename = filename; - xlssupport += 4; - catch - clear xlsin; - if (xlsinterfaces.POI) - printf ('... No luck with JXL either, unsupported file format.\n', filename); - endif - end_try_catch - endif + # Supported interfaces determined; Excel file type check moved to seperate interfaces. + chk1 = strcmpi (filename(end-3:end), '.xls'); # Regular (binary) BIFF + chk2 = strcmpi (filename(end-4:end-1), '.xls'); # Zipped XML / OOXML + + # Initialize file ptr struct + xls = struct ("xtype", 'NONE', "app", [], "filename", [], "workbook", [], "changed", 0, "limits", []); - if (xlsinterfaces.OXS && ~xlssupport) - if (~chk1) - error ("OXS can only read from .xls files") - endif - try - wb = javaObject ('com.extentech.ExtenXLS.WorkBookHandle', filename); - xls.xtype = 'OXS'; - xls.app = 'void - OpenXLS'; - xls.workbook = wb; - xls.filename = filename; - xlssupport += 8; - catch - printf ('Unsupported file format for OpenXLS - %s\n'); - end_try_catch - endif + # Keep track of which interface is selected + xlssupport = 0; - if (xlsinterfaces.UNO && ~xlssupport) - # First, the file name must be transformed into a URL - if (~isempty (strmatch ("file:///", filename)) || ~isempty (strmatch ("http:///", filename))... - || ~isempty (strmatch ("ftp:///", filename)) || ~isempty (strmatch ("www:///", filename))) - # Seems in proper shape for OOo (at first sight) - else - # Transform into URL form - fname = canonicalize_file_name (strsplit (filename, filesep){end}); - # On Windows, change backslash file separator into forward slash - if (strcmp (filesep, "\\")) - tmp = strsplit (fname, filesep); - flen = numel (tmp); - tmp(2:2:2*flen) = tmp; - tmp(1:2:2*flen) = '/'; - fname = [ tmp{:} ]; - endif - filename = [ 'file://' fname ]; - endif - try - xContext = java_invoke ("com.sun.star.comp.helper.Bootstrap", "bootstrap"); - xMCF = xContext.getServiceManager (); - oDesktop = xMCF.createInstanceWithContext ("com.sun.star.frame.Desktop", xContext); - # Workaround for <UNOruntime>.queryInterface(): - unotmp = java_new ('com.sun.star.uno.Type', 'com.sun.star.frame.XComponentLoader'); - aLoader = oDesktop.queryInterface (unotmp); - # Some trickery as Octave Java cannot create initialized arrays - lProps = javaArray ('com.sun.star.beans.PropertyValue', 1); - lProp = java_new ('com.sun.star.beans.PropertyValue', "Hidden", 0, true, []); - lProps(1) = lProp; - if (xwrite > 2) - xComp = aLoader.loadComponentFromURL ("private:factory/scalc", "_blank", 0, lProps); - else - xComp = aLoader.loadComponentFromURL (filename, "_blank", 0, lProps); - endif - # Workaround for <UNOruntime>.queryInterface(): - unotmp = java_new ('com.sun.star.uno.Type', 'com.sun.star.sheet.XSpreadsheetDocument'); - xSpdoc = xComp.queryInterface (unotmp); - # save in ods struct: - xls.xtype = 'UNO'; - xls.workbook = xSpdoc; # Needed to be able to close soffice in odsclose() - xls.filename = filename; - xls.app.xComp = xComp; # Needed to be able to close soffice in odsclose() - xls.app.aLoader = aLoader; # Needed to be able to close soffice in odsclose() - xls.odfvsn = 'UNO'; - xlssupport += 16; - catch - error ('Couldn''t open file %s using UNO', filename); - end_try_catch - endif + # Interface preference order is defined below: currently COM -> POI -> JXL -> OXS -> UNO + if (xlsinterfaces.COM && ~xlssupport) + # Excel functioning has been tested above & file exists, so we just invoke it + app = actxserver ("Excel.Application"); + try # Because Excel itself can still crash on file formats etc. + app.Application.DisplayAlerts = 0; + if (xwrite < 2) + # Open workbook + wb = app.Workbooks.Open (canonicalize_file_name (filename)); + elseif (xwrite > 2) + # Create a new workbook + wb = app.Workbooks.Add (); + ### Uncommenting the below statement can be useful in multi-user environments. + ### Be sure to uncomment correspondig stanza in xlsclose to avoid zombie Excels + # wb.SaveAs (canonicalize_file_name (filename)) + endif + xls.app = app; + xls.xtype = 'COM'; + xls.workbook = wb; + xls.filename = filename; + xlssupport += 1; + lastintf = 'COM'; + catch + warning ( sprintf ("ActiveX error trying to open or create file %s\n", filename)); + app.Application.DisplayAlerts = 1; + app.Quit (); + delete (app); + end_try_catch + endif + + if (xlsinterfaces.POI && ~xlssupport) + if ~(chk1 || chk2) + error ("Unsupported file format for Apache POI.") + endif + # Get handle to workbook + try + if (xwrite > 2) + if (chk1) + wb = java_new ('org.apache.poi.hssf.usermodel.HSSFWorkbook'); + elseif (chk2) + wb = java_new ('org.apache.poi.xssf.usermodel.XSSFWorkbook'); + endif + xls.app = 'new_POI'; + else + xlsin = java_new ('java.io.FileInputStream', filename); + wb = java_invoke ('org.apache.poi.ss.usermodel.WorkbookFactory', 'create', xlsin); + xls.app = xlsin; + endif + xls.xtype = 'POI'; + xls.workbook = wb; + xls.filename = filename; + xlssupport += 2; + lastintf = 'JXL'; + catch + clear xlsin; + if (xlsinterfaces.JXL) + printf ('Couldn''t open file %s using POI; trying Excel''95 format with JXL...\n', filename); + endif + end_try_catch + endif + + if (xlsinterfaces.JXL && ~xlssupport) + if (~chk1) + error ("JXL can only read reliably from .xls files") + endif + try + xlsin = java_new ('java.io.File', filename); + if (xwrite > 2) + # Get handle to new xls-file + wb = java_invoke ('jxl.Workbook', 'createWorkbook', xlsin); + else + # Open existing file + wb = java_invoke ('jxl.Workbook', 'getWorkbook', xlsin); + endif + xls.xtype = 'JXL'; + xls.app = xlsin; + xls.workbook = wb; + xls.filename = filename; + xlssupport += 4; + lastintf = 'POI'; + catch + clear xlsin; + if (xlsinterfaces.POI) + printf ('... No luck with JXL either, unsupported file format.\n', filename); + endif + end_try_catch + endif - # if - # ---- other interfaces - # endif + if (xlsinterfaces.OXS && ~xlssupport) + if (~chk1) + error ("OXS can only read from .xls files") + endif + try + wb = javaObject ('com.extentech.ExtenXLS.WorkBookHandle', filename); + xls.xtype = 'OXS'; + xls.app = 'void - OpenXLS'; + xls.workbook = wb; + xls.filename = filename; + xlssupport += 8; + lastintf = 'OXS'; + catch + printf ('Unsupported file format for OpenXLS - %s\n'); + end_try_catch + endif - if (~xlssupport) - if (isempty (reqinterface)) + if (xlsinterfaces.UNO && ~xlssupport) + # First, the file name must be transformed into a URL + if (~isempty (strmatch ("file:///", filename)) || ~isempty (strmatch ("http:///", filename))... + || ~isempty (strmatch ("ftp:///", filename)) || ~isempty (strmatch ("www:///", filename))) + # Seems in proper shape for OOo (at first sight) + else + # Transform into URL form + fname = canonicalize_file_name (strsplit (filename, filesep){end}); + # On Windows, change backslash file separator into forward slash + if (strcmp (filesep, "\\")) + tmp = strsplit (fname, filesep); + flen = numel (tmp); + tmp(2:2:2*flen) = tmp; + tmp(1:2:2*flen) = '/'; + fname = [ tmp{:} ]; + endif + filename = [ 'file://' fname ]; + endif + try + xContext = java_invoke ("com.sun.star.comp.helper.Bootstrap", "bootstrap"); + xMCF = xContext.getServiceManager (); + oDesktop = xMCF.createInstanceWithContext ("com.sun.star.frame.Desktop", xContext); + # Workaround for <UNOruntime>.queryInterface(): + unotmp = java_new ('com.sun.star.uno.Type', 'com.sun.star.frame.XComponentLoader'); + aLoader = oDesktop.queryInterface (unotmp); + # Some trickery as Octave Java cannot create initialized arrays + lProps = javaArray ('com.sun.star.beans.PropertyValue', 1); + lProp = java_new ('com.sun.star.beans.PropertyValue', "Hidden", 0, true, []); + lProps(1) = lProp; + if (xwrite > 2) + xComp = aLoader.loadComponentFromURL ("private:factory/scalc", "_blank", 0, lProps); + else + xComp = aLoader.loadComponentFromURL (filename, "_blank", 0, lProps); + endif + # Workaround for <UNOruntime>.queryInterface(): + unotmp = java_new ('com.sun.star.uno.Type', 'com.sun.star.sheet.XSpreadsheetDocument'); + xSpdoc = xComp.queryInterface (unotmp); + # save in ods struct: + xls.xtype = 'UNO'; + xls.workbook = xSpdoc; # Needed to be able to close soffice in odsclose() + xls.filename = filename; + xls.app.xComp = xComp; # Needed to be able to close soffice in odsclose() + xls.app.aLoader = aLoader; # Needed to be able to close soffice in odsclose() + xls.odfvsn = 'UNO'; + xlssupport += 16; + lastintf = 'UNO'; + catch + error ('Couldn''t open file %s using UNO', filename); + end_try_catch + endif + + # if + # ---- other interfaces + # endif + + # Rounding up. If none of the xlsinterfaces is supported we're out of luck. + if (~xlssupport) + if (isempty (reqinterface)) printf ("None.\n"); - warning ("No support for Excel .xls I/O"); - else - warning ("File type not supported by %s %s %s %s %s", reqinterface{:}); - endif - xls = []; - else - # From here on xwrite is tracked via xls.changed in the various lower - # level r/w routines and it is only used to determine if an informative - # message is to be given when saving a newly created xls file. - xls.changed = xwrite; + warning ("No support for Excel .xls I/O"); + else + warning ("File type not supported by %s %s %s %s %s", reqinterface{:}); + endif + xls = []; + # Reset found interfaces for re-testing in the next call. Add interfaces if needed. + chkintf = []; + else + # From here on xwrite is tracked via xls.changed in the various lower + # level r/w routines and it is only used to determine if an informative + # message is to be given when saving a newly created xls file. + xls.changed = xwrite; - # Until something was written to existing files we keep status "unchanged". - # xls.changed = 0 (existing/only read from), 1 (existing/data added), 2 (new, - # data added) or 3 (pristine, no data added). - if (xls.changed == 1), xls.changed = 0; endif - endif - -# Rounding up. If none of the xlsinterfaces is supported we're out of luck. - - if (~isempty (reqinterface) || ~xlssupport) - # Reset found interfaces for re-testing in the next call. Add interfaces if needed. - chkintf = []; - endif - + # Until something was written to existing files we keep status "unchanged". + # xls.changed = 0 (existing/only read from), 1 (existing/data added), 2 (new, + # data added) or 3 (pristine, no data added). + if (xls.changed == 1), xls.changed = 0; endif + endif + endfunction @@ -456,31 +466,32 @@ ## 2012-03-21 Print newline if COM found but no Java support ## '' Improved logic for finding out what interfaces to check ## '' Fixed bugs with Java interface checking (tmp1 initialization) +## 2012-06-06 Improved & simplified Java check code function [xlsinterfaces] = getxlsinterfaces (xlsinterfaces) # tmp1 = [] (not initialized), 0 (No Java detected), or 1 (Working Java found) - persistent tmp1 = []; persistent jcp; # Java class path + persistent tmp1 = []; persistent jcp; # Java class path persistent uno_1st_time = 0; - if (isempty (xlsinterfaces.COM) && isempty (xlsinterfaces.POI) && isempty (xlsinterfaces.JXL) + if (isempty (xlsinterfaces.COM) && isempty (xlsinterfaces.POI) && isempty (xlsinterfaces.JXL) && isempty (xlsinterfaces.OXS) && isempty (xlsinterfaces.UNO)) # Looks like first call to xlsopen. Check Java support - printf ("Detected XLS interfaces: "); + printf ("Detected XLS interfaces: "); tmp1 = []; - elseif (isempty (xlsinterfaces.POI) || isempty (xlsinterfaces.JXL) + elseif (isempty (xlsinterfaces.POI) || isempty (xlsinterfaces.JXL) || isempty (xlsinterfaces.OXS) || isempty (xlsinterfaces.UNO)) # Can't be first call. Here one of the Java interfaces is requested if (~tmp1) # Check Java support again tmp1 = []; endif - endif - deflt = 0; + endif + deflt = 0; - # Check if MS-Excel COM ActiveX server runs (only on Windows!) - if (isempty (xlsinterfaces.COM)) - xlsinterfaces.COM = 0; + # Check if MS-Excel COM ActiveX server runs (only on Windows!) + if (isempty (xlsinterfaces.COM)) + xlsinterfaces.COM = 0; if (ispc) try app = actxserver ("Excel.application"); @@ -500,130 +511,126 @@ endif end_try_catch endif - endif + endif - if (isempty (tmp1)) - # Check Java support. First try javaclasspath - try - jcp = javaclasspath ('-all'); # For java pkg > 1.2.7 - if (isempty (jcp)), jcp = javaclasspath; endif # For java pkg < 1.2.8 - # If we get here, at least Java works. Now check for proper version (>= 1.6) - jver = char (java_invoke ('java.lang.System', 'getProperty', 'java.version')); - cjver = strsplit (jver, '.'); - if (sscanf (cjver{2}, '%d') < 6) - warning ("\nJava version might be too old - you need at least Java 6 (v. 1.6.x.x)\n"); - return - endif - # Now check for proper entries in class path. Under *nix the classpath - # must first be split up. In java 1.2.8+ javaclasspath is already a cell array - if (isunix && ~iscell (jcp)); jcp = strsplit (char (jcp), ":"); endif - tmp1 = 1; - catch - # No Java support found - tmp1 = 0; - if ~(isempty (xlsinterfaces.POI) && isempty (xlsinterfaces.JXL)... - && isempty (xlsinterfaces.OXS) && isempty (xlsinterfaces.UNO)) - # Some Java-based interface explicitly requested but Java support is absent - xlsinterfaces.POI = 0; - xlsinterfaces.JXL = 0; - xlsinterfaces.OXS = 0; - xlsinterfaces.UNO = 0; - warning (' No Java support found (no Java JRE? no Java pkg installed AND loaded?'); - else - # No specific Java-based interface requested (first call?) - xlsinterfaces.POI = 0; - xlsinterfaces.JXL = 0; - xlsinterfaces.OXS = 0; - xlsinterfaces.UNO = 0; - printf ("\n"); - return; - endif - end_try_catch - endif + if (isempty (tmp1)) + # Check Java support. First try javaclasspath + try + jcp = javaclasspath ('-all'); # For java pkg > 1.2.7 + if (isempty (jcp)), jcp = javaclasspath; endif # For java pkg < 1.2.8 + # If we get here, at least Java works. Now check for proper version (>= 1.6) + jver = char (java_invoke ('java.lang.System', 'getProperty', 'java.version')); + cjver = strsplit (jver, '.'); + if (sscanf (cjver{2}, '%d') < 6) + warning ("\nJava version might be too old - you need at least Java 6 (v. 1.6.x.x)\n"); + return + endif + # Now check for proper entries in class path. Under *nix the classpath + # must first be split up. In java 1.2.8+ javaclasspath is already a cell array + if (isunix && ~iscell (jcp)); jcp = strsplit (char (jcp), ":"); endif + tmp1 = 1; + catch + # No Java support found + tmp1 = 0; + if (isempty (xlsinterfaces.POI) || isempty (xlsinterfaces.JXL)... + || isempty (xlsinterfaces.OXS) || isempty (xlsinterfaces.UNO)) + # Some or all Java-based interface(s) explicitly requested but no Java support + warning (' No Java support found (no Java JRE? no Java pkg installed AND loaded?'); + endif + # Set Java interfaces to 0 anyway as there's no Java support + xlsinterfaces.POI = 0; + xlsinterfaces.JXL = 0; + xlsinterfaces.OXS = 0; + xlsinterfaces.UNO = 0; + printf ("\n"); + # No more need to try any Java interface + return + end_try_catch + endif - # Try Java & Apache POI - if (isempty (xlsinterfaces.POI)) - xlsinterfaces.POI = 0; - # Check basic .xls (BIFF8) support - jpchk1 = 0; entries1 = {"poi-3", "poi-ooxml-3"}; - # Only under *nix we might use brute force: e.g., strfind (classname, classpath); - # under Windows we need the following more subtle, platform-independent approach: - for ii=1:length (jcp) - for jj=1:length (entries1) - if (isempty (strfind (tolower (jcp{ii}), entries1{jj}))), ++jpchk1; endif - endfor - endfor - if (jpchk1 > 1) - xlsinterfaces.POI = 1; - printf ("POI"); - endif - # Check OOXML support - jpchk2 = 0; entries2 = {"xbean", "poi-ooxml-schemas", "dom4j"}; - for ii=1:length (jcp) - for jj=1:length (entries2) - if (isempty (strfind (lower (jcp{ii}), entries2{jj}))), ++jpchk2; endif - endfor - endfor - if (jpchk2 > 2), printf (" (& OOXML)"); endif - if (xlsinterfaces.POI) - if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif - endif - endif + # Try Java & Apache POI + if (isempty (xlsinterfaces.POI)) + xlsinterfaces.POI = 0; + # Check basic .xls (BIFF8) support + jpchk1 = 0; entries1 = {"poi-3", "poi-ooxml-3"}; + # Only under *nix we might use brute force: e.g., strfind (classname, classpath); + # under Windows we need the following more subtle, platform-independent approach: + for ii=1:length (jcp) + for jj=1:length (entries1) + if (isempty (strfind (tolower (jcp{ii}), entries1{jj}))), ++jpchk1; endif + endfor + endfor + if (jpchk1 > 1) + xlsinterfaces.POI = 1; + printf ("POI"); + endif + # Check OOXML support + jpchk2 = 0; entries2 = {"xbean", "poi-ooxml-schemas", "dom4j"}; + for ii=1:length (jcp) + for jj=1:length (entries2) + if (isempty (strfind (lower (jcp{ii}), entries2{jj}))), ++jpchk2; endif + endfor + endfor + if (jpchk2 > 2), printf (" (& OOXML)"); endif + if (xlsinterfaces.POI) + if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif + endif + endif - # Try Java & JExcelAPI - if (isempty (xlsinterfaces.JXL)) - xlsinterfaces.JXL = 0; - jpchk = 0; entries = {"jxl"}; - for ii=1:length (jcp) - for jj=1:length (entries) - if (isempty (strfind (lower (jcp{ii}), entries{jj}))), ++jpchk; endif - endfor - endfor - if (jpchk > 0) - xlsinterfaces.JXL = 1; - printf ("JXL"); - if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif - endif - endif + # Try Java & JExcelAPI + if (isempty (xlsinterfaces.JXL)) + xlsinterfaces.JXL = 0; + jpchk = 0; entries = {"jxl"}; + for ii=1:length (jcp) + for jj=1:length (entries) + if (isempty (strfind (lower (jcp{ii}), entries{jj}))), ++jpchk; endif + endfor + endfor + if (jpchk > 0) + xlsinterfaces.JXL = 1; + printf ("JXL"); + if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif + endif + endif - # Try Java & OpenXLS - if (isempty (xlsinterfaces.OXS)) - xlsinterfaces.OXS = 0; - jpchk = 0; entries = {"openxls"}; - for ii=1:length (jcp) - for jj=1:length (entries) - if (isempty (strfind (lower (jcp{ii}), entries{jj}))), ++jpchk; endif - endfor - endfor - if (jpchk > 0) - xlsinterfaces.OXS = 1; - printf ("OXS"); - if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif - endif - endif + # Try Java & OpenXLS + if (isempty (xlsinterfaces.OXS)) + xlsinterfaces.OXS = 0; + jpchk = 0; entries = {"openxls"}; + for ii=1:length (jcp) + for jj=1:length (entries) + if (isempty (strfind (lower (jcp{ii}), entries{jj}))), ++jpchk; endif + endfor + endfor + if (jpchk > 0) + xlsinterfaces.OXS = 1; + printf ("OXS"); + if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; endif + endif + endif - # Try Java & UNO - if (isempty (xlsinterfaces.UNO)) - xlsinterfaces.UNO = 0; - # entries0(1) = not a jar but a directory (<00o_install_dir/program/>) - jpchk = 0; entries = {'program', 'unoil', 'jurt', 'juh', 'unoloader', 'ridl'}; - for jj=1:numel (entries) - for ii=1:numel (jcp) - jcplst = strsplit (jcp{ii}, filesep); - jcpentry = jcplst {end}; - if (~isempty (strfind (lower (jcpentry), lower (entries{jj})))); ++jpchk; endif - endfor - endfor - if (jpchk >= numel (entries)) - xlsinterfaces.UNO = 1; - printf ('UNO'); - if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; uno_1st_time = min (++uno_1st_time, 2); endif - endif - endif + # Try Java & UNO + if (isempty (xlsinterfaces.UNO)) + xlsinterfaces.UNO = 0; + # entries0(1) = not a jar but a directory (<00o_install_dir/program/>) + jpchk = 0; entries = {'program', 'unoil', 'jurt', 'juh', 'unoloader', 'ridl'}; + for jj=1:numel (entries) + for ii=1:numel (jcp) + jcplst = strsplit (jcp{ii}, filesep); + jcpentry = jcplst {end}; + if (~isempty (strfind (lower (jcpentry), lower (entries{jj})))); ++jpchk; endif + endfor + endfor + if (jpchk >= numel (entries)) + xlsinterfaces.UNO = 1; + printf ('UNO'); + if (deflt), printf ("; "); else, printf ("*; "); deflt = 1; uno_1st_time = min (++uno_1st_time, 2); endif + endif + endif - # ---- Other interfaces here, similar to the ones above + # ---- Other interfaces here, similar to the ones above - if (deflt), printf ("(* = active interface)\n"); endif + if (deflt), printf ("(* = active interface)\n"); endif ## FIXME the below stanza should be dropped once UNO is stable. # Echo a suitable warning about experimental status: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <prn...@us...> - 2012-06-06 21:08:00
|
Revision: 10572 http://octave.svn.sourceforge.net/octave/?rev=10572&view=rev Author: prnienhuis Date: 2012-06-06 21:07:54 +0000 (Wed, 06 Jun 2012) Log Message: ----------- Adapted to implementation of "formulas_as_text" option for COM Modified Paths: -------------- trunk/octave-forge/main/io/inst/io_xls_testscript.m trunk/octave-forge/main/io/inst/xls2oct.m Modified: trunk/octave-forge/main/io/inst/io_xls_testscript.m =================================================================== --- trunk/octave-forge/main/io/inst/io_xls_testscript.m 2012-06-06 17:15:51 UTC (rev 10571) +++ trunk/octave-forge/main/io/inst/io_xls_testscript.m 2012-06-06 21:07:54 UTC (rev 10572) @@ -20,6 +20,7 @@ ## Author: Philip Nienhuis ## Created: 2012-02-25 ## Updates: +## 2012-06-06 Adapted to COM implementation for "formulas_as_text" option printf ("\nTesting .xls interface %s ...\n", intf); @@ -93,14 +94,14 @@ endif end_try_catch -## Check if formulas_as_text works (doesn't with COM): +## Check if "formulas_as_text" option works: printf ("\n 8. Repeat reading, now return formulas as text\n"); opts.formulas_as_text = 1; xls = xlsopen ('io-test.xls', 0, intf); raw = xls2oct (xls, shnr, crange, opts); xls = xlsclose (xls); -## 9. Here come the tests, part 2. Fails on COM +## 9. Here come the tests, part 2. printf ("\n 9. Tests part 2 (read back formula):\n"); try Modified: trunk/octave-forge/main/io/inst/xls2oct.m =================================================================== --- trunk/octave-forge/main/io/inst/xls2oct.m 2012-06-06 17:15:51 UTC (rev 10571) +++ trunk/octave-forge/main/io/inst/xls2oct.m 2012-06-06 21:07:54 UTC (rev 10572) @@ -132,8 +132,9 @@ ## 2012-01-26 Fixed "seealso" help string ## 2012-02-25 Fixed missing quotes in struct check L.149-153 ## 2012-02-26 Updated texinfo header help text +## 2012-06-06 Implemented "formulas_as_text" option for COM ## -## Latest subfunc update: 2012-02-25 +## Latest subfunc update: 2012-06-06 function [ rawarr, xls, rstatus ] = xls2oct (xls, wsh=1, datrange='', spsh_opts=[]) @@ -168,9 +169,8 @@ # Select the proper interfaces if (strcmp (xls.xtype, 'COM')) - # Call Excel tru COM server. Excel/COM has no way of returning formulas - # as strings, so arg spsh_opts has no use (yet) - [rawarr, xls, rstatus] = xls2com2oct (xls, wsh, datrange); + # Call Excel tru COM / ActiveX server + [rawarr, xls, rstatus] = xls2com2oct (xls, wsh, datrange, spsh_opts); elseif (strcmp (xls.xtype, 'POI')) # Read xls file tru Java POI [rawarr, xls, rstatus] = xls2jpoi2oct (xls, wsh, datrange, spsh_opts); @@ -237,6 +237,7 @@ ## @deftypefn {Function File} [@var{obj}, @var{rstatus}, @var{xls} ] = xls2com2oct (@var{xls}) ## @deftypefnx {Function File} [@var{obj}, @var{rstatus}, @var{xls} ] = xls2com2oct (@var{xls}, @var{wsh}) ## @deftypefnx {Function File} [@var{obj}, @var{rstatus}, @var{xls} ] = xls2com2oct (@var{xls}, @var{wsh}, @var{range}) +## @deftypefnx {Function File} [@var{obj}, @var{rstatus}, @var{xls} ] = xls2com2oct (@var{xls}, @var{wsh}, @var{range}, @var{spsh_opts}) ## Get cell contents in @var{range} in worksheet @var{wsh} in an Excel ## file pointed to in struct @var{xls} into the cell array @var{obj}. ## @@ -265,8 +266,9 @@ ## 2010-11-12 Moved ptr struct check into main func ## 2010-11-13 Catch empty sheets when no range was specified ## 2012-01-26 Fixed "seealso" help string +## 2012-06-06 Implemented "formulas_as_text option" -function [rawarr, xls, rstatus ] = xls2com2oct (xls, wsh, crange) +function [rawarr, xls, rstatus ] = xls2com2oct (xls, wsh, crange, spsh_opts) rstatus = 0; rawarr = {}; @@ -339,7 +341,11 @@ if (nrows >= 1) # Get object from Excel sheet, starting at cell top_left_cell rr = sh.Range (crange); - rawarr = rr.Value; + if (spsh_opts.formulas_as_text) + rawarr = rr.Formula; + else + rawarr = rr.Value; + endif delete (rr); # Take care of actual singe cell range @@ -405,14 +411,14 @@ ## 2010-07-28 Added option to read formulas as text strings rather than evaluated value ## 2010-08-01 Some bug fixes for formula reading (cvalue rather than scell) ## 2010-10-10 Code cleanup: -getusedrange called; - fixed typo in formula evaluation msg; -## " moved cropping output array to calling function. +## '' moved cropping output array to calling function. ## 2010-11-12 Moved ptr struct check into main func ## 2010-11-13 Catch empty sheets when no range was specified ## 2010-11-14 Fixed sheet # index (was offset by -1) in call to getusedrange() in case ## of text sheet name arg ## 2012-01-26 Fixed "seealso" help string -function [ rawarr, xls, rstatus ] = xls2jpoi2oct (xls, wsh, cellrange=[], spsh_opts) +function [ rawarr, xls, rstatus ] = xls2jpoi2oct (xls, wsh, cellrange, spsh_opts) persistent ctype; if (isempty (ctype)) @@ -584,14 +590,14 @@ ## Added check for proper xls structure ## 2010-07-29 Added check for too latge requested data rectangle ## 2010-10-10 Code cleanup: -getusedrange(); moved cropping result array to -## " calling function +## '' calling function ## 2010-11-12 Moved ptr struct check into main func ## 2010-11-13 Catch empty sheets when no range was specified ## 2011-04-11 (Ron Goldman <ro...@oc...>) Fixed missing months var, wrong arg -## " order in strsplit, wrong isTime condition +## '' order in strsplit, wrong isTime condition ## 2012-01-26 Fixed "seealso" help string -function [ rawarr, xls, rstatus ] = xls2jxla2oct (xls, wsh, cellrange=[], spsh_opts) +function [ rawarr, xls, rstatus ] = xls2jxla2oct (xls, wsh, cellrange, spsh_opts) persistent ctype; persistent months; if (isempty (ctype)) @@ -788,7 +794,7 @@ ## Updates: ## 2012-02-25 Changed ctype into num array rather than cell array -function [ rawarr, xls, rstatus ] = xls2oxs2oct (xls, wsh, cellrange=[], spsh_opts) +function [ rawarr, xls, rstatus ] = xls2oxs2oct (xls, wsh, cellrange, spsh_opts) persistent ctype; if (isempty (ctype)) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-06 17:16:02
|
Revision: 10571 http://octave.svn.sourceforge.net/octave/?rev=10571&view=rev Author: jpicarbajal Date: 2012-06-06 17:15:51 +0000 (Wed, 06 Jun 2012) Log Message: ----------- statistics: adding regression based on gaussian processes. Modified Paths: -------------- trunk/octave-forge/main/statistics/INDEX trunk/octave-forge/main/statistics/NEWS Added Paths: ----------- trunk/octave-forge/main/statistics/inst/regress_gp.m Modified: trunk/octave-forge/main/statistics/INDEX =================================================================== --- trunk/octave-forge/main/statistics/INDEX 2012-06-05 14:20:07 UTC (rev 10570) +++ trunk/octave-forge/main/statistics/INDEX 2012-06-06 17:15:51 UTC (rev 10571) @@ -51,6 +51,7 @@ monotone_smooth princomp regress + regress_gp Plots boxplot normplot Modified: trunk/octave-forge/main/statistics/NEWS =================================================================== --- trunk/octave-forge/main/statistics/NEWS 2012-06-05 14:20:07 UTC (rev 10570) +++ trunk/octave-forge/main/statistics/NEWS 2012-06-06 17:15:51 UTC (rev 10571) @@ -1,3 +1,10 @@ +Summary of important user-visible changes for statistics 1.2.0: +------------------------------------------------------------------- + + ** The following functions are new in 1.2.0: + + regress_gp + Summary of important user-visible changes for statistics 1.1.3: ------------------------------------------------------------------- Added: trunk/octave-forge/main/statistics/inst/regress_gp.m =================================================================== --- trunk/octave-forge/main/statistics/inst/regress_gp.m (rev 0) +++ trunk/octave-forge/main/statistics/inst/regress_gp.m 2012-06-06 17:15:51 UTC (rev 10571) @@ -0,0 +1,123 @@ +## Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +## +## This program is free software; you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or +## (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/>. + +## -*- texinfo -*- +## @deftypefn {Function File} {[@var{m}, @var{K}] =} regress_gp (@var{x}, @var{y}, @var{Sp}) +## @deftypefnx {Function File} {[@dots{} @var{yi} @var{dy}] =} sqp (@dots{}, @var{xi}) +## Linear scalar regression using gaussian processes. +## +## It estimates the model @var{y} = @var{x}'*m for @var{x} R^D and @var{y} in R. +## The information about errors of the predictions (interpolation/extrapolation) is given +## by the covarianve matrix @var{K}. If D==1 the inputs must be column vectors, +## if D>1 then @var{x} is n-by-D, with n the number of data points. @var{Sp} defines +## the prior covariance of @var{m}, it should be a (D+1)-by-(D+1) positive definite matrix, +## if it is empty, the default is @code{Sp = 100*eye(size(x,2)+1)}. +## +## If @var{xi} inputs are provided, the model is evaluated and returned in @var{yi}. +## The estimation of the variation of @var{yi} are given in @var{dy}. +## +## Run @code{demo regress_gp} to see an examples. +## +## The function is a direc implementation of the formulae in pages 11-12 of +## Gaussian Processes for Machine Learning. Carl Edward Rasmussen and @ +## Christopher K. I. Williams. The MIT Press, 2006. ISBN 0-262-18253-X. +## available online at @url{http://gaussianprocess.org/gpml/}. +## +## @seealso{regress} +## @end deftypefn + +function [wm K yi dy] = regress_gp (x,y,Sp=[],xi=[]) + + if isempty(Sp) + Sp = 100*eye(size(x,2)+1); + end + + x = [ones(1,size(x,1)); x']; + + A = x*x' + inv (Sp); + K = inv (A); + wm = K*x*y; + + yi =[]; + dy =[]; + if !isempty (xi); + xi = [ones(size(xi,1),1) xi]; + yi = xi*wm; + dy = diag (xi*K*xi'); + end + +endfunction + +%!demo +%! % 1D Data +%! x = 2*rand (5,1)-1; +%! y = 2*x -1 + 0.3*randn (5,1); +%! +%! % Points for interpolation/extrapolation +%! xi = linspace (-2,2,10)'; +%! +%! [m K yi dy] = regress_gp (x,y,[],xi); +%! +%! plot (x,y,'xk',xi,yi,'r-',xi,yi+[-dy +dy],'b-'); + +%!demo +%! % 2D Data +%! x = 2*rand (4,2)-1; +%! y = 2*x(:,1)-3*x(:,2) -1 + 1*randn (4,1); +%! +%! % Mesh for interpolation/extrapolation +%! [xi yi] = meshgrid (linspace (-1,1,10)); +%! +%! [m K zi dz] = regress_gp (x,y,[],[xi(:) yi(:)]); +%! zi = reshape (zi, 10,10); +%! dz = reshape (dz,10,10); +%! +%! plot3 (x(:,1),x(:,2),y,'.g','markersize',8); +%! hold on; +%! h = mesh (xi,yi,zi,zeros(10,10)); +%! set(h,'facecolor','none'); +%! h = mesh (xi,yi,zi+dz,ones(10,10)); +%! set(h,'facecolor','none'); +%! h = mesh (xi,yi,zi-dz,ones(10,10)); +%! set(h,'facecolor','none'); +%! hold off +%! axis tight +%! view(80,25) + +%!demo +%! % Projection over basis function +%! pp = [2 2 0.3 1]; +%! n = 10; +%! x = 2*rand (n,1)-1; +%! y = polyval(pp,x) + 0.3*randn (n,1); +%! +%! % Powers +%! px = [sqrt(abs(x)) x x.^2 x.^3]; +%! +%! % Points for interpolation/extrapolation +%! xi = linspace (-1,1,100)'; +%! pxi = [sqrt(abs(xi)) xi xi.^2 xi.^3]; +%! +%! Sp = 100*eye(size(px,2)+1); +%! Sp(2,2) = 1; # We don't believe the sqrt is present +%! [m K yi dy] = regress_gp (px,y,Sp,pxi); +%! disp(m) +%! +%! plot (x,y,'xk;Data;',xi,yi,'r-;Estimation;',xi,polyval(pp,xi),'g-;True;'); +%! axis tight +%! axis manual +%! hold on +%! plot (xi,yi+[-dy +dy],'b-'); +%! hold off This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 15:27:30
|
Revision: 10570 http://octave.svn.sourceforge.net/octave/?rev=10570&view=rev Author: jpicarbajal Date: 2012-06-05 14:20:07 +0000 (Tue, 05 Jun 2012) Log Message: ----------- geometry: fixing demos in curve2polyline. Modified Paths: -------------- trunk/octave-forge/main/geometry/inst/shape2d/curve2polyline.m trunk/octave-forge/main/geometry/inst/shape2d/shapetransform.m Modified: trunk/octave-forge/main/geometry/inst/shape2d/curve2polyline.m =================================================================== --- trunk/octave-forge/main/geometry/inst/shape2d/curve2polyline.m 2012-06-05 13:52:10 UTC (rev 10569) +++ trunk/octave-forge/main/geometry/inst/shape2d/curve2polyline.m 2012-06-05 14:20:07 UTC (rev 10570) @@ -138,4 +138,4 @@ %! t = linspace(0,1,100)'; %! pc = curveval(curve,t); %! -%! plot(p(:,1),p(:,2),'-o',pc(:,1),pc(:,2),'-r') +%! plot(polyline(:,1),polyline(:,2),'-o',pc(:,1),pc(:,2),'-r') Modified: trunk/octave-forge/main/geometry/inst/shape2d/shapetransform.m =================================================================== --- trunk/octave-forge/main/geometry/inst/shape2d/shapetransform.m 2012-06-05 13:52:10 UTC (rev 10569) +++ trunk/octave-forge/main/geometry/inst/shape2d/shapetransform.m 2012-06-05 14:20:07 UTC (rev 10570) @@ -98,7 +98,7 @@ %! shape = shapetransform (shape,-T + [2; 0]); %! %! close -%! shapeplot (shape,'-r','linewidth',2) +%! shapeplot (shape,'-r','linewidth',2); %! hold on %! for i = 1:9 %! T = createRotation (i*pi/5)(1:2,1:2)/exp(0.3*i); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 13:52:19
|
Revision: 10569 http://octave.svn.sourceforge.net/octave/?rev=10569&view=rev Author: jpicarbajal Date: 2012-06-05 13:52:10 +0000 (Tue, 05 Jun 2012) Log Message: ----------- geometry: conversion form ellipse to covariance matrices. Modified Paths: -------------- trunk/octave-forge/main/geometry/DESCRIPTION trunk/octave-forge/main/geometry/NEWS trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m Modified: trunk/octave-forge/main/geometry/DESCRIPTION =================================================================== --- trunk/octave-forge/main/geometry/DESCRIPTION 2012-06-05 13:39:37 UTC (rev 10568) +++ trunk/octave-forge/main/geometry/DESCRIPTION 2012-06-05 13:52:10 UTC (rev 10569) @@ -1,6 +1,6 @@ Name: Geometry Version: 1.5.0 -Date: 2012-04-16 +Date: 2012-06-05 Author: David Legland <dav...@gr...>, José Luis García Pallero <jgp...@gm...>, Juan Pablo Carbajal <car...@if...>, Simeon Simeonov <ss...@ca...> Maintainer: Juan Pablo Carbajal <car...@if...> Title: Computational Geometry Modified: trunk/octave-forge/main/geometry/NEWS =================================================================== --- trunk/octave-forge/main/geometry/NEWS 2012-06-05 13:39:37 UTC (rev 10568) +++ trunk/octave-forge/main/geometry/NEWS 2012-06-05 13:52:10 UTC (rev 10569) @@ -1,7 +1,7 @@ Summary of important user-visible changes for releases of the geometry package =============================================================================== -geometry-1.5.0 Release Date: 2012-04-16 Release Manager: Juan Pablo Carbajal +geometry-1.5.0 Release Date: 2012-06-05 Release Manager: Juan Pablo Carbajal =============================================================================== * Added functions: Modified: trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m =================================================================== --- trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m 2012-06-05 13:39:37 UTC (rev 10568) +++ trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m 2012-06-05 13:52:10 UTC (rev 10569) @@ -20,7 +20,7 @@ %% Calculates covariance matrix from ellipse. %% %% If only one input is given, @var{elli} must define an ellipse as described in -%% @coomand{ellipses2d}. +%% @command{ellipses2d}. %% If two inputs are given, @var{ra} and @var{rb} define the half-lenght of the %% axes. %% If a third input is given, @var{theta} must be the angle of rotation of the This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 13:39:48
|
Revision: 10568 http://octave.svn.sourceforge.net/octave/?rev=10568&view=rev Author: jpicarbajal Date: 2012-06-05 13:39:37 +0000 (Tue, 05 Jun 2012) Log Message: ----------- geometry: conversion form ellipse to covariance matrices. Modified Paths: -------------- trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m Modified: trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m =================================================================== --- trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m 2012-06-05 12:26:47 UTC (rev 10567) +++ trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m 2012-06-05 13:39:37 UTC (rev 10568) @@ -19,13 +19,29 @@ %% @deftypefnx {Function File} {@dots{} = } cov2ellipse (@dots{}, @samp{tol},@var{tol}) %% Calculates ellipse parameters from covariance matrix. %% +%% @var{K} must be symmetric positive (semi)definite. The optional argument +%% @samp{tol} sets the tolerance for the verification of the +%% positive-(semi)definiteness of the matrix @var{K} (see @command{isdefinite}). +%% +%% If only one output argument is supplied a vector defining a ellipse is returned +%% as defined in @command{ellipses2d}. Otherwise the angle @var{theta} is given +%% in radians. +%% %% Run @code{demo cov2ellipse} to see an example. -%% +%% %% @seealso{ellipses2d, cov2ellipse, drawEllipse} %% @end deftypefn function varargout = cov2ellipse (K, varargin); + parser = inputParser (); + parser.FunctionName = "cov2ellipse"; + parser = addParamValue (parser,'Tol', 100*eps*norm (K, "fro"), @(x)x>0); + parser = parse(parser,varargin{:}); + + if isdefinite (K,parser.Results.Tol) == -1 + print_usage + end [R S W] = svd (K); theta = atan (R(1,1)/R(2,2)); v = sort (diag(S), 'ascend')'; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 12:26:58
|
Revision: 10567 http://octave.svn.sourceforge.net/octave/?rev=10567&view=rev Author: jpicarbajal Date: 2012-06-05 12:26:47 +0000 (Tue, 05 Jun 2012) Log Message: ----------- geometry: conversion form ellipse to covariance matrices. Modified Paths: -------------- trunk/octave-forge/main/geometry/INDEX trunk/octave-forge/main/geometry/NEWS Modified: trunk/octave-forge/main/geometry/INDEX =================================================================== --- trunk/octave-forge/main/geometry/INDEX 2012-06-05 12:23:02 UTC (rev 10566) +++ trunk/octave-forge/main/geometry/INDEX 2012-06-05 12:26:47 UTC (rev 10567) @@ -98,9 +98,11 @@ intersectLineCircle radicalAxis 2D Ellipses - ellipseAsPolygon + cov2ellipse drawEllipseArc drawEllipse + ellipseAsPolygon + ellipse2cov inertiaEllipse 2D Transformations transformLine Modified: trunk/octave-forge/main/geometry/NEWS =================================================================== --- trunk/octave-forge/main/geometry/NEWS 2012-06-05 12:23:02 UTC (rev 10566) +++ trunk/octave-forge/main/geometry/NEWS 2012-06-05 12:26:47 UTC (rev 10567) @@ -5,6 +5,8 @@ =============================================================================== * Added functions: + - cov2ellipse.m & ellipse2cov: transform between ellipses and covariances matrices. + - beltproblem.m : Finds the four lines tangent to two circles with given centers and radii. This is the solution to the belt problem in 2D. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 12:23:09
|
Revision: 10566 http://octave.svn.sourceforge.net/octave/?rev=10566&view=rev Author: jpicarbajal Date: 2012-06-05 12:23:02 +0000 (Tue, 05 Jun 2012) Log Message: ----------- geometry: conversion form ellipse to covariance matrices. Added Paths: ----------- trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m Added: trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m =================================================================== --- trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m (rev 0) +++ trunk/octave-forge/main/geometry/inst/geom2d/cov2ellipse.m 2012-06-05 12:23:02 UTC (rev 10566) @@ -0,0 +1,55 @@ +## Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +## +## This program is free software; you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or +## (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/>. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {@var{ellipse} = } cov2ellipse (@var{K}) +%% @deftypefnx {Function File} {[@var{ra} @var{rb} @var{theta}] = } cov2ellipse (@var{K}) +%% @deftypefnx {Function File} {@dots{} = } cov2ellipse (@dots{}, @samp{tol},@var{tol}) +%% Calculates ellipse parameters from covariance matrix. +%% +%% Run @code{demo cov2ellipse} to see an example. +%% +%% @seealso{ellipses2d, cov2ellipse, drawEllipse} +%% @end deftypefn + +function varargout = cov2ellipse (K, varargin); + + [R S W] = svd (K); + theta = atan (R(1,1)/R(2,2)); + v = sort (diag(S), 'ascend')'; + + if nargout == 1 + varargout{1} = [0 0 v theta*180/pi]; + elseif nargout == 3 + varargout{1} = v(1); + varargout{2} = v(2); + varargout{3} = theta; + end + +endfunction + +%!demo +%! K = [2 1; 1 2]; +%! L = chol(K,'lower'); +%! u = randn(1e3,2)*L'; +%! +%! elli = cov2ellipse (K) +%! +%! figure(1) +%! plot(u(:,1),u(:,2),'.r'); +%! hold on; +%! drawEllipse(elli,'linewidth',2); +%! hold off +%! axis tight Added: trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m =================================================================== --- trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m (rev 0) +++ trunk/octave-forge/main/geometry/inst/geom2d/ellipse2cov.m 2012-06-05 12:23:02 UTC (rev 10566) @@ -0,0 +1,92 @@ +## Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +## +## This program is free software; you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or +## (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/>. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {@var{K} = } ellipse2cov (@var{elli}) +%% @deftypefnx {Function File} {@var{K} = } ellipse2cov (@var{ra}, @var{rb}) +%% @deftypefnx {Function File} {@var{K} = } ellipse2cov (@dots{}, @var{theta}) +%% Calculates covariance matrix from ellipse. +%% +%% If only one input is given, @var{elli} must define an ellipse as described in +%% @coomand{ellipses2d}. +%% If two inputs are given, @var{ra} and @var{rb} define the half-lenght of the +%% axes. +%% If a third input is given, @var{theta} must be the angle of rotation of the +%% ellipse in radians, and in counter-clockwise direction. +%% +%% The output @var{K} contains the covariance matrix define by the ellipse. +%% +%% Run @code{demo ellipse2cov} to see an example. +%% +%% @seealso{ellipses2d, cov2ellipse, drawEllipse} +%% @end deftypefn + +function K = ellipse2cov (elli, varargin); + + ra = 1; + rb = 1; + theta = 0; + switch numel (varargin) + case 0 + %% ellipse format + if numel (elli) != 5 + print_usage (); + end + ra = elli(1,3); + rb = elli(1,4); + theta = elli(1,5)*pi/180; + + case 2 + %% ra,rb + if numel (elli) != 1 + print_usage (); + end + ra = elli; + rb = varargin{1}; + + case 3 + %% ra,rb, theta + if numel (elli) != 1 + print_usage (); + end + ra = elli; + rb = varargin{1}; + theta = varargin{2}; + + otherwise + print_usage (); + end + + T = createRotation (theta)(1:2,1:2); + K = T*diag([ra rb])*T'; + +endfunction + +%!demo +%! elli = [0 0 1 3 -45]; +%! +%! % Create 2D normal random variables with covarinace defined by elli. +%! K = ellipse2cov (elli) +%! L = chol(K,'lower'); +%! u = randn(1e3,2)*L'; +%! +%! Kn = cov (u) +%! +%! figure(1) +%! plot(u(:,1),u(:,2),'.r'); +%! hold on; +%! drawEllipse(elli,'linewidth',2); +%! hold off +%! axis tight This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 07:16:41
|
Revision: 10565 http://octave.svn.sourceforge.net/octave/?rev=10565&view=rev Author: jpicarbajal Date: 2012-06-05 07:16:31 +0000 (Tue, 05 Jun 2012) Log Message: ----------- system-identification: adding functions to TISEAN wrapper Modified Paths: -------------- trunk/octave-forge/main/system-identification/PKG_ADD trunk/octave-forge/main/system-identification/PKG_DEL trunk/octave-forge/main/system-identification/pre_install.m Added Paths: ----------- trunk/octave-forge/main/system-identification/inst/tisean_wrapper/ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m trunk/octave-forge/main/system-identification/inst/tisean_wrapper/falseneigh.m trunk/octave-forge/main/system-identification/inst/tisean_wrapper/mutualinf.m Removed Paths: ------------- trunk/octave-forge/main/system-identification/inst/tisean/ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m Modified: trunk/octave-forge/main/system-identification/PKG_ADD =================================================================== --- trunk/octave-forge/main/system-identification/PKG_ADD 2012-06-05 07:10:14 UTC (rev 10564) +++ trunk/octave-forge/main/system-identification/PKG_ADD 2012-06-05 07:16:31 UTC (rev 10565) @@ -1,5 +1,5 @@ %1 -dirlist = {"tisean"}; +dirlist = {"tisean_wrapper"}; dirname = fileparts (canonicalize_file_name (mfilename ("fullpath"))); pp = strsplit (dirname,filesep (), true); Modified: trunk/octave-forge/main/system-identification/PKG_DEL =================================================================== --- trunk/octave-forge/main/system-identification/PKG_DEL 2012-06-05 07:10:14 UTC (rev 10564) +++ trunk/octave-forge/main/system-identification/PKG_DEL 2012-06-05 07:16:31 UTC (rev 10565) @@ -1,5 +1,5 @@ %1 -dirlist = {"tisean"}; +dirlist = {"tisean_wrapper"}; dirname = fileparts (canonicalize_file_name (mfilename ("fullpath"))); pp = strsplit (dirname,filesep (), true); Deleted: trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m =================================================================== --- trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m 2012-06-04 11:41:32 UTC (rev 10561) +++ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m 2012-06-05 07:16:31 UTC (rev 10565) @@ -1,71 +0,0 @@ -%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> -%% -%% This program is free software: you can redistribute it and/or modify -%% it under the terms of the GNU General Public License as published by -%% the Free Software Foundation, either version 3 of the License, or -%% any later version. -%% -%% This program is distributed in the hope that it will be useful, -%% but WITHOUT ANY WARRANTY; without even the implied warranty of -%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -%% GNU General Public License for more details. -%% -%% You should have received a copy of the GNU General Public License -%% along with this program. If not, see <http://www.gnu.org/licenses/>. - -%% -*- texinfo -*- -%% @deftypefn {Function File} {[@var{c} @var{d} @var{h} @var{stat}] = } corrdim (@var{data}, @var{edim}) -%% Correlation dimension from @var{data}. This function calls @code{d2} @ -%% from the TISEAN package. -%% -%% @end deftypefn - -function [c d h] = corrdim (data, varargin) - - - [nT M] = size (data); - # --- Parse arguments --- # - parser = inputParser (); - parser.FunctionName = "corrdim"; - parser = addParamValue (parser,'Components', 1, @(x)x>0); - parser = addParamValue (parser,'MaxEdim', 10, @(x)x>0); - parser = addParamValue (parser,'Delay', 1, @(x)x>0); - parser = addParamValue (parser,'TheilerWindow', 0, @(x)x>=0); - parser = addParamValue (parser,'ScaleSpan', nT*[1e-3 1], @(x)all(x>0)); - parser = addParamValue (parser,'EpsilonCount', 100, @(x)x>0); - parser = addParamValue (parser,'PairCount', 1000, @(x)x>=0); - parser = addParamValue (parser,'Normalize', false); - parser = addParamValue (parser,'Verbose', false); - parser = parse(parser,varargin{:}); - - flag.E = ""; - if parser.Results.Normalize - flag.E = "-E"; - end - - infile = tmpnam (); - outfile = tmpnam (); - - %% Write data to file - save ('-ascii',infile, 'data'); - - %% Prepare format of the embedding vector - syscmd = sprintf ("d2 -d%d -M%d,%d -t%d -r%d -R%d -#%d -N%d %s -o%s -V0 %s", ... - parser.Results.Delay, ... - parser.Results.Components, ... - parser.Results.MaxEdim, ... - parser.Results.TheilerWindow, ... - parser.Results.ScaleSpan, ... - parser.Results.EpsilonCount, ... - parser.Results.PairCount, ... - flag.E, outfile, infile); - - %% Function call - system (syscmd); - - c = load ([outfile ".c2"]); - d = load ([outfile ".d2"]); - h = load ([outfile ".h2"]); - - -endfunction Copied: trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m (from rev 10564, trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m) =================================================================== --- trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m (rev 0) +++ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/corrdim.m 2012-06-05 07:16:31 UTC (rev 10565) @@ -0,0 +1,76 @@ +%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +%% +%% This program is free software: you can redistribute it and/or modify +%% it under the terms of the GNU General Public License as published by +%% the Free Software Foundation, either version 3 of the License, or +%% any later version. +%% +%% This program is distributed in the hope that it will be useful, +%% but WITHOUT ANY WARRANTY; without even the implied warranty of +%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +%% GNU General Public License for more details. +%% +%% You should have received a copy of the GNU General Public License +%% along with this program. If not, see <http://www.gnu.org/licenses/>. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {[@var{c} @var{d} @var{h}] = } corrdim (@var{data}) +%% @deftypefnx {Function File} {[@dots{}] = } corrdim (@dots{},@var{property},@var{value}) +%% Correlation dimension from @var{data}. This function calls @code{d2} @ +%% from the TISEAN package. +%% +%% @end deftypefn + +function [c d h] = corrdim (data, varargin) + + + [nT M] = size (data); + amp = max(max(data)-min(data)); + # --- Parse arguments --- # + parser = inputParser (); + parser.FunctionName = "corrdim"; + parser = addParamValue (parser,'Components', 1, @(x)x>0); + parser = addParamValue (parser,'MaxEdim', 10, @(x)x>0); + parser = addParamValue (parser,'Delay', 1, @(x)x>0); + parser = addParamValue (parser,'TheilerWindow', 0, @(x)x>=0); + parser = addParamValue (parser,'ScaleSpan', amp*[1e-3 1], @(x)all(x>0)); + parser = addParamValue (parser,'EpsilonCount', 100, @(x)x>0); + parser = addParamValue (parser,'PairCount', 1000, @(x)x>=0); + parser = addSwitch (parser,'Normalize'); + parser = addSwitch (parser,'Verbose'); + parser = parse(parser,varargin{:}); + + flag.E = ""; + if parser.Results.Normalize + flag.E = "-E"; + end + + infile = tmpnam (); + outfile = tmpnam (); + + %% Write data to file + save ('-ascii',infile, 'data'); + + %% Prepare format of system call + func = file_in_loadpath ("d2"); + + %% Prepare format of the embedding vector + syscmd = sprintf ("%s -d%d -M%d,%d -t%d -r%f -R%f -#%d -N%d %s -o%s -V0 %s", func, ... + parser.Results.Delay, ... + parser.Results.Components, ... + parser.Results.MaxEdim, ... + parser.Results.TheilerWindow, ... + parser.Results.ScaleSpan, ... + parser.Results.EpsilonCount, ... + parser.Results.PairCount, ... + flag.E, outfile, infile); + + %% Function call + system (syscmd); + + c = load ([outfile ".c2"]); + d = load ([outfile ".d2"]); + h = load ([outfile ".h2"]); + + +endfunction Added: trunk/octave-forge/main/system-identification/inst/tisean_wrapper/falseneigh.m =================================================================== --- trunk/octave-forge/main/system-identification/inst/tisean_wrapper/falseneigh.m (rev 0) +++ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/falseneigh.m 2012-06-05 07:16:31 UTC (rev 10565) @@ -0,0 +1,98 @@ +%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +%% +%% This program is free software: you can redistribute it and/or modify +%% it under the terms of the GNU General Public License as published by +%% the Free Software Foundation, either version 3 of the License, or +%% any later version. +%% +%% This program is distributed in the hope that it will be useful, +%% but WITHOUT ANY WARRANTY; without even the implied warranty of +%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +%% GNU General Public License for more details. +%% +%% You should have received a copy of the GNU General Public License +%% along with this program. If not, see <http://www.gnu.org/licenses/>. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {[@var{dim} @var{frac} @var{avgsize} @var{avg2size}] = } falseneigh (@var{data}) +%% @deftypefnx {Function File} {[@dots{}] = } pspec_mep (@dots{},@var{property},@var{value}) +%% Looks for the nearest neighbors of all data points in m dimensions and @ +%% iterates these neighbors one step (more precisely delay steps) into the future. @ +%% If the ratio of the distance of the iteration and that of the nearest neighbor @ +%% exceeds a given threshold the point is marked as a wrong neighbor. The output is @ +%% the fraction of false neighbors for the specified embedding dimensions. @ +%% This function calls @code{false_nearest} from the TISEAN package. +%% +%% @var{dim} is the dimension. +%% @var{frac} is the fraction of false nearest neighbors. +%% @var{avgsize} is the average size of the neighborhood. +%% @var{avg2size} is the average of the squared size of the neighborhood. +%% +%% Optional arguments are: +%% @table @samp +%% @item MinEdim +%% Minimal embedding dimensions of the vectors. Deafault 1. +%% @item MaxEdim +%% Maximal embedding dimensions of the vectors. Defualt 5. +%% @item Delay +%% Delay of the vectors. Default 1. +%% @item Ratio +%% Ratio factor. Default 2.0. +%% @item Theiler +%% Theiler window. Default 0. +%% @end table +%% +%% @end deftypefn + +function [dim frac avgsize avg2size] = falseneigh (data, varargin) + + # --- Parse arguments --- # + parser = inputParser (); + parser.FunctionName = "falseneigh"; + parser = addParamValue (parser,'MinEdim', 1, @(x)x>0); + parser = addParamValue (parser,'MaxEdim', 5, @(x)x>0); + parser = addParamValue (parser,'Delay', 1, @(x)x>0); + parser = addParamValue (parser,'Ratio', 2, @(x)x>0); + parser = addParamValue (parser,'Theiler', 0, @(x)x>=0); + parser = parse(parser,varargin{:}); + + medim = parser.Results.MinEdim; + Medim = parser.Results.MaxEdim; + delay = parser.Results.Delay; + ratio = parser.Results.Ratio; + thei = parser.Results.Theiler; + + clear parser + # ------ # + + [nT n] = size (data); + + flag.m = sprintf("-m%d", medim); + flag.M = sprintf("-M%d,%d", n,Medim); + flag.d = sprintf("-d%d", delay); + flag.f = sprintf("-f%f", ratio); + flag.t = sprintf("-t%d", thei); + + infile = tmpnam (); + outfile = tmpnam (); + + %% Write data to file + save ('-ascii',infile, 'data'); + + %% Prepare format of system call + func = file_in_loadpath ("false_nearest"); + syscmd = sprintf("%s %s %s %s %s %s -o%s %s", ... + func, flag.m, flag.M, flag.d, flag.f, flag.t, outfile, infile); + + %% Function call + system (syscmd); + + %% Read results + s = load (outfile); + + dim = s(:,1); % the dimension (counted like shown above) + frac = s(:,2); % the fraction of false nearest neighbors + avgsize = s(:,3); % the average size of the neighborhood + avg2size = s(:,4); % the average of the squared size of the neighborhood + +endfunction Added: trunk/octave-forge/main/system-identification/inst/tisean_wrapper/mutualinf.m =================================================================== --- trunk/octave-forge/main/system-identification/inst/tisean_wrapper/mutualinf.m (rev 0) +++ trunk/octave-forge/main/system-identification/inst/tisean_wrapper/mutualinf.m 2012-06-05 07:16:31 UTC (rev 10565) @@ -0,0 +1,81 @@ +%% Copyright (c) 2012 Juan Pablo Carbajal <car...@if...> +%% +%% This program is free software: you can redistribute it and/or modify +%% it under the terms of the GNU General Public License as published by +%% the Free Software Foundation, either version 3 of the License, or +%% any later version. +%% +%% This program is distributed in the hope that it will be useful, +%% but WITHOUT ANY WARRANTY; without even the implied warranty of +%% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +%% GNU General Public License for more details. +%% +%% You should have received a copy of the GNU General Public License +%% along with this program. If not, see <http://www.gnu.org/licenses/>. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {[@var{d} @var{minf} @var{shinf}] = } mutualinf (@var{data}) +%% @deftypefnx {Function File} {@dots{} = } nststz (@dots{},@var{property},@var{value}) +%% Estimates the time delayed mutual information of the data. It is the simplest +%% possible realization. It uses a fixed mesh of boxes. This function calls +%% @code{mutual} from the TISEAN package. +%% +%% @var{shinf} contains the Shannon entropy (normalized to the number of +%% occupied boxes). @var{d} contains the delays and @var{minf} the mutual information. +%% +%% Optional arguments are: +%% @table @samp +%% @item Nboxes +%% Number of boxes for the partition. Default 16. +%% @item MaxDelay +%% Maximal time delay. Default 20. +%% @end table +%% +%% @end deftypefn + +function [d minf shinf] = mutualinf (data, varargin) + + # --- Parse arguments --- # + parser = inputParser (); + parser.FunctionName = "nstatz"; + parser = addParamValue (parser,'Nboxes', 16, @(x)x>0); + parser = addParamValue (parser,'MaxDelay', 20, @(x)x>0); + parser = parse(parser,varargin{:}); + + nbox = parser.Results.Nboxes; + delay = parser.Results.MaxDelay; + clear parser + # ------ # + + flag.b = sprintf("-b%d", nbox); + flag.D = sprintf("-D%d", delay); + + infile = tmpnam (); + outfile = tmpnam (); + + [nT nc] = size (data); + shinf = zeros (1, nc); + minf = zeros (delay, nc); + + for i =1:nc + datac = data(:,i); + + %% Write data to file + save ('-ascii',infile, 'datac'); + + %% Prepare format of system call + func = file_in_loadpath ("mutual"); + syscmd = sprintf("%s %s %s -o%s -V0 %s", func, ... + flag.b, flag.D, outfile, infile); + + %% Function call + system (syscmd); + + [inftxt tmp] = textread (outfile, "%s%f", 1); + shinf(i) = tmp;%str2num (strsplit (inftxt, "="){2}) + + [d, tmp] = textread (outfile, "%f%f", delay, "headerlines", 1); + minf(:,i) = tmp; + end + +endfunction Modified: trunk/octave-forge/main/system-identification/pre_install.m =================================================================== --- trunk/octave-forge/main/system-identification/pre_install.m 2012-06-05 07:10:14 UTC (rev 10564) +++ trunk/octave-forge/main/system-identification/pre_install.m 2012-06-05 07:16:31 UTC (rev 10565) @@ -2,7 +2,7 @@ %% Prepares for installation a package that is organized in subfolders %% List of folders with src subfolder - subfld = {"tisean"}; + subfld = {"tisean_wrapper"}; %% Create correct strings subfld_ready = strcat ({[pwd() filesep() "inst" filesep()]}, This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-05 07:10:21
|
Revision: 10564 http://octave.svn.sourceforge.net/octave/?rev=10564&view=rev Author: jpicarbajal Date: 2012-06-05 07:10:14 +0000 (Tue, 05 Jun 2012) Log Message: ----------- system-identification: adding functions to TISEAN wrapper Modified Paths: -------------- trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m Modified: trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m =================================================================== --- trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m 2012-06-05 06:48:07 UTC (rev 10563) +++ trunk/octave-forge/main/system-identification/inst/tisean/corrdim.m 2012-06-05 07:10:14 UTC (rev 10564) @@ -14,7 +14,8 @@ %% along with this program. If not, see <http://www.gnu.org/licenses/>. %% -*- texinfo -*- -%% @deftypefn {Function File} {[@var{c} @var{d} @var{h} @var{stat}] = } corrdim (@var{data}, @var{edim}) +%% @deftypefn {Function File} {[@var{c} @var{d} @var{h}] = } corrdim (@var{data}) +%% @deftypefnx {Function File} {[@dots{}] = } corrdim (@dots{},@var{property},@var{value}) %% Correlation dimension from @var{data}. This function calls @code{d2} @ %% from the TISEAN package. %% @@ -24,6 +25,7 @@ [nT M] = size (data); + amp = max(max(data)-min(data)); # --- Parse arguments --- # parser = inputParser (); parser.FunctionName = "corrdim"; @@ -31,11 +33,11 @@ parser = addParamValue (parser,'MaxEdim', 10, @(x)x>0); parser = addParamValue (parser,'Delay', 1, @(x)x>0); parser = addParamValue (parser,'TheilerWindow', 0, @(x)x>=0); - parser = addParamValue (parser,'ScaleSpan', nT*[1e-3 1], @(x)all(x>0)); + parser = addParamValue (parser,'ScaleSpan', amp*[1e-3 1], @(x)all(x>0)); parser = addParamValue (parser,'EpsilonCount', 100, @(x)x>0); parser = addParamValue (parser,'PairCount', 1000, @(x)x>=0); - parser = addParamValue (parser,'Normalize', false); - parser = addParamValue (parser,'Verbose', false); + parser = addSwitch (parser,'Normalize'); + parser = addSwitch (parser,'Verbose'); parser = parse(parser,varargin{:}); flag.E = ""; @@ -49,8 +51,11 @@ %% Write data to file save ('-ascii',infile, 'data'); + %% Prepare format of system call + func = file_in_loadpath ("d2"); + %% Prepare format of the embedding vector - syscmd = sprintf ("d2 -d%d -M%d,%d -t%d -r%d -R%d -#%d -N%d %s -o%s -V0 %s", ... + syscmd = sprintf ("%s -d%d -M%d,%d -t%d -r%f -R%f -#%d -N%d %s -o%s -V0 %s", func, ... parser.Results.Delay, ... parser.Results.Components, ... parser.Results.MaxEdim, ... This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-05 06:48:13
|
Revision: 10563 http://octave.svn.sourceforge.net/octave/?rev=10563&view=rev Author: paramaniac Date: 2012-06-05 06:48:07 +0000 (Tue, 05 Jun 2012) Log Message: ----------- control: touch up lqe doc Modified Paths: -------------- trunk/octave-forge/main/control/devel/dlqe.m trunk/octave-forge/main/control/devel/lqe.m Modified: trunk/octave-forge/main/control/devel/dlqe.m =================================================================== --- trunk/octave-forge/main/control/devel/dlqe.m 2012-06-04 20:12:13 UTC (rev 10562) +++ trunk/octave-forge/main/control/devel/dlqe.m 2012-06-05 06:48:07 UTC (rev 10563) @@ -16,13 +16,11 @@ ## along with LTI Syncope. If not, see <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{sys}, @var{q}, @var{r}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{sys}, @var{q}, @var{r}, @var{s}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{a}, @var{b}, @var{q}, @var{r}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{a}, @var{b}, @var{q}, @var{r}, @var{s}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{a}, @var{b}, @var{q}, @var{r}, @var{[]}, @var{e}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqr (@var{a}, @var{b}, @var{q}, @var{r}, @var{s}, @var{e}) -## Linear-quadratic regulator. +## @deftypefn {Function File} {[@var{l}, @var{p}, @var{z}, @var{e}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{z}, @var{e}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}, @var{s}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{z}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{z}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}, @var{s}) +## Linear-quadratic estimator for discrete-time systems. ## ## @strong{Inputs} ## @table @var @@ -78,7 +76,9 @@ print_usage (); endif - if (isempty (s)) + if (isempty (g)) + [~, p, e] = dlqr (a.', c.', q, r, s); # dlqe (a, [], c, q, r, s), g=I + elseif (isempty (s)) [~, p, e] = dlqr (a.', c.', g*q*g.', r); else [~, p, e] = dlqr (a.', c.', g*q*g.', r, g*s); Modified: trunk/octave-forge/main/control/devel/lqe.m =================================================================== --- trunk/octave-forge/main/control/devel/lqe.m 2012-06-04 20:12:13 UTC (rev 10562) +++ trunk/octave-forge/main/control/devel/lqe.m 2012-06-05 06:48:07 UTC (rev 10563) @@ -16,12 +16,12 @@ ## along with LTI Syncope. If not, see <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{sys}, @var{q}, @var{r}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{sys}, @var{q}, @var{r}, @var{s}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}, @var{s}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}) -## @deftypefnx {Function File} {[@var{g}, @var{x}, @var{l}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}, @var{s}) +## @deftypefn {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{sys}, @var{q}, @var{r}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{sys}, @var{q}, @var{r}, @var{s}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{g}, @var{c}, @var{q}, @var{r}, @var{s}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}) +## @deftypefnx {Function File} {[@var{l}, @var{p}, @var{e}] =} lqe (@var{a}, @var{[]}, @var{c}, @var{q}, @var{r}, @var{s}) ## Linear-quadratic estimator. ## ## @strong{Inputs} @@ -79,7 +79,7 @@ endif if (isa (a, "lti")) - [l, p, e] = lqr (a.', g, c, q); # lqe (sys.', q, r, s), g=I, works like lqr (sys.', q, r, s).' + [l, p, e] = lqr (a.', g, c, q); # lqe (sys, q, r, s), g=I, works like lqr (sys.', q, r, s).' elseif (isempty (g)) [l, p, e] = lqr (a.', c.', q, r, s); # lqe (a, [], c, q, r, s), g=I, works like lqr (a.', c.', q, r, s).' elseif (isempty (s)) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <be...@us...> - 2012-06-04 20:12:20
|
Revision: 10562 http://octave.svn.sourceforge.net/octave/?rev=10562&view=rev Author: benjf5 Date: 2012-06-04 20:12:13 +0000 (Mon, 04 Jun 2012) Log Message: ----------- fastlscomplex compiles now, comments cleaned up. (Still not tested.) Modified Paths: -------------- trunk/octave-forge/extra/lssa/fastlscomplex.cpp Modified: trunk/octave-forge/extra/lssa/fastlscomplex.cpp =================================================================== --- trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-04 11:41:32 UTC (rev 10561) +++ trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-04 20:12:13 UTC (rev 10562) @@ -1,17 +1,7 @@ -/* Copyright (C) 2012 Benjamin Lewis <be...@gm...> +/* * fastlscomplex.cpp, compiles to fastlscomplex.oct - * This program is free software; you can redistribute it and/or modify it under - * the terms of the GNU General Public License as published by the Free Software - * Foundation; either version 3 of the License, or (at your option) any later - * version. - * - * This program is distributed in the hope that it will be useful, but WITHOUT - * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or - * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more - * details. - * - * You should have received a copy of the GNU General Public License along with - * this program; if not, see <http://www.gnu.org/licenses/>. + * Conversion to C++, with wrapper for Octave of the code from + * A. Mathias' nuspectral package. */ #include <octave/oct.h> @@ -20,21 +10,18 @@ #include <string> #include <iostream> -typedef struct -{ double x, t; -} XTElem; +ComplexRowVector fastlscomplex ( RowVector tvals, ComplexRowVector xvals, octave_idx_type n, + double length, int ncoeff, int noctaves, double omegamax ); inline double sqr(double x) { return x*x; } -#define PNUM 12; - #define SETXT(p_, op_, x_, t_) (p_)->x op_ x_; (p_++)->t op_ t_; #define SETT(p_, op_, x_, t_) *p_++ op_ t_; #define SETX(p_, op_, x_, t_) *p_++ op_ x_; /* h is a complex aux. variable; it is used for assignment times I everywhere */ -#define SETIX(p_, op_, x_, t_) h = x_; RE(*(p_)) op_ -IM(h); IM(*(p_)) op_ RE(h); p_++; +#define SETIX(p_, op_, x_, t_) h = x_; (*(p_)).real() op_ -(h.imag()); (*(p_)).imag() op_ h.real(); p_++; /* Macro that sums up the power series terms into the power series * element record pointed to by p_. @@ -46,7 +33,7 @@ */ // -10 points, comments don't match method. #define EXP_IOT_SERIES(p_, el_, t_, op_, setr_, seti_) \ -{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ +{ double t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ tt = -t; seti_(p_, op_, x*tt, tt) \ tt *= t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \ tt *= t*(-1.0/3.0); seti_(p_, op_, x*tt, tt) \ @@ -62,7 +49,7 @@ /* same as the above, but without alternating signs */ #define EXPIOT_SERIES(p_, el_, t_, op_, setr_, seti_) \ -{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ +{ double t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ seti_(p_, op_, x*t, t ) \ tt = t*t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \ tt *= t*(1.0/3.0); seti_(p_, op_, x*tt, tt) \ @@ -76,29 +63,21 @@ tt *= t*(1.0/11.0); seti_(p_, op_, x*tt, tt) \ } -# define SRCARG Real *tptr, Real *xptr, int n, double *lengthptr -# define SRCVAR int k; Real length = *lengthptr; -# define SRCT tptr[k] -# define SRCX xptr[k] +# define SRCARG double *tptr, double *xptr, int n, double *lengthptr +# define SRCVAR int k; double length = *lengthptr; + +//I'll remove these very shortly # define SRCFIRST k = 0 -# define SRCAVAIL (k<n) +# define SRCAVAIL (k<n) # define SRCNEXT k++ -#define FASTLSCOMPLEX_HELP \ - "Takes the fast complex least-squares transform of irregularly sampled data.\n\ +DEFUN_DLD(fastlscomplex, args, nargout, "Takes the fast complex least-squares transform of irregularly sampled data.\n\ When called, takes a time series with associated x-values, a number of octaves,\n\ a number of coefficients per octave, the maximum frequency, and a length term.\n\ -It returns a complex row vector of the transformed series." - -DEFUN_DLD(fastlscomplex, args, nargout, FASTLSCOMPLEX_HELP ){ - if ( args.length != 7 ) { - std::cout << "Improper number of arguments!" << std::endl; - std::cout << FASTLSCOMPLEX_HELP << std::endl; - return octave_value_list (); - } +It returns a complex row vector of the transformed series." ){ RowVector tvals = args(0).row_vector_value(); ComplexRowVector xvals = args(1).complex_row_vector_value(); - double length = ( tvals(tvals.numel()--) - tvals( octave_idx_type (0))); + double length = ( tvals((tvals.numel()-1)) - tvals( octave_idx_type (0))); int ncoeff = args(4).int_value(); int noctaves = args(5).int_value(); double omegamax = args(6).double_value(); @@ -121,7 +100,7 @@ } // possibly include error for when omegamax isn't right? ComplexRowVector results = fastlscomplex(tvals, xvals, n, length, ncoeff, noctaves, omegamax); - return results; + return octave_value_list ( (octave_value) results ); } ComplexRowVector fastlscomplex ( RowVector tvals, ComplexRowVector xvals, octave_idx_type n, @@ -134,11 +113,10 @@ struct SumVec { SumVec *next; std::complex<double> elements[12]; - int count; } - *shead, *stail, *sp, *sq; //Names have been kept from the C, may change if I want to. + int count; }; + SumVec *shead, *stail, *sp, *sq; double dtelems[12], /* power series elements of exp(-i dtau) */ *dte, *r, /* Pointers into dtelems */ - x, /* abscissa and ordinate value, p-th power of t */ tau, tau0, te, /* Precomputation range centers and range end */ tau0d, /* tau_h of first summand range at level d */ dtau = (0.5*M_PI)/omegamax,/* initial precomputation interval radius */ @@ -149,48 +127,45 @@ on_1, /* n_1*(omega/mu)^p, n_1*(2*omega/mu)^p */ mu = (0.5*M_PI)/length, /* Frequency shift: a quarter period of exp(i mu t) on length */ tmp; - std::complex<double> zeta, zz, // Accumulators for spectral coefficients to place in Complex<double> + std::complex<double> zeta, zz, // Accumulators for spectral coefficients to place in complex<double> e, emul, + x, h, eoelems[12], oeelems[12], *eop, *oep, *ep, *op, *p, *q, *pe; - ComplexRowVector results (ncoeff*noctaves); //This *should* be the right size. I'll check against later sections of code. + ComplexRowVector results (ncoeff*noctaves); int i , j; // Counters; coefficient and octave, respectively. - // probable doubles: zetar, zetai, zzr, zzi, er, ei, emulr, emuli, - // eoelemsr[12] eoelemsi[12], etc. (since I want to use Complex<double> - // as little as possible.) - + octave_idx_type k = 0; // Subdivision and precomputation, reinvisioned in an OOWorld. tau = tvals(k) + dtau; te = tau + dtau; tau0 = tau; - shead = stail = sp = alloca(sizeof(*sp)); + shead = stail = sp = (SumVec*) operator new (sizeof(SumVec)); sp->next = 0; { te = tvals(k) + ( 2 * dtau ); - while ( k < n ) { //THIS makes no sense, n is the wrong type. Will need to fix that. - EXP_IOT_SERIES(p, sp->elems, mu*(tvals(k)-tau), =, SETX, SETIX); + while ( k < n ) { + EXP_IOT_SERIES(p, sp->elements, mu*(tvals(k)-tau), =, SETX, SETIX); sp->count = 1; //Sum terms and show that there has been at least one precomputation. // I will probably replace the macro with a better explanation. for(SRCNEXT; SRCAVAIL && tvals(k) < te; SRCNEXT) { x = xvals(k); - EXP_IOT_SERIES(p,sp->elems,mu*(tvals(k)-tau), +=, SETX, SETIX); + EXP_IOT_SERIES(p,sp->elements,mu*(tvals(k)-tau), +=, SETX, SETIX); sp->count++; } if ( k >= n ) break; tau = te + dtau; te = tau + dtau; - sp = alloca(sizeof(*sp)); + sp = (SumVec*) operator new (sizeof(*sp)); stail->next = sp; stail = sp; sp->next = 0; sp->count = 0; } } - // Now isn't that just a nicer representation of much the same control structure as that ugly for-loop? // Defining starting values for the loops over octaves: ooct = omegamax / mu; omul = exp(-log(2)/ncoeff); @@ -199,7 +174,7 @@ dtaud = dtau; octave_idx_type iter ( 0 ); // Looping over octaves - for ( j = noctave ; ; ooct *= 0.5 , omegaoct *= 0.5 , tau0d += dtaud , dtaud *= 2 ) { + for ( j = noctaves ; ; ooct *= 0.5 , omegaoct *= 0.5 , tau0d += dtaud , dtaud *= 2 ) { // Looping over&results per frequency for ( i = ncoeff, o = ooct, omega = omegaoct; i-- ; o *= omul, omega *= omul ) { e.real() = cos( - ( omega * tau0d ) ); e.imag() = sin( - ( omega * tau0d ) ); @@ -208,7 +183,9 @@ // sets real, imag parts of emul for( zeta = 0 , sp = shead; sp; sp = sp->next, e *= emul ) { if ( sp->count ) { - for ( zz = std::complex<double>(0.0,0.0) , octave_idx_type i (0) , p = sp->elems , pe = p + 12 , on_1 = n_1 ; p < pe ; i++ ) { + zz = std::complex<double>(0.0,0.0); + octave_idx_type i ( 0 ); + for ( p = sp->elements , on_1 = n_1 ; i < (octave_idx_type) 12 ; i++ ) { zz += *p++ * on_1; on_1 *= 0; } @@ -222,41 +199,41 @@ EXPIOT_SERIES(r, dtelems, mu*dtaud, =, SETT, SETT); for(sp = shead; sp; sp = sp->next){ if(!(sq = sp->next) || !sq->count ) { - for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+12; p < pe; p++, dte++) + for(p = sp->elements, eop = eoelems, dte = dtelems+1, pe = p+12; p < pe; p++, dte++) { ep = eop; *eop++ = *p; for(r = dtelems, *p = *ep * *r; ; ) - { if(++r>=dte) break; --ep; h = *ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h); + { if(++r>=dte) break; --ep; h = *ep * *r; (*p).real() -= h.imag(); (*p).imag() += h.real(); if(++r>=dte) break; --ep; *p -= *ep * *r; - if(++r>=dte) break; --ep; h = -*ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h); + if(++r>=dte) break; --ep; h = -*ep * *r; (*p).real() -= h.imag(); (*p).imag() += h.real(); if(++r>=dte) break; --ep; *p += *ep * *r; } } if(!sq) break; /* reached the last precomputation range */ } else - if(sp->cnt) - for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+12; + if(sp->count) + for(p = sp->elements, q = sq->elements, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+12; p < pe; p++, q++, dte++) { ep = eop; *eop++ = *p+*q; *oep++ = *p-*q; op = oep; for(r = dtelems, *p = *ep * *r; ; ) - { if(++r>=dte) break; op -= 2; h = *op * *r; RE(*p) -= IM(h); IM(*p) += RE(h); + { if(++r>=dte) break; op -= 2; h = *op * *r; (*p).real() -= h.imag(); (*p).imag() += h.real(); if(++r>=dte) break; ep -= 2; *p -= *ep * *r; - if(++r>=dte) break; op -= 2; h = -*op * *r; RE(*p) -= IM(h); IM(*p) += RE(h); + if(++r>=dte) break; op -= 2; h = -*op * *r; (*p).real() -= h.imag(); (*p).imag() += h.real(); if(++r>=dte) break; ep -= 2; *p += *ep * *r; } } else - for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+12; q<pe; q++, dte++) + for(q = sq->elements, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+12; q<pe; q++, dte++) { ep = eop; *eop++ = *q; for(r = dtelems, *q = *ep * *r; ; ) - { if(++r>=dte) break; --ep; h = *ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h); + { if(++r>=dte) break; --ep; h = *ep * *r; (*q).real() -= h.imag(); (*q).imag() += h.real(); if(++r>=dte) break; --ep; *p -= *ep * *r; - if(++r>=dte) break; --ep; h = -*ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h); + if(++r>=dte) break; --ep; h = -*ep * *r; (*q).real() -= h.imag(); (*q).imag() += h.real(); if(++r>=dte) break; --ep; *q += *ep * *r; } } - sp->cnt += sq->cnt; sp->next = sq->next; /* free(sq) if malloc'ed */ + sp->count += sq->count; sp->next = sq->next; /* free(sq) if malloc'ed */ } } } This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <aba...@us...> - 2012-06-04 11:41:41
|
Revision: 10561 http://octave.svn.sourceforge.net/octave/?rev=10561&view=rev Author: abarth93 Date: 2012-06-04 11:41:32 +0000 (Mon, 04 Jun 2012) Log Message: ----------- ncread Added Paths: ----------- trunk/octave-forge/main/octcdf/inst/ncread.m Added: trunk/octave-forge/main/octcdf/inst/ncread.m =================================================================== --- trunk/octave-forge/main/octcdf/inst/ncread.m (rev 0) +++ trunk/octave-forge/main/octcdf/inst/ncread.m 2012-06-04 11:41:32 UTC (rev 10561) @@ -0,0 +1,68 @@ +% x = ncread(filename,varname) +% x = ncread(filename,varname,start,count,stride) +% read the variable varname from file filename. + +function x = ncread(filename,varname,start,count,stride) + +nc = netcdf(filename,'r'); +nv = nc{varname}; +sz = size(nv); sz = sz(end:-1:1); +nd = ndims(nv); + +if nargin == 2 + start = ones(1,nd); + count = inf*ones(1,nd); +end + +if nargin < 5 + stride = ones(1,nd); +end + +% end index + +endi = start + count.*stride; + +% replace inf in count + +i = endi == inf; +endi(i) = sz(i); + + +% load data + +% subsref structure +subsr.type = '()'; +subsr.subs = cell(1,nd); +for i=1:nd + subsr.subs{nd-i+1} = start(i):stride(i):endi(i); +end +start,endi + +x = subsref(nv,subsr); + +% apply attributes + +factor = nv.scale_factor(:); +offset = nv.add_offset(:); +fv = nv.FillValue_; + +if ~isempty(fv) + x(x == fv) = NaN; +else + fv = nv.missing_value; + + if ~isempty(fv) + x(x == fv) = NaN; + end +end + +if ~isempty(factor) + x = x * factor; +end + +if ~isempty(offset) + x = x + offset; +end + +x = permute(x,[ndims(x):-1:1]); + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-03 19:34:16
|
Revision: 10560 http://octave.svn.sourceforge.net/octave/?rev=10560&view=rev Author: paramaniac Date: 2012-06-03 19:34:09 +0000 (Sun, 03 Jun 2012) Log Message: ----------- control: add uplus operator for LTI models Modified Paths: -------------- trunk/octave-forge/main/control/INDEX Added Paths: ----------- trunk/octave-forge/main/control/inst/@lti/uplus.m Modified: trunk/octave-forge/main/control/INDEX =================================================================== --- trunk/octave-forge/main/control/INDEX 2012-06-03 14:04:04 UTC (rev 10559) +++ trunk/octave-forge/main/control/INDEX 2012-06-03 19:34:09 UTC (rev 10560) @@ -120,6 +120,7 @@ @lti/subsref @lti/transpose @lti/uminus + @lti/uplus @lti/vertcat Miscellaneous options Added: trunk/octave-forge/main/control/inst/@lti/uplus.m =================================================================== --- trunk/octave-forge/main/control/inst/@lti/uplus.m (rev 0) +++ trunk/octave-forge/main/control/inst/@lti/uplus.m 2012-06-03 19:34:09 UTC (rev 10560) @@ -0,0 +1,31 @@ +## Copyright (C) 2012 Lukas F. Reichlin +## +## This file is part of LTI Syncope. +## +## LTI Syncope is free software: you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation, either version 3 of the License, or +## (at your option) any later version. +## +## LTI Syncope is distributed in the hope that it will be useful, +## but WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +## GNU General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with LTI Syncope. If not, see <http://www.gnu.org/licenses/>. + +## -*- texinfo -*- +## Unary plus of LTI object. Used by Octave for "+sys". + +## Author: Lukas Reichlin <luk...@gm...> +## Created: June 2012 +## Version: 0.1 + +function sys = uplus (sys) + + if (nargin != 1) # prevent sys = uplus (sys1, sys2, sys3, ...) + error ("lti: uplus: this is an unary operator"); + endif + +endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-03 14:04:10
|
Revision: 10559 http://octave.svn.sourceforge.net/octave/?rev=10559&view=rev Author: paramaniac Date: 2012-06-03 14:04:04 +0000 (Sun, 03 Jun 2012) Log Message: ----------- control: remove whitespace between @ and { which caused problems with generate_html Modified Paths: -------------- trunk/octave-forge/main/control/doc/control.pdf trunk/octave-forge/main/control/inst/@lti/transpose.m Modified: trunk/octave-forge/main/control/doc/control.pdf =================================================================== --- trunk/octave-forge/main/control/doc/control.pdf 2012-06-03 13:51:44 UTC (rev 10558) +++ trunk/octave-forge/main/control/doc/control.pdf 2012-06-03 14:04:04 UTC (rev 10559) @@ -3279,8 +3279,8 @@ 966 0 obj << /Producer (pdfTeX-1.40.12) /Creator (TeX) -/CreationDate (D:20120603154532+02'00') -/ModDate (D:20120603154532+02'00') +/CreationDate (D:20120603160139+02'00') +/ModDate (D:20120603160139+02'00') /Trapped /False /PTEX.Fullbanner (This is pdfTeX, Version 3.1415926-2.3-1.40.12 (TeX Live 2011) kpathsea version 6.0.1) >> endobj @@ -3313,7 +3313,7 @@ /W [1 3 1] /Root 965 0 R /Info 966 0 R -/ID [<CA738F7C45E3205DA4F6C76549534E2C> <CA738F7C45E3205DA4F6C76549534E2C>] +/ID [<E6F8037627641F1152201403FA3918DF> <E6F8037627641F1152201403FA3918DF>] /Length 2241 /Filter /FlateDecode >> Modified: trunk/octave-forge/main/control/inst/@lti/transpose.m =================================================================== --- trunk/octave-forge/main/control/inst/@lti/transpose.m 2012-06-03 13:51:44 UTC (rev 10558) +++ trunk/octave-forge/main/control/inst/@lti/transpose.m 2012-06-03 14:04:04 UTC (rev 10559) @@ -18,7 +18,7 @@ ## -*- texinfo -*- ## Transpose of LTI objects. Used by Octave for "sys.'". ## Useful for dual problems, i.e. controllability and observability -## or designing estimator gains with @command{lqr} and @command {place}. +## or designing estimator gains with @command{lqr} and @command{place}. ## Author: Lukas Reichlin <luk...@gm...> ## Created: February 2010 This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-03 13:51:53
|
Revision: 10558 http://octave.svn.sourceforge.net/octave/?rev=10558&view=rev Author: paramaniac Date: 2012-06-03 13:51:44 +0000 (Sun, 03 Jun 2012) Log Message: ----------- control: prepare release of control-2.3.51 Modified Paths: -------------- trunk/octave-forge/main/control/DESCRIPTION trunk/octave-forge/main/control/NEWS trunk/octave-forge/main/control/devel/RELEASE_PACKAGE trunk/octave-forge/main/control/devel/lqe.m trunk/octave-forge/main/control/devel/makefile_control.m trunk/octave-forge/main/control/devel/pdfdoc/control.tex trunk/octave-forge/main/control/doc/control.pdf trunk/octave-forge/main/control/inst/filt.m trunk/octave-forge/main/control/inst/optiPID.m trunk/octave-forge/main/control/inst/test_control.m trunk/octave-forge/main/control/src/Makefile Added Paths: ----------- trunk/octave-forge/main/control/devel/lapack-3.4.1.tgz Removed Paths: ------------- trunk/octave-forge/main/control/src/lapack-3.4.1.tgz Modified: trunk/octave-forge/main/control/DESCRIPTION =================================================================== --- trunk/octave-forge/main/control/DESCRIPTION 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/DESCRIPTION 2012-06-03 13:51:44 UTC (rev 10558) @@ -1,6 +1,6 @@ Name: Control -Version: 2.3.50 -Date: 2012-03-06 +Version: 2.3.51 +Date: 2012-06-03 Author: Lukas Reichlin <luk...@gm...> Maintainer: Lukas Reichlin <luk...@gm...> Title: Control Systems Modified: trunk/octave-forge/main/control/NEWS =================================================================== --- trunk/octave-forge/main/control/NEWS 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/NEWS 2012-06-03 13:51:44 UTC (rev 10558) @@ -1,7 +1,7 @@ Summary of important user-visible changes for releases of the control package =============================================================================== -control-2.3.51 Release Date: 2012-xx-yy Release Manager: Lukas Reichlin +control-2.3.51 Release Date: 2012-06-03 Release Manager: Lukas Reichlin =============================================================================== ** filt, filtdata, tf @@ -16,15 +16,18 @@ ** ctranspose Conjugate transpose or pertransposition of LTI objects. - Used by Octave for "sys’". For a transfer-function matrix G, G’ denotes the - conjugate of G given by G.’(-s) for a continuous-time system or G.’(1/z) for + Used by Octave for "sys'". For a transfer-function matrix G, G' denotes the + conjugate of G given by G.'(-s) for a continuous-time system or G.'(1/z) for a discrete-time system. The frequency response of the pertransposition of G is the Hermitian (conjugate) transpose of G(jw), - i.e. freqresp (G’, w) = freqresp (G, w)’. - WARNING: Do NOT use this for dual problems, use the transpose "sys.’" + i.e. freqresp (G', w) = freqresp (G, w)'. + WARNING: Do NOT use this for dual problems, use the transpose "sys.'" (note the dot) instead. +** test_control + Add a few remarks to the help text regarding the severity of failing tests. + =============================================================================== control-2.3.50 Release Date: 2012-03-06 Release Manager: Lukas Reichlin =============================================================================== Modified: trunk/octave-forge/main/control/devel/RELEASE_PACKAGE =================================================================== --- trunk/octave-forge/main/control/devel/RELEASE_PACKAGE 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/devel/RELEASE_PACKAGE 2012-06-03 13:51:44 UTC (rev 10558) @@ -20,12 +20,12 @@ rm -R ~/octave/__TEMP__/control/devel cd ~/octave/__TEMP__ grep -i version control/DESCRIPTION -tar czf control-2.3.50.tar.gz control/ -md5 control-2.3.50.tar.gz -md5 control-2.3.50.tar.gz > md5_control_pkg.txt -uuencode control-2.3.50.tar.gz < control-2.3.50.tar.gz > control-2.3.50.tar.gz.uue +tar czf control-2.3.51.tar.gz control/ +md5 control-2.3.51.tar.gz +md5 control-2.3.51.tar.gz > md5_control_pkg.txt +uuencode control-2.3.51.tar.gz < control-2.3.51.tar.gz > control-2.3.51.tar.gz.uue octave -q --eval \ -"pkg install control-2.3.50.tar.gz" +"pkg install control-2.3.51.tar.gz" octave -q --eval \ "pkg load generate_html; generate_package_html ('control', 'control-html', 'octave-forge')" tar czf control-html.tar.gz control-html @@ -40,7 +40,7 @@ ===================================================================================== rm -R ~/octave/__TEMP__ -rm -R ~/octave/control-2.3.50 +rm -R ~/octave/control-2.3.51 ===================================================================================== Copied: trunk/octave-forge/main/control/devel/lapack-3.4.1.tgz (from rev 10554, trunk/octave-forge/main/control/src/lapack-3.4.1.tgz) =================================================================== (Binary files differ) Modified: trunk/octave-forge/main/control/devel/lqe.m =================================================================== --- trunk/octave-forge/main/control/devel/lqe.m 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/devel/lqe.m 2012-06-03 13:51:44 UTC (rev 10558) @@ -65,7 +65,7 @@ ## L = eig (A - B*G) ## @end group ## @end example -## @seealso{care, dare, dlqr} +## @seealso{lqr, care, dare, dlqe} ## @end deftypefn ## Author: Lukas Reichlin <luk...@gm...> Modified: trunk/octave-forge/main/control/devel/makefile_control.m =================================================================== --- trunk/octave-forge/main/control/devel/makefile_control.m 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/devel/makefile_control.m 2012-06-03 13:51:44 UTC (rev 10558) @@ -13,7 +13,7 @@ ## system ("make realclean"); # recompile slicotlibrary.a system ("make clean"); -system ("make -j4 all"); +system ("make -j1 all"); system ("rm *.o"); cd (homedir); Modified: trunk/octave-forge/main/control/devel/pdfdoc/control.tex =================================================================== --- trunk/octave-forge/main/control/devel/pdfdoc/control.tex 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/devel/pdfdoc/control.tex 2012-06-03 13:51:44 UTC (rev 10558) @@ -3,7 +3,7 @@ @setfilename control.info @settitle Octave Control Systems Package @afourpaper -@set VERSION 2.3.50 +@set VERSION 2.3.51 @finalout @c @afourwide @c %**end of header Modified: trunk/octave-forge/main/control/doc/control.pdf =================================================================== --- trunk/octave-forge/main/control/doc/control.pdf 2012-06-02 20:50:08 UTC (rev 10557) +++ trunk/octave-forge/main/control/doc/control.pdf 2012-06-03 13:51:44 UTC (rev 10558) @@ -42,7 +42,7 @@ >> stream xڅ\x90MK1\x86\xEF\xFB+rL3\x99\xCC&\xB9*\xB6 \xE2Gݞ\xC4ònk\xE9n\xAB\xE0\xBF7i\/\x8A\xC3̛y\x9F\xBC(L:(\xBC>\xF2\xA2\xDBV\xE6\xD4\x96\xA2\xB3i\x85_:\xC7\xEC(\xD5\xCA\xF4\xA8\xD3H\xDD?Kur\x8Eb\x9E7\xD5\xD9ij\xB0j,\x9AEƫ -\x817A4\xCF\xE2Qv\xFB\x9D"\x94\xC7a\xBFQ\x9A<K\xC9Ϩ\xA7檺l\xAA\xD7\xEF\xD5\xE9!GOaD\xC9W\xDA=\xF4bQ\xDDKf\x81\xA2\x89\x98-\xC9\xD4\xE8}\xF1\xBDH\xBE6\x8C\xBED\xF2\xE1\xE3p\xEC\xB7J\xA3<\x94\xC6]\x9E\xB7]\xBEY\xB6˾\xFB\xA1ӛy.\x9C\xBC\xED\x8Emֽ\xE7\xAB\xCF\xEC\xA7o\xA3r\xB5\xCD\xDA2\x9E \xB4\x8B\xA2-\xD7okU[\xD9fG\xEF\xE4J\xB3>c\xAC:EV\xBElV\xBB_\xC3\xF0\xEC\x8F(b\xA41\x8AO\x94\xC6s\x8C +\x817A4\xCF\xE2Qv\xFB\x9D"\x94\xC7a\xBFQ\x9A<K\xC9\xD5SsU]6\xD5\xEB\xF7\xEA\xF4\x90#\x83\xA70"\x8Cd\x89+\xEDz\xB1\xA8\xEE\x8B%\xB3@\xD1D̖d j\x8C\x82\x83\xF4\xBE\xF8^$_F_"\xF9\xF0q8\xF6[\xA5QJ\xE3.\xCF\xDB.\xDFk\xC5,\xDBe_\x8B\xFDP\x8A\xE9\xCD<N\xDEv\xC76\xEB\xDE\xF3\xD5g\xF6ӷ\xD1\xB9\xDAfmO\xDA\xC5 |
From: <mtm...@us...> - 2012-06-02 20:50:15
|
Revision: 10557 http://octave.svn.sourceforge.net/octave/?rev=10557&view=rev Author: mtmiller Date: 2012-06-02 20:50:08 +0000 (Sat, 02 Jun 2012) Log Message: ----------- marcumq: use tablify again now that it is in general The signal package now depends on general (>= 1.3.2) and no longer depends on image, no more circular dependency. Modified Paths: -------------- trunk/octave-forge/main/signal/DESCRIPTION trunk/octave-forge/main/signal/NEWS trunk/octave-forge/main/signal/inst/marcumq.m Modified: trunk/octave-forge/main/signal/DESCRIPTION =================================================================== --- trunk/octave-forge/main/signal/DESCRIPTION 2012-06-02 20:39:50 UTC (rev 10556) +++ trunk/octave-forge/main/signal/DESCRIPTION 2012-06-02 20:50:08 UTC (rev 10557) @@ -5,7 +5,7 @@ Maintainer: Octave-Forge community <oct...@li...> Title: Signal Processing. Description: Signal processing tools, including filtering, windowing and display functions. -Depends: octave (>= 3.6.0), optim (>= 1.0.0), specfun, control (>= 2.2.3), image +Depends: octave (>= 3.6.0), optim (>= 1.0.0), specfun, control (>= 2.2.3), general (>= 1.3.2) Autoload: no License: GPLv3+, public domain Url: http://octave.sf.net Modified: trunk/octave-forge/main/signal/NEWS =================================================================== --- trunk/octave-forge/main/signal/NEWS 2012-06-02 20:39:50 UTC (rev 10556) +++ trunk/octave-forge/main/signal/NEWS 2012-06-02 20:50:08 UTC (rev 10557) @@ -1,6 +1,14 @@ Summary of important user-visible changes for releases of the signal package =============================================================================== +signal-1.1.4 Release Date: 2012-XX-XX Release Manager: +=============================================================================== + + ** signal is no longer dependent on the image package. + + ** signal is now dependent on the general package. + +=============================================================================== signal-1.1.3 Release Date: 2012-05-12 Release Manager: Carnë Draug =============================================================================== Modified: trunk/octave-forge/main/signal/inst/marcumq.m =================================================================== --- trunk/octave-forge/main/signal/inst/marcumq.m 2012-06-02 20:39:50 UTC (rev 10556) +++ trunk/octave-forge/main/signal/inst/marcumq.m 2012-06-02 20:50:08 UTC (rev 10557) @@ -51,11 +51,7 @@ error("M must be a positive integer"); end - nr = max([size(a,1) size(b,1) size(M,1)]); - nc = max([size(a,2) size(b,2) size(M,2)]); - a = padarray(a, [nr - size(a,1) nc - size(a,2)], "replicate", "post"); - b = padarray(b, [nr - size(b,1) nc - size(b,2)], "replicate", "post"); - M = padarray(M, [nr - size(M,1) nc - size(M,2)], "replicate", "post"); + [a,b,M] = tablify(a,b,M); Q = arrayfun(@mq, a,b,M,tol); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <mtm...@us...> - 2012-06-02 20:39:56
|
Revision: 10556 http://octave.svn.sourceforge.net/octave/?rev=10556&view=rev Author: mtmiller Date: 2012-06-02 20:39:50 +0000 (Sat, 02 Jun 2012) Log Message: ----------- tablify: new function added to general Contributed by Robert T. Short. Function is used by marcumq in the signal package but could be useful to others. Closes artifact #3522120. Modified Paths: -------------- trunk/octave-forge/main/general/INDEX trunk/octave-forge/main/general/NEWS Added Paths: ----------- trunk/octave-forge/main/general/inst/tablify.m Modified: trunk/octave-forge/main/general/INDEX =================================================================== --- trunk/octave-forge/main/general/INDEX 2012-06-02 16:39:30 UTC (rev 10555) +++ trunk/octave-forge/main/general/INDEX 2012-06-02 20:39:50 UTC (rev 10556) @@ -27,6 +27,7 @@ packfields safeprod SHA1 + tablify unpackfields unresamp2 unvech Modified: trunk/octave-forge/main/general/NEWS =================================================================== --- trunk/octave-forge/main/general/NEWS 2012-06-02 16:39:30 UTC (rev 10555) +++ trunk/octave-forge/main/general/NEWS 2012-06-02 20:39:50 UTC (rev 10556) @@ -1,3 +1,9 @@ +Summary of important user-visible changes for general 1.3.2: +------------------------------------------------------------------- + + ** The following functions are new: + tablify + Summary of important user-visible changes for general 1.3.1: ------------------------------------------------------------------- Added: trunk/octave-forge/main/general/inst/tablify.m =================================================================== --- trunk/octave-forge/main/general/inst/tablify.m (rev 0) +++ trunk/octave-forge/main/general/inst/tablify.m 2012-06-02 20:39:50 UTC (rev 10556) @@ -0,0 +1,148 @@ +## Copyright (C) 2012 Robert T. Short +## +## This 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 software 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/>. + +## -*- texinfo -*- +## @deftypefn {Function File} {[@var{y1}, @dots{}] =} tablify (@var{x1}, @dots{}) +## +## @noindent +## Create a table out of the input arguments, if possible. The table is +## created by extending row and column vectors to like dimensions. If +## the dimensions of input vectors are not commensurate an error will +## occur. Dimensions are commensurate if they have the same number of +## rows and columns, a single row and the same number of columns, or the +## same number of rows and a single column. In other words, vectors +## will only be extended along singleton dimensions. +## +## @noindent +## For example: +## +## @example +## @group +## [a, b] = tablify ([1 2; 3 4], 5) +## @result{} a = [ 1, 2; 3, 4 ] +## @result{} b = [ 5, 5; 5, 5 ] +## @end group +## @end example +## @example +## @group +## [a, b, c] = tablify (1, [1 2 3 4], [5;6;7]) +## @result{} +## b = [ 1 1 1 1; 1 1 1 1; 1 1 1 1] +## @result{} b = [ 1 2 3 4; 1 2 3 4; 1 2 3 4] +## @result{} c = [ 5 5 5 5; 6 6 6 6; 7 7 7 7 ] +## @end group +## @end example +## +## @noindent +## The following example attempts to expand vectors that do not have +## commensurate dimensions and will produce an error. +## +## @example +## @group +## tablify([1 2],[3 4 5]) +## @end group +## @end example +## +## @noindent +## Note that use of array operations and broadcasting is more efficient +## for many situations. +## +## @seealso {common_size} +## +## @end deftypefn + +## Author: Robert T. Short +## Created: 3/6/2012 + +function [varargout] = tablify (varargin) + + if (nargin < 2) + varargout = varargin; + return; + endif + + empty = cellfun (@isempty, varargin); + + nrows = cellfun (@rows, varargin(~empty)); + ridx = nrows>1; + if (any(ridx)) + rdim = nrows(ridx)(1); + else + rdim = 1; + endif + + ncols = cellfun (@columns, varargin(~empty)); + cidx = ncols>1; + if (any(cidx)) + cdim = ncols(cidx)(1); + else + cdim = 1; + endif + + if ( any(nrows(ridx) != rdim ) || any(ncols(cidx) != cdim )) + error('tablify: incommensurate sizes'); + endif + + varargout = varargin; + varargout(~ridx) = cellindexmat(varargout(~ridx),ones(rdim,1),':'); + varargout(~cidx) = cellindexmat(varargout(~cidx),':',ones(1,cdim)); + +endfunction + +%!error tablify([1,2],[3,4,5]) +%!error tablify([1;2],[3;4;5]) + +%!test +%! a = 1; b = 1; c = 1; +%! assert(tablify(a),a); +%! [x,y,z]=tablify(a,b,c); +%! assert(x,a); +%! assert(y,b); +%! assert(z,c); + +%!test +%! a = 1; b = [1 2 3]; +%! [x,y]=tablify(a,b); +%! assert(x,[1 1 1]); +%! assert(y,[1 2 3]); + +%!test +%! a = 1; b = [1;2;3]; +%! [x,y]=tablify(a,b); +%! assert(x,[1;1;1]); +%! assert(y,[1;2;3]); + +%!test +%! a = [1 2]; b = [1;2;3]; c=1; +%! [x,y,z]=tablify(a,b,c); +%! assert(x,[1 2; 1 2; 1 2]); +%! assert(y,[1 1; 2 2; 3 3]); +%! assert(z,[1 1; 1 1; 1 1]); + +%!test +%! a = [1 2]; b = [1;2;3]; c=[2 3]; +%! [x,y,z]=tablify(a,b,c); +%! assert(x,[1 2; 1 2; 1 2]); +%! assert(y,[1 1; 2 2; 3 3]); +%! assert(z,[2 3; 2 3; 2 3]); + +%!test +%! a = [1 2]; b = [1;2;3]; c=[]; +%! [x,y,z]=tablify(a,b,c); +%! assert(x,[1 2; 1 2; 1 2]); +%! assert(y,[1 1; 2 2; 3 3]); +%! assert(z,[]); + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <be...@us...> - 2012-06-02 16:39:37
|
Revision: 10555 http://octave.svn.sourceforge.net/octave/?rev=10555&view=rev Author: benjf5 Date: 2012-06-02 16:39:30 +0000 (Sat, 02 Jun 2012) Log Message: ----------- A few fixes to fastlscomplex, it's getting closer to compiling. Modified Paths: -------------- trunk/octave-forge/extra/lssa/fastlscomplex.cpp Modified: trunk/octave-forge/extra/lssa/fastlscomplex.cpp =================================================================== --- trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-02 07:31:35 UTC (rev 10554) +++ trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-02 16:39:30 UTC (rev 10555) @@ -15,6 +15,7 @@ */ #include <octave/oct.h> +#include <math.h> #include <complex> #include <string> #include <iostream> @@ -28,6 +29,7 @@ #define PNUM 12; + #define SETXT(p_, op_, x_, t_) (p_)->x op_ x_; (p_++)->t op_ t_; #define SETT(p_, op_, x_, t_) *p_++ op_ t_; #define SETX(p_, op_, x_, t_) *p_++ op_ x_; @@ -130,11 +132,11 @@ * unused entries. */ struct SumVec { - struct SumVec *next; - std::complex<double> elements[PNUM]; + SumVec *next; + std::complex<double> elements[12]; int count; } *shead, *stail, *sp, *sq; //Names have been kept from the C, may change if I want to. - double dtelems[PNUM], /* power series elements of exp(-i dtau) */ + double dtelems[12], /* power series elements of exp(-i dtau) */ *dte, *r, /* Pointers into dtelems */ x, /* abscissa and ordinate value, p-th power of t */ tau, tau0, te, /* Precomputation range centers and range end */ @@ -149,7 +151,7 @@ tmp; std::complex<double> zeta, zz, // Accumulators for spectral coefficients to place in Complex<double> e, emul, - h, eoelems[PNUM], oeelems[PNUM], + h, eoelems[12], oeelems[12], *eop, *oep, *ep, *op, *p, *q, *pe; @@ -157,7 +159,7 @@ int i , j; // Counters; coefficient and octave, respectively. // probable doubles: zetar, zetai, zzr, zzi, er, ei, emulr, emuli, - // eoelemsr[PNUM] eoelemsi[PNUM], etc. (since I want to use Complex<double> + // eoelemsr[12] eoelemsi[12], etc. (since I want to use Complex<double> // as little as possible.) octave_idx_type k = 0; @@ -191,7 +193,7 @@ // Now isn't that just a nicer representation of much the same control structure as that ugly for-loop? // Defining starting values for the loops over octaves: ooct = omegamax / mu; - omul = exp(-M_LN2/ncoeff); + omul = exp(-log(2)/ncoeff); omegaoct = omegamax; tau0d = tau0; dtaud = dtau; @@ -206,7 +208,7 @@ // sets real, imag parts of emul for( zeta = 0 , sp = shead; sp; sp = sp->next, e *= emul ) { if ( sp->count ) { - for ( zz = std::complex<double>(0.0,0.0) , octave_idx_type i (0) , p = sp->elems , pe = p + PNUM , on_1 = n_1 ; p < pe ; i++ ) { + for ( zz = std::complex<double>(0.0,0.0) , octave_idx_type i (0) , p = sp->elems , pe = p + 12 , on_1 = n_1 ; p < pe ; i++ ) { zz += *p++ * on_1; on_1 *= 0; } @@ -220,7 +222,7 @@ EXPIOT_SERIES(r, dtelems, mu*dtaud, =, SETT, SETT); for(sp = shead; sp; sp = sp->next){ if(!(sq = sp->next) || !sq->count ) { - for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+PNUM; p < pe; p++, dte++) + for(p = sp->elems, eop = eoelems, dte = dtelems+1, pe = p+12; p < pe; p++, dte++) { ep = eop; *eop++ = *p; for(r = dtelems, *p = *ep * *r; ; ) { if(++r>=dte) break; --ep; h = *ep * *r; RE(*p) -= IM(h); IM(*p) += RE(h); @@ -233,7 +235,7 @@ } else if(sp->cnt) - for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+PNUM; + for(p = sp->elems, q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = p+12; p < pe; p++, q++, dte++) { ep = eop; *eop++ = *p+*q; *oep++ = *p-*q; op = oep; for(r = dtelems, *p = *ep * *r; ; ) @@ -244,7 +246,7 @@ } } else - for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+PNUM; q<pe; q++, dte++) + for(q = sq->elems, eop = eoelems, oep = oeelems, dte = dtelems+1, pe = q+12; q<pe; q++, dte++) { ep = eop; *eop++ = *q; for(r = dtelems, *q = *ep * *r; ; ) { if(++r>=dte) break; --ep; h = *ep * *r; RE(*q) -= IM(h); IM(*q) += RE(h); This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-02 07:31:41
|
Revision: 10554 http://octave.svn.sourceforge.net/octave/?rev=10554&view=rev Author: paramaniac Date: 2012-06-02 07:31:35 +0000 (Sat, 02 Jun 2012) Log Message: ----------- optim: use the same argument names as matlab Modified Paths: -------------- trunk/octave-forge/main/optim/inst/fminsearch.m Modified: trunk/octave-forge/main/optim/inst/fminsearch.m =================================================================== --- trunk/octave-forge/main/optim/inst/fminsearch.m 2012-06-02 07:16:53 UTC (rev 10553) +++ trunk/octave-forge/main/optim/inst/fminsearch.m 2012-06-02 07:31:35 UTC (rev 10554) @@ -14,8 +14,8 @@ ## this program; if not, see <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn {Function File} {@var{x} =} fminsearch (@var{f}, @var{x0}) -## @deftypefnx {Function File} {[@var{x}, @var{fval}] =} fminsearch (@var{f}, @var{x0}, @var{options}, @var{grad}, @var{P1}, @var{P2}, @dots{}) +## @deftypefn {Function File} {@var{x} =} fminsearch (@var{fun}, @var{x0}) +## @deftypefnx {Function File} {[@var{x}, @var{fval}] =} fminsearch (@var{fun}, @var{x0}, @var{options}, @var{grad}, @var{P1}, @var{P2}, @dots{}) ## ## Find the minimum of a funtion of several variables. ## By default the method used is the Nelder&Mead Simplex algorithm. @@ -24,13 +24,13 @@ ## TODO: describe arguments in texinfo help string -function [x, fval] = fminsearch (funfun, x0, options = [], grad = [], varargin) +function [x, fval] = fminsearch (fun, x0, options = [], grad = [], varargin) if (nargin < 2) print_usage (); endif - x = fmins (funfun, x0, options, grad, varargin{:}); - fval = feval (funfun, x, varargin{:}); + x = fmins (fun, x0, options, grad, varargin{:}); + fval = feval (fun, x, varargin{:}); -endfunction; +endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-06-02 07:23:31
|
Revision: 10553 http://octave.svn.sourceforge.net/octave/?rev=10553&view=rev Author: paramaniac Date: 2012-06-02 07:16:53 +0000 (Sat, 02 Jun 2012) Log Message: ----------- optim: fix and simplify fminsearch Modified Paths: -------------- trunk/octave-forge/main/optim/inst/fminsearch.m Modified: trunk/octave-forge/main/optim/inst/fminsearch.m =================================================================== --- trunk/octave-forge/main/optim/inst/fminsearch.m 2012-06-01 22:37:39 UTC (rev 10552) +++ trunk/octave-forge/main/optim/inst/fminsearch.m 2012-06-02 07:16:53 UTC (rev 10553) @@ -14,18 +14,23 @@ ## this program; if not, see <http://www.gnu.org/licenses/>. ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{x}] =} fminsearch(@var{f},@var{X0},@var{options},@var{grad},@var{P1},@var{P2}, @dots{}) +## @deftypefn {Function File} {@var{x} =} fminsearch (@var{f}, @var{x0}) +## @deftypefnx {Function File} {[@var{x}, @var{fval}] =} fminsearch (@var{f}, @var{x0}, @var{options}, @var{grad}, @var{P1}, @var{P2}, @dots{}) ## ## Find the minimum of a funtion of several variables. -## By default the method used is the Nelder&Mead Simplex algorithm +## By default the method used is the Nelder&Mead Simplex algorithm. ## @seealso{fmin,fmins,nmsmax} ## @end deftypefn -function [x fval] = fminsearch(funfun, X0, options, grad, varargin) - if (nargin == 0); usage('[x fval] = fminsearch(funfun, X0, options, grad, varargin)'); end - if (nargin < 3); options=[]; end - if (nargin < 4); grad=[]; end - if (nargin < 5); varargin={}; end - x = fmins(funfun, X0, options, grad, varargin{:}); - fval = feval(funfun, x, varargin{:}); +## TODO: describe arguments in texinfo help string + +function [x, fval] = fminsearch (funfun, x0, options = [], grad = [], varargin) + + if (nargin < 2) + print_usage (); + endif + + x = fmins (funfun, x0, options, grad, varargin{:}); + fval = feval (funfun, x, varargin{:}); + endfunction; This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <be...@us...> - 2012-06-01 22:37:46
|
Revision: 10552 http://octave.svn.sourceforge.net/octave/?rev=10552&view=rev Author: benjf5 Date: 2012-06-01 22:37:39 +0000 (Fri, 01 Jun 2012) Log Message: ----------- Wrote the Octave wrapper for fastlscomplex, doesn't compile yet. Modified Paths: -------------- trunk/octave-forge/extra/lssa/fastlscomplex.cpp Modified: trunk/octave-forge/extra/lssa/fastlscomplex.cpp =================================================================== --- trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-01 18:34:17 UTC (rev 10551) +++ trunk/octave-forge/extra/lssa/fastlscomplex.cpp 2012-06-01 22:37:39 UTC (rev 10552) @@ -17,15 +17,113 @@ #include <octave/oct.h> #include <complex> #include <string> -#include <stdio> +#include <iostream> +typedef struct +{ double x, t; +} XTElem; + +inline double sqr(double x) +{ return x*x; } + #define PNUM 12; -// In order to reduce memory overhead, fastlscomplex will use a reference as the last input. +#define SETXT(p_, op_, x_, t_) (p_)->x op_ x_; (p_++)->t op_ t_; +#define SETT(p_, op_, x_, t_) *p_++ op_ t_; +#define SETX(p_, op_, x_, t_) *p_++ op_ x_; +/* h is a complex aux. variable; it is used for assignment times I everywhere */ +#define SETIX(p_, op_, x_, t_) h = x_; RE(*(p_)) op_ -IM(h); IM(*(p_)) op_ RE(h); p_++; -void fastlscomplex ( RowVector tvals, ComplexRowVector xvals, int n, - double length, int ncoeff, int noctaves, double omegamax, - ComplexRowVector& result ) { + /* Macro that sums up the power series terms into the power series + * element record pointed to by p_. + * By using = and += for op_, initial setting and accumulation can be selected. + * t_ is the expression specifying the abscissa value. set_ can be either + * SETXT to set the x and t fields of an XTElem record, or SETT/SETX to set + * the elements of a Real array representing alternately real and imaginary + * values. + */ + // -10 points, comments don't match method. +#define EXP_IOT_SERIES(p_, el_, t_, op_, setr_, seti_) \ +{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ + tt = -t; seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(-1.0/3.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(-1.0/5.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(-1.0/7.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(-1.0/9.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(-1.0/11.0); seti_(p_, op_, x*tt, tt) \ +} + +/* same as the above, but without alternating signs */ +#define EXPIOT_SERIES(p_, el_, t_, op_, setr_, seti_) \ +{ Real t = t_, tt; p_ = el_; setr_(p_, op_, x, 1) \ + seti_(p_, op_, x*t, t ) \ + tt = t*t*(1.0/2.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/3.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/4.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/5.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/6.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/7.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/8.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/9.0); seti_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/10.0); setr_(p_, op_, x*tt, tt) \ + tt *= t*(1.0/11.0); seti_(p_, op_, x*tt, tt) \ +} + +# define SRCARG Real *tptr, Real *xptr, int n, double *lengthptr +# define SRCVAR int k; Real length = *lengthptr; +# define SRCT tptr[k] +# define SRCX xptr[k] +# define SRCFIRST k = 0 +# define SRCAVAIL (k<n) +# define SRCNEXT k++ + +#define FASTLSCOMPLEX_HELP \ + "Takes the fast complex least-squares transform of irregularly sampled data.\n\ +When called, takes a time series with associated x-values, a number of octaves,\n\ +a number of coefficients per octave, the maximum frequency, and a length term.\n\ +It returns a complex row vector of the transformed series." + +DEFUN_DLD(fastlscomplex, args, nargout, FASTLSCOMPLEX_HELP ){ + if ( args.length != 7 ) { + std::cout << "Improper number of arguments!" << std::endl; + std::cout << FASTLSCOMPLEX_HELP << std::endl; + return octave_value_list (); + } + RowVector tvals = args(0).row_vector_value(); + ComplexRowVector xvals = args(1).complex_row_vector_value(); + double length = ( tvals(tvals.numel()--) - tvals( octave_idx_type (0))); + int ncoeff = args(4).int_value(); + int noctaves = args(5).int_value(); + double omegamax = args(6).double_value(); + + if ( error_state ) { + //print an error statement, see if error_state can be printed or interpreted. + return octave_value_list (); + } + if ( tvals.length() != xvals.length() ) { + std::cout << "Unmatched time series elements." << std::endl; + // return error state if possible? + return octave_value_list (); + } + octave_idx_type n = tvals.numel(); /* I think this will actually return octave_idx_type. + * In that case I'll change the signature of fastlscomplex. */ + if ( ( ncoeff == 0 ) || ( noctaves == 0 ) ) { + std::cout << "Cannot operate over either zero octaves or zero coefficients." << std::endl; + // return error state if possible + return octave_value_list (); + } + // possibly include error for when omegamax isn't right? + ComplexRowVector results = fastlscomplex(tvals, xvals, n, length, ncoeff, noctaves, omegamax); + return results; +} + +ComplexRowVector fastlscomplex ( RowVector tvals, ComplexRowVector xvals, octave_idx_type n, + double length, int ncoeff, int noctaves, double omegamax ) { /* Singly-linked list which contains each precomputation record * count stores the number of samples for which power series elements * were added. This will be useful for accelerating computation by ignoring @@ -33,7 +131,7 @@ */ struct SumVec { struct SumVec *next; - Complex<double> elements[PNUM]; + std::complex<double> elements[PNUM]; int count; } *shead, *stail, *sp, *sq; //Names have been kept from the C, may change if I want to. double dtelems[PNUM], /* power series elements of exp(-i dtau) */ @@ -43,7 +141,7 @@ tau0d, /* tau_h of first summand range at level d */ dtau = (0.5*M_PI)/omegamax,/* initial precomputation interval radius */ dtaud, /* precomputation interval radius at d'th merging step */ - n_1 = 1.0/n, /* reciprocal of sample count */ + n_1 = 1.0/n, /* reciprocal of sample count */ // n is implicitly cast to double. ooct, o, omul, /* omega/mu for octave's top omega and per band, mult. factor */ omegaoct, omega, /* Max. frequency of octave and current frequency */ on_1, /* n_1*(omega/mu)^p, n_1*(2*omega/mu)^p */ @@ -55,7 +153,7 @@ *eop, *oep, *ep, *op, *p, *q, *pe; - + ComplexRowVector results (ncoeff*noctaves); //This *should* be the right size. I'll check against later sections of code. int i , j; // Counters; coefficient and octave, respectively. // probable doubles: zetar, zetai, zzr, zzi, er, ei, emulr, emuli, @@ -97,6 +195,7 @@ omegaoct = omegamax; tau0d = tau0; dtaud = dtau; + octave_idx_type iter ( 0 ); // Looping over octaves for ( j = noctave ; ; ooct *= 0.5 , omegaoct *= 0.5 , tau0d += dtaud , dtaud *= 2 ) { // Looping over&results per frequency @@ -107,13 +206,14 @@ // sets real, imag parts of emul for( zeta = 0 , sp = shead; sp; sp = sp->next, e *= emul ) { if ( sp->count ) { - for ( zz = std::complex<double>(0.0,0.0) , p = sp->elems , pe = p + PNUM , on_1 = n_1 ; p < pe ; ) { + for ( zz = std::complex<double>(0.0,0.0) , octave_idx_type i (0) , p = sp->elems , pe = p + PNUM , on_1 = n_1 ; p < pe ; i++ ) { zz += *p++ * on_1; on_1 *= 0; } zeta += e * zz; } - *rp++ = zeta; + results(iter) = std::complex<double>(zeta.real(), zeta.imag()); + iter++; } if ( --j <= 0 ) break; //Avoids excess merging This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-01 18:34:24
|
Revision: 10551 http://octave.svn.sourceforge.net/octave/?rev=10551&view=rev Author: jpicarbajal Date: 2012-06-01 18:34:17 +0000 (Fri, 01 Jun 2012) Log Message: ----------- statistics: fixed bug in histfit. It was calling an old function. Modified Paths: -------------- trunk/octave-forge/main/statistics/inst/histfit.m Modified: trunk/octave-forge/main/statistics/inst/histfit.m =================================================================== --- trunk/octave-forge/main/statistics/inst/histfit.m 2012-06-01 12:42:55 UTC (rev 10550) +++ trunk/octave-forge/main/statistics/inst/histfit.m 2012-06-01 18:34:17 UTC (rev 10551) @@ -21,7 +21,7 @@ ## @code{histfit (@var{data}, @var{nbins})} plots a histogram of the values in ## the vector @var{data} using @var{nbins} bars in the histogram. With one input ## argument, @var{nbins} is set to the square root of the number of elements in -## data. +## data. ## ## Example ## @@ -61,10 +61,9 @@ sr = nanstd(data); ## Estimates the parameter, SIGMA, of the normal distribution. x=(-3*sr+mr:0.1*sr:3*sr+mr)';## Evenly spaced samples of the expected data range. [xb,yb] = bar(xbin,n); - y = normal_pdf(x,mr,sr.^2); + y = normpdf(x,mr,sr); binwidth = xbin(2)-xbin(1); y = row*y*binwidth; ## Normalization necessary to overplot the histogram. plot(xb,yb,";;b",x,y,";;r-"); ## Plots density line over histogram. endfunction - This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <jpi...@us...> - 2012-06-01 12:43:06
|
Revision: 10550 http://octave.svn.sourceforge.net/octave/?rev=10550&view=rev Author: jpicarbajal Date: 2012-06-01 12:42:55 +0000 (Fri, 01 Jun 2012) Log Message: ----------- signal: edit to xcorr.m and xcov.m texinfo Modified Paths: -------------- trunk/octave-forge/main/signal/inst/xcorr.m trunk/octave-forge/main/signal/inst/xcov.m Modified: trunk/octave-forge/main/signal/inst/xcorr.m =================================================================== --- trunk/octave-forge/main/signal/inst/xcorr.m 2012-05-31 13:51:53 UTC (rev 10549) +++ trunk/octave-forge/main/signal/inst/xcorr.m 2012-06-01 12:42:55 UTC (rev 10550) @@ -15,60 +15,106 @@ ## You should have received a copy of the GNU General Public License along with ## this program; if not, see <http://www.gnu.org/licenses/>. -## usage: [R, lag] = xcorr (X [, Y] [, maxlag] [, scale]) +## -*- texinfo -*- +## @deftypefn {Function File} {[@var{R}, @var{lag}] =} xcorr ( @var{X} ) +## @deftypefnx {Function File} {@dots{} =} xcorr ( @var{X}, @var{Y} ) +## @deftypefnx {Function File} {@dots{} =} xcorr ( @dots{}, @var{maxlag}) +## @deftypefnx {Function File} {@dots{} =} xcorr ( @dots{}, @var{scale}) +## Estimates the cross-correlation. ## -## Estimate the cross correlation R_xy(k) of vector arguments X and Y -## or, if Y is omitted, estimate autocorrelation R_xx(k) of vector X, -## for a range of lags k specified by argument "maxlag". If X is a -## matrix, each column of X is correlated with itself and every other +## Estimate the cross correlation R_xy(k) of vector arguments @var{X} and @var{Y} +## or, if @var{Y} is omitted, estimate autocorrelation R_xx(k) of vector @var{X}, +## for a range of lags k specified by argument "maxlag". If @var{X} is a +## matrix, each column of @var{X} is correlated with itself and every other ## column. ## ## The cross-correlation estimate between vectors "x" and "y" (of -## length N) for lag "k" is given by -## R_xy(k) = sum_{i=1}^{N}{x_{i+k} conj(y_i), -## where data not provided (for example x(-1), y(N+1)) is zero. +## length N) for lag "k" is given by ## -## ARGUMENTS -## X [non-empty; real or complex; vector or matrix] data +## @iftex +## @tex +## $$ R_{xy}(k) = \sum_{i=1}^{N} x_{i+k} \conj(y_i), +## @end tex +## @end iftex +## @ifnottex +## @example +## N +## R_xy(k) = sum x_@{i+k@} conj(y_i), +## i=1 +## @end example +## @end ifnottex ## -## Y [real or complex vector] data -## If X is a matrix (not a vector), Y must be omitted. -## Y may be omitted if X is a vector; in this case xcorr -## estimates the autocorrelation of X. +## where data not provided (for example x(-1), y(N+1)) is zero. Note the +## definition of cross-correlation given above. To compute a +## cross-correlation consistent with the field of statistics, see @command{xcov}. ## -## maxlag [integer scalar] maximum correlation lag -## If omitted, the default value is N-1, where N is the -## greater of the lengths of X and Y or, if X is a matrix, -## the number of rows in X. +## @strong{ARGUMENTS} +## @table @var +## @item X +## [non-empty; real or complex; vector or matrix] data ## -## scale [character string] specifies the type of scaling applied -## to the correlation vector (or matrix). is one of: -## 'none' return the unscaled correlation, R, -## 'biased' return the biased average, R/N, -## 'unbiased' return the unbiassed average, R(k)/(N-|k|), -## 'coeff' return the correlation coefficient, R/(rms(x).rms(y)), -## where "k" is the lag, and "N" is the length of X. -## If omitted, the default value is "none". -## If Y is supplied but does not have the ame length as X, -## scale must be "none". -## -## RETURNED VARIABLES -## R array of correlation estimates -## lag row vector of correlation lags [-maxlag:maxlag] +## @item Y +## [real or complex vector] data ## -## The array of correlation estimates has one of the following forms. -## (1) Cross-correlation estimate if X and Y are vectors. -## (2) Autocorrelation estimate if is a vector and Y is omitted, -## (3) If X is a matrix, R is an matrix containing the cross- -## correlation estimate of each column with every other column. -## Lag varies with the first index so that R has 2*maxlag+1 -## rows and P^2 columns where P is the number of columns in X. -## If Rij(k) is the correlation between columns i and j of X -## R(k+maxlag+1,P*(i-1)+j) == Rij(k) -## for lag k in [-maxlag:maxlag], or -## R(:,P*(i-1)+j) == xcorr(X(:,i),X(:,j)). -## "reshape(R(k,:),P,P)" is the cross-correlation matrix for X(k,:). +## If @var{X} is a matrix (not a vector), @var{Y} must be omitted. +## @var{Y} may be omitted if @var{X} is a vector; in this case xcorr +## estimates the autocorrelation of @var{X}. ## +## @item maxlag +## [integer scalar] maximum correlation lag +## If omitted, the default value is N-1, where N is the +## greater of the lengths of @var{X} and @var{Y} or, if @var{X} is a matrix, +## the number of rows in @var{X}. +## +## @item scale +## [character string] specifies the type of scaling applied +## to the correlation vector (or matrix). is one of: +## @table @samp +## @item none +## return the unscaled correlation, R, +## @item biased +## return the biased average, R/N, +## @item unbiased +## return the unbiassed average, R(k)/(N-|k|), +## @item coeff +## return the correlation coefficient, R/(rms(x).rms(y)), +## where "k" is the lag, and "N" is the length of @var{X}. +## If omitted, the default value is "none". +## If @var{Y} is supplied but does not have the same length as @var{X}, +## scale must be "none". +## @end table +## @end table +## +## @strong{RETURNED VARIABLES} +## @table @var +## @item R +## array of correlation estimates +## @item lag +## row vector of correlation lags [-maxlag:maxlag] +## @end table +## +## The array of correlation estimates has one of the following forms: +## (1) Cross-correlation estimate if @var{X} and @var{Y} are vectors. +## +## (2) Autocorrelation estimate if is a vector and @var{Y} is omitted. +## +## (3) If @var{X} is a matrix, R is an matrix containing the cross-correlation +## estimate of each column with every other column. Lag varies with the first +## index so that R has 2*maxlag+1 rows and P^2 columns where P is the number of +## columns in @var{X}. +## +## If Rij(k) is the correlation between columns i and j of @var{X} +## +## @code{R(k+maxlag+1,P*(i-1)+j) == Rij(k)} +## +## for lag k in [-maxlag:maxlag], or +## +## @code{R(:,P*(i-1)+j) == xcorr(X(:,i),X(:,j))}. +## +## @code{reshape(R(k,:),P,P)} is the cross-correlation matrix for @code{X(k,:)}. +## +## @seealso{xcov} +## @end deftypefn ## The cross-correlation estimate is calculated by a "spectral" method ## in which the FFT of the first vector is multiplied element-by-element @@ -89,7 +135,7 @@ ## ( hankel(x(1:k),x(k:N-k)) * y ) ./ N function [R, lags] = xcorr (X, Y, maxlag, scale) - + if (nargin < 1 || nargin > 4) print_usage; endif @@ -110,7 +156,7 @@ endif ## assign defaults to missing arguments - if isvector(X) + if isvector(X) ## if isempty(Y), Y=X; endif ## this line disables code for autocorr'n N = max(length(X),length(Y)); else @@ -121,7 +167,7 @@ ## check argument values if isempty(X) || isscalar(X) || ischar(Y) || ! ismatrix(X) - error("xcorr: X must be a vector or matrix"); + error("xcorr: X must be a vector or matrix"); endif if isscalar(Y) || ischar(Y) || (!isempty(Y) && !isvector(Y)) error("xcorr: Y must be a vector"); @@ -130,7 +176,7 @@ error("xcorr: X must be a vector if Y is specified"); endif if !isscalar(maxlag) || !isreal(maxlag) || maxlag<0 || fix(maxlag)!=maxlag - error("xcorr: maxlag must be a single non-negative integer"); + error("xcorr: maxlag must be a single non-negative integer"); endif ## ## sanity check on number of requested lags @@ -141,7 +187,7 @@ if (maxlag > N-1) pad_result = maxlag - (N - 1); maxlag = N - 1; - %error("xcorr: maxlag must be less than length(X)"); + %error("xcorr: maxlag must be less than length(X)"); else pad_result = 0; endif @@ -149,15 +195,15 @@ !strcmp(scale,'none') error("xcorr: scale must be 'none' if length(X) != length(Y)") endif - + P = columns(X); M = 2^nextpow2(N + maxlag); - if !isvector(X) + if !isvector(X) ## For matrix X, correlate each column "i" with all other "j" columns R = zeros(2*maxlag+1,P^2); ## do FFTs of padded column vectors - pre = fft (postpad (prepad (X, N+maxlag), M) ); + pre = fft (postpad (prepad (X, N+maxlag), M) ); post = conj (fft (postpad (X, M))); ## do autocorrelations (each column with itself) @@ -181,7 +227,7 @@ post = fft( postpad(X(:),M) ); cor = ifft( post .* conj(post) ); R = [ conj(cor(maxlag+1:-1:2)) ; cor(1:maxlag+1) ]; - else + else ## compute cross-correlation of X and Y ## If one of X & Y is a row vector, the other can be a column vector. pre = fft( postpad( prepad( X(:), length(X)+maxlag ), M) ); @@ -193,7 +239,7 @@ ## if inputs are real, outputs should be real, so ignore the ## insignificant complex portion left over from the FFT if isreal(X) && (isempty(Y) || isreal(Y)) - R=real(R); + R=real(R); endif ## correct for bias @@ -218,7 +264,7 @@ elseif !strcmp(scale, 'none') error("xcorr: scale must be 'biased', 'unbiased', 'coeff' or 'none'"); endif - + ## Pad result if necessary ## (most likely is not required, use "if" to avoid uncessary code) ## At this point, lag varies with the first index in R; @@ -229,7 +275,7 @@ endif ## Correct the shape (transpose) so it is the same as the first input vector if isvector(X) && P > 1 - R = R.'; + R = R.'; endif ## return the lag indices if desired @@ -241,7 +287,7 @@ endfunction ##------------ Use brute force to compute the correlation ------- -##if !isvector(X) +##if !isvector(X) ## ## For matrix X, compute cross-correlation of all columns ## R = zeros(2*maxlag+1,P^2); ## for i=1:P @@ -258,7 +304,7 @@ ##elseif isempty(Y) ## ## reshape X so that dot product comes out right ## X = reshape(X, 1, N); -## +## ## ## compute autocorrelation for 0:maxlag ## R = zeros (2*maxlag + 1, 1); ## for k=0:maxlag @@ -271,7 +317,7 @@ ## ## reshape and pad so X and Y are the same length ## X = reshape(postpad(X,N), 1, N); ## Y = reshape(postpad(Y,N), 1, N)'; -## +## ## ## compute cross-correlation ## R = zeros (2*maxlag + 1, 1); ## R(maxlag+1) = X*Y; Modified: trunk/octave-forge/main/signal/inst/xcov.m =================================================================== --- trunk/octave-forge/main/signal/inst/xcov.m 2012-05-31 13:51:53 UTC (rev 10549) +++ trunk/octave-forge/main/signal/inst/xcov.m 2012-06-01 12:42:55 UTC (rev 10550) @@ -13,24 +13,42 @@ ## You should have received a copy of the GNU General Public License along with ## this program; if not, see <http://www.gnu.org/licenses/>. -## usage: [c, lag] = xcov (X [, Y] [, maxlag] [, scale]) -## +## -*- texinfo -*- +## @deftypefn {Function File} {[@var{R}, @var{lag}] =} xcov ( @var{X} ) +## @deftypefnx {Function File} {@dots{} =} xcov ( @var{X}, @var{Y} ) +## @deftypefnx {Function File} {@dots{} =} xcov ( @dots{}, @var{maxlag}) +## @deftypefnx {Function File} {@dots{} =} xcov ( @dots{}, @var{scale}) ## Compute covariance at various lags [=correlation(x-mean(x),y-mean(y))]. ## -## X: input vector -## Y: if specified, compute cross-covariance between X and Y, +## @table @var +## @item X +## input vector +## @item Y +## if specified, compute cross-covariance between X and Y, ## otherwise compute autocovariance of X. -## maxlag: is specified, use lag range [-maxlag:maxlag], +## @item maxlag +## is specified, use lag range [-maxlag:maxlag], ## otherwise use range [-n+1:n-1]. -## Scale: -## 'biased' for covariance=raw/N, -## 'unbiased' for covariance=raw/(N-|lag|), -## 'coeff' for covariance=raw/(covariance at lag 0), -## 'none' for covariance=raw -## 'none' is the default. +## @item scale: +## @table @samp +## @item biased +## for covariance=raw/N, +## @item unbiased +## for covariance=raw/(N-|lag|), +## @item coeff +## for covariance=raw/(covariance at lag 0), +## @item none +## for covariance=raw +## @item none +## is the default. +## @end table +## @end table ## -## Returns the covariance for each lag in the range, plus an +## Returns the covariance for each lag in the range, plus an ## optional vector of lags. +## +## @seealso{xcorr} +## @end deftypefn function [retval, lags] = xcov (X, Y, maxlag, scale) @@ -58,5 +76,5 @@ else [retval, lags] = xcorr(center(X), maxlag, scale); endif - + endfunction This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <par...@us...> - 2012-05-31 13:52:04
|
Revision: 10549 http://octave.svn.sourceforge.net/octave/?rev=10549&view=rev Author: paramaniac Date: 2012-05-31 13:51:53 +0000 (Thu, 31 May 2012) Log Message: ----------- control-devel: revert example Modified Paths: -------------- trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m Modified: trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m =================================================================== --- trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-05-31 13:49:40 UTC (rev 10548) +++ trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-05-31 13:51:53 UTC (rev 10549) @@ -48,12 +48,9 @@ dat = iddata (Y, U) -[sys, x0, info] = n4sid (dat, 's', 10, 'n', 5) % s=10, n=5 +[sys, x0, info] = moen4 (dat, 's', 10, 'n', 5) % s=10, n=5 -L = chol (info.Ry, 'lower') -Be = info.K * L - [y, t] = lsim (sys, U, [], x0); err = norm (Y - y, 1) / norm (Y, 1) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |