forked from JiadingGai/Impatient-MRI
-
Notifications
You must be signed in to change notification settings - Fork 0
SimpleTheoryOfTypes/Impatient-MRI
Folders and files
| Name | Name | Last commit message | Last commit date | |
|---|---|---|---|---|
Repository files navigation
(C) Copyright 2010 The Board of Trustees of the University of Illinois.
All rights reserved.
Developed by:
IMPACT & MRFIL Research Groups
University of Illinois, Urbana Champaign
Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to
deal with the Software without restriction, including without limitation the
rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
sell copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:
Redistributions of source code must retain the above copyright notice, this
list of conditions and the following disclaimers.
Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimers in the documentation
and/or other materials provided with the distribution.
Neither the names of the IMPACT Research Group, MRFIL Research Group, the
University of Illinois, nor the names of its contributors may be used to
endorse or promote products derived from this Software without specific
prior written permission.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS WITH
THE SOFTWARE.
------------------------------ Chapter separator ------------------------------
Prolog
------
This is the readme file of the MRI package. This package is contributed by the
IMPACT and MRFIL research groups in University of Illinois at Urbana-Champaign.
------------------------------ Chapter separator ------------------------------
Table of contents
-----------------
1. Introduction
2. Package revision history
3. Environment setup
4. How to compile and run
5. How to integrate it with your project
6. Todo list and future work
------------------------------ Chapter separator ------------------------------
1. Introduction
The MRI toolset is an implementation in CUDA for iterative MR image reconstruction
using Graphics Processing Units (GPUs). It is used in the research of medical ima-
ging, especially in the area of image reconstruction for magnetic resonance imaging
(MRI). In our experiments, this implementation significantly reduces the computation
times by two orders of magnitude (compared with the non-GPU implementation) while
accurately compensating for field inhomogeneity, parallel imaging(SENSE).
For more information about this work, please refer to the following publications.
* Jiading Gai, Joseph L. Holtrop, Xiao-Long Wu, Fan Lam, Maojing Fu, Justin P. Haldar,
Wen-mei W. Hwu, Zhi-Pei Liang, Bradley P. Sutton, More IMPATIENT: A Gridding-Accelerated
Toeplitz-Based Strategy for Non-Cartesian High-Resolution 3D MRI on GPU, Proceedings of
the International Society for Magnetic Resonance in Medicine (ISMRM), May 2012.
* Xiao-Long Wu, Jiading Gai, Fan Lam, Maojing Fu, Justin Haldar, Yue Zhuo, Zhi-Pei
Liang, Wen-Mei Hwu, Bradley Sutton, IMPATIENT MRI: Illinois Massively Parallel Accel-
eration Toolkit for Image reconstruction with ENhanced Throughput in MRI, Proceedings
of the International Society for Magnetic Resonance in Medicine (ISMRM), May 2011.
* Xiao-Long Wu, Jiading Gai, Fan Lam, Maojing Fu, Justin Haldar, Yue Zhuo, Zhi-Pei
Liang, Wen-Mei Hwu, Bradley Sutton, IMPATIENT MRI: Illinois Massively Parallel Accele-
ration Toolkit for Image reconstruction with ENhanced Throughput in MRI, Proceedings
of the IEEE International Symposium on Biomedical Imaging (ISBI), March 2011.
* Xiao-Long Wu, Yue Zhuo, Jiading Gai, Fan Lam, Maojing Fu, Justin P. Haldar, Wen-mei
Hwu, Zhi-Pei Liang, Bradley P. Sutton, "Advanced MRI Reconstruction Toolbox with Accel-
erating on GPU," SPIE 2011.
* Yue Zhuo, Xiao-Long Wu, Justin P. Haldar, Thibault Marin, Wen-mei Hwu, Zhi-Pei Liang,
Bradley P. Sutton, "Using GPUs to Accelerate Advanced MRI Reconstruction with Field
Inhomogeneity Compensation," GPU Computing Gems, Book Chapter 44, Elsevier Inc., 2011.
* Yue Zhuo, Xiao-Long Wu, Justin P. Haldar, Wen-mei Hwu, Zhi-Pei Liang, Bradley P. Sutton,
Multi-GPU Implementation for Iterative MR Image Reconstruction with Field Correction,
Proceeding of International Society for Magnetic Resonance in Medicine (ISMRM) 2010.
* Yue Zhuo, Xiao-Long Wu, Justin P. Haldar, Wen-mei Hwu, Zhi-Pei Liang, Bradley P. Sutton,
Accelerating Iterative Field-Compensated MR Image Reconstruction on GPUs, IEEE Internatio-
nal Symposium on Biomedical Imaging (ISBI) 2010.
------------------------------ Chapter separator ------------------------------
2. Package revision history
Package version : 3.1 alpha
Release date : 03/21/2012
Author list : Jiading Gai, Beckman Institute UIUC ([email protected])
Joseph L. Holtrop, Bioengineering UIUC ([email protected])
Xiao-Long Wu, ECE UIUC ([email protected])
Fan Lam, ECE UIUC ([email protected])
Maojing Fu, ECE UIUC ([email protected])
Synopsis : The third major release.
New Features :
- Two flags (-writeQ and -reuseQ <path>) have been added to facilitate the
customization of IMPATIENT to specific MRI application, such as fMRI. Recall
that Toeplitz method is perfect for applications such as fmri, where each
slice in every data point shares the same scan trajectory and image s-
ize. A significant portion of the computation in Toeplitz method is the
calculation of Q matrices, which facilitates multiplication operations
involving F^HF and depends only on the scan trajectory (not the scan da-
ta) and the size of the image. Hence, in a task like fmri, Q can be com-
puted once prior to acquiring an image's scan data and re-used on each
slice. In this way, the critical path for a given reconstruction consists
only of compputing F^H(d) and executing the CG solver. Please refer to section 4.3
for more details on the usage of -writeQ and -reuseQ.
Package version : 3.0 alpha
Release date : 03/16/2012
Author list : Jiading Gai, Beckman Institute UIUC ([email protected])
Joseph L. Holtrop, Bioengineering UIUC ([email protected])
Xiao-Long Wu, ECE UIUC ([email protected])
Fan Lam, ECE UIUC ([email protected])
Maojing Fu, ECE UIUC ([email protected])
Synopsis : The third major release.
New Features :
- In the previous release (2.0 alpha), we implemented a fast image reconstruction routine
for reconstruction of arbitrary 3D imaging trajectories using the Toeplitz computation
structure. In this release (3.0 alpha), we introduce a new IMPATIENT which removes bottlenecks
to the computation by implementing the formation of the Toeplitz operator and formation
of the data input structure in a fast, gridding implementation. Further, we adapt
the algorithm to perform field correction and parallel imaging for large, arbitrary
3D sampling trajectories. Speed-ups of up to a factor of 1000 are provided in the
improved GPU implementation compared to the previous accelerated GPU code.
- 3.0 alpha provides users with three different types of MRI reconstruction strategies.
In reverse chronological order, they are:
2) Toeplitz with gridding:
. Support both 2D and 3D, arbitrary trajectory.
. Fast evaluations of Q (step 1) and FHD (step 2) using gridding algorithm.
. Quickest among the three strategies, but approximate.
-> With a decent gridding oversampling factor, the average
percentage error introduced by gridding is around 1%
. Arbitrary single precision gridding oversampling ratio between [1.0,2.0].
-> It is found that significant reductions in computation memory and
time can be obtained while maintaining high accuracy by using a
minimal oversampling ratio, from 1.125 to 1.375, instead of the
typically employed grid oversampling ratio of two.
. Equipped with various regularization forms:
-> spatial roughness penalty,
-> quadratic regularization
-> total variation.
IMPATIENT 3.0 alpha provides two equivalent ways to implement the regularization
function: 1) sparse matrix ; 2) explicit finite difference calculations.
1) Toeplitz with direct evaluations of $Q$ and $F^HD$:
. Support both 2D and 3D, arbitrary trajectory.
. Direct evaluations of Q (step 1) and FHD (step 2) using their equations.
. slower than Toeplitz with gridding, but as accurate as Brute Force.
-> If Q can be re-used in an MRI experiment on multiple slice, it
can also be much faster than Brute Force
. Equipped with various regularization forms:
-> spatial roughness penalty,
-> quadratic regularization
-> total variation.
IMPATIENT 3.0 alpha provides two equivalent ways to implement the regularization
function: 1) sparse matrix ; 2) explicit finite difference calculations.
0) Brute Force:
. Support only 2D data, arbitrary trajectory.
. Exhaustive evaluation of the entire system equation (both forward and adjoint operators).
. Equipped with various regularization forms:
-> spatial roughness penalty,
-> quadratic regularization
-> total variation.
IMPATIENT 3.0 alpha provides two equivalent ways to implement the regularization
function: 1) sparse matrix ; 2) explicit finite difference calculations.
Package version : 2.0 alpha
Release date : 05/06/2011
Author list : Xiao-Long Wu, ECE UIUC ([email protected])
Jiading Gai, Beckman Institute UIUC ([email protected])
Fan Lam, ECE UIUC ([email protected])
Maojing Fu, ECE UIUC ([email protected])
Yue Zhuo, BIOE UIUC ([email protected])
Synopsis : The second major release.
New Features :
- Toeplitz reconstruction strategy is implemented to provide users with
a choice of the tradeoff between accuracy (Brute Force) and speed (To-
eplitz).
. Takes advantages of Toeplitz formulation for fast computation of sys-
tem matrix computations.
. Field inhomogeneity performed via time segmentation.
. Toeplitz involves approximation, but is also faster. Tradeoff between
number of time segments for accuracy and computation time.
. Toeplitz method is perfect for applications such as fmri, where each
slice in every data point shares the same scan trajectory and image s-
ize. A significant portion of the computation in Toeplitz method is the
calculation of Q matrices, which facilitates multiplication operations
involving F^HF and depends only on the scan trajectory (not the scan da-
ta) and the size of the image. Hence, in a task like fmri, Q can be com-
puted once prior to acquiring an image's scan data and re-used on each
slice. In this way, the critical path for a given reconstruction consists
only of compputing F^H(d) and executing the CG solver.
. More details on Toeplitz: JA Fessler, S Lee, VT Olafsson, HR Shi, DC Noll.
IEEE Trans Sign Proc, 53(9): 3393-3402.
- Toeplitz strategy supports sensitivity Encoding for Fast MRI recon-
struction (SENSE) and total variation (TV) based regularization in our SENSE
framework (solved using an alternating minimization algorithm called half-
quadratic minimization). We provide two ways to encode the regularization matrix
. finite difference (directly calculated on-the-fly).
. compressed row sparse matrix (pre-calculated prior to the recon task).
- Toeplitz method is designed to handle larger image size from 256x256 to
512x512.
- No CPU code is available.
Package version : 1.0 alpha
Release date : 11/19/2010
Author list : Xiao-Long Wu, ECE UIUC ([email protected])
Jiading Gai, Beckman Institute UIUC ([email protected])
Fan Lam, ECE UIUC ([email protected])
Maojing Fu, ECE UIUC ([email protected])
Yue Zhuo, BIOE UIUC ([email protected])
Synopsis : The first major release.
New Features :
- Supports Sensitivity Encoding for Fast MRI reconstruction (SENSE).
- Incorporates a total variation based regularization in our SENSE
framework. (solved using an alternating minimization algorithm called
half-quadratic minimization). We provide two ways to encode the
regularization matrix
. sparse matrix (pre-calculated prior to the recon task).
. finite difference (directly calculated on-the-fly).
- Enable SMVM to increase image quality but some bugs exist.
- FT and IFT are extended to handle larger image size from 256x256 to
512x512. (We have some bugs on 1024x1024 images.)
Package version : 0.1
Release date : 04/16/2010
Author list : Xiao-Long Wu, ECE UIUC ([email protected])
Yue Zhuo, BIOE UIUC ([email protected])
Synopsis : The first release.
Description : Here lists the key functions of this package.
- Use the algorithms mainly referred in the papers, ISMRM 2010 and ISBI
2010.
- Multi-GPU platforms are supported.
- GPU code is modified to accommodate GPU platforms with/out double-
precision hardware.
- CPU code is optimized with OpenMP.
- C++ language is used to facilitate the cooperation with other
projects.
------------------------------ Chapter separator ------------------------------
3. Environment setup
To compile this package you must have the following package installed.
3.1 A Linux-Like Environment
We have successfully tested the package on the following environment.
However, your system shall work if that's not listed below because
we didn't use platform-specific libraries.
- Ubuntu: 9.10, 10.04
- Debian: 4.3.5
- Fedora: 12 13 14
- Mac OS X: 10.06
We haven't tested the package under Windows platforms. But since we
don't use specific libraries for Linux-like platforms, you should be
able to use Windows-based compilers to build the package with/out
minor modifications.
3.2 NVIDIA Graphics Card
We have successfully tested the package with the following cards. Yet,
your graphics card shall work as long as it's supported by the NVIDIA
CUDA.
- Tesla C2050
- GeForce GTX 480
- GeForce GTX 280
- GeForce 8800 GTX
- GeForce 8600M GT
- Quadro FX 5600
3.3 NVIDIA CUDA Developer Driver
After you installed the operating system, the NVIDIA graphics card
device driver should have been installed. If not, you can download the
driver from one of the following links.
Driver version: 260.19.21
http://developer.nvidia.com/object/cuda_3_2_downloads.html
Driver version: 195.36.15
http://developer.nvidia.com/object/cuda_3_0_downloads.html
Driver version: 190.53
http://developer.nvidia.com/object/cuda_2_3_downloads.html
3.4 NVIDIA CUDA Toolkit
CUDA toolkit contains NVCC compiler and the relevant libraries.
You can download the toolkit from one of the following links.
If no special issues, we recommend you to use the latest one.
Toolkit version: 3.2
http://developer.nvidia.com/object/cuda_3_2_downloads.html
Toolkit version: 3.0
http://developer.nvidia.com/object/cuda_3_0_downloads.html
Toolkit version: 2.3
http://developer.nvidia.com/object/cuda_2_3_downloads.html
3.5 Compilers
We have successfully tested the package with the following compiler
versions. It is noted that GCC/++ 4.4 are not supported by NVCC
compilers yet. Besides, if you want to apply OpenMP to parallelize the
CPU code, you must have GCC/++ with version 4.2 or up.
- G++: 3.4, 4.1, 4.3, 4.4
- GCC: 3.4, 4.1, 4.3, 4.4
- CUDA Toolkit: 2.3, 3.0, 3.2, 4.0 RC2
------------------------------ Chapter separator ------------------------------
4. How to compile and run
4.1 Environment setting of the Makefile in the package
After you unzipped the package, please modify the Makefile file to fit
your environment setting and specific requirements for your need.
If you want to have a faster CPU implementation, you can enable OpenMP
to take advantage of the underlying multi-core platform. To enable
OpenMP, you must have GCC/++-4.2 or up and set OPENMP to be "true". The
default is false.
export OPENMP := "false"
To use the double precision computation hardware of your graphics card,
you can set DOUBLE_PRECISION to be "true". GT200 based cards support
double-precision, such as GTX260/275/280/285/295, Telsa C1060, Telsa
S1070, Quadro FX5800 and the new Fermi based cards. The default value
is "false". Note, enabling this variable on the platform not supporting
double precision hardware may cause execution failure.
export DOUBLE_PRECISION := "false"
Besides, please make sure your NVCC compiler is installed under
"/usr/local/cuda". If not, please modify the following variable.
export CUDA_INSTALL_PATH ?= /usr/local/cuda
Finally, you must make sure you provide the right compiler version for
your environment. Note, if your system has two compilers, for example,
GCC/++-4.3 and GCC/++-4.4, you must make sure /usr/bin/gcc and g++
point to the right compiler, GCC/++ version 4.3 in this example. NVCC
seems to use the library of /usr/bin/gcc as the default one during the
compilation without honoring the following variable definitions.
export CXX ?= g++
export CC ?= gcc
export LINKER ?= g++ -fPIC
export AR ?= ar
export RANLIB ?= ranlib
4.2 Compile and make
At the first time, to build the package, please execute the following
command. This shall build all files, including necessary libraries.
$> make clean_all; make CUDA_ARCH=-arch\ sm_xx
Later, if you want to modify some of the source files without
rebuilding all necessary libraries, please execute the following
command.
$> make clean; make CUDA_ARCH=-arch\ sm_xx
To build with more information, you can add "verbose=1" in the command
line, as listed below.
$> make clean verbose=1; make verbose=1 CUDA_ARCH=-arch\ sm_xx
If you modified some macro definitions in any header files, .cuh or .h,
you must "make clean; make" to make sure the changes are propagated to
all source files.
For advanced users who would like to work under debug mode, please use
the following make command.
$> make clean_all debug=1 debug_xcpplib=1;make debug=1 debug_xcpplib=1 CUDA_ARCH=-arch\ sm_xx
Note that, to compile correctly, it is necessary to specify the compute capability
of your GPU card in the Makefile command line (e.g., for a GTX 580, it should be
CUDA_ARCH=-arch\ sm_20; wile for a GTX8800, it should be CUDA_ARCH=-arch\ sm_10).
This is because gridding implementation prefers atomic operation on GPU if available.
4.3 Run the program
This chapter describes the basic usage of the IMPATIENT toolkit. For more
information, please see http://impact.crhc.illinois.edu/mri.php. Before
running the program, you must have data prepared. We provided a handful of
datasets for your reference under the "data" directory.The basic usage of
the IMPATIENT package to reconstruct an image is simple, and it typically
looks something like the following commands:
Example 1: Use the brute force strategy to reconstruct a dataset called
64x64x1 located under data directory.
(a) The number of CG iterations is 8 (by default). Use quadratic
regularization with the weighting parameter lambda set to be
10^-2.
$> ./mriSolver -idir mriData/64x64x1 -fdp 1e-2
(b) The number of CG iterations is 21. Use total variation regularization
with 1) the weighting parameter lambda set to be 10.0 and 2) the num-
ber of TV iteration set to be 16.
$> ./mriSolver -idir mriData/64x64x1 -cg_num 21 -tv -tv_num 16 -fdp 10.0
Example 2: Use the toeplitz strategies to reconstruct a dataset called 32x32x4:
(12 time segments are used to approximate field inhomogeneity)
(a) Toeplitz with gridding is used. The number of CG iterations is 10 (-cg_num). Use total
variation regularization (-tv) with the weighting parameter lambda set to be 10^4 (-tv_num).
The gridding oversampling factor for computing Q is set to 1.125 and the gridding
oversampling factor for computing FHD is set to 1.375 (-gridOS_Q and -gridOS_FH).
Sparse matrix encoded with the finite difference regularization matrix is used (-reg).
to perform the regularization.
$> ./mriSolver -idir mriData/32x32x4 -reg -cg_num 10 -toeplitzGridding \
-tv -tv_num 10 -gpu_id 2 -ntime_segs 8 -fdp 1e4 \
-gridOS_Q 1.125 -gridOS_FH 1.375
(b) To replicate all the parameters and experimental settings in (a).
Only this time, Toeplitz with direct evaluation is used and regularization
is done by explicit finite difference calculations (i.e., without -reg flag).
$> ./mriSolver -idir mriData/32x32x4 -cg_num 10 -toeplitzDirect \
-tv -tv_num 10 -gpu_id 2 -ntime_segs 8 -fdp 1e4
Example 3: To run Toeplitz recon on all the sample data sets and then display
all the sample recons via Matlab. Use the following Bash script:
#!/bin/bash
for file in `ls mriData`; do
export dataDir=mriData/$file;
./mriSolver -idir $dataDir -cg_num 10 -toeplitzGridding -tv -tv_num 10 \
-fdp 1e-6 -ntime_segs 8 -gridOS_Q 1.125 -gridOS_FH 1.375;
done
Also, there are more command options as listed below :
$> ./mriSolver <Options>
Mandatory options:
-idir <dir> MRI input data directory.
Image enhancement options:
-toeplitzDirect Use the toeplitz reconstruction strategy with direct evaluation.
The default reconstruction strategy is brute force.
-toeplitzGridding Use the toeplitz reconstruction strategy with gridding.
The default gridding oversampling factor for Q is 1.125.
The default gridding oversampling factor for FHD is 1.5.
The default reconstruction strategy is brute force.
-gridOS_Q <num> Specify the gridding oversampling factor for Q. <num> has
to be in the range [1.0,2.0]
-gridOS_FH <num> Specify the gridding oversampling factor for FHD. <num> has
to be in the range [1.0,2.0]
-ntime_segs <num> Specify how many time segments to be used by the Toeplitz
reconstruction strategies.
-writeQ This flag will cause the computed Q data structure to be written
to hard disk. This can be useful when a user wants to re-use the saved
Q for reconstruction of later slices, so as to save the time spent
to re-compute the same Q again.
-reuseQ <dir-to-Q> This flag tells IMPATIENT to skip the computation of Q (i.e., step 1).
Instead, a previously saved Q from hard disk will be re-used for
reconstruction of the current slice in the CG step (i.e., step 3).
Note that <dir-to-Q> has to be the path name to the directory that contains
the to-be-reused Q file. The to-be-reused Q file has to bear the filename
'Q_stone.file'
-cg_num <num> Number of iterations for the Conjugate Gradient method.
Default number of iteration is 8.
-reg Enable roughness regularization that incorporates a priori
information to stablize reconstructions when data-consistency
constraints alone do not result in a sufficiently well-posed
inverse problem.
-fd Enable finite difference penalizer. This penalizer denotes the
finite difference of every pixel pair along the horizontal and
vertical directions of the image respectively.
-fdp <num> Specify finite difference penalizer value. WARNING to the users
who would like to use regularization: this weight parameter lambda
is expected to significantly impact your reconstruction quality.
You are forced to choose the value of lambda properly based on your
data and application. IMPATIENT HAS NO DEFAULT VALUE FOR lambda and
NO AUTOMATIC WAY TO SELECT AN OPTIMAL lambda FOR YOU.
-tv Enable total variation regularization that penalizes the total amount
of gradient while preserving the edge information in the image.
-tv_num <num> Number of iterations for the total variation regularization. Default
number of iteration is 10.
-tv_update Enable updates of the outputs for each TV iteration.
Miscellaneous options:
-odir <dir> MRI output data directory. If the output directory is not provided,
"<idir>/output" is used.
-cpu Enable image reconstruction on CPUs as a comparison to the GPU counter-
parts. Default is not enabled. This flag is not applicable to toeplitz
strategy as there is no corresponding cpu code availabe for toeplitz.
-mgpu Enable multiple GPU computatoin when possible. When this option is not
specified, only one GPU is used.
-nogpu Disable the GPU part computation. This enables the CPU computation automatically.
-version Version of this release and the license.
-help View this message.
4.4 Input/Output data and the corresponding formats
Here we list the input/output files and their corresponding information
so that you can feed our program with your data. If you'd like to
modify the program to fit your needs, you can start from the function,
loadInputData(), in "structures.cpp".
MRI conventions in IMPATIENT:
1) k-space is unitless, and should be contained within the window -N/2
to N/2. Take actual k-space trajectory (with units of 1/m) and multiply
by the field of view (in m) to get unitless k-space. Gridding specially
requires this and is part of the preferred Toeplitz pathway. However,
if only brute force or Toeplitz with direct evaluation is used, then the
limit will not be strictly enforced.
2) Field map has a minus sign convention, field map phase term in the forward
signal model is exp(-i*\omega*t).
3) Field map is in units of radians/second. Time is in units of seconds
4.4.1 Inputs data
Each of our input MRI file consists of two parts: a header (in text format) and the
data (in binary format). The header consists of lines in text format of the following
form:
"<keywords> = <value>\n"
A complete list of header keywords, with brief explanations:
- "version":
A version identifier. The newest version number is equal to 0.20000.
- "xDimension"
An integer number denoting the number of columns in the output image.
- "yDimension"
An integer number denoting the number of rows in the output image.
- "zDimension"
An integer number denoting the number of depths in the output image.
- "coil_number"
The number of coils.
- "slice_number"
The number of slices.
- "file_size"
The total number of single precision floating numbers in the data section.
- "Binary_Size" and "Binary:"
The data in the data section may be divided into meaningful pieces, each
of which is labelled by "<Binary_size>=<value>\n" to denote the size (in
floats), then immediately followed by "Binary:" and the actual binary data
for the piece.
The data part of the input file are also written in binary format using the following mode:
- single precision floats with the vector contents.
The size are N pixels in the image space, M k-space samples per coil and P coils.
fm.dat Field map in units of Hertz. Positive sign
convention - (N-by-1)
idata_r.dat Real part of the initial image estiamte (typically
zeros) - (N-by-1)
idata_i.dat Imaginary part of the initial image estiamte
(typically zeros) - (N-by-1)
ix.dat X-coordinates in image space (pixel grid) - (N-by-1)
X-coordinates in terms of field of view, going from
-0.5 to 0.5-delta_x, indicating the proportion of the
field of view for the current spatial position.
iy.dat Y-coordinates in image space (pixel grid) - (N-by-1)
iz.dat Z-coordinates in image space (pixel grid) - (N-by-1)
kdata_i.dat Real part of the k-space data (M*P-by-1)
kdata_r.dat Imaginary part of the k-space data (M*P-by-1)
kx.dat X-coordinate of the k-space trajectory - (M-by-1)
Unitless k-space data format, where k-space typically
runs from -matrix_size/2 to matrix_size/2-1. This is
the k-space locations times the field of view.
ky.dat Y-coordinate of the k-space trajectory - (M-by-1)
kz.dat Z-coordinate of the k-space trajectory - (M-by-1)
mask.dat Image space computation mask - (N-by-1)
t.dat Acquisition time vector in seconds - (M-by-1)
sensi_r.dat Real part of the Coil sensitivity - (N*P-by-1)
sensi_i.dat Imaginary part of the Coil sensitivity - (N*P-by-1)
Note that in the 1.0 alpha release the iz and kz vectors are not used
(only the first element of the vector is used (should be zero). They
will be integrated in a further release.
Note that toeplitz part of the code also generates five intermediate
data files to facilites the computation of Q matrices etc. The five
intermediate files are ixQ.dat, iyQ.dat, izQ.dat and phiR.dat, phiI.dat.
These are computed based on the input files given above and are saved
in the same directory as the input files given above. The five inter-
mediate files are managed by the internal implementation and should be
irrelevant to end users.
4.4.2 Output data
(a) Brute Force strategy:
Files are written in binary format using the following mode:
- single precision floats with the vector contents.
output_cpu_i.dat Real part of the reconstructed image using CPU -
(N-by-1)
output_cpu_r.dat Imaginary part of the reconstructed image using
CPU - (N-by-1)
output_gpu_i.dat Real part of the reconstructed image using GPU -
(N-by-1)
output_gpu_r.dat Imaginary part of the reconstructed image using
GPU - (N-by-1)
(b) Both Toeplitz strategies:
The single output file named out.file is written in binary format
in the input data directory using the following mode:
- single precision floats
- out.file starts with real part of the reconstructed image and
concludes with imaginary part the reconstructed image.
- one possible Matlab script to display the reconstructed image:
fp=fopen('out.file','r');
recon = fread(fp,'single');
len = length(recon);
recon_r = recon(1:len/2);
recon_i = recon(len/2+1:end);
N = (len/2) ^ 0.5;
gpu_output = reshape(recon_r+1i*recon_i,[N,N]);
figure;
imagesc(abs(gpu_output));
colormap(gray);
Note that two sample Matlab scripts have been provided to display the
two output formats. They reside in the mriSolver directory:
- display_recon.m (Toeplitz)
e.g., MATLAB >> display_recon('mriData/64x64x1/out.file')
- display_recon2.m (Brute Force)
e.g., MATLAB >> display_recon2('mriData/64x64x1/')
4.5 Two more tutorials including code samples in the code package
The IMPATIENT package provides two code demonstrations under
matlab_utils/
to help you get started on the path of reconstructing your own data
with IMPATIENT.
4.5.1 matlab_utils/create_sparse_matrix
This folder contains the Matlab scripts that are used to generate
the regularization sparse matrices (i.e., those c.rmtx files in each
data sub-directory, residing along with other input *.dat files).
The sparse matrices are first generated in memory using Matlab script
createDWithPeriodicBoundary.m, then written to hard disk via mmwrite.m.
c.rmtx is a format invented by Matrix Market. For more details, please
refer to the README file under that directory.
4.5.2 matlab_utils/format_conversion
This folder demonstrates an example of converting into IMPATIENT data format,
an MRI data stored in Matlab's mat file (i.e., recon_data_32x32x4.mat).
Hopefully, this example will save our users' the pain to figure out the hard
way how to convert their favorite MRI data format into the IMPATIENT's format.
For more details, please refer to the README file under that directory.
------------------------------ Chapter separator ------------------------------
5. How to integrate it with your project
The main function in the file, main.cpp, is designed as an example for your
reference. The following code segment is the launch point.
result = mainMriSolver(input_dir, output_dir, iter_num, ...);
if (result == false) {
cerr<< "Error: Failed to execute the "<< prog<< ".\n";
goto _ERROR_mri_main;
}
------------------------------ Chapter separator ------------------------------
6. Todo list and future work
This section lists our ongoing tasks and future work.
Ongoing tasks:
- Further optimize the performance (e.g., thread coarsening).
- Improve the multi-gpu job dispatch queue and distribute the per-coil
computation across many GPUs.
- Extend the gridding kernels to a tiled implementation so that gridding
can be done on larger data sets, whose size would not fit in GPU at once.
- Thorough comparison with existing methods (e.g., NUFFT)
Future work:
- Automatically detect the underlying GPU hardware capability and
choose the right parameters.
------------------------------ End of README File -----------------------------
About
Fast MRI reconstruction on CUDA GPUs
Resources
Stars
Watchers
Forks
Releases
No releases published
Packages 0
No packages published
Languages
- C 43.4%
- C++ 30.1%
- Cuda 22.3%
- Makefile 3.6%
- Other 0.6%