Как преобразовать 2D-точки в 3D?

Ответ 1

Хорошо, я пришел сюда, чтобы найти ответ, и не нашел ничего простого и простого, поэтому я пошел дальше и сделал тупым, но эффективным (и относительно простым): оптимизация в Монте-Карло.

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

Вот фото из Thomas the Tank Engine:

Скажем, мы используем GIMP, чтобы найти 2D-координаты того, что мы считаем квадратом на плоскости земли (независимо от того, действительно ли это квадрат, зависит от вашего суждения о глубине):

Я получаю четыре точки в 2D-изображении: (318, 247), (326, 312), (418, 241) и (452, 303).

По соглашению мы говорим, что эти точки должны соответствовать трехмерным точкам: (0, 0, 0), (0, 0, 1), (1, 0, 0) и (1, 0, 1). Другими словами, единичный квадрат в плоскости y = 0.

Проецирование каждой из этих трехмерных координат в 2D выполняется путем умножения 4D-вектора [x, y, z, 1] на матрицу проецирования 4x4, а затем разделение x и y на z на фактическую коррекцию перспективы. Это более или менее то, что делает gluProject(), за исключением того, что gluProject() также учитывает текущий видовой экран и учитывает отдельную матрицу просмотра модели ( мы можем просто предположить, что матрица вида модели является единичной матрицей). Очень удобно смотреть документацию gluProject(), потому что я действительно хочу решение, которое работает для OpenGL, но остерегайтесь того, что в документации отсутствует разделение на z в формуле.

Помните, что алгоритм должен начинаться с некоторой матрицы проекций и произвольно возмущать его, пока он не даст желаемую проекцию. Итак, что мы собираемся делать, это спроектировать каждый из четырех трехмерных точек и посмотреть, как близко мы добираемся до 2D-точек, которые мы хотели. Если наши случайные возмущения заставляют проецируемые 2D-точки приблизиться к тем, которые были отмечены выше, то мы сохраняем эту матрицу как улучшение по сравнению с нашей начальной (или предыдущей) догадкой.

Определим наши точки:

# Known 2D coordinates of our rectangle
i0 = Point2(318, 247)
i1 = Point2(326, 312)
i2 = Point2(418, 241)
i3 = Point2(452, 303)

# 3D coordinates corresponding to i0, i1, i2, i3
r0 = Point3(0, 0, 0)
r1 = Point3(0, 0, 1)
r2 = Point3(1, 0, 0)
r3 = Point3(1, 0, 1)

Нам нужно начать с некоторой матрицы, тождественная матрица кажется естественным выбором:

mat = [
    [1, 0, 0, 0],
    [0, 1, 0, 0],
    [0, 0, 1, 0],
    [0, 0, 0, 1],
]

Нам нужно фактически реализовать проекцию (которая в основном является матричным умножением):

def project(p, mat):
    x = mat[0][0] * p.x + mat[0][1] * p.y + mat[0][2] * p.z + mat[0][3] * 1
    y = mat[1][0] * p.x + mat[1][1] * p.y + mat[1][2] * p.z + mat[1][3] * 1
    w = mat[3][0] * p.x + mat[3][1] * p.y + mat[3][2] * p.z + mat[3][3] * 1
    return Point(720 * (x / w + 1) / 2., 576 - 576 * (y / w + 1) / 2.)

Это в основном то, что gluProject(), 720 и 576 - ширина и высота изображения соответственно (т.е. окно просмотра), и мы вычитаем из 576, чтобы подсчитать, что мы подсчитали y координат сверху OpenGL обычно учитывает их снизу. Вы заметите, что мы не вычисляем z, потому что нам это действительно не нужно (хотя было бы удобно, чтобы он попадал в диапазон, который OpenGL использует для буфера глубины).

Теперь нам нужна функция для оценки того, насколько мы близки к правильному решению. Значение, возвращаемое этой функцией, - это то, что мы будем использовать, чтобы проверить, лучше ли одна матрица, чем другая. Я решил пойти по сумме квадратов расстояний, то есть:

# The squared distance between two points a and b
def norm2(a, b):
    dx = b.x - a.x
    dy = b.y - a.y
    return dx * dx + dy * dy

def evaluate(mat): 
    c0 = project(r0, mat)
    c1 = project(r1, mat)
    c2 = project(r2, mat)
    c3 = project(r3, mat)
    return norm2(i0, c0) + norm2(i1, c1) + norm2(i2, c2) + norm2(i3, c3)

Чтобы возмутить матрицу, мы просто выбираем элемент для возмущения случайным количеством в некотором диапазоне:

def perturb(amount):
    from copy import deepcopy
    from random import randrange, uniform
    mat2 = deepcopy(mat)
    mat2[randrange(4)][randrange(4)] += uniform(-amount, amount)

(Стоит отметить, что наша функция project() фактически не использует mat[2] вообще, так как мы не вычисляем z, и поскольку все наши координаты y равны 0, значения mat[*][1] также не имеют значения. Мы могли бы использовать этот факт и никогда не пытаться нарушать эти значения, что дало бы небольшое ускорение, но это оставлено как упражнение...)

Для удобства добавим функцию, которая выполняет основную часть аппроксимации, вызывая perturb() снова и снова на том, какая лучшая матрица мы нашли до сих пор:

def approximate(mat, amount, n=100000):
    est = evaluate(mat)

    for i in xrange(n):
        mat2 = perturb(mat, amount)
        est2 = evaluate(mat2)
        if est2 < est:
            mat = mat2
            est = est2

    return mat, est

Теперь все, что осталось сделать, это запустить его...:

for i in xrange(100):
    mat = approximate(mat, 1)
    mat = approximate(mat, .1)

Я нахожу, что это уже дает довольно точный ответ. После некоторого времени работы найденная мной матрица была:

[
    [1.0836000765696232,  0,  0.16272110011060575, -0.44811064935115597],
    [0.09339193527789781, 1, -0.7990570384334473,   0.539087345090207  ],
    [0,                   0,  1,                    0                  ],
    [0.06700844759602216, 0, -0.8333379578853196,   3.875290562060915  ],
]

с ошибкой около 2.6e-5. (Обратите внимание, что элементы, которые мы сказали, не использовались при вычислении, фактически не были изменены из нашей исходной матрицы, потому что изменение этих элементов не изменило бы результат оценки, и поэтому изменение никогда не будет перенесено вдоль.)

Мы можем передать матрицу в OpenGL с помощью glLoadMatrix() (но не забудьте сначала перенести его и не забудьте загрузить вашу модельную матрицу с идентификационной матрицей):

def transpose(m):
    return [
        [m[0][0], m[1][0], m[2][0], m[3][0]],
        [m[0][1], m[1][1], m[2][1], m[3][1]],
        [m[0][2], m[1][2], m[2][2], m[3][2]],
        [m[0][3], m[1][3], m[2][3], m[3][3]],
    ]

glLoadMatrixf(transpose(mat))

Теперь мы можем, например, перевести вдоль оси z, чтобы получить разные позиции вдоль дорожек:

glTranslate(0, 0, frame)
frame = frame + 1

glBegin(GL_QUADS)
glVertex3f(0, 0, 0)
glVertex3f(0, 0, 1)
glVertex3f(1, 0, 1)
glVertex3f(1, 0, 0)
glEnd()

Конечно, это не очень элегантно с математической точки зрения; вы не получаете уравнения закрытой формы, которое вы можете просто подключить к своим номерам и получить прямой (и точный) ответ. ОДНАКО, это позволяет вам добавлять дополнительные ограничения, не беспокоясь о усложнении своих уравнений; например, если бы мы хотели включить высоту, мы могли бы использовать этот уголок дома и сказать (в нашей оценочной функции), что расстояние от земли до крыши должно быть таким-то и повторить алгоритм. Так что да, это грубая сила, но работает и работает хорошо.

Ответ 2

Д. DeMenthon разработал алгоритм вычисления позы объекта (его положение и ориентацию в пространстве) из точек функции в 2D-изображении при знании модели объекта - это ваша точная проблема:

Мы описываем метод нахождения позы объекта из одного изображения. Мы предполагаем, что мы можем обнаруживать и сопоставлять в изображении четыре или более некомпланарных точек объекта, и что мы знаем их относительную геометрию на объекте.

Алгоритм известен как Позит и описан в его классической статье "Постановка объекта на основе модели в 25 строк кода" (доступно на его веб-сайт, раздел 4).

Прямая ссылка на статью: http://www.cfar.umd.edu/~daniel/daniel_papersfordownload/Pose25Lines.pdf Реализация OpenCV: http://opencv.willowgarage.com/wiki/Posit

Идея состоит в том, чтобы неоднократно аппроксимировать перспективную проекцию масштабированной орфографической проекцией, пока она не сходится к точной позе.

Ответ 3

Это классическая проблема для расширенной реальности на основе маркеров.

У вас есть квадратный маркер (2D-штрих-код), и вы хотите найти его Позу (перевод и поворот по отношению к камере), после нахождения четырех краев маркера. Обзор-Картинка

Я не знаю последних вкладов в эту область, но, по крайней мере, до одной точки (2009) RPP должен был превзойти POSIT, о котором упоминалось выше (и это действительно классический подход для этого) См. Ссылки, они также предоставляют источник.

• http://www.emt.tugraz.at/~vmg/schweighofer

• http://www.emt.tugraz.at/publications/EMT_TR/TR-EMT-2005-01.pdf

• http://www.emt.tugraz.at/system/files/rpp_MATLAB_ref_implementation.tar.gz

(PS - Я знаю это немного старая тема, но в любом случае сообщение может быть полезно кому-то)

Ответ 4

В двухмерном пространстве будут два правильных прямоугольника, которые можно построить. Не зная исходную матричную проекцию, вы не будете знать, какой из них правильный. Это то же самое, что проблема с "ящиком": вы видите два квадрата, один внутри другого, с четырьмя внутренними вершинами, связанными с 4 соответствующими внешними вершинами. Вы смотрите на коробку сверху вниз или снизу вверх?

При этом вы ищете матричное преобразование T, где...

{{x1, y1, z1}, {x2, y2, z2}, {x3, y3, z3}, {x4, y4, z4}} x T = {{x1, y1}, {x2, y2 }, {x3, y3}, {x4, y4}}

(4 x 3) x T = (4 x 2)

Итак, T должна быть матрицей (3 x 2). Итак, у нас есть 6 неизвестных.

Теперь построим систему ограничений на T и решим с помощью Simplex. Чтобы построить ограничения, вы знаете, что строка, проходящая через первые две точки, должна быть параллельна линии, проходящей во вторую две точки. Вы знаете, что линия, проходящая через точки 1 и 3, должна быть параллельна линиям, проходящим через точки 2 и 4. Вы знаете, что линия, проходящая через 1 и 2, должна быть ортогональной линии, проходящей через точки 2 и 3. Вы знаете, что длина строки из 1 и 2 должны равняться длине линии от 3 и 4. Вы знаете, что длина линии от 1 до 3 должна быть равна длине линии от 2 до 4.

Чтобы сделать это еще проще, вы знаете о прямоугольнике, поэтому вы знаете длину всех сторон.

Это должно дать вам множество ограничений для решения этой проблемы.

Конечно, чтобы вернуться, вы можете найти T-обратное.

@Rob: Да, существует бесконечное количество проекций, но не бесконечное количество проектов, где точки должны удовлетворять требованиям прямоугольника.

@nlucaroni: Да, это разрешимо только в том случае, если у вас есть четыре точки в проекции. Если прямоугольник проецируется всего на 2 точки (т.е. Плоскость прямоугольника ортогональна поверхности проекции), то это не может быть разрешено.

Хммм... Я должен пойти домой и написать этот маленький камень. Это звучит весело.

Обновление:

• Существует бесконечное количество прогнозов, если вы не исправите одну из точек. Если вы фиксируете точки исходного прямоугольника, тогда возможны два возможных исходных прямоугольника.

Ответ 5

Для моего движка OpenGL следующий снипп преобразует координаты мыши/экрана в 3D-координаты мира. Прочитайте комментарии для фактического описания того, что происходит.

/*   FUNCTION:        YCamera :: CalculateWorldCoordinates
     ARGUMENTS:       x         mouse x coordinate
                      y         mouse y coordinate
                      vec       where to store coordinates
     RETURN:          n/a
     DESCRIPTION:     Convert mouse coordinates into world coordinates
*/

void YCamera :: CalculateWorldCoordinates(float x, float y, YVector3 *vec) { // START GLint viewport[4]; GLdouble mvmatrix[16], projmatrix[16];

GLint real_y;
GLdouble mx, my, mz;

glGetIntegerv(GL_VIEWPORT, viewport);
glGetDoublev(GL_MODELVIEW_MATRIX, mvmatrix);
glGetDoublev(GL_PROJECTION_MATRIX, projmatrix);

real_y = viewport[3] - (GLint) y - 1;   // viewport[3] is height of window in pixels
gluUnProject((GLdouble) x, (GLdouble) real_y, 1.0, mvmatrix, projmatrix, viewport, &mx, &my, &mz);

/*  'mouse' is the point where mouse projection reaches FAR_PLANE.
    World coordinates is intersection of line(camera->mouse) with plane(z=0) (see LaMothe 306)

    Equation of line in 3D:
        (x-x0)/a = (y-y0)/b = (z-z0)/c      

    Intersection of line with plane:
        z = 0
        x-x0 = a(z-z0)/c  <=> x = x0+a(0-z0)/c  <=> x = x0 -a*z0/c
        y = y0 - b*z0/c

*/
double lx = fPosition.x - mx;
double ly = fPosition.y - my;
double lz = fPosition.z - mz;
double sum = lx*lx + ly*ly + lz*lz;
double normal = sqrt(sum);
double z0_c = fPosition.z / (lz/normal);

vec->x = (float) (fPosition.x - (lx/normal)*z0_c);
vec->y = (float) (fPosition.y - (ly/normal)*z0_c);
vec->z = 0.0f;

} Код >

Ответ 6

Предполагая, что точки действительно являются частью прямоугольника, я даю общую идею:

Найдите две точки с максимальным расстоянием между ними: они, скорее всего, определяют диагональ (исключение: специальные случаи, когда прямоугольник почти паралелирован к плоскости YZ, слева для ученика). Назовите их A, C. Вычислите углы BAD, BCD. Эти, по сравнению с прямыми углами, дают вам ориентацию в 3d-пространстве. Чтобы узнать о расстоянии z, вам необходимо сопоставить проецируемые стороны с известными сторонами, а затем, основываясь на методе 3d-проекции (это 1/z?), Вы на правильном пути, чтобы знать расстояния.

Ответ 7

Чтобы отслеживать подход Rons: вы можете найти свои z-значения, если знаете, как вы повернули свой прямоугольник.

Фокус в том, чтобы найти проективную матрицу, которая сделала проекцию. К счастью, это возможно и даже дешево. Соответствующую математику можно найти в статье "Проективные сопоставления для искажения изображения" Пола Хекберта.

http://pages.cs.wisc.edu/~dyer/cs766/readings/heckbert-proj.pdf

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

Теперь у вас все еще осталось четыре строки вместо точек (как объяснил Рон). Поскольку вы знаете размер исходного прямоугольника, но ничего не потеряно. Теперь вы можете подключить данные из метода Рона и от 2D-подхода к линейному решателю уравнений и решить для z. Таким образом вы получаете точные значения z каждой вершины.

Примечание. Это просто работает, потому что:

• Первоначальная форма была прямоугольником
• Вы знаете точный размер прямоугольника в трехмерном пространстве.

Это особый случай.

Надеюсь, что это поможет, Nils

Ответ 8

Я получу свою линейную книгу Алгебры, когда вернусь домой, если никто не ответит. Но @D G, не все матрицы обратимы. Сингулярные матрицы не обратимы (когда определитель = 0). Это будет происходить все время, так как проекционная матрица должна иметь собственные значения 0 и 1 и быть квадратной (так как она идемпотентна, поэтому p ^ 2 = p).

Простым примером является [[0 1] [0 1]], так как определитель = 0, и это проекция на прямую x = y!

Ответ 9

Проецирование, которое вы имеете на 2D-поверхность, имеет бесконечно много 3D-прямоугольников, которые будут проецироваться в одну и ту же 2D-форму.

Подумайте об этом так: у вас есть четыре 3D-точки, которые составляют 3D-прямоугольник. Вызовите их (x0, y0, z0), (x1, y1, z1), (x2, y2, z2) и (x3, y3, z3). Когда вы проецируете эти точки на плоскость x-y, вы отбрасываете координаты z: (x0, y0), (x1, y1), (x2, y2), (x3, y3).

Теперь вы хотите проецировать обратно в 3D-пространство, вам нужно перепроектировать то, что было z0,.., z3. Но любой набор координат z, который a) сохраняет одно и то же расстояние x-y между точками, и b) сохраняет форму прямоугольником. Итак, любой член этого (бесконечного) набора будет делать: {(z0 + i, z1 + i, z2 + i, z3 + i) | я < - R}.

Изменить @Jarrett: Представьте, что вы решили это и закончили с прямоугольником в 3D-пространстве. Теперь представьте скольжение этого прямоугольника вверх и вниз по оси z. Эти бесконечные количества переведенных прямоугольников имеют одну и ту же проекцию x-y. Откуда вы знаете, что нашли "правильный"?

Редактировать # 2: Хорошо, это из комментария, который я сделал по этому вопросу, - более интуитивный подход к рассуждению об этом.

Представьте, что вы держите лист бумаги над своим столом. Притворяйтесь, что в каждом углу бумаги прикреплена невесомая лазерная указка, указывающая вниз на стол. Бумага представляет собой трехмерный объект, а точки лазерного указателя на столе - это 2D-проекция.

Теперь, как вы можете сказать, насколько высок на столе бумага, глядя на точки лазерной указки?

Вы не можете. Переместите бумагу вверх и вниз. Лазерные указатели по-прежнему будут сиять на тех же местах на столе, независимо от высоты бумаги.

Поиск z-координат в обратном проецировании аналогично попытке найти высоту бумаги на основе точек лазерной указки только на столе.

Ответ 10

При проецировании с 3D на 2D вы теряете информацию.

В простом случае одной точки обратная проекция даст вам бесконечный луч через трехмерное пространство.

Стереоскопическая реконструкция, как правило, начинается с двух изображений 2d и проецирует их обратно в 3D. Затем найдите пересечение двух полученных трехмерных лучей.

Проецирование может иметь разные формы. Ортогональный или перспективный. Я предполагаю, что вы предполагаете ортогональную проекцию?

В вашем случае, если у вас есть оригинальная матрица, у вас будет 4 луча в 3D-пространстве. Тогда вы сможете ограничить проблему своими измерениями 3d-прямоугольника и попытаться решить эту проблему.

Решение не будет уникальным, так как вращение вокруг любой оси, параллельной плоскости проекции 2d, будет неоднозначным по направлению. Другими словами, если 2d-изображение перпендикулярно оси z, то вращение 3d-прямоугольника по часовой стрелке или против часовой стрелки вокруг оси x даст одно и то же изображение. Аналогично для оси y.

В случае, когда плоскость прямоугольника параллельна оси z, у вас еще больше решений.

Поскольку у вас нет исходной матрицы прогноза, дополнительная двусмысленность вводится с помощью произвольного масштабного коэффициента, который существует в любой проекции. Вы не можете различать масштабирование в проекции и перевод в 3d в направлении оси z. Это не проблема, если вас интересуют только относительные позиции 4-х точек в трехмерном пространстве, когда они связаны друг с другом, а не с плоскостью 2-й проекции.

В перспективной проекции все становится сложнее...

Ответ 11

Если вы знаете, что форма представляет собой прямоугольник в плоскости, вы можете значительно ограничить проблему. Вы, конечно же, не можете вычислить "какую" плоскость, поэтому вы можете выбрать, что она лежит на плоскости, где z = 0, а один из углов - при x = y = 0, а ребра параллельны оси x/y.

Таким образом, точки в 3d являются {0,0,0}, {w, 0,0}, {w, h, 0} и {0, h, 0}. Я почти уверен, что абсолютный размер не будет найден, поэтому только отношение w/h освобождается, поэтому это неизвестно.

Относительно этой плоскости камера должна находиться в некоторой точке cx, cy, cz в пространстве, должна указывать направление nx, ny, nz (вектор длины один, так что один из них является избыточным) и имеют focal_length/image_width фактор w. Эти числа превращаются в матрицу проекций 3x3.

Это дает в общей сложности 7 неизвестных: w/h, cx, cy, cz, nx, ny и w.

У вас есть всего 8 известных: пары 4 x + y.

Итак, это можно решить.

Следующий шаг - использовать Matlab или Mathmatica.

Ответ 12

http://library.wolfram.com/infocenter/Articles/2794/

http://campar.in.tum.de/Students/SepPoseEstimation

http://opencvlibrary.sourceforge.net/Posit будет сортировать работу, но это может расходиться или циклироваться.

Ответ 13

Да, Монте-Карло работает, но я нашел лучшее решение для этой проблемы. Этот код работает отлично (и использует OpenCV):

Cv2.CalibrateCamera(new List<List<Point3f>>() { points3d }, new List<List<Point2f>>() { points2d }, new Size(height, width), cameraMatrix, distCoefs, out rvecs, out tvecs, CalibrationFlags.ZeroTangentDist | CalibrationFlags.FixK1 | CalibrationFlags.FixK2 | CalibrationFlags.FixK3);

Эта функция принимает известные 3d и 2d точки, размер экрана и возвращает вращение (rvecs [0]), перевод (tvecs [0]) и матрицу внутренних значений камеры. Это все, что вам нужно.

Ответ 14

Благодаря @Vegard за отличный ответ. Я немного очистил код:

import pandas as pd
import numpy as np

class Point2:
    def __init__(self,x,y):
        self.x = x
        self.y = y

class Point3:
    def __init__(self,x,y,z):
        self.x = x
        self.y = y
        self.z = z

# Known 2D coordinates of our rectangle
i0 = Point2(318, 247)
i1 = Point2(326, 312)
i2 = Point2(418, 241)
i3 = Point2(452, 303)

# 3D coordinates corresponding to i0, i1, i2, i3
r0 = Point3(0, 0, 0)
r1 = Point3(0, 0, 1)
r2 = Point3(1, 0, 0)
r3 = Point3(1, 0, 1)

mat = [
    [1, 0, 0, 0],
    [0, 1, 0, 0],
    [0, 0, 1, 0],
    [0, 0, 0, 1],
]

def project(p, mat):
    #print mat
    x = mat[0][0] * p.x + mat[0][1] * p.y + mat[0][2] * p.z + mat[0][3] * 1
    y = mat[1][0] * p.x + mat[1][1] * p.y + mat[1][2] * p.z + mat[1][3] * 1
    w = mat[3][0] * p.x + mat[3][1] * p.y + mat[3][2] * p.z + mat[3][3] * 1
    return Point2(720 * (x / w + 1) / 2., 576 - 576 * (y / w + 1) / 2.)

# The squared distance between two points a and b
def norm2(a, b):
    dx = b.x - a.x
    dy = b.y - a.y
    return dx * dx + dy * dy

def evaluate(mat): 
    c0 = project(r0, mat)
    c1 = project(r1, mat)
    c2 = project(r2, mat)
    c3 = project(r3, mat)
    return norm2(i0, c0) + norm2(i1, c1) + norm2(i2, c2) + norm2(i3, c3)    

def perturb(mat, amount):
    from copy import deepcopy
    from random import randrange, uniform
    mat2 = deepcopy(mat)
    mat2[randrange(4)][randrange(4)] += uniform(-amount, amount)
    return mat2

def approximate(mat, amount, n=1000):
    est = evaluate(mat)
    for i in xrange(n):
        mat2 = perturb(mat, amount)
        est2 = evaluate(mat2)
        if est2 < est:
            mat = mat2
            est = est2

    return mat, est

for i in xrange(1000):
    mat,est = approximate(mat, 1)
    print mat
    print est

Приблизительный вызов с .1 не сработал у меня, поэтому я взял его. Я тоже запустил его, и, наконец, я проверил его на

[[0.7576315397559887, 0, 0.11439449272592839, -0.314856490473439], 
[0.06440497208710227, 1, -0.5607502645413118, 0.38338196981556827], 
[0, 0, 1, 0], 
[0.05421620936883742, 0, -0.5673977598434641, 2.693116299312736]]

с ошибкой около 0,02.