Limits of pylidar spatial processing mode

Original report by Neil Flood (Bitbucket: neilflood, GitHub: neilflood).

(Note - also posted this to the discussion group, but I don't know if anyone sees that)

Now that few users have gained some experience in using pylidar for real tasks (it has taken a while), it is starting to become apparent that there are some real limitations with the spatial processing mode of pylidar. This is when each block of data is constructed as an array of binned data, with the first two dimensions representing X and Y position, and the third dimension being a ragged array giving each point which falls within that bin. The ragged array is implemented using a numpy masked array, so that for each bin, the third dimension only extends to as many points actually occur in that bin, and the remaining elements are masked.

The main problems arise with the time taken to construct this array, and the size it can occupy in memory. It would seem that these are both tied to the fact the point density can be quite variable across a scan. This means that for a particular block, there can be a single bin which has many points, and this is used as the maximum size of the third dimension. However, if this number is considerably larger than the average across the block, then it makes the whole array very large in that direction, resulting in considerable "empty" space in the rest of the array. Since the block size (in terms of X and Y) is fixed across a whole call to doProcessing(), it must be small enough that the block with the largest point count in a single bin can be accommodated in memory, however, for some scans (TLS in particular), this can impose a very small block size, in order to fit within available memory.

I am not sure, but I think this also ties in to the run time, as the construction of this array takes considerable time, and with a small block size, I think that much of the data would have to be handled multiple times to select which points apply to a given block. Not sure about that part.

In any case, on TLS scans (where point density is highest in the centre), the memory requirements can go into many Gb, for files which are not that big in the first place. The run times can be many hours to process across such a file.

Jasmine and I put together a non-spatial processing framework, in which the non-spatial mode is used, and the user function does what it requires in spatial terms, just on-the-fly. In this case, the operation was a simple min(Z) operation, and was implemented using a numba jit function to do the per-point section. The memory requirement was tiny, and the run time went from on the order of 12 hours down to around 5 minutes.

Obviously this is a vast improvement.

My feeling is that the spatial processing mode, with the masked arrays, was worth testing out, but I suspect it is not a feasible model for doing this sort of thing, and we are better off using the non-spatial mode and handling the spatial aspects in the particular case, for a given task. Using numba will help greatly in this.

With this in mind, I would like to suggest that we change the default control setting to non-spatial, and shelve all further plans for development on the spatial-mode parts of the code.

It also raises the question (which Nick has raised before) about what purpose the indexing of a file is actually serving for us. The TLS example mentioned above made no use of the index, and runs vastly faster, so it is not clear that indexing files is a useful thing to do.

These are just my impressions, and I recognise that there could be a lot else I don't know about. As always, I emphasize that I am just helping out here, and don't really know anything about lidar, but these results suggest that there are some important questions to consider.

Any thoughts are welcome.

Neil

ubarsc / pylidar

Limits of pylidar spatial processing mode #11