Структура от движения от двух представлений

Структура от движения (SfM) является процессом оценки 3-D структуры сцены из набора 2D изображений. Этот пример показывает вам, как оценить положения калиброванной камеры от двух изображений, восстановить 3-D структуру сцены до неизвестного масштабного коэффициента, и затем восстановить фактический масштабный коэффициент путем обнаружения объекта известного размера.

Обзор

Этот пример показывает, как восстановить 3-D сцену из пары 2D изображений, взятых с камерой, калиброванной с помощью приложения Camera Calibrator. Алгоритм состоит из следующих шагов:

  1. Совпадайте с разреженным набором точек между двумя изображениями. Существует несколько способов найти соответствия точки между двумя изображениями. Этот пример обнаруживает углы в первом изображении с помощью функции detectMinEigenFeatures и отслеживает их во второе изображение с помощью vision.PointTracker. Кроме того, можно использовать extractFeatures, сопровождаемый matchFeatures.

  2. Оцените основную матрицу с помощью estimateFundamentalMatrix.

  3. Вычислите движение камеры с помощью функции relativeCameraPose.

  4. Совпадайте с плотным набором точек между двумя изображениями. Повторно обнаружьте точку с помощью detectMinEigenFeatures с уменьшаемым 'MinQuality', чтобы понять больше мыслей. Затем отследите плотные точки во второе изображение с помощью vision.PointTracker.

  5. Определите 3-D местоположения совпадающих точек с помощью triangulate.

  6. Обнаружьте объект известного размера. В этой сцене существует земной шар, радиус которого, как известно, составляет 10 см. Используйте pcfitsphere, чтобы найти земной шар в облаке точек.

  7. Восстановите фактическую шкалу, приводящую к метрической реконструкции.

Считайте пару изображений

Загрузите пару изображений в рабочую область.

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 и включают камеру intrinsics и коэффициенты искажения объектива.

% 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');

Оцените существенную матрицу

Используйте функцию estimateEssentialMatrix, чтобы вычислить существенную матрицу и найти точки inlier, которые соответствуют epipolar ограничению.

% Estimate the fundamental matrix
[E, epipolarInliers] = estimateEssentialMatrix(...
    matchedPoints1, matchedPoints2, cameraParams, 'Confidence', 99.99);

% Find epipolar inliers
inlierPoints1 = matchedPoints1(epipolarInliers, :);
inlierPoints2 = matchedPoints2(epipolarInliers, :);

% Display inlier matches
figure
showMatchedFeatures(I1, I2, inlierPoints1, inlierPoints2);
title('Epipolar Inliers');

Вычислите положение камеры

Вычислите местоположение и ориентацию второй камеры относительно первой. Обратите внимание на то, что t является единичным вектором, потому что перевод может только быть вычислен, чтобы масштабироваться.

[orient, loc] = relativeCameraPose(E, cameraParams, inlierPoints1, inlierPoints2);

Восстановите 3-D местоположения совпадающих точек

Повторно обнаружьте точки в первом изображении с помощью более низкого 'MinQuality', чтобы понять больше мыслей. Отследите новые точки во второе изображение. Оцените 3-D местоположения, соответствующие совпадающим точкам с помощью функции triangulate, которая реализует алгоритм Прямого линейного преобразования (DLT) [1]. Поместите источник в оптическом центре камеры, соответствующей первому изображению.

% Detect dense feature points. Use an ROI to exclude points close to the
% image edges.
roi = [30, 30, size(I1, 2) - 30, size(I1, 1) - 30];
imagePoints1 = detectMinEigenFeatures(rgb2gray(I1), 'ROI', roi, ...
    '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 Z-axis. Thus, its
% rotation matrix is identity, and its translation vector is 0.
camMatrix1 = cameraMatrix(cameraParams, eye(3), [0 0 0]);

% Compute extrinsics of the second camera
[R, t] = cameraPoseToExtrinsics(orient, loc);
camMatrix2 = cameraMatrix(cameraParams, R, t);

% 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);

Отобразите 3-D облако точек

Используйте функцию 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', loc, 'Orientation', orient, '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);
loc = loc * 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', loc, 'Orientation', orient, '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.