Tạp chí Khoa học - Công nghệ Thủy sản<br />
<br />
Số 4/2014<br />
<br />
THOÂNG BAÙO KHOA HOÏC<br />
<br />
TÍNH ĐỘ CONG BỀ MẶT CHO PHÂN VÙNG BỀ MẶT TỰ DO<br />
DỰA TRÊN PHẦN MỀM MATLAB<br />
SURFACE CURVATURE COMPUTATION FOR FREE-FORM SURFACE<br />
PARTITIONING BASED ON MATLAB PROGRAM<br />
Nguyễn Văn Tường1<br />
Ngày nhận bài: 05/8/2014; Ngày phản biện thông qua: 11/8/2014; Ngày duyệt đăng: 01/12/2014<br />
<br />
TÓM TẮT<br />
Có thể chia bề mặt tự do thành các vùng lồi, lõm, yên ngựa khác nhau nhờ vào đặc điểm độ cong bề mặt tại các điểm<br />
trên bề mặt. Bài báo này trình bày việc tính độ cong của bề mặt tự do cho mục đích phân vùng. Các thông số độ cong bề<br />
mặt được sử dụng làm dữ liệu cho quá trình phân vùng và xác định biên các vùng khi sử dụng các đặc tính hình học của<br />
bề mặt và kỹ thuật mã xích trong lĩnh vực xử lý ảnh. Quá trình tính toán được thực hiện nhờ chương trình được viết bằng<br />
phần mềm Matlab. Dữ liệu đầu vào cho chương trình là phương trình toán học của bề mặt tự do ở dạng tường minh hoặc<br />
bề mặt Bspline. Tọa độ các điểm trên bề mặt cũng như trên biên của các vùng trong dữ liệu đầu ra được sử dụng cho việc<br />
mô hình hóa bề mặt với các vùng riêng biệt trong môi trường CAD (Computer Aided Design).<br />
Từ khóa: Bề mặt tự do, độ cong bề mặt, phân vùng bề mặt<br />
<br />
ABSTRACT<br />
A free-form surface can be partitioned into different convex, concave and saddle regions thanks to the characteristics<br />
of surface curvatures at points on the surface. This paper presents the work of surface curvature computation for free-form<br />
surface partitioning. The surface curvatures are used as data for partitioning and defining the boundaries of regions on<br />
the surface when using the characteristics of surface geometry and the chain code technique in image processing field.<br />
The computation process is performed by a Matlab program. The input data of the program are mathematical equations of<br />
free-form surfaces in explicit form or Bspline surfaces. The coordinates of points on the surface and on the region<br />
boundaries in the output data are used for modelling the surface with separate regions in CAD environment.<br />
Keywords: Free-form surface, surface curvatures, surface partitioning<br />
I. ĐẶT VẤN ĐỀ<br />
Bề mặt tự do là mặt đều, trơn, thường được<br />
sử dụng trong các ngành thiết kế và chế tạo khuôn<br />
mẫu, thiết kế thân ô tô, tàu thủy, máy bay và trong<br />
các tác phẩm nghệ thuật. Quá trình gia công bề mặt<br />
tự do trên máy CNC (Computer Numerical Control)<br />
thường tốn nhiều thời gian do đường kính dao bị<br />
hạn chế bởi bán kính cong nhỏ nhất của bề mặt<br />
cần gia công. Một trong những phương pháp nâng<br />
cao năng suất gia công bề mặt tự do là chia bề mặt<br />
thành các vùng khác nhau và mỗi vùng có thể được<br />
gia công bằng các dao có đường kính khác nhau.<br />
Chen và cs [2] đã tính các tính chất hình<br />
học của bề mặt tự do như độ cong Gauss, độ<br />
cong trung bình, độ cong cực đại và cực tiểu và<br />
1<br />
<br />
pháp véc tơ bề mặt để thành lập một véc tơ đa chiều<br />
cho quá trình chia vùng bề mặt tự do. Phương pháp<br />
nhóm cụm mờ (fuzzy clustering) và phương pháp<br />
C-means mờ (fuzzy C-means) đã được các tác giả<br />
sử dụng để chia bề mặt tự do thành các vùng lồi,<br />
lõm và yên ngựa. Các tính chất hình học nói trên<br />
cũng được Roman và cs [7, 8] để xác định biên và<br />
phân vùng bề mặt tự do.<br />
Bey và cộng sự [1] đã xấp xỉ bề mặt tự do thành<br />
các tam giác. Các thông số pháp véc tơ và độ cong<br />
bề mặt được tính để xác định hình dạng cục bộ vùng<br />
bề mặt tại các đỉnh tam giác. Từ đó, các đỉnh này<br />
được nhóm thành các vùng có đồ hình khác nhau.<br />
Để nâng cao hiệu suất gia công bằng cách sử<br />
dụng dao lớn nhất có thể, Li và Zhang [4] cũng chia<br />
<br />
TS. Nguyễn Văn Tường: Khoa Cơ khí – Trường Đại học Nha Trang<br />
<br />
TRƯỜNG ĐẠI HỌC NHA TRANG • 65<br />
<br />
Tạp chí Khoa học - Công nghệ Thủy sản<br />
<br />
Số 4/2014<br />
<br />
bề mặt tự do thành các vùng lồi, lõm và yên ngựa<br />
nhờ độ cong bề mặt. Họ gọi các vùng lõm và yên<br />
ngựa là những vùng tới hạn mà ở đó có thể xảy ra<br />
việc cắt lẹm. Các tác giả đã xây dựng chương trình<br />
tính toán chia vùng và kiểm tra cắt lẹm bằng ngôn<br />
ngữ C++.<br />
Elber và Cohen [3] đã tiến hành phân tích độ<br />
cong của bề mặt tự do để nghiên cứu đặc tính đồ<br />
hình bề mặt. Họ đã phát triển một phương pháp lai<br />
sử dụng các toán tử ký hiệu và toán tử số để tính<br />
các độ cong bề mặt, từ đó xác định biên các đường<br />
biên của các vùng riêng biệt trên bề mặt. Tuy nhiên,<br />
các phương trình đường biên là những đa thức bậc<br />
cao, rất khó giải.<br />
Nói tóm lại, cho đến nay, có nhiều công trình<br />
nghiên cứu phân vùng bề mặt tự do dựa trên cách<br />
tiếp cận dùng độ cong bề mặt. Tuy nhiên đa số các<br />
phương pháp đã đưa ra là khá phức tạp, khó áp<br />
dụng. Tường và Pokorny [9] đã đưa ra một phương<br />
pháp đơn giản nhưng hiệu quả để phân vùng bề mặt<br />
tự do. Ở đây, bề mặt tự do được chia thành các vùng<br />
lồi, lõm và yên ngựa dựa trên độ cong bề mặt. Đường<br />
biên của các vùng được xác định nhờ áp dụng kỹ<br />
thuật mã xích dùng trong xử lý ảnh. Quá trình tính<br />
toán phân vùng và xác định biên bề mặt được thực<br />
hiện nhờ một chương trình Matlab được viết cho bề<br />
mặt tự do ở dạng tường minh hoặc bề mặt Bspline.<br />
Bài báo này tập trung giới thiệu phương pháp tính<br />
toán độ cong bề mặt tự do cho mục đích phân vùng<br />
với sự hỗ trợ của phần mềm Matlab.<br />
II. ĐỐI TƯỢNG NGHIÊN CỨU VÀ PHƯƠNG PHÁP<br />
NGHIÊN CỨU<br />
1. Thông số hình học của bề mặt tự do<br />
Bề mặt tự do có thể được biểu diễn theo:<br />
- Dạng ẩn: f(x, y, z) = 0<br />
(1)<br />
- Dạng tường minh: z = f(x, y)<br />
(2)<br />
- Dạng tham số:<br />
S(u,v) = {Sx(u,v), Sy(u,v), Sz(u,v)}<br />
(3)<br />
Một số thông số hình học chủ yếu của bề mặt<br />
tự do là:<br />
a. Pháp véc tơ tại một điểm<br />
Cho một bề mặt tự do S(u, v) và một điểm bất<br />
kỳ trên bề mặt. Tại điểm này, Su và Sv là hai véc tơ<br />
tiếp tuyến theo hai phương tham số u và. Su và Sv<br />
không song song nhau, và một véc tơ vuông góc với<br />
cả hai véc tơ này là véc tơ đơn vị, được xác định bởi<br />
công thức (4) và được biểu diễn như trên hình 1 [6].<br />
(4)<br />
<br />
Véc tơ t bất kỳ vuông góc với nS được gọi là<br />
véc tơ tiếp tuyến với S(u,v) tại điểm p. Mặt phẳng<br />
chứa tất cả các véc tơ tiếp tuyến với mặt S tại điểm<br />
p được gọi là mặt phẳng tiếp tuyến tại điểm p, và<br />
được ký hiệu là Tp(S) (hình 1).<br />
<br />
Hình 1. Pháp véc tơ đơn vị và mặt phẳng tiếp tuyến<br />
tại một điểm<br />
<br />
b. Dạng toàn phương thứ nhất, F1<br />
Dạng toàn phương thứ nhất của một bề mặt tự<br />
do S biểu diễn các tính chất khối của bề mặt, được<br />
xác định bởi [6]:<br />
F1= dS.dS = Edu2 + 2Fdudv + Gdv2<br />
(5)<br />
trong đó<br />
(6)<br />
là các hệ số của dạng toàn phương thứ nhất.<br />
c. Dạng toàn phương thứ hai, F2:<br />
Dạng toàn phương thứ hai mô tả độ cong của<br />
một bề mặt tự do, được xác định bởi [35]:<br />
F2= −dnS .dS = Ldu2 + 2Mdudv + Ndv2 (7)<br />
trong đó<br />
(8)<br />
là các hệ số của dạng toàn phương thứ hai.<br />
d. Độ cong Gauss (K) và độ cong trung bình (H)<br />
Cho một bề mặt tự do S(u,v) và p là một điểm<br />
bất kỳ trên nó. Gọi (Q) là mặt phẳng chứa pháp véc<br />
tơ của mặt S tại điểm p. Giao tuyến của Q là S là<br />
một đường cong có độ cong nhất định (hình 2). Khi<br />
mặt (Q) quay xung quanh pháp véc tơ nói trên thì<br />
độ cong của đường cong thay đổi. Ơ-le đã chính<br />
minh rằng tồn tại các hướng mà ở đó độ cong của<br />
đường cong đạt đạt giá trị cực tiểu và cực đại [6].<br />
Các độ cong ở các hướng này được gọi là các độ<br />
cong chính tắc và các hướng độ cong chính tắc<br />
vuông góc nhau.<br />
<br />
66 • TRƯỜNG ĐẠI HỌC NHA TRANG<br />
<br />
Tạp chí Khoa học - Công nghệ Thủy sản<br />
<br />
Số 4/2014<br />
<br />
Hình 2. Độ cong của bề mặt tự do<br />
<br />
Độ cong Gaussian (K), độ cong trung bình (H) và<br />
các độ cong chính tắc Kmin và Kmax của bề mặt S(u, v),<br />
tại điểm p, được tính bằng các công thức [6]:<br />
(9)<br />
<br />
(10)<br />
(11)<br />
(12)<br />
Dựa vào các giá trị của độ cong Gauss, độ cong<br />
trung bình và các độ cong chính tắc, các điểm trên<br />
một bề mặt tự do có thể được chia thành sáu loại<br />
như sau [5]:<br />
- Điểm eliptic lõm:<br />
Nếu K > 0 và H > 0.<br />
- Điểm eliptic lồi:<br />
Nếu K > 0 và H < 0.<br />
- Điểm hyperbolic:<br />
Nếu K 0.<br />
- Điểm parabolic lồi: Nếu K = 0 và H < 0.<br />
- Điểm rốn phẳng:<br />
Nếu K = 0 và H = 0.<br />
Để chia một bề mặt tự do thành các vùng lồi (kể<br />
cả vùng phẳng), vùng lõm và vùng yên ngựa, hình<br />
dạng bề mặt cục bộ quanh một điểm có thể được<br />
chia thành ba loại vùng khác nhau như sau [5]:<br />
* K ≥ 0 và H £ 0: hình dạng bề mặt cục bộ lồi.<br />
* K ≥ 0 và H > 0: hình dạng bề mặt cục bộ lõm.<br />
* K < 0 và H ¹ 0: hình dạng bề mặt cục bộ yên ngựa.<br />
2. Tính độ cong của bề mặt tự do<br />
Trong nghiên cứu này, thuật toán phân vùng bề<br />
mặt tự do như sau:<br />
(a) Tạo tập hợp lưới điểm bề mặt {p} từ mô hình<br />
toán học của bề mặt S và lưu tất cả các điểm vào<br />
một ma trận chung.<br />
(b) Tính các thông số K and H tại mỗi điểm pi,j.<br />
(c) Xét mỗi điểm pi,j thuộc tập {p}:<br />
- Nếu K ≥ 0 và H £ 0: lưu điểm vào ma trận<br />
các điểm vùng lồi, mã hóa điểm lưu thành số 1.<br />
<br />
Các điểm không có tính chất này được mã hóa<br />
thành số 0.<br />
- Nếu K ≥ 0 và H > 0: lưu điểm vào ma trận các<br />
điểm vùng lõm, mã hóa điểm lưu thành số 2. Các<br />
điểm không có tính chất này được mã hóa thành số 0.<br />
- Nếu K < 0: lưu điểm vào ma trận các điểm<br />
vùng yên ngựa, mã hóa điểm lưu thành số 3. Các<br />
điểm không có tính chất này được mã hóa thành<br />
số 0.<br />
Cấu trúc ma trận các vùng riêng biệt nêu trên<br />
tương tự như cấu trúc của ma trận của ảnh nhị<br />
phân. Điều này tạo điều kiện dễ dàng cho việc áp<br />
dụng kỹ thuật mã xích trong xử lý ảnh để xác định<br />
biên các vùng đã phân. Các điểm biên sẽ được<br />
dùng để xây dựng các đường cong không gian ba<br />
chiều dùng cho việc chia bề mặt tự do thành các<br />
vùng khác nhau trong môi trường CAD.<br />
Để thực hiện tính toán theo thuật toán nói trên<br />
và căn cứ các phương trình (9) đến (12), các bước<br />
để tính độ cong của một bề mặt tự do được thực<br />
hiện như sau:<br />
- Tạo tập hợp các điểm trên bề mặt tự do<br />
đã cho.<br />
- Tính toán các pháp véc tơ đơn vị tại tất cả các<br />
điểm trên bề mặt.<br />
- Tính các hệ số dạng toàn phương thứ nhất<br />
và thứ hai.<br />
- Tính độ cong Gaussian, độ cong trung bình và<br />
các độ cong chính tắc.<br />
Trong nghiên cứu này, việc tính toán phân vùng<br />
và xác định biên các vùng được thực hiện bằng một<br />
chương trình Matlab. Chương trình gồm các tập tin<br />
M-function và M-script để tạo mô hình toán của bề<br />
mặt, tính độ cong bề mặt, phân vùng và xác định<br />
biên các vùng. Hàm tính toán độ cong có nội dung<br />
cơ bản như sau:<br />
- Tính đạo hàm bậc nhất và đạo hàm bậc hai<br />
theo các biến u và v bằng cách sử dụng hàm tiêu<br />
chuẩn gradient.<br />
[Xu,Xv] = gradient(X);<br />
[Yu,Yv] = gradient(Y);<br />
[Zu,Zv] = gradient(Z);<br />
[Xuu,Xuv] = gradient(Xu);<br />
[Yuu,Yuv] = gradient(Yu);<br />
[Zuu,Zuv]= gradient(Zu);<br />
[Xuv,Xvv] = gradient(Xv);<br />
[Yuv,Yvv] = gradient(Yv);<br />
[Zuv,Zvv] = gradient(Zv);<br />
X, Y và Z ở đây là những mảng 2 chiều của<br />
các điểm trên bề mặt. Những mảng này phải được<br />
chuyển thành các véc tơ để thực hiện các phép tính<br />
về sau.<br />
<br />
TRƯỜNG ĐẠI HỌC NHA TRANG • 67<br />
<br />
Tạp chí Khoa học - Công nghệ Thủy sản<br />
<br />
Số 4/2014<br />
<br />
- Tính pháp véc tơ tại các điểm trên bề mặt:<br />
m = cross(Xu,Xv,2);<br />
q = sqrt(dot(m,m,2));<br />
- Tính các hệ số dạng toàn phương thứ nhất:<br />
E = dot(Xu,Xu,2);<br />
F = dot(Xu,Xv,2);<br />
- Tính các hệ số dạng toàn phương thứ hai:<br />
L = dot(Xuu,n,2);<br />
M = dot(Xuv,n,2);<br />
- Tính độ cong Gauss:<br />
K = (L.*N - M.^2)./(E.*G - F.^2);<br />
- Tính độ cong trung bình:<br />
H = (E.*N + G.*L - 2.*F.*M)./(2*(E.*G - F.^2));<br />
- Tính độ cong chính tắc:<br />
Kmax = H + sqrt(H.^2 - K);<br />
Kmin = H - sqrt(H.^2 - K);<br />
III. KẾT QUẢ NGHIÊN CỨU VÀ THẢO LUẬN<br />
Trong nghiên cứu này, chương trình Matlab<br />
được viết cho tường minh hoặc bề mặt Bspline. Bài<br />
báo này sử dụng bề mặt tự do ở dạng tường minh<br />
để minh họa cho việc tính toán độ cong bề mặt.<br />
Ví dụ: Cho bề mặt tự do được biểu diễn bởi<br />
phương trình<br />
<br />
, trong đó x và y có<br />
<br />
G = dot(Xv,Xv,2);<br />
N = dot(Xvv,n,2);<br />
<br />
giá trị trong đoạn [-1,3]. Giả sử ma trận điểm lưới<br />
cần tạo trên bề mặt có cỡ là 41´41 theo hai phương<br />
x và y.<br />
Trong nghiên cứu này, chương trình Matlab<br />
được chạy trên máy tính xách tay (Intel Core i5,<br />
1,80GHz, RAM 4 GB) cài đặt hệ điều hành Windows<br />
7. Hình 3 và hình 4 trình bày kết quả tính K, H, Kmax<br />
và Kmin tại một số điểm lưới trên bề mặt.<br />
<br />
Hình 3. Minh họa giá trị K, H tại một số điểm lưới<br />
<br />
Hình 4. Minh họa giá trị Kmax, Kmin tại một số điểm lưới<br />
<br />
68 • TRƯỜNG ĐẠI HỌC NHA TRANG<br />
<br />
n = m./[q q q];<br />
<br />
Tạp chí Khoa học - Công nghệ Thủy sản<br />
Kết quả tính toán cho thấy bề mặt đã cho có<br />
2 vùng lõm, 1 vùng lồi và 3 vùng yên ngựa. Chương<br />
trình cũng cho kết quả tọa độ các điểm biên của<br />
6 vùng này. Tuy nhiên, đối với ví dụ này, chỉ cần sử<br />
dụng các điểm biên của các vùng lõm và vùng lồi để<br />
<br />
Số 4/2014<br />
xây dựng hai đường cong không gian cho mục đích<br />
chia vùng bề mặt trong môi trường CAD. Trên hình 5<br />
là mô hình bề mặt đã cho với các điểm biên của<br />
các vùng lõm và vùng lồi, được hiển thị trong môi<br />
trường Matlab.<br />
<br />
Hình 5. Mô hình bề mặt với các điểm biên của vùng lồi, lõm<br />
<br />
Các điểm biên của vùng lồi được chứa trong<br />
ma trận PB24, gồm 71 điểm (hình 6a). Ma trận<br />
PB13cell{1,1} trên hình 6b là ma trận chứa 27 điểm<br />
biên của vùng lõm thứ nhất và ma trận PB13cell{2,1}<br />
trên hình 6c chứa 60 điểm biên của vùng lõm thứ<br />
hai. Trong các ma trận này, các cột 1, 2 và 3 tương<br />
ứng với tọa độ x, y và z của các điểm.<br />
Trong chương trình tính toán này, các giá trị<br />
tọa độ của các điểm trên bề mặt (ở một ma trận<br />
riêng) cũng như trên các biên có thể dễ dàng được<br />
lưu ở dạng file Excel. Điều này tạo thuận lợi cho<br />
việc nhập dữ liệu xây dựng bề mặt đã cho và các<br />
đường cong biểu diễn biên của các vùng trên bề<br />
<br />
(a)<br />
<br />
mặt bằng các phần mềm CAD như Catia, Creo, NX,<br />
SolidWorks,… Trong môi trường CAD, các đường<br />
cong biên sẽ được sử dụng làm công cụ xén bề mặt<br />
để chia bề mặt nguyên thành các vùng riêng biệt.<br />
Để tăng độ chính xác của các đường cong biên,<br />
có thể tạo ma trận điểm lưới bề mặt với mức độ<br />
điểm dày hơn. Khi đó số lượng các điểm trên các<br />
biên sẽ tăng. Tuy nhiên, lúc này máy tính sẽ tính<br />
toán lâu hơn. Ở ví dụ này, với ma trận điểm 41´41<br />
thì thời gian chạy chương trình là 0,35 giây. Thời<br />
gian chạy chương trình khi ma trận điểm 201´201<br />
là 5,7 giây. Bề mặt có kích thước càng lớn thì thời<br />
gian chạy chương trình của máy tính sẽ càng lâu.<br />
<br />
(b)<br />
(c)<br />
Hình 6. Minh họa tọa độ các điểm biên trên các vùng lồi (a) và lõm (b, c)<br />
<br />
TRƯỜNG ĐẠI HỌC NHA TRANG • 69<br />
<br />