矢量转栅格 栅格转矢量

2019-11-21  本文已影响0人  不玩了啊

#!/usr/bin/env python

# -*- coding: utf-8 -*-

from osgeo import ogr

from osgeo import gdal

# set pixel size

pixel_size = 0.00002

no_data_value = -9999

# Shapefile input name

# input projection must be in cartesian system in meters

# input wgs 84 or EPSG: 4326 will NOT work!!!

input_shp = r'../geodata/data/ply_golfcourse-strasslach3857.shp'

# TIF Raster file to be created

output_raster = r'../geodata/data/ply_golfcourse-strasslach.img'

# Open the data source get the layer object

# assign extent coordinates

#shp_driver = ogr.GetDriverByName('ESRI Shapefile')

open_shp = ogr.Open(input_shp)

shp_layer = open_shp.GetLayer()

x_min, x_max, y_min, y_max = shp_layer.GetExtent()

print(x_min)

print(x_max)

# calculate raster resolution

x_res = int((x_max - x_min) / pixel_size)

y_res = int((y_max - y_min) / pixel_size)

print(x_res)

print(y_res)

# set the image type for export

image_type = 'GTiff'

image_type = 'HFA'

driver = gdal.GetDriverByName(image_type)

driver.Register()

# create a new raster takes Parameters

# Filename    the name of the dataset to create. UTF-8 encoded.

# nXSize    width of created raster in pixels.

# nYSize    height of created raster in pixels.

# nBands    number of bands.

# eType    type of raster.

new_raster = driver.Create(output_raster, x_res, y_res, 1, gdal.GDT_Byte)

new_raster.SetGeoTransform((x_min, pixel_size, 0, y_max, 0, -pixel_size))

# get the raster band we want to export too

raster_band = new_raster.GetRasterBand(1)

# assign the no data value to empty cells

raster_band.SetNoDataValue(no_data_value)

# run vector to raster on new raster with input Shapefile

gdal.RasterizeLayer(new_raster, [1], shp_layer, burn_values=[255])

#!/usr/bin/env python

# -*- coding: utf-8 -*-

from osgeo import ogr

from osgeo import gdal

#  get raster datasource

open_image = gdal.Open( "../geodata/cadaster_borders-2tone-black-white.png" )

input_band = open_image.GetRasterBand(3)

#  create output datasource

output_shp = "../geodata/data/cadaster_raster"

shp_driver = ogr.GetDriverByName("ESRI Shapefile")

# create output file name

output_shapefile = shp_driver.CreateDataSource( output_shp + ".shp" )

new_shapefile = output_shapefile.CreateLayer(output_shp, srs = None )

gdal.Polygonize(input_band, None, new_shapefile, -1, [], callback=None)

new_shapefile.SyncToDisk()

上一篇下一篇

猜你喜欢

热点阅读