[Numpy-discussion] summarizing blocks of an array using a moving window
Robin Kraft
rkraft4 at gmail.com
Thu Jul 22 13:30:27 EDT 2010
Vincent, Pauli,
> From: Vincent Schut <schut at sarvision.nl>
> - an other option would be some smart reshaping, which finally gives you
> a [y//2, x//2, 2, 2] array, which you could then reduce to calculate
> stats (mean, std, etc) on the last two axes. I *think* you'd have to
> first reshape both x and y axes, and then reposition the axes. An
> example:
> a = gdal_array.BandReadAsArray(blabla)
> y,x = a.shape #y and x need be divideable by 2!
> b = a.reshape(y/2, 2, x/2, x).transpose(0,2,1,3).reshape(y/2, x/2, 4)
> bMean = b.mean(axis=-1)
> bMax = ......etc.
You and Pauli agree on this - seems like a good option.
> - a third option would be to create an index array, which has a unique
> value per 2x2 square, and then use histogram2d. This would, if you use
> its 'weight' functionality, at least enable you to get efficient counts
> and sums/means. Other stats might be hard, though.
Hmmmm I don't get this, but I'll experiment. I've seen some really useful things done with histogram2d so it seems worthwhile to figure out.
> Message: 6
> From: Pauli Virtanen <pav at iki.fi>
> > Let's say the image looks like this: np.random.randint(0,2,
> > 16).reshape(4,4)
> >
> > array([[0, 0, 0, 1],
> > [0, 0, 1, 1],
> > [1, 1, 0, 0],
> > [0, 0, 0, 0]])
> >
> > I want to use a square, non-overlapping moving window for resampling, so
> > that I get a count of all the 1's in each 2x2 window.
> >
> > 0, 0, 0, 1
> > 0, 0, 1, 1 0 3
> > => 2 0
> > 1, 1, 0, 0
> > 0, 0, 0, 0
> >
> > In another situation with similar data I'll need the average, or the
> > maximum value, etc..
>
> Non-overlapping windows can be done by reshaping:
>
> x = np.array([[0, 0, 0, 1, 1, 1],
> [0, 0, 1, 1, 0, 0],
> [1, 1, 0, 0, 1, 1],
> [0, 0, 0, 0, 1, 1],
> [1, 0, 1, 0, 1, 1],
> [0, 0, 1, 0, 0, 0]])
>
> y = x.reshape(3,2,3,2)
> y2 = y.sum(axis=3).sum(axis=1)
This is perfect - and so fast! Thanks! Now I just have to understand why it works ...
Can anyone recommend a tutorial on working with (slicing, reshaping, etc.) multi-dimensional arrays? Stefan's slides are beautiful but my brain starts to hurt if I try to follow them line by line.
-Robin
More information about the NumPy-Discussion
mailing list