Topics

Transform Conversion to Lat/Long


j.carnes2553@...
 

Hello everyone,

In the picture below I am using rasterio to open a .tif file and then transform the those bounds to get the upper left and lower right points of the dataset. My question is around how to convert a point that is not from the bounds and is just an x,y point in meters in a specific epsg. I want to be able to convert multiple points that are in the x,y meter format to a lat/long format all within this epsg.

Another question I have is around how applicable this is to creating new distances. For example is I am at the upper left spot and increase the x by 5 meters and convert that coordinate back to lat/long, does that mean that the two lat/long coordinates are now exactly 5 meters apart? 

I am quite new to the GIS world so any help would be greatly appreciated, thanks!


Guillaume Lostis
 

Hi,

From the screenshot you have provided, we can see that the CRS of your image (EPSG:32617) is a UTM zone, so its unit is indeed in meters. That means that if you take a point with coordinates expressed in the CRS of your image, that you add 5 meters to the x coordinates, and that you convert both points to lat/lon, your two points will indeed be separated by 5 meters.

Here is a code snippet showing what you want to achieve, that uses the pyproj package to do the coordinate conversion:

import rasterio
import pyproj

with rasterio.open("S2A_MSIL1C_20190701T174911_N0207_R141_T13UFQ_20190701T212910__B11.jp2") as src:
    bounds = src.bounds
    crs = src.crs

left, bottom, right, top = bounds

in_proj = pyproj.Proj(init="epsg:{}".format(crs.to_epsg()))
out_proj = pyproj.Proj(init="epsg:4326")

lons, lats = pyproj.transform(in_proj, out_proj, [left, left + 5], [top, top])

for lat, lon in zip(lats, lons):
    print(lat, lon)

This snippet uses a Sentinel 2 band with a different CRS than yours, but it should work the same way.

A few things to note:

  • you can directly access the left, bottom, right and top bounds from the bounds object, there is no need to use the transform for that.
  • EPSG:4326 is the CRS for coordinates expressed in lat/lon
  • pyproj.transform can take either single x, y coordinates, or iterables of x, y coordinates (as in the example above where I pass a list of two x's and two y's). If you want to convert a lot of points between two CRS's, it is much more efficient to pass them as a list in a single call to pyproj.transform because the transformation is vectorized.

Guillaume Lostis


j.carnes2553@...
 

Guillaume, thanks for the swift reply, this sounds exactly like what I'm looking for. I will work on implementing it into my code now, thank again :)