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
|
function [mask weight] = mesh2vol(node, elem, xi, yi, zi)
%
% [mask weight]=mesh2vol(node,face,Nxyz)
% [mask weight]=mesh2vol(node,face,[Nx,Ny,Nz])
% [mask weight]=mesh2vol(node,face,xi,yi,zi,hf)
% or
% newval=mesh2vol(node_val,face,...)
%
% fast rasterization of a 3D mesh to a volume with tetrahedron index labels
%
% author: Qianqian Fang <q.fang at neu.edu>
% date for initial version: Feb 10,2014
%
% input:
% node: node coordinates, dimension N by 2 or N by 3 array
% nodeval: a 4-column array, the first 3 columns are the node coordinates,
% the last column denotes the values associated with each node
% face: a triangle surface, N by 3 or N by 4 array
% Nx,Ny,Nxy: output image in x/y dimensions, or both
% xi,yi: linear vectors for the output pixel center positions in x/y
% hf: the handle of a pre-created figure window for faster rendering
%
% output:
% mask: a 3D image, the value of each pixel is the index of the
% enclosing triangle, if the pixel is outside of the mesh, NaN
% weight: (optional) a 3 by Nx by Ny array, where Nx/Ny are the dimensions
% for the mask
% newval: when the node has 4 columns, the last column represents the
% values (color) at each node, the output newval is the rasterized
% mesh value map over the specified grid.
%
% note: This function only works for matlab
%
% example:
%
% [no,el]=meshgrid6(0:5,0:5,0:3);
% mask=mesh2vol(no,el(:,1:4),0:0.1:5,0:0.1:5,0:0.1:4);
% imagesc(mask(:,:,3))
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%
nodeval = [];
if (size(node, 2) == 4)
nodeval = node(:, 4);
node(:, 4) = [];
end
if (nargin == 3 && length(xi) == 1 && xi > 0)
mn = min(node);
mx = max(node);
df = (mx(1:min(3, size(node, 2))) - mn(1:min(3, size(node, 2)))) / xi;
elseif (nargin == 3 && length(xi) == 3 && all(xi > 0))
mn = min(node);
mx = max(node);
df = (mx(1:min(3, size(node, 2))) - mn(1:min(3, size(node, 2)))) ./ (xi(:)');
elseif (nargin == 5)
mx = [max(xi) max(yi) max(zi)];
mn = [min(xi) min(yi) min(zi)];
df = [min(diff(xi(:))) min(diff(yi(:))) min(diff(zi(:)))];
else
error('you must give at least xi input');
end
xi = mn(1):df(1):mx(1);
yi = mn(2):df(2):mx(2);
zi = mn(3):df(3):mx(3);
if (size(node, 2) ~= 3 || size(elem, 2) <= 3)
error('node must have 3 columns; face can not have less than 3 columns');
end
nz = length(zi);
mask = zeros(length(xi) - 1, length(yi) - 1, length(zi) - 1);
if (nargout > 1 || ~isempty(nodeval))
weight = zeros(4, length(xi) - 1, length(yi) - 1, length(zi) - 1);
end
hf = figure('visible', 'on');
for i = 1:nz
if (~isempty(nodeval))
[cutpos, cutvalue, facedata, elemid] = qmeshcut(elem, node, nodeval, ['z=' num2str(zi(i))]);
else
[cutpos, cutvalue, facedata, elemid] = qmeshcut(elem, node, node(:, 1), ['z=' num2str(zi(i))]);
end
if (isempty(cutpos))
continue
end
if (nargout > 1 || ~isempty(nodeval))
[maskz, weightz] = mesh2mask(cutpos, facedata, xi, yi, hf);
weight(:, :, :, i) = weightz;
else
maskz = mesh2mask(cutpos, facedata, xi, yi, hf);
end
idx = find(~isnan(maskz));
if (~isempty(nodeval))
eid = facedata(maskz(idx), :);
maskz(idx) = cutvalue(eid(:, 1)) .* weightz(1, idx)' + cutvalue(eid(:, 2)) .* weightz(2, idx)' + ...
cutvalue(eid(:, 3)) .* weightz(3, idx)' + cutvalue(eid(:, 4)) .* weightz(4, idx)';
else
maskz(idx) = elemid(maskz(idx));
end
mask(:, :, i) = maskz;
end
close(hf);
|
