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")
Taras
  • 32,823
  • 4
  • 66
  • 137
vigna purohit
  • 441
  • 2
  • 12
  • 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 06 '17 at 05: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 06 '17 at 06:11

1 Answers1

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'
das-g
  • 1,445
  • 2
  • 11
  • 33
whyzar
  • 12,053
  • 23
  • 38
  • 72