diff --git a/zUtil_Deformation/gtGetODFFromGvdm.m b/zUtil_Deformation/gtGetODFFromGvdm.m index 6b78a0585c522052176747be7b2f62ca8ba38d74..f64e89092306b16a583261d7a050ffe9e9d803f3 100644 --- a/zUtil_Deformation/gtGetODFFromGvdm.m +++ b/zUtil_Deformation/gtGetODFFromGvdm.m @@ -4,26 +4,42 @@ function odf = gtGetODFFromGvdm(gvdm, grid_points, gvint) % ddm is a <n_voxels x 3> vector, where each row is the Rodriguez vector of % each voxel - size_odf = size(grid_points); - odf = zeros(size_odf); - - if (size(gvdm, 1) == 3) - gvdm = gvdm'; - end if (~exist('gvint', 'var')) gvint = ones(size(gvdm, 1), 1); else gvint = reshape(gvint, [], 1); end - voxel_size = grid_points{2, 2, 2}.R_vector - grid_points{1, 1, 1}.R_vector; + if (size(gvdm, 1) == 3) + gvdm = gvdm'; + end - z_edges = linspace(grid_points{1, 1, 1}.R_vector(3)-voxel_size(3)/2, ... - grid_points{1, 1, end}.R_vector(3)+voxel_size(3)/2, size_odf(3)+1); - y_edges = linspace(grid_points{1, 1, 1}.R_vector(2)-voxel_size(2)/2, ... - grid_points{1, end, 1}.R_vector(2)+voxel_size(2)/2, size_odf(2)+1); - x_edges = linspace(grid_points{1, 1, 1}.R_vector(1)-voxel_size(1)/2, ... - grid_points{end, 1, 1}.R_vector(1)+voxel_size(1)/2, size_odf(1)+1); + if (iscell(grid_points)) + voxel_size = grid_points{2, 2, 2}.R_vector - grid_points{1, 1, 1}.R_vector; + + size_odf = size(grid_points); + + z_edges = linspace(grid_points{1, 1, 1}.R_vector(3)-voxel_size(3)/2, ... + grid_points{1, 1, end}.R_vector(3)+voxel_size(3)/2, size_odf(3)+1); + y_edges = linspace(grid_points{1, 1, 1}.R_vector(2)-voxel_size(2)/2, ... + grid_points{1, end, 1}.R_vector(2)+voxel_size(2)/2, size_odf(2)+1); + x_edges = linspace(grid_points{1, 1, 1}.R_vector(1)-voxel_size(1)/2, ... + grid_points{end, 1, 1}.R_vector(1)+voxel_size(1)/2, size_odf(1)+1); + else + % We expect a <1x9> vector with <1x6> o-space bb, plus <1x3> the + % resolution + voxel_size = grid_points(7:9) + + size_odf = (grid_points(4:6) - grid_points(1:3)) ./ voxel_size + + z_edges = linspace(grid_points(3) - voxel_size(3) / 2, ... + grid_points(6) + voxel_size(3) / 2, size_odf(3) + 1); + y_edges = linspace(grid_points(2) - voxel_size(2) / 2, ... + grid_points(5) + voxel_size(2) / 2, size_odf(2) + 1); + x_edges = linspace(grid_points(1) - voxel_size(1) / 2, ... + grid_points(4) + voxel_size(1) / 2, size_odf(1) + 1); + end + odf = zeros(size_odf); fprintf('Producing ODF..'); c = tic();