Skip to content
Snippets Groups Projects
Commit 0b00e11f authored by Yoann Guilhem's avatar Yoann Guilhem Committed by Nicola Vigano
Browse files

Update and add some VTK utility functions 

git-svn-id: https://svn.code.sf.net/p/dct/code/trunk@947 4c865b51-4357-4376-afb4-474e03ccb993
parent 7ff9ef36
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
zUtil_VTK/gtVTKConvertTypeVTK2Matlab.m 0 → 100644
+ 50
− 0
View file @ 0b00e11f
function [mType, inFMT, inType] = gtVTKConvertTypeVTK2Matlab(vtkType)
% GTVTKCONVERTTYPEVTK2MATLAB Convert VTK datatype to Matlab datatype.
switch vtkType
case 'unsigned_char'
mType = 'uint8';
inFMT = '%d';
inType = 'uchar';
case 'unsigned_short'
mType = 'uint16';
inFMT = '%d';
inType = 'ushort';
case 'unsigned_int'
mType = 'uint32';
inFMT = '%d';
inType = 'uint';
case 'unsigned_long'
mType = 'uint64';
inFMT = '%d';
inType = 'uint64';
case 'char'
mType = 'int8';
inFMT = '%d';
inType = 'schar';
case 'short'
mType = 'int16';
inFMT = '%d';
inType = 'short';
case 'int'
mType = 'int32';
inFMT = '%d';
inType = 'int';
case 'long'
mType = 'int64';
inFMT = '%d';
inType = 'int64';
case 'float'
mType = 'single';
inFMT = '%f';
inType = 'float32';
case 'double'
mType = 'double';
inFMT = '%f';
inType = 'double';
otherwise
gtError('gtVTKConvertTypeVTK2Matlab:wrong_datatype', ...
['Unimplemented datatype "', vtkType, '", sorry...']);
end
end % end of function
zUtil_VTK/gtVTKMeshReader.m
+ 34
− 18
View file @ 0b00e11f
......@@ -14,6 +14,18 @@ function mesh = gtVTKMeshReader(filename, varargin)
%
% Version 001 20-10-2012 by YGuilhem yoann.guilhem@esrf.fr
%% Guessing output directory, filename and extension
[fpath, fname, fext] = fileparts(filename);
if isempty(fext)
fext = '.vtk';
end
% Trick to handle ~ path in unix
if ~isempty(fpath) && isunix && fpath(1)=='~'
fpath(1) = '';
fpath = fullfile(getenv('HOME'), fpath);
end
filename = fullfile(fpath, [fname fext]);
fid = fopen(filename, 'r');
if fid ==-1
gtError('gtVTKMeshReader:bad_input_file', 'Can''t open the file.');
......@@ -24,7 +36,7 @@ if ~strcmp(file.desc(3:5), 'vtk')
gtError('gtVTKMeshReader:bad_input_file', 'Not a valid VTK file.');
end
% Read header, mode and dataset
%% Read header, mode and dataset
file.header = fgets(fid);
str = fgets(fid);
file.mode = sscanf(str, '%s', 1);
......@@ -39,7 +51,7 @@ if ~strcmp(file.dataset, 'POLYDATA')
'This vtk file does not contain PolyData dataset!');
end
% Read number of vertices and coordinate datatype
%% Read number of vertices and coordinate datatype
str = fgets(fid);
nVert = sscanf(str, 'POINTS %d', 1);
file.nVertices = nVert;
......@@ -50,7 +62,8 @@ if binary
case {'float', 'double'}
coordFMT = coordType;
otherwise
true;
gtError('gtVTKMeshReader:wrong_coord_type', ...
['cannot handle coordinate with type: ' coordType]);
end
else
switch coordType
......@@ -63,41 +76,42 @@ else
end
end
% Read vertices
%% Read vertices
nCoord = 3*nVert;
if binary
coord = fread(fid, nCoord, ['3*' coordFMT '=>' coordFMT], 0 , 'ieee-be');
else
[coord, sz] = fscanf(fid, '%f %f %f', nCoord);
if sz ~= nCoord
warning('Problem in reading vertices.');
warning('gtVTKMeshReader:vertice_warning', ...
'Problem in reading vertices.');
end
end
mesh.vertices = reshape(coord, 3, nVert)';
% Read polygons
%% Read polygons
str = fgets(fid);
while str(1) ~= 'P' %&& str(1) ~= 'V'
str = fgets(fid);
end
if str(1) == 'P'
numbers = sscanf(str,'%*s %d %d', 2);
numbers = sscanf(str, '%*s %d %d', 2);
nFace = numbers(1);
nValue = numbers(2);
nNode = 3*nFace;
if binary
nbytes = 4;
type = 'uint32';
fseek(fid, nbytes, 0);
connectivity = fread(fid, nValue, '3*int32=>uint32', nbytes, 'ieee-be');
else
connectivity = uint32(fscanf(fid, '%*d %d %d %d\n', nNode));
%[connectivity, cnt] = fscanf(fid, '%*d %d %d %d\n', nNode);
%if cnt ~= nNode
% warning('Problem in reading faces.');
%end
%connectivity = uint32(connectivity);
[connectivity, cnt] = fscanf(fid, '%*d %d %d %d\n', nNode);
if cnt ~= nNode
warning('gtVTKMeshReader:face_warning', ...
'Problem in reading faces.');
end
connectivity = uint32(connectivity);
class(connectivity)
end
mesh.faces = reshape(connectivity+1, 3, nFace)';
else
......@@ -106,21 +120,23 @@ end
% read vertex ids
%if info == 'V'
% nv = sscanf(str,'%*s %d %*s', 1);
% nv = sscanf(str, '%*s %d %*s', 1);
% [A,cnt] = fscanf(fid,'%d %d \n', 2*nv);
% [A, cnt] = fscanf(fid, '%d %d \n', 2*nv);
% if cnt~=2*nv
% warning('Problem in reading faces.');
% warning('gtVTKMeshReader:vertex_warning', ...
% 'Problem in reading vertex.');
% end
% A = reshape(A, 2, cnt/2);
% faces = repmat(A(2,:)+1, 3, 1);
% faces = repmat(A(2, :)+1, 3, 1);
%end
%if str(1) ~= 'P' %&& str(1) ~= 'V'
% faces = 0;
%end
%% Close the file end terminate
fclose(fid);
end % end of function
zUtil_VTK/gtVTKMeshWriter.m
+ 19
− 16
View file @ 0b00e11f
......@@ -34,7 +34,7 @@ par.facefield = [];
par.facecolor = [];
par = parse_pv_pairs(par, varargin);
% Check coordinate data type
%% Check coordinate data type
switch par.coordType
case {'float', 'double'}
coordFMT = '%f %f %f\n';
......@@ -43,23 +43,23 @@ otherwise
['Unsupported coordinates datatype: ' par.coordType]);
end
% Get file parts and check extension
%% Get file parts and check extension
[fpath, fname, fext] = fileparts(filename);
if isempty(fext)
fext = '.vtk';
filename = fullfile(fpath, [fname fext]);
elseif ~strcmp(fext, '.vtk')
gtError('gtVTKMeshWriter:wrong_file_extension', ...
'Output file extension should be ''.vtk''!');
fext = '.vtk';
end
% Open file for writing
%% Open file for writing
fid = fopen(filename, 'w');
if fid ==-1
gtError('gtVTKMeshWriter:bad_output_file', 'Can''t open the file.');
end
% Header
%% Header
fprintf(fid, '# vtk DataFile Version 3.0\n');
fprintf(fid, '%s\n', par.header);
if par.binary
......@@ -69,20 +69,20 @@ else
end
fprintf(fid, 'DATASET POLYDATA\n');
% Points / Vertices
%% Points / Vertices
fprintf(fid, 'POINTS %d %s\n', size(mesh.vertices, 1), par.coordType);
if par.binary
count = fwrite(fid, mesh.vertices', par.coordType, 0, 'ieee-be');
%if count ~= size(mesh.vertices, 1)
if count ~= numel(mesh.vertices)
warning('Problem in writing vertices.');
warning('gtVTKMeshWriter:vertice_warning', ...
'Problem in writing vertices.');
end
fprintf(fid, '\n');
else
fprintf(fid, coordFMT, mesh.vertices');
end
% Polygons / Faces / Cells
%% Polygons / Faces / Cells
nFaces = size(mesh.faces, 1);
nValues = 4*nFaces;
fprintf(fid, 'POLYGONS %d %d\n', nFaces, nValues);
......@@ -90,14 +90,14 @@ output = [3*ones(1, nFaces) ; mesh.faces'-1];
if par.binary
count = fwrite(fid, output, 'uint32', 0, 'ieee-be');
if count ~= nValues
warning('Problem in writing faces.');
warning('gtVTKMeshWriter:face_warning', 'Problem in writing faces.');
end
fprintf(fid, '\n');
else
fprintf(fid, '%d %d %d %d\n', output);
end
% Face field
%% Face field
if ~isempty(par.facefield)
for ifield=1:length(par.facefield)
fField = par.facefield(ifield);
......@@ -127,7 +127,8 @@ if ~isempty(par.facefield)
if par.binary
count = fwrite(fid, fField.data, outType, 0, 'ieee-be');
if count ~= numel(fField.data)
warning(['Problem in writing field of face field ''' fField.name '']);
warning('gtVTKMeshWriter:face_warning', ...
['Problem in writing field of face field ''' fField.name '']);
end
fprintf(fid, '\n');
else
......@@ -140,7 +141,8 @@ if ~isempty(par.facefield)
if par.binary
count = fwrite(fid, uint8(255*output), 'uchar', 0, 'ieee-be');
if count ~= numel(output)
warning(['Problem in writing lookup table of face field ''' fField.name '']);
warning('gtVTKMeshWriter:lookup_table_warning', ...
['Problem in writing lookup table of face field ''' fField.name '']);
end
fprintf(fid, '\n');
else
......@@ -150,7 +152,7 @@ if ~isempty(par.facefield)
end
end
% Face color
%% Face color
if ~isempty(par.facecolor)
for ifield=1:length(par.facecolor)
fField = par.facecolor(ifield);
......@@ -171,7 +173,8 @@ if ~isempty(par.facecolor)
if par.binary
count = fwrite(fid, uint8(255*output), 'uchar', 0, 'ieee-be');
if count ~= numel(output)
warning(['Problem in writing field of face field ''' fField.name '']);
warning('gtVTKMeshWriter:face_field_warning', ...
['Problem in writing field of face field ''' fField.name '']);
end
fprintf(fid, '\n');
else
......@@ -180,7 +183,7 @@ if ~isempty(par.facecolor)
end
end
% Close the file end terminate
%% Close the file end terminate
fclose(fid);
end % end of function
zUtil_VTK/gtVTKStructGridWriter.m 0 → 100644
+ 131
− 0
View file @ 0b00e11f
function gtVTKStructGridWriter(vol, filename, varargin)
% GTVTKSTRUCTGRIDWRITER Write a volume in VTK file.
%
% gtVTKStructGridWriter(vol, filename, varargin)
% --------------------------------------------------------------------------
%
% INPUT:
% vol =
% filename = <string> Path to the output VTK file
%
% OPTIONAL INPUT:
% 'binary' = <bool> Write output in binary format (true by default)
% 'spacing' = <double> Vector defining spacing in each space dimension
% ([1 1 1] by default)
% 'origin' = <double> Vector defining volume origin
% ([0 0 0] by default)
% 'header' = <string> VTK file header description
% 'datatype' = <string> Data type to consider for input volume
% 'cmap' = <colormap>
%
%
% Version 001 20-10-2012 by YGuilhem yoann.guilhem@esrf.fr
par.spacing = ones(1, 3);
par.origin = zeros(1, 3);
par.dataname = 'Values';
par.header = 'VTK structured volume created by matlabDCT';
par.datatype = '';
par.cmap = [];
par.opacity = 1;
par.verbose = false;
par.binary = true;
par = parse_pv_pairs(par, varargin);
if iscell(vol)
error('Not implemented: only one volume can be handled by now!');
else
[nx, ny, nz] = size(vol);
end
ntot = nx*ny*nz;
if ~isempty(par.cmap) && size(par.cmap, 1) < max(max(max(vol)))
warning('gtVTKStructGridWriter:cmap_warning', ...
'You may have not enough values in your cmap...');
end
%% Guess input & output datatype
if isempty(par.datatype)
par.datatype = class(vol);
end
[outputTypeVTK, outputType, outputFMT] = gtVTKConvertTypeMatlab2VTK(par.datatype);
if par.verbose
disp(['Spacing : ' num2str(par.spacing)]);
disp(['Origin : ' num2str(par.origin)]);
disp(['Color map : ' num2str(par.cmap)]);
disp(['Dataname : ' par.dataname]);
disp(['Header : ' par.header]);
disp(['OutputType : ' outputType]);
disp(['OutputTypeVTK : ' outputTypeVTK]);
end
%% Open the file
fid = fopen(filename, 'w');
if fid ==-1
gtError('gtVTKStructGridWriter:bad_output_file', 'Can''t open the file.');
end
%% Write ASCII file header
fprintf(fid, '# vtk DataFile Version 3.0\n');
fprintf(fid, '%s\n', par.header);
if par.binary
fprintf(fid, 'BINARY\n');
else
fprintf(fid, 'ASCII\n');
end
fprintf(fid, 'DATASET STRUCTURED_POINTS\n');
fprintf(fid, 'DIMENSIONS %d %d %d\n', nx, ny, nz);
fprintf(fid, 'ORIGIN %4.3f %4.3f %4.3f\n', par.origin);
fprintf(fid, 'SPACING %4.3f %4.3f %4.3f\n', par.spacing);
fprintf(fid, 'POINT_DATA %d\n', ntot);
fprintf(fid, 'SCALARS %s %s\n', par.dataname, outputTypeVTK);
if ~isempty(par.cmap)
fprintf(fid, 'LOOKUP_TABLE map\n');
else
fprintf(fid, 'LOOKUP_TABLE default\n');
end
%% Write volume data
if par.binary
count = fwrite(fid, reshape(vol, 1, ntot), par.datatype, 0, 'ieee-be');
if count ~= numel(vol)
warning('gtVTKStructGridWriter:structgrid_volume_warning', ...
'Problem in writing structured grid volume.');
end
fprintf(fid, '\n');
else
fprintf(fid, [outputFMT '\n'], reshape(vol, 1, ntot));
end
%% Write lookup table if needed
if ~isempty(par.cmap)
mapSize = size(par.cmap);
if mapSize(2) == 3
output = [par.cmap'; par.opacity*ones(1, mapSize(1))];
elseif mapSize(2) == 4
output = par.cmap';
else
gtError('gtVTKStructGridWriter:wrong_cmap_size', ...
'Color map sould be N row vector of size 3 or 4!');
end
fprintf(fid, 'LOOKUP_TABLE map %d\n', mapSize(1));
if par.binary
count = fwrite(fid, uint8(225*output), 'uchar', 0, 'ieee-be');
if count ~= numel(output)
warning('gtVTKMeshWriter:cmap_warning', ...
'Problem in writing color map!');
end
fprintf(fid, '\n');
else
fprintf(fid, '%f %f %f %f\n', output);
end
end
%% Close file and terminate
fclose(fid);
end % end of function
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment