It works perfectly well, when I rasterize all the features at once, into one 2d numpy array;
However, what I want now is to produce a "cube", having separate raster for each feature separately, with a corresponding value (say, number of category). Plan is then to perform aggregation - say, take moda, from this cube.
However, that takes a lot of time on my machine (~3 hours vs 30 sec for one raster map), so I wonder if there are any suggestions/hacks on how I can speed it up.