## Using a given colour table with Image?

11

6

Context

I would like to represent large images with a given colour table.

Now, if I use Image

dat = RandomReal[{0, 1}, {1024, 1024}];
dat // Image; // Timing

(* ==> {0.000027, Null} *)


its fast, but in grayscales; on the other hand if I use, say MatrixPlot

dat // MatrixPlot[#, ColorFunction -> "Temperature"] &
dat // MatrixPlot[#, ColorFunction -> "Temperature"] &; // Timing

(* ==> {1.5748, Null} *)


its in colour, but its slow.

Question

Is there a method to get the best of both worlds? (i.e. Speed and chosen colour table).

2For completeness, I'll point out that the built-in function for applying a colour map to a grayscale image is Colorize. This is about twice as fast as ArrayPlot, though of course not nearly as fast as Mr.Wizard's renderImage. – None – 2015-01-16T17:05:17.580

I updated my answer. I wonder if there is a faster method available. I think perhaps a compiled function working on the image data directly would do it. – Mr.Wizard – 2013-03-16T15:38:17.223

I finally remembered why I thought Raster was faster: this comment by Vitaliy Kaurov.

– Mr.Wizard – 2013-03-16T16:04:03.283

14

I think I finally succeeded in creating something faster.

Edit: now ~40X faster than ArrayPlot.

renderImage[
array_?MatrixQ,
cf_,
q_Integer: 2048,
opts : OptionsPattern[Image]
] :=
Module[{tbl},
tbl = List @@@ Array[cf, q, {0, 1}] // N // DeveloperToPackedArray;
Image[tbl[[# + 1]] & /@ Round[(q - 1) array], opts]
]


A test of function:

dat = Map[Mean, ImageData[Import["ExampleData/lena.tif"]], {2}];

ArrayPlot[dat, ColorFunction -> "Rainbow"]

renderImage[dat, ColorData["Rainbow"], ImageSize -> 300]


A test of speed:

big = RandomReal[1, {1500, 1500}];

ArrayPlot[big, ColorFunction -> "Rainbow"] // Timing // First

renderImage[big, ColorData["Rainbow"], ImageSize -> 300] // Timing // First


2.325

0.0624

And this time that's correct timing data.

### Update

I have added a parameter q to control the number of quantization steps used. It arbitrarily defaults to 2048 which appears to be visually sufficient for most schemes and images. Examples of effect on quality and timing:

renderImage[dat, ColorData["Rainbow"], #, ImageSize -> 300] & /@ {7, 10000}


Needs["GeneralUtilities"]

BenchmarkPlot[
{renderImage[big, ColorData["Rainbow"], #] &},
Identity,
5^Range[9]
]


Is it fair to say it is still significantly slower than Image? – chris – 2013-03-16T16:49:37.940

@chris Certainly, but anything is likely to be as the color look-up has to to take some amount of time. Still I believe a compiled function would be faster but that's not my strength so I'll leave it for someone else. This at least lays the foundation. – Mr.Wizard – 2013-03-16T16:52:40.763

@chris I upgraded my function and it is now about 40 times faster than ArrayPlot. Please take a look. – Mr.Wizard – 2013-03-16T17:15:36.153

That's a significant improvement indeed. – chris – 2013-03-16T17:21:46.547

The thing is: by creating a lookup table for the colors you do a quantization which is not done by ArrayPlot. Although, the visual outcome may look the same, we compare two different things here. – halirutan – 2013-03-16T17:47:55.630

The performance drops significantly if a customized colour table is used cf this question

– chris – 2013-03-16T17:50:26.843

@Mr.Wizard The steps are clearly visible when somewhat customized data is used: renderImage[ConstantArray[Append[Table[x, {x, 0, 0.01, 0.0001}], 1], 100],ColorData["Rainbow"], ImageSize -> 300] compare this by using the same array with ArrayPlot. – halirutan – 2013-03-16T17:52:51.417

@chris I'll look into that. – Mr.Wizard – 2013-03-16T17:55:49.193

@halirutan you too; I increased the sampling rate 8X -- I hope that solves the problem. – Mr.Wizard – 2013-03-16T18:12:48.417

@Mr.Wizard perfect! – chris – 2013-03-16T18:14:13.483

@chris Thanks for the Accept, but this almost certainly deserves some additional tuning. That will have to wait for tomorrow. – Mr.Wizard – 2013-03-16T18:15:01.087

This is brilliant. I have been using Raster and then rasterizing to create images with a ColorFunction, this is so much faster. This is going straight into the toolbox (with the addition of a Clip to handle data outside the range 0-1). – Simon Woods – 2013-03-16T18:18:42.733

@Simon I'm sure you can refine this code; please let me know what you come up with. I intend to explore the sampling a bit more. Maybe use Rescale if Max[data] > 1? Also flattening the data and partitioning after should help with images with many rows and fewer columns. Like I said, tuning. – Mr.Wizard – 2013-03-16T18:30:15.503

@Mr.Wizard, I have nothing to add :-) It seems to me that the overhead of flattening and partitioning won't pay for itself unless the image is so tall and thin that you probably don't want to display it as an image anyway. Personally I prefer not to automatically Rescale, as I often want a fixed mapping of value to colour across multiple images. – Simon Woods – 2013-03-17T12:37:08.623

As for the sampling, I don't think there's a single "correct" answer unless you can make some assumptions about how smoothly the ColorFunction varies. For the nice smooth built-in gradients like "Rainbow" and "SolarColors" your 2048 samples gets you about 90% of the total number of unique colours that the function can create, which is probably good enough for most purposes. – Simon Woods – 2013-03-17T12:37:41.807

1

This answer was posted in error. Nevertheless I think the information below is helpful.

I believe the fastest general method is Raster, like this:

Graphics[Raster[RandomReal[1, {10, 20}], ColorFunction -> "Rainbow"]]


Actually, this isn't any faster than MatrixPlot, it's just different. With MatrixPlot the time is spent when the graphic is created, and with Raster it is spent when it is displayed:

Timing[g1 = MatrixPlot[dat, ColorFunction -> "Temperature"];]
Timing[g2 = Graphics[Raster[dat, ColorFunction -> "Temperature"]];]

{0.639, Null}

{0., Null}


To see the rendering time set:

SetOptions[\$FrontEndSession, EvaluationCompletionAction->"ShowTiming"]


Then:

g1

g2


and you will see that g1 displays immediately, whereas g2 takes about as long to render as it did to create g1`.