Разложение матрицы аффинного преобразования

Не так давно в процессе разработки редактора 2D-графики возникла задача разложить матрицу аффинного преобразования на плоскости, на произведение матриц простых преобразований с тем, чтобы отобразить их пользователю и предложить какую-то более-менее адекватную интерпретацию того, что произошло с объектом на канвасе. Честно говоря, эта задача вызвала у меня определенные трудности. Университет я закончил уже давно, и мне было непонятно, а возможно ли это сделать в принципе, учитывая, что исходная матрица могла быть результатом произвольной последовательности сдвигов, масштабов, поворотов, и переносов, причем каждое преобразование могло иметь свой произвольный центр. И, во-вторых, непонятно было, как найти семь параметров, имея всего шесть коэффициентов матрицы. Ключом к решению этой задачи оказалась статья "Разложение матрицы центроаффинного преобразования для нормализации изображения"¹, в которой рассматривается такая же задача, но без учета преобразования переноса и для преобразований относительно центра координат. Далее я фактически просто адаптирую результаты этой статьи с учетом переноса и для произвольного центра преобразований.

Итак, пусть матрица, задающая произвольное преобразование на плоскости, имеет вид:

    ┌           ┐
    │ a11 a12 0 │
M = │ a21 a22 0 │.
    │ a31 a32 1 │
    └           ┘

Ее определитель

det = a11⋅a22 - a12⋅a21, det ≠ 0.

Пусть точка на плоскости задается вектор-строкой вида (x, y, 1), а ее преобразование – умножением справа на матрицу преобразования:

                            ┌           ┐
              ┌         ┐   │ a11 a12 0 │
p1 = p0 • M = │ x0 y0 1 │ • │ a21 a22 0 │
              └         ┘   │ a31 a32 1 │
                            └           ┘

В упомятнутой выше статье¹ утверждается, что произвольная матрица M центроаффинного преобразования может быть представлена как произведение матриц R поворота, матрицы Hx сдвига вдоль оси X и матрицы S масштабирования:

Mc = R • Hx • S

Здесь

    ┌                  ┐
    │  cos(α) sin(α) 0 │
R = │ -sin(α) cos(α) 0 │,
    │     0      0   1 │
    └                  ┘
     ┌         ┐
     │ 1  hx 0 │
Hx = │ 0  1  0 │,
     │ 0  0  1 │
     └         ┘
    ┌         ┐
    │ sx 0  0 │
S = │ 0  sy 0 │,
    │ 0  0  1 │
    └         ┘

При разложении произвольной матрица аффинного преобразования необходимо привести центр трансформации к центру координат, а также учесть преобразование переноса:

M = T0 • S • H • R • T0⁻¹ • T,
     ┌             ┐
     │ 1    0    0 │
T0 = │ 0    1    0 │,
     │ -tx0 -ty0 1 │
     └             ┘
       ┌           ┐ 
       │ 1   0   0 │
T0⁻¹ = │ 0   1   0 │,
       │ tx0 ty0 1 │
       └           ┘
    ┌         ┐ 
    │ 1  0  0 │
T = │ 0  1  0 │.
    │ tx ty 1 │
    └         ┘

Подставив выражения для матриц простых преобразований получим следующие выражения для коэффициентов M:

a11 = sx⋅cos(α) - sx⋅hx⋅sin(α)
a12 = sx⋅hx⋅cos(α) + sx⋅sin(α)
a21 = -sy⋅sin(α)
a22 = sy⋅cos(α)
a31 = tx + tx0⋅(1 - (sx⋅cos(α) - sx⋅hx⋅sin(α))) + ty0⋅sy⋅sin(α) 
    = tx + tx0⋅(1 - a11) - ty0⋅a21
a32 = ty + ty0⋅(1 - cos(α)⋅sy) - tx0⋅(sx⋅hx⋅cos(α) + sx⋅sin(α))
    = ty + ty0⋅(1 - a22) - tx0⋅a12

Решая данные уравнения относительно sx, sy, α, hx, tx, ty, после некоторых упрощений, получим выражения для искомых параметров:

if a22=0
   α = π/2,
   sy = -a21
else
   α = atan(-a21/a22),
   sy = a22/cos(α),

sx = det(M)/sy,
hx = (a11⋅a21 + a12⋅a22)/det,

tx = a31 + ty0⋅a21 + tx0⋅(a11 - 1),
ty = a32 + tx0⋅a12 + ty0⋅(a22 - 1).

Выражения для α, sx, sy, hx сооветствуют аналогичным в статье¹, хотя и несколько отличаются от них по форме. Кроме того мы получили формулы расчета параметров преобразования переноса tx и ty. Хочется также заметить, что даже если в оригинальной последовательности присутствовали сдвиги вдоль обеих осей, в разложении достаточно лишь сдвига вдоль одной из осей (здесь – вдоль оси X). Кроме того, поскольку угол поворота определен как результат функиции арктангенса, он принципиально ограничен значениями от -90˚ до +90˚. Учитывая также, что угол поворота на 180˚ соответсвует sx=-1 и sx=-1, мы имеем здесь некоторую неоднозначность. Например, изначально имея поворот на 120˚ при разложении по данному алгоритму мы получим -60˚ и sx=sy=-1.



¹) Путятин Е.П., Яковлева Е.В., Любченко В.А. «Разложение матрицы центроаффинного преобразования для нормализации изображения»
Поделиться публикацией
AdBlock похитил этот баннер, но баннеры не зубы — отрастут

Подробнее
Реклама
Комментарии 16
  • +1
    Пожалейте неучей, покажите картинки! Вы же редактор 2D графики делаете, должны же быть иллюстрации?
    • 0
      Я его делал. :) С задачей этой я разобрался, но в редакторе не реализовал. Поэтому показывать вроде и нечего.
      • +10
        Можно вопрос, зачем вы всё транспонировали? Традиционно геометрический вектор обозначается вектором-столбцом, а у матрицы преобразования строка (0, 0, 1) — внизу, а не справа.
        Вектор-строка же обозначает обычно ковектор, который принципиально другой объект.

        Что же касается формул, держите и больше так не делайте (  для отступов я уж не стал делать, можете сами):

        Картинки
        M = \begin{pmatrix}\ a_{11} & a_{12} & 0 \\ a_{21} & a_{22} & 0 \\ a_{31} & a_{32} & 1\end{pmatrix}
        \operatorname{det} M = a_{11}\,a_{22} - a_{12}\,a_{21},\qquad \operatorname{det} M \ne 0
        p_1 = p_0 \cdot M = (x_0, y_0, 1) \begin{pmatrix}\ a_{11} & a_{12} & 0 \\ a_{21} & a_{22} & 0 \\ a_{31} & a_{32} & 1\end{pmatrix}
        M_c = R \cdot H_x \cdot S
        R = \begin{pmatrix}\ \cos \alpha & \sin \alpha & 0 \\ -\sin \alpha & \cos \alpha & 0 \\ 0 & 0 & 1\end{pmatrix}
        H_x = \begin{pmatrix}\ 1 & hx & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1\end{pmatrix}
        S = \begin{pmatrix}\ sx & 0 & 0 \\ 0 & sy & 0 \\ 0 & 0 & 1\end{pmatrix}
        M = T_0 \cdot S \cdot H \cdot R \cdot T_0^{-1} \cdot T
        T_0 = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ -tx_0 & -ty_0 & 1\end{pmatrix}
        T_0^{-1} = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ tx_0 & ty_0 & 1\end{pmatrix}
        T = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ tx & ty & 1\end{pmatrix}
        a_{11} = sx \cdot \cos \alpha - sx \cdot hx \cdot \sin \alpha
        a_{12} = sx \cdot hx \cdot \cos \alpha + sx \cdot \sin \alpha
        a_{21} = -sy \cdot \sin \alpha
        a_{22} = sy \cdot \cos \alpha
        a_{31} = tx + tx_0 \cdot (1 - (sx \cdot \cos \alpha - sx \cdot hx \cdot \sin \alpha)) + ty_0 \cdot sy \cdot \sin \alpha
        = tx + tx_0 \cdot (1 - a_{11}) - ty_0 \cdot a_{21}
        a_{32} = ty + ty_0 \cdot (1 - \cos \alpha \cdot sy) - tx_0 \cdot (sx \cdot hx \cdot \cos \alpha + sx \cdot \sin \alpha)
        = ty + ty_0 \cdot (1 - a_{22}) - tx_0 \cdot a_{12}
        \operatorname{if} a_{22}=0
        \alpha = \frac{\pi}{2},
        sy = -a_{21}
        \operatorname{else}
        \alpha = \operatorname{atan}\left(-\frac{a_{21}}{a_{22}}\right),
        sy = \frac{a_{22}}{\cos \alpha}
        sx = \frac{\operatorname{det}(M)}{sy}},
        hx = \frac{a_{11} \cdot a_{21} + a_{12} \cdot a_{22}}{\operatorname{det}(M)},
        tx = a_{31} + ty_0 \cdot a_{21} + tx_0 \cdot (a_{11} - 1),
        ty = a_{32} + tx_0 \cdot a{12} + ty_0 \cdot (a_{22} - 1).

        TeX
        $$ M = \begin{pmatrix}\ a_{11} & a_{12} & 0 \\ a_{21} & a_{22} & 0 \\ a_{31} & a_{32} & 1\end{pmatrix}$$
        
        $$ \operatorname{det} M = a_{11}\,a_{22} - a_{12}\,a_{21},\qquad \operatorname{det} M \ne 0$$
        
        $$ p_1 = p_0 \cdot M = (x_0, y_0, 1) \begin{pmatrix}\ a_{11} & a_{12} & 0 \\ a_{21} & a_{22} & 0 \\ a_{31} & a_{32} & 1\end{pmatrix} $$
        
        $$ M_c = R \cdot H_x \cdot S $$
        
        $$ R = \begin{pmatrix}\ \cos \alpha & \sin \alpha & 0 \\ -\sin \alpha & \cos \alpha & 0 \\ 0 & 0 & 1\end{pmatrix}$$
        
        $$ H_x = \begin{pmatrix}\ 1 & hx & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1\end{pmatrix}$$
        
        $$ S = \begin{pmatrix}\ sx & 0 & 0 \\ 0 & sy & 0 \\ 0 & 0 & 1\end{pmatrix}$$
        
        $$ M = T_0 \cdot S \cdot H \cdot R \cdot T_0^{-1} \cdot T $$
        
        $$ T_0 = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ -tx_0 & -ty_0 & 1\end{pmatrix}$$
        
        $$ T_0^{-1} = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ tx_0 & ty_0 & 1\end{pmatrix}$$
        
        $$ T = \begin{pmatrix}\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ tx & ty & 1\end{pmatrix}$$
        
        $$a_{11} = sx \cdot \cos \alpha - sx \cdot hx \cdot \sin \alpha$$
        
        $$a_{12} = sx \cdot hx \cdot \cos \alpha + sx \cdot \sin \alpha$$
        
        $$a_{21} = -sy \cdot \sin \alpha$$
        
        $$a_{22} = sy \cdot \cos \alpha$$
        
        $$a_{31} = tx + tx_0 \cdot (1 - (sx \cdot \cos \alpha - sx \cdot hx \cdot \sin \alpha)) + ty_0 \cdot sy \cdot \sin \alpha $$
        
        $$= tx + tx_0 \cdot (1 - a_{11}) - ty_0 \cdot a_{21}$$
        
        $$a_{32} = ty + ty_0 \cdot (1 - \cos \alpha \cdot sy) - tx_0 \cdot (sx \cdot hx \cdot \cos \alpha + sx \cdot \sin \alpha) $$
        
        $$ = ty + ty_0 \cdot (1 - a_{22}) - tx_0 \cdot a_{12}$$
        
        $$\operatorname{if} a_{22}=0 $$
        
        $$ \alpha = \frac{\pi}{2},$$
        
        $$ sy = -a_{21} $$
        
        $$ \operatorname{else} $$
        
        $$ \alpha = \operatorname{atan}\left(-\frac{a_{21}}{a_{22}}\right), $$
        
        $$ sy = \frac{a_{22}}{\cos \alpha} $$
        
        $$ sx = \frac{\operatorname{det}(M)}{sy}},$$
        
        $$ hx = \frac{a_{11} \cdot a_{21} + a_{12} \cdot a_{22}}{\operatorname{det}(M)},$$
        
        $$ tx = a_{31} + ty_0 \cdot a_{21} + tx_0 \cdot (a_{11} - 1),$$
        
        $$ ty = a_{32} + tx_0 \cdot a{12} + ty_0 \cdot (a_{22} - 1).$$


        Сконвертить TeX в картинки можно тут.
        • 0
          Если рассматривать это как систему уравнений в которой мы умножаем перпеменные координаты точки на коэффициеты матриы преобразования, то, наверно, вариант с вектором-столбцом действительно смотрится логичнее. По поводу транспонирования это сюда: Quartz 2D Programming Guide: The Math Behind the Matrices. Поскольку это все делалось под Мак. Лично я ничего криминального в умножении справа не вижу.

          По поводу формул, спасибо вам за картинки
          • 0
            Кроме того, в статье на которую я ссылаюсь, судя по виду матрицы поворота, также умножение справа.
            По поводу формул, спасибо вам за картинки. Пусть они будут здесь. Но, мне все-таки мой вариант нравится больше, за исключением того, что из-за увеличенного интервала между строк разбился "ASCII art" для скобок вокруг матриц. Интересно, можно ли этот интервал локально убрать? А нравится мне мой вариант больше, потому что визуально большой разницы нет (по-крайней мере не для математика), а поскольку это не картинки а текст, его можно просто скопировать и куда-то вставить, например себе в комментарии в код.
      • НЛО прилетело и опубликовало эту надпись здесь
        • +1
          Тут прекрасно всё — и преобразования имени славного греческого города Афины, и программирование алгоритма.
          • 0
            Про Афины – спасибо, поправил. А вот про "программирование алгоритма" не понял.
            • +1
              Программирования в материале не видно: как самой программы, так и результатов работы ее алгоритма. Написание нескольких строк кода по формуле — не выглядит серьезным достижением.
              • 0
                А, вы имеет в виду, что я поместил в хаб Программирование и не привел примера ее реализации? Понятно. Ну я в этот хаб поместил, потому, что эта задача имеет очевидное отношение к программированию, а не потому что, я много и упорно программировал, пытаясь ее решить. Очевидно, что реализация этого алгоритма элементарна. По этой причине я и не стал приводить ее.
                На серьезное достижение я не претендую, что вы. Достижение не у меня, а в статье, на которую я ссылаюсь. Я всего лишь вытащил из нее то, что непосредственно требуется для решения поставленной задачи, добавил отсутствующие компоненты и опубликовал здесь в надежде на то, что это может пригодиться еще кому-нибудь. Потому что, это теперь статья Путятина у меня в Гугле выскакивает на первой странице, а когда я искал соответствующую информацию то ли год, то ли полгода назад, как-то ничего не попадалось. То ли я запрос по другому составил, то ли в Гугле что-то подкрутили, или может этой статьи там тогда и не было добавлено. Не знаю.
        • 0
          Вот тут за матрицы вам расскажут таки всё: http://www.macaronikazoo.com/?page_id=413

          Втч и про декомпозицию.
          • 0
            За ссылку спасибо. Почитал. Есть интересные моменты, но про то, что там "за матрицы таки все", это вы погорячились. :)

            Technically its impossible to extract scale and rotation from a general transformation matrix, but in reality its generally easy. What does this mean? Well basically a general transformation matrix can have shear as well as scale as well as rotation. And all three of these pieces of data live in the upper 3×3 quadrant of the matrix. But if you know you’re just working with rotations, translations and scale (which in my professional life has always been the case – but perhaps others have different stories?) then its possible and in fact pretty easy.

            Человек пишет, фактически, что в общем случае задача декомпозиции неразрешима и рассматривает один отдельный случай, в котором отсутствуют преобразования сдвига и порядок применения матриц известен.
            • 0
              Просто после прочтения этих 6 статей, мозги встают на место и начинаешь понимать что же такое матрица. По крайней мере я после этого смог довольно быстро сбацать свой собственный деформер для Autodesk Maya как раз на основе матриц.

              Всё же поражаюсь иногда, насколько неэффективно и непонятно могут обьяснять простые вещи в учебных заведениях.
            • 0
              (Продолжение, опять не на ту клавишу нажал.)
              По поводу неразрешимости задачи декомпозиции в общем случае для 3D сказать сразу ничего не могу. Но как мы видим для случая 2D он не прав.
              В статьях его еще меня смутило то, что он вроде бы использует умножение справа, как и у меня, но в случае иерархии объектов сначала применяет матрицу трансформации родительского объекта (глобальную), а не локальную. Это странно. Например, у меня есть квадрат, растянутый в прямоугольник преобразованием масштабирования по оси Y. Он является частью вращающегося объекта. Кажется очевидным, что для того чтобы получить вращающийся прямоугольник, который сохраняет свою геометрию, надо сначала растянуть квадрат и только потом вращать, а не наоборот.
              • 0
                Это может быть связано с конкретной майской имплементацией, автор — майский TD.
                Там есть разные порядки применения трансформов — все эти TRS, SRT и проч и проч.
                • 0
                  Да, возможно это "особенности" майи. То есть, они могут формировать локальную и глобальную матрицы (независимо от того какой у них там порядок: TSR, SRT) как для умножения справа, но результирующую матрицу формировать для умножения слева, транспонируя глобальную и локальную матрицы.

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