Plot histogram with rasterio

Question

This is not something I do often, but I recall it being fairly easy. Today however I can't make it work:

Python 3.8.5 (default, Jul 28 2020, 12:59:40) 
[GCC 9.3.0] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import rasterio
>>> from matplotlib import pyplot
>>> src = rasterio.open("ProtectedAreas_152160_mask_5x5.gtiff")
>>> from rasterio.plot import show_hist
>>> show_hist(src, bins=50, lw=0.0, stacked=False, alpha=0.3,histtype='stepfilled', title="Histogram")
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/lib/python3/dist-packages/rasterio/plot.py", line 250, in show_hist
    arr = source.read(masked=masked)
  File "rasterio/_io.pyx", line 336, in rasterio._io.DatasetReaderBase.read
MemoryError: Unable to allocate array with shape (1, 87425, 89057) and data type float64
>>> 
>>> show_hist(src, title="Histogram")
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/lib/python3/dist-packages/rasterio/plot.py", line 250, in show_hist
    arr = source.read(masked=masked)
  File "rasterio/_io.pyx", line 336, in rasterio._io.DatasetReaderBase.read
MemoryError: Unable to allocate array with shape (1, 87425, 89057) and data type float64

Is there something missing here? Some type conversion perhaps?

Answer 1

You are passing a rasterio.dataset object, i.e. src to the show_hist function which is fine.

But then you hit a MemoryError. This tells you that the array you are trying to load in memory (this is done under the hood by the show_hist function, see arr = source.read(masked=masked) in your traceback) is just too big in size to be allocated. You could think to split your raster in subsets and use the show_hist on those.

EDIT 1

One way I can think of working with the whole array, i.e. to avoid to split it into subsets would be to cast it into a different dtype that would take less memory. On my computer trying to initialize an array of your size (and your data type) would result in the same Memory error:

>>> arr = np.ones((1, 87425, 89057)) # By default dtype is numpy.float64
>>> MemoryError: Unable to allocate 58.0 GiB for an array with shape (1, 87425, 89057) and data type float64

Trying with a lower precision float data type, i.e. np.float32 would result in the same error:

>>> arr = np.ones((1, 87425, 89057), dtype=np.float32)
>>> MemoryError: Unable to allocate 29.0 GiB for an array with shape (1, 87425, 89057) and data type float32

As you can see the array that it is trying to allocate is now smaller in size (29 GB compare to the previous 58 GB) but still too big for my available memory.

Only by using an int8 dtype I am able to load the array in memory:

>>> arr = np.ones((1, 87425, 89057), dtype=np.int8)
>>> arr.shape
>>> (1, 87425, 89057)

But I guess there is a reason why you are using a float dtype.