forked from fangq/iso2mesh
-
Notifications
You must be signed in to change notification settings - Fork 0
/
mesh2mask.m
137 lines (130 loc) · 4.38 KB
/
mesh2mask.m
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
function [mask, weight]=mesh2mask(node,face,xi,yi,hf)
%
% [mask weight]=mesh2mask(node,face,Nxy)
% or
% [mask weight]=mesh2mask(node,face,[Nx,Ny])
% or
% [mask weight]=mesh2mask(node,face,xi,yi,hf)
%
% fast rasterization of a 2D mesh to an image with triangle index labels
%
% author: Qianqian Fang <q.fang at neu.edu>
% date for initial version: July 18,2013
%
% input:
% node: node coordinates, dimension N by 2 or N by 3 array
% 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: (optional) the handle of a pre-created figure window, for faster
% rendering
%
% output:
% mask: a 2D 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
%
% note: This function only works in MATLAB when the DISPLAY is not
% disabled. The maximum size of the mask output is limited by the
% screen size.
%
% example:
%
% [no,fc]=meshgrid6(0:5,0:5);
% [mask weight]=mesh2mask(no,fc,-1:0.1:5,0:0.1:5);
% imagesc(mask);
%
% -- this function is part of iso2mesh toolbox (http://iso2mesh.sf.net)
%
if(nargin==3 && length(xi)==1 && xi>0)
mn=min(node);
mx=max(node);
df=(mx(1:2)-mn(1:2))/xi;
elseif(nargin==3 && length(xi)==2 && all(xi>0))
mn=min(node);
mx=max(node);
df=(mx(1:2)-mn(1:2))./xi;
elseif(nargin==4 || nargin==5)
mx=[max(xi) max(yi)];
mn=[min(xi) min(yi)];
df=[min(diff(xi(:))) min(diff(yi(:)))];
else
error('you must give at least xi input');
end
if(size(node,2)<=1 || size(face,2)<=2)
error('node must have 2 or 3 columns; face can not have less than 2 columns');
end
if(nargin<5)
hf=figure('visible','on');
else
clf(hf);
end
patch('Vertices',node,'Faces',face,'linestyle','none','FaceColor','flat',...
'FaceVertexCData',(1:size(face,1))','CDataMapping', 'scaled');
set(gca, 'Position', [0 0 1 1]);
cm=jet(size(face,1));
colormap(cm);
axis off
set(gca,'xlim',[mn(1) mx(1)]);
set(gca,'ylim',[mn(2) mx(2)]);
set(gca,'clim',[1 size(face,1)]);
output_size = round((mx(1:2)-mn(1:2))./df);%Size in pixels
if(isoctavemesh || isempty(getenv('DISPLAY')))
resolution = 300; %Resolution in DPI
set(gcf,'PaperPositionMode','manual')
set(gcf,'paperunits','inches','paperposition',[0 0 output_size/resolution]);
deletemeshfile(mwpath('post_mesh2mask.png'));
print(mwpath('post_mesh2mask.png'),'-dpng',['-r' num2str(resolution)]);
mask=imread(mwpath('post_mesh2mask.png'));
else
pos=get(hf,'position');
pos(3:4)=max(pos(3:4),output_size+20);
set(hf,'position',pos);
set(gca, 'Units','pixels','position',[1, 1, output_size(1), output_size(2)]);
mask=getframe(gca);
if(any(size(mask.cdata)<[output_size([2 1]) 3]))
error('the requested rasterization grid is larger than the screen resolution');
end
mask=mask.cdata(1:output_size(2),1:output_size(1),:);
end
if(nargin<5)
close(hf);
end
mask=int32(reshape(mask,[size(mask,1)*size(mask,2) size(mask,3)]));
[isfound,locb]=ismember(mask,floor(cm*255),'rows');
locb(isfound==0)=nan;
mask=rot90(reshape(locb,output_size([2 1]))');
if(nargout>=2)
xi=mn(1)+df(1)/2:df(1):mx(1);
yi=mn(2)+df(2)/2:df(2):mx(2);
weight=barycentricgrid(node,face,xi,yi,mask);
if(size(face,2)>=4)
badidx=find(weight(1,:,:)<0 | weight(2,:,:)<0 | weight(3,:,:)<0);
badidx=badidx(face(mask(badidx),3)~=face(mask(badidx),4));
weight2=barycentricgrid(node,face(:,[1 3 4]),xi,yi,mask);
weight(:,badidx)=0;
weight([1 3 4],badidx)=weight2(:,badidx);
end
end
function weight=barycentricgrid(node,face,xi,yi,mask)
[xx,yy]=meshgrid(xi,yi);
idx=find(~isnan(mask));
eid=mask(idx);
t1=node(face(:,1),:);
t2=node(face(:,2),:);
t3=node(face(:,3),:);
tt=(t2(:,2)-t3(:,2)).*(t1(:,1)-t3(:,1))+(t3(:,1)-t2(:,1)).*(t1(:,2)-t3(:,2));
w(:,1)=(t2(eid,2)-t3(eid,2)).*(xx(idx)-t3(eid,1))+(t3(eid,1)-t2(eid,1)).*(yy(idx)-t3(eid,2));
w(:,2)=(t3(eid,2)-t1(eid,2)).*(xx(idx)-t3(eid,1))+(t1(eid,1)-t3(eid,1)).*(yy(idx)-t3(eid,2));
w(:,1)=w(:,1)./tt(eid);
w(:,2)=w(:,2)./tt(eid);
w(:,3)=1-w(:,1)-w(:,2);
weight=zeros(3,size(mask,1),size(mask,2));
ww=zeros(size(mask));
ww(idx)=w(:,1);
weight(1,:,:)=ww;
ww(idx)=w(:,2);
weight(2,:,:)=ww;
ww(idx)=w(:,3);
weight(3,:,:)=ww;