Difference between revisions of "3D Reslicing using COMKAT image tool (basic)"

From COMKAT wiki
Jump to navigation Jump to search
Line 2: Line 2:
  
 
===Overview===
 
===Overview===
Reslicing a 3D (or 3D vs time) image dataset can be accomplished using the COMKAT image tool and sliceVolume().  
+
Reslicing a 3D (or 3D vs time) image dataset can be accomplished using various components of COMKAT including the function sliceVolume().  
This example explains how to create image slices from a volume in at a position, plane orientation, and magnification.
+
This example explains how to create an image be slicing from a volume at a position, plane orientation, and magnification specified by the user.
The approach is to load the image volume dataset into an instance if an ImageVolumeData (abbreviated IVD) object and to use the sliceVolume() method.
+
The approach is to load the image volume dataset into an instance oif an ImageVolumeData (abbreviated IVD) object and to use the sliceVolume() method.
  
 
====Background====
 
====Background====
Line 11: Line 11:
 
----
 
----
 
====Approach I. Demonstrate the method for coordinate transformations====
 
====Approach I. Demonstrate the method for coordinate transformations====
Load a file using COMKAT image tool and use the functions to translate and rotate the image volume as you desire.  
+
Create an instance of an IVD loaded with an image volume.
  
  ivd = ImageVolumeData();
 
  
 +
  ivd = ImageVolumeData(); % create an instance of an ImageVolumeData object
  
Read data to ivd,  e.g. read_DICOM()
+
 
 +
Load an image volume data into ivd,  e.g. by reading DICOM offline storage files.
  
 
   ivd = read_DICOM(ivd, pathName, fileName); % load volume into an instance of IVD object;
 
   ivd = read_DICOM(ivd, pathName, fileName); % load volume into an instance of IVD object;
  
  
Create a list of indicies of all pixels in slice that we are creating
+
Create lists of indices of all pixels in the 2D slice (rectangular grid) that we are creating
  
   [i, j] = meshgrid(0 : Nc-1, 0 : Nr-1);
+
   [i, j] = meshgrid(0 : Nc-1, 0 : Nr-1); % i and j will be 2D arrays, meshgrid is a function built into MATLAB
 
     ij = [c(:)’ ; r(:)’];  % make matrix, each column corresponding to a single pixel in the slice we are creating
 
     ij = [c(:)’ ; r(:)’];  % make matrix, each column corresponding to a single pixel in the slice we are creating
  
  
Specify the transformation matrix based on desired pixel spacing (zoom), position (location), orientation. ref [http://medical.nema.org/dicom/2004/04_03PU.PDF] p. 275
+
Compute the physical x,y,z locations, in mm, from pixel indices according to the DICOM coordinate system ref [http://medical.nema.org/dicom/2004/04_03PU.PDF] p. 275.
 +
This is a two-step process.  The first step is to compute the coordinate transformation matrix.  Note that pixel spacing/zoom, orientation, and position for the slice are specified in the transformation matrix.
  
 
   M = ( Insert the method for generating the transformation matrix );
 
   M = ( Insert the method for generating the transformation matrix );
  
  
Calculate physical (mm) location of each pixel in the desired slice
+
The second step is to use the transformation matrix to calculate the physical (mm) location of each pixel in the desired slice
  
 
   xyz = M * ij;
 
   xyz = M * ij;
  
 +
 +
These xyz locations are the same as those in the image volume that is being sliced to make the 2D image.  From these xyz locations, we find the corresponding 3D indices, (u,v,w), into the volume.  This uses the transformation matrix for the volume, Mhat, that relates the indices to the xyz physical location.  This is analogous to M used for the desired slice but here the pixel spacing, orientation, and position indicate how the volume data are stored.
  
 
Specify the reverse mapping matrix ( xyz --> index space of the original image volume )
 
Specify the reverse mapping matrix ( xyz --> index space of the original image volume )
Line 44: Line 48:
 
Calculate voxel indcies into the volume corresponding to xyz physical location
 
Calculate voxel indcies into the volume corresponding to xyz physical location
  
   uvw =  Mhat * xyz;
+
   uvw =  Mhat * xyz;   NO, you have to invert Mhat
 +
 
 +
      comment on the size of uvw
 +
 
 +
Separate uvw into the components
 +
    u =
 +
    v =
 +
    w =
  
  
Line 50: Line 61:
  
 
   slice = sliceVolume(idv, v, u, w, backgroundPixelValue, ‘linear);
 
   slice = sliceVolume(idv, v, u, w, backgroundPixelValue, ‘linear);
 
+
EXPLAIN WHAT IS backgroundPixelValue
  
 
Display new slice
 
Display new slice

Revision as of 01:56, 8 August 2012

Reslicing 3D image volume using COMKAT image tool (basic)

Overview

Reslicing a 3D (or 3D vs time) image dataset can be accomplished using various components of COMKAT including the function sliceVolume(). This example explains how to create an image be slicing from a volume at a position, plane orientation, and magnification specified by the user. The approach is to load the image volume dataset into an instance oif an ImageVolumeData (abbreviated IVD) object and to use the sliceVolume() method.

Background

sliceVolume() is a mex-file written in c with an interface to MATLAB that makes the operation particularly efficient. COMKATImageTool uses sliceVolume() and you can use it too.


Approach I. Demonstrate the method for coordinate transformations

Create an instance of an IVD loaded with an image volume.


  ivd = ImageVolumeData(); % create an instance of an ImageVolumeData object


Load an image volume data into ivd, e.g. by reading DICOM offline storage files.

  ivd = read_DICOM(ivd, pathName, fileName); % load volume into an instance of IVD object;


Create lists of indices of all pixels in the 2D slice (rectangular grid) that we are creating

  [i, j] = meshgrid(0 : Nc-1, 0 : Nr-1); % i and j will be 2D arrays, meshgrid is a function built into MATLAB
   ij = [c(:)’ ; r(:)’];  % make matrix, each column corresponding to a single pixel in the slice we are creating


Compute the physical x,y,z locations, in mm, from pixel indices according to the DICOM coordinate system ref [1] p. 275. This is a two-step process. The first step is to compute the coordinate transformation matrix. Note that pixel spacing/zoom, orientation, and position for the slice are specified in the transformation matrix.

  M = ( Insert the method for generating the transformation matrix );


The second step is to use the transformation matrix to calculate the physical (mm) location of each pixel in the desired slice

  xyz = M * ij;


These xyz locations are the same as those in the image volume that is being sliced to make the 2D image. From these xyz locations, we find the corresponding 3D indices, (u,v,w), into the volume. This uses the transformation matrix for the volume, Mhat, that relates the indices to the xyz physical location. This is analogous to M used for the desired slice but here the pixel spacing, orientation, and position indicate how the volume data are stored.

Specify the reverse mapping matrix ( xyz --> index space of the original image volume )

  Mhat = ( Insert the method for generate the mapping);


Calculate voxel indcies into the volume corresponding to xyz physical location

  uvw =  Mhat * xyz;    NO, you have to invert Mhat
      comment on the size of uvw

Separate uvw into the components

    u = 
   v =
   w =


Use sliceVolume() to interpolate the slice

  slice = sliceVolume(idv, v, u, w, backgroundPixelValue, ‘linear);

EXPLAIN WHAT IS backgroundPixelValue

Display new slice

  figure, imagesc(slice); axis image % isotropic



Approach II. Use coordinateGen() to do the coordinate transformation

  • This should create same result as approach I but require fewer lines of coding since coordinateGen() does most things that are needed.


Read the image volume into an ImageVolumeData object

  ivd = ( Insert the method for reading data );


Use coordinateGen() to generate uvw

  [u, v, w] = coordinateGen(ivd, Nc, Nr, pixelSpacing, planePos, orientation); % Input the desired pixelSpacing, planePos and orientation matrices


Use sliceVolume() to interpolate

  slice = sliceVolume(idv, v, u, w, backgroundPixelValue, ‘linear);


Display slice

  figure, imagesc(slice); axis image % isotropic