Artifacts in reconstructed images

[Kevin-Mattheus-Moerman]

I added the following to your code to visualize errors in the reconstructed images.

%Erors in reconstructed images
figure

% reference image
subplot( 2, 4, 1 )
imagesc( IRef(:,:,midSliceIndex).' )
axis image
colormap( gray ); colorbar; 
title( 'reference image' )

for i=1:1:numel(G)    
    IRef_inverse = dvf.imdeform( IDfm, G{i} );    
    % error of undeformed (study) image
    subplot( 2, 4, i+1)
    imagesc( (IRef_inverse(:,:,midSliceIndex)-IRef(:,:,midSliceIndex)).' )
    axis image
    colorbar; caxis([-1 1]);
    title( ['Recovered undeformed, scheme: ',num2str(i)] );
    drawnow
end

For the 2D demo I then get: And for the 3D demo I get:

1. Can you explain these patterns?
2. Why does method 4 contain hole like artifacts for the 2D method?
3. Why does method 4-7 contain hole like artifacts for the 3D method? Why is performance so different from the 2D method?
4. Why does method 7 create an additional blue rim around the image? Are NaN's introduced around the border for this method perhaps?

[ailiop]

The reviewer correctly pointed out certain artifacts in the inverse DVF results. Most of the artifacts emerged in displacement space (i.e., the displacement and IC residual vector fields) during the inversion iterations. I have located and fixed the issues that led to these artifacts.

The reviewer also warped the demo study image with the inverse DVF estimates, per the code posted above, and found grid-like artifacts in image space (i.e., the scalar intensity fields). These artifacts affected all images and were the combined effect of several factors. Some were translated from the displacement-space artifacts. The others pertain to warping the reference image, which entails interpolation.

I have amended the code to resolve these issues and remove the artifacts that the reviewer pointed out. Below are described the specific sources of the artifacts and the solutions in the updated code.

Removal of erroneous constraints in the source code

There were two inexcusable mistakes/bugs in the source code which were responsible for the major artifacts that arose in displacement space (questions 2 and 3 by the reviewer).

One mistake was in the calculation of locally optimal feedback control values with 2D DVFs in function dvf.feedbackControlVal. Specifically, the case of complex eigenvalues in the deformation Jacobian was missing. Like in the 3D demos, the forward deformation in our 2D demo features many complex eigenvalues. I have corrected the calculation with complex eigenvalue pairs.

The other mistake was the setting of the variable windowControl in function dvf.inversion. It was implicitly assumed that the size of the spatial domain in every dimension would be sufficiently larger than the maximum forward displacement. The data in our 3D demo were too small along the z-axis, violating this assumption. This led to degenerate neighborhoods for local search in spatially variant control schemes and spurious errors in the 3D demo results. I have removed the unnecessary restriction.

Mitigation of boundary effects

The "rim" artifact is resolved by padding the DVF domain before initiating the inversion iteration. I have updated function dvf.inversion to do that automatically in the case of multi-scale iterations. This functionality may be disabled via the new flagPrepad parameter.

Inversion assessment and image-space visualization

Our means of assessing the inverse DVF estimates are the IC residual vector fields. The residuals are directly related to the DVF inversion error, which is unknown. In addition to the percentile curves of IC residual magnitudes throughout the iteration, I have added visualizations of IC residual magnitude maps over the spatial domain at the end of the iteration, to facilitate visual inspection in displacement space.

It is often desirable to visualize the effect of the inverse DVF in image space, although the inversion process does not involve the image intensity field. The reviewer suggested inspecting the difference between the reference image and its recovery via backwards-warping of the deformed study image. This approach gave rise to very prominent grid-like artifacts, even where the IC residuals were effectively eliminated. The reason was interpolation errors in the forward-warped study image; these were exacerbated by the coarse intensity quantization in our synthetic image. I have replaced the image with a smoother one, and reformulated the reference image recovery calculations to circumvent forward warping. The formulation is equivalent to forward-backward warping in free space, but greatly reduces interpolation errors that are unrelated to the quality of the inverse DVF.

[Kevin-Mattheus-Moerman]

Thanks for fixing this and for the detailed reply.