-
Notifications
You must be signed in to change notification settings - Fork 13
/
Copy pathdrawBsplineSurf.m
111 lines (102 loc) · 3.21 KB
/
drawBsplineSurf.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
%
% Project: Approximation and Finite Elements in Isogeometric Problems
% Author: Luca Carlon
% Date: 2009.11.06
%
% Copyright (c) 2009-2021 Luca Carlon. All rights reserved.
%
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program. If not, see <http://www.gnu.org/licenses/>.
%
function [] = drawBsplineSurf(n, p, Xi, m, q, Eta, P, varargin)
% Draw the physical domain.
hold on;
box on;
grid on;
colormap('jet');
xlabel('x');
ylabel('y');
zlabel('z');
if nargin == 8
Q = P;
Q(:, :, 3) = varargin{1};
end
% Define the sample values.
spacexi = linspace(Xi(1), Xi(end), 100);
spaceeta = linspace(Eta(1), Eta(end), 100);
% Define the points to be drawn.
C = zeros(length(spacexi), length(spaceeta));
X = C;
Y = C;
Z = C;
for i = 1:length(spacexi)
for j = 1:length(spaceeta)
S = computeBsplineSurfPoint(n, p, Xi,...
m, q, Eta, P, spacexi(i), spaceeta(j));
X(i, j) = S(1);
Y(i, j) = S(2);
Z(i, j) = S(3);
if nargin == 8
T = computeBsplineSurfPoint(n, p, Xi,...
m, q, Eta, Q, spacexi(i), spaceeta(j));
C(i, j) = T(3);
end
end
end
if nargin == 8
surf(X, Y, Z, C);
else
surf(X, Y, Z);
end
shading interp;
clear S;
% Draw the control points.
for i = 1:length(P(:, 1, 1))
for j = 1:length(P(1, :, 1))
plot3(P(i, j, 1), P(i, j, 2),...
P(i, j, 3),...
'.', 'MarkerSize', 10, 'Color', 'red');
end
end
% Draw the lines.
for i = 1:length(P(:, 1, 1))
for j = 1:length(P(1, :, 1))
if i+1<=length(P(:, 1, 1))
plot3([P(i, j, 1), P(i+1, j, 1)],...
[P(i, j, 2), P(i+1, j, 2)],...
[P(i, j, 3), P(i+1, j, 3)],...
'--');
end
if j+1<=length(P(1, :, 1))
plot3([P(i, j, 1), P(i, j+1, 1)],...
[P(i, j, 2), P(i, j+1, 2)],...
[P(i, j, 3), P(i, j+1, 3)],...
'--');
end
end
end
% Draw the isocurves on Eta.
for j = 1:length(Eta)
for i = 1:length(spacexi)
S(i, :) = computeBsplineSurfPoint(n, p, Xi, m, q, Eta, P, spacexi(i), Eta(j));
end
plot3(S(:, 1), S(:, 2), S(:, 3), 'k');
end
% Draw the isocurves on Xi.
for i = 1:length(Xi)
for j = 1:length(spaceeta)
S(j, :) = computeBsplineSurfPoint(n, p, Xi, m, q, Eta, P, Xi(i), spaceeta(j));
end
plot3(S(:, 1), S(:, 2), S(:, 3), 'k');
end
endfunction