2

Continuing what is written in this post, I have another question:

for each point (X, Y) of my GEOTIFF file I extract the coordinates as follows:

Envelope2D pixelEnvelop = geometry.gridToWorld(new GridEnvelope2D(x, y, 1, 1));

double[] coords = new double[]{pixelEnvelop.getCenterX(), pixelEnvelop.getCenterY()};

For example I obtain

-577460.1192321777, 640293.6331481934
-577460.1192321777, 639293.9553527832
... (and so on)

However, I would like to have latitude and longitude available. Is it possible to convert them?

For information, I print also the content of pixelEnvelop.getCoordinateReferenceSystem() method

PROJCS["unnamed", 
  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 longitude", EAST], 
    AXIS["Geodetic latitude", NORTH], 
    AUTHORITY["EPSG","4326"]], 
  PROJECTION["Transverse_Mercator"], 
  PARAMETER["central_meridian", 12.5], 
  PARAMETER["latitude_of_origin", 42.0], 
  PARAMETER["scale_factor", �], 
  PARAMETER["false_easting", 0.0], 
  PARAMETER["false_northing", 0.0], 
  UNIT["m", 1.0], 
  AXIS["Easting", EAST], 
  AXIS["Northing", NORTH]]

As you can see, scale_factor parameter does not contain a readable value.

1

It is fairly simple to reproject points in GeoTools, as the following shows. If GeoTools is unable to work out your projection from the coverage then you can use an EPSG code that represents your data's projection.

AbstractGridFormat format = GridFormatFinder.findFormat(args[0]);
Hints hints = null;
if (format instanceof GeoTiffFormat) {
  hints = new Hints(Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER, Boolean.TRUE);
}
AbstractGridCoverage2DReader reader = format.getReader(args[0],hints);

GridCoverage2D coverage = null;
try {
  coverage = reader.read(null);
} catch (IOException giveUp) {
  throw new RuntimeException(giveUp);
}
reader.dispose();
CoordinateReferenceSystem wgs84= CRS.decode("EPSG:4326", true);
CoordinateReferenceSystem target = coverage.getCoordinateReferenceSystem();//CRS.decode("EPSG:3997", true);
GeometryFactory gf = new GeometryFactory();
MathTransform targetToWgs = CRS.findMathTransform( target,wgs84);
GridEnvelope gridRange2D = coverage.getGridGeometry().getGridRange();
int count = 0;
for (int i = gridRange2D.getLow(0); i < gridRange2D.getHigh(0); i++) {
  for (int j = gridRange2D.getLow(1); j < gridRange2D.getHigh(1); j++) {
    if(count++>10)System.exit(1);
    GridCoordinates2D coord = new GridCoordinates2D(i, j);
    DirectPosition p = coverage.getGridGeometry().gridToWorld(coord);
    Point point = gf.createPoint(new Coordinate(p.getOrdinate(0),p.getOrdinate(1)));
    Geometry wgsP = JTS.transform(point, targetToWgs);
    System.out.println(i+","+j+" -> "+point+" -> "+wgsP);
  }
}

You might want to check what the storage order of your coverage is using gdalinfo to decide if you want to loop through in row or column first order.

  • When I execute the instruction Geometry wgsP = JTS.transform(point, targetToWgs); I get the following error: org.geotools.referencing.operation.projection.ProjectionException: Transformation doesn't convergence. – Alessio Frabotta Apr 6 '18 at 8:08
  • then there is probably a problem understanding your coverage's projection, try setting it from the epsg code – Ian Turton Apr 6 '18 at 8:29
  • How can I do this? – Alessio Frabotta Apr 6 '18 at 8:33
  • use something like CRS.decode("EPSG:3997", true); epag.io and prj2epsg.org will help you work out the correct code, based on knowing where your image is – Ian Turton Apr 6 '18 at 8:36
  • for my purpose I need to find a Coordinate reference system valid for Italy, but I do not know which to choose – Alessio Frabotta Apr 6 '18 at 10:41

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.