Структура из движения (SfM) - процесс оценки 3-D структуры сцены по набору 2-D изображений. В этом примере показано, как оценить позы калиброванной камеры по двум изображениям, восстановить 3-D структуру сцены до неизвестного масштабного коэффициента, а затем восстановить действительный масштабный коэффициент путем обнаружения объекта известного размера.
В этом примере показано, как восстановить сцену 3-D из пары изображений 2-D сделанных с помощью камеры, откалиброванной с помощью калибратора камеры. Алгоритм состоит из следующих этапов:
Сопоставление разреженного набора точек между двумя изображениями. Существует несколько способов нахождения точечных соответствий между двумя изображениями. Этот пример определяет углы на первом изображении с помощью 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');
В этом примере используются параметры камеры, рассчитанные с помощью приложения «Калибратор камеры». Параметры сохраняются в cameraParams объект, и включают характеристику камеры и коэффициенты искажения объектива. helpview (fullfile (docroot, 'toolbox', 'matlab _ prog', 'matlab _ prog.map'), 'nested _ functions')
% 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(im2gray(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 функция для вычисления фундаментальной матрицы и нахождения внутренних точек, удовлетворяющих эпиполярному ограничению.
% 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(im2gray(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.