I wanted to use geotools (15.3) to reproject data from EPSG:4326 to a local projection. To see if everything works as expected I printed the length of two test linestrings. The expected result is around 200 meters.
However I noticed that the result varies greatly depending on how I define the source CoordinateReferenceSystem. I thought these two CRS should be identical, however they are not:
- DefaultGeographicCRS.WGS84 -> works as expected
- CRS.decode("EPSG:4326") -> completely wrong results
For DefaultGeographicCRS.WGS84 the defined CRS and the (correct) lengths are:
GEOGCS["WGS84(DD)",
DATUM["WGS84",
SPHEROID["WGS84", 6378137.0, 298.257223563]],
PRIMEM["Greenwich", 0.0],
UNIT["degree", 0.017453292519943295],
AXIS["Geodetic longitude", EAST],
AXIS["Geodetic latitude", NORTH]]
northSouth: 199.70672773686056m
eastWest: 200.33865746495295m
For CRS.decode("EPSG:4326") the defined CRS and the (completely wrong) lengths are:
GEOGCS["WGS 84",
DATUM["World Geodetic System 1984",
SPHEROID["WGS 84", 6378137.0, 298.257223563, 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"]]
northSouth: 222.76855610197927m
eastWest: 346.2445236833851m
How is this possible or what am I doing wrong? The result is the same with geotools 14.0.
import org.geotools.factory.Hints;
import org.geotools.geometry.jts.JTS;
import org.geotools.geometry.jts.JTSFactoryFinder;
import org.geotools.geometry.jts.WKTReader2;
import org.geotools.referencing.CRS;
import org.geotools.referencing.crs.DefaultGeographicCRS;
import org.opengis.referencing.crs.CoordinateReferenceSystem;
import org.opengis.referencing.operation.MathTransform;
import com.vividsolutions.jts.geom.Geometry;
import com.vividsolutions.jts.geom.GeometryFactory;
public class TransformFailMain {
public static void main(String[] args) throws Exception {
CoordinateReferenceSystem sourceCRS = CRS.decode("EPSG:4326");
// CoordinateReferenceSystem sourceCRS = DefaultGeographicCRS.WGS84;
CoordinateReferenceSystem targetCRS = CRS.decode("EPSG:31256");
Hints hints = new Hints(Hints.CRS, sourceCRS);
GeometryFactory geometryFactory = JTSFactoryFinder.getGeometryFactory(hints);
WKTReader2 reader = new WKTReader2(geometryFactory);
// two lineStrings in Vienna, both around 200m long
Geometry linestringNorthSouth = reader.read("LINESTRING(16.368831 48.201796,16.368831 48.2)");
Geometry linestringEastWest = reader.read("LINESTRING(16.368831 48.2,16.371526 48.2)");
MathTransform transform = CRS.findMathTransform(sourceCRS, targetCRS);
System.out.println(sourceCRS);
System.out.println("northSouth: " + JTS.transform(linestringNorthSouth, transform).getLength() + "m");
System.out.println("eastWest: " + JTS.transform(linestringEastWest, transform).getLength() + "m");
}
}
Best Answer
I stumbled over a common pitfall.
So you should choose wisely which CRS you use when e.g. parsing external WKT strings as I did.
More in-depth explanations of why this is the case are found in https://stackoverflow.com/a/13579921/1648538 and http://docs.geotools.org/latest/userguide/library/referencing/order.html. (Thank you for the hint, iant!)
Good to know: the order of lat/lon can be changed in geotools: