Gigapixel Images in Mathematica
April 29, 2013 — Piotr Wendykier, Mathematica Algorithm R&D
Professional cameras offer a resolution of 50 megapixels and more. In addition, projects like GigaPan allow one to create gigapixel panoramas with billions of pixels. How can we process these images on a desktop computer with 8 GB of RAM?
One of Mathematica 9′s new and exciting features is out-of-core image processing. What does the out-of-core term really mean? It is a way to process very large images that are too big to fit into main memory. Let’s say we have a machine with 8 GB of RAM, and let’s assume that Mathematica can use up to 7.2 GB of that memory (the remaining 0.8 GB will be used by the operating system). Freshly started, Mathematica 9 on Windows 8 takes up about 200 MB of memory, so the kernel can use about 7 GB of RAM. What is the maximal size of the image that we can load into the kernel (we don’t want to visualize it at this point)? If we assume that the image is in the RGB color space and a single byte encoding, then the following formula gives a maximal width (and height) of an image that can be loaded at once into the memory:
50 k by 50 k is a very large image, but remember that we haven’t left any memory for processing and visualization. If we want to perform some operations on that image and visualize the results, then the realistic memory limit would be something around 0.5 GB. Why only 0.5 GB? We should leave at least 2 GB for the front end, and we have to keep in mind that any nontrivial image-processing algorithm requires additional storage for temporary objects. So, when we restrict a single image storage to 0.5 GB, then the maximal dimensions would be:
The above value is a good approximation provided that we only use image processing functions that preserve byte encoding and do not require internal conversion to Real32 or Real type. One example of such a function is ColorNegate. What if we want to call ImageConvolve? Then, not only the internal computations have to be performed using floating points, but also the resulting image is stored in the Real type (8 bytes per value). Our memory limit decreases to something like 0.1 GB, and the maximal width and resolution (in megapixels) become:
Thus, to do image processing on images larger than about 36 megapixels, we need out-of-core processing. Before we show how to do out-of-core processing on an image, let’s write some code that will allow us to visualize large images. In this post we will use an image from the NASA Visible Earth collection. It is a 222 megapixel JPEG file:
Let’s start with inspecting the data and visualizing a thumbnailed version of our large image. The function LargeImageInspect allows us to specify the range of rows and columns of the image that we want to see. Notice that we use the “FirstRow” and “LastRow” options of the Import function that were introduced in Mathematica 9:
Now, let’s create a function that shows a thumbnail of the large image. The function LargeImageThumbnail takes two arguments: a file path and optionally the width of a thumbnail image. It uses Import with the “FirstRow” and “LastRow” options to read a strip of a large image, re-sizes it according to the value of the second argument, and places the result in a list. At the end, ImageAssemble is called on that list to create a final thumbnail.
To explore large images, we need one more utility—an out-of-core image viewer that will allow us to visualize actual pixel values stored in a file. Here I implement the function LargeImageViewer that shows the thumbnail of a large image in the top pane, and a part of the actual image data in the bottom pane. To see the implementation and all the code, please download this post as a Computable Document Format (CDF) file.
Now that we know how to explore large images, let’s take a closer look at the out-of-core image processing functions. Mathematica 9 offers three functions (ImageFileApply, ImageFileFilter, and ImageFileScan) that allow processing of large images in chunks and provide file-to-file workflow. All three functions work with TIFF, JPEG, and PNG files.
ImageFileApply is an out-of-core equivalent of the ImageApply function introduced in Mathematica 7. It applies a function specified as a first argument to the list of channel values for each pixel of the image.
The above function call is equivalent to Export[outputFile, ColorNegate[Import[inputFile]]], however, ImageFileApply reads, processes, and writes successive blocks of image data so it requires significantly less memory than the in-core workflow. Here is the proof. The following code shows that ImageFileApply can be successfully executed when memory is constrained to 100 MB. Notice that in the call to ImageFileApply we used a non-default value of the “MaxBlockSize” option that defines the maximal size of image blocks in pixels. The default value of 8000000 is too large to run computations with the 100 MB memory constraint.
Trying to do the same thing with the in-core image processing functionality fails because it runs out of memory.
ImageFileApply‘s function parameter can be arbitrarily complex as long as it yields a number or a list of numbers. Here, we create a “blue snow” effect:
ImageFileApply is a pixel operation and cannot be used for neighborhood operations such as filters. For this purpose, Mathematica 9 has a function called ImageFileFilter. We will use that function on an image from the ExampleData paclet:
A simple example of ImageFileFilter is to replace every value by the minimum in its range-r neighborhood (MinFilter). In the example below we have not specified the output file. In that case, both ImageFileApply and ImageFileFilter create an output file in the user directory with a name that is formed from the input file name by appending a time stamp.
A popular Gaussian filtering is not much harder than the MinFilter. We create a kernel using the GaussianMatrix function, and then we sum up a product of that kernel with a neighborhood of each pixel, that is, Total[ker #,2]&:
The last function that I want to discuss here is ImageFileScan. It is an out-of-core equivalent of ImageScan (also new in Mathematica 9) that evaluates a given function applied to each pixel of an image. Since ImageFileScan discards the results and does not create any file, it is useful in carrying out operations that have a “side effect.” For example, we can use this function to obtain different image statistics such as minimum and maximum values or image histograms:
All three out-of-core functions have one more interesting option called “ImageList”. It allows us to specify which frames in a multiframe TIFF file should be processed. In this example we color negate only the second and third frame of a multiframe image:
The three functions that I have described here—ImageFileApply, ImageFileFilter, and ImageFileScan—are just the beginning of our work on out-of-core image processing. Ultimately, in future versions of Mathematica, the file-to-file workflow will be possible for all built-in image processing functions.