[Proj] Grown error if re-projecting from 4269 to LCC (2285) and backward multiple times
Thomas Knudsen
knudsen.thomas at gmail.com
Tue Oct 31 07:59:02 EST 2017
Thanks for an admirably complete question, with sufficient information
to meaningfully attempt to answer it.
Unfortunately, my answer will probably disappoint you.
After 7 iterations your deviation from the origin is less than a micrometre
in X and less than one tenth of a millimetre in Y.
Since PROJ.4 uses numerical, not analytical, mathematics, and since the
resolution is limited by the mantissa size of IEEE-754 64 bit doubles, the
roundtrips will never (for a rather all-encompassing value of never) be
exact.
Typically, they will converge after a few hundred iterations, after which
the
deviation will cease to grow overall, but they will still not be exact.
To me, it is slightly surprising that you get exact roundtrips for mercator,
although the ellipsoidal mercator projection is somewhat simpler than
ellipsoidal LCC, so it will probably be more probable.
So my overall message is that you should not expect exact roundtrips.
Another thing is that it surprises me you need exact roundtrips after
an operation that involves direct user interaction. Which human
interaction is it that touches up data at the micrometer level, and needs
exact replication? If you could tell a bit about this it might be easier to
see potential workarounds.
2017-10-31 7:50 GMT+01:00 Dmitry Mefed <dmitrymefed at gmail.com>:
> Hi,
>
>
> Our software stores the data with EPSG 4269 and present to the User in a
> local projected CRS, EPSG 2285 in examples described below. In case if an
> entity is modified on the UI we need to push the changes back. In that
> scenario we have to do
>
> 1. 4269 -> 2285 projection,
>
> 2. then User modifies the entity (in 2285 CRS),
>
> 3. In order to store the changes, we do 2285 -> 4269 projection
>
> 4. and next time the entity is requested we do 4269 -> 2285
> projection again.
>
> The problem is that proj4 does not produce same results in step 2 and 4.
> It would be expected behavior as 2285 provides more precision due to used
> units US-ft, in comparison to 4269 which units are degrees. But the real
> issue for us is that subsequent transformations between 2285 -> 4269 ->
> 2285… produces growing error. Please consider the following output of cs2cs
> program, (in parallel I did inverse projection to obtain coordinates in
> EPSG 4269):
>
> C:\PROJ\bin>cs2cs.exe -v +init=epsg:2285 +to +init=epsg:4269 -f %.16f
>
> pj_open_lib(epsg): call fopen(C:\proj\share\epsg) - succeeded
>
>
>
> pj_open_lib(epsg): call fopen(C:\proj\share\epsg) - succeeded
>
>
>
> # ---- From Coordinate System ----
>
> #Lambert Conformal Conic
>
> # Conic, Sph&Ell
>
> # lat_1= and lat_2= or lat_0
>
> # +init=epsg:2285 +proj=lcc +lat_1=48.73333333333333 +lat_2=47.5 +lat_0=47
>
> # +lon_0=-120.8333333333333 +x_0=500000.0001016001 +y_0=0 +datum=NAD83
>
> # +units=us-ft +no_defs +ellps=GRS80 +towgs84=0,0,0
>
> # ---- To Coordinate System ----
>
> #Lat/long (Geodetic alias)
>
> #
>
> # +init=epsg:4269 +proj=longlat +datum=NAD83 +no_defs +ellps=GRS80
>
> # +towgs84=0,0,0
>
> 2196293.3184066000 <(318)%20406-6000> 643350.39301622100
>
> -118.5293900000000300 48.7408309999860520 0.0000000000000000
>
> 2196293.3184067495000000 643350.3930111383100000
>
> -118.5293900000000300 48.7408309999721180 0.0000000000000000
>
> 2196293.3184069018000000 643350.3930060574800000
>
> -118.5293900000000300 48.7408309999581850 0.0000000000000000
>
> 2196293.3184070541000000 643350.3930009742000000
>
> -118.5293900000000300 48.7408309999442440 0.0000000000000000
>
> 2196293.3184072063000000 643350.3929958911600000
>
> -118.5293900000000300 48.7408309999302960 0.0000000000000000
>
> 2196293.3184073586000000 643350.3929908055600000
>
> -118.5293900000000300 48.7408309999163550 0.0000000000000000
>
> 2196293.3184075109000000 643350.3929857177400000
>
> -118.5293900000000300 48.7408309999023930 0.0000000000000000
>
> 2196293.3184076636000000 643350.3929806276000000
>
> -118.5293900000000300 48.7408309998884450 0.0000000000000000
>
>
>
> *Simple Excel based analysis shows the following deltas on each step:*
>
> X
>
> dX
>
> Y
>
> dY
>
> 0,31840660000
>
>
>
> 0,39301622100
>
>
>
> 0,31840674950
>
> 0,00000014950
>
> 0,39301113831
>
> -0,00000508269
>
> 0,31840690180
>
> 0,00000015230
>
> 0,39300605748
>
> -0,00000508083
>
> 0,31840705410
>
> 0,00000015230
>
> 0,39300097420
>
> -0,00000508328
>
> 0,31840720630
>
> 0,00000015220
>
> 0,39299589116
>
> -0,00000508304
>
> 0,31840735860
>
> 0,00000015230
>
> 0,39299080556
>
> -0,00000508560
>
> 0,31840751090
>
> 0,00000015230
>
> 0,39298571774
>
> -0,00000508782
>
> 0,31840766360
>
> 0,00000015270
>
> 0,39298062760
>
> -0,00000509014
>
>
>
> I also tried other LCC CRSs the picture is similar. And I also tried
> Mercator projection EPSG 3857 it does not seem to have such issue, after
> full cycle of projections, it produces constant results in both directions.
>
> Please let me know if it is known issue or it is WAD (Works As Designed)
> or if I can do something with this.
>
>
> PS
>
> this is the second window output:
>
> C:\PROJ\bin>cs2cs.exe -v +init=epsg:4269 +to +init=epsg:2285 -f %.16f
> # ---- From Coordinate System ----
> #Lat/long (Geodetic alias)
> #
> # +init=epsg:4269 +proj=longlat +datum=NAD83 +no_defs +ellps=GRS80
> # +towgs84=0,0,0
> # ---- To Coordinate System ----
> #Lambert Conformal Conic
> # Conic, Sph&Ell
> # lat_1= and lat_2= or lat_0
> # +init=epsg:2285 +proj=lcc +lat_1=48.73333333333333 +lat_2=47.5 +lat_0=47
> # +lon_0=-120.8333333333333 +x_0=500000.0001016001 +y_0=0 +datum=NAD83
> # +units=us-ft +no_defs +ellps=GRS80 +towgs84=0,0,0
> -118.5293900000000300 48.7408309999860520
> 2196293.3184067495000000 643350.3930111383100000 0.0000000000000000
> -118.5293900000000300 48.7408309999721180
> 2196293.3184069018000000 643350.3930060574800000 0.0000000000000000
> -118.5293900000000300 48.7408309999581850
> 2196293.3184070541000000 643350.3930009742000000 0.0000000000000000
> -118.5293900000000300 48.7408309999442440
> 2196293.3184072063000000 643350.3929958911600000 0.0000000000000000
> -118.5293900000000300 48.7408309999302960
> 2196293.3184073586000000 643350.3929908055600000 0.0000000000000000
> -118.5293900000000300 48.7408309999163550
> 2196293.3184075109000000 643350.3929857177400000 0.0000000000000000
> -118.5293900000000300 48.7408309999023930
> 2196293.3184076636000000 643350.3929806276000000 0.0000000000000000
>
>
> Thanks,
>
> -Dmitry Mefed
>
> _______________________________________________
> Proj mailing list
> Proj at lists.maptools.org
> http://lists.maptools.org/mailman/listinfo/proj
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.maptools.org/pipermail/proj/attachments/20171031/da30f417/attachment-0001.htm
More information about the Proj
mailing list