Thuật toán 8 điểm và vì sao chuẩn hoá toạ độ lại quan trọng
1. Bối cảnh
Trong bài trước, ta giả định hai camera đã biết và khôi phục một điểm 3D từ hai ảnh của nó. Trong thực tế, camera hiếm khi được cho sẵn. Thứ bạn có là một tập các cặp điểm tương ứng: các cặp pixel trông như cùng một điểm vật lý trong hai bức ảnh. Bước đầu tiên của hầu hết pipeline structure-from-motion là rút ra từ các cặp này một ma trận 3×3 duy nhất, ma trận cơ bản \(F\), mã hoá quan hệ hình học tương đối giữa hai góc nhìn.
Cách kinh điển để tính nó từ tám cặp điểm trở lên là thuật toán 8 điểm. Longuet-Higgins (1981) đề xuất thuật toán này cho camera đã hiệu chuẩn; sau đó trong một thời gian dài nó bị coi là quá nhạy với nhiễu để dùng được. Hartley (1997) chỉ ra rằng một bước tiền xử lý rất đơn giản — tịnh tiến và co giãn toạ độ pixel — giúp chính thuật toán tuyến tính đó cạnh tranh được với các phương pháp lặp. Bài viết này dẫn xuất thuật toán, cài đặt cả hai phiên bản và đo sự khác biệt.
2. Hình học epipolar và ma trận cơ bản
Xét điểm 3D \(\mathbf X\) được hai camera nhìn thấy tại các điểm ảnh \(\mathbf x=(u,v,1)^\top\) và \(\mathbf x'=(u',v',1)^\top\). Hai tâm camera \(C,C'\) và \(\mathbf X\) xác định một mặt phẳng, gọi là mặt phẳng epipolar. Cả hai điểm ảnh nằm trên giao tuyến của mặt phẳng này với mỗi mặt phẳng ảnh, nên \(\mathbf x'\) phải nằm trên một đường thẳng trong ảnh thứ hai, đường này chỉ phụ thuộc vào \(\mathbf x\) và hình học camera. Đường thẳng nối \(C\) và \(C'\) (baseline) cắt các mặt phẳng ảnh tại các epipole \(\mathbf e,\mathbf e'\).
Sự thật hình học này được biểu diễn đại số bằng ràng buộc epipolar
$$\mathbf x'^\top F\,\mathbf x = 0 ,$$trong đó, với các camera \(P=K[I\mid\mathbf 0]\) và \(P'=K'[R\mid\mathbf t]\), ma trận cơ bản là \(F=K'^{-\top}[\mathbf t]_\times R\,K^{-1}\), với \([\mathbf t]_\times\) là ma trận phản đối xứng của tích có hướng (Hartley & Zisserman, Ch. 9). Hai tính chất quan trọng khi ước lượng: \(F\) chỉ xác định sai khác một hệ số tỉ lệ, và nó có hạng 2 (vì \([\mathbf t]_\times\) có hạng 2), do đó còn 7 bậc tự do. Nếu camera đã hiệu chuẩn, \(E=K'^\top F K\) là ma trận cốt yếu (essential matrix) mà Longuet-Higgins làm việc với.
3. Thuật toán tuyến tính
Viết \(\mathbf x'^\top F\mathbf x=0\) cho một cặp điểm. Với \(F=(f_{ij})\), khai triển cho
$$u'u\,f_{11}+u'v\,f_{12}+u'\,f_{13}+v'u\,f_{21}+v'v\,f_{22}+v'\,f_{23}+u\,f_{31}+v\,f_{32}+f_{33}=0 ,$$tuyến tính theo chín ẩn \(\mathbf f=(f_{11},\dots,f_{33})^\top\). Do đó mỗi cặp điểm cho một hàng \((u'u,\;u'v,\;u',\;v'u,\;v'v,\;v',\;u,\;v,\;1)\), và xếp \(N\ge 8\) cặp điểm ta được
$$A\,\mathbf f=\mathbf 0,\qquad A\in\mathbb R^{N\times 9}.$$Với tám cặp điểm chính xác ở vị trí tổng quát, \(A\) có hạng 8 và nghiệm là vector không gian null của nó. Với dữ liệu nhiễu và \(N\ge 8\), ta giải
$$\hat{\mathbf f}=\arg\min_{\|\mathbf f\|=1}\ \|A\mathbf f\|^2 ,$$nghiệm là vector kỳ dị phải của \(A\) ứng với giá trị kỳ dị nhỏ nhất, giống hệt như trong bài tam giác hoá. Ma trận \(\hat F\) sau khi sắp xếp lại nhìn chung có hạng 3, vi phạm hình học epipolar (các đường epipolar sẽ không cùng đi qua một điểm). Vì vậy ta thay nó bằng ma trận hạng 2 gần nhất theo chuẩn Frobenius: nếu \(\hat F=U\,\mathrm{diag}(\sigma_1,\sigma_2,\sigma_3)\,V^\top\), đặt
$$F=U\,\mathrm{diag}(\sigma_1,\sigma_2,0)\,V^\top .$$4. Vì sao toạ độ pixel thô gây hại: điều kiện số
Hãy nhìn các phần tử của một hàng trong \(A\). Với toạ độ pixel cỡ hàng trăm, \(u'u\) cỡ \(10^5\), các số hạng như \(u\) hay \(u'\) cỡ \(10^2\), còn phần tử cuối là \(1\). Các cột của \(A\) chênh nhau về thang đo tới năm bậc độ lớn. SVD cực tiểu hoá \(\|A\mathbf f\|\) theo nghĩa Euclid, nên thực chất nó gán trọng số lớn hơn nhiều cho các cột lớn so với các cột nhỏ; bài toán điều kiện kém.
Chuẩn hoá Hartley sửa thang đo trước khi dựng \(A\). Với từng ảnh riêng biệt, tìm phép biến đổi đồng dạng \(T\) đưa trọng tâm các điểm về gốc toạ độ và co giãn sao cho khoảng cách trung bình tới gốc bằng \(\sqrt 2\):
$$\tilde{\mathbf x}_i=T\mathbf x_i,\qquad \tilde{\mathbf x}'_i=T'\mathbf x'_i .$$Chạy thuật toán ở trên với các điểm đã chuẩn hoá để có \(\tilde F\), rồi hoàn tác chuẩn hoá. Vì \(\tilde{\mathbf x}'^\top\tilde F\tilde{\mathbf x}=\mathbf x'^\top (T'^\top\tilde F\,T)\,\mathbf x\), kết quả trong toạ độ pixel là
$$F=T'^\top\,\tilde F\,T .$$(Hartley & Zisserman gọi đây là normalized 8-point algorithm, Ch. 11.)
5. Cài đặt
Phần lõi của bản cài đặt (NumPy). Toàn bộ script, gồm cả thí nghiệm và hai hình, có thể tải về.
import numpy as np
def normalise(pts):
"""Hartley normalisation: centroid -> origin, mean distance -> sqrt(2)."""
c = pts.mean(axis=0)
d = np.sqrt(((pts - c) ** 2).sum(axis=1)).mean()
s = np.sqrt(2) / d
T = np.array([[s, 0, -s * c[0]], [0, s, -s * c[1]], [0, 0, 1.0]])
return np.column_stack([pts, np.ones(len(pts))]) @ T.T, T
def eight_point(p1, p2, normalised):
"""p1, p2: (N,2) pixel coordinates, N>=8. Convention: x2^T F x1 = 0."""
if normalised:
h1, T1 = normalise(p1); h2, T2 = normalise(p2)
else:
h1 = np.column_stack([p1, np.ones(len(p1))])
h2 = np.column_stack([p2, np.ones(len(p2))])
T1 = T2 = np.eye(3)
A = np.stack([np.kron(b, a) for a, b in zip(h1, h2)]) # N x 9
_, S, Vt = np.linalg.svd(A)
F = Vt[-1].reshape(3, 3)
U, D, Vt2 = np.linalg.svd(F) # enforce rank 2
F = U @ np.diag([D[0], D[1], 0.0]) @ Vt2
F = T2.T @ F @ T1 # denormalise
return F / np.linalg.norm(F), S
Kiểm tra tính đúng đắn
Tôi sinh 20 cặp điểm không nhiễu từ một \(F_{\text{true}}\) đã biết. Cả hai phiên bản đều khôi phục được nó tới mức sai số làm tròn: \(\|F-F_{\text{true}}\|=4.7\times10^{-15}\) (thô) và \(5.9\times10^{-15}\) (chuẩn hoá), sau khi cố định dấu. Khi không có nhiễu, hai thuật toán đúng như nhau; mọi khác biệt phải đến từ điều kiện số.
6. Thí nghiệm: thô so với chuẩn hoá
Thiết lập (dữ liệu tổng hợp, tái lập hoàn toàn). Tham số nội \(f=800\) px, điểm chính \((320,240)\), giống nhau cho hai camera. Camera 1 là \([I\mid\mathbf 0]\); camera 2 quay \(-8^\circ\) quanh trục y với \(\mathbf t=(-1,\,0.05,\,0.05)\). Các điểm 3D ngẫu nhiên với \(x\in[-1,1]\), \(y\in[-0.75,0.75]\), \(z\in[4,8]\). Mỗi lần thử dùng \(N=20\) cặp điểm, cộng nhiễu Gauss độc lập với độ lệch chuẩn \(\sigma\) vào mọi toạ độ pixel. \(F\) ước lượng được chấm điểm trên 1000 cặp điểm không nhiễu riêng biệt không dùng khi ước lượng, bằng khoảng cách trung bình đối xứng từ điểm đến đường epipolar (đơn vị pixel). Tôi báo cáo trung vị trên 300 lần thử (NumPy seed 11).
Điều kiện số. Với dữ liệu không nhiễu, trung vị tỉ số giữa giá trị kỳ dị lớn nhất và giá trị kỳ dị thứ tám của \(A\) là \(4.6\times10^{5}\) với toạ độ pixel thô và \(18.0\) sau khi chuẩn hoá (seed 5, 300 lần thử).
| Nhiễu σ (px) | 8 điểm thô (px) | 8 điểm chuẩn hoá (px) |
|---|---|---|
| 0.00 | 0.0000 | 0.0000 |
| 0.10 | 0.1943 | 0.0917 |
| 0.50 | 1.8377 | 0.4428 |
| 1.00 | 8.6689 | 0.9326 |
| 2.00 | 15.1767 | 1.8583 |
Các con số cho thấy điều gì
- Chuẩn hoá giảm số điều kiện bốn bậc độ lớn (từ khoảng \(4.6\times10^{5}\) xuống 18) trong thiết lập này, đây là lý do số học khiến hai phiên bản cư xử khác nhau khi có nhiễu.
- Ước lượng chuẩn hoá suy giảm êm. Sai số trên điểm giữ lại gần với mức nhiễu được thêm vào và tăng xấp xỉ tỉ lệ với \(\sigma\): 0.092, 0.443, 0.933 và 1.858 px với \(\sigma=0.1, 0.5, 1, 2\) px.
- Ước lượng thô tệ hơn nhiều và tăng nhanh hơn nhiễu. Tại \(\sigma=1\) px sai số của nó là 8.67 px so với 0.93 px (lớn hơn 9.3 lần); tại \(\sigma=0.5\) lớn hơn 4.2 lần và tại \(\sigma=0.1\) lớn hơn 2.1 lần.
- Không nhiễu thì không có khác biệt, như phép kiểm tra đã cho thấy: chuẩn hoá là về ổn định số học, không phải về một nghiệm khác.
Hạn chế. Dữ liệu là tổng hợp, nhiễu là Gauss đẳng hướng, không có cặp điểm ngoại lai, tôi dùng một hình học camera cố định với \(N=20\) cặp điểm, và chỉ đo khoảng cách epipolar trên các điểm giữ lại. Tôi không so sánh với các bộ ước lượng lặp (hợp lý cực đại); Hartley (1997) báo cáo các so sánh như vậy trên ảnh thật và thấy thuật toán chuẩn hoá tương đương với các phương pháp lặp tốt nhất. Các cặp điểm thật còn chứa ngoại lai, đòi hỏi ước lượng bền vững (RANSAC) bao quanh bộ giải 8 điểm.
7. Bài học thực hành
- Luôn chuẩn hoá toạ độ trước khi chạy một bộ ước lượng tuyến tính (kiểu DLT); chỉ tốn vài dòng code mà loại bỏ một nguồn sai số số học lớn.
- Ép hạng 2 sau SVD và hoàn tác chuẩn hoá ở cuối.
- Kiểm thử bằng dữ liệu tổng hợp không nhiễu trước, rồi đánh giá trên các điểm giữ lại, không phải trên chính các điểm đã khớp.
- Dùng kết quả 8 điểm làm ước lượng khởi tạo; với dữ liệu thật, bọc nó trong RANSAC và tinh chỉnh bằng cách cực tiểu sai số hình học.
- Khi đã biết \(F\) (hoặc \(E\)) và khôi phục được các camera, bạn có thể tam giác hoá các điểm như trong bài trước.
Tôi nhận làm các việc này với tư cách freelancer. Hãy mô tả dự án trong vài câu, tôi sẽ phản hồi bằng các câu hỏi làm rõ và báo giá theo phạm vi cố định.
Nhận báo giá Xem dịch vụ
Tài liệu tham khảo
- H. C. Longuet-Higgins. “A computer algorithm for reconstructing a scene from two projections.” Nature 293:133–135, 10 September 1981. doi:10.1038/293133a0.
- R. I. Hartley. “In defense of the eight-point algorithm.” IEEE Transactions on Pattern Analysis and Machine Intelligence 19(6):580–593, June 1997. doi:10.1109/34.601246.
- R. Hartley and A. Zisserman. Multiple View Geometry in Computer Vision, 2nd ed. Cambridge University Press, 2004. ISBN 978-0-521-54051-3. (Chương 9, hình học epipolar và ma trận cơ bản; Chương 11, tính ma trận cơ bản.) Trang nhà xuất bản.