How to generate a color coded vector field

Hello everyone,

I am working on Diffusion Tensor Imaging and I am trying to generate a color coded FA map. I have generated the diffusion tensor matrix and extracted the eigenvalues and eigenvectors for each voxel of a 256 x 256 image. For each voxel i have an intensity (between 0 and 1) and a 3x1 array representing its direction in space.

My question is how can I color code each voxel, acording to the direction of the main eigenvector? For example, let's say that if the vector has a direction of x axis (1,0,0) the color of the voxel would be red; y axis (0,1,0) green; and z axis (0,0,1) blue. And if the vector has a direction in between two axis (example x an y) the color of the voxel would be a mix of red and green, weightned accordingly if it is closer to x or y axis. Example: http://imagebank.osa.org/getImage.xqy?img=LmxhcmdlLG9lLTIzLTEyLTE1MTM0LWcwMDg The color of each voxel would be a color on the surface of the sphere, according to the direction of the vector corresponding to that voxel (the origin of the vector is in the center).

I know there are other sofware (e.g. FSL) that compute color coded FA maps but I am designing a matlab code from scratch for a project and haven't figured out how to generate it.

Any help would be appreciated! Thanks

[pina ricciuti]

Hi Bruno, I have the same problem. How did you solve it in matlab?

Answer 1

[Star Strider]

This is one possibility:

data = rand(100,3);                                 % Original Data
datascl = bsxfun(@times, data, [100 50 1]);         % ‘Weigh’ Data
datadisp = reshape(sum(datascl,2), 10, 10);         % Sum & Reshape
colormap(jet);                                      % Set Appropriate Colormap

figure(1)
image(datadisp)                                     % Display Result

The ‘data’ matrix creates a three-element vector for each row on the interval (0,1). The bsxfun call creates an ‘artificial’ value matrix, weighing the x-direction by 100, the y-direction by 50, and keeps the z-direction unscaled. The idea is that in the jet colormap, the highest values are red, the centre values are green, and the lowest values are blue, so the vector created in the datadisp assignment is now weighted to colour the vectors in the x-direction red, y green, and z blue. It then sets the colormap to jet, since it most closely corresponds to the colours you want. The [100 50 1] approach may not be the best weighting scheme for your data, so you will likely have to experiment with it to get the result you want. This is simply an approach to setting the direction and colour vectors to the appropriate values. There may be other (and more efficient) approaches.

[bruno robalo]

Hello Star Trider,

Your answer is almost what I was looking for but I think it gave me an idea for the solution. If I understood well, in your example each component of the vector will give a color on the jet map, accordingly if it as the higher or lower value (0 to 1). Bu this means that if x component is higher, it will give red but in the next pixel, if x component is the middle value it is no longer red-ish. What I would like is that x component would be always red-ish (with intensity varying according to x value), y component always green-ish, and z component allways blue-ish. And in the end I would mix the three colors to form an RGB image...

[Star Strider]

The only way I can imagine to code the directions separately is to set the individual R, G, and B matrices, so that they will display each individually and all when combined. I don’t know how your original data are organised, so I’m totally guessing with this:

M = randi([0 1], 10, 10, 3);                            % Create Data (10x10 Binary ‘Image’ Of (X,Y,Z) Vectors, One In Each ‘Page’ Of The Image)

Rc = uint8(255*~M(:,:,1));                              % Set ‘Red’ If M(:,:,1) == 1
Gc = uint8(255*~M(:,:,2));                              % Set ‘Green’ If M(:,:,2) == 1
Bc = uint8(255*~M(:,:,3));                              % Set ‘Blue’ If M(:,:,3) == 1
Ac = cat(3, Rc, Gc, Bc);                                % Concatenate To Create Synthesis

figure(1)
subplot(4,1,1)
imshow(Rc)
title('X')
subplot(4,1,2)
imshow(Gc)
title('Y')
subplot(4,1,3)
imshow(Bc)
title('Z')
subplot(4,1,4)
imshow(Ac)
title('All')

The negation operator (~) is necessary because of the way imshow displays data, so that 0 is black, and 255 is white.