I have a GeoTIFF raster file in the EPSG 4326 (WGS 84) geographic coordinate system.

I want to crop a 20 km x 20 km square around a particular (lat, long) coordinate, and save that as a new raster.

But since the units for my image are in lat long degrees, but I want a box in terms of kilometers, is the best way to do this to use the Haversine formula to convert my (lat, long) degrees to coordinates in meters, compute the bounding box in meters, and then invert the Haversine computation to get back the box in degrees?

OR do I just project the entire image into a projected coordinate system, so that I can then easily find the bounding box in meters, and then re-project all the bounding box coords back into degrees in my original geographic coordinate system?

Is there any difference in accuracy when doing this? I understand that using Haversine may not be that accurate, but it seems to me that projecting would also lose some accuracy?

(I'm doing this in Python with rasterio or GDAL if that makes a difference.)

My data is around lon: -120.011597902, lat: 39.54719334000001
Resolution (cell size): 0.0035 deg/pixel
- I think this is around 300 meters/pixel (not sure).

| improve this question | | | | |
  • 2
    Given the small area of 20x20 km your best bet is to reproject/transform the grid to the locally commonly used coordinate system, your distance errors will be acceptable for most usage. Can you elaborate what the future use of the cropped 20x20 km raster will be, what is your cell size eg? Is there a special reason why you need to work with geographic coordinates in such a small area? – Hans Erren Jun 1 '19 at 22:57
  • I want to use a statistical model to return some estimates of properties of each 20x20km area. I have a much larger area, and I want to grid the entire area into 20x20 boxes around specific points of interest. Each box will then be individually fed into my model. – rasen58 Jun 3 '19 at 15:04
  • Also, 20x20 km is actually still too large, since I want to get more precise estimates of the properties around each of my points of interest, so I'm planning on making these boxes perhaps even 5x5 km. – rasen58 Jun 3 '19 at 15:05
  • So I tried to do the projection, but is the distance too small already? Since when I projected nearby points, the get turned into the same projected coordinates, so the distance between them is 0: pastebin.com/uekGa3QM – rasen58 Jun 3 '19 at 15:29
  • 1
    Does no-one else see a preoblem with the statement "a GeoTIFF raster file in the EPSG 4326 (WGS 84) geographic coordinate system"? A raster file -- of square pixels -- is still a projection, no matter what the units are. See gis.stackexchange.com/questions/146321/… – Martin F Jun 3 '19 at 19:34

Your Answer

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

Browse other questions tagged or ask your own question.