Access the array of a loaded layer in QGIS3 with pyQGIS

Question

The code below is a clunky attempt to get a raster numpy array loaded up for basic machine learning from the QGIS Python console. It usually works but has stability issues with the GetRasterBand() and I have asked a separate question. As an improvement I am trying to figure how to get the array values directly from the loaded rasters in QGIS. I.e. you open the raster layers in QGIS, do whatever is needed (e.g. training polygons) and then run the code.

#Part 1. Importing Python modules
import numpy as np
import gdal
from sklearn.cluster import KMeans
import os



#Part 2. Specify the image file to be classified and the class raster with the training data.

MyFolder = 'D:\\RemoteSensing\\Practical_7\\' #Change this to match the folder on your drive.  Note that we use \\ instead of \
MyImageRoot = 'S2_mini' #Specify the root image name. Do not nclude the extension, e.g. 'glacier_'.  If a color image, just write the file name without extension.
Nclasses = 4 #how many classes are in the image
MyClassRaster = 'KmeansTEST_4cls'
BandList = ['B02'] #Indicate the bands you want to load.  B


#Part 3 Load the image 
fileName = os.path.join(MyFolder + MyImageRoot + BandList[0] + '.tif')
ImageFile = gdal.Open(fileName)
band1 = ImageFile.GetRasterBand(1)
array1 = band1.ReadAsArray()

I am trying to get array1 by direct querry to a loaded band in the QGIS layers. I have found QgisProject:

layer = QgsProject.instance().mapLayersByName(BandList[0])

But I can't find how to access the raster values for the layer and have them as a numpy array for the next processing steps. I'm working with QGIS 3.2 and 3.4.

Answer 1

For accessing raster values for the layer in PyQGIS you can use QgsRasterBlock objects. It is exemplified in the following code for a raster named 'aleatorio':

layer = QgsProject.instance().mapLayersByName('aleatorio')

provider = layer[0].dataProvider()

extent = provider.extent()

rows = layer[0].height()
cols = layer[0].width()

# 1 is referring to band 1
block = provider.block(1, extent, cols, rows)

#accessing one particular value; e.g. value[4,3]
print (block.value(4,3))

After running above code, its value for (4,3) indexes is printed at Python Console (equal to 64). With Value Tool plugin it can be corroborated (following image).