11/16/2006

如何畫圖

圓形是畫圖之重要項目之一,但是在MATLAB並無直接指令可進行畫圓。畫一個圓需要圓心及半徑,亦需要顏色的參數。問題是圓在繪製時,必須先清礎座標的比例,否則會造成楕圓或比例不對稱的圓,此點可以在畫圓前先作等軸的宣告,至於終端機本身之軸向若有不等比例,則必須由終端機之控制鈕調整。即:


>> axis equal;



>>axis image

後者之指令功能與前者大略相同,但它會將所繪的圖自動調整至中央部份,可以一覽全圖。無論如何,經由上述任一個宣告之後,兩軸之比例亦會相等,若不宣告,則必須自行調整視窗之長與寬,使其近似等比的情況。

基本上圓周之構成可用三角函數的關係式計算:
x=rcosθ
y=rsinθ
    1.1
其中角度θ則應自零至360度。而其區間應為θ/ Nb。

繪圓之指令(circle.m)則如下述,其輸入參數包括半徑、圓心位置及繪製之點數。本指令本身並不自行繪圖,主要提供圓周各點之座標,供其他座標操作用途,故必須配合一小程式執行繪圖之功能:


function [coords] = circle(r,x0,y0,nn)
% This function draws a circle in a radius r,
% at the center (x0,y0)
% The inputs:
% r = radius of circle
% x0, y0= coordinates of the circular center
% nn = number of drawing points
% coords(nn,1-2)= vectors to store the coordinates
% Example: circle(10,0,0,10)
jj=0:2*pi/nn:2*pi;
coords=[x0+r*cos(jj);y0+r*sin(jj)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


實例:繪一半徑為5,原點處之圓。

>> coord=circle(5,0,0,20)
coord =
Columns 1 through 8
5.0000 4.7553 4.0451 2.9389 1.5451 0.0000 -1.5451 -2.9389
0 1.5451 2.9389 4.0451 4.7553 5.0000 4.7553 4.0451
Columns 9 through 16
-4.0451 -4.7553 -5.0000 -4.7553 -4.0451-2.9389 -1.5451 -0.0000
2.9389 1.5451 0.0000 -1.5451 -2.9389-4.0451 -4.7553 -5.0000
Columns 17 through 21
1.5451 2.9389 4.0451 4.7553 5.0000
-4.7553 -4.0451 -2.9389 -1.5451 -0.0000
>> plot(coord(1,:),coord(2,:))
>> axis equal
>> grid on




上述指令中,plot()是常用的指令,但也可以用line()這個指令取代(讀者可以自行在MATLAB上測試)。兩個間之不同處是,可以同時繪多圖,亦可在其後面之參數添加顏色值,後則是僅能繪一組線資料,其顏色則必須利用set()之指令進行設定。一般認知上,應為圓滑之曲線,但就circle()指令之操作,實際上圓仍由線段組成。因此段數或點數nn愈多,所繪製之曲線會愈為平滑;而其平滑度仍然會與半徑大小有關。故nn值之設定在本例中也有其重要性。

上述circle()之繪圓指令,只是先獲得圓周點之資料,供其他程式計算之用途。指令若能直接繪圖,則其功能將更為強大。函數circle1()是另一種多功能的繪圓指令,其內容如下:

function h=circle1(r,x0,y0,Nb,C)
% function h=circle1(r,x0,y0,Nb,C)
% Adding circles to the current plot
% Variables:
% r:radi of circle, a scalar or row matrix.
% x0,y0: Centers of circles, a scalar or row matrix
% C:line colors, a string ('r','b'...),
% or RGB values in column
% Nb:No. of drawing points, a scalar or
% row matrix(default=300)
% the size of Nb must be one
% or equal to the size of r.
% h:handles to the circles
% Rules:
% 1. r=a matrix, (x0,y0)=a scalar:Multiple
% co-centered circles
% 2. r=a scalar, (x0,y0)=row matrix: circle with
% r at each center
% 3. r,(x0,y0)=same length row matrix: circle
% with coresponding r at cooresponding center
% 4. r,(x0,y0)=different-length row matrix: Mutiple
% circles with different r at each center.
% Examples: (execute the commands "clf;" first)
%% Example 0:> circle1
%% Example 1
% circle1([1 2 3],[2 3 5],[1 2 1],20);
%
%% Example 2
% the=linspace(0,pi,200);
% r=cos(5*the);
% circle1(0.1,r.*sin(the),r.*cos(the),20,hsv(40));
%
%% Example 3
% [x,y]=meshgrid(1:10,1:10);
% circle1([0.5,0.3,0.1],x,y,20,['r';'y']);
%
%% Example 4
% circle1(1:10,0,0,3:12,[]);
%
%% Example 5
% circle1((1:10),[0,0,20,20],[0,20,20,0]);
%
% rewritten by Din-sue Fon. BIME, NTU. Date:Nov 18,2004.
% Check the number of input arguments
axis equal;
if nargin <5,C=get(gca,'colororder');end
if nargin <4,Nb=300;end
if nargin <3,y0=0;end
if nargin <2,x0=0;end
if nargin <1,r=1;end
% Change matrices into row-wise ones
x0=x0(:);y0=y0(:);r=r(:);Nb=Nb(:);nx=length(x0);
if length(y0)>nx&nx==1,x0=ones(length(y0),1)*x0;end;
if nx>length(y0)&length(y0)==1,y0=ones(nx,1)*y0;end;
nr=length(r);nx=length(x0);nnb=length(Nb);

if length(y0)~=nx,
error('The lengths of x0 and y0 must be identical');
return;
end;
% plotting
for k=1:nx,
if nx==nr,
coords=circle(r(k),x0(k),y0(k),Nb(rem(k-1,nnb)+1)+1);
h(k)=line(coords(1,:),coords(2,:));
set(h(k),'color',C(rem(k-1,size(C,1))+1,:));
else
for m=1:nr
coords=circle(r(m),x0(k),y0(k),...
Nb(rem(m-1,nnb)+1)+1);
h(k,m)=line(coords(1,:),coords(2,:));
set(h(k,m),'color',C(rem(k*m-1,size(C,1))+1,:));
end
end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [coords] = circle(r,x0,y0,nn)
% function [coords] = circle(r,x0,y0,nn)
% This function draws a circle in a radius r,
% at the center (x0,y0)
% The inputs:
% r = radius of circle
% x0, y0= coordinates of the circular center
% nn = number of drawing points
% coords(nn,1-2)= vectors to store the coordinates
% Example: circle(10,0,0,10)
jj=0:2*pi/nn:2*pi;
coords=[x0+r*cos(jj);y0+r*sin(jj)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

1.1.2畫圓程式之應用

在MATLAB中繪線,可用line的功能,這個繪線指令係根據座標點矩陣逐點連線繪製。但是連線時係以直線表示,因此必須在適當的點數下所繪製的圖才能近似圓形,故點數也相當重要。在circle之函數當中,其輸入之變數分別如下:

  •  r::圓之半徑,可為列矩陣,若為列矩陣時,代表可同時繪製許多圓。
  •  x0, y0: 圓心之座標,可為矩陣,若為矩陣時,代表可同時繪製許多不同圓心位置之圓。
  •  C:圓之顏色
  •  Nb: 繪圓時所用之點數,其預設值為300點。

函數circle1()除可繪圖外,可以繪製多種型式的圓,並附上顏色。程式內容看來較為複雜,實際執行僅是其中之下半部,其餘均在核對輸入參數之正確性,使程式更能處理不同的需求。這個程式即使僅打入circle1亦可運作,如:

>> circle1
ans =
152.0021




圖2所示即為其繪圖之結果,由於輸入參數均未設定,故程式僅採用預設值,即半徑為1,圓心為原點。指令窗中所顯示之ans=152.0021為此圖之握把值,該值可以作為程式中指定該圖之參考指標。 

實例:

>> circle1([1:10],0,0)
ans =
Columns 1 through 8
152.0050 156.0038 157.0037 158.0033 159.0033 160.0032 161.0032 162.0032
Columns 9 through 10
163.0032 164.0032



圖3為同心圓,只要半徑為列矩陣,即可繪出不同的同心圓,而其顏色則依預設之順序分配,若欲固定顏色,則可在最後一參數中輸入其RGB值。

同心圓之方式也可稍作改變,這是將半徑10分成1,2,3,…,10之同心圓。若將中心點之位置移動,而半徑值固定,其程式如下:

實例:

>> clf;circle1(20,1:2:10,1:2:10)
ans =
152.0054
156.0042
157.0040
158.0035
159.0035



在此實例中,開始下clf是先將前圖洗淨的意思,若不下此指令,則後圖會與前面的圖重疊,不然就要先以configure(h)設定另一個圖窗。若要不同半徑之圓處於不相干之位置,亦可利用下述之方法求得:

實例:

>clf;circle1([1 2 3],[2 3 5],[1 2 1],20);



在此指令下,可以看到半徑r=[1 2 3]的列向量,而圓心x,y也是同元素之列向量,故繪出其對應的圓。其顏色則依內定的排列方式自動顯出。在此每圖之點數均為20點。這種變化有時可以應用在其他有規則之圓心軌跡上,如下面之實例:
實例:

>> the=linspace(0,pi,200);
>> r=cos(5*the);
>> circle1(0.1,r.*sin(the),r.*cos(the),20,hsv(40));
>> circle1(0.1,r.*sin(the),r.*cos(the),20,hsv(40));




在第一行指令中,有一個新指令linspace(),它是分點的函數,亦即在一數列上自開始點至終點間分成若干線段點,其最後之參數即為段數。請注意它不是linespace()。

在這個例子中,另外加顏色碼,使其變化依繪圓的過程而轉變顏色。hsv(40)是一個色彩度之顏色值,分成40個點,但每點有3項RGB值,故其形成之資料大小應為40x3。顏色亦可用文字代表,即'r'(紅色);'y'(黃色); 'k'(黑色);'b'(藍色)等,但必須以行矩陣表示。

若半徑與圓心座標之元素個數不同,則會在不同的圓心上作同心圓,如下實例:

實例:

%test1
clf;
axis equal;
[x,y]=meshgrid(1:5,1:5);
circle1([0.5,0.3,0.1],x,y,['r';'b']);




網格之製作可利用meshgrid()函數,其參數包含X軸與Y軸部份,以方格進行切割,得到不同之座標位置。利用circle1()函數即可以網格點為圓心,繪製不同之同心圓。

除此以外,亦可利用圓周之點繪製多邊形,其實例如下:

實例:

>> circle1(1:10,0,0,3:12)
ans =
Columns 1 through 8
152.0056 156.0044 157.0043 158.0038 159.0038 160.0034 161.0034 162.0034
Columns 9 through 10
163.0034 164.0034



比較複雜之運用則是除點數外,可以運用參數之可變特性繪製特殊之圓形圖形。

實例:

%test2
%
clf;
x=zeros(1,200);
y=cos(linspace(0,1,200)*4*pi);
rad=linspace(1,0,200);
cmap=hot(50);
circle1(rad,x,y,300,[flipud(cmap);cmap]);




此程式中,所用之函數大部份均已用過,其中zeros()函數較為特殊,通常一個矩陣中各元素要令其為零時,可用此函數,其括號內為矩陣之大小。此函數亦可作為某一矩陣大小之宣告,事先預約空間,使電腦在執行時加快速度。其相對應之函數為ones(),其功能相同,只是此函數將所有元素均設定為一。

在顏色之參數中,有一flipud()函數,這是將矩陣上下鏡射交換,使色調成為對稱狀態。

1.2如何繪製橢圓

某一特定點離兩固定點間之距離和為一定時之軌跡為ellipse 圓。垂直之兩軸是由半長軸Ra與半短軸Rb,構成,其與橢圓周上任意點P之座標關係如下:

x²/Ra²+y²/Rb²=1
x = Racosθ, y=Rb sinθ      1.2



若橢圓之長軸旋轉一個角度φ,則,則其新座標將變為:

x'=x cosφ-y sinφ
y'=x sinφ+y cosφ    1.3

將公式1.2代入1.3則關係式變為:

x'=Racosθ cosφ-Rb sinθ sinφ
y'=Racosθ sinφ+Rb sinθcosφ        1.4

事實上圓的圖形亦可作為ellipse 圓之特殊形,因為只要ellipse 圖之長軸與短軸相等時,即可以作成圓。與前面作圓之情況相同,在執行此指令時,若要避免發生重疊,可以先清除圖窗(即clf;)。

橢圓繪圖程式之應用

依據公式1.4,橢圓之座標點可用ellipsexy()這項函數求得,其輸入參數包括半長短軸Ra、Rb,中心座標(x0,y0),迴轉角度ang及點數 nn等。其座標值則儲存在coords中,可以直接用來繪圖。其程式之內容如下:


function [coords] = ellipsexy(ra,rb,x0,y0,ang,nn)
% function [coords] = ellipsexy(ra,rb,x0,y0,ang,nn)
% This function draws a ellipse with a long radius ra,
% and short radius rb at the center (x0,y0)
% The inputs:
% ra, rb = long & short radii of the ellipse
% x0, y0= coordinates of the ellipse center
% nn = number of drawing points
% ang = angle of the long axis, in radians
% coords(nn,1-2)= vectors to store the coordinates
% Example: ellixy(10,0,0,10)
jj=0:2*pi/nn:2*pi;
coords=[ra*cos(ang)*cos(jj)-rb*sin(ang)*sin(jj)+x0;...
ra*sin(ang)*cos(jj)+rb*cos(ang)*sin(jj)+y0];

本程式ellipsexy()並未有繪線之功能,可以配合plot()或linke()等繪線函數為之。仍採用line 的功能,逐點連線繪製。由於連線係以直線表示,因此必須在適當的點數下才能繪出近似橢圓形,故點數也相當重要。下面之實例中,點數僅10點,故其橢圓未能成形,可以比較其所繪出之圖形。

實例:半長短軸分別為10、5,中心座標為原點,長軸迴轉30度,試繪出其外形。

>> coord=ellipsexy(10,5,0,0,pi/6,10)
coord =
Columns 1 through 8
8.6603 5.5368 0.2985 -5.0538 -8.4758 -8.6603 -5.5368 -0.2985
5.0000 6.5903 5.6633 2.5731 -1.4999 -5.0000 -6.5903 -5.6633
Columns 9 through 11
5.0538 8.4758 8.6603
-2.5731 1.4999 5.0000
>> plot(coord(1,:),coord(2,:))
>> axis equal;grid on



具有繪製橢圓功能之函數為ellipse0()。可以多重繪製各種橢圓。其程式內容如後,主程式並呼叫副程式ellipsexy(),後者可以另行建檔,或附於呼叫程式之後。其相關參數如後:

  • ra, rb:橢圓之半長短軸,可為列矩陣,若為列矩陣時,代表可同時繪製許多橢圓,原則上兩陣列之大小應相同,若其中一個為常數(即為1x1陣列),程式會自動調整與其他一項目同。
  • ang:長軸之傾斜角,以度度表示。可為列矩陣,若與ra或rb之大小相同,則會在各中心點處繪製對應之單一橢圓。
  • x0,y0:橢圓之中心座標,可為矩陣,若為矩陣時,代表可同時繪製許多不同圓心位置之橢圓。
  • C:橢圓線之顏色,可為RGB之三組數值,是為行矩陣。亦可用代碼表示。如['r';'b';'k'...]。注意為使其形成行矩陣,中間必須用分號隔開。
  • Nb: 繪ellipse 圓時所用之點數。預設值為300。

其他參數配合使用原則:

  • 若ra 為單一項,x0,y0 為向量矩陣,則會繪製向量矩陣數之橢 圓。若x0,y0 為單一項,則會繪製同心橢圓。
  • ra 為向量矩陣,則會繪製同一橢圓心之不同半徑向量矩陣數之橢圓。
  • ra 均為同大小之矩陣向量,則會繪製數目相同之橢圓。 若x0,y0 與ra 均為不同大小之矩陣,則會繪製總數為兩矩陣數目之乘積。
  • 實例:繪製一個傾斜某30角度之橢圓。

>> ellipse0(2,1,0,0,30,'r')
ans = 152.0020


實例:在相同半長短軸及中心點下,可以變化角度,產生旋轉之橢圓群:

>> ellipse0(2,1,0,0,0:20:360)
ans =
Columns 1 through 8
160.0020 161.0020 162.0020 163.0020 164.0020 165.0020 166.0020 167.0020
Columns 9 through 16
168.0020 169.0020 170.0020 171.0020 172.0020 173.0020 174.0020 175.0020
Columns 17 through 19
176.0020 177.0020 178.0020



實例:若半短軸維持不變(可以輸入單一值),配合半長軸進行變化。

>> ellipse0(1:10,4,0,0)
ans =
Columns 1 through 8
152.0024 156.0024 157.0024 158.0024 159.0024 160.0024 161.0024 162.0024
Columns 9 through 10
163.0024 164.0024



實例:同時變化半長短軸,其他維持不變。

>> ellipse0(1:6,0.5:0.1:1.0)
ans = 152.0026 156.0026 157.0026 158.0026 159.0026 160.0026




實例:變化半長短軸,也變化Y軸之高度。

>> ellipse0([1:10]*2,1:10,0,1:10)
ans =
152.0037 156.0037 157.0037 158.0037 159.0037 160.0037
161.0037 162.0037 163.0037 164.0037



實例:曲線變化半長短軸,也變化Y軸之高度。

>>ellipse0([1:10].^2,1:10,0,[1:10]*10);




程式內容:

function h=ellipse0(ra,rb,x0,y0,ang,C,Nb)
% h=ellipse0(ra,rb,ang,x0,y0,C,Nb)
% Adding ellipse to the current plot
% Variables:
% ra,rb:longitudinal & horizontal axes of a ellipse, scalar or matrix.
% ang: angle the ellipse inclines, in deg.
% x0,y0: Centers of the ellipses, a scalar or row matrix
% C:line colors, a string ('r','b','k'...),or RGB values in column
% Nb:No. of drawing points, a scalar or row matrix(default=300)
% h:handles to the ellipse
% Rules:
% 1. ra,rb= matrix, (x0,y0)=a scalar:Multiple co-centered ellipes
% 2. ra,rb= scalar, (x0,y0)=row matrix: ellipse with ra & rb at each center
% 3. ra,rb,(x0,y0)=same length row matrix: ellipse with coresponding r at
% cooresponding center
% 4. ra,rb & (x0,y0)=different-length row matrix: multiple ellipses with
% different pairs of ra & rb at each center
% 5. ra, rb should have same size
% Examples: (execute the commands "clf;" first)
%% Example 0: ellipse0
% ellipse0(2,1,0,0,0:10:360);
% clf;ellipse0(3,1,[0 2 5],[0 1 4],0:10:360)
% ellipse0(ra,rb,ang,x0,y0,C,Nb), Nb specifies the number of points
% ellipse(1,2,pi/8,1,1,'r')
% author: Din-sue Fon. Bime, NTU, Date:November 18, 2004
axis equal;
if nargin <7,Nb=300;end
if nargin <6,C=get(gca,'colororder');end
if nargin <5,ang=0;end
if nargin <4,y0=0;end
if nargin <3,x0=0;end
if nargin <2,rb=1;end
if nargin <1,ra=2;end
% change the matrices into one row
x0=x0(:);y0=y0(:);ra=ra(:);rb=rb(:);ang=ang(:)*pi/180;
nr=length(ra);nx=length(x0);
if length(rb)>nr&nr==1,ra=ones(length(rb),1)*ra;end;
if nr>length(rb)&length(rb)==1,rb=ones(nr,1)*rb;end;
if length(y0)>nx&nx==1,x0=ones(length(y0),1)*x0;end;
if nx>length(y0)&length(y0)==1,y0=ones(nx,1)*y0;end;
nr=length(ra);nx=length(x0);nnb=length(Nb);nang=length(ang);
if length(y0)~=nx|length(rb)~=nr,
error('The lengths of ra & rb should have same size!');
return;
end
if isstr(C),C=C(:);end;
for k=1:nx,
if nx==nr,
for n=1:nang
coords=ellipsexy(ra(k),rb(k),x0(k),y0(k),ang(n),...
Nb(rem(k-1,nnb)+1)+1);
h(k,n)=line(coords(1,:),coords(2,:));
set(h(k,n),'color',C(rem(k*n-1,size(C,1))+1,:));
end
else
for m=1:nr
for n=1:nang
coords=ellipsexy(ra(m),rb(m),x0(k),y0(k),ang(n),...
Nb(rem(m-1,nnb)+1)+1);
h(k,m)=line(coords(1,:),coords(2,:));
set(h(k,m),'color',C(rem(k*m*n-1,size(C,1))+1,:));
end
end
end
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [coords] = ellipsexy(ra,rb,x0,y0,ang,nn)
% function [coords] = ellipsexy(ra,rb,x0,y0,ang,nn)
% This function draws a ellipse with a long radius ra,
% and short radius rb at the center (x0,y0)
% The inputs:
% ra, rb = long & short radii of the ellipse
% x0, y0= coordinates of the ellipse center
% nn = number of drawing points
% ang = angle of the long axis, in radians
% coords(nn,1-2)= vectors to store the coordinates
% Example: ellixy(10,0,0,10)
jj=0:2*pi/nn:2*pi;
coords=[ra*cos(ang)*cos(jj)-rb*sin(ang)*sin(jj)+x0;...
ra*sin(ang)*cos(jj)+rb*cos(ang)*sin(jj)+y0];