Difference between revisions of "Matlab based temperature solver"

From MGH/MIT Parallel Transmission Resources
Jump to navigation Jump to search
Line 3: Line 3:
 
 
 
=== What this code does ===
 
=== What this code does ===
This Matlab code (no MEX, no GPU code) computes least-squares and magnitude least-squares pTx spoke pulses. The MLS problem can be solved using the phase adoption approach of Setsompop et al. or a full optimization strategy (see the file spokes_3/run_spokes.m and options therein). Spokes can be optimized for several frequencies (spatio-spectral design) in order to build in robustness to off-resonance effects.
+
This is a Crank-Nicholson Matlab-based temperature solver. The code mixed spatial and temporal finite differences using the Crank-Nicholson scheme. Non-homogeneous tissues properties (including tissue perfusion, thermal conductivity and heat capacity) are supported. The code is fully validated by comparison with the SEMCAD temperature solver. Two solve options are available: Transient temperature and steady-state calculations.
  
 
=== Download files ===
 
=== Download files ===
The code can be downloaded [https://ptx.martinos.org/images/f/f2/SpokesMinExcErrLocGlobSARPuPow_fmincon_v12.zip here]. A test dataset can be downloaded [https://ptx.martinos.org/images/e/e2/Spokes_3.zip here]. If you wonder how I obtained the field maps, ROI and SODA (=slice information) files, download [https://ptx.martinos.org/images/c/c0/Prescan_data.zip this file].
+
The code can be downloaded [https://ptx.martinos.org/index.php/File:PennesBioHeatSolver_LITE.zip here].
  
 
=== Instructions ===
 
=== Instructions ===
==== Pre-scan data format ====
+
==== Transient simulation ====
The format of the pre-scan data (e.g., B0 map, B1+ maps, slice information in the SODA file) matches the Siemens pTx "Step 2" data format. To take a closer look at this, download the prescan_data.zip file,  
+
To run a simple example, go to the sub-folder SEMCAD_validation/Simu_1_Results and run the script "semcad_vs_pbhs.m". This runs a comparison of the transient solver with SEMCAD. This script should make it clear how the code is run. The input to the code are, in this order:
* An ROI file in NIFTI format.
+
*T0: The initial temperature map (can be non-uniform) in K
* A B0 map in NIFTI format.
+
*DUR: The total duration of the simulation in s
* An FLD file containing the B1+ maps for the slice considered as obtained by pushing "save dataset" in the B1 mapping tab of the Siemens pTx "Step 2" interface.
+
*[dx dy dz]: A 1x3 Matlab vector with the x, y and z spatial resolutions in m
To get the ROI and B0 maps, I use the Siemens product gre_fieldmapping sequence (make sure to chose "save magnitude/phase" in order to get both the B0 and ROI maps in the same scan). The B0, B1 and ROI datasets must have the same number of pixels and same dimensions.
+
*RHO: The density map (can be non-uniform) in kg/m^3
==== Spokes optimization ====
+
*C: The heat capacity map (can be non-uniform) in J/kg/K
Unzip the code folder and add it to your path. Unzip the data folder and run the script run_spokes.m. There are several options in it that are straightforward to understand (e.g., LS/MLS, spatio-spectral design). The frequencies of the spatio-spectral design are expressed in Hz and are offset with respect to the Larmor frequency (i..e, +50 means Larmor frequency+50 Hz). Most of the spoke options (such as slice-thickness, slice-selection gradient, time-bandwidth product etc...) are specified in the text file spokes_def.txt and can be changed there.
+
*K: The thermal conductivity map (can be non-uniform) in W/m/K
 +
*Wb: The blood perfusion map (can be non-uniform) in kg/s/m^3
 +
*Cb: The blood heat capacity in J/kg/K
 +
*Tb: The blood temperature in K
 +
*Q: The heat rate generation map (i.e., SAR*RHO which can be non-uniform) in W/m^3
 +
 
 +
==== Steady state simulation ====
 +
To run a simple example of the steady-state simulation, go to the TEST_SteadyState/ folder and run "TEST_SS.m". Again, the way to run the code should be straightforward. The input to the steady-state function are similar to the input to the transient solve function, except obviously the duration is not needed.

Revision as of 14:10, 15 April 2016

Citation

Guérin B, Gebhardt M, Cauley S, Adalsteinsson E, Wald LL (2014). "Local specific absorption rate (SAR), global SAR, transmitter power, and excitation accuracy trade‐offs in low flip‐angle parallel transmit pulse design." Magnetic Resonance in Medicine 71(4): 1446-1457.

What this code does

This is a Crank-Nicholson Matlab-based temperature solver. The code mixed spatial and temporal finite differences using the Crank-Nicholson scheme. Non-homogeneous tissues properties (including tissue perfusion, thermal conductivity and heat capacity) are supported. The code is fully validated by comparison with the SEMCAD temperature solver. Two solve options are available: Transient temperature and steady-state calculations.

Download files

The code can be downloaded here.

Instructions

Transient simulation

To run a simple example, go to the sub-folder SEMCAD_validation/Simu_1_Results and run the script "semcad_vs_pbhs.m". This runs a comparison of the transient solver with SEMCAD. This script should make it clear how the code is run. The input to the code are, in this order:

  • T0: The initial temperature map (can be non-uniform) in K
  • DUR: The total duration of the simulation in s
  • [dx dy dz]: A 1x3 Matlab vector with the x, y and z spatial resolutions in m
  • RHO: The density map (can be non-uniform) in kg/m^3
  • C: The heat capacity map (can be non-uniform) in J/kg/K
  • K: The thermal conductivity map (can be non-uniform) in W/m/K
  • Wb: The blood perfusion map (can be non-uniform) in kg/s/m^3
  • Cb: The blood heat capacity in J/kg/K
  • Tb: The blood temperature in K
  • Q: The heat rate generation map (i.e., SAR*RHO which can be non-uniform) in W/m^3

Steady state simulation

To run a simple example of the steady-state simulation, go to the TEST_SteadyState/ folder and run "TEST_SS.m". Again, the way to run the code should be straightforward. The input to the steady-state function are similar to the input to the transient solve function, except obviously the duration is not needed.