2

I am trying to create a VRT file using Python (GDAL). I am getting the following error:

"ValueError: array larger than output file, or offset off edge".

How to solve this error?

The code written for the same is as follows.

import osr
import gdal
import numpy as np
from gdalconst import *

ds = gdal.Open('/home/vigna/Downloads/2740115314.GMODO-SVM12_npp_d20170824_t0721423_e0727209_b30176_c20170824140610010359_nobc_ops.h5')

subds = ds.GetSubDatasets()
#print(subds)
bt = gdal.Open('HDF5:"/home/vigna/Downloads/2740115314.GMODO-SVM12_npp_d20170824_t0721423_e0727209_b30176_c20170824140610010359_nobc_ops.h5"://All_Data/VIIRS-M12-SDR_All/BrightnessTemperature')
bt_array = bt.ReadAsArray()
#print(bt_array)
long = gdal.Open('HDF5:"/home/vigna/Downloads/2740115314.GMODO-SVM12_npp_d20170824_t0721423_e0727209_b30176_c20170824140610010359_nobc_ops.h5"://All_Data/VIIRS-MOD-GEO_All/Longitude')
long_array = long.ReadAsArray()
lat = gdal.Open('HDF5:"/home/vigna/Downloads/2740115314.GMODO-SVM12_npp_d20170824_t0721423_e0727209_b30176_c20170824140610010359_nobc_ops.h5"://All_Data/VIIRS-MOD-GEO_All/Latitude')
lat_array = lat.ReadAsArray()

def raster(ds):
    xsize = ds.RasterXSize
    ysize = ds.RasterYSize
    GeoT = ds.GetGeoTransform()
    Projection = osr.SpatialReference()
    Projection.ImportFromWkt(ds.GetProjectionRef())
    DataType = 'UInt16'
    X_dataset=long_array
    Y_dataset=lat_array
    return xsize, ysize, GeoT, Projection,DataType,X_dataset,Y_dataset

def Create_vrt(Name, Array, driver,xsize, ysize, GeoT, Projection,X_dataset,Y_dataset):
#for DataType:
    DataType = gdal.GDT_Float32
    NewFileName = 'bt3'+'.vrt'
    DataSet = driver.Create( NewFileName, xsize,ysize, 1, DataType )
        # the '1' is for band 1.
    DataSet.SetGeoTransform(GeoT)
    DataSet.SetProjection( Projection.ExportToWkt() )
    # Write the array
    DataSet.GetRasterBand(1).WriteArray( Array )
    return NewFileName

DataSet = ds
Band = DataSet.GetRasterBand(1)
# Open as an array.
Array = bt_array
NewArray = Array
print(NewArray)
xsize, ysize, GeoT, Projection, DataType, X_dataset, Y_dataset = raster(ds)
driver = gdal.GetDriverByName('VRT')
NewFileName = Create_vrt('bt3', NewArray, driver, xsize, ysize, GeoT, Projection, X_dataset, Y_dataset)
#Create_vrt()
print("ok")
2
  • 2
    You can't write an array to a VRT. A VRT is a virtual raster (just a small XML text file). Write to something like a GeoTIFF instead.
    – user2856
    Nov 6, 2017 at 5:20
  • 2
    If you don't want to write a GeoTIFF to the disk, you could use the MEM (in-memory raster) driver (driver = gdal.GetDriverByName('MEM')), or write a GeoTIFF to the /vsimem virtual in-memory filesystem (driver = gdal.GetDriverByName('GTIFF'); driver.Create('/vsimem/some.tiff', etc...))
    – user2856
    Nov 6, 2017 at 6:11

1 Answer 1

1

I would recommend as suggested here: Python GDAL: Save array as raster with projection from other file

write the array as a raster

def array_to_raster(array):
    """Array > Raster
    Save a raster from a C order array.

    :param array: ndarray
    """
    dst_filename = '/a_file/name.tiff'

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.