Funciones con PyQGIS: archivos ráster

En el post anterior comenzamos a considerar la sistematización de todo el código de prueba que hasta ahora se ha realizado con PyQGIS y que, por tanto, pueda ser reutilizado fácilmente. Comenzamos con archivos vectoriales y ahora le tocará el turno a los ráster. El post sobre el cual se basa la sistematización es:

Estadísticas de un raster “Single Band” usando PyQGIS

e incluye una serie de procedimientos para obtener las estadísticas al partir del provider de un ráster ‘singlebandgray’.

El código que adapté a continuación como la función statsRaster(), si no existe capa activa o ésta es de tipo vectorial lo señala con una salida. Por otra parte, cuando la capa es de tipo ráster, la salida incluye el número de bandas, el tipo de renderizado del ráster, las unidades del mapa (metros, grados, etc), la resolución (X, Y) y finalmente, la media, la desviación estándar, el valor mínimo, el valor máximo y el rango de valores para cada banda usada (determinadas con renderer.usesBands).

#Funciones para usar con archivos ráster

def unitsTypes():
    
    UnitType = {0 : 'Meters ', 
                1 : 'Feet ', 
                2 : 'Degrees', 
                3 : 'UnknownUnit', 
                4 : 'DecimalDegrees', 
                5 : 'DegreesMinutesSeconds', 
                6 : 'DegreesDecimalMinutes', 
                7 : 'NauticalMiles'}
        
    return UnitType
    
def statsRaster():
    
    layer = iface.activeLayer()
    
    if layer is not None:
        if layer.type() != layer.RasterLayer:
            print "La capa no es de tipo raster"
            return
    else:
        print "No hay active layer"
        return
        
    number_bands = layer.bandCount()
       
    resol_x = layer.rasterUnitsPerPixelX()
    resol_y = layer.rasterUnitsPerPixelY()
    
    layer_extent = layer.extent()
    provider = layer.dataProvider()
    renderer = layer.renderer()
    layer_crs = layer.crs() 
    
    units_types = unitsTypes() 
    map_units = layer_crs.mapUnits()
    uses_band = renderer.usesBands()
    
    print "El raster tiene " + str(number_bands) + " banda(s)"
    print "Es un raster de tipo " + renderer.type()
    print "Las unidades del mapa son " + units_types[map_units]
    
    print "rasterUnitsPerPixelX = ", resol_x
    print "rasterUnitsPerPixelY = ", resol_y
    
    print
    
    new_stats = QgsRasterBandStats()

    for band in uses_band:
        
        if provider.hasStatistics(band, QgsRasterBandStats.All) is not True:
            provider.initStatistics(new_stats,
                                    band, 
                                    QgsRasterBandStats.All,
                                    layer_extent,
                                    0)
        
        stats = provider.bandStatistics(band, 
                                        QgsRasterBandStats.All, 
                                        layer_extent, 
                                        0)
                                                     
        print "media banda" + str(band) + ": " + str(stats.mean)
        print "stdv banda" + str(band) + ": " + str(stats.stdDev)
        print "Min Value banda" + str(band) + ": " + str(stats.minimumValue)
        print "Max Value banda" + str(band) + ": " + str(stats.maximumValue)
        print "Range banda" + str(band) + ": " + str(stats.range)
        print
    

A continuación, tenemos la salida respectiva en la Python Console para un ráster ‘multibancolor’ (red, green, blue y transparency) y un ráster ‘singlebandgray’.

>>>statsRaster()
El raster tiene 4 banda(s)
Es un raster de tipo multibandcolor
Las unidades del mapa son Meters 
rasterUnitsPerPixelX =  15.0
rasterUnitsPerPixelY =  15.0

media banda1: 246.901782243
stdv banda1: 29.0155763882
Min Value banda1: 0.0
Max Value banda1: 255.0
Range banda1: 255.0

media banda2: 91.1657625487
stdv banda2: 62.6236229709
Min Value banda2: 0.0
Max Value banda2: 255.0
Range banda2: 255.0

media banda3: 7.99217407492
stdv banda3: 28.7666416747
Min Value banda3: 0.0
Max Value banda3: 255.0
Range banda3: 255.0

>>>statsRaster()
El raster tiene 1 banda(s)
Es un raster de tipo singlebandgray
Las unidades del mapa son Meters 
rasterUnitsPerPixelX =  73.9887232281
rasterUnitsPerPixelY =  73.9887232281

media banda1: 1824.71800614
stdv banda1: 373.273662546
Min Value banda1: 1323.0
Max Value banda1: 3556.0
Range banda1: 2233.0

En la imagen siguiente tenemos la visualización de los resultados cuando se ejecuta statsRaster() en la Python Console y parte del código (utils.py) en la ventana del editor.

raster

Esta entrada fue publicada en PyQGIS, SIG, Software Libre. Guarda el enlace permanente.

Una respuesta a Funciones con PyQGIS: archivos ráster

  1. Anónimo dijo:

    “Hay un libro abierto siempre para todos los ojos: La Naturaleza”.

    Flor Rocío Espinosa Jiménez
    Lic. Ciencias de la Tierra

Responder

Por favor, inicia sesión con uno de estos métodos para publicar tu comentario:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s