It is quite easy to get raster width and height using snappy:
width = firstBand.getRasterWidth()
height = firstBand.getRasterHeight()
but I can’t http://step.esa.int/docs/v5.0/apidoc/engine/org/esa/snap/core/datamodel/AbstractBand.html a method to get a particular pixel size, so that I will know the spatial resolution of my raster file.
I’ve tried to implement this as follows:
but it seems that my function returns different resolution/pixel size than that which can see gdalinfo (file which I used in testing is contained within the repository):
Since gdalinfo do not handle all ESA SNAP formats, I can’t just call gdalinfo in my script.
I thought about converting each product to GeoTIFF firstly and then call gdalinfo, but it will be very time consuming.
In most cases it is sufficient to compute the resolution just at the centre pixel.
As far as I understand your calculations correctly, at least part of the problems lies in the way your trying to extract the corners coodinates. By taking the min/max lat/lon values and calculating the difference to get the extent, you obtain the extent along a particular longitude or latitude without considering the fact of the rotation of the scene relative to the lat/lon grid - i.e. two corner points have no lat or lon value in common.
Maybe try the following and see if this works out:
import numpy as np
from geographiclib.geodesic import Geodesic
corner_top_left = prod.getSceneGeoCoding().getGeoPos(snappy.PixelPos(0,0), None)
corner_bottom_left = prod.getSceneGeoCoding().getGeoPos(snappy.PixelPos(0, prod.getSceneRasterHeight()), None)
corner_top_right = prod.getSceneGeoCoding().getGeoPos(snappy.PixelPos(prod.getSceneRasterWidth(), 0), None)
height_range_m = Geodesic.WGS84.Inverse(corner_top_left.getLat(), corner_top_left.getLon(), corner_bottom_left.getLat(), corner_bottom_left.getLon())["s12"]
width_range_m = Geodesic.WGS84.Inverse(corner_top_left.getLat(), corner_top_left.getLon(), corner_top_right.getLat(), corner_top_right.getLon())["s12"]
height_n_pixels = prod.getSceneRasterHeight()
width_n_pixels = prod.getSceneRasterWidth()
mean_pixel_size = (round(height_range_m/height_n_pixels,2), round(width_range_m/width_n_pixels,2))