RDP gradient can generate NaNs if epsilon is zero

[KrisThielemans]

Reconstructing the mMR NEMA data gives me NaNs in the background of the image, then quickly ruining everything of course. Obviously, when epsilon=0, the function is not differentiable at x,y=0. We try to set the gradient to zero, but apparently still miss some cases. Relevant code is

STIR/src/recon_buildblock/RelativeDifferencePrior.cxx

Lines 403 to 424 in 8ced2d7

	elemT current;
	if (this->epsilon == 0.0 && current_image_estimate[z][y][x] == 0.0
	&& current_image_estimate[z + dz][y + dy][x + dx] == 0.0)
	{
	// handle the undefined nature of the gradient
	current = 0.0;
	}
	else
	{
	current
	= weights[dz][dy][dx]
	* (((current_image_estimate[z][y][x] - current_image_estimate[z + dz][y + dy][x + dx])
	* (this->gamma
	* abs(current_image_estimate[z][y][x] - current_image_estimate[z + dz][y + dy][x + dx])
	+ current_image_estimate[z][y][x] + 3 * current_image_estimate[z + dz][y + dy][x + dx]
	+ 2 * this->epsilon))
	/ (square((current_image_estimate[z][y][x] + current_image_estimate[z + dz][y + dy][x + dx])
	+ this->gamma
	* abs(current_image_estimate[z][y][x]
	- current_image_estimate[z + dz][y + dy][x + dx])
	+ this->epsilon)));
	}

@Imraj-Singh have you seen this?

[robbietuk]

Could be because of small x/y values?

[KrisThielemans]

This was fixed by #1410