Perform band calc in numpy array

Question

I'd like to perform a band calc in some raster bands. What I am trying to do is import the raster (ENVI format) like array with gdal, then I do the calc but when I try to save the raster I get an error:

'NoneType' object has no attribute 'SetMetadata'

This is my code that I take from http://geoexamples.blogspot.com.es/2012/12/raster-calculations-with-gdal-and-numpy.html

from osgeo import gdal  
    from osgeo.gdalnumeric import *  
    from osgeo.gdalconst import *
    import numpy as np

    slope = 1.35
    intercept = -26.82

    fileName = r'mypath\b1.img'
    outFile = r'mypath\b1_out.img' 

    #Open the dataset  
    ds1 = gdal.Open(fileName, GA_ReadOnly)
    band1 = ds1.GetRasterBand(1)
    data1 = BandReadAsArray(band1)

    #The actual calculation  
    dataOut = data1*slope + intercept  

    #Write the out file  
    driver = gdal.GetDriverByName("ENVI")  
    dsOut = driver.Create(outFile, ds1.RasterXSize, ds1.RasterYSize, 1, band1.DataType)  
    CopyDatasetInfo(ds1,dsOut)  
    bandOut=dsOut.GetRasterBand(1)  
    BandWriteArray(bandOut, dataOut) 

Antoher big doubt is that the original data are bytes (0-255) but when I perform the calc I get data values biggers and lowers than 0 and 255. So I am not very sure if I could save the result like byte of should I save it in ther format and then run a gdal_translate ?

Answer 1

Presumably one of the parameters to CopyDatasetInfo is None, so verify that dsOut is not None.

You could save your results as Float32, or you could scale them back to 0-255 to save as Byte, or scale them to some other integer range. It depends on what you want to do with the result, and how much information you can afford to lose.

Answer 2

Eventuall I use this code from rasterio examples:

import numpy as np
import rasterio
import subprocess


    # Register GDAL format drivers and configuration options with a
    # context manager.
    with rasterio.drivers():

    # Read raster bands directly to Numpy arrays.
    #
    with rasterio.open(r'mypath\band1.img') as src:
        rs = src.read()

        NoData_rs = np.ma.masked_where(rs==255,rs) 
        total = NoData_rs*slope + intercept

        min_msk = (total<0)
        max_msk = (total>255)
        total[min_msk] = 0
        total[max_msk] = 254


        profile = src.profile
        profile.update(type=rasterio.uint8, NoData = 255)


        with rasterio.open(r'my_path\ex_nd.img', 'w', **profile) as dst:
            dst.write(total.astype(rasterio.uint8))

I find this a kind of weird, because I thought that I should use 'gdal_translate' to get my raster in byte again. When I perform the calc, I get values above and below 0 and 255, and if I runned this script without the min and max masks, It send these values to NoData (255) which isn't correct. Anyway, now with this scripts ans min/max masks, I get the same values that doing the gdal_translate step and that the ones that I got it running the normalization in Idrisi (which is my goal here, to do the same that Idrisi Regress in python). So I believe that I am done with this!

Thanks for your answer