Opened 10 years ago

Closed 10 years ago

#5277 closed defect (invalid)

Something is wrong with ESGP:31370

Reported by: caller82 Owned by: warmerdam
Priority: normal Milestone:
Component: default Version: 1.10.1
Severity: normal Keywords:
Cc:

Description

hi,

I am using your java bingings. If I try to convert coordinates between WGS84 latlon (4326) and EPSG::31370 I got huge errors (50-100m) if I select the SpacialReference by the EPSG code. If I instead use the Proj4 definition string which is stored in your "epsg" file, everything is fine. Do you may have a bug in the "pcs.csv" file for this code?

   SpatialReference src = new SpatialReference("");
	        src.ImportFromEPSGA(4326);
	        SpatialReference dst = new SpatialReference("");
//	        dst.ImportFromEPSGA(31370);//Belgian projection, using this gives awful results!
	        dst.ImportFromProj4("<31370> +proj=lcc +lat_1=51.16666723333333 +lat_2=49.8333339 +lat_0=90 +lon_0=4.367486666666666 +x_0=150000.013 +y_0=5400088.438 +ellps=intl +towgs84=106.869,-52.2978,103.724,-0.33657,0.456955,-1.84218,1 +units=m +no_defs <>");

CoordinateTransformation ct = new CoordinateTransformation(src, dst);
	        double[] p = new double[3];
	        p[0] = 3.5;
	        p[1] = 50.5;
	        p[2] = 100;
	        System.out.println("x-org  :" + p[0] + " y:" + p[1] + " z:" + p[2]);
	        ct.TransformPoint(p);
	        System.out.println("x-trans:" + p[0] + " y:" + p[1] + " z:" + p[2]);

thx

-marco

Change History (4)

comment:1 by Jukka Rahkonen, 10 years ago

Hi,

I checked from a couple of Spatialite versions the content of spatial_ref_sys as well as from a few EPSG files from Mapserver MS4W installations. All of those have:

+towgs84=-106.869,52.2978,-103.724,0.3366,-0.457,1.8422,-1.2747

Are you sure about your own parameters with reversed signs? If they are correct then Proj.4 has been wrong for a long time.

-Jukka Rahkonen-

comment:2 by caller82, 10 years ago

Hi,

maybe you got me wrong. I have extracted the proj4 values from the "epsg" file in the current GDAL source tree. If you copy this values, and initialize the projection using this string by: dst.ImportFromProj4("... everything is fine!!! But if I using this initialisation dst.ImportFromEPSGA(31370); I got wrong results, even they should be equal. the "ImportFromEPSGA" takes its values from "pcs.csv" as far as I know... so I guess something is wrong with this file...

THX

-marco

comment:3 by Even Rouault, 10 years ago

Frank,

I'm deeply confused by this issue... Isn't there an issue with the sign of towgs84 parameters depending on the method 9606 vs 9607 ? My understanding of http://trac.osgeo.org/proj/wiki/GenParms ( "Note that EPSG method 9607 (coordinate frame rotation) coefficients can be converted to EPSG method 9606 (position vector 7-parameter) supported by PROJ.4 by reversing the sign of the rotation vectors. The methods are otherwise the same. ") would be that in the Proj.4 EPSG file the parameters should always be normalized to 9606, so for GCS with method 9607, there should be a sign inversion between PROJ epsg file and GDAL gcs.csv ? Right ? This is not what I can see :

  • gcs.csv (since GDAL 1.8, no TOWGS84 before):

4313,Belge 1972,6313,Reseau National Belge 1972,6313,9122,7022,8901,1,0,6422,15929,1,9607,-106.8686,52.2978,-103.7239,-0.3366,0.457,-1.8422,-1.2747

  • epsg (PROJ4.8 and trunk)

# Belge 1972 <4313> +proj=longlat +ellps=intl +towgs84=-106.869,52.2978,-103.724,0.3366,-0.457,1.8422,-1.2747 +no_defs <>

  • epsg (PROJ4.7, no TOWGS84 before ):

# Belge 1972 <4313> +proj=longlat +ellps=intl +towgs84=106.869,-52.2978,103.724,-0.33657,0.456955,-1.84218,1 +no_defs <>

I can see that in ogr_fromepsg.cpp, we reverse the sign when reading from gcs.csv and the method is 9607 (which is consistant with the above Proj wiki page). This logic was added 11 years ago ( http://trac.osgeo.org/gdal/changeset/3959/trunk/ogr/ogr_fromepsg.cpp ), so I don't understand why PROJ epsg file has changed between PROJ 4.7 and PROJ 4.8 wheras GDAL gcs.csv hasn't changed since it has TOWGS84 parameters for that SRS...

comment:4 by Even Rouault, 10 years ago

Resolution: invalid
Status: newclosed

Hum, ok actually the story is a bit more complicated than that and I've found the explanation for the sign changes. There has been a related past GDAL ticket about that ( http://trac.osgeo.org/gdal/ticket/3362 ) that lead to a gcs.override.csv with

4313,Belge 1972,6313,Reseau National Belge 1972,6313,9122,7022,8901,1,0,6422,9607,-106.8686,52.2978,-103.7239,-0.3366,0.457,-1.8422,-1.2747

The override has disappeared in a later update of the EPSG database, but the normal gcs.csv seems to contain that line now.

caller82, my belief is that the old values in the epsg file you are using (proj 4.7) were wrong and that has been corrected by #3362.

So I'm closing the ticket. Please reopen if you have evidences (based on reference points for example) that the new result is incorrect, but the analysis done in #3362 seems to be quite deeps. And upgrade your epsg PROJ file to proj 4.8 to have consistant results between proj and GDAL.

Note: See TracTickets for help on using tickets.