[Merkaartor] Inaccuracy in Projection

[Chris Browet]

I already did some rounding in my current code. I'm stuck with a problem,
that's why I didn' commit.

Thanks anyway.

- Chris -

2009/4/14 Timo Schlüßler <timo at schluessler.org>

> Hi,
>
> I noticed a small (or little?) inaccuracy when e.g. creating a node. The
> node doesn't appear there where I clicked but some pixels away from that
> point. When zooming in to a high zoom value (1meter scale) this effect is
> intensified.
>
> The appended patch changes some (int) - casts to qRound() s but the problem
> still lasts but much lesser.
> Could it be that the lasting inaccuracy results through converting every
> coordinate to an integer? The inserted debug outputs show that the
> projection (the boost projection library) itself doesn't contain rounding
> mistakes.
>
> Any ideas how to solve this completely? Or has nobody else noticed this
> variances?
>
> regards,
> timo
>
> Index: src/Interaction/CreateNodeInteraction.cpp
> ===================================================================
> --- src/Interaction/CreateNodeInteraction.cpp   (revision 14474)
> +++ src/Interaction/CreateNodeInteraction.cpp   (working copy)
> @@ -24,6 +24,7 @@
>        if (ev->button() == Qt::LeftButton)
>        {
>                Coord P(projection().inverse(ev->pos()));
> +               qDebug() << "inverting " << ev->pos() << " to " <<
> projection().project(P);
>                if (aRoad)
>                {
>                        main()->properties()->setSelection(0);
> Index: src/Map/Coord.h
> ===================================================================
> --- src/Map/Coord.h     (revision 14479)
> +++ src/Map/Coord.h     (working copy)
> @@ -11,7 +11,7 @@
>
>  inline double angToRad(double a)
>  {
> -       return a*M_PI/180;
> +       return a*M_PI/180.0;
>  }
>  //
>  //inline double radToAng(double a)
> @@ -21,12 +21,13 @@
>
>  inline int angToInt(double a)
>  {
> -       return static_cast<int>(a*INT_MAX/180.0);
> +       //return static_cast<int>(a*INT_MAX/180.0);
> +       return qRound(a/180.0*INT_MAX);
>  }
>
>  inline double intToAng(int a)
>  {
> -       return double(a)*180/INT_MAX;
> +       return double(a)*180.0/INT_MAX;
>  }
>
>  inline double intToRad(int a)
> @@ -36,7 +37,7 @@
>
>  inline int radToInt(double x)
>  {
> -       return static_cast<int>(x*INT_MAX/M_PI);
> +       return qRound(x/M_PI*INT_MAX);
>  }
>
>  class Coord
> Index: src/Map/Projection.cpp
> ===================================================================
> --- src/Map/Projection.cpp      (revision 14479)
> +++ src/Map/Projection.cpp      (working copy)
> @@ -19,6 +19,9 @@
>
>  #include <math.h>
>
> +#include <QDebug>
> +#include <iomanip>
> +
>  // from wikipedia
>  #define EQUATORIALRADIUS 6378137.0
>  #define POLARRADIUS      6356752.0
> @@ -66,13 +69,13 @@
>
>  double Projection::latAnglePerM() const
>  {
> -       double LengthOfOneDegreeLat = EQUATORIALRADIUS * M_PI / 180;
> +       double LengthOfOneDegreeLat = EQUATORIALRADIUS * M_PI / 180.0;
>        return 1 / LengthOfOneDegreeLat;
>  }
>
>  double Projection::lonAnglePerM(double Lat) const
>  {
> -       double LengthOfOneDegreeLat = EQUATORIALRADIUS * M_PI / 180;
> +       double LengthOfOneDegreeLat = EQUATORIALRADIUS * M_PI / 180.0;
>        double LengthOfOneDegreeLon = LengthOfOneDegreeLat * fabs(cos(Lat));
>        return 1 / LengthOfOneDegreeLon;
>  }
> @@ -106,6 +109,7 @@
>        point_2d out;
>
>        theProj->forward(in, out);
> +       std::cout << "noverse; in: " << in << ", out: " << out <<
> std::endl;
>
>        return QPoint(out.x(), out.y());
>  }
> @@ -116,6 +120,7 @@
>        point_ll_rad out;
>
>        theProj->inverse(in, out);
> +       std::cout << "inverse; in: " << in << ", out: " << out <<
> std::endl;
>
>        return Coord(radToInt(out.lat()), radToInt(out.lon()));
>  }
> @@ -123,26 +128,26 @@
>  QPoint Projection::project(const Coord & Map) const
>  {
>        QPoint p = projProject(Map);
> -       return QPoint(int(p.x() * ScaleLon + DeltaLon), int(-p.y() *
> ScaleLat + DeltaLat));
> +       return QPoint(qRound(p.x() * ScaleLon + (double)DeltaLon),
> qRound(-p.y() * ScaleLat + (double)DeltaLat));
>  }
>
>  QPoint Projection::project(TrackPoint* aNode) const
>  {
>        if (aNode && aNode->projectionType() == theProjectionType &&
> !aNode->projection().isNull())
> -               return QPoint(int(aNode->projection().x() * ScaleLon +
> DeltaLon),
> -                                          int(-aNode->projection().y() *
> ScaleLat + DeltaLat));
> +               return QPoint(qRound(aNode->projection().x() * ScaleLon +
> (double)DeltaLon),
> +                                          qRound(-aNode->projection().y()
> * ScaleLat + (double)DeltaLat));
>
>        QPoint p = projProject(aNode->position());
>
>        aNode->setProjectionType(theProjectionType);
>        aNode->setProjection(p);
>
> -       return QPoint(int(p.x() * ScaleLon + DeltaLon), int(-p.y() *
> ScaleLat + DeltaLat));
> +       return QPoint(qRound(p.x() * ScaleLon + (double)DeltaLon),
> qRound(-p.y() * ScaleLat + (double)DeltaLat));
>  }
>
>  Coord Projection::inverse(const QPointF & Screen) const
>  {
> -       Coord c = projInverse(QPoint((Screen.x() - DeltaLon ) / ScaleLon,
> -(Screen.y() - DeltaLat) / ScaleLat));
> +       Coord c = projInverse(QPoint(qRound((Screen.x() - (double)DeltaLon
> ) / ScaleLon), qRound(-(Screen.y() - (double)DeltaLat) / ScaleLat)));
>        return c;
>  }
>
>
> _______________________________________________
> Merkaartor mailing list
> Merkaartor at openstreetmap.org
> http://lists.openstreetmap.org/listinfo/merkaartor
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.openstreetmap.org/pipermail/merkaartor/attachments/20090414/1e999114/attachment.html>