Python GDAL: Write new raster using projection from old

http://gis.stackexchange.com/questions/57005/python-gdal-write-new-raster-using-projection-from-old
This is the routine I developed to convertISIS3 cubes to GTiffs. I expect a similar approach should work between anytypes of drivers (though I think the driver.Create() method might limit thechoice of output file).
 
import numpy as np
import gdal
from gdalconst import *
from osgeo import osr
 
# Function to read the original file's projection:
def GetGeoInfo(FileName):
    SourceDS =gdal.Open(FileName, GA_ReadOnly)
    NDV =SourceDS.GetRasterBand(1).GetNoDataValue()
    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 NDV,xsize, ysize, GeoT, Projection, DataType
 
# Function to write a new file.
def CreateGeoTiff(Name, Array, driver, NDV, 
                 xsize, ysize, GeoT, Projection, DataType):
    if DataType== 'Float32':
       DataType = gdal.GDT_Float32
    NewFileName= Name+'.tif'
    # Set nansto the original No Data Value
   Array[np.isnan(Array)] = NDV
    # Set upthe dataset
    DataSet =driver.Create( NewFileName, xsize, ysize, 1, DataType )
            #the '1' is for band 1.
   DataSet.SetGeoTransform(GeoT)
    DataSet.SetProjection(Projection.ExportToWkt() )
    # Write thearray
   DataSet.GetRasterBand(1).WriteArray( Array )
   DataSet.GetRasterBand(1).SetNoDataValue(NDV)
    returnNewFileName
 
# Open the original file
FileName = 'I29955002trim.cub'    # This is the ISIS3 cube file
                                  # It's aninfra-red photograph
                                  # taken bythe 2001 Mars Odyssey orbiter.
DataSet = gdal.Open(FileName, GA_ReadOnly)
# Get the first (and only) band.
Band = DataSet.GetRasterBand(1)
# Open as an array.
Array = Band.ReadAsArray()
# Get the No Data Value
NDV = Band.GetNoDataValue()
# Convert No Data Points to nans
Array[Array == NDV] = np.nan
 
# Now I do some processing on Array, it's prettycomplex 
# but for this example I'll just add 20 to eachpixel.
NewArray = Array + 20 # If only it were that easy
 
# Now I'm ready to save the new file, in the meantimeI have 
# closed the original, so I reopen it to get theprojection
# information...
NDV, xsize, ysize, GeoT, Projection, DataType =GetGeoInfo(FileName)
 
# Set up the GTiff driver
driver = gdal.GetDriverByName('GTiff')
 
# Now turn the array into a GTiff.
NewFileName = CreateGeoTiff('I29955002trim',NewArray, driver, NDV, 
                            xsize, ysize, GeoT, Projection, DataType)
 
And that's it. I can open both images inQGIS. And gdalinfo on either file shows that I have the same projectioninformation and georeferencing.