Re: [Geotools-devel] Axis problems again with CRS.lookupEpsgCode

SourceForge Headquarters 225 Broadway Suite 1600 San Diego, CA 92101 +1 (858) 454-5900

Hey Andrea

Thanks a lot for the feedback. I wish to apologize for not sending a new 
e-mail yesterday with what I was finding out.

Turns out that, yeah, when I came back the next morning, this exact WKT 
was working (but not the real one I'm getting in my system here). I must 
have messed something up when I tested it before sending the e-mail, 
really sorry for that, guess I was too tired.

But don't be completely mad at me, there are real issues going on. I was 
waiting to send an e-mail with those issues sorted out, but, anyway, 
let's talk about what I think at this moment.

1. The main issue seems to be how some referencing object comparisons 
are made. I noticed that AbstractDatum.equals() compares Datums using 
the Datum name if "compareMetadata" is false. That's kind of OK in 
principle, but it doesn't even check if there are Identifiers present. 
So, try running the WKT below, where I just changed the line 
"DATUM[\"South American Datum 1969\"," to "DATUM[\"SAD69\"," (and please 
shoot me if it works, I just tested *exactly* that one). Note that the 
Datum does have an identifier stating it to be EPSG:6618.

String wkt = "PROJCS[\"SAD69 / UTM zone 23S\"," +
"GEOGCS[\"SAD69\"," +
"DATUM[\"SAD69\"," +
"SPHEROID[\"GRS 1967 (SAD69)\", 6378160.0, 298.25, 
AUTHORITY[\"EPSG\",\"7050\"]]," +
"TOWGS84[-66.87, 4.37, -38.52, 0.0, 0.0, 0.0, 0.0]," +
"AUTHORITY[\"EPSG\",\"6618\"]]," +
"PRIMEM[\"Greenwich\", 0.0, AUTHORITY[\"EPSG\",\"8901\"]]," +
"UNIT[\"degree\", 0.017453292519943295]," +
"AXIS[\"Geodetic longitude\", EAST]," +
"AXIS[\"Geodetic latitude\", NORTH]," +
"AUTHORITY[\"EPSG\",\"4618\"]]," +
"PROJECTION[\"Transverse Mercator\", AUTHORITY[\"EPSG\",\"9807\"]]," +
"PARAMETER[\"central_meridian\", -45.0]," +
"PARAMETER[\"latitude_of_origin\", 0.0]," +
"PARAMETER[\"scale_factor\", 0.9996]," +
"PARAMETER[\"false_easting\", 500000.0]," +
"PARAMETER[\"false_northing\", 10000000.0]," +
"UNIT[\"m\", 1.0]," +
"AXIS[\"Easting\", EAST]," +
"AXIS[\"Northing\", NORTH]]";

So, this one is not found, which sounds odd to me because it would mean 
that WKTs generated by 3rd parties must know the exact Datum name used 
in this particular epsg database. If that's really how it works, well, 
then I'll have to fix it by hand before making the lookup.

1b. This issue has other consequences. Since many such object 
comparisons seem to ignore the EPSG code, I have the following problem: 
a WKT was generated where the numbers were printed with only a certain 
level of precision, making the WGS84 SPHEROID's inverse flattening turn 
from 298.257223563 to 298.2572. Again, I agree that if no EPSG code is 
given, then it is a different SPHEROID. But when the code IS given, the 
comparison still ignores it and says they are different. So, if I get 
the WKT for 32723 (WGS84 UTM 23S) and simply exchange the spheroid's 
flattening as above, it does not work either:

String wkt = "PROJCS[\"WGS 84 / UTM zone 23S\"," +
"GEOGCS[\"WGS 84\"," +
"DATUM[\"World Geodetic System 1984\"," +
"SPHEROID[\"WGS 84\", 6378137.0, 298.2572, AUTHORITY[\"EPSG\",\"7030\"]]," +
"AUTHORITY[\"EPSG\",\"6326\"]]," +
"PRIMEM[\"Greenwich\", 0.0, AUTHORITY[\"EPSG\",\"8901\"]]," +
"UNIT[\"degree\", 0.017453292519943295]," +
"AXIS[\"Geodetic latitude\", NORTH]," +
"AXIS[\"Geodetic longitude\", EAST]," +
"AUTHORITY[\"EPSG\",\"4326\"]]," +
"PROJECTION[\"Transverse Mercator\", AUTHORITY[\"EPSG\",\"9807\"]]," +
"PARAMETER[\"central_meridian\", -45.0]," +
"PARAMETER[\"latitude_of_origin\", 0.0]," +
"PARAMETER[\"scale_factor\", 0.9996]," +
"PARAMETER[\"false_easting\", 500000.0]," +
"PARAMETER[\"false_northing\", 10000000.0]," +
"UNIT[\"m\", 1.0]," +
"AXIS[\"Easting\", EAST]," +
"AXIS[\"Northing\", NORTH]]";

What do you think? Should it work this way? My personal opinion is that 
the identifiers should be checked always when available..

2. A smaller different issue may be related to that forceXY property you 
have in your code. I guess it turns out that now a ProjectedCRS such as 
these which specify "AXIS[\"x\", EAST]" and "AXIS[\"y\", NORTH]" instead 
of "Easting" and "Northing" don't work? It looked like that in the axis 
comparison code in DefaultCoordinateSystemAxis.equals(). But I tried to 
set the property and it still didn't work (see WKT below). Am I missing 
something? By the way, why is this not a hint, but a system property? 
Oh, just to let you know: this isn't so much of an issue for me right 
now, it's just to know what works and what doesn't.

String wkt = "PROJCS[\"SAD69 / UTM zone 23S\"," +
"GEOGCS[\"SAD69\"," +
"DATUM[\"South American Datum 1969\"," +
"SPHEROID[\"GRS 1967 (SAD69)\", 6378160.0, 298.25, 
AUTHORITY[\"EPSG\",\"7050\"]]," +
"TOWGS84[-66.87, 4.37, -38.52, 0.0, 0.0, 0.0, 0.0]," +
"AUTHORITY[\"EPSG\",\"6618\"]]," +
"PRIMEM[\"Greenwich\", 0.0, AUTHORITY[\"EPSG\",\"8901\"]]," +
"UNIT[\"degree\", 0.017453292519943295]," +
"AXIS[\"Geodetic longitude\", EAST]," +
"AXIS[\"Geodetic latitude\", NORTH]," +
"AUTHORITY[\"EPSG\",\"4618\"]]," +
"PROJECTION[\"Transverse Mercator\", AUTHORITY[\"EPSG\",\"9807\"]]," +
"PARAMETER[\"central_meridian\", -45.0]," +
"PARAMETER[\"latitude_of_origin\", 0.0]," +
"PARAMETER[\"scale_factor\", 0.9996]," +
"PARAMETER[\"false_easting\", 500000.0]," +
"PARAMETER[\"false_northing\", 10000000.0]," +
"UNIT[\"m\", 1.0]," +
"AXIS[\"x\", EAST]," +
"AXIS[\"y\", NORTH]]";

*uff*, that's all for now then. Sorry for the long e-mail!

Cheers
Milton

Andrea Aime wrote:
> Milton Jonathan ha scritto:
>> Hi there
>>
>> Once again, the axis issue attacked me :P. I guess this may have 
>> something to do with changes from 2.5.x to 2.6.x (has LONGITUDE_FIRST 
>> become the default now?).
>>
>> The situation is like this: I have a (valid) WKT generated for a 
>> Projected CRS with longitude first (it is a SouthAmericanDatum1969 
>> with a UTM projection). When I try to lookup its epsg code, it tries, 
>> tries, and fails. But it fails due to the axis order of the DATUM 
>> (i.e., not due to the projected part).
>>
>> To illustrate, here goes the code that fails:
>>
>> String wkt = "PROJCS[\"SAD69 / UTM zone 23S\"," +
>> "GEOGCS[\"SAD69\"," +
>> "DATUM[\"South American Datum 1969\"," +
>> "SPHEROID[\"GRS 1967 (SAD69)\", 6378160.0, 298.25, 
>> AUTHORITY[\"EPSG\",\"7050\"]]," +
>> "TOWGS84[-66.87, 4.37, -38.52, 0.0, 0.0, 0.0, 0.0]," +
>> "AUTHORITY[\"EPSG\",\"6618\"]]," +
>> "PRIMEM[\"Greenwich\", 0.0, AUTHORITY[\"EPSG\",\"8901\"]]," +
>> "UNIT[\"degree\", 0.017453292519943295]," +
>> "AXIS[\"Geodetic longitude\", EAST]," +
>> "AXIS[\"Geodetic latitude\", NORTH]," +
>> "AUTHORITY[\"EPSG\",\"4618\"]]," +
>> "PROJECTION[\"Transverse Mercator\", AUTHORITY[\"EPSG\",\"9807\"]]," +
>> "PARAMETER[\"central_meridian\", -45.0]," +
>> "PARAMETER[\"latitude_of_origin\", 0.0]," +
>> "PARAMETER[\"scale_factor\", 0.9996]," +
>> "PARAMETER[\"false_easting\", 500000.0]," +
>> "PARAMETER[\"false_northing\", 10000000.0]," +
>> "UNIT[\"m\", 1.0]," +
>> "AXIS[\"Easting\", EAST]," +
>> "AXIS[\"Northing\", NORTH]]";
>>
>> CoordinateReferenceSystem crs = CRS.parseWKT(wkt);
>> Integer srid = CRS.lookupEpsgCode(crs, true);
> 
> Milton, I just tried your code on my computer and the lookup
> returns back 29193 for me, no matter how I set the
> axis order hints (or even if I set them at all).
> 
> Full program attached. What am I missing?
> 
> Cheers
> Andrea
> 
> 

-- 

Milton Jonathan
Grupo GIS e Meio Ambiente
Tecgraf/PUC-Rio
Tel: +55-21-3527-2502

Re: [Geotools-devel] Axis problems again with CRS.lookupEpsgCode

Toolkit for working with and mapping geospatial data

Re: [Geotools-devel] Axis problems again with CRS.lookupEpsgCode