Keep track of density by hash rather than in an array

[geowurster]

Run each point through the affine transform and generate a dictionary of indices and counts. Stitch everything together at the end. Will likely be much faster with a lower memory footprint since most of the "pixels" are empty.

density = defaultdict(int)

pos = (x, y) * ~affine
density[pos] += 1

out = np.array((rows, cols))

for pos, count in density.items():
    r, c = pos
    out[r][c] += count