volshow

Load MRI data into the workspace and remove the singleton dimension.

load mri
V = squeeze(D);

Generate a colormap and transparency (alpha) map suitable for MRI images.

intensity = [0 20 40 120 220 1024];
alpha = [0 0 0.15 0.3 0.38 0.5];
color = [0 0 0; 43 0 0; 103 37 20; 199 155 97; 216 213 201; 255 255 255]/255;
queryPoints = linspace(min(intensity),max(intensity),256);
alphamap = interp1(intensity,alpha,queryPoints)';
colormap = interp1(intensity,color,queryPoints);

This MRI scan has a non-uniform, or anisotropic, voxel size of 1-by-1-by-2.5 mm. Specify the transformation matrix that scales the image to the correct voxel dimensions.

sx = 1;
sy= 1;
sz = 2.5;
A = [sx 0 0 0; 0 sy 0 0; 0 0 sz 0; 0 0 0 1];

Create an affinetform3d object that performs the scaling.

tform = affinetform3d(A);

View the volume with the custom colormap, transparency map, and transformation. Drag the mouse to rotate the volume. Use the scroll wheel to zoom in and out of the volume.

vol = volshow(V,Colormap=colormap,Alphamap=alphamap,Transformation=tform);

This example uses a subset of the Medical Segmentation Decathlon data set [1]. The subset of data includes two CT chest volumes and corresponding label images, stored in the NIfTI file format.

Run this code to download the MedicalVolumNIfTIData.zip file from the MathWorks® website, then unzip the file. The size of the data file is approximately 76 MB.

zipFile = matlab.internal.examples.downloadSupportFile("medical", ...
    "MedicalVolumeNIfTIData.zip");
filepath = fileparts(zipFile);
unzip(zipFile,filepath)

The folder dataFolder contains the downloaded and unzipped data.

dataFolder = fullfile(filepath,"MedicalVolumeNIfTIData");

Specify the filenames of the volume and label image used in this example.

dataFile = fullfile(dataFolder,"lung_043.nii.gz");
labelDataFile = fullfile(dataFolder,"LabelData","lung_043.nii.gz");

Read the image data and the metadata from the image file.

V = niftiread(dataFile);
info = niftiinfo(dataFile);

Define a transparency map and color map for this volume. The values used in this example were determined using manual trial and error.

alpha = [0 0 0.7 1.0];
color = [0 0 0; 200 140 75; 231 208 141; 255 255 255] ./ 255;
intensity = [-3024 -700 -400 3071];
queryPoints = linspace(min(intensity),max(intensity),256);
alphamap = interp1(intensity,alpha,queryPoints)';
colormap = interp1(intensity,color,queryPoints);

This MRI scan has a nonuniform, or anisotropic, voxel size. Extract the voxel spacing from the file metadata, and define the transformation to display the volume with correct dimensions.

voxelSize = info.PixelDimensions;
sx = voxelSize(2);
sy= voxelSize(1);
sz = voxelSize(3);
A = [sx 0 0 0; 0 sy 0 0; 0 0 sz 0; 0 0 0 1];

Create an affinetform3d object that performs the scaling.

tform = affinetform3d(A);

View the volume as a 3-D object. Specify the rendering style as "CinematicRendering". The cinematic rendering style displays the volume based on the specified color and transparency for each voxel,with iterative postprocessing that produces photorealistic shadows and lighting.

viewer = viewer3d;
vol = volshow(V,Parent=viewer, ...
    RenderingStyle="CinematicRendering", ...
    Colormap=colormap, ...
    Alphamap=alphamap, ...
    Transformation=tform);

Pause to apply all postprocessing iterations before updating the display in Live Editor.

pause(3)
drawnow

Load a grayscale volume into the workspace and display the volume using volshow.

load("spiralVol.mat")
h = volshow(spiralVol);
viewer = h.Parent;
hFig = viewer.Parent;
drawnow

Specify the name of the GIF file in which to save the animation.

filename = "animatedSpiral.gif";

Aim the camera at the center of the volume.

sz = size(spiralVol);
center = sz/2 + 0.5;
viewer.CameraTarget = center;

Specify the number of frames in the animation, then create an array of camera positions in a circle around the center of the volume.

numberOfFrames = 12;
vec = linspace(0,2*pi,numberOfFrames)';
dist = sqrt(sz(1)^2 + sz(2)^2 + sz(3)^2);
myPosition = center + ([cos(vec) sin(vec) ones(size(vec))]*dist);

At each camera position, update the display and write the frame to the GIF file. You can play the file in a video viewer.

for idx = 1:length(vec)
    % Update the current view
    viewer.CameraPosition = myPosition(idx,:);
    % Capture the image using the getframe function
    I = getframe(hFig);
    [indI,cm] = rgb2ind(I.cdata,256);
    % Write the frame to the GIF file
    if idx==1
        % Do nothing. The first frame displays only the viewer, not the
        % volume.
    elseif idx == 2
        imwrite(indI,cm,filename,"gif",Loopcount=inf,DelayTime=0)
    else
        imwrite(indI,cm,filename,"gif",WriteMode="append",DelayTime=0)
    end
end

This example creates a large 500-by-500-by-2500 image volume. If your machine does not have enough memory to create and store the 2.5 GB volume, decrease imSize before running this example.

imSize = [500,500,2500];

Create a simulated 3-D image of bubbles, V. This can take several minutes.

V = rand(imSize,"single");
BW = false(size(V));
BW(V < 0.000001) = true;
V = bwdist(BW);
V(V <= 20) = 1;
V(V > 20) = 0;

If you try to display V directly, volshow returns an error that the volume is too large. Instead, create a blockedImage object that points to V and has a block size of 500-by-500-by-500 voxels.

bim = blockedImage(V,BlockSize=[500,500,500]);

Display the blockedImage using volshow. The volshow function reads blocks into memory one at a time and stitches individual block renderings to produce the final volume.

bVol = volshow(bim);