tessl/cmake-simpleitk
SimpleITK is a simplified interface to the Insight Toolkit (ITK) for image registration and segmentation
—
Pending
Quality
Pending
Does it follow best practices?
Impact
Pending
No eval scenarios have been run
registration.mddocs/
Registration
Complete image registration framework for aligning images through optimization of similarity metrics. SimpleITK provides a comprehensive registration system that can handle rigid, affine, and deformable transformations with various similarity metrics and optimization methods.
Capabilities
ImageRegistrationMethod
The main registration framework that coordinates all registration components.
class ImageRegistrationMethod:
def __init__(self):
"""Initialize image registration method"""
def SetFixedImage(self, image: Image):
"""Set the fixed (reference) image"""
def GetFixedImage(self) -> Image:
"""Get the fixed image"""
def SetMovingImage(self, image: Image):
"""Set the moving (image to be registered) image"""
def GetMovingImage(self) -> Image:
"""Get the moving image"""
def SetInitialTransform(self, transform: Transform):
"""
Set initial transformation.
Args:
transform: Initial transform (will be optimized)
"""
def GetInitialTransform(self) -> Transform:
"""Get initial transform"""
def SetMovingInitialTransform(self, transform: Transform):
"""Set initial transform for moving image"""
def GetMovingInitialTransform(self) -> Transform:
"""Get moving initial transform"""
def SetFixedInitialTransform(self, transform: Transform):
"""Set initial transform for fixed image"""
def GetFixedInitialTransform(self) -> Transform:
"""Get fixed initial transform"""
def Execute(self) -> Transform:
"""
Execute registration and return optimized transform.
Returns:
Optimized transformation
"""
# Metric configuration
def SetMetricAsMeanSquares(self):
"""Use mean squares similarity metric"""
def SetMetricAsCorrelation(self):
"""Use normalized correlation metric"""
def SetMetricAsANTSNeighborhoodCorrelation(self, radius: int = 5):
"""
Use ANTS neighborhood correlation metric.
Args:
radius: Neighborhood radius
"""
def SetMetricAsJointHistogramMutualInformation(self, numberOfHistogramBins: int = 20,
varianceForJointPDFSmoothing: float = 1.5):
"""
Use joint histogram mutual information metric.
Args:
numberOfHistogramBins: Histogram bins
varianceForJointPDFSmoothing: Smoothing variance
"""
def SetMetricAsMattesMutualInformation(self, numberOfHistogramBins: int = 50):
"""
Use Mattes mutual information metric.
Args:
numberOfHistogramBins: Histogram bins for estimation
"""
def SetMetricFixedMask(self, mask: Image):
"""Set mask for fixed image metric computation"""
def SetMetricMovingMask(self, mask: Image):
"""Set mask for moving image metric computation"""
def SetMetricSamplingStrategy(self, strategy: str):
"""
Set sampling strategy for metric computation.
Args:
strategy: 'NONE', 'REGULAR', or 'RANDOM'
"""
def SetMetricSamplingPercentage(self, percentage: float):
"""
Set percentage of pixels to sample for metric computation.
Args:
percentage: Sampling percentage (0.0 to 1.0)
"""
# Optimizer configuration
def SetOptimizerAsRegularStepGradientDescent(self, learningRate: float = 1.0,
minStep: float = 0.001,
numberOfIterations: int = 100,
gradientMagnitudeTolerance: float = 1e-6):
"""
Use regular step gradient descent optimizer.
Args:
learningRate: Initial learning rate
minStep: Minimum step size for convergence
numberOfIterations: Maximum iterations
gradientMagnitudeTolerance: Gradient tolerance for convergence
"""
def SetOptimizerAsGradientDescent(self, learningRate: float = 1.0,
numberOfIterations: int = 100,
convergenceMinimumValue: float = 1e-6,
convergenceWindowSize: int = 10):
"""Use gradient descent optimizer"""
def SetOptimizerAsGradientDescentLineSearch(self, learningRate: float = 1.0,
numberOfIterations: int = 100,
convergenceMinimumValue: float = 1e-6,
convergenceWindowSize: int = 10):
"""Use gradient descent with line search optimizer"""
def SetOptimizerAsConjugateGradientLineSearch(self, learningRate: float = 1.0,
numberOfIterations: int = 100,
convergenceMinimumValue: float = 1e-6,
convergenceWindowSize: int = 10):
"""Use conjugate gradient line search optimizer"""
def SetOptimizerAsLBFGSB(self, gradientConvergenceTolerance: float = 1e-5,
numberOfIterations: int = 500,
maximumNumberOfCorrections: int = 5,
maximumNumberOfFunctionEvaluations: int = 2000,
costFunctionConvergenceFactor: float = 1e7):
"""Use L-BFGS-B optimizer for bounded optimization"""
def SetOptimizerAsOnePlusOneEvolutionary(self, numberOfIterations: int = 100,
epsilon: float = 1.5e-4,
initialRadius: float = 1.01,
growthFactor: float = -1.0,
shrinkFactor: float = -1.0):
"""Use (1+1) evolutionary strategy optimizer"""
def SetOptimizerAsPowell(self, maximumNumberOfIterations: int = 100,
maximumNumberOfLineIterations: int = 100,
stepLength: float = 1.0, stepTolerance: float = 1e-6,
valueTolerance: float = 1e-6):
"""Use Powell optimizer (derivative-free)"""
def SetOptimizerScales(self, scales: list):
"""
Set parameter scaling for optimizer.
Args:
scales: Scaling factors for each parameter
"""
def SetOptimizerWeights(self, weights: list):
"""Set parameter weights for optimizer"""
# Interpolator configuration
def SetInterpolator(self, interpolator: int):
"""
Set interpolation method.
Args:
interpolator: Interpolation type (e.g., sitk.sitkLinear, sitk.sitkBSpline)
"""
def GetInterpolator(self) -> int:
"""Get interpolation method"""
# Multi-resolution configuration
def SetShrinkFactorsPerLevel(self, factors: list):
"""
Set shrink factors for multi-resolution registration.
Args:
factors: List of shrink factors per resolution level
"""
def GetShrinkFactorsPerLevel(self) -> list:
"""Get shrink factors per level"""
def SetSmoothingSigmasPerLevel(self, sigmas: list):
"""
Set smoothing sigmas for multi-resolution registration.
Args:
sigmas: List of smoothing sigmas per level
"""
def GetSmoothingSigmasPerLevel(self) -> list:
"""Get smoothing sigmas per level"""
def SetSmoothingSigmasAreSpecifiedInPhysicalUnits(self, inPhysicalUnits: bool):
"""Specify if smoothing sigmas are in physical units"""
def GetSmoothingSigmasAreSpecifiedInPhysicalUnits(self) -> bool:
"""Check if smoothing sigmas are in physical units"""
# Results and monitoring
def GetOptimizerIteration(self) -> int:
"""Get current optimizer iteration"""
def GetOptimizerPosition(self) -> tuple:
"""Get current parameter values"""
def GetOptimizerLearningRate(self) -> float:
"""Get current learning rate"""
def GetOptimizerConvergenceValue(self) -> float:
"""Get optimizer convergence value"""
def GetMetricValue(self) -> float:
"""Get current metric value"""
def GetOptimizerStopConditionDescription(self) -> str:
"""Get description of why optimizer stopped"""
# Callbacks for monitoring
def AddCommand(self, event: int, command: object):
"""
Add command callback for registration events.
Args:
event: Event type (e.g., sitk.sitkIterationEvent)
command: Callback function or command object
"""
def RemoveAllCommands(self):
"""Remove all command callbacks"""Transform Initialization
Utilities for initializing transforms based on image geometry and landmarks.
def CenteredTransformInitializer(fixedImage: Image, movingImage: Image,
transform: Transform, operationMode: str = "GEOMETRY") -> Transform:
"""
Initialize transform based on image centers.
Args:
fixedImage: Fixed image
movingImage: Moving image
transform: Transform to initialize
operationMode: 'GEOMETRY' or 'MOMENTS'
Returns:
Initialized transform
"""
def CenteredVersorTransformInitializer(fixedImage: Image, movingImage: Image,
transform: VersorRigid3DTransform) -> Transform:
"""Initialize versor transform using image centers"""
def LandmarkBasedTransformInitializer(fixedLandmarks: list, movingLandmarks: list,
transform: Transform) -> Transform:
"""
Initialize transform using corresponding landmark points.
Args:
fixedLandmarks: List of points in fixed image
movingLandmarks: List of corresponding points in moving image
transform: Transform to initialize
Returns:
Initialized transform
"""
class BSplineTransformInitializer:
def __init__(self):
"""Initialize B-spline transform initializer"""
def SetTransformDomainMeshSize(self, meshSize: list):
"""Set mesh size for B-spline grid"""
def SetImage(self, image: Image):
"""Set reference image for transform domain"""
def Execute(self) -> BSplineTransform:
"""Create initialized B-spline transform"""Resampling
Apply transforms to images for registration results or coordinate transformation.
def Resample(image: Image, transform: Transform, interpolator: int = None,
defaultPixelValue: float = 0.0, outputPixelType: int = None,
outputSpacing: tuple = None, outputOrigin: tuple = None,
outputDirection: tuple = None, outputSize: tuple = None,
referenceImage: Image = None) -> Image:
"""
Resample image using transform.
Args:
image: Input image to resample
transform: Transformation to apply
interpolator: Interpolation method
defaultPixelValue: Value for pixels outside image
outputPixelType: Output pixel type
outputSpacing: Output image spacing
outputOrigin: Output image origin
outputDirection: Output image direction
outputSize: Output image size
referenceImage: Use geometry from reference image
Returns:
Resampled image
"""
class ResampleImageFilter:
def __init__(self):
"""Initialize resampling filter"""
def SetTransform(self, transform: Transform):
"""Set transformation"""
def GetTransform(self) -> Transform:
"""Get transformation"""
def SetInterpolator(self, interpolator: int):
"""Set interpolation method"""
def GetInterpolator(self) -> int:
"""Get interpolation method"""
def SetDefaultPixelValue(self, value: float):
"""Set default pixel value for outside regions"""
def GetDefaultPixelValue(self) -> float:
"""Get default pixel value"""
def SetSize(self, size: tuple):
"""Set output image size"""
def GetSize(self) -> tuple:
"""Get output image size"""
def SetOutputSpacing(self, spacing: tuple):
"""Set output image spacing"""
def GetOutputSpacing(self) -> tuple:
"""Get output image spacing"""
def SetOutputOrigin(self, origin: tuple):
"""Set output image origin"""
def GetOutputOrigin(self) -> tuple:
"""Get output image origin"""
def SetOutputDirection(self, direction: tuple):
"""Set output image direction"""
def GetOutputDirection(self) -> tuple:
"""Get output image direction"""
def SetReferenceImage(self, image: Image):
"""Set reference image for output geometry"""
def Execute(self, image: Image) -> Image:
"""Execute resampling"""Registration Evaluation
Metrics and utilities for evaluating registration quality.
def LabelOverlapMeasures(sourceImage: Image, targetImage: Image) -> dict:
"""
Compute overlap measures between labeled images.
Args:
sourceImage: Source labeled image
targetImage: Target labeled image
Returns:
Dictionary with overlap measures per label
"""
def HausdorffDistance(image1: Image, image2: Image) -> dict:
"""
Compute Hausdorff distance between binary images.
Args:
image1: First binary image
image2: Second binary image
Returns:
Dictionary with distance measures
"""
def SimilarityIndex(image1: Image, image2: Image) -> float:
"""
Compute similarity index (Dice coefficient) between binary images.
Args:
image1: First binary image
image2: Second binary image
Returns:
Similarity index (0 to 1)
"""
class LabelOverlapMeasuresImageFilter:
def __init__(self):
"""Initialize label overlap measures filter"""
def Execute(self, sourceImage: Image, targetImage: Image):
"""Compute overlap measures"""
def GetJaccardCoefficient(self, label: int = 1) -> float:
"""Get Jaccard coefficient for label"""
def GetDiceCoefficient(self, label: int = 1) -> float:
"""Get Dice coefficient for label"""
def GetVolumeSimilarity(self, label: int = 1) -> float:
"""Get volume similarity for label"""
def GetFalseNegativeError(self, label: int = 1) -> float:
"""Get false negative error for label"""
def GetFalsePositiveError(self, label: int = 1) -> float:
"""Get false positive error for label"""
def GetTargetOverlap(self, label: int = 1) -> float:
"""Get target overlap for label"""
def GetLabels(self) -> list:
"""Get list of labels processed"""Advanced Registration Components
Specialized registration components for complex scenarios.
# Demons registration (deformable)
def DemonsRegistration(fixedImage: Image, movingImage: Image,
numberOfIterations: int = 10, standardDeviations: float = 1.0) -> Image:
"""
Demons deformable registration.
Args:
fixedImage: Fixed image
movingImage: Moving image
numberOfIterations: Number of iterations
standardDeviations: Smoothing standard deviation
Returns:
Displacement field image
"""
def DiffeomorphicDemonsRegistration(fixedImage: Image, movingImage: Image,
numberOfIterations: int = 10,
standardDeviations: float = 1.0) -> Image:
"""Diffeomorphic demons registration for topology preservation"""
def FastSymmetricForcesDemonsRegistration(fixedImage: Image, movingImage: Image,
numberOfIterations: int = 10,
standardDeviations: float = 1.0) -> Image:
"""Fast symmetric forces demons registration"""
def SymmetricForcesDemonsRegistration(fixedImage: Image, movingImage: Image,
numberOfIterations: int = 10,
standardDeviations: float = 1.0) -> Image:
"""Symmetric forces demons registration"""
# Level set motion registration
def LevelSetMotionRegistration(fixedImage: Image, movingImage: Image,
numberOfIterations: int = 10,
standardDeviations: float = 1.0) -> Image:
"""Level set motion registration"""
# Displacement field utilities
def InverseDisplacementField(displacementField: Image,
maximumNumberOfIterations: int = 10,
maxErrorToleranceThreshold: float = 0.1,
meanErrorToleranceThreshold: float = 0.001,
enforceBoundaryCondition: bool = True) -> Image:
"""
Compute inverse of displacement field.
Args:
displacementField: Forward displacement field
maximumNumberOfIterations: Maximum iterations for inversion
maxErrorToleranceThreshold: Maximum error tolerance
meanErrorToleranceThreshold: Mean error tolerance
enforceBoundaryCondition: Enforce boundary conditions
Returns:
Inverse displacement field
"""
def IterativeInverseDisplacementField(displacementField: Image,
numberOfIterations: int = 5,
stopValue: float = 0.0) -> Image:
"""Iterative computation of inverse displacement field"""
def TransformToDisplacementField(transform: Transform, outputSize: tuple,
outputSpacing: tuple, outputOrigin: tuple,
outputDirection: tuple = None) -> Image:
"""
Convert transform to displacement field.
Args:
transform: Input transformation
outputSize: Size of output displacement field
outputSpacing: Spacing of output field
outputOrigin: Origin of output field
outputDirection: Direction of output field
Returns:
Displacement field image
"""
def DisplacementFieldJacobianDeterminant(displacementField: Image,
useImageSpacing: bool = True,
derivativeWeights: tuple = None) -> Image:
"""
Compute Jacobian determinant of displacement field.
Args:
displacementField: Input displacement field
useImageSpacing: Account for image spacing in derivatives
derivativeWeights: Weights for derivative computation
Returns:
Jacobian determinant image
"""Usage Examples
Basic Rigid Registration
import SimpleITK as sitk
# Load images
fixed_image = sitk.ReadImage('fixed.nii')
moving_image = sitk.ReadImage('moving.nii')
# Create registration method
registration = sitk.ImageRegistrationMethod()
# Set up similarity metric
registration.SetMetricAsMeanSquares()
# Set up optimizer
registration.SetOptimizerAsRegularStepGradientDescent(
learningRate=1.0,
minStep=0.001,
numberOfIterations=200
)
# Set interpolator
registration.SetInterpolator(sitk.sitkLinear)
# Create and initialize transform
initial_transform = sitk.CenteredTransformInitializer(
fixed_image, moving_image,
sitk.Euler3DTransform(),
sitk.CenteredTransformInitializerFilter.GEOMETRY
)
registration.SetInitialTransform(initial_transform)
# Execute registration
final_transform = registration.Execute()
print(f"Final metric value: {registration.GetMetricValue()}")
print(f"Optimizer stop condition: {registration.GetOptimizerStopConditionDescription()}")
# Apply transform to moving image
resampled = sitk.Resample(moving_image, final_transform, sitk.sitkLinear, 0.0,
moving_image.GetPixelID())
# Save result
sitk.WriteImage(resampled, 'registered_result.nii')Multi-modal Registration with Mutual Information
import SimpleITK as sitk
# Load multi-modal images (e.g., T1 and T2 MRI)
fixed_image = sitk.ReadImage('t1_image.nii')
moving_image = sitk.ReadImage('t2_image.nii')
# Registration setup
registration = sitk.ImageRegistrationMethod()
# Use mutual information for multi-modal registration
registration.SetMetricAsMattesMutualInformation(numberOfHistogramBins=50)
registration.SetMetricSamplingStrategy(registration.RANDOM)
registration.SetMetricSamplingPercentage(0.01) # Sample 1% of pixels
# Use more robust optimizer
registration.SetOptimizerAsGradientDescentLineSearch(
learningRate=1.0,
numberOfIterations=300,
convergenceMinimumValue=1e-6,
convergenceWindowSize=10
)
# Multi-resolution registration
registration.SetShrinkFactorsPerLevel([4, 2, 1])
registration.SetSmoothingSigmasPerLevel([2, 1, 0])
registration.SetSmoothingSigmasAreSpecifiedInPhysicalUnits(True)
# Initialize with similarity transform
initial_transform = sitk.CenteredTransformInitializer(
fixed_image, moving_image,
sitk.Similarity3DTransform(),
sitk.CenteredTransformInitializerFilter.GEOMETRY
)
registration.SetInitialTransform(initial_transform, inPlace=False)
# Execute registration
final_transform = registration.Execute()
print(f"Final metric value: {registration.GetMetricValue()}")
print(f"Final parameters: {final_transform.GetParameters()}")Registration with Monitoring
import SimpleITK as sitk
def command_iteration(method):
"""Callback function to monitor registration progress"""
print(f"Iteration {method.GetOptimizerIteration():3d} = "
f"{method.GetMetricValue():10.5f} : "
f"{method.GetOptimizerPosition()}")
# Set up registration
registration = sitk.ImageRegistrationMethod()
registration.SetMetricAsMeanSquares()
registration.SetOptimizerAsRegularStepGradientDescent(1.0, 0.001, 200)
registration.SetInterpolator(sitk.sitkLinear)
# Add monitoring callback
registration.AddCommand(sitk.sitkIterationEvent, lambda: command_iteration(registration))
# Set initial transform
initial_transform = sitk.Euler3DTransform()
registration.SetInitialTransform(initial_transform)
# Execute with monitoring
final_transform = registration.Execute(fixed_image, moving_image)Deformable Registration with B-splines
import SimpleITK as sitk
# Load images
fixed_image = sitk.ReadImage('fixed.nii')
moving_image = sitk.ReadImage('moving.nii')
# First perform rigid registration
rigid_registration = sitk.ImageRegistrationMethod()
rigid_registration.SetMetricAsMeanSquares()
rigid_registration.SetOptimizerAsRegularStepGradientDescent(1.0, 0.001, 100)
rigid_registration.SetInterpolator(sitk.sitkLinear)
rigid_transform = sitk.CenteredTransformInitializer(
fixed_image, moving_image, sitk.Euler3DTransform(),
sitk.CenteredTransformInitializerFilter.GEOMETRY
)
rigid_registration.SetInitialTransform(rigid_transform)
rigid_result = rigid_registration.Execute()
# Apply rigid transform to moving image
rigidly_registered = sitk.Resample(moving_image, rigid_result, sitk.sitkLinear, 0.0)
# Set up B-spline deformable registration
bspline_registration = sitk.ImageRegistrationMethod()
bspline_registration.SetMetricAsMeanSquares()
bspline_registration.SetOptimizerAsLBFGSB(
gradientConvergenceTolerance=1e-5,
numberOfIterations=100
)
bspline_registration.SetInterpolator(sitk.sitkLinear)
# Create B-spline transform
mesh_size = [8, 8, 8] # Control point grid size
bspline_transform = sitk.BSplineTransformInitializer(fixed_image, mesh_size)
bspline_registration.SetInitialTransform(bspline_transform, inPlace=False)
# Multi-resolution deformable registration
bspline_registration.SetShrinkFactorsPerLevel([4, 2, 1])
bspline_registration.SetSmoothingSigmasPerLevel([2, 1, 0])
# Execute deformable registration
final_bspline = bspline_registration.Execute(fixed_image, rigidly_registered)
# Apply final transform
final_result = sitk.Resample(moving_image,
sitk.CompositeTransform([rigid_result, final_bspline]),
sitk.sitkLinear, 0.0)
sitk.WriteImage(final_result, 'deformable_registration_result.nii')Registration Evaluation
import SimpleITK as sitk
# Load ground truth and registered segmentations
ground_truth = sitk.ReadImage('ground_truth_segmentation.nii')
registered_seg = sitk.ReadImage('registered_segmentation.nii')
# Compute overlap measures
overlap_filter = sitk.LabelOverlapMeasuresImageFilter()
overlap_filter.Execute(ground_truth, registered_seg)
# Get overlap statistics for each label
labels = overlap_filter.GetLabels()
print("Label\tDice\tJaccard\tVolume Similarity")
for label in labels:
dice = overlap_filter.GetDiceCoefficient(label)
jaccard = overlap_filter.GetJaccardCoefficient(label)
vol_sim = overlap_filter.GetVolumeSimilarity(label)
print(f"{label}\t{dice:.3f}\t{jaccard:.3f}\t{vol_sim:.3f}")
# Compute Hausdorff distance for binary images
if len(labels) == 2: # Binary case
hausdorff_filter = sitk.HausdorffDistanceImageFilter()
hausdorff_filter.Execute(ground_truth, registered_seg)
print(f"Hausdorff distance: {hausdorff_filter.GetHausdorffDistance():.2f}")
print(f"Average Hausdorff distance: {hausdorff_filter.GetAverageHausdorffDistance():.2f}")Displacement Field Analysis
import SimpleITK as sitk
# Load displacement field from registration
displacement_field = sitk.ReadImage('displacement_field.nii')
# Compute Jacobian determinant
jacobian = sitk.DisplacementFieldJacobianDeterminant(displacement_field)
# Analyze deformation
stats = sitk.StatisticsImageFilter()
stats.Execute(jacobian)
print(f"Jacobian statistics:")
print(f" Mean: {stats.GetMean():.3f}")
print(f" Std: {stats.GetSigma():.3f}")
print(f" Min: {stats.GetMinimum():.3f}")
print(f" Max: {stats.GetMaximum():.3f}")
# Check for folding (negative Jacobian)
negative_jacobian = sitk.BinaryThreshold(jacobian, -1e10, 0, 1, 0)
negative_count = sitk.StatisticsImageFilter()
negative_count.Execute(negative_jacobian)
print(f"Pixels with folding: {int(negative_count.GetSum())}")
# Save Jacobian for visualization
sitk.WriteImage(jacobian, 'jacobian_determinant.nii')Install with Tessl CLI
npx tessl i tessl/cmake-simpleitk