Chap 5. 特殊圖形 方煒 台大生機系
二維 特殊圖形
特殊繪圖函數 指令 說明 bar, barh, bar3, bar3h 長條圖 Area 面積圖 pie, pie3 扇形圖 stem, stem3 針頭圖 stairs 階梯圖 fill, fill3 實心圖 quiver, quiver3 向量場圖 contour, contourf, contour3 等高線圖
長條圖之繪製 長條圖 (Bar Graphs)特別適用於少量且離散的資料。欲畫出垂直長條圖,可用 bar 指令。 範例:bar01.m x = [1 3 4 5 2]; bar(x);
長條圖之繪製(cont.) bar 指令也可接受矩陣輸入,它會將同一橫列的資料聚集在一起。 範例:bar02.m x = [2 3 4 5 7; 1 2 3 2 1]; bar(x);
長條圖之繪製(cont.) bar 及 barh 指令還有一項特異功能,就是可以將同一橫列的資料以堆疊(Stack)方式來顯示。 範例:bar03.m x = [2 3 4 5 7; 1 2 3 2 1]; bar(x,'stack')
長條圖之繪製(cont.) 除了平面長條圖之外,MATLAB 亦可使用 bar3 指令來畫出立體長條圖。 範例:bar04.m x = [2 3 4 5 7; 1 2 3 2 1]; bar3(x)
長條圖之繪製(cont.) bar3 指令還可以使用群組(Group)方式來呈現長條圖 範例:bar05.m x = [2 3 4 5 7; 1 2 3 2 1]; bar3(x, 'group')
長條圖之繪製(cont.) 長條圖的指令和類別 : 垂直長條圖 水平長條圖 平面 bar barh 立體 bar3 bar3h
長條圖之繪製(cont.) 若要指定長條圖的 x 座標,可使用兩個輸入向量給 bar 指令。假設新竹的月平均溫度如下: 範例:bar06.m y = 35*rand(1, 6); % 溫度值(假設是介於 0~35 的亂數) bar(x, y); xlabel('月份'); % x 軸的說明文字 ylabel('平均溫度 (^{o}c)'); % y 軸的說明文字 % 下列指令將 x 軸的數字改成月數 set(gca, 'xticklabel', {'一月','二月','三月', '四月', '五月', '六月'});
長條圖之繪製(cont.)
面積圖之繪製 面積圖(Area Graphs)和以堆疊方式呈現的長條圖很類似,特別適用於具有疊加關係的資料。舉例來說,若要顯示清華大學在過去 10 年來的人數(含大學部,研究生,及教職員)變化情況,可用面積圖顯示。 範例:area01.m y = rand(10,3)*100; x = 1:10; area(x, y); xlabel('Year'); ylabel('Count')
扇形圖之繪製 使用 pie 指令,可畫出平面扇形圖(Pie Charts),並可加上說明。 範例:pie01.m x = [2 3 5 4]; label={'東','南','西','北'}; pie(x, label);
扇形圖之繪製(cont.) pie 指令直接將 x 元素視為面積百分比,因此可畫出不完全的扇形圖。 範例:pie02.m
扇形圖之繪製(cont.) pie 指令還有一特異功能,可將某個或數個扇形圖向外拖出,以強調部份資料。 範例:pie03.m x = [2 3 5 4]; explode = [1 1 0 0]; pie(x, explode); 其中指令 explode 中非零的 元素即代表要向外拖出的扇形。
扇形圖之繪製(cont.) 欲畫出立體扇形圖,可用 pie3 指令。 範例 : pie301.m x = [2 3 5 4]; explode = [1 1 0 0]; label = {'春','夏','秋','冬'}; pie3(x, explode, label);
針頭圖之繪製 顧名思義,針頭圖(Stem Plots)就是以一個大頭針來表示某一點資料,其指令為 stem。 範例:stem01.m t = 0:0.2:4*pi; y = cos(t).*exp(-t/5); stem(t, y)
針頭圖之繪製(cont.) 針頭圖特別適用於表示「數位訊號處理」(DSP,Digital Signal Processing)中的數位訊號。若要畫出實心的針頭圖,可加“fill”選項。 範例:stem02.m t = 0:0.2:4*pi; y = cos(t).*exp(-t/5); stem(t, y, 'fill');
5-4 針頭圖之繪製(cont.) 欲畫出立體的針頭圖, 可用 stem3 指令。 範例5-14:stem301.m theta = -pi:0.05:pi; x = cos(theta); y = sin(theta); z = abs(cos(3*theta)).*exp(-abs(theta/3)); stem3(x, y, z); Fig. 5-14
階梯圖之繪製 使用 stairs 指令,可畫出階梯圖(Stairstep Plots),其精神和針頭圖很相近,只是將目前資料點的高度向右水平畫至下一點為止。(在數位訊號處理,此種作法稱為 Zero-order Hold。) 範例:stairs01.m t = 0:0.4:4*pi; y = cos(t).*exp(-t/5); stairs(t, y);
階梯圖之繪製(cont.) 若再加上針頭圖,則可見兩 者相似之處。 範例:stairs02.m t = 0:0.4:4*pi; y = cos(t).*exp(-t/5); stairs(t, y); hold on % 保留舊圖形 stem(t, y); % 疊上針頭圖 hold off
實心圖之繪製 MATLAB 指令 fill 將資料點視為多邊形頂點,並將此多邊形塗上顏色,呈現出實心圖(Filled Plots)的結果。 範例 : fill01.m t = 0:0.4:4*pi; y = sin(t).*exp(-t/5); fill(t, y, 'b'); % 'b'為藍色
實心圖之繪製(cont.) 若與 stem 合用,則可創造出 一些不同的視覺效果。 範例 : fill02.m t = 0:0.4:4*pi; y = sin(t).*exp(-t/5); fill(t, y, 'y'); % 'y' 為黃色 hold on % 保留舊圖形 stem(t, y, 'b'); % 疊上藍色針頭圖 hold off
實心圖之繪製(cont.) fill3 可用於三維的實心圖。 範例: fill301.m X = [0 0 1 1]; Y = [0 1 1 0]; Z = [0 1 1 0]; C = [0 0.3 0.6 0.9]'; fill3(X, Y, Z, C);
實心圖之繪製(cont.) 使用 fill3 指令,我們亦可以 畫出各種酷酷的圖形。 範例 : fill302.m t = (1/16:1/8:1)'*2*pi; x = sin(t); y = cos(t); c = linspace(0, 1, length(t)); fill3(x, y/sqrt(2), y/sqrt(2), c, x/sqrt(2), y, x/sqrt(2), c); axis tight
常用三維圖形 quiver(x, y, z) meshc contour3 contourf quiver3
向量場圖之繪製 使用 quiver 指令可畫出平面 上的向量場圖(Quiver Plots) ,特別適用於表示分布於平面 的向量場 (Vector Fields), 例如平面上的電場分布,或是流 速分布 。 範例:quiver01.m [x, y, z] = peaks(20); [u, v] = gradient(z); contour(x, y, z, 10); hold on, quiver(x,y,u,v); hold off; axis image
向量場圖之繪製(cont.) 欲畫出空間中的向量場圖, 可用 quiver3 指令。 範例:quiver301.m [x, y] = meshgrid(-2:0.2:2, -1:0.1:1); z = x.*exp(-x.^2-y.^2); [u, v, w] = surfnorm(x, y, z); quiver3(x, y, z, u, v, w); hold on, surf(x, y, z); hold off axis equal
等高線圖之繪製 我們可用 contour 指令來畫出「等高線圖」(Contour Plots)。 範例:contour01.m z = peaks; contour(z, 30); % 畫出 30 條等高線
等高線圖之繪製(cont.) 若要畫出特定高度的等高線,可執行如下: 範例:contour02.m z = peaks; contour(z, [0 2 5]);
等高線圖之繪製(cont.) 欲標明等高線的高度, 可用 clabel 指令。 範例:contour03.m z = peaks; [c,handle] = contour(z, 10); clabel(c, handle);
等高線圖之繪製(cont.) 若欲在等高線之間填入顏色,可用 contourf 指令。 範例:contour04.m z = peaks; contourf(z);
等高線圖之繪製(cont.) 若要使畫出的等高線對 應至正確的 x 及 y 座標,則可執行如下: 範例:contour05.m [x,y,z] = peaks; contour(x, y, z); % 使用三個輸入
等高線圖之繪製(cont.) contourf 亦可接受 x、y、z 輸入引數。若要將等高線畫在曲面的正下方,可用 surfc 或 meshc 指令。 範例:contour06.m [x, y, z] = peaks; meshc(x, y, z); axis tight
等高線圖之繪製(cont.) 若要畫出三度空間中的等高線,可用 contour3 指令。 範例:contour301.m [x, y, z] = peaks; contour3(x, y, z, 30); axis tigh
等高線圖之繪製(cont.) 使用 contour 指令亦可畫出極座標中的等高線,但過程較為複雜,以下列複數函數為例: 其中 z 代表複數平面中的任一點複數,如果我們要畫出此函數的等高線,可見下列範例: 範例:contour07.m t = linspace(0, 2*pi, 61); % 角度的格子點 r = 0:0.05:1; % 長度的格子點 [tt, rr] = meshgrid(t, r); % 產生二維的格子點 [xx, yy] = pol2cart(tt, rr); % 將極座標轉換至直角座標 zz = xx + sqrt(-1)*yy; % 複數表示 ff = abs(zz.^3-1); % 曲面的函數 contour(xx, yy, ff, 50); % 畫出等高線 axis image
等高線圖之繪製(cont.)
等高線圖之繪製(cont.) 上例中,座標的標示仍為直 角座標。 欲將等高線顯示於極座標上,需先用 polar 指令產生一個極座標圖,再移除圖形,留下圖軸,然後再進行作圖。 範例:contour08.m h = polar([0 2*pi], [0 1]); delete(h); hold on contour(x, y, abs(f), 30); hold off % 產生在極座標上的一條直線 % 移除上述圖形,但留下極座標圖軸 *要注意的是,你必須先執行 contour07.m, 然後再執行 contour08.m,才能得到上述 的極座標等高線的效果。
等高線圖之繪製(cont.) 我們也可以同時畫出複數函數的曲面和等高線圖,例如,下列範例可以畫出複數函數: 範例:contour09.m t = linspace(0, 2*pi, 61); % 角度的格子點 r = 0:0.05:1; % 長度的格子點 [tt, rr] = meshgrid(t, r); % 產生二維的格子點 [xx, yy] = pol2cart(tt, rr); % 將極座標轉換至直角座標 zz = xx + sqrt(-1)*yy; % 複數表示 ff = abs(zz.^3-1); % 曲面的函數值 h = polar([0 2*pi], [0 1]); % 產生在極座標上的一條直線 delete(h); % 移除上述圖形,但留下極座標圖軸 hold on contour(xx, yy, ff, 20); % 等高線 surf(xx, yy, ff); % 曲面圖 hold off view(-19, 22); % 設定觀測角度
等高線圖之繪製(cont.)
Ex5_1a meshgrid %********************************************************* % Define the arrays x and y. Warning: don’t make the step size too small clear; close all; x=-1:.1:1; y=0:.1:1.5; % Use meshgrid to convert these 1-d arrays into 2-d matrices of % x and y values over the plane [X,Y]=meshgrid(x,y); % Get f(x,y) by using f(X,Y). Note the use of .* with X and Y rather than with x and y f=(2-cos(pi*X)).*exp(Y); % Note that this function is uphill in y between x=-1 and x=1 and has a valley at x=0 surf(X,Y,f); length(x), length(y), size(f)
Ex5_1b ndgrid [X,Y]=ndgrid(x,y); f=(2-cos(pi*X)).*exp(Y); surf(X,Y,f); length(x) length(y) size(f)
title(’Contour Plot’); xlabel(’x’); ylabel(’y’); pause Ex5_2a Contour Plots %***************************************************** % make a contour plot by asking Matlab to evenly space N contours % between the minimum of f(x,y) and the maximum f(x,y) (the default) N=40; contour(X, Y, f, N); title(’Contour Plot’); xlabel(’x’); ylabel(’y’); pause
Ex5_2b Contour Plots %**************************************************** % You can also tell Matlab which contour levels you want to plot. close all; top=max(max(f)); % find the max and min of f bottom=min(min(f)); dv=(top-bottom)/20; % interval for 21 equally spaced contours V=bottom:dv:top; cs=contour(X,Y,f,V); clabel(cs,V(1:2:21)) % give clabel the name of the plot and % an array of the contours to label title(’Contour Plot’) xlabel(’x’); ylabel(’y’) pause
Ex5_2c Surface Plots close all; %********************************************************* % Now make a surface plot of the function with the viewing % point rotated by AZ degrees from the x-axis and % elevated by EL degrees above the xy plane close all; surf(X,Y,f); % or you can use mesh(X,Y,f) to make a wire plot AZ=30;EL=45; view(AZ,EL); title(’Surface Plot’); xlabel(’x’); ylabel(’y’); pause
Ex5_2d Surface Plots close all; surf(X,Y,f); title(’Surface Plot’) %********************************************************************* % Here’s a piece of code that lets you fly around the % surface plot by continually changing the viewing angles % and using the pause command. %******************************************************************** close all; surf(X,Y,f); title(’Surface Plot’) xlabel(’x’); ylabel(’y’); EL=45; for m=1:100 AZ=30+m/100*360; view(AZ,EL); pause(.1); % pause units are in seconds end pause
Ex5_2e Surface Plots dt=.1; for m=1:100 %******************************************************************************* % This same trick will let you make animations of both xy and surface plots. To make this surface % oscillate up and down like a manta ray you could do this. dt=.1; for m=1:100 t=m*dt; g=f*cos(t); surf(X,Y,g); AZ=30;EL=45; view(AZ,EL); title(’Surface Plot’); xlabel(’x’); ylabel(’y’) axis([-1 1 -1 1 min(min(f)) max(max(f))]); pause(.1) end help graph3d
Ex5_3 Evaluating Fourier Series clear; close all; % set some constants a=2;b=1;Vo=1; Nx=80;Ny=40; % build the x and y grids dx=2*b/Nx;dy=a/Ny; x=-b:dx:b; y=0:dy:a; [X,Y]=meshgrid(x,y); % build the 2-d grids for plotting % set the number of terms to keep and do the sum Nterms=20; % zero V out so we can add into it V=zeros(Ny+1,Nx+1); % add the terms of the sum into V for m=0:Nterms V=V+cosh((2*m+1)*pi*X/a)/cosh((2*m+1)*pi*b/a).*sin((2*m+1)*pi*Y/a)/(2*m+1); end % put on the multiplicative constant V=4*Vo/pi*V; % surface plot the result surf(X,Y,V); xlabel(’x’); ylabel(’y’); zlabel(’V(x,y)’);
Ex5_4a Vector Field Plots clear;close x=-5.25:.5:5.25;y=x; % define the x and y grids (avoid (0,0)) [X,Y]=meshgrid(x,y); % Electric field of a long charged wire Ex=X./(X.^2+Y.^2); Ey=Y./(X.^2+Y.^2); quiver(X,Y,Ex,Ey); % make the field arrow plot title(’E of a long charged wire’) axis equal % make the x and y axes be equally scaled pause % Magnetic field of a long current-carrying wire Bx=-Y./( X.^2+Y.^2); By=X./(X.^2+Y.^2); % make the field arrow plot quiver(X,Y,Bx,By) axis equal title(’B of a long current-carrying wire’)
Ex5_4b Vector Field Plots %******************************************************** % The big magnitude difference across the region makes most arrows too small % to see. This can be fixed by plotting unit vectors instead (losing all % magnitude information B=sqrt(Bx.^2+By.^2); Ux=Bx./B; Uy=By./B; quiver(X,Y,Ux,Uy); axis equal title(’B(wire): unit vectors’) pause
Ex5_4c Vector Field Plots %********************************************************* % You can still see qualitative size information % without such a big variation in arrow size by % having the arrow length be logarithmic. If s is % the desired ratio between the longest arrow and % the shortest one, this code will make the appropriate % field plot. Bmin=min(min(B)); Bmax=max(max(B)); s=2; % choose an arrow length ratio k=(Bmax/Bmin)^(1/(s-1)); logsize=log(k*B/Bmin); Lx=Ux.*logsize; Ly=Uy.*logsize; quiver(X,Y,Lx,Ly); axis equal title(’B(wire): logarithmic arrows’)
其他進階繪圖功能 MATLAB 在 5.3 版後,開始支援「容積目視法」(Volume Visualization) 能夠畫出在三度空間中的 流線圖 向量場圖 等高面圖(Isosurfaces) 切面圖(Slices)等
容積目視法 相關指令 v5.3以上版本 說明 coneplot 以圓錐瓶畫出三度空間的向量場圖 contourslice 在三度空間的切面上畫出等高線 isosurface 從容積資料中算出等高面資料 isocaps 計算等高面在端點切片的等高資訊 isonormals 計算等高面的法向量 slice 在三度空間的切片 streamline 從 2D 或 3D 的流線資料來畫流線圖
容積目視法 相關指令 v6.0以上版本 streamtube streamribbon streamslice isocolors 計算等高區面頂點的顏色 divergence 計算3-D向量場的亂度(Divergence) curl 計算3-D向量場的curl 及垂直方向的角速度 streamtube 由向量資料畫出流線管(Stream Tubes) streamribbon 由向量資料畫出流線緞帶(Stream Ribbons) streamslice 由向量資料畫出間隔分明的流線 streamparticles 由向量資料畫出流線粒子(Stream Particles) interpstreamspeed 由速度對流線頂點做內差(Interpolation) volumebounds 傳回容積資料的座標及顏色極限值
容積目視法 相關示範 在 MATLAB 指令視窗下輸入volvec,以開啟展示視窗 共有六個示範程式 容積目視法 相關示範 在 MATLAB 指令視窗下輸入volvec,以開啟展示視窗 共有六個示範程式 在點選「Multiple」之後,產生如右圖形
植物工廠內栽培床架上 人工光源配置與 給光均勻度分析 Contour map 應用例 植物工廠內栽培床架上 人工光源配置與 給光均勻度分析
植物工廠示意圖 保溫層讓內外部熱傳導降至最小 人工光源 Artificial light
植物工廠
植物工廠
光分佈量測與改善 燈具高度 10 公分(量測距離),栽培床為120 x 40 公分保麗龍板,於長邊取13點,短邊取6點,共78個測量點 循序改善:
光分佈改善 (1=>2)
光分佈改善 (2=>3)
光分布測量結果 條件 光量 ( μmol/m2/s ) 消耗功率 ( W ) 1.等距離排列 246 85 348 168 W 六支長燈管 平均值 標準差 最大值 1.等距離排列 246 85 348 168 W 六支長燈管 2.等距排列,中間間隔增加 215 63 284 3.中間間隔增加且兩側補光 236 45 286 184 W 六支長燈管+兩支短燈管
加反光紙 (6支燈管減少為4支燈管) 減少燈管至4支 雖然降低光量,光分佈能夠較均勻 較節省能源 側視圖 6 支燈管 +兩側補光+兩側反光 4 支燈管+ 兩側補光 +兩側反光 10 cm 10 cm 5cm 5cm 反光紙
五種給光條件之彙整 條件 光量 ( μmol/m2/s ) 消耗功率 ( W ) 1.等距離排列 246 85 348 168 六支長燈管 平均值 標準差 最大值 1.等距離排列 246 85 348 168 六支長燈管 2.等距排列但中間間隔增加 215 63 284 3.中間間隔增加且兩側補光 236 45 286 184 六支長燈管+兩支短燈管 4.於條件3加反光紙 337 40 405 5.於條件4減兩支長燈管 230 30 282 128 四支長燈管+兩支短燈管
改善前後光分佈立體圖 改善前 (條件1) 改善後(條件5)
改善前後光分佈等高線圖
結論 透過燈管排列的調整與增加反光紙 可達到 目的 減少燈管數量 (6長2短 4 長2短) 減少燈管功率 (184 W降至128 W) 固定成本降低接近 1/3 減少燈管功率 (184 W降至128 W) 操作成本減少30 % 提高均勻度 (標準差從85降至30 μmol/m2/s) 目的