I'm iterating through thousands of polygons. Each polygon clips a tif raster producing a numpy array. Then the numpy array values are manipulated and the tif raster is updated with those new values. The script works for a few thousand polygons but at some point it throws an error.
import os, sys, datetime, time
import geopandas as gpd
import rasterio
from rasterio.mask import mask
from fiona.crs import from_epsg
import numpy as np
import pycrs
from rasterio.windows import Window
from shutil import copy
import petl as etl
import shapely.wkt
alr_path = os.path.join(working_folder, 'INPUT', 'NI_COMBINED_ALR.tif')
alr = rasterio.open(alr_path)
out_meta = alr.meta.copy()
alr_avg_path = os.path.join(working_folder, 'OUTPUT', 'alr_avg.tif')
copy(alr_path, alr_avg_path)
feature_gdf = gpd.read_file(gdb_polygs, driver='FileGDB', layer=fc_vml_polygs)
for index, row in feature_gdf.iterrows():
window, out_transform = mask(alr, row, all_touched=True, crop=True)
x = alr.meta['transform'][2]
y = alr.meta['transform'][5]
width = window.shape[2]
height = window.shape[1]
col_off = out_transform[2]
row_off = out_transform[5]
window_orig = alr.read(1, window=Window(int((col_off - x)/5), int((y - row_off)/5), width, height))
if np.all([window < 0]):
avg = -1
else:
avg = window[window != -1].mean().item()
window_avg = (np.where(window!=-1, int(round(avg,0)), window_orig))
with rasterio.open(alr_avg_path, "r+", **out_meta) as dest:
dest.write(window_avg, window=Window(int((col_off - x)/5), int((y - row_off)/5), width, height))
This is the error that the last line of the code above throws:
ERROR 1: LZWDecode:Corrupted LZW table at scanline 10808
ERROR 1: TIFFReadEncodedStrip() failed.
ERROR 1: .\OUTPUT\alr_avg.tif, band 1: IReadBlock failed at X offset 0, Y offset 10808: TIFFReadEncodedStrip() failed.
Traceback (most recent call last):
File "footprint.py", line 277, in <module>
process_raster(con, pg_table_name, path_tif)
File "footprint.py", line 206, in process_raster
dest.write(numpy_array, window=Window(x, y, width, height))
File "rasterio\_io.pyx", line 1409, in rasterio._io.DatasetWriterBase.write
OSError: Read or write failed
Any ideas about what this error means and how it could be overcome?