Hyperspectral data is big!

Development trade-offs in hyperspectral software

Compared to remote sensing, hyperspectral core scanning produces fairly manageable file sizes.

For a single core box (SWIR data in this example, as that is what I mostly work with):

Dataset	Size	Data type	Format
Raw data	~400 MB	float32	ENVI
Cropped reflectance	~170 MB	float32	npy/ENVI
Derived data
Smoothed	~170 MB	float32	npy/ENVI
Normalised (continuum removed)	~300 MB	float64	npy/ENVI
Mask, bands, metadata, interpretation images	~10 MB	various	npy, ENVI, jpg, json, xml

But there are a lot of them!

A typical project is 10–12 holes of ~80 boxes each, so call it 1000 boxes for a small to medium project. At ~1 GB per box, that is ~1 TB, and the derived data accounts for almost half of it.

But it is derived data. Why not just derive it when you need it?

If you have used CoreSpecViewer, you know that I do store it. All of it, for every box. Here is why.

Three strategies

A system working with hyperspectral data has three resources to spend:

• Disk
• Time
• RAM

And there are three ways to handle derived products:

1. Derive on demand. Compute it when needed, use it, throw it away. Costs time.
2. Cache in RAM. Compute everything once at startup and hold it. Costs RAM.
3. Cache to disk. Compute once, store it, read it back. Costs disk.

Let’s go through them, worst first.

Cache in RAM

Many of my processes need both the smoothed and the normalised data. We could compute it for a whole hole at startup: a one-time slowdown, then fast iterative work and no wasted disk.

For an 80-box hole at ~300 MB per normalised dataset, that is 24 GB of RAM, before any algorithm allocates a single intermediate array, and before you meet a hole longer than average.

I do not see a >64 GB machine in my professional near future, never mind >128 GB. An 8 GB laptop is probably standard for the average geo-about-town. So this one is out.

Derive on demand

This is the option that deserves a proper look. It costs no disk, and its RAM footprint is small, because you only hold one box at a time:

Dataset	Size
Cropped reflectance	~170 MB
Smoothed	~170 MB
Normalised	~300 MB
Inputs held in RAM	~640 MB

Add the intermediates of a feature extraction and you are still under 1 GB. That is nothing, so RAM is not the problem here. Time is.

Derived product	Time
Smoothing	~1–2 s
Continuum removal	~10–30 s

Here is the call signature for feature extraction in CoreSpecViewer:

def Combined_MWL(savgol, savgol_cr, mask, bands, feature, technique='QUAD',
                 use_width=False, cached_arrays=None):

It takes savgol (the smoothed reflectance) and savgol_cr (the normalised reflectance). Mask and bands have to be stored regardless. So, starting from cropped reflectance:

Step	Derive on demand	Cached to disk
Smoothing	~1 s	–
Continuum removal	~10 s	–
Read from disk	–	~1 s
Feature extraction	~15 s	~15 s
Per box	~26 s	~16 s
Per 80-box hole	~35 min	~21 min
Per 10-hole project	~5.8 h	~3.5 h

In this straightforward path, caching saves about 40%, which is real, but a slow feature extraction is still slow.

You will never run a feature extraction once. You pick parameters, look at the result, change a threshold, try a different fitting method, look again. You flick between boxes in the viewer, and you build downhole plots. Every one of those is an access.

Derive on demand pays the derivation cost on every access. Caching pays it once, and every access after that costs a disk read.

Ten iterations of tuning across one hole means 800 box accesses. At ~11 s of re-derivation each, that is nearly 2.5 hours spent recomputing data that has not changed. The cached version spends about 13 minutes reading it back. In exploratory, iterative work, the work this software exists for, the number of accesses explodes and the gap explodes with it.

So, disk

Disk caching vastly increases the access speed for derived data, but what pushes it over the edge into the absolute preferred option?
Memmapping.

Because everything lives on disk as arrays, CoreSpecViewer can memory-map it. That lets me have a 1 km hole “open” (over 150 GB of data) on an 8 GB RAM machine. Arrays are paged in as they are read and dropped back to the memmap when I am done with them.

You cannot memmap an array that does not exist yet. A whole-hole view built on derive-on-demand would mean deriving the whole hole, which drops us straight back into the RAM problem.

Caching everything looked wasteful when I started, since it nearly doubles the storage. But of the three resources, disk is by far the most abundant and the cheapest. A terabyte costs less than an afternoon of a geologist’s time.

What looked extremely wasteful is actually spending the abundant resource to protect the scarce ones.

“Conclusions”…

This conclusion depends on my stack.

“You fool!”, I hear you cry, “you wrote it in Python! Of course it is slow!”

You would not be wrong. Derivation could be pushed into compiled kernels, rewritten in another language, made blazingly fast. If continuum removal took one second instead of ten, the access-count argument would shrink a lot and I would revisit this. Maybe the commercial software world has already solved it. But Python gives me too much (NumPy, SciPy, hylite) for me to move away from it, and with this stack, disk wins.

In a nod to the fact disk caching doubles the storage CoreSpecViewer has an archive file type for long-term storage. This is the minimum dataset needed to deterministically reproduce a fully analysed box. Cache everything while you are working, and archive the minimum when you are done. The rehydrate process is the time penalty up-front when you open an archive, then you are back to working iteratively on memmaps.

The Original Linkedin post had some good points made by commenters that are worth reading, especially if you made it all the way through this!