Registering two volumes along with mask is giving expection for sampler type "RandomSparseMask"

[yackmurali]

We have tried to perform the elastix registration for two 3D volumes and used the mask along with them. It is giving issue in C++ API but working file when passing through exe with sampler "RandomSparseMask" set it parameter map.

The following is the output we get:

Running elastix with parameter map 0
Current time: Fri Oct 18 12:38:29 2024.

ELASTIX version: 5.2.0

Command line options from ElastixBase:

-out D:\Python_projects\new_WMH_code\WMH2.0\275_100002450\20170417\cpp_results\Elastix_results/

-priority unspecified, so NORMAL process priority

-threads unspecified, so all available threads are used

Command line options from TransformBase:

-t0 unspecified, so no initial transform used

Reading images...

Reading images took 0 ms.

WARNING: the fixed pyramid schedule is not fully specified!

A default pyramid schedule is used.

WARNING: the moving pyramid schedule is not fully specified!

A default pyramid schedule is used.

WARNING: The parameter "CenterOfGravityUsesLowerThreshold", requested at entry number 0, does not exist at all.

The default value "false" is used instead.

WARNING: The parameter "NumberOfSamplesForCenteredTransformInitialization", requested at entry number 0, does not exist at all.

The default value "10000" is used instead.

Transform parameters are initialized as: [1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0]

InitializeTransform took 0.00s

Scales are estimated automatically.

Scales for transform parameters are: [2693.969331634363, 4794.691475386296, 6352.22389482894, 2693.969331634363, 4794.691475386296, 6352.22389482894, 2693.969331634363, 4794.691475386296, 6352.22389482894, 1, 1, 1]

Initialization of all components (before registration) took: 9 ms.

Preparation of the image pyramids took: 13 ms.

Resolution: 0

WARNING: The parameter "CheckNumberOfSamples", requested at entry number 0, does not exist at all.

The default value "true" is used instead.

WARNING: The parameter "UseMultiThreadingForMetrics", requested at entry number 0, does not exist at all.

The default value "true" is used instead.

Setting the fixed masks took: 0 ms.

Setting the moving masks took: 0 ms.

WARNING: The parameter "NumberOfFixedHistogramBins", requested at entry number 0, does not exist at all.

The default value "32" is used instead.

WARNING: The parameter "NumberOfMovingHistogramBins", requested at entry number 0, does not exist at all.

The default value "32" is used instead.

WARNING: The parameter "FixedLimitRangeRatio", requested at entry number 0, does not exist at all.

The default value "0.01" is used instead.

WARNING: The parameter "MovingLimitRangeRatio", requested at entry number 0, does not exist at all.

The default value "0.01" is used instead.

WARNING: The parameter "UseFastAndLowMemoryVersion", requested at entry number 0, does not exist at all.

The default value "true" is used instead.

WARNING: The parameter "UseJacobianPreconditioning", requested at entry number 0, does not exist at all.

The default value "false" is used instead.

WARNING: The parameter "FiniteDifferenceDerivative", requested at entry number 0, does not exist at all.

The default value "false" is used instead.

WARNING: The parameter "SP_A", requested at entry number 0, does not exist at all.

The default value "20" is used instead.

WARNING: The parameter "SigmoidInitialTime", requested at entry number 0, does not exist at all.

The default value "0" is used instead.

WARNING: The parameter "MaxBandCovSize", requested at entry number 0, does not exist at all.

The default value "192" is used instead.

WARNING: The parameter "NumberOfBandStructureSamples", requested at entry number 0, does not exist at all.

The default value "10" is used instead.

WARNING: The parameter "UseConstantStep", requested at entry number 0, does not exist at all.

The default value "false" is used instead.

WARNING: The parameter "MaximumStepLengthRatio", requested at entry number 0, does not exist at all.

The default value "1" is used instead.

WARNING: The parameter "MaximumStepLength", requested at entry number 0, does not exist at all.

The default value "1" is used instead.

WARNING: The parameter "NumberOfGradientMeasurements", requested at entry number 0, does not exist at all.

The default value "0" is used instead.

WARNING: The parameter "NumberOfJacobianMeasurements", requested at entry number 0, does not exist at all.

The default value "1000" is used instead.

WARNING: The parameter "NumberOfSamplesForExactGradient", requested at entry number 0, does not exist at all.

The default value "100000" is used instead.

WARNING: The parameter "SigmoidScaleFactor", requested at entry number 0, does not exist at all.

The default value "0.1" is used instead.

Elastix initialization of all components (for this resolution) took: 50 ms.

Initialization of AdvancedMattesMutualInformation metric took: 195 ms.

Starting automatic parameter estimation for AdaptiveStochasticGradientDescent ...

WARNING: The parameter "ASGDParameterEstimationMethod", requested at entry number 0, does not exist at all.

The default value "Original" is used instead.

Computing JacobianTerms ...

Computing the Jacobian terms took 0.001074s

NumberOfGradientMeasurements to estimate sigma_i: 3

Sampling gradients ...

Time spent in resolution 0 (ITK initialization and iterating): 0.216

Stopping condition: Error in metric.

Settings of AdaptiveStochasticGradientDescent in resolution 0:

( SP_a 1.000000 )

( SP_A 20.000000 )

( SP_alpha 0.602000 )

( SigmoidMax 1.000000 )

( SigmoidMin -0.800000 )

( SigmoidScale 0.000000 )

itk::ExceptionObject (00000088A7AFD7D0)

Location: "ElastixTemplate - Run()"

File: D:\Cpp_projects\elastix-5.2.0\Common\ImageSamplers\itkImageRandomSamplerSparseMask.hxx

Line: 44

Description: ITK ERROR: RandomSamplerSparseMask(000001C61599AC90): ERROR: do not call this function when no mask is supplied. When using the ImageRandomSamplerSparseMask sampler, a mask is required. Otherwise you may consider using a sampler that does not require a mask, for example, ImageRandomSampler.

Error occurred during actual registration.

Errors occurred!

This is the way we create the elastix mask from ITK Image Pointer.
elastixFixedMask = elastix::ELASTIX::ImagePointer(fixedMask.GetPointer());

[N-Dekker]

@yackmurali Thank you for your question! Can you possibly share a few more lines of your code? How do you pass the mask image to elastix?

Did you locally build elastix? If so, maybe you can make debug build, and put a break-point in the debugger at itkImageRandomSamplerSparseMask.hxx line 44, to see where it goes wrong!

[yackmurali]

Yes we have locally build elastix let me share the details of what we are doing and some observations we found after debugging so basically we are doing two step registration using the following api:
RegisterImages(ImagePointer fixedImage,
ImagePointer movingImage,
const std::vector & parameterMaps,
const std::string & outputPath,
bool performLogging,
bool performCout,
ImagePointer fixedMask = nullptr,
ImagePointer movingMask = nullptr,
ObjectPointer transform = nullptr);

Task and observations: In the first step we are doing the course registration and later using the results of the first we are futher refining the results in second registration. Using the above api we are passing both the parameterMaps in the parameterMap vector, everything works fine except in the second registration the parameterMap has RandomSparseMask as the imagesampler but first uses RandomCoordinate, as the RandomSparseMask needs masks we are passing masks in the api as arguments but after the first registration the masks are overwritten from the results of first registration to NULL (which may because we are not using masks for the first registration) and we get error in the second registration that we didn't supply masks.

One solution we found is to do both the registrations separately and pass the resultant parameterMap of first to the second as initial transform (ObjectPointer transform) but the resultant parameterMap after first registration (using ParameterMapType GetTransformParameterMap()) and the transform have different datatypes and we couldn't find a way to convert the ParameterMapType to ObjectPointer, are we missing something here or is there a solution for this issue?

When we tried passing the -t0 with first registration transform results (as Euler transform). using command prompt it is working.
But when we use the RegisterImages API we are not sure on how to create the intialTransform(of type objectpointer) from result of first registration result which is of type ParameterMapType .

We tried converting the ParameterMapType to InitialTransform as follows:
using ImageType = itk::Image<float, 3>;
using ParameterValuesType = std::vectorstd::string;
using ParameterMapType = std::map<std::string, ParameterValuesType>;
using ObjectPointer = itk::Object::Pointer;

ObjectPointer CreateInitialTransform(const ParameterMapType& parameterMap) {
using TransformType = itk::EulerTransform<double, 3>;
auto transform = TransformType::New();

// Parse TransformParameters
if (parameterMap.count("TransformParameters")) {
	std::vector<double> transformParams;
	const ParameterValuesType& paramsVector = parameterMap.at("TransformParameters");

	for (const auto& paramStr : paramsVector) {
		try {
			transformParams.push_back(std::stod(paramStr));
		}
		catch (const std::invalid_argument& e) {
			std::cerr << "Invalid TransformParameters value: " << paramStr << std::endl;
			return nullptr;
		}
	}

	if (transformParams.size() >= 6) {
		// Set translation using itk::Vector
		TransformType::OutputVectorType translation;
		translation[0] = transformParams[0];
		translation[1] = transformParams[1];
		translation[2] = transformParams[2];
		transform->SetTranslation(translation);
	}

	// Rotation parameters are omitted as AffineTransform does not directly use Euler angles
}

// Parse CenterOfRotationPoint if available
if (parameterMap.count("CenterOfRotationPoint")) {
	const ParameterValuesType& centerVector = parameterMap.at("CenterOfRotationPoint");
	std::vector<double> centerOfRotation;

	for (const auto& centerStr : centerVector) {
		try {
			centerOfRotation.push_back(std::stod(centerStr));
		}
		catch (const std::invalid_argument& e) {
			std::cerr << "Invalid CenterOfRotationPoint value: " << centerStr << std::endl;
			return nullptr;
		}
	}

	if (centerOfRotation.size() == 3) {
		// Set center using itk::Point
		TransformType::InputPointType center;
		center[0] = centerOfRotation[0];
		center[1] = centerOfRotation[1];
		center[2] = centerOfRotation[2];
		transform->SetCenter(center);
	}
}

return transform;

}
We are getting now following error:
elastix is started at Tue Nov 5 11:57:43 2024.

---

Running elastix with parameter map 0
Current time: Tue Nov 5 12:01:21 2024.
ELASTIX version: 5.2.0
Command line options from ElastixBase:
-out D:\Python_projects\new\cpp_results\Elastix_results/
-priority unspecified, so NORMAL process priority
-threads unspecified, so all available threads are used
Command line options from TransformBase:
-t0 unspecified, so no initial transform used

Reading images...
Reading images took 15 ms.

WARNING: the fixed pyramid schedule is not fully specified!
A default pyramid schedule is used.
WARNING: the moving pyramid schedule is not fully specified!
A default pyramid schedule is used.
WARNING: The parameter "CenterOfGravityUsesLowerThreshold", requested at entry number 0, does not exist at all.
The default value "false" is used instead.
WARNING: The parameter "NumberOfSamplesForCenteredTransformInitialization", requested at entry number 0, does not exist at all.
The default value "10000" is used instead.
Transform parameters are initialized as: [1, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0]
InitializeTransform took 0.01s
Scales are estimated automatically.
Scales for transform parameters are: [2693.969331634363, 4794.691475386296, 6352.22389482894, 2693.969331634363, 4794.691475386296, 6352.22389482894, 2693.969331634363, 4794.691475386296, 6352.22389482894, 1, 1, 1]
Initialization of all components (before registration) took: 75 ms.
Preparation of the image pyramids took: 27 ms.

Resolution: 0
WARNING: The parameter "UseMultiThreadingForMetrics", requested at entry number 0, does not exist at all.
The default value "true" is used instead.
Setting the fixed masks took: 0 ms.
Setting the moving masks took: 0 ms.
WARNING: The parameter "NumberOfFixedHistogramBins", requested at entry number 0, does not exist at all.
The default value "32" is used instead.
WARNING: The parameter "NumberOfMovingHistogramBins", requested at entry number 0, does not exist at all.
The default value "32" is used instead.
WARNING: The parameter "FixedLimitRangeRatio", requested at entry number 0, does not exist at all.
The default value "0.01" is used instead.
WARNING: The parameter "MovingLimitRangeRatio", requested at entry number 0, does not exist at all.
The default value "0.01" is used instead.
WARNING: The parameter "UseFastAndLowMemoryVersion", requested at entry number 0, does not exist at all.
The default value "true" is used instead.
WARNING: The parameter "UseJacobianPreconditioning", requested at entry number 0, does not exist at all.
The default value "false" is used instead.
WARNING: The parameter "FiniteDifferenceDerivative", requested at entry number 0, does not exist at all.
The default value "false" is used instead.
WARNING: The parameter "SP_A", requested at entry number 0, does not exist at all.
The default value "20" is used instead.
WARNING: The parameter "SigmoidInitialTime", requested at entry number 0, does not exist at all.
The default value "0" is used instead.
WARNING: The parameter "MaxBandCovSize", requested at entry number 0, does not exist at all.
The default value "192" is used instead.
WARNING: The parameter "NumberOfBandStructureSamples", requested at entry number 0, does not exist at all.
The default value "10" is used instead.
WARNING: The parameter "UseConstantStep", requested at entry number 0, does not exist at all.
The default value "false" is used instead.
WARNING: The parameter "MaximumStepLengthRatio", requested at entry number 0, does not exist at all.
The default value "1" is used instead.
WARNING: The parameter "MaximumStepLength", requested at entry number 0, does not exist at all.
The default value "1" is used instead.
WARNING: The parameter "NumberOfGradientMeasurements", requested at entry number 0, does not exist at all.
The default value "0" is used instead.
WARNING: The parameter "NumberOfJacobianMeasurements", requested at entry number 0, does not exist at all.
The default value "1000" is used instead.
WARNING: The parameter "NumberOfSamplesForExactGradient", requested at entry number 0, does not exist at all.
The default value "100000" is used instead.
WARNING: The parameter "SigmoidScaleFactor", requested at entry number 0, does not exist at all.
The default value "0.1" is used instead.
Elastix initialization of all components (for this resolution) took: 30 ms.
Initialization of AdvancedMattesMutualInformation metric took: 2048 ms.
Starting automatic parameter estimation for AdaptiveStochasticGradientDescent ...
WARNING: The parameter "ASGDParameterEstimationMethod", requested at entry number 0, does not exist at all.
The default value "Original" is used instead.
Computing JacobianTerms ...
Computing the Jacobian terms took 0.013143s
NumberOfGradientMeasurements to estimate sigma_i: 3
Sampling gradients ...
Time spent in resolution 0 (ITK initialization and iterating): 3
Stopping condition: Error in metric.
Settings of AdaptiveStochasticGradientDescent in resolution 0:
( SP_a 1.000000 )
( SP_A 20.000000 )
( SP_alpha 0.602000 )
( SigmoidMax 1.000000 )
( SigmoidMin -0.800000 )
( SigmoidScale 0.000000 )

[yackmurali]

@yackmurali Thank you for your question! Can you possibly share a few more lines of your code? How do you pass the mask image to elastix?

Did you locally build elastix? If so, maybe you can make debug build, and put a break-point in the debugger at itkImageRandomSamplerSparseMask.hxx line 44, to see where it goes wrong!

We build elastix locally. While debugging we were able to find out the crash is due to type mismatch of Mask buffer type. We have given the type as Float which was supposed to be of type unsigned char type. There is no crash observed after this change.

Now the issue left is on InitialTransformation parameter.

In our code we are doing registration of 3D MR sequence in two stages.

Stage1 : We perform the registration with Fixed and Moving 3D images with Euler as transform.
For this we use the first registration API.

Stage 2: For the second registration we pass the same Fixed and Moving 3D images that were used in first registration along with FixedMask and Intial Tranformation (which Euler). For the second registration we chosse Affine Transformation.
We were not clear on how to create the InitialTranformation from the result of first registration (ParameterMapType). We don't know if the convertion of ParameterType (of euler transform) to InitialTransform which is of type ObjectPointer.

The following code is used for converting the ParameterMapType (First registration result of type Euler Transform) to Object Pointer.
"
ObjectPointer CreateInitialTransform(const ParameterMapType& parameterMap) {
using TransformType = itk::EulerTransform<double, 3>;
auto transform = TransformType::New();

// Parse TransformParameters
if (parameterMap.count("TransformParameters")) {
std::vector transformParams;
const ParameterValuesType& paramsVector = parameterMap.at("TransformParameters");

for (const auto& paramStr : paramsVector) {
	try {
		transformParams.push_back(std::stod(paramStr));
	}
	catch (const std::invalid_argument& e) {
		std::cerr << "Invalid TransformParameters value: " << paramStr << std::endl;
		return nullptr;
	}
}

if (transformParams.size() >= 6) {
	// Set translation using itk::Vector
	TransformType::OutputVectorType translation;
	translation[0] = transformParams[0];
	translation[1] = transformParams[1];
	translation[2] = transformParams[2];
	transform->SetTranslation(translation);
}

// Rotation parameters are omitted as AffineTransform does not directly use Euler angles

}

// Parse CenterOfRotationPoint if available
if (parameterMap.count("CenterOfRotationPoint")) {
const ParameterValuesType& centerVector = parameterMap.at("CenterOfRotationPoint");
std::vector centerOfRotation;

for (const auto& centerStr : centerVector) {
	try {
		centerOfRotation.push_back(std::stod(centerStr));
	}
	catch (const std::invalid_argument& e) {
		std::cerr << "Invalid CenterOfRotationPoint value: " << centerStr << std::endl;
		return nullptr;
	}
}

if (centerOfRotation.size() == 3) {
	// Set center using itk::Point
	TransformType::InputPointType center;
	center[0] = centerOfRotation[0];
	center[1] = centerOfRotation[1];
	center[2] = centerOfRotation[2];
	transform->SetCenter(center);
}

}

return transform;
}
"

We felt this transformation results after second registration using Intial Transformation and Mask is same as First registration.

The same is working fine if we pass the command line parameters for elastix.exe with Fixed ,Moving , Fixed mask and t0 InitialTransform(from first registration) as files.

Can you please provide a sample code to convert the ParameterMapType to InitialTranform.