Reconstruction fails when supplied initial image is zero

[smajee]

A) When performing a reconstruction with a supplied initial image (using -t flag) that is all zeros

1. Reconstruction fails (output image remains zero)
2. It does not produce any error messages
3. The printed non-homogeneous iterations is wrong (Equivalent iterations = 0.0, (non-homogeneous iterations = 300))

B) If the InitImageValue parameter in reconparams is set to zero for the initial condition, it shows a warning: 'Warning in ../params/sv-mbirct/object.reconparams: InitImageValue should be positive. Reverting to default.'

I understand that using a small positive scalar as an initial value is beneficial for speeding up convergence in most cases. However, in my opinion, the user should be free to choose any initial value/image that is suitable for the application.

Here is an excerpt of the command line output for case (A):

SUPER-VOXEL MBIR RECONSTRUCTION FOR 3D PARALLEL-BEAM CT build time: Jul 4 2020, 15:14:06

Parsing command line... -> will perform reconstruction (QGGMRF) -> will read system matrix from file -> will read initial condition from file(s) -> will compute forward project of initial condition NOTE you may save run time by pre-computing the initial projection See help (-e option) -> will save projection of output image state to file(s)

Filenames provided: Sino params = /depot/bouman/users/smajee/codedExp_project/workDir/demo/sv-mbirct_params/object.sinoparams Image params = /depot/bouman/users/smajee/codedExp_project/workDir/demo/sv-mbirct_params/object.imgparams Recon params = /depot/bouman/users/smajee/codedExp_project/workDir/demo/sv-mbirct_params/object.reconparams Sinogram data = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/sino/object_sliceNNN.2Dsinodata Weight data = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/weight/object_sliceNNN.2Dweightdata Output images = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/recon/object_sliceNNN.2Dimgdata Initial image = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/init/object_sliceNNN.2Dimgdata System matrix = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/sysmatrix/object.2Dsvmatrix Output projection = /depot/bouman/users/smajee/codedExp_project/workDir/data/sv-mbirct/proj/object_sliceNNN.2Dprojection SINOGRAM PARAMETERS:

• Number of sinogram views per slice = 2053
• Number of detector channels per slice = 91
• Number of slices = 1
• First slice index (wrt file names) = 0
• Spacing between detector channels = 1.0000000 (mm)
• Center of rotation offset = 0.0000000 (channels)
• Spacing between slices = 1.0000000 (mm) IMAGE PARAMETERS:
• Number of pixels per slice in (X,Y)-directions = (91,91)
• Number of slices (to reconstruct if output param) = 1
• First slice index (wrt sino/img file names) = 0
• Pixel width in XY plane = 1.0000000 (mm)
• Spacing between slices = 1.0000000 (mm)
• ROIRadius = 45.5000000 (mm) RECONSTRUCTION/PRIOR PARAMETERS:
• Prior Type = QGGMRF
• Q-GGMRF Prior Parameter, q = 1.200000
• Q-GGMRF Prior Parameter, p = 2.000000
• Q-GGMRF Prior Parameter, T = 1.000000
• Prior Regularization parameter, SigmaX = 0.0030000 (mm^-1)
• Scaling for sino weights, SigmaY (W=exp(-y)/SigmaY^2) = 1.0000000 (mm^-1)
• Prior weight for nearest neighbors within slice = 0.1132759
• Prior weight for diagonal neighbors within slice = 0.0800861
• Prior weight for nearest neighbors in adjacent slices = 0.1132759
• Inital image value = 0.020253 (mm-1)
• Stop threshold for convergence = 0.0000000 %
• Maximum number of ICD iterations = 3
• Positivity constraint flag = 1

Reading system matrix... Reading initial image... Projecting initial image... Projection time = 358 ms Reconstructing... No stopping condition--running fixed iterations Equivalent iterations = 0.0, (non-homogeneous iterations = 300) Average update in last iteration (relative) = 0.000000 % Average update in last iteration (magnitude) = 0.000000 mm^-1 Reconstruction time = 6350 ms Writing image files... Writing projection of image state... Done.

[sjkisner]

Result of "zero-skipping". Pixels surrounded by zeros are skipped from updating, so it skips everything and with no image change it looks like it's converged. A couple possible solutions: 1) Disable zero-skipping for the first iteration 2) Remove zero-skipping altogether. In the SV code the updates are relatively fast, so the additional branches and reads/writes with zero-skip checking for every pixel might offset the benefit of skipping some of the updates.

[sjkisner]

Added flag to disable zero-skipping. Currently set to disable zero-skipping only for 1st iteration.