don't save thousands of images in a directory

[RebeccaWPerry]

SEAS IT gave us the feedback that we should not be using our current method of organizing data series by using directories full of thousands of .tif images. Instead, we can use HDF5 which is an extremely general data storage framework.

I propose:

• Still save individual images that you capture one at a time as .tif's
• For slow and fast time series, save all images into a single HDF5 file
• For slow and fast time series, also save the first image of the time series as a .tif alongside the HDF5 as a thumbnail to show what is inside the HDF5.

Multipage tif's max out at 4GB unless you use bigtiff. I'm not sure how to write out a multipage tiff without having all the images in memory. All our images are in the frame grabber's memory, and then we like to transfer the images one at a time through the operating system's memory and to the hard disk.

Writing out 100 images to an HDF5 file could look like this:

import h5py

f = h5py.File("mytestfile.hdf5","w")

for i in range(100):
fakeim = around(10*rand(1024,1024))
dset = f.create_dataset(str(i),data=fakeim)

f.close()

[vnmanoharan]

I like the proposal and am coming around to HDF5. I would like us to think about how we can store the data inside the HDF5 in order to be able to compress efficiently. Our usual procedure is to store data as a 1D array of 2D frames, but it may be better for compression to store in a different way. For example, because each individual pixel doesn't change very often in a high-speed data set, we might get very efficient compression by storing a 2D array of 1D arrays, where each 1D array corresponds to the value of a particular pixel as a function of time.

Anyway, that's just one idea. Since we are moving to a new format, now might be the time to think about how we can make our data storage more efficient.

[tdimiduk]

The simplest way to compress will probably be to use the built in hdf5 compression:
http://docs.h5py.org/en/2.3/high/dataset.html

dset = f.create_dataset("zipped", (100, 100), compression="gzip")

we might be able to get what @vnmanoharan suggests transparently by how h5py chunks the data, for example chunking as:

dset = f.create_dataset("zipped", (1024, 1024, 1000), compression="gzip", chunks=(1, 1, 1000))

I haven't tested it yet, but it looks like that would instruct hdf5 to store each t stack of pixels as a chunk of data and compress it. There may be some performance hit to storing it like that, but hdf5 may be enough faster than what we were using before for it to still be fast

[tdimiduk]

Did some quick testing, and at least when we have all of the data already in memory (I just loaded 100 holograms from some old data), just doing h5py's default compression is both faster and smaller than telling it what to do with junks:
The data set is 100 images, (1024, 1024, 100)
Uncompressed: 801MB
r.create_dataset("data", data=data)

Default h5py compression: 94MB
g.create_dataset("data", data=data, compression='gzip')

Manually chunked compression: 165MB
f.create_dataset("zipped", data=data, compression='gzip', chunks=(1, 1, 100))

[tdimiduk]

I finished some basic benchmarks of h5py compression. See below for the full script.

raw: 16.6s, 419.0MB
gzipped: 44.9s, 78.0MB
tlines: 2927.9s, 155.0MB
tlines_x16: 451.5s, 61.0MB
tlines_x64: 451.9s, 62.0MB

Iterating through writing 1 image at a time seems to be a very poor performance case for the 2d array of 1d arrays that vinny suggested, I would guess it is recompressing every time it adds a new image.

import h5py
import glob
import os.path
import time
import holopy as hp
import tempfile

ims = [hp.load(f) for f in glob.glob('/group/manoharan2/tom/data/2014-09-22/10/*.tif')]

def setup_hdf5_file(shape, name, **options):
    f = h5py.File(name+'.hdf5', 'w')
    f.create_dataset('data', shape, **options)
    return f

def save_images(f):
    for i, im in enumerate(ims):
        f['data'][..., i] = im
    f.close()

def benchmark(name, **options):
    tstart = time.time()
    f = setup_hdf5_file((1024, 1024, 100), name, **options)
    save_images(f)
    tstop = time.time()
    size = os.path.getsize(name+'.hdf5')
    print("{0}: {1:.1f}s, {2}MB".format(name, tstop-tstart, size//1e6))

benchmark('raw')
benchmark('gzipped', compression='gzip')
benchmark('tlines', compression='gzip', chunks=(1, 1, 100))
benchmark('tlines_x16', compression='gzip', chunks=(16, 16, 100))
benchmark('tlines_x64', compression='gzip', chunks=(64, 64, 100))

[tdimiduk]

If instead we write images it blocks (say 100 images at a time shouldn't strain the computer's memory too badly), things perform much better:

raw: 0.3s, 838.0MB
autochunk: 8.8s, 98.0MB
chunk_16x16x100: 7.2s, 71.0MB

The compression is still a little slow for doing online, compared to the raw writes, but we can get a 10x size reduction, so not too shabby.

import h5py
import glob
import os.path
import time
import holopy as hp
import tempfile
import numpy as np

ims = [hp.load(f) for f in glob.glob('/group/manoharan2/tom/data/2014-09-22/10/*.tif')]
data = np.dstack(ims)

def print_stats(name, tstart, tstop):
    size = os.path.getsize(name+'.hdf5')
    print("{0}: {1:.1f}s, {2}MB".format(name, tstop-tstart, size//1e6))

def write_all_images(name, **options):
    f = h5py.File(name+'.hdf5', 'w')
    f.create_dataset('data', data=data, **options)

def benchmark(name, operation):
    tstart = time.time()
    operation(name)
    tstop = time.time()
    print_stats(name, tstart, tstop)

benchmark('raw', lambda x: write_all_images(x))
benchmark('autochunk', lambda x: write_all_images(x, compression='gzip'))
benchmark('chunk_16x16x100', lambda x: write_all_images(x, compression='gzip', chunks=(16,16,100)))

[vnmanoharan]

I don't quite understand what "chunks=(16,16,100)" means, but this seems like a great way of storing data efficiently.

One thing to try for the future (though not necessarily for resolving this issue) is to store differences of subsequent frames. If the data is high-speed and there is little noise, most values in the difference frames will be zero, so one can use compression routines designed for sparse matrices. For example, we could do run-length encoding.

Anyway I think we should implement Becca's proposed protocol using Tom's compression scheme by default. Also I think we should have a way to store the relevant metadata inside the hdf5, so that it is self-contained.

[tdimiduk]

Chunks=(16,16,100) means that on disk the data is broken into chunks which are 16 by 16 by 100 pixels (so a little rectangle of the image, across 100 frames). This allows gzip to find the compressibility along the time axis.

In this manner, I think we are effectively getting what @vnmanoharan suggests (gzip is essentially a very sophisticated run-length coding). I don't think storing actual numbers vs deltas matters for compression: 6, 6, 6, ... should compress just as well as 6, 0, 0, ... I would guess writing our own compressor could get us a noticeable improvement in write speed (since we could specialize it to our data), but will probably not get as good of compression ratios as gzip (since it has just had a lot more engineering put into it). Additionally, it would make the files less portable.

One way we could get around the write speed problem would be to write out uncompressed data when the user hits save, then spawn a thread which would compress the data into a new hdf5 file. Since we have multicore machines, this would let you get back to taking data quickly, and as long as you are not taking a ton of high frame rate data sets right in a row, the thread would hopefully finish before you are saving the next data.

HDF5 has an official standard for storing metadata:
http://docs.h5py.org/en/2.3/high/attr.html
attrs are a set of key value pairs associated with a data set, and h5py provides them as a dictalike in dataset.attrs. So we should put all of our metadata into those attrs when we write images, and then we don't need the seperate metadata yaml.

[RebeccaWPerry]

I can try implementing this tomorrow so that we can discuss around group
meeting if anything comes up. Just to be clear-- do you all want to save
the images as separate datasets inside the hdf5file or as a single 3D
dataset in the hdf5file? As separate images lets you have separate
metadata for each one I think, but there might be advantages to having a
3D dataset to slice through.

On 10/22/2014 07:14 PM, Tom Dimiduk wrote:

Chunks=(16,16,100) means that on disk the data is broken into chunks
which are 16 by 16 by 100 pixels (so a little rectangle of the image,
across 100 frames). This allows gzip to find the compressibility along
the time axis.

In this manner, I think we are effectively getting what @vnmanoharan
https://github.com/vnmanoharan suggests (gzip is essentially a very
sophisticated run-length coding). I don't think storing actual numbers
vs deltas matters for compression: 6, 6, 6, ... should compress just
as well as 6, 0, 0, ... I would guess writing our own compressor could
get us a noticeable improvement in write speed (since we could
specialize it to our data), but will probably not get as good of
compression ratios as gzip (since it has just had a lot more
engineering put into it). Additionally, it would make the files less
portable.

One way we could get around the write speed problem would be to write
out uncompressed data when the user hits save, then spawn a thread
which would compress the data into a new hdf5 file. Since we have
multicore machines, this would let you get back to taking data
quickly, and as long as you are not taking a ton of high frame rate
data sets right in a row, the thread would hopefully finish before you
are saving the next data.

HDF5 has an official standard for storing metadata:
http://docs.h5py.org/en/2.3/high/attr.html
attrs are a set of key value pairs associated with a data set, and
h5py provides them as a dictalike in dataset.attrs. So we should put
all of our metadata into those attrs when we write images, and then we
don't need the seperate metadata yaml.

—
Reply to this email directly or view it on GitHub
#13 (comment).

[tdimiduk]

I think it needs to be one 3d dataset to compress well. Since most of the compressability is along the t axis, splitting each frame into a different data set would hide that from the compression algorithm.

For fast timeseries, I think what you will want to do is create a 3d data set that is the size of the whole movie (so say 1024, 1024, 1000), and then write into it in chunks of say 100 images (see what plays well with memory on the computer). Not writing individual images one at a time will speed things up a lot, but by not trying to write the whole movie at once will stop it from having memory problems if you try save a movie with tons of frames.

For slow time series I think it will be fine to just write each frame as we get it, though we will need to check, since I think it may have to re-compress the whole chunks each time you add a frame. This is where @vnmanoharan suggestions of a more specialized encoder could help a lot. But lets see how things go with the simple solution.

[RebeccaWPerry]

Tom, your last benchmark set isn't quite realistic because you create the dataset with 100 images and then never append to it. But, if you want to load ten sets of 100 images into the same dataset, I think it gets slow. We might have to save 100 image datasets, or save each image as a separate set and then compress offline like you suggested.