Структура от движения (SfM) является процессом оценки 3-D структуры сцены из набора 2D изображений. Этот пример показывает вам, как оценить положения калиброванной камеры от двух изображений, восстановить 3-D структуру сцены до неизвестного масштабного коэффициента, и затем восстановить фактический масштабный коэффициент путем обнаружения объекта известного размера.
В этом примере показано, как восстановить 3-D сцену из пары 2D изображения, взятые с камерой, калиброванной с помощью приложения Camera Calibrator. Алгоритм состоит из следующих шагов:
Совпадайте с разреженным набором точек между двумя изображениями. Существует несколько способов найти соответствия точки между двумя изображениями. Этот пример обнаруживает углы в первом изображении с помощью detectMinEigenFeatures функция, и отслеживает их во второе изображение с помощью vision.PointTracker. Кроме того, можно использовать extractFeatures сопровождаемый matchFeatures.
Оцените основную матрицу с помощью estimateFundamentalMatrix.
Вычислите движение камеры с помощью cameraPose функция.
Совпадайте с плотным набором точек между двумя изображениями. Повторно обнаружьте точку с помощью detectMinEigenFeatures с уменьшаемым 'MinQuality' понимать больше мыслей. Затем отследите плотные точки во второе изображение с помощью vision.PointTracker.
Определите 3-D местоположения совпадающих точек с помощью triangulate.
Обнаружьте объект известного размера. В этой сцене существует земной шар, радиус которого, как известно, составляет 10 см. Используйте pcfitsphere найти земной шар в облаке точек.
Восстановите фактическую шкалу, приводящую к метрической реконструкции.
Загрузите пару изображений в рабочую область.
imageDir = fullfile(toolboxdir('vision'), 'visiondata','upToScaleReconstructionImages'); images = imageDatastore(imageDir); I1 = readimage(images, 1); I2 = readimage(images, 2); figure imshowpair(I1, I2, 'montage'); title('Original Images');
Этот пример использует параметры камеры, вычисленные cameraCalibrator приложением. Параметры хранятся в cameraParams объект, и включает внутренние параметры камеры и коэффициенты искажения объектива.
% Load precomputed camera parameters load upToScaleReconstructionCameraParameters.mat
Искажение объектива может влиять на точность итоговой реконструкции. Можно удалить искажение из каждого из изображений с помощью undistortImage функция. Этот процесс выправляет линии, которые изгибаются радиальным искажением линзы.
I1 = undistortImage(I1, cameraParams); I2 = undistortImage(I2, cameraParams); figure imshowpair(I1, I2, 'montage'); title('Undistorted Images');
Обнаружьте хорошие функции, чтобы отследить. Уменьшайте 'MinQuality' обнаружить меньше точек, которые были бы более равномерно распределены в изображении. Если движение камеры не является очень большим, то отслеживание использования алгоритма KLT является хорошим способом установить соответствия точки.
% Detect feature points imagePoints1 = detectMinEigenFeatures(rgb2gray(I1), 'MinQuality', 0.1); % Visualize detected points figure imshow(I1, 'InitialMagnification', 50); title('150 Strongest Corners from the First Image'); hold on plot(selectStrongest(imagePoints1, 150));
% Create the point tracker tracker = vision.PointTracker('MaxBidirectionalError', 1, 'NumPyramidLevels', 5); % Initialize the point tracker imagePoints1 = imagePoints1.Location; initialize(tracker, imagePoints1, I1); % Track the points [imagePoints2, validIdx] = step(tracker, I2); matchedPoints1 = imagePoints1(validIdx, :); matchedPoints2 = imagePoints2(validIdx, :); % Visualize correspondences figure showMatchedFeatures(I1, I2, matchedPoints1, matchedPoints2); title('Tracked Features');
Используйте estimateFundamentalMatrix функция, чтобы вычислить основную матрицу и найти точки inlier, которые соответствуют epipolar ограничению.
% Estimate the fundamental matrix [fMatrix, epipolarInliers] = estimateFundamentalMatrix(... matchedPoints1, matchedPoints2, 'Method', 'MSAC', 'NumTrials', 10000); % Find epipolar inliers inlierPoints1 = matchedPoints1(epipolarInliers, :); inlierPoints2 = matchedPoints2(epipolarInliers, :); % Display inlier matches figure showMatchedFeatures(I1, I2, inlierPoints1, inlierPoints2); title('Epipolar Inliers');
Вычислите вращение и перевод между положениями камеры, соответствующими двум изображениям. Обратите внимание на то, что t единичный вектор, потому что перевод может только быть вычислен, чтобы масштабироваться.
[R, t] = cameraPose(fMatrix, cameraParams, inlierPoints1, inlierPoints2);
Повторно обнаружьте точки в первом изображении с помощью более низкого 'MinQuality' понимать больше мыслей. Отследите новые точки во второе изображение. Оцените 3-D местоположения, соответствующие совпадающим точкам с помощью triangulate функция, которая реализует алгоритм Прямого линейного преобразования (DLT) [1]. Поместите источник в оптическом центре камеры, соответствующей первому изображению.
% Detect dense feature points imagePoints1 = detectMinEigenFeatures(rgb2gray(I1), 'MinQuality', 0.001); % Create the point tracker tracker = vision.PointTracker('MaxBidirectionalError', 1, 'NumPyramidLevels', 5); % Initialize the point tracker imagePoints1 = imagePoints1.Location; initialize(tracker, imagePoints1, I1); % Track the points [imagePoints2, validIdx] = step(tracker, I2); matchedPoints1 = imagePoints1(validIdx, :); matchedPoints2 = imagePoints2(validIdx, :); % Compute the camera matrices for each position of the camera % The first camera is at the origin looking along the X-axis. Thus, its % rotation matrix is identity, and its translation vector is 0. camMatrix1 = cameraMatrix(cameraParams, eye(3), [0 0 0]); camMatrix2 = cameraMatrix(cameraParams, R', -t*R'); % Compute the 3-D points points3D = triangulate(matchedPoints1, matchedPoints2, camMatrix1, camMatrix2); % Get the color of each reconstructed point numPixels = size(I1, 1) * size(I1, 2); allColors = reshape(I1, [numPixels, 3]); colorIdx = sub2ind([size(I1, 1), size(I1, 2)], round(matchedPoints1(:,2)), ... round(matchedPoints1(:, 1))); color = allColors(colorIdx, :); % Create the point cloud ptCloud = pointCloud(points3D, 'Color', color);
Используйте plotCamera функция, чтобы визуализировать местоположения и ориентации камеры и pcshow функция, чтобы визуализировать облако точек.
% Visualize the camera locations and orientations cameraSize = 0.3; figure plotCamera('Size', cameraSize, 'Color', 'r', 'Label', '1', 'Opacity', 0); hold on grid on plotCamera('Location', t, 'Orientation', R, 'Size', cameraSize, ... 'Color', 'b', 'Label', '2', 'Opacity', 0); % Visualize the point cloud pcshow(ptCloud, 'VerticalAxis', 'y', 'VerticalAxisDir', 'down', ... 'MarkerSize', 45); % Rotate and zoom the plot camorbit(0, -30); camzoom(1.5); % Label the axes xlabel('x-axis'); ylabel('y-axis'); zlabel('z-axis') title('Up to Scale Reconstruction of the Scene');
Найдите земной шар в облаке точек путем подбора кривой сфере к 3-D точкам с помощью pcfitsphere функция.
% Detect the globe globe = pcfitsphere(ptCloud, 0.1); % Display the surface of the globe plot(globe); title('Estimated Location and Size of the Globe'); hold off
Фактический радиус земного шара составляет 10 см. Можно теперь определить координаты 3-D точек в сантиметрах.
% Determine the scale factor scaleFactor = 10 / globe.Radius; % Scale the point cloud ptCloud = pointCloud(points3D * scaleFactor, 'Color', color); t = t * scaleFactor; % Visualize the point cloud in centimeters cameraSize = 2; figure plotCamera('Size', cameraSize, 'Color', 'r', 'Label', '1', 'Opacity', 0); hold on grid on plotCamera('Location', t, 'Orientation', R, 'Size', cameraSize, ... 'Color', 'b', 'Label', '2', 'Opacity', 0); % Visualize the point cloud pcshow(ptCloud, 'VerticalAxis', 'y', 'VerticalAxisDir', 'down', ... 'MarkerSize', 45); camorbit(0, -30); camzoom(1.5); % Label the axes xlabel('x-axis (cm)'); ylabel('y-axis (cm)'); zlabel('z-axis (cm)') title('Metric Reconstruction of the Scene');
Этот пример показал, как восстановить движение камеры и восстановить 3-D структуру сцены из двух изображений, взятых с калиброванной камерой.
[1] Хартли, Ричард и Эндрю Зиссермен. Несколько просматривают геометрию в компьютерном зрении. Второй выпуск. Кембридж, 2000.
