2

I have lat, long, and data arrays

data = np.array [[1,3],[2,3],[3,4]]

lat = np.array [[44,43],[46,44],[43,45]]

long = np.array [[33,32],[31,33],[36,33]]

The coordinates here are random but my question is how to get the rasterio transform from such objects, my original image is a slanted rectangle so the rasterio.transform.from_bounds(west, south, east, north, width, height) didn't work out.

My final goal is to write the data array as a GeoTIFF.

0

1 Answer 1

3

You need rasterio.transform.from_gcps:

  1. Get the upper-left, upper-right, lower-left & lower-right corner coordinates from your arrays in lon, lat / x, y order.
  2. Create GCPs for each corner
  3. Pass a list of GCPs to the from_gcps function

In the example below I've used an actual dataset that I know the coordinates of.

import rasterio as rio
from rasterio.transform import from_gcps
from rasterio.control import GroundControlPoint as GCP

ul = (10.0000000, 30.0000000)  # in lon, lat / x, y order
ll = (80.7106781, -40.7106781)
ur = (80.7110000, 100.711)
lr = (151.4213562, 30.0000000)
cols, rows = 20, 20

gcps = [
    GCP(0, 0, *ul),
    GCP(0, cols, *ur),
    GCP(rows, 0, *ll),
    GCP(rows, cols, *lr)
]

transform = from_gcps(gcps)

# Check derived transform angainst dataset transform
with rio.open('byte_rotated.tif') as tif:
    print(tif.transform)
    print(transform)

Output:

| 3.54, 3.54, 10.00|
| 3.54,-3.54, 30.00|
| 0.00, 0.00, 1.00|

| 3.54, 3.54, 10.00|
| 3.54,-3.54, 30.00|
| 0.00, 0.00, 1.00|

enter image description here

However your arrays of coordinates likely represent the centre of each pixel, so your GCPs should look like:

gcps = [
    GCP(0.5, 0.5, *ul),  # Centre of top-left pixel
    GCP(0.5, cols, *ur),
    GCP(rows-0.5, 0.5, *ll),
    GCP(rows-0.5, cols-0.5, *lr)
]

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.