Обсуждение

Чтобы задать аффинное преобразование, достаточно выбрать в качестве опорных точек три любые точки, не лежащие на одной прямой, и указать старые и новые координаты этих точек. Удобно принять в качестве опорных точек начало координат и две точки, лежащие соответственно на осях X и Y.

- Аффинное преобразование с вычислением коэффициентов методом наименьших квадратов предназначено для случая, когда число заданных опорных точек больше трех. Это означает, что решается задача уточнения по ориентирам. Если число ориентиров превышает 3 и некоторые из значений их координат (старых или новых) получены с некоторой погрешностью, то нельзя рассчитывать на построение аффинного преобразования, точно переводящего все ориентиры в новое положение. Можно обеспечить точное преобразование любых трех точек, но чем хуже остальные опорные точки? Чтобы не прибегать к произвольному выбору трех точек для построения преобразования, можно использовать известный в математике метод наименьших квадратов. При построении преобразования по этому методу может оказаться, что ни для одной из опорных точек расчетное положение не совпадет с измеренными новыми координатами, однако средняя ошибка (точнее, средний квадрат расстояния между измеренным и рассчитанным положением каждой опорной точки) будет минимальной среди всевозможных аффинных преобразований. Таким образом, данный тип преобразования не может обеспечить точность преобразования опорных точек, но зато сохраняет все прямые линии прямыми."

Теоретические выкладки для аффинного преобразования с подбором коэффициентов методом наименьших квадратов здесь нет возможности привести. Приведем только фрагмент кода на C расчета коэффициентов, который используется в трансформаторе ГИС ObjectLand:

/******************** AFFI_LSCOEFS ********************/
U32 WINAPI AFFI_LSCOEFS ( ORD_COLLECTION_OOP pSrcPointList,
                          ORD_COLLECTION_OOP pDestPointList,
                                   FLOAT_OOP *coefs)
/*----------------------------------------------------------------
 * Вычисляет коэффициенты аффинного преобразования, переводящего
 * набор точек srcPoints в точки, наиболее близкие к заданным
 * destPoints по критерию наименьшей суммы квадратов отклонений.
 * Всегда возвращает 0.
 *
 * ЗАМЕЧАНИЕ. Используется простейший подход: построение системы нормальных
 * уравнений и ее решение по компактной схеме метода Гаусса (Корн&Корн, с.662).
 * Преобразование задается формулами:
 *   xd = C00 + C01 * xs + C02 * ys
 *   yd = C10 + C11 * xs + C12 * ys
 * Коэффициенты помещаются в массив coef по строкам: 
 *   C00,C01,C02,C10,C11,C12.
 * ЗАМЕЧАНИЕ. В случае вырожденных данных (когда исходные точки лежат на
 * одной прямой) возвращает неверные результаты. Это не проверяется.
 */
{
    INT32  N;             /* Количество точек */
    INT32  k, m;          /* Индексы циклов */
    double a11, a12, a13, a21, a22, a23, a31, a32, a33;    /* Матрица 3*3 */
    double b1, b2, b3;    /* Свободный член */
    double xx, yy;        /* Координаты текущей точки */
    double zz;            /* Либо x, либо y */
    double x2;            /* Рабочая величина при вычислении коэффициентов */
    INT32  srcStart, destStart;
    ARRAY_OOP srcContents, destContents;

    srcStart = fromSmallInteger(pSrcPointList->startPosition) - 1;
    N = fromSmallInteger(pSrcPointList->endPosition) - srcStart;
    srcContents = (pSrcPointList->contents) + srcStart;
    destStart = fromSmallInteger(pDestPointList->startPosition) - 1;
    destContents = (pDestPointList->contents) + destStart;

    /* Нормальная система состоит из 2-х подсистем по 3 уравнения, отличающихся
       свободными членами. */
    /* Подсчет общих коэффициентов системы */
    a11 = N;
    a12 = a13 = a22 = a23 = a33 = 0.0L;
    for (k = 0; k < N; k++)    {
        xx = *(((FL_POINT_OOP)(srcContents[k]))->x);
        yy = *(((FL_POINT_OOP)(srcContents[k]))->y);
        a12 += xx;
        a13 += yy;
        a22 += xx * xx;
        a23 += xx * yy;
        a33 += yy * yy;
    }    /*  for k */
    /*  Преобразование коэффициентов (треугольное разложение матрицы) */
    a21 = a12;
    a31 = a13;
    a12 /= a11;
    a13 /= a11;
    a22 -= a21 * a12;
    a32 = a23 - a31 * a12;
    a23 = a32 / a22;
    a33 -= a31 * a13 + a32 * a23;
    /*  Теперь вещи, различные для подсистем X и Y */
    for (m = 0; m < 2; m++) {
        /*  Подсчет свободных членов подсистемы */
        b1 = b2 = b3 = 0.0L;
        for (k = 0; k < N; k++)    {
            zz = (!m 
                ? *(((FL_POINT_OOP)(destContents[k]))->x)
                : *(((FL_POINT_OOP)(destContents[k]))->y));
            b1  += zz;
            b2  += *(((FL_POINT_OOP)(srcContents[k]))->x) * zz;
            b3  += *(((FL_POINT_OOP)(srcContents[k]))->y) * zz;
        }    /*  for k */
        /*  Преобразование свободных членов */
        b1 /= a11;
        b2 = (b2 - a21 * b1) / a22;
        b3 = (b3 - a31 * b1 - a32 * b2) / a33;
        /*  Решение подсистемы */
        x2 = b2 - a23 * b3;
        **coefs++ = b1 - a12 * x2 - a13 * b3;
        **coefs++ = x2;
        **coefs++ = b3;
    }    /*  for m */
    return 0L;
}/* AFFI_LSCOEFS */