flatten.py

#!/usr/bin/env python

import rasterio as rio
from sys import argv
from skimage import io
import numpy as np

# Reads an image with a bunch of bands (i.e., layers on the z axis)
# and returns an image with one band, each value in which is the mean
# along the z axis at that x,y.

with rio.open(argv[1]) as src:
  lidar = src.read()

flat = np.mean(lidar, axis=0)

io.imsave('flat.tif', flat)