def subarray(rast, filtered_value=-9999, indices=False):
'''
Given a (p x m x n) raster (or array), returns a (p x z) subarray where
z is the number of cases (pixels) that do not contain the filtered value
(in any band, in the case of a multi-band image). Arguments:
rast The input gdal.Dataset or a NumPy array
filtered_value The value to remove from the raster array
indices If True, return a tuple: (indices, subarray)
'''
# Can accept either a gdal.Dataset or numpy.array instance
if not isinstance(rast, np.ndarray):
rastr = rast.ReadAsArray()
else:
rastr = rast.copy()
shp = rastr.shape
if len(shp) == 1:
# If already raveled
return rastr[rastr != filtered_value]
if len(shp) == 2 or shp[0] == 1:
# If a "single-band" image
arr = rastr.reshape(1, shp[-2]*shp[-1])
return arr[arr != filtered_value]
# For multi-band images
arr = rastr.reshape(shp[0], shp[1]*shp[2])
idx = (arr != filtered_value).any(axis=0)
if indices:
# Return the indices as well
rast_shp = (shp[-2], shp[-1])
return (np.indices(rast_shp)[:,idx.reshape(rast_shp)], arr[:,idx])
return arr[:,idx]
评论列表
文章目录
