Lending you a hand with image processing – basic techniques 1

Literally!

This is a PA ulnar deviation x-ray of my left wrist from last month, which gives a good view of the scaphoid bone from above.

The bone is chipped in the area pointed by the arrow, due to a fall that occurred 20 or so years ago. Somewhere in there, there’s also a tiny detached fragment of cartilage that calcified (as seen in a CT scan at the time). I was lucky, because typically the result of a fall with outstretched hand for people aged 17-40 is the scaphoid fracture, which are known to have unpredictable healing. Lately, however, due to a tendonitis, the fragment too is acting out. I’m left handed so this is causing some trouble, and that’s why the recent x-rays.

Family doctor, radiologist, and specialist all agreed there’s nothing I can do about the fragment, but taking care of the tendonitis with physiotherapy should improve the situation. The long-term outcome will depend on whether at the time of the fall the gliding surface of the cartilage was damaged in which case I’d get arthritis and there’s no way to know that unless I went in for an endoscopic procedure with a microcamera.

Nothing I can do. That’s when I decided I WAS going to do something at least with these x-ray. The idea stems from some good examples of image processing using a hand x-ray I saw some time ago on The Image Processing Cookbook by John C. Russ. I plan a whole (semi-regular) series taking examples from that book and from The Scientist and Engineer’s Guide to Digital Signal Processing. I will start with some basic image display, processing, and analysis techniques using my hand x-ray and other images, and then work my way up to a few more advanced techniques. I will post my Matlab code, which was written in Matlab version 2007a so you will need this or later version to run the code, but not the Image Processing Toolbox.

So yes, I am literally lending you a hand with image processing.

In this first post I will start with importing one of the x-ray images, then try a number of displays to see if they enhance details in the image. I will not be using the image above but rather one that shows more details of the hand and fingers. You can download it in here if you want to replicate exactly what I did. To import into Matlab just use the import command in the file menu.

Let’s start with displaying the original in grayscale. Here’s the code below. The first few lines will convert the data from RGB truecolor to intensity which are more compact (single matrix), easy to manipulate, and can still be displayed with different colormaps (for more on displaying data as images in Matlab read this excellent blog series).

%% convert inported jpeg x-ray to intensity matrix
hand=zeros(1044,797);
n=1:1:1044; k=1:1:797;
hand(n,k)=0.2989 * xray(n,k,1) + 0.5870* xray(n,k,2) ...
         + 0.1140 * xray(n,k,3); % creates intensity matrix
  % this is essentially what your printer does when you tell it to
  % print a color figure in black and white
%% display original for comparison
figure1=figure;
imagesc(hand);
colormap(bone);
colorbar;
set(gca,'YDir','reverse');
daspect([0.5 0.5 1]);
view(0,90);
lighting phong;
material shiny;
h=lightangle(45,60);

%% some figure and house cleaning
set(figure1,'Position',[640 300 950 708]);
set(figure1,'Color',[0.9 0.9 0.9]);
set(figure1,'OuterPosition',[636 396 958 790]);
title('Original X-ray - left wrist','FontSize',14);
axis off
axis tight

Here’s the resulting figure. This is how radiologists and doctors in emergency look at x-rays (for more on that read this post), the reason being that it requires virtually no processing, it is simple to view them on bright light panels, and fractures will stick out nicely.

But we are not interested in this image for diagnostic purposes only, so let’s see if we can enhance this display. Very often the first thing it is tried is to replace pseudo-colors for the intensity values. The rationale is that human eye in normal conditions can detect far more color hues (from many hundreds to up to millions) than shades of gray levels, and also gradual variations in color more easily. I wil ltry this next and use a rainbow colormap. To call the rainbow colormap in Matlab (called jet) and change the title we modify two lines in the code above to:

colormap(jet);
title('Original X-ray - left wrist - jet (rainbow) colormap','FontSize',14);

And here’s the new figure:

Looks a bit weird? Not surprisingly. It is unfortunate that many of the standard colormaps, like rainbow and spectrum do more harm than the good above mentioned: they often obscure some details and confuse the image by introducing artifacts. I will tackle this topic in depth in a series of future posts future post , and I will give some alternatives to the standard colormaps. For now, if you’re interested in the subject, please read this note by the IBM Research Centre, particularly the first two examples, which highlight how the rainbow colormap creates these perceptual artifacts. In addition, rainbow is very confusing for people with color deficiencies, which occur in only 0.4% of women, but up to an 8% in caucasian men. There’s an excellent simulation of these confusing effects in Fig.1 of this article of the American Geophysical Union’s EOS publication.

As suggested by Dr. Russ in The Image Processing Cookbook “there is no general rule for how people interpret different colors and therefore the successful use of pseudo-color displays is difficult”. A more successful approach is to create a pseudo-relief image by replacing elevation for intensity, and then combining it with the use of illumination, because surfaces are part of everyone’s everyday experience and we all intuitively recognize and interpret the effects of direct, reflected, and scattered light on surfaces and the effect of changes in its direction. It is the same ability that allows us to navigate around objects and tell the difference between a bump and a hole in the ground (we may decide to walk on the bump but we will avoid the hole). Steve Lynch provides in this paper a compelling empirical demonstration (through a survey) of why pseudo-surface is more successful than color, followed by sound scientific explanation based on the theory of human trivariate color vision.

Now let’s create this pseudo-elevation version of the image:

%% display result as 3D surface with phong lighting
figure2=figure;
surf(hand,'FaceColor','texturemap','EdgeColor','none');
% used hand for both the relief and the color intensity to be mapped
colormap(bone);
colorbar;
set(gca,'YDir','reverse');
daspect([0.5 0.5 15]);
view(0,90);
lighting phong;
material shiny;
h=lightangle(45,60);  % sets the position of the specified light
%% some figure and house cleaning
set(figure2,'Position',[640 300 950 708]);
set(figure2,'Color',[0.9 0.9 0.9]);
set(figure2,'OuterPosition',[636 396 958 790]);
title('Render of X-ray - left wrist','FontSize',14);
axis off
axis tight

Here is the result. This is to me a far superior image than both the color and the original. Notice how details reallu jum out of the page.

Can we do even better? I sometimes like to create a different matrix for the elevation to increase its dynamic range and enhance the relief effect. This is easier than just using the daspect command, and it gave good results with topographic or geophysical data. Let’s try it and see what happens with the hand:

%% creates elevation matrix
elev=log2(hand+1)*10; % creates arbitrary elevation matrix (I did this
   % through trial and error).
   % N.B. The function surf used below, would have
   % taken the variable hand as input for both relief
   % and color, but I like more control that what is
   % afforded with the command daspect alone
%% creates elevation matrix
elev=log2(hand+1)*10; % creates arbitrary elevation matrix (I did this
   % through trial and error).
   % N.B. The function surf used below, would have
   % taken the variable hand as input for both relief
   % and color, but I like more control that what is
   % afforded with the command daspect alone
%% display result as 3D surface with phong lighting
figure3=figure;
surf(elev,hand,'FaceColor','texturemap','EdgeColor','none');
  % used elev for relief and hand for color intensity to be mapped
  % to heated body colormap
colormap(bone);
colorbar;
set(gca,'YDir','reverse');
daspect([0.5 0.5 1]);
view(0,90);
lighting phong;
material shiny;
h=lightangle(45,60);  % sets the position of the specified light
%% some figure and house cleaning
set(figure3,'Position',[640 300 950 708]);
set(figure3,'Color',[0.9 0.9 0.9]);
set(figure3,'OuterPosition',[636 396 958 790]);
title('Render of X-ray - left wrist','FontSize',14);
axis off
axis tight
clear n k xray figure1 figure2 figure3;

And here’s the result.

To me this really stands our even more. But it also brought up the low intensity areas corresponding to the fleshy part of the hand, which now looks like the surface of the sea as seen from above when hit by sunlight. You know, all that glare? Is this perhaps distracting from the task of trying to interpret what you are looking at, i.e. the bones.

Let me know what is your opinion, which one do you like best?