|
| 1 | +""" |
| 2 | + Change directory to provide relative paths for doctests |
| 3 | + >>> import os |
| 4 | + >>> filepath = os.path.dirname( os.path.realpath( __file__ ) ) |
| 5 | + >>> datadir = os.path.realpath(os.path.join(filepath, '../../testing/data')) |
| 6 | + >>> os.chdir(datadir) |
| 7 | +
|
| 8 | +""" |
| 9 | +import os |
| 10 | + |
| 11 | +from nipype.interfaces.base import (CommandLineInputSpec, CommandLine, traits, |
| 12 | + TraitedSpec, File, StdOutCommandLine, |
| 13 | + StdOutCommandLineInputSpec, isdefined) |
| 14 | +from nipype.utils.filemanip import split_filename |
| 15 | + |
| 16 | +class SFPICOCalibDataInputSpec(StdOutCommandLineInputSpec): |
| 17 | + snr = traits.Float(argstr='-snr %f', units='NA', |
| 18 | + desc=('Specifies the signal-to-noise ratio of the ' |
| 19 | + 'non-diffusion-weighted measurements to use in simulations.')) |
| 20 | + scheme_file = File(exists=True, argstr='-schemefile %s', mandatory=True, |
| 21 | + desc='Specifies the scheme file for the diffusion MRI data') |
| 22 | + info_file = File(desc='The name to be given to the information output filename.', |
| 23 | + argstr='-infooutputfile %s', mandatory=True, genfile=True, |
| 24 | + hash_files=False) # Genfile and hash_files? |
| 25 | + trace = traits.Float(argstr='-trace %f', units='NA', |
| 26 | + desc='Trace of the diffusion tensor(s) used in the test function.') |
| 27 | + onedtfarange = traits.List(traits.Float, argstr='-onedtfarange %s', |
| 28 | + minlen=2, maxlen=2, units='NA', |
| 29 | + desc=('Minimum and maximum FA for the single tensor ' |
| 30 | + 'synthetic data.')) |
| 31 | + onedtfastep = traits.Float(argstr='-onedtfastep %f', units='NA', |
| 32 | + desc=('FA step size controlling how many steps there are ' |
| 33 | + 'between the minimum and maximum FA settings.')) |
| 34 | + twodtfarange = traits.List(traits.Float, argstr='-twodtfarange %s', |
| 35 | + minlen=2, maxlen=2, units='NA', |
| 36 | + desc=('Minimum and maximum FA for the two tensor ' |
| 37 | + 'synthetic data. FA is varied for both tensors ' |
| 38 | + 'to give all the different permutations.')) |
| 39 | + twodtfastep = traits.Float(argstr='-twodtfastep %f', units='NA', |
| 40 | + desc=('FA step size controlling how many steps there are ' |
| 41 | + 'between the minimum and maximum FA settings ' |
| 42 | + 'for the two tensor cases.')) |
| 43 | + twodtanglerange = traits.List(traits.Float, argstr='-twodtanglerange %s', |
| 44 | + minlen=2, maxlen=2, units='NA', |
| 45 | + desc=('Minimum and maximum crossing angles ' |
| 46 | + 'between the two fibres.')) |
| 47 | + twodtanglestep = traits.Float(argstr='-twodtanglestep %f', units='NA', |
| 48 | + desc=('Angle step size controlling how many steps there are ' |
| 49 | + 'between the minimum and maximum crossing angles for ' |
| 50 | + 'the two tensor cases.')) |
| 51 | + twodtmixmax = traits.Float(argstr='-twodtmixmax %f', units='NA', |
| 52 | + desc=('Mixing parameter controlling the proportion of one fibre population ' |
| 53 | + 'to the other. The minimum mixing parameter is (1 - twodtmixmax).')) |
| 54 | + twodtmixstep = traits.Float(argstr='-twodtmixstep %f', units='NA', |
| 55 | + desc=('Mixing parameter step size for the two tensor cases. ' |
| 56 | + 'Specify how many mixing parameter increments to use.')) |
| 57 | + seed = traits.Float(argstr='-seed %f', units='NA', |
| 58 | + desc='Specifies the random seed to use for noise generation in simulation trials.') |
| 59 | + |
| 60 | +class SFPICOCalibDataOutputSpec(TraitedSpec): |
| 61 | + PICOCalib = File(exists=True, desc='Calibration dataset') |
| 62 | + calib_info = File(exists=True, desc='Calibration dataset') |
| 63 | + |
| 64 | +class SFPICOCalibData(StdOutCommandLine): |
| 65 | + """ |
| 66 | + Generates Spherical Function PICo Calibration Data. |
| 67 | + |
| 68 | + SFPICOCalibData creates synthetic data for use with SFLUTGen. The |
| 69 | + synthetic data is generated using a mixture of gaussians, in the |
| 70 | + same way datasynth generates data. Each voxel of data models a |
| 71 | + slightly different fibre configuration (varying FA and fibre- |
| 72 | + crossings) and undergoes a random rotation to help account for any |
| 73 | + directional bias in the chosen acquisition scheme. A second file, |
| 74 | + which stores information about the datafile, is generated along with |
| 75 | + the datafile. |
| 76 | +
|
| 77 | + Example 1 |
| 78 | + --------- |
| 79 | + To create a calibration dataset using the default settings |
| 80 | + |
| 81 | + >>> import nipype.interfaces.camino as cam |
| 82 | + >>> calib = cam.SFPICOCalibData() |
| 83 | + >>> calib.inputs.scheme_file = 'A.scheme' |
| 84 | + >>> calib.inputs.snr = 20 |
| 85 | + >>> calib.inputs.info_file = 'PICO_calib.info' |
| 86 | + >>> calib.run() # doctest: +SKIP |
| 87 | + |
| 88 | + The default settings create a large dataset (249,231 voxels), of |
| 89 | + which 3401 voxels contain a single fibre population per voxel and |
| 90 | + the rest of the voxels contain two fibre-populations. The amount of |
| 91 | + data produced can be varied by specifying the ranges and steps of |
| 92 | + the parameters for both the one and two fibre datasets used. |
| 93 | + |
| 94 | + Example 2 |
| 95 | + --------- |
| 96 | + To create a custom calibration dataset |
| 97 | + |
| 98 | + >>> import nipype.interfaces.camino as cam |
| 99 | + >>> calib = cam.SFPICOCalibData() |
| 100 | + >>> calib.inputs.scheme_file = 'A.scheme' |
| 101 | + >>> calib.inputs.snr = 20 |
| 102 | + >>> calib.inputs.info_file = 'PICO_calib.info' |
| 103 | + >>> calib.inputs.twodtfarange = [0.3, 0.9] |
| 104 | + >>> calib.inputs.twodtfastep = 0.02 |
| 105 | + >>> calib.inputs.twodtanglerange = [0, 0.785] |
| 106 | + >>> calib.inputs.twodtanglestep = 0.03925 |
| 107 | + >>> calib.inputs.twodtmixmax = 0.8 |
| 108 | + >>> calib.inputs.twodtmixstep = 0.1 |
| 109 | + >>> calib.run() # doctest: +SKIP |
| 110 | + |
| 111 | + This would provide 76,313 voxels of synthetic data, where 3401 voxels |
| 112 | + simulate the one fibre cases and 72,912 voxels simulate the various |
| 113 | + two fibre cases. However, care should be taken to ensure that enough |
| 114 | + data is generated for calculating the LUT. # doctest: +SKIP |
| 115 | + """ |
| 116 | + _cmd = 'sfpicocalibdata' |
| 117 | + input_spec=SFPICOCalibDataInputSpec |
| 118 | + output_spec=SFPICOCalibDataOutputSpec |
| 119 | + |
| 120 | + def _list_outputs(self): |
| 121 | + outputs = self.output_spec().get() |
| 122 | + outputs['PICOCalib'] = os.path.abspath(self._gen_outfilename()) |
| 123 | + outputs['calib_info'] = os.path.abspath(self.inputs.info_file) |
| 124 | + return outputs |
| 125 | + |
| 126 | + def _gen_outfilename(self): |
| 127 | + _, name , _ = split_filename(self.inputs.scheme_file) |
| 128 | + return name + '_PICOCalib.Bfloat' |
| 129 | + |
| 130 | + |
| 131 | + |
| 132 | +class SFLUTGenInputSpec(StdOutCommandLineInputSpec): |
| 133 | + in_file = File(exists=True, argstr='-inputfile %s', mandatory=True, |
| 134 | + desc='Voxel-order data of the spherical functions peaks.') |
| 135 | + info_file = File(argstr='-infofile %s', mandatory=True, |
| 136 | + desc=('The Info file that corresponds to the calibration ' |
| 137 | + 'datafile used in the reconstruction.')) |
| 138 | + outputstem = traits.Str('LUT', argstr='-outputstem %s', |
| 139 | + desc=('Define the name of the generated luts. The form of the filenames will be ' |
| 140 | + '[outputstem]_oneFibreSurfaceCoeffs.Bdouble and ' |
| 141 | + '[outputstem]_twoFibreSurfaceCoeffs.Bdouble'), |
| 142 | + usedefault=True) |
| 143 | + pdf = traits.Enum('bingham', 'watson', argstr='-pdf %s', |
| 144 | + desc=('Sets the distribution to use for the calibration. The default is the Bingham ' |
| 145 | + 'distribution, which allows elliptical probability density contours. ' |
| 146 | + 'Currently supported options are: ' |
| 147 | + ' bingham - The Bingham distribution, which allows elliptical probability ' |
| 148 | + ' density contours. ' |
| 149 | + ' watson - The Watson distribution. This distribution is rotationally symmetric.'), |
| 150 | + usedefault=True) |
| 151 | + binincsize = traits.Int(argstr='-binincsize %d', units='NA', |
| 152 | + desc=('Sets the size of the bins. In the case of 2D histograms such as the ' |
| 153 | + 'Bingham, the bins are always square. Default is 1.')) |
| 154 | + minvectsperbin = traits.Int(argstr='-minvectsperbin %d', units='NA', |
| 155 | + desc=('Specifies the minimum number of fibre-orientation estimates a bin ' |
| 156 | + 'must contain before it is used in the lut line/surface generation. ' |
| 157 | + 'Default is 50. If you get the error "no fibre-orientation estimates ' |
| 158 | + 'in histogram!", the calibration data set is too small to get enough ' |
| 159 | + 'samples in any of the histogram bins. You can decrease the minimum ' |
| 160 | + 'number per bin to get things running in quick tests, but the sta- ' |
| 161 | + 'tistics will not be reliable and for serious applications, you need ' |
| 162 | + 'to increase the size of the calibration data set until the error goes.')) |
| 163 | + directmap = traits.Bool(argstr='-directmap', |
| 164 | + desc=('Use direct mapping between the eigenvalues and the distribution parameters ' |
| 165 | + 'instead of the log of the eigenvalues.')) |
| 166 | + order = traits.Int(argstr='-order %d', units='NA', |
| 167 | + desc=('The order of the polynomial fitting the surface. Order 1 is linear. ' |
| 168 | + 'Order 2 (default) is quadratic.')) |
| 169 | + |
| 170 | +class SFLUTGenOutputSpec(TraitedSpec): |
| 171 | + lut_one_fibre = File(exists=True, desc='PICo lut for one-fibre model') |
| 172 | + lut_two_fibres = File(exists=True, desc='PICo lut for two-fibre model') |
| 173 | + |
| 174 | +class SFLUTGen(StdOutCommandLine): |
| 175 | + """ |
| 176 | + Generates PICo lookup tables (LUT) for multi-fibre methods such as |
| 177 | + PASMRI and Q-Ball. |
| 178 | + |
| 179 | + SFLUTGen creates the lookup tables for the generalized multi-fibre |
| 180 | + implementation of the PICo tractography algorithm. The outputs of |
| 181 | + this utility are either surface or line coefficients up to a given |
| 182 | + order. The calibration can be performed for different distributions, |
| 183 | + such as the Bingham and Watson distributions. |
| 184 | + |
| 185 | + This utility uses calibration data generated from SFPICOCalibData |
| 186 | + and peak information created by SFPeaks. |
| 187 | + |
| 188 | + The utility outputs two lut's, *_oneFibreSurfaceCoeffs.Bdouble and |
| 189 | + *_twoFibreSurfaceCoeffs.Bdouble. Each of these files contains big- |
| 190 | + endian doubles as standard. The format of the output is: |
| 191 | + dimensions (1 for Watson, 2 for Bingham) |
| 192 | + order (the order of the polynomial) |
| 193 | + coefficient_1 |
| 194 | + coefficient_2 |
| 195 | + ... |
| 196 | + coefficient_N |
| 197 | + In the case of the Watson, there is a single set of coefficients, |
| 198 | + which are ordered: |
| 199 | + constant, x, x^2, ..., x^order. |
| 200 | + In the case of the Bingham, there are two sets of coefficients (one |
| 201 | + for each surface), ordered so that: |
| 202 | + for j = 1 to order |
| 203 | + for k = 1 to order |
| 204 | + coeff_i = x^j * y^k |
| 205 | + where j+k < order |
| 206 | +
|
| 207 | + Example |
| 208 | + --------- |
| 209 | + To create a calibration dataset using the default settings |
| 210 | + |
| 211 | + >>> import nipype.interfaces.camino as cam |
| 212 | + >>> lutgen = cam.SFLUTGen() |
| 213 | + >>> lutgen.inputs.in_file = 'QSH_peaks.Bdouble' |
| 214 | + >>> lutgen.inputs.info_file = 'PICO_calib.info' |
| 215 | + >>> lutgen.run() # doctest: +SKIP |
| 216 | + """ |
| 217 | + _cmd = 'sflutgen' |
| 218 | + input_spec=SFLUTGenInputSpec |
| 219 | + output_spec=SFLUTGenOutputSpec |
| 220 | + |
| 221 | + def _list_outputs(self): |
| 222 | + outputs = self.output_spec().get() |
| 223 | + outputs['lut_one_fibre'] = self.inputs.outputstem + '_oneFibreSurfaceCoeffs.Bdouble' |
| 224 | + outputs['lut_two_fibres'] = self.inputs.outputstem + '_twoFibreSurfaceCoeffs.Bdouble' |
| 225 | + return outputs |
| 226 | + |
| 227 | + def _gen_outfilename(self): |
| 228 | + return '/dev/null' |
0 commit comments