Note
Click here to download the full example code
rasterMath with custom block size, mask, and in 3 dimensions¶
Tips to use rasterMath by defining its block size and to receive a full block (not a array with one pixel per row.)
Tips : A function readBlockPerBlock() yields each block, without saving results to a new raster.
Import librairies¶
from museotoolbox.processing import RasterMath,image_mask_from_vector
from museotoolbox import datasets
from matplotlib import pyplot as plt
import numpy as np
Load HistoricalMap dataset¶
raster,vector = datasets.load_historical_data()
Initialize rasterMath with raster¶
# Set return_3d to True to have full block size (not one pixel per row)
# Create raster mask to only keep pixel inside polygons.
image_mask_from_vector(vector,raster,'/tmp/mask.tif',invert=False)
rM = RasterMath(raster,in_image_mask='/tmp/mask.tif',return_3d=True)
#rM.addInputRaster('/tmp/mask.tif')
print(rM.get_random_block().shape)
Out:
Total number of blocks : 15
(256, 256, 3)
Plot blocks
x = rM.get_random_block()
rM.add_function(np.mean,'/tmp/mean.tif',axis=2,out_np_dt=np.int16)
rM.run()
from osgeo import gdal
dst = gdal.Open('/tmp/mean.tif')
arr = dst.GetRasterBand(1).ReadAsArray()
plt.imshow(np.ma.masked_where(arr == rM._outputs[0]['nodata'], arr))
Out:
Detected 1 band for function mean.
No data is set to : -32768.
Batch processing (15 blocks using 11Mo of ram)
rasterMath... [##......................................]6%
rasterMath... [#####...................................]13%
rasterMath... [########................................]20%
rasterMath... [##########..............................]26%
rasterMath... [#############...........................]33%
rasterMath... [################........................]40%
rasterMath... [##################......................]46%
rasterMath... [#####################...................]53%
rasterMath... [########################................]60%
rasterMath... [##########################..............]66%
rasterMath... [#############################...........]73%
rasterMath... [################################........]80%
rasterMath... [##################################......]86%
rasterMath... [#####################################...]93%
rasterMath... [########################################]100%
<matplotlib.image.AxesImage object at 0x7f7fdf87d690>
Total running time of the script: ( 0 minutes 0.250 seconds)