To imread Parula image in Matlab without losing resolution

Question

There is no bijection between RGB and Parula, discussed here. I am thinking how to do well the image processing of files in Parula. This challenge has been developed from this thread about removing black color from ECG images by extending the case to a generalized problem with Parula colors.

Data:

enter image description here

which is generated by

[X,Y,Z] = peaks(25);
imgParula = surf(X,Y,Z);
view(2);
axis off;

It is not the point of this thread to use this code in your solution to read the second image.

Code:

[imgParula, map, alpha] = imread('http://i.stack.imgur.com/tVMO2.png');

where map is [] and alpha is a completely white image. Doing imshow(imgParula) gives

enter image description here

where you see a lot of interference and lost of resolution because Matlab reads images as RGB, although the actual colormap is Parula. Resizing this picture does not improve resolution.

How can you read image into corresponding colormap in Matlab? I did not find any parameter to specify the colormap in reading.

OK. If you don't know what was the setting of caxis (i.e. the colormap's range) there is an indetermination. For instance, the blue could be any value depending on this setting ... A solution would give you an array ranged between 0 and 1 for instance, which may not be the initial data range. — Ratbert
An image file is not meant to represent arbitrary data, but graphics, which may or may not be representing data. If you want to store data generated in Matlab, use a mat file or maybe another binary or text-based file format. — A. Donda
I think I see how to do it and will try a solution in a few hours / tomorrow when I'm back home, if no solution has been posted in between. — Ratbert
I think it can be done in the following way: convert the black lines into white, apply regionprops on the non-white pixels, use the barycenters to map the regions on a regular grid and, finally, average the color on the region and map it to the parula colormap by looking at the closest color. Unfortunately, the first step depends too much on how the initial image has been generated (format, resolution), and the fact that the black grid is not completely black. — Ratbert
@Ratbert: see my attempt below. I'm not sure I covered every aspect of the problem – at least as you detailed – but I did implement a solution to find the "closest" color. — horchler

horchler horchler · Accepted Answer · 2015-05-05T19:52:37

The Problem

There is a one-to-one mapping from indexed colors in the parula colormap to RGB triplets. However, no such one-to-one mapping exists to reverse this process to convert a parula indexed color back to RGB (indeed there are an infinite number ways to do so). Thus, there is no one-to-one correspondence or bijection between the two spaces. The plot below, which shows the R, G, and B values for each parula index, makes this clearer.

Parula to RGB plot

This is the case for most indexed colors. Any solution to this problem will be non-unique.

A Built-in Solution

I after playing around with this a bit, I realized that there's already a built-in function that may be sufficient: rgb2ind, which converts RGB image data to indexed image data. This function uses dither (which in turn calls the mex function ditherc) to perform the inverse colormap transformation.

Here's a demonstration that uses JPEG compression to add noise and distort the colors in the original parula index data:

img0 = peaks(32);                     % Generate sample data
img0 = img0-min(img0(:));
img0 = floor(255*img0./max(img0(:))); % Convert to 0-255
fname = [tempname '.jpg'];            % Save file in temp directory
map = parula(256);                    % Parula colormap
imwrite(img0,map,fname,'Quality',50); % Write data to compressed JPEG
img1 = imread(fname);                 % Read RGB JPEG file data

img2 = rgb2ind(img1,map,'nodither');  % Convert RGB data to parula colormap

figure;
image(img0);                          % Original indexed data
colormap(map);
axis image;

figure;
image(img1);                          % RGB JPEG file data
axis image;

figure;
image(img2);                          % rgb2ind indexed image data
colormap(map);
axis image;

This should produce images similar to the first three below.

Example original data and converted images

Alternative Solution: Color Difference

Another way to accomplish this task is by comparing the difference between the colors in the RGB image with the RGB values that correspond to each colormap index. The standard way to do this is by calculating ΔE in the CIE L*a*b* color space. I've implemented a form of this in a general function called rgb2map that can be downloaded from my GitHub. This code relies on makecform and applycform in the Image Processing Toolbox to convert from RGB to the 1976 CIE L*a*b* color space.

The following code will produce an image like the one on the right above:

img3 = rgb2map(img1,map);

figure;
image(img3);                          % rgb2map indexed image data
colormap(map);
axis image;

For each RGB pixel in an input image, rgb2map calculates the color difference between it and every RGB triplet in the input colormap using the CIE 1976 standard. The min function is used to find the index of the minimum ΔE (if more than one minimum value exists, the index of the first is returned). More sophisticated means can be used to select the "best" color in the case of multiple ΔE minima, but they will be more costly.

Conclusions

As a final example, I used an image of the namesake Parula bird to compare the two methods in the figure below. The two results are quite different for this image. If you manually adjust rgb2map to use the more complex CIE 1994 color difference standard, you'll get yet another rendering. However, for images that more closely match the original parula colormap (as above) both should return more similar results. Importantly, rgb2ind benefits from calling mex functions and is almost 100 times faster than rgb2map despite several optimizations in my code (if the CIE 1994 standard is used, it's about 700 times faster).

RGB image of bird converted using two methods

Lastly, those who want to learn more about colormaps in Matlab, should read this four-part MathWorks blog post by Steve Eddins on the new parula colormap.

Update 6-20-2015: rgb2map code described above updated to use different color space transforms, which improves it's speed by almost a factor of two.

To imread Parula image in Matlab without losing resolution

1 Answers