常见离散信号的MATLAB产生和图形显示

实验一常见离散信号的MATLAB产生和图形显示

实验目的：加深对常用离散信号的理解;

实验原理：

1.单位抽样序列

“、Un=0

o(n)=«

[0

在MATLAB中可以利用zeros()函数实现。

x=zeros。,N);

XD=1；

如果6(〃)在时间轴上延迟了k个单位，得至1」6(〃-女)即:

1n=k

d(n-k)="

0〃w0

2.单位阶跃序列

1n>()

〃(〃)〈

0n<0

在MATLAB中可以利用ones()函数实现。

x=ones(l,N);

3.正弦序列

x(n)=Asin(2加/Fs+(p)

在MATLAB中

n=0:N-l

x=X*sin(2*pi*广nlFs+fid)

4.复正弦序列

x(n)=ejan

在MATLAB中

n=0:N-\

x=exp(j*

5.指数序列

x(n)=an

在MATLAB中

〃=O:N-1

x=a.^n

实验内容：编制程序产生上述5种信号(长度可输入确定)，并绘出其图形。

实验要求：讨论复指数序列的性质。

实验过程：

L单位冲击序列：

»n=0:10;

»xl=[lzeros(l,10)];

»x2=[zeros(l,8)1zeros(lz8)];

»subplot(l,2,l);

»stem(n,xl);

»xlabel('时间序列n');

»ylabe(幅度)；

单位冲激序列8(n),)；

»subplot(lf2,2);

»stem(x2);

»xlabel('时间序列n');

»ylabel('幅度');

>>title(,延时了8个单位的冲激序列6(n-8),);

单位冲激序列b(n)延时了8个单位的冲激序列6(n-8)

1

0.90.9

0808

0707

0.60.6

05

04

0.30.3

0202

01

0L0_____00.00000

051005101520

时间序列n时间序列n

2.单位阶跃序列:

>>n=0:10;

»u=(ones(lzll)］；

»stem(nzu);

>>xlabel('时间序列n');

»ylabe(信号幅度,)；

>>title(，单位阶跃序列u(n),)；

3.正弦序列：

»n=l:30;

»x=2*sin(pi*n/6+pi/4);

»stem(n,x);

»xlabel(时间序列n');

»ylabe(振幅)

»titled正弦函数序列x=2*sin(pi*n/6+pi/4),);

»

4.复指数序列:

>>n=l:30;

»x=5*exp(j*3*n);

»stem(nzx);

>>xlabel('时间序列n');

»ylabel('振幅。；

>>tit®，复指数序列x=5*exp(j*3*n)');

5.指数序列：

»n=l:30;

»x=1.8.An;

»stem(n,x);

»xlabel('时间序歹ijn')；

»ylabel。振幅，);

»title('指数序歹Ux=1.8.An');

复指数序列的周期性讨论:

为了研究复指数序列的周期性质，我们分别作了正弦函数xl=1.5sin(03n)

和x2=sin(0.6n);的幅度特性图像。由下图看出：xl=1.5sin(0.3冗n)的周期是20,

而x2=sin(0.6n)是非周期的。理论计算中对第一个，N=2n/(0.3n)=20/3,

为有理数。而第二个0.6不是兀的倍数，所以不是周期的。因此可以看出，实

验结果和理论相符。

对于离散复指数函数X=只有当z是纯虚数，且纯虚数的系数是兀的

倍数时，才是周期的。其它情况下均不是。

MATLAB实现过程：

»n=0:80;

»xl=1.5*sin(0.3*pi*n);

»x2=sin(0.6*n);

»subplotfl,2,1);

»stem(n,xl);

»xlabel('时间序列n');

»ylabel。振幅,)；

»title('正弦序列xl=1.5*sin(0.3*pi*n)');

»subplot(l,2,2);

»stem(n,x2);

>>xlabelC时间序列n')；

»ylabel('振幅’)；

»title('正弦序列xl=1.5*sin(0.3*pi*n)');

实验2离散系统的差分方程、冲激响应和卷积分析

实验目的：加深对离散系统的差分方程、冲激响应和卷积分析方法的理解。

实验原理：离散系统

x[n]--------------[----乂川

--------->Discrete-time

system

其输入、输出关系可用以下差分方程描述:

NM

zdky[n-k]=zpkx[n-k]

k=0Z=0

00

M川二―加--加]

输入信号分解为冲激信号，帕。o记系统单位冲激响应

况刈f[网，则系统响应为如下的卷积计算式：

y[n]=*h[n]=Yx[m]h[n-m]

当dk'k=1'2，…N时，h[n]是有限长度的(n：[0,M]),称系

统为FIR系统；反之，称系统为IIR系统。

在MATLAB中，可以用函数y-Filter(p,d,x)求解差分方程，也可以用函数

y=Conv(x,h)计算卷积。

实验内容：编制程序求解下列两个系统的单位冲激响应和阶跃响应，并绘出其

图形。

y[n]+0.6兄〃-1]+0.08乂〃-2]=——x[n-1]

y[n]=0.2{x[n-1]+-2]+x[n-3]+x[n-4]+-5]+-6]}

实验要求：给出理论计算结果和程序计算结果并讨论。

实验过程：

+0.6况〃一1]+0.08M〃-2J=x[ri]-x[n—1]

(i)单位冲激响应：

»a=[l,0.6,0.08];

»b=[l,-l];

»N=20;

»n=0:l:N;

»x=[l/zeros(l/N)];

»y=filter(b,a,x);

»stem(n,y);

>>xlabel(时间序列n');

>>ylabe(信号幅度。;

>>title(，单位冲激响应h(n)，);

»

单位冲激响应h(n)

IIIIII

Q

5

o@.

0Looo。。c-e—oo。oc

5o

o

2468101214161820

时间序列n

(2)单位阶跃响应:

»a=[l,0.6,0.08];

»b=[lrl]；

»N=20;

»n=0:l:N-l;

»x=[ones{l,N)];

»y=filter(b,a,x);

»stem(n,y);

»xlabel('时间序号)；

:*〉ylabel('信号巾鬲度')；

>>title(，单位阶跃响应h(n),);

单位阶跃响应h(n)

Q

OOQQO。。。。。。。。。

68101214161820

时间序号

理论分析：

-,

l-z

由差分方程得系统函数为：H(z)=------------：--------------

l+0.6z-,+0.08z-2

利用分部分式法可得：H(Z)=---r-----------r,Z反变换得：

l+0.4z-1l+0.2z-1

h(n)=[7•(-0.4)"-6-(—0.2〃)]w(n)

h(n)即为单位冲击响应。

对于阶跃响应，X(z)=，所以阶跃响应K(z)=X(z)•H(z)

1—z

_________1

-l+0.6z-l+0.08z-2

2_________1

l+0.4z-,-l+0.2z-1

所以，y(〃)=[2•(-0.4)”一(一0.2")]u(n)

可见，实验结果与理论分析是一致的。

=0.2(x[n-1J+M〃-2]+x[n-3]+x[n-4J+x[n-5]+-6]}

(1)单位冲激响应

»a=[l];

»b=[0.2/0.2/0.2/0.2/0.2/0.2];

»N=20;

»n=0:l:N;

»x=[l,zeros。,N)];

»y=filter(b,a,x);

»stem(n,y);

>>xlabel('时间序号n');

»ylabe(信号幅度,)；

“titled单位冲激响应h(n)');

单位冲激响应h(n)

o.

18

O._

-

O5.6

11

-4

1

51

2

^_

超1

婴1

0._

5-

Dlr8

四

^1

0.06

0.1

04

0_

02

0一

Jo

时间序号n

(2)单位阶跃响应

»a=l;

»b=[0.2/0.2,0.2,0.2/0.2/0.2];

»N=20;

»n=0:l:N-l:

»x=[ones(l,N)];

»y=filter(b,a,x);

〉〉stem(n,y);

»xlabel('时间序号,)；

»ylabelC信号幅度)；

>>titled单位阶跃响应h(n),)；

单位阶跃响应h(n)

68101214161820

时间序号

理论分析：

由差分方程得系统函数为：H(z)=0.2•(z~+z-2+z-3+zY+z-5+z«)

z反变换得；

h(n)=0.2•(b(〃-1)+S(n-2)+S(n-3)+S(n-4)+6(n-5)+8{n-6))

h(n)即为单位冲击响应。

对于阶跃响应，X(z)=」■7r，所以阶跃响应：

1—Z

Y(z)=X(z)H(z)

0.2-(z-1+Z-2+Z-3+z"+Z-5+Z-6)

=

所以，y(n)=0.2•(〃(〃-1)4-u(n-2)+u(n-3)+u(n-4)+u(n-5)4-u[n-6))

可见，实验结果与理论分析是一致的。

实验3离散系统的频率响应分析和零、极点分布

实验目的：加深对离散系统的频率响应分析和零、极点分布的概念理解。

实验原理：离散系统的时域方程为

£dky(n-=£Pkx{n-k)

k=0k=O

其变换域分析方法如下：

y[h]=*h[n]=f耳团]川〃-m]oY(ejM)=加)

频域

H(屋)-P(/)_〃o+i

系统的频率响应为Dd)do+V+...+

8

M〃]=x[n]*h[n]=-m]<=>Y(z)=X(z)H(z)

Z域2f

/、-1-M

H(z)=型1=Po+*+…+&z

系统的转移函数为D(z)do+4zT+...+"/'

MM,

W(z)=^——=K^---------

1

Ydkz-na-B22.

分解因式i=。i=l,其中'和"称为零、极

点。

在MATLAB中，可以用函数[z,p,K]-tf2zp(num,den)求得有理分式形式的

系统转移函数的零、极点，用函数zplane(z,p)绘出零、极点分布图；也可

以用函数zplane(num,den)直接绘出有理分式形式的系统转移函数的零、极

点分布图。

另外，在MATLAB中，可以用函数[r,p,k]=residuez(num,den)完成部

分分式展开计算；可以用函数sos=zp2sos(z,p,K)完成将高阶系统分解为2

阶系统的串联。

实验内容：求系统

0.0528+0.797z-,+0.1295z"+0.1295z-3+0.797Z-4+0.0528z-5

(l-1.8107z-1+2.4947z-2-1.8801Z-3+0.9537Z-4-0.2336z-5

的零、极点和幅度频率响应。

实验要求：编程实现系统参数输入，绘出幅度频率响应曲线刃零、极点分布图。

实验过程：

画出零极点图：

»num=[O.0528,0.797,0.1295.0.1295,0.797,0.0528];

»den=[l,-l.8107,2.4947,-1.8801,0.9537,-0.2336];

>>[z,p,K]=tf2zp(num,den);

>>zplane(b,a);

»title('零极点图')；

»

零极点图

绘制幅度频率响应曲线：

»num=[0.0528,0.797,0.1295,0.1295,0.797,0.0528];

»den=[l,-l.8107,2.4947,T.8801,0.9537,-0.2336];

>>[H,W]=freqz(num,den);

»M=abs(H);

»subplot(2,1,1);

»xlabel(*MormalizedFrequency(xnrad/sample),);

»ylabel(,Magritude(dB),);

»title('幅度频率响应曲线')；

»grid;

>>subplot(2,1,2);

»plot(W);

>>xlabel(,MormalizedFrequency(xnrad/sample),);

»ylabel(*Phase(degree)f);

»grid;

»

幅度频率响应曲线

m

p

)

a

p4

n

K

昼2

0

o

100200300400500600

MormalizedFrequency(xrirad/sample)

(

①

出

6

a

p

)

a

s

w

d

100200300400500600

MormalizedFrequency(xrirad/sample)

实验4离散信号的DTFT和DFT

实验目的：加深对离散信号的DTFT和DFT的及其相互关系的理解。

_.OO

实验原理：序列x[n]的DTFT定义：X(*)=^x[n]ejnco

〃=—co

N点序列x[n]的DFT定义：

.祖kNT-kn

x伙]=X(JM

〃=o

N-\

n=0

在MATLAB中，对形式为二P(〃")%+〃，力+・・・+。叱一心的DTDFT可以

D(ejto)4+4"a+...+公<刎

用函数H=Freqz(num,den,w)计算；可以用函数U=fft(u,N)和u=ifft(U,

N)计算N点序列的DFT正、反变换。

....x(n)=cos—«,0<n<\5

实验内容：分别计算16点序列16的16点和32点DFT,

绘出幅度谱图形，并绘出该序列的DTFT图形。

实验要求：讨论DTFT和DFT之间的相互关系。说明实验产生的现象的原因。

实验步骤：

①16点序列X(n)的16点及32点DFT:

»N=16;

»n=l:16;

»x=sin(5*pi*n/16);

»Xl=fft(x,16);

»Xll=abs(XI);

»subplot(2,1,1);

»stem(Xll);

»xlabel('频率');

»ylabclC幅度')；

»title(T6点序列x(n)的16点DFT);

»X2=fft(x,32);

»X22=abs(X2);

»subplot(2,1,2);

»stem(X22);

»xlabelC频率');

»ylabel('幅度');

»title(T6点序列x(n)的32点DFT');

»

②序列x(n)的DTFT:

»x=sin(5*pi*n/16);

>>Xl=fft(x);

»Xll=abs(Xl);

»stem(Xll);

»xlabelC频率')；

»ylabel('幅度');

»title('序列x(n)的DTFT');

»

DTFT和DFT之间的区别和关系：

1、DTFT是离散时间傅里叶变换，DFT是离散傅里叶变换。

2.DTFT变换后的图形中的频率是一般连续的(cos(wn)等这样的特殊函数除外，

其变换后是冲击串)，而DFT是DTFT的等间隔抽样，是离散的点，其函数表示为

X(k),而DTFT的函数表示为(DFT是DTFT的等间隔抽样，DTFT变化后的

频率响应一般是连续的，DFT变换后的频率响应是离散的)。

3、DTFT是以2万为周期的。而DFT的序列X(k)是有限长的,

4、DTFT是以复指数序列｛X(e-M)｝的加权和来表示的，而DFT是等间隔抽样,

抽样间隔为二(N为离散序列的长度)。

5、DTFT和DFT都能表征原序列的信息。由于现在计算主要使用计算机，必需要

是离散的值才能参与运算，因此在工程中DFT应用比较广泛，FFT是DFT的快速

算法。

实验5FFT算法的应用

实验目的：加深对离散信号的DFT的理解及其FFT算法的运用。

实验原理：N点序列的DFT和IDFT变换定义式如下：

N-l1N—1

x[k]=£同川叱＞河川=[7sx伙]叼如

n=0Nk=Q

kn

利用旋转因子W皿，_一e。N具有周期性，可以得到快速算法(FFT)o

在MATLAB中,可以用函数X=fft(x,N)和x=ifft(X,N)计算N点序列的

DFT正、反变换。

实验内容：

2N点实数序列

2〃12万

cos(十7/i)+-cos(十19几)/=0,1,2,…2N—1

x{n}=v

-0,其它〃

N=64。用一个64点的复数FFT程序，一次算出乂㈤二不如⑺以，并绘

#(6|

LU1o

(2)已知某序列M九)在单位圆上的N=64等分样点的Z变换为

X@)=X(Q=匚小^攵=。』2...63。用N点IFFT程序计算

x(n)=IDFT[X(k)]绘出和尤⑺

实验要求：利用MATLAB编程完成计算，绘出相应图形。并与理论计算相比较,

说明实验结果的原因。

实验步骤：

2N点实数序列x(n)的64点的复数FFT：

»k=O:N-l;

»nl=2*n;

»n2=2*n+l;

>>xl=cos(2*pi*7*nl/N)+l/2*cos(2*pi*19*nl/N);

»x2=cos(2*pi*7*n2/N)+l/2*cos(2*pi*19*n2/N);

»XKl=fft(xl,64);

»XK2=fft(x2,64);

»XKll=XKl+(exp(-j*pi*k/N)).*XK2;

»XK22=XKl-(exp(-j*pi*k/N)).*XK2;

»XKll=[zeros(l,N),XKll];

»XK22=[XK22,zeros(1,N)];

»Xk=XKll+XK22;

»k=0:l:2*N;

»mag=abs(Xk);

>>stem(mag);

»xlabel(，kJ);

»ylabelC幅度|X(k)|');

»title,128点序列x(n)的64点的复数FFT');

»

128点序列x(n)的64点的复数FFT

Q

理论计算:

.,,八5I12冗_、12TTIC.

由欧拉公式：品〃]=cos(^x7〃)+5cos^~xl9〃)

17—x7n/—x(A/-7n)17—xl9n1J-x(N-l9n)

N

=-(e+小+_LJN+—eN)

222

j—Kkn

设p[m=eN,则其2N点的DFT变换为:

2N-1・2加2N-1・2m/”、

-Jn^―./(2K-W)

p[m]=^p[n]e2N=2N

w=0n=0

京X2N3-M

当22。机时，p[m]=-------------=0

1-e2N

当2无时,即〃[2月=2N

由此可得M灯当k=14,38,90.114时有值其余为0(04ZW2N—1)

414]=4114]=64.438]=490]=32

与程序计算有相同的结论,可见，理论与实验结果相符合。

(2)X(k)0<jNAIFFT:

»N=64;

»n=0:N-l;

»Xk=l./(l-0.8.*exp(-2*pi*j*n/N));

»xn=ifft(Xk,N);

»stem(n,xn);

»xlabel('时间序列n');

»ylabelC序列x(n)');

»title('X(k)的ifft变换')；

»

X(k)的ifft变换

Q

)

x

宏

微

理论计算：

由M川="勿[川，其Z变换为X(z)二—二可得

i-az

x[n]=0.8"〃[川,其N=64点的DFT为：

63.2成1-0.8641

X[灯=设乂川J'64”下x下

n=0=1-0.8/百1-0.8/百

得

X(zJ=X(k)=工0J/2加N，2=0/2…63

其IDFT为：乂〃]=—LpxO.8"B0.8n(几=0,1,2…63)

1-0.864

则x[O]=l,x[l]=O.8,x[2]=0.64,x[3]=0.512…

与图像所示的x[n]相差很小，可以认为程序计算值与理论计算值相符。

实验6基于MATLAB的数字滤波器设计

实验目的：加深对数字滤波器的常用指标和设计过程的理解。

实验原理：低通滤波器的常用指标：

1-^P<|G(/)归1+%,foi\a\<①p

|G(^)|<8s,forcosw阿W乃

通带边缘频率：@P,阻带边缘频率：-s,

PasshandStnphand

通带起伏：",通带峰，直起伏：Transition

band

Fig7.1Typicalmagnitude

3=-20*(1-“阳，阻带起伏：e,

specificationforadigitalLPF

最小阻带衰减：&s=-201ogo(a)["B]。

数字滤波器有IIR和FIR两种类型，它们的特点和设计方法不同。

在MATLAB中，可以用[b,a]=butter(N,Wn)等函数辅助设计HR数字滤波

器，也可以用b二firl(N,Wn,'ftype')等函数辅助设计FIR数字滤波器。

实验内容：利用MATLAB编程设计一个数字带通滤波器，指标要求如下：

通带边缘频率：%二0・44，%2=0・6万，通带峰值起伏：叫

阻带边缘频率：&用=0,3乃，%2=。・7》,最小阻带衰减,as>40[dB]o

分别用IIR和FIR两种数字滤波器类型进行设计。

实验要求：给出HR数字滤波器参数和FIR数字滤波器的冲激响应，绘出它们

的幅度和相位频响曲线，讨论它们各自的实现形式和特点。

实验步骤：

①Butterworth滤波器的设计（HR）

»wp=[0.4*pi,0.6*pi];

»wr=[0.3*piz0.7*pi];

»Ap=l;

»Ar=40;

»[N,Wn]=buttord(wp/pizwr/pi,Ap,Ar)

N=

7

Wn=

0.38540.6146

,

»[b,a]=butter(N,Wn/bandpass)

b=

Columns1through12

0.00020-0.001400.00420-0.00710

0.00710-0.00420

Columns13through15

0.00140-0.0002

a=

Columns1through12

1.00000.00003.77380.00006.56140.00006.65180.0000

4.20300.00001.64370.0000

Columns13through15

0.36660.00000.0359

»[H,w]=freqz(b,a);

»mag=abs(H);

»plot(w/pi,mag);

>>xlabelC角频率(\Omega)，)；

»ylabel('幅度|Ha(j\Omega)|

»titled数字butterworth带通滤波器幅度响应|Ha(j\Omega)|')；

»phase=angle(H);

»plot(w/pizphase);

>>xlabel('角频率(\0mega)。；

»ylabel(湘位，)；

>>title(，数字butterworth带通滤波器相位响应曲线')；

»

②FIR滤波器的设计：

»wpl=0.4*pi;wp2=0.6*pi;

»wsl=0.3*pi;ws2=0.7*pi;

»tr_width=min((wpl-wsl)/(ws2-wp2))

tr_width=

0.3142

»M=ceil(6.2*pi/tr_width)+1

M=

63

»n=[0:l:M-l];

»wcl=(wsl+wpl)/2;wc2=(wp2+ws2)/2;

»wc=[wcl/pi/wc2/pi];

〉〉wlndow=hannlng(M);

»[hl/w]=freqz(window/l)；

»figure(l);

»subplot(2,l,l)

»stem(window);

»axis([O6001.2]);

»grid;

»xlabel('n');

»title('Hanning窗函数')；

»subplot(2,l,2)

»plot(w/pi,20*log(abs(hl)/abs(hl|l))));

»axis((O1-3500]);

»grid;

»xlabel('wApi')；

»ylabelC幅度(dB))

»title(，Hanning窗函数的频谱,);

»hn=firl(M-lzwc,banning(M));

»[h2,w]=freqz(hn,1/512);

»figure(2);

»subplot(2,l,l)

»stem(n,hn);

»axis([O60-0.250.25));

»grid;

»xlabel('n');

»ylabel('h(n)');

»titlef'Hanning窗函数的单位脉冲响应');

»subplot(2,l,2)

»plot(w/pi,20*log(abs(h2)/abs(h2ll))));

»grid;

»xlabel('wApi')；

»ylabelC幅度(dB)，)；

»figure(3);

»phase=angle(hl);

»plot(phase);

»axis([lpi-10]);

»xlabel('wApi')；

>>ylabe(线性相位,)；

»title('Hanning窗函数相位特性曲线,);

»

实现形式及特点分析：

1.在本例中，相同的技术指标下，HR滤波器实现的阶数为N=7,

而FIR滤波器的阶数N=63。因此，相同的技术指标，用HR滤

波器实现的阶数远远小于用FIR滤波器的阶数。这是由于IIR

滤波器存在着输出对输入的反馈。

2.从实验中绘制的相位特性曲线可以看出，FIR滤波器可以得到

严格的线性相位，而HR滤波器做不到这一点。HR滤波器的

选择性愈好，其相位的非线性愈严重。因而，如果HR滤波器

要得到线性相位，又要满足滤波的技术要求，必须加全通网络

进行相位校正，这同样会大大增加滤波器的阶数。因此，从相

位特性考虑，FIR滤波器又优于IIR滤波器。

附录资料：matlab画二次曲面

一、螺旋线

1.静态螺旋线

a=0:0.1:20*pi;

,

h=plot3(a.*cos(a),a.*sin(a)/2.*a/b'/linewidth',2);

axis([-50,50,-50,50,0,150]);

gridon

set(h/erasemode7none7markersize,z22);

xlabel('x轴');ylabel('yft');zlabel('z轴')；

title('静态螺旋线')；

2.动态螺旋线

t=0:0.1:10*pi;

i=l；

,

h=plot3(sin(t(i))/cos(t(i)),t(i)/*'/'erasemode7none);

gridon

axis([-22-22035])

fori=2:length(t)

setO'xdata'sinmiD/ydatalcos"⑴),'zdata',t(i));

drawnow

pause(O.Ol)

end

title（'动态螺旋线'）;

（图略）

3.圆柱螺旋线

t=0:0.1:10*pi;

x=r.*cos(t);

y=r.*sin(t);

z=t;

plot3(x,y,z,'h','linewidth',2);

gridon

axis('square')

xlabel('x轴');ylabel('y轴');zlabel('z轴');

title('圆柱螺旋线')

二、旋转抛物面

b=0:0.2:2*pi;

[X,Y]=meshgrid(-6:0.1:6);

Z=(X.A2+Y.A2)./4;

meshc(X,XZ);

axis('square')

xlabel('x轴');ylabel('y轴');zlabel('z轴');

shadingflat;

title('旋转抛物面')

或直接用：ezsurfc('(X.A2+Y.A2)./4')

三、椭圆柱面

loadclown

,,

ezsurf((2*cos(u))/4*sin(u)'/'v'/[0/2*piz0/2*pi])

view(-105,40)%视角处理

shadinginterp%灯光处理

colormap(map)%颜色处理

gridon%添加网格线

axisequal%使x,y轴比例一致

xlabel('x轴');ylabel('y轴');zlabel('z轴');

shadingflat;

title（'椭圆柱面'）%添加标题

於回住蜀

四、椭圆抛物面

b=0:0.2:2*pi;

[X,Y]=meshgrid(-6:0.1:6);

Z=X.人2./9+Y.人2./4;

meshc(X,Y,Z);

axis('square')

xlabel('x轴');ylabel('y轴');zlabel('z轴');

shadingflat;

title(椭圆抛物面')

或直接用：ezsurfc('X.A2./9+Y.A2J4')

梢IB食物直

五、'双叶双曲面

,,

ezsurf('8*tan(u)*cos(v)'z8.*tan(u)*sin(v)z'2.*sec(u)'/[-pi./2/3*pi./2J0z2*pi])

axisequal

gridon

axissquare

xlabel('x轴');ylabel('y轴');zlabel('z轴');

shadingflat;

title('双叶双曲面')

六、双曲柱面

loadclown

,,,

ezsurf(2*sec(u)；2*tan(u)；v\[-pi/2,pi/2/-3*pi,3*pi])

holdon%在原来的图上继续作图

,

ezsurf(2*sec(u)72*tan(u)'；v\[pi/2,3*pi/2/-3*piz3*pi])

colormap(map)

shadinginterp

view(-15,30)

axisequal

gridon

axisequal

xlabel('x$ih');ylabel('y轴');zlabel('z轴')；

shadingflat;

title('双曲柱面')