Anatomy Processing

This recipe gets you from DICOM files to mesh head models and mesh ROIs.

The examples will be for subject skeri9999. The functions in the following pipeline assume ID codes in the form of skeri####

NOTE: I'm using linux-style path strings in all these examples. If you're running on another platform you'll need to use appropriate path strings for the Matlab commands e.g. on a PC '/raid/MRI/data/...' should be replaced with 'X:\data\...' or something similar.

Getting ready

Linux commands will in green.

Set up Freesurfer version 4.0 or later. There are various versions on mono in /raid/MRI/toolbox/MGH/ if you don't have it locally. Choose one appropriate for your OS whether 32- or 64-bit.

    • setenv FREESURFER_HOME /raid/MRI/toolbox/MGH/fs440_x64

    • source $FREESURFER_HOME/FreeSurferEnv.csh

Set up FSL version 4.1 or later. Right now the only FSL 4.1 on mono is in /raid/MRI/toolbox/FSL/fsl_centos5_x64/fsl/

    • setenv FSLDIR /raid/MRI/toolbox/FSL/fsl_centos5_x64/fsl

    • source $FSLDIR/etc/fslconf/fsl.csh

    • setenv PATH $FSLDIR/bin:$PATH

Set up MNE Suite

    • setenv MNE_ROOT /raid/MRI/toolbox/MNESuite/mne_linux_x86_64

    • source $MNE_ROOT/bin/mne_setup

Matlab commands will be in red

Vistasoft, the MNE toolbox, and Spero's "HeadModels" toolbox need to be on your Matlab path.

Assuming /raid/MRI/toolbox/VISTASOFT_DEVEL/matlab_scripts is already on your path, you can run


    • addMNE

    • addSperoToolbox

Othwise run

    • addpath(genpath('/raid/MRI/toolbox/VISTASOFT_DEVEL/trunk'),0)

    • addpath('/raid/MRI/toolbox/MNESuite/mne/matlab',0)

    • addpath('/raid/MRI/toolbox/SperoToolbox/HeadModels',0)

Convert DICOM files to NIfTI images

Create anatomy directories. Put skeri9999's anatomical scan data in the RawDicom directory.

    • mkdir /raid/MRI/anatomy/skeri9999

    • mkdir /raid/MRI/anatomy/skeri9999/nifti

    • mkdir /raid/MRI/anatomy/skeri9999/RawDicom

Convert Dicom files directly to nii.gz with Freesurfer's mri_convert. This example assumes two T1 scans and one T2 scan were collected on 1/1/09 as we typically collect at UCSF's NIC.

    • cd /raid/MRI/anatomy/skeri9999/nifti

    • mri_convert ../RawDicom/t1_mprage08iso_3/IM-0002-0001.dcm skeri9999_010109_T1-3.nii.gz

    • mri_convert ../RawDicom/t1_mprage08iso_4/IM-0003-0001.dcm skeri9999_010109_T1-4.nii.gz

    • mri_convert ../RawDicom/t2_spc_ns_sag_p2_iso_5/IM-0004-0001.dcm skeri9999_010109_T2-5.nii.gz

Reorient NIfTI images LAS RADIOLOGICAL for FSL compatibility

FSL's routines require a left, anterior, superior (LAS) dimension ordering. Check the orientation in fslview and fslorient, and adjust as needed with fslswapdim and/or fslorient...

    • fslswapdim skeri9999_010109_T1-3.nii.gz ±? ±? ±? skeri9999_010109_T1-3.nii.gz

    • fslorient -swaporient...

When correctly oriented:

    • fslorient -getorient skeri9999_010109_T1-3.nii.gz

will return RADIOLOGICAL, and the volume will be oriented like this in fslview .

Our current anatomy prescriptions at NIC have posterior,inferior,left dimension ordering, thus the following commands convert to LAS. This will not necessarily be the same for data from different scanners, or if the prescription is changed!

    • fslswapdim skeri9999_010109_T1-3.nii.gz z -x -y skeri9999_010109_T1-3.nii.gz

    • fslswapdim skeri9999_010109_T1-4.nii.gz z -x -y skeri9999_010109_T1-4.nii.gz

    • fslswapdim skeri9999_010109_T2-5.nii.gz z -x -y skeri9999_010109_T2-5.nii.gz

Test for radiological orientation with

    • fslorient -getorient skeri9999_010109_T1-3.nii.gz

Average NIfTI images [Optional]

You can input multiple T1 images to Freesurfer, and it will do between-scan motion corrections and averaging prior to starting its segmentation pipeline. However, this happens after resampling to 1x1x1mm voxels. To average 2 volumes at the native T1 resolution, e.g. 0.8mm cubic voxels at NIC, use

    • cd /raid/MRI/anatomy/skeri9999/nifti

    • /MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_010109_T1-3 skeri9999_010109_T1-4 skeri9999_T1-avg [angle] [cost]

The optional inputs angle (default = 15), and cost (default = mutualinfo) specify the search ranges for rotations, and the cost algorithm for FSL's flirt. See flirt help for more info.

Import T1 image to Freesurfer4 and run autorecon

    • recon-all -i skeri9999_010109_T1-avg.nii -subjid skeri9999_fs4

    • recon-all -all -subjid skeri9999_fs4

If you're in a hurry to get a mrGray segmentation as fast as possible, Freesurfer divides the autorecon process in 3 stages, and the files you'll need for mrGray are generated in the 2nd stage. The recon-all -all... step above can be replaced with

    • recon-all -autorecon1 -autorecon2 -subjid skeri9999_fs4

    • recon-all -autorecon3 -subjid skeri9999_fs4

and you can proceed to segmentation when the first 2 stages are complete and the 3rd autorecon stage is still running

See the Freesurfer Wiki for more info.

Convert Freesurfer4 volumes to .nii.gz & reorient to LAS for FSL compatibility

Freesurfer's volumes are oriented so that the +X, +Y, and +Z dimensions correspond to left, inferior, and anterior. FSL's routines require LAS dimension ordering. The script converts Freesurfer mgz-files to .nii.gz format using mri_convert, then reorients to LAS using FSL's fslswapdim. It is called as follows:

    • cd /MRI/anatomy/skeri9999/nifti

    • /MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_fs4 skeri9999

It will the following files in your current working directory





To verify that +X=Left, +Y=Anterior, and +Z=Superior in the *_FS4_*.nii.gz files, open in fslview and move the cursor around and look at the numbers in the X,Y,Z boxes. Other MRI volume viewing software such as MRIcron might change the ordering and/or directions of the axes, I wouldn't use it for this purpose.

To verify that the NIfTI headers all indicate a RADIOLOGICAL format, use fslorient e.g.

    • fslorient -getorient skeri9999_FS4_orig.nii.gz

If there are any problems see the usage info on fslswapdim and fslorient to get your NIfTI files coverted to RADIOLOGICAL LAS format.

Note: Recently (FS5? New version of FS? VISTASOFT_GIT? ) we have found that the resulting anatomies are not correctly oriented. To fix this we do....


fslswapdim R3074_FS4_brain x z -y out3

The other thing you might have to do is delete all the .nii files (but leave the .nii.gz). FSL can't operate if both the .nii and .nii.gz versions of the file are in the same directory.

You can get around this by setting your default output type in .cshrc to .nii.gz

Create a vAnatomy.dat and class-files for mrGray

createVolAnatSKERI.m and fs4_ribbon2class.m can both be called without inputs allowing choice of input and output files via GUI dialogs. Both of these m-files can also detect a volume's orientation so this step can be run before or after the reorientation to LAS described below.

    • cd /raid/MRI/anatomy/skeri9999/nifti

    • createVolAnatSKERI('skeri9999_FS4_nu.nii.gz')

    • fs4_ribbon2class('skeri9999_FS4_ribbon.nii.gz')

Create a skeri9999_fiducials.txt file for aligning Polhemus elp-files. After choosing the nasion and left/right auricular locations, save the result in the Freesurfer bem directory

    • vAnatomyFiducials('/raid/MRI/anatomy/skeri9999/vAnatomy.dat')

Segment in mrGray

You'll need to create a new project for each hemisphere, and load the vAnatomy.dat and class-files created in the previous step. The gray matter in the class file should be deleted and regrown in mrGray. Freesurfer's segmentation should have few if any satellite volumes, cavities, and/or handles. It does usually have a bit of the optic nerve that I like to edit out by hand.

See Alex's video tutorial

Once segmentation is complete you can install it in Vistasoft and build Gray and Flat views, and define ROIs etc.

save /raid/MRI/anatomy/FREESURFER_SUBS/skeri9999_fs4/bem/skeri9999_fiducials.txt

save /raid/MRI/anatomy/skeri9999/left/left.Class

and /raid/MRI/anatomy/skeri9999/right/right.Class

save output as /raid/MRI/anatomy/skeri9999/vAnatomy.dat

    • make stuctural and white matter classification nifti files in the subject's main level anatomy directory

$ mri_convert $SUBJECTS_DIR/subjid/mri/T1.mgz t1.nii.gz $ fslswapdim t1 -x z -y t1 $ fslorient -swaporient t1 $ mri_convert $SUBJECTS_DIR/subjid/mri/ribbon.mgz t1_class.nii.gz $ fslswapdim t1_class -x z -y t1_class $ fslorient -swaporient t1_class $ fslmaths t1_class -thr 2.5 -min 1 -mul 38 t1_temp $ fslmaths t1_class -sub t1_temp -add 1 -thr 2.5 -uthr 4.5 t1_class $ rm t1_temp.nii.gz

    • Fix satellite volumes, cavities, and handles in itkgray and save the corrected t1_class.nii.gz

or, try newer ITKGray segmentation pathway...

NOTE: mrVista alignments may crash when using t1.nii.gz instead of a vAnatomy.dat, and making unfolds might crash mex files when using t1_class.nii.gz instead of *.MrM files from mrGray, so you might have to replicate the whole mrGray pipeline anyway, but the above works and will yield a segmentation in the same voxel coords as Freesurfer. Itkgray can save VTK, STL, or BYU meshes, but I don't see how to use these in mrVista.

    • Convert t1.nii.gz into vAnatomy.dat.

    • createVolAnat.m is supposed to do this, but it has a bug.

>> anat = mrLoad('t1.nii.gz','nifti'); >> = mrAnatRotateAnalyze(; >> mrSave(anat,'vAnatomy.dat','vanat')

    • Convert t1_class.nii.gz into *.Class files

>> itkGray2mrGrayClass

Align the T2 volume with the T1

This is a multi-step process making use of various FSL tools that has been combined into a single script. The steps are:

    1. extract the brain from the T2 volume

    2. use the brain mask to get a non-bias-corrected T2 brain from the original volume

    3. estimate the bias correction field from the T2 brain

    4. correct the whole T2 volume using this estimate

    5. estimate the rigid-body transform to align the T2 brain with Freesurfer's T1 brain

    6. transform the T2 volume into alignment.

The script performs this alignment and is called as follows:

    • /raid/MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_010109_T2-5.nii.gz skeri9999_FS4_brain.nii.gz 0.3

The 3rd input argument is the brain extraction level. If you get poor T1-T2 alignment (you can overlay these in fslview to check) look at your T2 brain extraction and experiment with this parameter (see bet help).

In detail, does the following

    • bet skeri9999_010109_T2-5 T2_brain -B -f 0.3

    • fslmaths skeri9999_010109_T2-5 -mas T2_brain_mask T2_brain

    • fast -t 2 -n 4 -l 15 -b --nopve T2_brain

    • fslmaths -dt float skeri9999_010109_T2-5 -div T2_brain_bias T2_unbias

    • flirt -in T2_brain -ref skeri9999_FS4_brain -dof 6 -usesqform -omat T2_align.txt

    • flirt -in T2_unbias -ref skeri9999_FS4_brain -out T2_FINAL -applyxfm -init T2_align.txt -interp sinc

Make non-cortical meshes

In preparation open the skeri9999_FS4_orig.nii.gz volume in fslview and get the coordinates of the brain center: Cx, Cy, Cz (voxel coordinates). Choose something near the location shown below.

Using the script

    1. get a T2 brain mesh with bet2

    2. get the alignment of T1 to a standard brain

    3. make meshes with betsurf

    • /raid/MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_FS4_orig T2_YYYYMMDD_HHMMPM_FINAL Cx Cy Cz

The individual FSL commands combined in getBEMsurfs are:

    • bet2 T2_YYYYMMDD_HHMMPM_FINAL headmodel_T2brain -e -m -c Cx Cy Cz

    • flirt -usesqform -ref $FSLDIR/data/standard/MNI152_T1_1mm -in skeri9999_FS4_orig -out headmodel_T1ref -omat headmodel_T1ref.txt

    • betsurf -m -s skeri9999_FS4_orig T2_YYYYMMDD_HHMMPM_FINAL headmodel_T2brain_mesh.vtk headmodel_T1ref.txt headmodel

Convert,, and to tri-files in the Freesurfer bem directory for MNE suite

    • FSLoff2tri('/raid/MRI/anatomy/skeri9999/nifti/headmodel','skeri9999_fs4')

Setup MNE source space and BEM models

Make lh and rh midgray (halfway between white and pial) Freesurfer surfaces

    • makeFreesurferMidgraySurface('skeri9999_fs4')

Setup an MNE source space file for the midgray layer, then set up the BEM model

    • cd /raid/MRI/anatomy/FREESURFER_SUBS/skeri9999_fs4/bem

    • /raid/MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_fs4

    • /raid/MRI/toolbox/SperoToolbox/HeadModels/ skeri9999_fs4

Make a high-res scalp (/surf/lh.smseghead) in Freesurfer

    • mkheadsurf -subjid skeri9999_fs4 -srcvol orig.mgz -thresh1 30 -thresh2 1

Check it out and make sure it looks like a head you could register surface data to. Any of the following will do, and each have their pluses and minuses. If there are problems look at the help and try adjusting the thresh1 and thresh2 options to mkheadsurf. 20 is the default for both thresholds, using tkmedit below you can read voxel values and pick a thresh1 high enough to not get fooled by junk outside the head.

    • tkmedit skeri9999_fs4 orig.mgz lh.smseghead

    • scuba -s skeri9999_fs4 -v orig.mgz -f lh.smseghead

    • freeview -v ../mri/orig.mgz -f ../surf/lh.smseghead

Convert the lh.smseghead scalp to MNE format. If you have problems, try recreating the scalp with a low thresh2, like 1.

    • mne_surf2bem --surf ../surf/lh.smseghead --id 4 --check --fif skeri9999_fs4-head.fif

NOTE: In at least one installation of Freesurfer 4.5.0, mkheadsurf has a bug where it still does smoothing but overwrites lh.seghead instead of creating lh.smseghead. Use the former if the latter doesn't exist.

Create decimated mrMesh-format cortex

Use FS4toDefaultCortex.m. This m-file makes use of mrMesh server, and isn't currently running on wocket. Try a windows machine. We're switching to MNE decimations which have more uniform faces. Set 2nd flag input to true (default), otherwise it'll do a reducepatch.mdecimation of the pial surface.

Syntax: FS4toDefaultCortex(FS4subjid,MNEdecimationFlag)

    • FS4toDefaultCortex('skeri9999_fs4',true)

    • Before saving the result, click on the mesh window and rotate the cortex checking for any obvious problems with the midgray surface.

    • Save as /raid/MRI/anatomy/skeri9999/Standard/meshes/defaultCortex.mat. If you don't have permissions to write in the Standard folder, save it somewhere else under the default name in the GUI dialog, and get somebody with permissions to rename and move/copy it to the right location.

    • After saving, close the mrMesh server window.

Create mesh ROIs

Convert the Desikan ROIs from Freesurfer 4's parcellations to mesh-index format

    • FS4parc2cortex('skeri9999_fs4')

Once you have defined functional gray ROIs in VistaSoft, and saved them to the shared directory for a given subject, you can then proceed to map them to the cortex mesh. These are necessary for source-localized EEG projects.

The matlab command to generate mesh ROIs is in Spero's HeadModels toolbox

    • MLRrois2mesh('/raid/MRI/anatomy/skeri9999/vAnatomy.dat')

The rest is rather user friendly. This command will pop up a window asking you to specify which Gray ROIs you want to map onto the mesh. The program should automatically find the gray ROIs associated with the vAnatomy file, asusming you have saved the gray ROIs to the appropriate folder (MRI/anatomy/skeri####/Standard/Gray/ROIs). Last but not least, it will next pop open a window asking you to chose the cortex mesh mat file. Again, it automatically looks for the corresponding file for the given subject in MRI/anatomy/skeri####/Standard/meshes

Two figures appear, one showing the location of the specified ROIs onto the mesh, and the second with histograms of the distances from VistaSoft Gray nodes in each ROI to the nearest decimated cortex vertex. Look to see if you are happy with the anatomical location of the ROI, and that the histogram shows most ROI nodes less than 1 mm or so from their mapped vertex (the former is more of an issue these days, and latter not really so much). If you like what you see, you can save these ROIs to the mesh by choosing SAVE ROIs under the HemiViews menu option on the anatomical figure.