(Software recommendations?) Converting a 32bit ome.tiff to an 8bit ome.tiff

Hi all!

Overall question:

I will have a set of fifty 32bit, ~2-10GB pyramidal whole slide images in ome.tiff format that I would like to convert to 8bit pyramidal whole slide images (including multiplying the original 32 bit values by 100). Does anyone have a software recommendation to accomplish this?

Specifics:

Each image (unfortunately cannot be shared, but I can probably mock up something from several open source images) was written out using QuPath from a set of 8 deconvolved images. The final images have 9 channels and are generally around 7GB. The default bit depth there is 32, as it handles negative numbers and all of the decimal places needed to cover the optical densities from the deconvolved channels. I checked, and there does not seem to currently be a way to export the image “with some math.”

Goal
Pyramidal ome.tiff, 32bit, 9channel, values ~0-2.5 -> ?? 8bit, 9channel, values 0-255 (truncated on either end)

All values in all channels tend to be in the range of 0.01 to 2.00. There are some values below 0 and above 2.55, but I am fine with clipping those as they are not generally biologically relevant (black smudges, or below significant thresholds). So my hope was to multiply by 100 and then save as a 0-255 8bit image (per channel). And hopefully this will reduce the file size a bit.

My first thought was FIJI. The one image I have tested did open in FIJI (after downsampling they are only 22,000 by 25,000 or so pixels), but it took a very, very long time to open. I am not sure how long as I ended up going to bed and it was finished in the morning. So a macro there might be the only way to go if I can get everything written out well. The image math is easy, but I am not sure how to convert to 8bit while clipping values above and below 0 and 255. I think the default 8bit command interpolates based on max and min values, which would be bad for consistency.

Does anyone have any other suggestions? I looked at the bfconvert commands but was not able to find anything about performing pixel math or changing bit depth. I am also considered Imaris, but the image math there would be manual per image, and the output would be a .ims file (which would be an additional step to convert back to something like a pyramidal ome.tiff).

Anyway, I will cut off my rambling there, but please ask if you have any more questions or I can specify anything else useful about my inputs or desired outputs.

*None of this is for quantification, as all measurements were made from the brightfield images. I am not extremely worried about the exact accuracy of the pixel values as long as they are generally correct. This is purely for visualization of the different stains in combination with cluster analysis/gating.

Thanks,
Mike

2 Likes

Hmm, come to think of it, if I do stick with FIJI, I can probably use something like this for the moment.

run("Multiply...", "value=100 stack");
run("Min...", "value=0 stack");
run("Max...", "value=255 stack");
run("8-bit");

The primary issue is getting the image open in the first place! :slight_smile:

I didn’t think it was possible in QuPath yet, although it kind of is:

/**
 * Transform an ImageServer to 8-bit in QuPath, with scaling.
 *
 * Note that this *kind of* works (I think) in v0.2... but is quite hack-y,
 * and may not be entirely reliable (i.e. I've seen some errors).
 *
 * @author Pete Bankhead
 */

import qupath.lib.images.servers.PixelType
import qupath.opencv.ops.ImageOps
import qupath.lib.images.servers.ImageServers
import qupath.lib.images.writers.ome.OMEPyramidWriter

import static qupath.lib.gui.scripting.QPEx.*

def imageData = getCurrentImageData()

// Output downsample (may be 1 for full-resolution)
double outputDownsample = 4

// Output server path
def path = buildFilePath(PROJECT_BASE_DIR, 'converted-8bit.ome.tif')

// Create an op to handle all transforms
def op = ImageOps.buildImageDataOp()
    .appendOps(
            ImageOps.Core.multiply(100),
            ImageOps.Core.ensureType(PixelType.UINT8)
    )

// Create a scaling & bit-depth-clipping server
def server = ImageOps.buildServer(
        imageData,
        op,
        imageData.getServer().getPixelCalibration().createScaledInstance(outputDownsample, outputDownsample))

// Ensure the server is pyramidal
List<Double> downsamples = [outputDownsample]
double d = outputDownsample * 4
while (Math.min(server.getWidth()/d, server.getHeight()/d) >= 512) {
    downsamples << d
    d *= 4
}
server = ImageServers.pyramidalize(server, downsamples as double[])

// Write the pyramid
new OMEPyramidWriter.Builder(server)
        .parallelize()
        .downsamples(downsamples as double[])
        .bigTiff()
        .channelsInterleaved()
        .build()
        .writePyramid(path)

My attempt at running this was neither a complete success nor a complete failure.

Note that the ImageOps server won’t be pyramidal, therefore the pyramid is effectively having to be recalculated from the full resolution. This is not exactly very efficient.

1 Like

libvips recently added support for reading and writing OME-TIFF pyramids. I would expect (but have not tried) that it can efficiently convert your files.

Also tifffile might work. This is an example for conversion in RAM without compresion or transferring metadata:

# convert a multi-channel, tiled, pyramidal 32-bit OME-TIFF to 8-bit
from tifffile import TiffFile, TiffWriter

with TiffFile('pyramidal_32bit.ome.tif') as tif:
    with TiffWriter('pyramidal_8bit.ome.tif', bigtiff=True) as tif8:
        for series in tif.series:
            for i, level in enumerate(series.levels):
                image = level.asarray()
                image *= 100
                numpy.clip(image, 0, 255, out=image)
                image = image.astype(numpy.uint8)
                tile = level.keyframe.tilelength, level.keyframe.tilewidth
                photometric = level.keyframe.photometric
                subifds = len(series.levels) - 1 if i == 0 else None
                tif8.save(
                    image, subifds=subifds, tile=tile, photometric=photometric
                )
                del image
1 Like

As long as I can manage adjust it so that I can “Run for project…” and walk away for the weekend, it would suddenly be very efficient for me indeed :slight_smile:

@petebankhead
Do you know offhand what happens to negative and 256+ values at the ensureType step?

I’m gonna say clipped? And also that I didn’t check that to make entirely sure, just in case I’m wrong. It’s ultimately calling OpenCV methods, which I believe clip (unlike numpy).

I should admit that when I used it, I once saw a weird mixture of resolution in a single plane of an image (i.e. unexpected blockiness). I shall hope that was just a rare unlucky event.

If tifffile can do the job you need, I suspect it could be more efficient. Ideally, it would be nice to be able to convert each layer and each tile as it goes. This QuPath method won’t do that (upon reflection, a bit more substantial scripting that subclasses AbstractTileableImageServer probably could…).

A QuPath RenderedImageServer could be an alternative if you’re happy with an RGB output. But if you just want a scaled version of the original multichannel image then it wouldn’t be enough.

1 Like

Yep, 9 channels would be a tad too many to visualize easily in RGB, there is a lot of overlap. I’ll give this a shot first though, as it would be nice to keep the whole workflow in QuPath. If it terribly slow, I’ll explore some of the other options.

And it is very nice to have options!

Thank you both :slight_smile:

@petebankhead I ran into a problem trying to run the script:
I added two print statements to describe the inputs to the line that is causing the error (shown below), edited the file name, and changed the downsample to 1 (the image had already been downsampled at this point), but otherwise left the script the same. The line causing the error is:

server = ImageServers.pyramidalize(server, downsamples)

INFO: server ---ImageOp server: ImageOp server: ImageData: Fluorescence, 6454-02B_IF1.ome.tif {"op":{"type":"op.core.sequential","ops":[{"type":"op.core.multiply","values":[100.0]},{"type":"op.core.convert","pixelType":"UINT8"}]}} (6454-02B_IF1.ome.tif)
INFO: downsamples ---[1.0, 4.0, 16.0]
ERROR: MissingMethodException at line 47: No signature of method: static qupath.lib.images.servers.ImageServers.pyramidalize() is applicable for argument types: (qupath.opencv.ops.ImageOpServer, ArrayList) values: [ImageOp server: ImageOp server: ImageData: Fluorescence, 6454-02B_IF1.ome.tif {"op":{"type":"op.core.sequential","ops":[{"type":"op.core.multiply","values":[100.0]},{"type":"op.core.convert","pixelType":"UINT8"}]}} (6454-02B_IF1.ome.tif), ...]
Possible solutions: pyramidalize(qupath.lib.images.servers.ImageServer, [D)

ERROR
ERROR: groovy.lang.MetaClassImpl.invokeStaticMissingMethod(MetaClassImpl.java:1568)
    groovy.lang.MetaClassImpl.invokeStaticMethod(MetaClassImpl.java:1554)
    org.codehaus.groovy.runtime.callsite.StaticMetaClassSite.call(StaticMetaClassSite.java:50)
    org.codehaus.groovy.runtime.callsite.CallSiteArray.defaultCall(CallSiteArray.java:47)
    org.codehaus.groovy.runtime.callsite.AbstractCallSite.call(AbstractCallSite.java:125)
    org.codehaus.groovy.runtime.callsite.AbstractCallSite.call(AbstractCallSite.java:148)
    Script4.run(Script4.groovy:48)
    org.codehaus.groovy.jsr223.GroovyScriptEngineImpl.eval(GroovyScriptEngineImpl.java:317)
    org.codehaus.groovy.jsr223.GroovyScriptEngineImpl.eval(GroovyScriptEngineImpl.java:155)
    qupath.lib.gui.scripting.DefaultScriptEditor.executeScript(DefaultScriptEditor.java:926)
    qupath.lib.gui.scripting.DefaultScriptEditor.executeScript(DefaultScriptEditor.java:859)
    qupath.lib.gui.scripting.DefaultScriptEditor.executeScript(DefaultScriptEditor.java:782)
    qupath.lib.gui.scripting.DefaultScriptEditor$2.run(DefaultScriptEditor.java:1271)
    java.base/java.util.concurrent.Executors$RunnableAdapter.call(Unknown Source)
    java.base/java.util.concurrent.FutureTask.run(Unknown Source)
    java.base/java.util.concurrent.ThreadPoolExecutor.runWorker(Unknown Source)
    java.base/java.util.concurrent.ThreadPoolExecutor$Worker.run(Unknown Source)
    java.base/java.lang.Thread.run(Unknown Source)

downsamples seems to be an aray of doubles, and the server at that point is an imageOps server, but maybe not an ImageServer? Any ideas off the top of your head? I’ll try and troubleshoot what is happening on my end more over the weekend if not :slight_smile:

This is 0.2.3.

SCRIPT
/**
 * Transform an ImageServer to 8-bit in QuPath, with scaling.
 *
 * Note that this *kind of* works (I think) in v0.2... but is quite hack-y,
 * and may not be entirely reliable (i.e. I've seen some errors).
 *
 * @author Pete Bankhead
 */

import qupath.lib.images.servers.PixelType
import qupath.opencv.ops.ImageOps
import qupath.lib.images.servers.ImageServers
import qupath.lib.images.writers.ome.OMEPyramidWriter

import static qupath.lib.gui.scripting.QPEx.*

def imageData = getCurrentImageData()

// Output downsample (may be 1 for full-resolution)
double outputDownsample = 1

// Output server path
def path = buildFilePath(PROJECT_BASE_DIR, '6454-02B-8bit.ome.tif')

// Create an op to handle all transforms
def op = ImageOps.buildImageDataOp()
    .appendOps(
            ImageOps.Core.multiply(100),
            ImageOps.Core.ensureType(PixelType.UINT8)
    )

// Create a scaling & bit-depth-clipping server
def server = ImageOps.buildServer(
        imageData,
        op,
        imageData.getServer().getPixelCalibration().createScaledInstance(outputDownsample, outputDownsample))

// Ensure the server is pyramidal
List<Double> downsamples = [outputDownsample]
double d = outputDownsample * 4
while (Math.min(server.getWidth()/d, server.getHeight()/d) >= 512) {
    downsamples << d
    d *= 4
}
print "server ---" +server
print "downsamples ---" +downsamples
server = ImageServers.pyramidalize(server, downsamples)

// Write the pyramid
new OMEPyramidWriter.Builder(server)
        .parallelize()
        .downsamples(downsamples as double[])
        .bigTiff()
        .channelsInterleaved()
        .build()
        .writePyramid(path)

Just recently saw the tip about the groovy addition to the scripting formatting. Nice!

@Mike_Nelson oops, mustn’t has posted the final version – you need to include as double[], which is the groovy trick to convert a list of Doubles (what you have) to an array of doubles (what you need).

Meanwhile, here’s a completely different approach that you might try:

/**
 * Transform an ImageServer to 8-bit in QuPath, with scaling.
 *
 * Note that this *kind of* works (I think) in v0.2... but is quite hack-y,
 * and may not be entirely reliable (i.e. I've seen some errors).
 *
 * @author Pete Bankhead
 */


import qupath.lib.color.ColorModelFactory
import qupath.lib.images.servers.ImageServer
import qupath.lib.images.servers.ImageServerBuilder
import qupath.lib.images.servers.ImageServerMetadata
import qupath.lib.images.servers.PixelType
import qupath.lib.images.servers.TransformingImageServer
import qupath.lib.images.writers.ome.OMEPyramidWriter
import qupath.lib.regions.RegionRequest

import java.awt.image.BufferedImage
import java.awt.image.DataBuffer
import java.awt.image.WritableRaster

import static qupath.lib.gui.scripting.QPEx.*

def imageData = getCurrentImageData()

// Output server path
def path = buildFilePath(PROJECT_BASE_DIR, 'converted-8bit-other.ome.tif')

// Create a scaling & bit-depth-clipping server
def server = new TypeConvertServer(imageData.getServer(), 100f, 0f)

// Write the pyramid
new OMEPyramidWriter.Builder(server)
        .parallelize()
        .downsamples(server.getPreferredDownsamples())
        .bigTiff()
        .channelsInterleaved()
        .build()
        .writePyramid(path)
print 'Done!'

class TypeConvertServer extends TransformingImageServer<BufferedImage> {

    private float scale = 1f
    private float offset = 0
    private ImageServerMetadata originalMetadata
    def cm = ColorModelFactory.getDummyColorModel(8)

    protected TypeConvertServer(ImageServer<BufferedImage> server, float scale, float offset) {
        super(server)
        this.scale = scale
        this.offset = offset
        this.originalMetadata = new ImageServerMetadata.Builder(server.getMetadata())
            .pixelType(PixelType.UINT8)
            .build()
    }

    public ImageServerMetadata getOriginalMetadata() {
        return originalMetadata
    }

    @Override
    protected ImageServerBuilder.ServerBuilder<BufferedImage> createServerBuilder() {
        throw new UnsupportedOperationException()
    }

    @Override
    protected String createID() {
        return TypeConvertServer.class.getName() + ": " + getWrappedServer().getPath() + " scale=" + scale + ", offset=" + offset
    }

    @Override
    String getServerType() {
        return "Type converting image server"
    }

    public BufferedImage readBufferedImage(RegionRequest request) throws IOException {
        def img = getWrappedServer().readBufferedImage(request);
        def raster = img.getRaster()
        int nBands = raster.getNumBands()
        int w = img.getWidth()
        int h = img.getHeight()
        def raster2 = WritableRaster.createInterleavedRaster(DataBuffer.TYPE_BYTE, w, h, nBands, null)
        float[] pixels = null
        for (int b = 0; b < nBands; b++) {
            pixels = raster.getSamples(0, 0, w, h, b, (float[])pixels)
            for (int i = 0; i < w*h; i++) {
                pixels[i] = (float)GeneralTools.clipValue(pixels[i] * scale + offset, 0, 255)
            }
            raster2.setSamples(0, 0, w, h, b, pixels)
        }
        return new BufferedImage(cm, raster2, false, null)
    }

}

Main difference: it retains the same resolution, doing it with the ‘tile at a time’ method. Avoids ops altogether.

Edit: Fixed scale and offset…

3 Likes

Haha, I noticed that as I was trying to understand the code while the first run was going. That also writes out incredibly fast, I think it about two minutes per image. That will do perfectly!

1 Like

@cgohlke Is it possible to write channel by channel or more generally, plane by plane with pyramidal ome.tiff in tifffile? I have large planes and for purposes of RAM, I’d like to write one in at a time as they are processed.

Here is an example of what happens using a modified example from the docs:

import numpy
import tifffile
import zarr

data = numpy.random.randint(0, 255, (30, 301, 219), 'uint8')

with tifffile.TiffWriter('frame_temp.tif') as tif:
     for channel in range(data.shape[0]):
         tif.write(data[channel,:,:], contiguous=True)

imdata = tifffile.TiffFile("frame_temp.tif").series[0]
# each channel is saved as a page
print(tifffile.TiffFile("frame_temp.tif").series[0].pages)
>>> 30

data = numpy.arange(3*1024*1024, dtype='uint8').reshape((3, 1024, 1024))
with tifffile.TiffWriter('pyr_temp.ome.tif') as tif:
    for channel in range(data.shape[0]):
        options = dict(tile=(256, 256))
        tif.write(data[channel,:,:], subifds=2, **options)
        # save pyramid levels to the two subifds
        # in production use resampling to generate sub-resolutions!
        tif.write(data[channel,::2,::2], subfiletype=1, **options)
        tif.write(data[channel,::4,::4], subfiletype=1, **options)

tfzarrstore = tifffile.imread("pyr_temp.ome.tif", aszarr=True, series=0)
zarr_pyr = zarr.open(tfzarrstore)

# shows channels are separated as series
[print(zarr_pyr[z].shape) for z in zarr_pyr.array_keys()]
>>> (1024, 1024)
>>> (512, 512)
>>> (256, 256)

print(len(tifffile.TiffFile("pyr_temp.ome.tif").series))
>>> 3

Another thing from the docs… is this intended to generate 900 pages in a single series?

data = numpy.random.randint(0, 255, (30, 301, 219), 'uint8')
with tifffile.TiffWriter('temp.tif') as tif:
    for frame in data:
        tif.write(data, contiguous=True)

print(len(tifffile.TiffFile("temp.tif").series[0].pages))
>>> 900

Yes. You probably meant to write frame by frame, not write 30 times the whole data array:

    for frame in data:
        tif.write(frame, contiguous=True)

I will reply later re. writing OME-TIFF channel by channel.

2 Likes

Thanks, just FYI that is what is in the docs here:

Thanks. I’ll correct that.

To understand why the code is creating three series: TiffWriter.write is meant and optimized for writing numpy nd-arrays as a series of TIFF pages (IFDs) along with metadata (OME, ImageJ, or JSON) in the first page to later retrieve the series as a numpy ndarray. An exception is contiguous=True, which is meant for streaming uncompressed images such that they can later be memory mapped directly as a numpy ndarray. This mode does not work with tiling, compression, subifds, or OME metadata.

Although not that convenient or efficient for large number of pages, TiffWriter.write can be used to create a TIFF file page by page. Disable the automatic processing of metadata (metadata=None) and provide your own metadata:

import numpy
import tifffile

data = numpy.arange(3 * 1024 * 1024, dtype='uint8').reshape((3, 1024, 1024))

# create OME-XML metadata for the image
omexml = tifffile.OmeXml()
omexml.addimage(
    dtype=data.dtype,
    shape=data.shape,
    # specify how the image is stored in the TIFF file
    storedshape=(data.shape[0], 1, 1, data.shape[1], data.shape[2], 1)
)
omexml = omexml.tostring()

# write the image to a pyramidal OME-TIFF page by page, subifd by subifd
with tifffile.TiffWriter('pyramidal.ome.tif') as tif:
    for i, channel in enumerate(range(data.shape[0])):
        options = dict(
            tile=(256, 256), photometric='miniswhite', metadata=None
        )
        # write OME-XML to the ImageDescription tag of the first page
        description = omexml if i == 0 else None
        # write channel data
        tif.write(data[channel], subifds=2, description=description, **options)
        # save pyramid levels to the two subifds
        # in production use resampling to generate sub-resolutions!
        tif.write(data[channel,::2,::2], subfiletype=1, **options)
        tif.write(data[channel,::4,::4], subfiletype=1, **options)

# access the levels in the pyramidal OME-TIFF file through zarr
with tifffile.imread('pyramidal.ome.tif', aszarr=True) as store:
    zarr_group = zarr.open(store, mode='r')
    print(zarr_group)
    for level in zarr_group.array_keys():
        zarr_array = zarr_group[level]
        print(zarr_array)

Out:

<zarr.hierarchy.Group '/' read-only>
<zarr.core.Array '/0' (3, 1024, 1024) uint8 read-only>
<zarr.core.Array '/1' (3, 512, 512) uint8 read-only>
<zarr.core.Array '/2' (3, 256, 256) uint8 read-only>
2 Likes