1

I am using following code to create bt2.vrt file using python.

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

# Function to read the original file's projection:

def GetGeoInfo(FileName):
    SourceDS = gdal.Open(FileName, GA_ReadOnly)
    xsize = SourceDS.RasterXSize
    ysize = SourceDS.RasterYSize
    GeoT = SourceDS.GetGeoTransform()
    Projection = osr.SpatialReference()
    Projection.ImportFromWkt(SourceDS.GetProjectionRef())
    DataType = SourceDS.GetRasterBand(1).DataType
    #DataType = gdal.GetDataTypeName(DataType)
    return xsize, ysize, GeoT, Projection, DataType

# Function to write a new file.


def Create_vrt(Name, Array, driver,xsize, ysize, GeoT, Projection, DataType):
    if DataType == 'UInt16':
        DataType = gdal.GDT_Float32
    NewFileName = 'bt2'+'.vrt'
    # Set nans to the original No Data Value
    # Set up the dataset
    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

# Open the original file
FileName='/home/vigna/TIFFs/bt.tif'    

DataSet = gdal.Open(FileName, GA_ReadOnly)
# Get the first (and only) band.
Band = DataSet.GetRasterBand(1)
# Open as an array.
Array = Band.ReadAsArray()

# Now I do some processing on Array, it's pretty complex 
NewArray = Array
print(NewArray)

#save the new file, in the meantime I have 
xsize, ysize, GeoT, Projection, DataType = GetGeoInfo(FileName)
# Set up the vrt driver
driver = gdal.GetDriverByName('VRT')

# Now turn the array into a vrt.
NewFileName = Create_vrt('bt3', NewArray, driver,  
                        xsize, ysize, GeoT, Projection, DataType)
#Create_vrt()
 print("ok")

bt2.vrt file created using above code is : [1]: https://i.stack.imgur.com/JYqVf.png

According to [1]: Unable to warp HDF5 files I need to include lat.vrt and long.vrt as y_dataset and X_dataset respectively to perform "gdal_warp" on bt.vrt file. How to include it using python code?

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.