Обрезать растр в массивы на основе шейп-файла или полигонов - PullRequest
Новая
       13

Обрезать растр в массивы на основе шейп-файла или полигонов

0 голосов
/ 26 апреля 2018

У меня есть большой растровый файл, который я хочу обрезать до 6 массивов на основе полигонов или шейп-файлов.У меня есть шейп-файл, а также полигоны в качестве кадра данных геопанды.Может кто-нибудь, пожалуйста, помогите мне, как это сделать с помощью Python (без arcpy)

1 Ответ

0 голосов
/ 26 апреля 2018

Я создал маленький генератор, который должен делать то, что вам нужно.Я выбрал генератор, а не перебирал объекты напрямую, потому что удобнее, если вы хотите проверить массивы.Если вы хотите, вы все равно можете перебрать генератор. 

import gdal
import ogr, osr

# converts coordinates to index

def bbox2ix(bbox,gt):
    xo = int(round((bbox[0] - gt[0])/gt[1]))
    yo = int(round((gt[3] - bbox[3])/gt[1]))
    xd = int(round((bbox[1] - bbox[0])/gt[1]))
    yd = int(round((bbox[3] - bbox[2])/gt[1]))
    return(xo,yo,xd,yd)

def rasclip(ras,shp):
    ds = gdal.Open(ras)
    gt = ds.GetGeoTransform()

    driver = ogr.GetDriverByName("ESRI Shapefile")
    dataSource = driver.Open(shp, 0)
    layer = dataSource.GetLayer()

    for feature in layer:

        xo,yo,xd,yd = bbox2ix(feature.GetGeometryRef().GetEnvelope(),gt)
        arr = ds.ReadAsArray(xo,yo,xd,yd)
        yield arr

    layer.ResetReading()
    ds = None
    dataSource = None

Если ваш шейп-файл называется shapefile.shp, а ваш растр big_raster.tif, вы можете использовать его следующим образом: 

gen = rasclip('big_raster.tif','shapefile.shp')

# manually with next

clip = next(gen)

## some processing or inspection here

# clip with next feature

clip = next(gen)

# or with iteration

for clip in gen:

    ## apply stuff to clip
    pass # remove
...