opening vtis in matlab/octave

plus a little example

post-processing/base64decode.m

@@ -0,0 +1,81 @@
+function y = base64decode(x)
+%BASE64DECODE Perform base64 decoding on a string.
+%
+%   BASE64DECODE(STR) decodes the given base64 string STR.
+%
+%   Any character not part of the 65-character base64 subset set is silently
+%   ignored.
+%
+%   This function is used to decode strings from the Base64 encoding specified
+%   in RFC 2045 - MIME (Multipurpose Internet Mail Extensions).  The Base64
+%   encoding is designed to represent arbitrary sequences of octets in a form
+%   that need not be humanly readable.  A 65-character subset ([A-Za-z0-9+/=])
+%   of US-ASCII is used, enabling 6 bits to be represented per printable
+%   character.
+%
+%   See also BASE64ENCODE.
+
+%   Author:      Peter J. Acklam
+%   Time-stamp:  2004-09-20 08:20:50 +0200
+%   E-mail:      [email protected]
+%   URL:         http://home.online.no/~pjacklam
+
+%   Modified by Guillaume Flandin, May 2008
+
+
+% Perform the following mapping
+%--------------------------------------------------------------------------
+%   A-Z  ->  0  - 25         a-z  ->  26 - 51         0-9  ->  52 - 61
+%   +    ->  62              /    ->  63              =    ->  64
+%   anything else -> NaN
+
+base64chars = NaN(1,256);
+base64chars('ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/=') = 0:64;
+x = base64chars(x);
+
+% Remove/ignore any characters not in the base64 characters list or '='
+%--------------------------------------------------------------------------
+
+x = x(~isnan(x));
+
+% Replace any incoming padding ('=' -> 64) with a zero pad
+%--------------------------------------------------------------------------
+
+if     x(end-1) == 64, p = 2; x(end-1:end) = 0;
+elseif x(end)   == 64, p = 1; x(end) = 0;
+else                   p = 0;
+end
+
+% Allocate decoded data array
+%--------------------------------------------------------------------------
+
+n = length(x) / 4;                               % number of groups
+x = reshape(uint8(x), 4, n);                     % input data
+y = zeros(3, n, 'uint8');                        % decoded data
+
+% Rearrange every 4 bytes into 3 bytes
+%--------------------------------------------------------------------------
+%    00aaaaaa 00bbbbbb 00cccccc 00dddddd
+%
+% to form
+%
+%    aaaaaabb bbbbcccc ccdddddd
+
+y(1,:) = bitshift(x(1,:), 2);                    % 6 highest bits of y(1,:)
+y(1,:) = bitor(y(1,:), bitshift(x(2,:), -4));    % 2 lowest bits of y(1,:)
+
+y(2,:) = bitshift(x(2,:), 4);                    % 4 highest bits of y(2,:)
+y(2,:) = bitor(y(2,:), bitshift(x(3,:), -2));    % 4 lowest bits of y(2,:)
+
+y(3,:) = bitshift(x(3,:), 6);                    % 2 highest bits of y(3,:)
+y(3,:) = bitor(y(3,:), x(4,:));                  % 6 lowest bits of y(3,:)
+
+% Remove any zero pad that was added to make this a multiple of 24 bits
+%--------------------------------------------------------------------------
+
+if p, y(end-p+1:end) = []; end
+
+% Reshape to a row vector
+%--------------------------------------------------------------------------
+
+y = reshape(y, 1, []);

post-processing/opening_vti_example.m

@@ -0,0 +1,14 @@
+% Script to decode a .vti (palabos output) into a matlab matrix
+% I gratefully acknowledge the help of Guillaume Flandin 
+
+vti_struct=xml2struct('rho_f1_001_00291000.vti'); %read output file
+vti_str=base64decode(vti_struct.VTKFile.ImageData.Piece.PointData.DataArray.Text);
+vti_no=typecast([0 0 vti_str],'single');
+vti_size=str2num(vti_struct.VTKFile.ImageData.Attributes.WholeExtent); %#ok<ST2NM>
+vti_x=vti_size(1)+vti_size(2)+1;
+vti_y=vti_size(3)+vti_size(4)+1;
+vti_z=vti_size(5)+vti_size(6)+1;
+
+
+%vti_matrix=reshape(vti_no(2:end),[3 vti_x vti_y vti_z]);  %for velocity
+vti_matrix=reshape(vti_no(2:end),[vti_x vti_y vti_z]);  %for density

post-processing/xml2struct.m

@@ -0,0 +1,184 @@
+function [ s ] = xml2struct( file )
+%Convert xml file into a MATLAB structure
+% [ s ] = xml2struct( file )
+%
+% A file containing:
+% <XMLname attrib1="Some value">
+%   <Element>Some text</Element>
+%   <DifferentElement attrib2="2">Some more text</Element>
+%   <DifferentElement attrib3="2" attrib4="1">Even more text</DifferentElement>
+% </XMLname>
+%
+% Will produce:
+% s.XMLname.Attributes.attrib1 = "Some value";
+% s.XMLname.Element.Text = "Some text";
+% s.XMLname.DifferentElement{1}.Attributes.attrib2 = "2";
+% s.XMLname.DifferentElement{1}.Text = "Some more text";
+% s.XMLname.DifferentElement{2}.Attributes.attrib3 = "2";
+% s.XMLname.DifferentElement{2}.Attributes.attrib4 = "1";
+% s.XMLname.DifferentElement{2}.Text = "Even more text";
+%
+% Please note that the following characters are substituted
+% '-' by '_dash_', ':' by '_colon_' and '.' by '_dot_'
+%
+% Written by W. Falkena, ASTI, TUDelft, 21-08-2010
+% Attribute parsing speed increased by 40% by A. Wanner, 14-6-2011
+% Added CDATA support by I. Smirnov, 20-3-2012
+%
+% Modified by X. Mo, University of Wisconsin, 12-5-2012
+
+    if (nargin < 1)
+        clc;
+        help xml2struct
+        return
+    end
+
+    if isa(file, 'org.apache.xerces.dom.DeferredDocumentImpl') || isa(file, 'org.apache.xerces.dom.DeferredElementImpl')
+        % input is a java xml object
+        xDoc = file;
+    else
+        %check for existance
+        if (exist(file,'file') == 0)
+            %Perhaps the xml extension was omitted from the file name. Add the
+            %extension and try again.
+            if (isempty(strfind(file,'.xml')))
+                file = [file '.xml'];
+            end
+
+            if (exist(file,'file') == 0)
+                error(['The file ' file ' could not be found']);
+            end
+        end
+        %read the xml file
+        xDoc = xmlread(file);
+    end
+
+    %parse xDoc into a MATLAB structure
+    s = parseChildNodes(xDoc);
+
+end
+
+% ----- Subfunction parseChildNodes -----
+function [children,ptext,textflag] = parseChildNodes(theNode)
+    % Recurse over node children.
+    children = struct;
+    ptext = struct; textflag = 'Text';
+    if hasChildNodes(theNode)
+        childNodes = getChildNodes(theNode);
+        numChildNodes = getLength(childNodes);
+
+        for count = 1:numChildNodes
+            theChild = item(childNodes,count-1);
+            [text,name,attr,childs,textflag] = getNodeData(theChild);
+
+            if (~strcmp(name,'#text') && ~strcmp(name,'#comment') && ~strcmp(name,'#cdata_dash_section'))
+                %XML allows the same elements to be defined multiple times,
+                %put each in a different cell
+                if (isfield(children,name))
+                    if (~iscell(children.(name)))
+                        %put existsing element into cell format
+                        children.(name) = {children.(name)};
+                    end
+                    index = length(children.(name))+1;
+                    %add new element
+                    children.(name){index} = childs;
+                    if(~isempty(fieldnames(text)))
+                        children.(name){index} = text; 
+                    end
+                    if(~isempty(attr)) 
+                        children.(name){index}.('Attributes') = attr; 
+                    end
+                else
+                    %add previously unknown (new) element to the structure
+                    children.(name) = childs;
+                    if(~isempty(text) && ~isempty(fieldnames(text)))
+                        children.(name) = text; 
+                    end
+                    if(~isempty(attr)) 
+                        children.(name).('Attributes') = attr; 
+                    end
+                end
+            else
+                ptextflag = 'Text';
+                if (strcmp(name, '#cdata_dash_section'))
+                    ptextflag = 'CDATA';
+                elseif (strcmp(name, '#comment'))
+                    ptextflag = 'Comment';
+                end
+
+                %this is the text in an element (i.e., the parentNode) 
+                if (~isempty(regexprep(text.(textflag),'[\s]*','')))
+                    if (~isfield(ptext,ptextflag) || isempty(ptext.(ptextflag)))
+                        ptext.(ptextflag) = text.(textflag);
+                    else
+                        %what to do when element data is as follows:
+                        %<element>Text <!--Comment--> More text</element>
+
+                        %put the text in different cells:
+                        % if (~iscell(ptext)) ptext = {ptext}; end
+                        % ptext{length(ptext)+1} = text;
+
+                        %just append the text
+                        ptext.(ptextflag) = [ptext.(ptextflag) text.(textflag)];
+                    end
+                end
+            end
+
+        end
+    end
+end
+
+% ----- Subfunction getNodeData -----
+function [text,name,attr,childs,textflag] = getNodeData(theNode)
+    % Create structure of node info.
+
+    %make sure name is allowed as structure name
+    name = toCharArray(getNodeName(theNode))';
+    name = strrep(name, '-', '_dash_');
+    name = strrep(name, ':', '_colon_');
+    name = strrep(name, '.', '_dot_');
+
+    attr = parseAttributes(theNode);
+    if (isempty(fieldnames(attr))) 
+        attr = []; 
+    end
+
+    %parse child nodes
+    [childs,text,textflag] = parseChildNodes(theNode);
+
+    if (isempty(fieldnames(childs)) && isempty(fieldnames(text)))
+        %get the data of any childless nodes
+        % faster than if any(strcmp(methods(theNode), 'getData'))
+        % no need to try-catch (?)
+        % faster than text = char(getData(theNode));
+        text.(textflag) = toCharArray(getTextContent(theNode))';
+
+    end
+
+end
+
+% ----- Subfunction parseAttributes -----
+function attributes = parseAttributes(theNode)
+    % Create attributes structure.
+
+    attributes = struct;
+    if hasAttributes(theNode)
+       theAttributes = getAttributes(theNode);
+       numAttributes = getLength(theAttributes);
+
+       for count = 1:numAttributes
+            %attrib = item(theAttributes,count-1);
+            %attr_name = regexprep(char(getName(attrib)),'[-:.]','_');
+            %attributes.(attr_name) = char(getValue(attrib));
+
+            %Suggestion of Adrian Wanner
+            str = toCharArray(toString(item(theAttributes,count-1)))';
+            k = strfind(str,'='); 
+            attr_name = str(1:(k(1)-1));
+            attr_name = strrep(attr_name, '-', '_dash_');
+            attr_name = strrep(attr_name, ':', '_colon_');
+            attr_name = strrep(attr_name, '.', '_dot_');
+            attributes.(attr_name) = str((k(1)+2):(end-1));
+       end
+    end
+end