I am working on a plugin for Qgis to calculate spatial Kernel density maps. I have all the calculations working, all I am missing is a way to turn a Numpy Array, with density values into a multiband raster layer.
Do I have to create a geotiff on a temp file using Gdal and then load it?
Or is there a direct way to create the layer from data in memory?
if so, how to do it?
Answer
Here is the code that I use to convert an array to gdal raster saving it to the disk, "param" is a dicitionary containing gdal parameters (check the gdal documentation) and "array" is a numpy array. Than you can instantiate a QgsMapLayer with your file as source. You have to create the geotiff in the disk.
from osgeo import gdal as osgdal # Adapt the import to fit yor environement.
driver = osgdal.GetDriverByName(param['out_format'])
dataset = driver.Create(
param['dst_filename'],
param['x_pixels'],
param['y_pixels'],
1,
osgdal.GDT_Float32,
)
dataset.SetGeoTransform((
param['xmin'], #0
param['pixel_size'], #1
0, #2
param['ymin'], #3
0, #4
param['pixel_size'])) #5
out_srs = osr.SpatialReference()
out_srs.ImportFromEPSG(param['SRID'])
dataset.SetProjection(out_srs.ExportToWkt())
dataset.GetRasterBand(1).WriteArray(array.T) # Remove "T" if it's inverted.
dataset = None
No comments:
Post a Comment