Failure when using custom gradient deviations file

[akhanf]

Thanks for the quick reply to the last (non)-issue!

This one is stumping me a bit more -- I am running MDT on some 7T data where we have performed gradient non-linearity correction, and created a gradient deviations file in the same format as HCP data (9-dim in dim4). I have been using these datasets (with gradient deviations) in FSL dtifit/bedpost with no issues, but running into an error when trying it on MDT.

I've shared the dataset here as well (250MB tar): https://www.dropbox.com/s/veqn9vsi5fhtyas/test_dwi_grad_dev.tar?dl=0

It runs fine without the --gradient-deviation option.

But if I run the following commands,

mdt-create-protocol sub-CT01_dwi_space-T1wGC_preproc.bvec sub-CT01_dwi_space-T1wGC_preproc.bval -o sub-CT01_dwi_space-T1wGC_preproc.prtcl
mdt-model-fit NODDI sub-CT01_dwi_space-T1wGC_preproc.nii.gz sub-CT01_dwi_space-T1wGC_preproc.prtcl sub-CT01_dwi_space-T1wGC_brainmask.nii.gz --gradient-deviations sub-CT01_dwi_space-T1wGC_preproc.grad_dev.nii.gz

I get this output:

[2019-04-01 13:33:57,755] [INFO] [mdt.lib.model_fitting] [get_model_fit] - Starting intermediate optimization for generating initialization point.
[2019-04-01 13:33:57,960] [INFO] [mdt.lib.model_fitting] [_apply_user_provided_initialization_data] - Preparing model BallStick_r1 with the user provided initialization data.
[2019-04-01 13:33:57,967] [INFO] [mdt.lib.model_fitting] [fit_composite_model] - Using MDT version 0.20.3
[2019-04-01 13:33:57,967] [INFO] [mdt.lib.model_fitting] [fit_composite_model] - Preparing for model BallStick_r1
[2019-04-01 13:33:57,968] [INFO] [mdt.lib.model_fitting] [fit_composite_model] - Current cascade: ['BallStick_r1']
[2019-04-01 13:33:58,991] [INFO] [mdt.models.composite] [_prepare_input_data] - No volume options to apply, using all 103 volumes.
[2019-04-01 13:33:59,041] [INFO] [mdt.models.composite] [set_input_data] - Using the gradient deviations in the model optimization.
[2019-04-01 13:33:59,041] [INFO] [mdt.utils] [estimate_noise_std] - Trying to estimate a noise std.
[2019-04-01 13:33:59,107] [INFO] [mdt.utils] [estimate_noise_std] - Estimated global noise std 39.15896224975586.
[2019-04-01 13:33:59,108] [INFO] [mdt.lib.model_fitting] [_model_fit_logging] - Fitting BallStick_r1 model
[2019-04-01 13:33:59,108] [INFO] [mdt.lib.model_fitting] [_model_fit_logging] - The 4 parameters we will fit are: ['S0.s0', 'w_stick0.w', 'Stick0.theta', 'Stick0.phi']
[2019-04-01 13:33:59,108] [INFO] [mdt.lib.model_fitting] [fit_composite_model] - Saving temporary results in /scratch/akhanf/test_7T_mdt_grad_dev_enabled/sub-CT01/BallStick_r1/tmp_results.
[2019-04-01 13:33:59,139] [INFO] [mdt.lib.processing_strategies] [_process_chunk] - Computations are at 0.00%, processing next 100000 voxels (309970 voxels in total, 0 processed). Time spent: 0:00:00:00, time left: ? (d:h:m:s).
[2019-04-01 13:33:59,168] [INFO] [mdt.lib.processing_strategies] [_process] - Starting optimization
[2019-04-01 13:33:59,168] [INFO] [mdt.lib.processing_strategies] [_process] - Using MOT version 0.9.1
[2019-04-01 13:33:59,168] [INFO] [mdt.lib.processing_strategies] [_process] - We will use a single precision float type for the calculations.
[2019-04-01 13:33:59,168] [INFO] [mdt.lib.processing_strategies] [_process] - Using device 'CPU - Intel(R) Xeon(R) CPU E5-2683 v4 @ 2.10GHz (Intel(R) OpenCL)'.
[2019-04-01 13:33:59,169] [INFO] [mdt.lib.processing_strategies] [_process] - Using compile flags: ['-cl-denorms-are-zero', '-cl-mad-enable', '-cl-no-signed-zeros']
[2019-04-01 13:33:59,169] [INFO] [mdt.lib.processing_strategies] [_process] - We will use the optimizer Powell with default settings.
[2019-04-01 13:34:17,811] [INFO] [mdt.lib.processing_strategies] [_process] - Finished optimization
[2019-04-01 13:34:17,813] [INFO] [mdt.lib.processing_strategies] [_process] - Starting post-processing
2 errors generated.
/usr/lib/python3/dist-packages/pyopencl/__init__.py:63: CompilerWarning: Non-empty compiler output encountered. Set the environment variable PYOPENCL_COMPILER_OUTPUT=1 to see more.
  "to see more.", CompilerWarning)
Traceback (most recent call last):
  File "/usr/bin/mdt-model-fit", line 9, in <module>
    load_entry_point('mdt==0.20.3', 'console_scripts', 'mdt-model-fit')()
  File "/usr/lib/python3/dist-packages/mdt/lib/shell_utils.py", line 47, in console_script
    cls().start(sys.argv[1:])
  File "/usr/lib/python3/dist-packages/mdt/lib/shell_utils.py", line 66, in start
    self.run(args, {})
  File "/usr/lib/python3/dist-packages/mdt/cli_scripts/mdt_model_fit.py", line 168, in run
    fit_model()
  File "/usr/lib/python3/dist-packages/mdt/cli_scripts/mdt_model_fit.py", line 162, in fit_model
    use_cascaded_inits=args.use_cascaded_inits)
  File "/usr/lib/python3/dist-packages/mdt/__init__.py", line 197, in fit_model
    inits = get_optimization_inits(model_name, input_data, output_folder, cl_device_ind=cl_device_ind)
  File "/usr/lib/python3/dist-packages/mdt/__init__.py", line 89, in get_optimization_inits
    return get_optimization_inits(model_name, input_data, output_folder, cl_device_ind=cl_device_ind)
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 162, in get_optimization_inits
    return get_init_data(model_name)
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 94, in get_init_data
    fit_results = get_model_fit('BallStick_r1')
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 70, in get_model_fit
    initialization_data={'inits': get_init_data(model_name)}).run()
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 336, in run
    _, maps = self._run(self._model, self._recalculate, self._only_recalculate_last)
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 382, in _run
    apply_user_provided_initialization=not _in_recursion)
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 393, in _run_composite_model
    optimizer_options=self._optimizer_options)
  File "/usr/lib/python3/dist-packages/mdt/lib/model_fitting.py", line 465, in fit_composite_model
    return processing_strategy.process(worker)
  File "/usr/lib/python3/dist-packages/mdt/lib/processing_strategies.py", line 75, in process
    self._process_chunk(processor, chunks)
  File "/usr/lib/python3/dist-packages/mdt/lib/processing_strategies.py", line 120, in _process_chunk
    process()
  File "/usr/lib/python3/dist-packages/mdt/lib/processing_strategies.py", line 117, in process
    processor.process(chunk, next_indices=next_chunk)
  File "/usr/lib/python3/dist-packages/mdt/lib/processing_strategies.py", line 293, in process
    self._process(roi_indices, next_indices=next_indices)
  File "/usr/lib/python3/dist-packages/mdt/lib/processing_strategies.py", line 490, in _process
    results = self._model.get_post_optimization_output(x_final, results['status'])
  File "/usr/lib/python3/dist-packages/mdt/models/composite.py", line 257, in get_post_optimization_output
    volume_maps = self.post_process_optimization_maps(volume_maps, results_array=optimized_parameters)
  File "/usr/lib/python3/dist-packages/mdt/models/composite.py", line 348, in post_process_optimization_maps
    fim = self._compute_fisher_information_matrix(results_array)
  File "/usr/lib/python3/dist-packages/mdt/models/composite.py", line 589, in _compute_fisher_information_matrix
    cl_runtime_info=CLRuntimeInfo(double_precision=True)
  File "/usr/lib/python3/dist-packages/mot/cl_routines/numerical_differentiation.py", line 133, in estimate_hessian
    hessian_kernel.evaluate(kernel_data, nmr_voxels, use_local_reduction=True, cl_runtime_info=cl_runtime_info)
  File "/usr/lib/python3/dist-packages/mot/lib/cl_function.py", line 248, in evaluate
    use_local_reduction=use_local_reduction, cl_runtime_info=cl_runtime_info)
  File "/usr/lib/python3/dist-packages/mot/lib/cl_function.py", line 608, in apply_cl_function
    cl_function, kernel_data, cl_runtime_info.double_precision, use_local_reduction)
  File "/usr/lib/python3/dist-packages/mot/lib/cl_function.py", line 653, in __init__
    self._kernel = self._build_kernel(self._get_kernel_source(), compile_flags)
  File "/usr/lib/python3/dist-packages/mot/lib/cl_function.py", line 712, in _build_kernel
    return cl.Program(self._cl_context, kernel_source).build(' '.join(compile_flags))
  File "/usr/lib/python3/dist-packages/pyopencl/__init__.py", line 213, in build
    options=options, source=self._source)
  File "/usr/lib/python3/dist-packages/pyopencl/__init__.py", line 253, in _build_and_catch_errors
    raise err
pyopencl.RuntimeError: clBuildProgram failed: build program failure -

Build on <pyopencl.Device 'Intel(R) Xeon(R) CPU E5-2683 v4 @ 2.10GHz' on 'Intel(R) OpenCL' at 0x2751f48>:

Compilation started
6:529:35: error: scalar operand type has greater rank than the type of the vector element. ('mot_float_type' (aka 'double') and 'float4' (vector of 4 'float' values))
6:529:35: error: can't convert between vector values of different size ('float4' (vector of 4 'float' values) and 'mot_float_type' (aka 'double'))
Compilation failed

(options: -cl-denorms-are-zero -cl-mad-enable -cl-no-signed-zeros -I /usr/lib/python3/dist-packages/pyopencl/cl)
(source saved as /tmp/tmp_koj1nk4.cl)

[robbert-harms]

Hi Ali,

Thanks for reporting this, this was a bug in the handling of the gradient deviations. When I programmed this it was valid code. Maybe the OpenCL specifications changed or the compiler became more strict in handling the vector types.

Anyway, I fixed the issue and released a new version fixing this issue. If you would use v0.21.0 you should be able to run your data with the gradient deviations.

Best wishes,

Robbert

[akhanf]

Thanks Robbert, works great now!