1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
|
function Dump2VTK(filename, fields, mesh, fieldname, varargin)
% Dump2VTK(filename, fields, mesh, fieldname, varargin)
%
% Dump fields extraced from an hdf5 file to a vtk file format
%
% possible arguments:
% 'NativeDump': 0 (default) / 1, dump in native coordinate system
% 'CloseAlpha': 0 (default) / 1, repeat first/last line in
% alpha-direction for a full cylindrical mesh
%
% example:
%
%
% openEMS matlab interface
% -----------------------
% author: Thorsten Liebig
%
% See also ReadHDF5FieldData ReadHDF5Mesh GetField_TD2FD GetField_Interpolation
NativeDump = 0;
CloseAlpha = 0;
for n=1:2:numel(varargin)
if (strcmp(varargin{n},'NativeDump')==1);
NativeDump = varargin{n+1};
elseif (strcmp(varargin{n},'CloseAlpha')==1);
CloseAlpha = varargin{n+1};
end
end
x = mesh.lines{1};
y = mesh.lines{2};
z = mesh.lines{3};
fid = fopen(filename,'w+');
% set nan values to zero
ind = find(isnan(fields));
if (~isempty(ind))
warning('openEMS:Dump2VTK','field contains nan, setting to zero');
fields(ind)=0;
end
% set inf values to zero
ind = find(isinf(fields));
if (~isempty(ind))
warning('openEMS:Dump2VTK','field contains inf, setting to zero');
fields(ind)=0;
end
if ((CloseAlpha~=0) && (mesh.type==1) && (range(y)<2*pi))
y(end+1) = y(1)+2*pi;
fields(:,end+1,:,:) = fields(:,1,:,:);
end
if (mesh.type==0) %write cartesian mesh to vtk
fprintf(fid,'# vtk DataFile Version 2.0\n');
fprintf(fid,'Rectilinear Grid by matlab-interface of openEMS\n');
fprintf(fid,'ASCII\n');
fprintf(fid,'DATASET RECTILINEAR_GRID\n');
fprintf(fid,'DIMENSIONS %d %d %d\n',numel(x),numel(y),numel(z));
fprintf(fid,'X_COORDINATES %d double\n',numel(x));
fprintf(fid,'%e',x(1));
for n=2:numel(x)
fprintf(fid,' %e',x(n));
end
fprintf(fid,'\n');
fprintf(fid,'Y_COORDINATES %d double\n',numel(y));
fprintf(fid,'%e',y(1));
for n=2:numel(y)
fprintf(fid,' %e',y(n));
end
fprintf(fid,'\n');
fprintf(fid,'Z_COORDINATES %d double\n',numel(z));
fprintf(fid,'%e',z(1));
for n=2:numel(z)
fprintf(fid,' %e',z(n));
end
elseif (mesh.type==1) %write cylindrical mesh to vtk
fprintf(fid,'# vtk DataFile Version 3.0\n');
fprintf(fid,'Structured Grid by matlab-interface of openEMS\n');
fprintf(fid,'ASCII\n');
fprintf(fid,'DATASET STRUCTURED_GRID\n');
fprintf(fid,'DIMENSIONS %d %d %d\n',numel(x),numel(y),numel(z));
fprintf(fid,'POINTS %d double\n',numel(x)*numel(y)*numel(z));
for nz=1:numel(z)
for ny=1:numel(y)
for nx=1:numel(x)
fprintf(fid,'%e %e %e\n',x(nx)*cos(y(ny)),x(nx)*sin(y(ny)),z(nz));
end
end
end
if ((ndims(fields)==4) && (NativeDump==0))
[R A Z] = ndgrid(x,y,z);
sinA = sin(A);
cosA = cos(A);
field_CC(:,:,:,1) = fields(:,:,:,1) .* cosA - fields(:,:,:,2) .* sinA;
field_CC(:,:,:,2) = fields(:,:,:,1) .* sinA + fields(:,:,:,2) .* cosA;
field_CC(:,:,:,3) = fields(:,:,:,3);
fields = field_CC;
clear R A Z sinA cosA field_CC
end
elseif (mesh.type==2) %write spherical mesh to vtk
fprintf(fid,'# vtk DataFile Version 3.0\n');
fprintf(fid,'Structured Grid by matlab-interface of openEMS\n');
fprintf(fid,'ASCII\n');
fprintf(fid,'DATASET STRUCTURED_GRID\n');
fprintf(fid,'DIMENSIONS %d %d %d\n',numel(x),numel(y),numel(z));
fprintf(fid,'POINTS %d double\n',numel(x)*numel(y)*numel(z));
for nz=1:numel(z)
for ny=1:numel(y)
for nx=1:numel(x)
fprintf(fid,'%e %e %e\n',...
x(nx)*sin(y(ny))*cos(z(nz)),...
x(nx)*sin(y(ny))*sin(z(nz)),...
x(nx)*cos(y(ny)));
end
end
end
if ((ndims(fields)==4) && (NativeDump==0))
[R T A] = ndgrid(x,y,z);
sinA = sin(A);
cosA = cos(A);
sinT = sin(T);
cosT = cos(T);
field_CC(:,:,:,1) = fields(:,:,:,1) .* sinT .* cosA + fields(:,:,:,2) .*cosT .* cosA - fields(:,:,:,3) .* sinA;
field_CC(:,:,:,2) = fields(:,:,:,1) .* sinT .* cosA + fields(:,:,:,2) .*cosT .* sinA + fields(:,:,:,3) .* cosA;
field_CC(:,:,:,3) = fields(:,:,:,1) .* cosT - fields(:,:,:,2) .*sinT;
fields = field_CC;
clear R A T sinA cosA sinT cosT field_CC
end
end
fprintf(fid,'\n\n');
fprintf(fid,'POINT_DATA %d\n',numel(x)*numel(y)*numel(z));
% dump vector field data
if (size(fields,4)>1)
if (nargin>3)
fprintf(fid,['VECTORS ' fieldname ' double\n']);
else
fprintf(fid,'VECTORS field double\n');
end
fclose(fid);
field_x = fields(:,:,:,1);
field_y = fields(:,:,:,2);
field_z = fields(:,:,:,3);
clear fields
dumpField(:,1) = field_x(:);
dumpField(:,2) = field_y(:);
dumpField(:,3) = field_z(:);
save('-ascii','-append',filename,'dumpField')
return
elseif (size(fields,4)==1) % scalar field
if (nargin>3)
fprintf(fid,['SCALARS ' fieldname ' double 1\nLOOKUP_TABLE default\n']);
else
fprintf(fid,'SCALARS field double 1\nLOOKUP_TABLE default\n');
end
fclose(fid);
dumpField = fields(:);
save('-ascii','-append',filename,'dumpField')
return
end
fclose(fid);
|
