Cómo crear máscara de bandas en GDAL y Python?

Question

He pequeño programa en Python que utiliza gdal y numpy para realizar una escala logarítmica de algunos datos en el RGB rango 0-255.

Ahora esta de entrada de datos tiene un valor nodata conjunto para enmascarar los valores válidos, y me gustaría tener la salida a conservar la máscara.

Obviamente no puedo usar el nodata para mi salida como soy de ussing el total de bytes rango de valores de color.

De lo que puedo decir de la GDAL API necesito para crear una máscara de la banda en el archivo de salida.

He añadido el código para hacer esto, pero de repente el programa empezó a funcionar muy lento y el tamaño del archivo de salida llegó a ser enorme.

¿Cuál es la mejor manera de hacer esto?

Aquí está mi código:

#!/usr/bin/env python

import gdal
from gdalconst import *
import numpy
import sys

infile = sys.argv[1]
outfile = sys.argv[2]
format = 'GTiff'
type = GDT_Byte

indataset = gdal.Open(infile, GA_ReadOnly)

gdal.SetConfigOption('GDAL_TIFF_INTERNAL_MASK', 'YES')
out_driver = gdal.GetDriverByName(format)
outdataset = out_driver.Create(outfile, indataset.RasterXSize, indataset.RasterYSize, indataset.RasterCount, type)
outdataset.CreateMaskBand(gdal.GMF_PER_DATASET)

gt = indataset.GetGeoTransform()
outdataset.SetGeoTransform(gt)

prj = indataset.GetProjectionRef()
outdataset.SetProjection(prj)

for iBand in range(1, indataset.RasterCount + 1):
    inband = indataset.GetRasterBand(iBand)
    inmaskband = inband.GetMaskBand()

    outband = outdataset.GetRasterBand(iBand)
    outmaskband = outband.GetMaskBand()

    for i in range(inband.YSize - 1, -1, -1):
        scanline = inband.ReadAsArray(0, i, inband.XSize, 1, inband.XSize, 1)
        scanline = numpy.multiply(numpy.divide(numpy.log1p(numpy.minimum(10000, numpy.maximum(0, scanline))), numpy.log1p(10000)), 255)
        outband.WriteArray(scanline, 0, i)

        scanline = inmaskband.ReadAsArray(0, i, inband.XSize, 1, inband.XSize, 1)
        outmaskband.WriteArray(scanline, 0, i)

Answer 1

para una banda en particular, he aquí cómo puedo crear un nodata máscara.

import gdal
import numpy as np

nodata = -9999#or set it based on the raster nodata value.
r = gdal.Open(rastername)
raster_arr = np.array(r.GetRasterBand(band).ReadAsArray())

nodatamask = raster_arr == nodata

#do your thing and end up with 
#a result_raster that you'd like to save

result_raster[nodatamask] = nodata
#then save.

Por supuesto, si usted está haciendo logaritmos y tiene valores negativos, sería necesario hacer un temporal valor nodata.

nodatamask = raster_arr == nodata

raster_arr[nodatamask] = 100000000

#do your thing

result_raster[nodatamask] = nodata
#then save.