NAMIC Wiki:DTI:Tensor format
THIS DOCUMENT IS FOR HISTORICAL PURPOSE ONLY
Initial suggestions by Martin Styner (please feel free to comment/add etc)
- Simple Tensor Datastructure and IO, each voxel simply contains a 3x3 matrix.
- Only support for 3D. Is 2D or 4D support necessary?
- Often useful to store a "mask" or "confidence" value alongside the matrix coefficients.
Datastructure for processing within VTK
- vtkTensor class: http://www.vtk.org/doc/nightly/html/classvtkTensor.html
Datastructure for processing within ITK - [Done - itk::SymmetricSecondRankTensor, DiffusionTensor3D]
- 'PixelType' datastructure itk::TensorPixel<T> would be needed (current special types within itk are RGBPixel, RGBAPixel, and BloxPixel)
- The TensorPixel datastructure contains a 3x3 matrix: itk::Matrix<T,3,3>
- Image datastructure itk::Image<itk::TensorPixel<double>,3> . Should the dimensionality of the image be fixed to 3 and a special type be declared? E.g. itkTensorImage<T> = itk::Image<itk::TensorPixel<T>,3>
There seems to be a need for converting a vkTensor into an itk based Tensor Image as part of the vtkItk package.
- Should we support multiple formats (one for vtk and for itk IO)? I think so.
- VTK format is simply ASCII readable format using the [vtkDataObjectWriter] / [vtkDataObjectReader]
NOTE: Those readers/writer might disappear in future version of VTK in favor of the VTK-XML format. So I'd suggest to directly use vtkImageData / vtkPointData / SetTensors and save the result using those classes instead: [vtkXMLImageDataWriter] / [vtkXMLImageDataReader]
- Image format readable with ITK: As far as I know, none of the current Imaging formats (DICOM, Analyze, GIPL, Meta etc) support Tensor data (correct?). Now the [MetaImageIO] class supports Tensordata.
- Within ITK there is support for NRRD format, that allows the I/O of tensor data (BWH group uses this format)
Using NRRD for DWI and DTI File I/O (from Gordon Kindlmann)
The NRRD format is one way of conveying all the important information about diffusion-weighted MRI and diffusion tensor datasets. An up-to-date specification of the NRRD format [is available here]. There is a reader in Insight/Utilities/NrrdIO, which is being improved to better map NRRD fields to the corresponding ITK Image fields. The NRRD format was chosen for the Harvard datasets in NAMIC, [as described here]. To summarize that information, the following information is handled through existing fields in the NRRD header:
- Image resolution and basic geometrical properties handled by [dimension], [sizes], [thicknesses], and [centerings] fields.
- Image orientation is handled by [space], [space origin], [space directions], and [space units] fields.
- Numerical data representation and file info handled by [type], [encoding], [endian], and [data file] fields.
- Basic meta-information about value attributes is handled by [content] and [kinds] fields.
There is a recent addition to the NRRD format, the [[http://teem.sourceforge.net/nrrd/format.html#measurementframe measurement frame]] field, which clarifies what has always been a confusing aspect of working with DWI data- the relationship between the coordinate frame in which the coefficients of the diffusion gradients (or the diffusion tensor) are measured, and the coordinate frame of the image itself. The "measurement frame" field identifies the mapping from coordinates associated with the diffusion gradients/tensors to the world space coordinates associated with the image orientation. Pending finalization of the key/value pair proposal below, all the Harvard NAMIC datasets will be updated to include the new "measurement frame" field.
The fact that NRRD ultimately has an "everything is a scalar" mentality means that a DWI scan is represented as a 4-dimensional array. All the per-axis fields in the NRRD header reflect this. In particular:
- "sizes: sizeX sizeY sizeZ numImg" : This gives the sizes ("dimensions") of the array, and indicates exactly how many MR images there are (numImg). Note: this number includes the one or more non-diffusion-wieghted images that are part of the acquisition. If there are two non-diffusion-weighted image, and 12 DWIs, then numImg is 14. Example: "sizes: 256 256 80 14"
- "kinds: space space space list" : The [kinds] field in the NRRD header describes what kind of information lies along each axis of the array, so there is one kind for each of the 4 axes. Like all per-axis fields, the ordering of per-axis information goes from fastest to slowest axis. This means that the first three axes ("space space space") correspond to the spatial domain, and the last/slowest axis ("list") gives all the different DWI values. The DWI axis could hypothetically be on any axis, but it seems simplest to standardize on either the fastest or the slowest. Using the slowest axis allows the raw datafile to play well with the Analyze and NIFTY-1 formats.
The following proposes some clarifications to how the NRRD header is used to represent the information specific to a DWI acquisition, through the use of certain key/value pairs:
- "modality:=DWMRI" : This key/value pair explicitly indicates that the image is a diffusion-weighted MRI scan, and it implies that all of the following key/value pairs are also set in the header:
- "DWMRI_b-value:=b" : This key/value pair gives the (scalar) diffusion-weighting value, in units of s/mm^2. Example: "DWMRI_b-value:=1000". As described below, the actual magnitude of diffusion-weighting for each DWI is determined with some simple calculations based on the individual per-DWI gradient directions or B-matrices.
- For every index position NNNN along the slowest (DWI) axis, "DWMRI_gradient_NNNN:=x y z" or "DWMRI_B-matrix_NNNN:=xx xy xz yy yz zz" should record either the diffusion-sensitizing gradient direction, or the full B-matrix, respectively, as follows:
- For simplicity, it is probably best to have all DWMRI_gradient or all DWMRI_B-matrix keys, rather than some intermixing of both
- NNNN should start at 0000, and go up to numImg-1, always using 4 digits to represent the index (as if by the "%04d" format string to sprintf()).
- The gradient "direction" corresponding to a truly non-diffusion-weighted reference image is "0 0 0", and the B-matrix is "0 0 0 0 0 0". Example: "DWMRI_gradient_000:=0 0 0". However, there is no reason not to record the actual diffusion weighting associated with the nominally non-diffusion-weighted image, if it is known.
- Giving the B-matrix allows recording diffusion weighting which results from both diffusion-sensitizing gradients as well as imaging gradients. Ignoring the effects of imaging gradients, the B-matrix associated with gradient (x,y,z) is (x*x, x*y, x*z, y*y, y*z, z*z).
- Note that the B-matrix off-diagonal entries are not pre-multiplied by 2, so the actual DWI value for B-matrix "Bxx Bxy Bxz Byy Byz Bzz" and diffusion tensor D is modeled by: I_0*exp(-b*(Bxx*Dxx + 2*Bxy*Dxy + 2*Bxz*Dxz + Byy*Dyy + 2*Byz*Dyz + Bzz*Dzz)).
Note that the NNNN values inside the keys correspond exactly to a zero-based numbering of the samples along the DWI axis of the 4-D volume, and that either gradient or B-matrix information has to somehow identified for each and every axis position.
There are two simple conventions that increase the expressivity of the scheme described above:
- To facilitate representing the situation where multiple images have been acquired with the same diffusion weighting, and these images are contiguous along the slowest axis, the "DWMRI_gradient_NNNN" or "DWMRI_B-matrix_NNNN" key can be omitted if the value is identical to the preceding one. The index value (NNNN) along the slowest axis still has to be correct for the next (different) diffusion weighting.
- To facilitate representing acquisitions where there are a range of different diffusion weightings, the magnitudes of all the gradient direction vectors (or B-matrices) are calculated, and implicitly scaled so that the largest magnitude is 1.0. The magnitude of a gradient direction vector is determined with the usual L^2 norm, and the magnitude of a B-matrix is via the [Frobenius Norm]. It is after this magnitude rescaling that the nominal b-value (given via "DWMRI_b-value:=b") applies. Thus, gradient "directions" of different and non-unit lengths may be given, to signify different amounts of diffusion weighting. Also, gradient directions that are traditionally given by non-unit vectors (such as (1,1,0)) can retain their expected appearance, without having to be normalized to unit length.
Here is a complete though very contrived example to demonstrate the scheme described above. The volume is 256x256x50, and for every voxel there are 38 values, 2 non-diffusion-weighted values, and then multiple diffusion-weighted values. There are 2 DW values for each of 6 gradient directions at b=1000 s/mm^2, and then 4 DW value for each of (the same) 6 gradient directions at b=2000 s/mm^2.
NRRD0005 content: SomeIDNumber42 type: short dimension: 4 space: right-anterior-superior sizes: 256 256 50 38 thicknesses: NaN NaN 4 NaN space directions: (-0.859375,0,0) (0,0,-0.859375) (0,5,0) none centerings: cell cell cell none kinds: space space space list endian: big encoding: raw space units: "mm" "mm" "mm" space origin: (117.5,-93,119) data file: rawdatafile.img measurement frame: (1,0,0) (0,0,1) (0,1,0) modality:=DWMRI DWMRI_b-value:=2000 DWMRI_gradient_0000:= 0 0 0 DWMRI_gradient_0002:= 1 0 1 DWMRI_gradient_0004:= 1 1 0 DWMRI_gradient_0006:= 0 1 1 DWMRI_gradient_0008:= 1 -1 0 DWMRI_gradient_0010:=-1 0 1 DWMRI_gradient_0012:= 0 1 -1 DWMRI_gradient_0014:= 2 0 2 DWMRI_gradient_0018:= 2 2 0 DWMRI_gradient_0022:= 0 2 2 DWMRI_gradient_0026:= 2 -2 0 DWMRI_gradient_0030:=-2 0 2 DWMRI_gradient_0034:= 0 2 -2
Note that the [[http://teem.sourceforge.net/nrrd/format.html#measurementframe measurement frame]] field identifies the coordinate space in which the coefficients of the diffusion gradients are measured, which will in turn determine the coordinate frame of the estimated diffusion tensors.
Addendum (From Gordon and Raul)
This addendum proposes a change regarding the way the number of excitations (NEX) is encoded. In the current proposal, you can explicitly list each and every repeated gradient, or, you can represent the NEX implicitly by the absence of the DWMRI_gradient_### keys for the repeated gradients. Thus, after parsing the key/value pairs, the number of repetition has to be figured out by some non-trivial logic on the presence and absence of keys. This goes counter to a convention in key/value formats, in which every relevant value should be encoded with a corresponding key. That is, NEX should be explicitly encoded.
This addendum proposes that the implicit representation of NEX in the above proposal suffers from "too clever is dumb", and so it should be disallowed. Instead, there are two options for identifying repetitions of gradient directions:
- Explicitly indicating NEX by repeating gradient directions. This is consistent with the description above; it just says you can't omit repeated gradients. Consequently, the example header shown above would be instead:
.... DWMRI_gradient_0000:= 0 0 0 DWMRI_gradient_0001:= 0 0 0 DWMRI_gradient_0002:= 1 0 1 DWMRI_gradient_0003:= 1 0 1 DWMRI_gradient_0004:= 1 1 0 DWMRI_gradient_0005:= 1 1 0 ....
- Indicating NEX by adding a new DWMRI_NEX_#### key with a an integral value giving the number of excitation for the given gradient. The header would read then:
.... DWMRI_gradient_0000:= 0 0 0 DWMRI_NEX_0000:= 2 DWMRI_gradient_0002:= 1 0 1 DWMRI_NEX_0002:= 2 DWMRI_gradient_0004:= 1 1 0 DWMRI_NEX_0004:= 2 ....
Not having a "NEX" key implies the default value of 1.
The upshot of all this is that there is more opportunities for error checking. The absence of a DWMRI_gradient_#### key no longer has semantic significance: it is either an error, or it was expected from the presence of an earlier DWMRI_nex_#### key.
What about DICOM?
- What is the potential for using DICOM for the DT format?
- What is the potential for using DICOM for processed DT format?
Discussion was moved to NAMIC_Wiki:DTI:Nrrd_format#DICOM:_Current_state