src/map/transform.cpp

   1 #include "projection.h"
   2 #include "matrix.h"
   3 #include "transform.h"
   4
   5
   6 #define NULL_QTRANSFORM 0,0,0,0,0,0,0,0,0
   7
   8 void Transform::simple(const ReferencePoint &p1, const ReferencePoint &p2)
   9 {
  10         if (p1.xy().x() == p2.xy().x() || p1.xy().y() == p2.xy().y()) {
  11                 _errorString = "Invalid reference points tuple";
  12                 return;
  13         }
  14
  15         double sX = (p1.xy().x() - p2.xy().x()) / (p1.pp().x() - p2.pp().x());
  16         double sY = (p2.xy().y() - p1.xy().y()) / (p2.pp().y() - p1.pp().y());
  17         double dX = p2.xy().x() - p2.pp().x() * sX;
  18         double dY = p1.xy().y() - p1.pp().y() * sY;
  19
  20         _proj2img = QTransform(sX, 0, 0, sY, dX, dY);
  21         _img2proj = _proj2img.inverted();
  22 }
  23
  24 void Transform::affine(const QList<ReferencePoint> &points)
  25 {
  26         Matrix c(3, 2);
  27         c.zeroize();
  28         for (size_t i = 0; i < c.h(); i++) {
  29                 for (size_t j = 0; j < c.w(); j++) {
  30                         for (int k = 0; k < points.size(); k++) {
  31                                 double f[3], t[2];
  32
  33                                 f[0] = points.at(k).pp().x();
  34                                 f[1] = points.at(k).pp().y();
  35                                 f[2] = 1.0;
  36                                 t[0] = points.at(k).xy().x();
  37                                 t[1] = points.at(k).xy().y();
  38                                 c.m(i,j) += f[i] * t[j];
  39                         }
  40                 }
  41         }
  42
  43         Matrix Q(3, 3);
  44         Q.zeroize();
  45         for (int qi = 0; qi < points.size(); qi++) {
  46                 double v[3];
  47
  48                 v[0] = points.at(qi).pp().x();
  49                 v[1] = points.at(qi).pp().y();
  50                 v[2] = 1.0;
  51                 for (size_t i = 0; i < Q.h(); i++)
  52                         for (size_t j = 0; j < Q.w(); j++)
  53                                 Q.m(i,j) += v[i] * v[j];
  54         }
  55
  56         Matrix M = Q.augemented(c);
  57         if (!M.eliminate()) {
  58                 _errorString = "Singular transformation matrix";
  59                 return;
  60         }
  61
  62         _proj2img = QTransform(M.m(0,3), M.m(0,4), M.m(1,3), M.m(1,4), M.m(2,3),
  63           M.m(2,4));
  64         _img2proj = _proj2img.inverted();
  65 }
  66
  67
  68 Transform::Transform()
  69   : _proj2img(NULL_QTRANSFORM), _img2proj(NULL_QTRANSFORM)
  70 {
  71 }
  72
  73 Transform::Transform(const QList<ReferencePoint> &points)
  74   : _proj2img(NULL_QTRANSFORM), _img2proj(NULL_QTRANSFORM)
  75 {
  76         if (points.count() < 2)
  77                 _errorString = "Insufficient number of reference points";
  78         else if (points.size() == 2)
  79                 simple(points.at(0), points.at(1));
  80         else
  81                 affine(points);
  82 }
  83
  84 Transform::Transform(const ReferencePoint &p1, const ReferencePoint &p2)
  85   : _proj2img(NULL_QTRANSFORM), _img2proj(NULL_QTRANSFORM)
  86 {
  87         simple(p1, p2);
  88 }
  89
  90 Transform::Transform(const ReferencePoint &p, const PointD &scale)
  91   : _proj2img(NULL_QTRANSFORM), _img2proj(NULL_QTRANSFORM)
  92 {
  93         if (scale.x() == 0.0 || scale.y() == 0.0) {
  94                 _errorString = "Invalid scale factor";
  95                 return;
  96         }
  97
  98         _img2proj = QTransform(scale.x(), 0, 0, -scale.y(), p.pp().x() - p.xy().x()
  99           / scale.x(), p.pp().y() + p.xy().x() / scale.y());
 100         _proj2img = _img2proj.inverted();
 101 }
 102
 103 Transform::Transform(double matrix[16])
 104   : _proj2img(NULL_QTRANSFORM), _img2proj(NULL_QTRANSFORM)
 105 {
 106         _img2proj = QTransform(matrix[0], matrix[1], matrix[4], matrix[5],
 107           matrix[3], matrix[7]);
 108         if (!_img2proj.isInvertible())
 109                 _errorString = "Singular transformation matrix";
 110         else
 111                 _proj2img = _img2proj.inverted();
 112 }
 113
 114 #ifndef QT_NO_DEBUG
 115 QDebug operator<<(QDebug dbg, const ReferencePoint &p)
 116 {
 117         dbg.nospace() << "ReferencePoint(" << p.xy() << ", " << p.pp() << ")";
 118         return dbg.space();
 119 }
 120 #endif // QT_NO_DEBUG