MATLAB非矩形区域上曲面的绘制

一、引言

Matlab提供了mesh和surf等绘制曲面图形的命令，在使用这些命令绘图的过程中，经常需要用meshgrid生成网格坐标矩阵，该命令要求横坐标和纵坐标的取值范围都是常数，这就导致了利用mesh和surf等命令绘制的曲面都是矩形区域内的图形。实际应用中经常需要绘制非矩形区域内的曲面。Matlab提供fmesh、fsurf、ezmesh和ezsurf等命令来绘制非矩形区域内的图形，但是这要求函数表达式必须已知。事实上，也可以利用mesh和surf等绘制非矩形区域内的曲面，只是需要配合给出非矩形区域内的坐标网格矩阵。
本文以球面、椭球面、降落伞、封闭的圆柱面等非矩形区域内的几何曲面为例，给出了一些比较复杂的曲面图形的绘制。

二、绘制球面

（一）绘制球面的方法有很多，常见的有：
1）使用matlab内嵌的sphere命令直接绘制球面；
2）使用fsurf绘制曲面图，或者使用fmesh绘制网格图。
3）使用mesh或者surf等绘制球面。
说明：sphere命令是由系统自动处理横坐标和纵坐标的网格点，然后绘制曲面或者网格面。fsurf和fmesh必须事先已知曲面的函数表达式。
而mesh和surf是由自己手工定义横坐标和纵坐标的网格点，为了方便定义非矩形区域内的网格坐标，使用了球面坐标系。
关于fsurf，fmesh，surf和mesh等命令的用法已经在其它博文中给出了详解，此处不再赘叙，仅给出sphere的常见用法。
（二）sphere的基本语法：
1） sphere
功能是绘制20x20带有网格单位球面图形
2） sphere( n )
功能是绘制由n×n曲面片构成的单位球面
3）[ X, Y, Z ] = sphere( n )
功能是获取单位球面的横坐标、纵坐标和竖坐标的网格坐标矩阵，之后利用mesh或者surf等绘图命令绘制相应的曲面。
3）绘制球面的例子
（三）下面逐一给出上述三类绘制球面的例子。
Demo1： 利用sphere命令绘制半径为1的单位球面。

%syntax 1
figure; 
subplot( 1, 2, 1 )
sphere;%绘制20x20网格单位球面
axis equal;
title( 'sphere' )
%syntax 2
subplot( 1, 2, 2 )
sphere( 40 ); %绘制40x40网格单位球面
axis equal;
title( 'sphere(40)' )
 %syntax 3
[ X, Y, Z ] = sphere( 30 );%获取单位球面的网格坐标
figure;
subplot( 1, 3, 1 )
mesh( X, Y, Z ) %利用mesh绘制X,Y,Z生成的网格球面
axis equal;
hidden off
title( 'mesh' )
subplot( 1, 3, 2 )
surf( X, Y, Z ) %利用surf绘制X,Y,Z生成的网格球面
axis equal;
title( 'surf' )
subplot( 1, 3, 3 )
plot3( X, Y, Z, 'b' ) %利用surf绘制X,Y,Z生成的网格球面
axis equal;
title( 'plot3' )
text( 0, 9, '[ X, Y, Z ] = sphere( 30 )' )
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

运行结果：


Demo2： 利用fsurf和fmesh命令分别绘制半径为1的单位球面。

fun = @(x,y)sqrt( 1 - x.^2 - y.^2 );
figure; 
subplot( 1, 2, 1 )
fmesh( fun )
hold on
fun1 = @(x,y)-sqrt( 1 - x.^2 - y.^2 );
fmesh( fun1 )
axis equal
title( 'fmesh:单位球面' )
%3)fsurf绘制球面
subplot( 1, 2, 2 )
fsurf( fun )
hold on
fsurf( fun1 )
axis equal
title( 'fsurf:单位球面' )
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16

运行结果：

Demo3： 利用surf和mesh命令分别绘制半径为1的单位球面。

%利用球面坐标
r = 1;
phi = [ 0 : 0.02 : 1 ] * pi;
theta = [ -1 : 0.05 : 1 ] * pi;
[P, T] = meshgrid( phi, theta );
X = r * sin( P ) .* cos( T );
Y = r * sin( P ) .* sin( T );
Z = r * cos( P );
figure;
subplot( 1, 2, 1 )
mesh( X, Y, Z )
axis equal
title( 'mesh : 球面坐标绘制球面' )
subplot( 1, 2, 2 )
surf( X, Y, Z )%surf绘制球面
axis equal
title( 'surf : 球面坐标绘制球面' )
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

运行结果：

三、绘制椭球面

（一）绘制椭球面常见的方法：
1）使用matlab内嵌的ellipsoid命令直接绘制椭球面；
2）使用fsurf绘制曲面图，或者使用fmesh绘制网格图。
3）使用mesh或者surf等绘制球面。
（二）ellipsoid的语法
[X,Y,Z] = ellipsoid(xc,yc,zc,a,b,c,n)
其中xc、xy和xz表示椭球面的中心坐标，a、b和c分别表示椭球面在三个坐标轴方向的半径，n表示网格经纬线的条数。
X,Y和Z分别表示椭球面网格坐标矩阵，后续可以使用mesh(X,Y,Z)或者surf(X,Y,Z)来绘制椭球面。
（三）绘制椭球面的示例
Demo1：利用ellipsoid绘制椭球面

clear all
clc
cx = 0;
cy = 0;
cz = 0;
a = 3;
b = 2;
c = 1;
subplot( 1, 2, 1 )
ellipsoid( cx, cy, cz, a, b, c );
axis equal
title( 'ellipsoid : ellipse' )
xlabel( 'x' )
ylabel( 'y' )
zlabel( 'z' )
subplot( 1, 2, 2 )
ellipsoid( cx, cy, cz, a, b, c, 40 );
axis equal
title( 'ellipsoid : ellipse' )
xlabel( 'x' )
ylabel( 'y' )
zlabel( 'z' )
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22

运行结果

Demo2：利用fmesh和fsurf分别绘制椭球面

clear all
clc
syms phi theta
a = 4;
b = 3;
c = 2;
X = a * sin( phi ) .* cos( theta );
Y = b * sin( phi ) .* sin( theta );
Z = c * cos( phi );
subplot( 1, 2, 1 );
fmesh( X, Y, Z, [0, pi, -pi, pi], 'EdgeColor', [ 0, 0, 1 ] );
axis equal
title( 'fmesh : Ellipse' )
subplot( 1, 2, 2 );
fsurf( X, Y, Z, [0, pi, -pi, pi] )
axis equal
title( 'fsurf : Ellipse' )
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

运行结果：

Demo3：利用mesh和surf分别绘制椭球面

clear all
clc
a = 4;
b = 3;
c = 2;
phi = [ 0 : 0.02 : 1 ] * pi;
theta = [ -1 : 0.05 : 1 ] * pi;
[P, T] = meshgrid( phi, theta );
X = a * sin( P ) .* cos( T );
Y = b * sin( P ) .* sin( T );
Z = c * cos( P );
subplot( 1, 2, 1 );
mesh( X, Y, Z, 'EdgeColor', [0, 1, 1] );
axis equal
subplot( 1, 2, 2 );
surf( X, Y, Z, 'FaceColor', [ 1, 0.2, 0 ], 'EdgeColor', [0, 1, 1] )
axis equal
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

运行结果：

四、绘制降落伞图——demo1（降落伞下悬挂六面体）

这属于多张去组合而成的曲面，需要分别绘制各自的曲面片，然后拼接在一起即可。需要注意的必须确保不同的曲面片之间要无缝拼接，因此对坐标网格矩阵要求就比较严格了。
1.绘制的思路：
1）绘制伞衣，利用surf和球面坐标，同时调整Z的取值范围；
2）绘制伞绳，利用plot3画直线即可；
3）伞绳悬挂的储物筐，利用surf绘制一个六面体。
2.绘制降落伞的代码：

clear
clc
r = 4;
%绘制降落伞伞衣
phi = [ 0 : 1/80 : 1/4 ] * pi;
theta = [ 0 : 1/20 : 2 ] * pi;
[ P, T ] = meshgrid( phi, theta );
X = r .* sin(P) .* cos( T );
Y = r .* sin(P) .* sin( T );
Z = 2 * r .* cos(P);
surf( X, Y, Z, 'EdgeColor', 'y', 'Facecolor', 'r' )
 
%绘制降落伞伞绳
hold on
lines = 30;
Xx = r * sin( pi/4 ) * cos( linspace( 0, 2*pi, lines ) );
Yy = r * sin( pi/4 ) * sin( linspace( 0, 2*pi, lines ) );
Zz = 2 * r * cos( pi/4 ) * ones( lines );
for i = 1 : lines
    plot3( [Xx(i), 0], [Yy(i),0], [Zz(i), 0], 'b' );
    hold on
end
axis equal
axis( [ -r, r, -r, r, -1, 2*r ] )
view( [2, 2, -1] )
axis off
 
 
%绘制降落伞下的悬挂物——储物筐（六面体）
hold on
[ X, Y ] = meshgrid( -0.5 : 0.5  );
Z1 = ones( size(X) )-1;
Z2 = zeros( size(X) )-1;
loc = 0.5;
surf( X, Y, Z1, 'FaceColor', [1,0,0] ); %上面
surf( X, Y, Z2, 'FaceColor', [0,0,0] ); %下面
 
surf( Y, Z1+loc, X-loc, 'FaceColor', [0.5, 1, 1] );%右面
surf( Y, Z2+loc, X-loc, 'FaceColor', [0.5, 1, 1] );%左面
surf( Z1+loc, X, Y-loc, 'FaceColor', [0.5, 1, 1] );%前面
surf( Z2+loc, X, Y-loc, 'FaceColor', [0.5, 1, 1] );%后面
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

3.运行结果：

五、绘制降落伞图——demo2（降落伞下悬挂类圆柱体）

clear
clc
r = 60;
%绘制降落伞伞衣
phi = [ 0 : 1/80 : 1/4 ] * pi;
theta = [ 0 : 1/20 : 2 ] * pi;
[ P, T ] = meshgrid( phi, theta );
X = r .* sin(P) .* cos( T );
Y = r .* sin(P) .* sin( T );
Z = 2 * r .* cos(P);
surf( X, Y, Z, 'EdgeColor', 'y', 'Facecolor', 'r' )

%绘制降落伞伞绳
hold on
lines = 30;
Xx = r * sin( pi/4 ) * cos( linspace( 0, 2*pi, lines ) );
Yy = r * sin( pi/4 ) * sin( linspace( 0, 2*pi, lines ) );
Zz = 2 * r * cos( pi/4 ) * ones( lines );
for i = 1 : lines
    plot3( [Xx(i), 0], [Yy(i),0], [Zz(i), 0], 'b' );
    hold on
end
axis equal
axis( [ -r, r, -r, r, -1, 2*r ] )
axis off

%绘制降落伞下的悬挂物——储物罐（类圆柱体）
%绘制上窄下宽的圆柱体
theta = [ 0 : 0.05 : 2 ] * pi;
z = ones( size( theta ) );
X = [];
Y = [];
Z = [];
for i = 19 : -1 : 0
    r = 10 - i / 10;
    X = [ X; r*cos(theta) ];
    Y = [ Y; r*sin(theta) ];
    Z = [ Z; (i-19)*z ];
end
surf( X, Y, Z, 'EdgeColor', 'None', 'FaceColor', [ 0.3, 0.6, 0.3 ] )

%绘制储物罐的大底(球冠)
r = 9.9;
syms the phi
fsurf( r*sin(phi)*cos(the), r*sin(phi)*sin(the), r*cos(phi)/5-19, [ pi/2, pi, 0, 2*pi ], 'FaceColor', [ 0.3, 0.3, 0.9 ], 'EdgeColor', 'None' );
axis equal
view( [ 1, 1, -0.2 ] )
%绘制储物罐的上盖(球冠)
R = 8;
fsurf( R*sin(phi)*cos(the), R*sin(phi)*sin(the), R*cos(phi)/5, [ 0, pi/2, 0, 2*pi ], 'FaceColor', [ 0.9, 0.3, 0.2 ], 'EdgeColor', 'None' );
axis equal
view( [2, 2, -1] ) %修改视角
%view( [2, 2, 1] ) %修改视角


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

运行结果：