Export VRT to GeoTIFF in Python

Question

I'm trying to export a Virtual raster (VRT) to GeoTIFF in Python 2.7 (Windows 10 System), but code runs with some problems.

When i try export the following image (.vrt file) to GeoTIFF i had a problem that i don't know how to resolve and this is the image when i export the vrt file

The code used is the following:

from osgeo import gdal
import os, glob, osr, numpy as np
path = 'D:\Tesis\Imagenes\L8 OLI TIRS'
VRT = "C:/DUMP_PYTHON"
pathout = "C:\Users\pp\Documents\LANDSAT\NUEVAS IMGS LANDSAT\Input ETa"

def save_raster ( output_name, dataset, raster_data, driver ,NaN_Value):
    """
    A function to save a 1-band raster using GDAL to the file indicated
    by ``output_name``. It requires a GDAL-accesible dataset to collect 
    the projection and geotransform.
    """
    # Open the reference dataset
    g = ( dataset )
    # Get the Geotransform vector
    if raster_data is False:
        raster_data=g.ReadAsArray()
    geo_transform = g.GetGeoTransform ()
    x_size = g.RasterXSize # Raster xsize
    y_size = g.RasterYSize # Raster ysize
    srs = g.GetProjectionRef () # Projection
    NaN_rast = NaN_Value
    # raster_data[raster_data == NaN_rast] = 'NaN'
    raster_data[raster_data == NaN_rast] = np.NaN
    # Need a driver object. By default, we use GeoTIFF
    driver = gdal.GetDriverByName ( driver )
    dataset_out = driver.Create ( output_name, x_size, y_size, 1, \
            gdal.GDT_Float32 )
    dataset_out.SetGeoTransform ( geo_transform )
    dataset_out.SetProjection ( srs )
    dataset_out.GetRasterBand ( 1 ).WriteArray ( \
            raster_data.astype(np.float32) )


#READ DATA
i=672
files= glob.glob(os.path.join(path,'*.tif'))
albedo = gdal.Open(VRT+"/dummy_alb"+str(i)+'.vrt')
save_raster(pathout+'/Albedo_L8_'+files[i][39:46]+'.tif',albedo1,False,"GTiff",99999)

i also tried with the following code, but the result is the same:

i=672    
albedo = gdal.Open(VRT+"/dummy_alb"+str(i)+'.vrt')
        src_ds = gdal.Open(VRT+"/dummy_alb"+str(i)+'.vrt')
#Open output format driver, see gdal_translate --formats for list
    format = "GTiff"
    driver = gdal.GetDriverByName( format )
#Output to new format
    dst_ds = driver.CreateCopy( pathout+'/Albedo_L8.tif', src_ds, -10 )

#Properly close the datasets to flush to disk
    dst_ds = None
    src_ds = None

So.. where is the problem?

Answer 1

GDAL does not write everything to file directly but keeps it in cache until the limit is reached. If you're done flush the cache to file and make sure to destroy all GDAL-related objects:

import numpy as np
from osgeo import gdal
def save_raster(output_name, dataset, raster_data, driver, NaN_Value):
    """
    A function to save a 1-band raster using GDAL to the file indicated
    by ``output_name``. It requires a GDAL-accesible dataset to collect
    the projection and geotransform.
    """
    # Open the reference dataset
    g = gdal.Open(dataset)
    # Get the Geotransform vector
    if raster_data is False:
        raster_data = g.ReadAsArray()
    geo_transform = g.GetGeoTransform()
    x_size = g.RasterXSize  # Raster xsize
    y_size = g.RasterYSize  # Raster ysize
    srs = g.GetProjectionRef()  # Projection
    g = None
    # Need a driver object. By default, we use GeoTIFF
    driver = gdal.GetDriverByName(driver)
    dataset_out = driver.Create(output_name, x_size, y_size, 1, gdal.GDT_Float32)
    driver = None
    dataset_out.SetGeoTransform(geo_transform)
    dataset_out.SetProjection(srs)
    band_out = dataset_out.GetRasterBand(1)
    band_out.SetNoDataValue(NaN_Value)
    band_out.WriteArray(raster_data.astype(np.float32))
    band_out.FlushCache()
    band_out = None
    dataset_out = None

Or alternatively use my package spatialist:

from spatialist import Raster

i=672
files= glob.glob(os.path.join(path,'*.tif'))
# always concatenate paths like this:
albedo = os.path.join(VRT, 'dummy_alb{}.vrt'.format(i))
outname = os.path.join(pathout, 'Albedo_L8_{}.tif'.format(files[i][39:46]))
with Raster(albedo) as ras:
    # ras.assign(<somearray>) # use this to assign a different array to the Raster object
    ras.write(outname=outname, format='GTiff', nodata=99999, dtype='float32')