This program calculates the depth of ROIs, masks and 'background', using the fun Euclidean Distance Transform (EDT). Basically, this means calculating the Euclidean distance of each voxel's centroid to the nearest boundary with a separate ROI (well, to be brutally technical, to centroid of the nearest voxel in a neighboring ROI. The input dataset should be a map of ROIs (so, integer-valued). The EDT values are calculated throughout the entire FOV by default, even in the zero/background regions (there is an option to control this). written by: PA Taylor and P Lauren (SSCC, NIMH, NIH)
This code calculates the Euclidean Distance Transform (EDT) for 3D volumes following this nice, efficient algorithm, by Felzenszwalb and Huttenlocher (2012; FH2012): Felzenszwalb PF, Huttenlocher DP (2012). Distance Transforms of Sampled Functions. Theory of Computing 8:415-428. https://cs.brown.edu/people/pfelzens/papers/dt-final.pdf Thanks to C. Rorden for pointing this paper out and discussing it. The current code here extends/tweaks the FH2012 algorithm to a more general case of having several different ROIs present, for running in 3D (trivial extension), and for having voxels of non-unity and non-isotropic lengths. It does this by utilizing the fact that at its very heart, the FH2012 algorithm works line by line and can even be thought of as working boundary-by-boundary. Here, the zero-valued 'background' is also just treated like an ROI, with one difference. At a FOV boundary, the zero-valued ROI/backgroud is treated as open, so that the EDT value at each 'zero' voxel is always to one of the shapes within the FOV. For nonzero ROIs, one can treat the FOV boundary *either* as an ROI edge (EDT value there will be 1 edge length) *or* as being open.
3dDepthMap [options] -prefix PREF -input DSET where: -input DSET :(req) input dataset -prefix PREF :(req) output prefix name -mask MASK :mask dataset. NB: this mask is only applied *after* the EDT has been calculated. Therefore, the boundaries of this mask have no affect on the calculated distance values, except for potentially zeroing some out at the end. -dist_sq :by default, the output EDT volume contains distance values. By using this option, the output values are distance**2. -ignore_voxdims :this EDT algorithm works in terms of physical distance and uses the voxel dimension info in each direction, by default. However, using this option will ignore voxel size, producing outputs as if each voxel dimension was unity. -rimify RIM :instead of outputting a depthmap for each ROI, output a map of each ROI's 'rim' voxels---that is, the boundary layer or periphery up to thickness RIM---if RIM>0. + Note that RIM is applied to whatever kind of depth information you are calculating: if you use '-dist_sq' then the voxel's distance-squared value to the ROI edge is compared with RIM; if using '-ignore_voxdims', then the number-of-voxels to the edge is compared with RIM. The depthmap thresholding is applied as: abs(DEPTH)<=RIM. + When using this opt, any labeltable/atlastable from the original should be passed along, as well. + A negative RIM value inverts the check, and the output is kept if the depth info is: abs(DEPTH)>=abs(RIM). NB: with a negative RIM value, it is possible an ROI could disappear! -zeros_are_zero :by default, EDT values are output for the full FOV, even zero-valued regions. If this option is used, EDT values are only reported within the nonzero locations of the input dataset. -zeros_are_neg :if this option is used, EDT in the zero/background of the input will be negative (def: they are positive). This opt cannot be used if '-zeros_are_zero' is. -nz_are_neg :if this option is used, EDT in the nonzero ROI regions of the input will be negative (def: they are positive). -bounds_are_not_zero :this flag affects how FOV boundaries are treated for nonzero ROIs: by default, they are viewed as ROI boundaries (so the FOV is a closed boundary for an ROI, as if the FOV were padded by an extra layer of zeros); but when this option is used, the ROI behaves as if it continued 'infinitely' at the FOV boundary (so it is an open boundary). Zero-valued ROIs (= background) are not affected by this option. -only2D SLI :instead of running full 3D EDT, run just in 2D, per. plane. Provide the slice plane you want to run along as the single argument SLI: "axi" -> for axial slice "cor" -> for coronal slice "sag" -> for sagittal slice -binary_only :if the input is a binary mask or should be treated as one (all nonzero voxels -> 1; all zeros stay 0), then using this option will speed up the calculation. See Notes below for more explanation of this. NOT ON YET! -verb V :manage verbosity when running code (def: 1). Providing a V of 0 means to run quietly.
The original EDT algorithm of FH2012 was developed for a simple binary mask input (and actually for homogeneous data grids of spacing=1). This program, however, was built to handle more generalized cases of inputs, namely ROI maps (and arbitrary voxel dimensions). The tradeoff of the expansion to handling ROI maps is an increase in processing time---the original binary-mask algorithm is *very* efficient, and the generalized one is still pretty quick but less so. So, if you know that your input should be treated as a binary mask, then you can use the '-binary_only' option to utilize the more efficient (and less generalized) algorithm. The output dataset should be the same in either case---this option flag is purely about speed of computation. All other options about outputting dist**2 or negative values/etc. can be used in conjunction with the '-binary_only', too.
1) Basic case: 3dDepthMap \ -input roi_map.nii.gz \ -prefix roi_map_EDT.nii.gz 2) Same as above, but only output distances within nonzero regions/ROIs: 3dDepthMap \ -zeros_are_zero \ -input roi_map.nii.gz \ -prefix roi_map_EDT_ZZ.nii.gz 3) Output distance-squared at each voxel: 3dDepthMap \ -dist_sq \ -input mask.nii.gz \ -prefix mask_EDT_SQ.nii.gz 4) Distinguish ROIs from nonzero background by making the former have negative distance values in output: 3dDepthMap \ -nz_are_neg \ -input roi_map.nii.gz \ -prefix roi_map_EDT_NZNEG.nii.gz 5) Have output voxel values represent (number of vox)**2 from a boundary; voxel dimensions are ignored here: 3dDepthMap \ -ignore_voxdims \ -dist_sq \ -input roi_map.nii.gz \ -prefix roi_map_EDT_SQ_VOX.nii.gz 6) Basic case, with option for speed-up because the input is a binary mask (i.e., only ones and zeros); any of the other above options can be combined with this, too: 3dDepthMap \ -binary_only \ -input roi_mask.nii.gz \ -prefix roi_mask_EDT.nii.gz 7) Instead of outputting ROI depth, output a map of the ROI rims, keeping: the part of the ROI up to where depth is >=1.6mm 3dDepthMap \ -input roi_map.nii.gz \ -rimify 1.6 \ -prefix roi_map_rim.nii.gz