SimpleITK 配准与重抽样预告

等忙过这段时间会写sitk的配准与重抽样。先写一个配准时遇到的小错误与解决方案。
错误如下：

Initial Number of Parameters: 375
Traceback (most recent call last):
File “”, line 2, in
File “”, line 18, in only_regis
File “”, line 20, in registration
File “C:\Users\sparrow\Anaconda3\lib\site-packages\SimpleITK\SimpleITK.py”, line 10387, in Execute
return _SimpleITK.ImageRegistrationMethod_Execute(self, fixed, moving)
RuntimeError: Exception thrown in SimpleITK ImageRegistrationMethod_Execute: D:\a\1\sitk\Code\Registration\src\sitkImageRegistrationMethod.cxx:791:
sitk::ERROR: Filter does not support fixed image type: 16-bit signed integer

如果你也遇到了这个错误，可以在加载数据的数据的时候指定类型为sitk.sitkFloat32，像这样：

fix_path = './data/struc_mri.nii.gz'
mov_path = './data/func_mri.nii.gz'
fix_img = sitk.ReadImage(fix_path, outputPixelType=sitk.sitkFloat32)
mov_img = sitk.ReadImage(mov_path, outputPixelType=sitk.sitkFloat32)

附一个简单的配准代码：

def command_iteration(method, bspline_transform):
    if method.GetOptimizerIteration() == 0:
        # The BSpline is resized before the first optimizer
        # iteration is completed per level. Print the transform object
        # to show the adapted BSpline transform.
        print(bspline_transform)

    print(
        f"{method.GetOptimizerIteration():3} "
        + f"= {method.GetMetricValue():10.5f}"
    )


def command_multi_iteration(method):
    # The sitkMultiResolutionIterationEvent occurs before the
    # resolution of the transform. This event is used here to print
    # the status of the optimizer from the previous registration level.

    if method.GetCurrentLevel() > 0:
        print(
            "Optimizer stop condition: "
            + f"{method.GetOptimizerStopConditionDescription()}"
        )
        print(f" Iteration: {method.GetOptimizerIteration()}")
        print(f" Metric value: {method.GetMetricValue()}")

    print("--------- Resolution Changing ---------")



def registration(refimg, movimg, out_path):
    fixed = refimg

    moving = movimg

    transformDomainMeshSize = [2] * fixed.GetDimension()
    tx = sitk.BSplineTransformInitializer(fixed, transformDomainMeshSize)

    print(f"Initial Number of Parameters: {tx.GetNumberOfParameters()}")

    R = sitk.ImageRegistrationMethod()
    R.SetMetricAsJointHistogramMutualInformation()

    R.SetOptimizerAsGradientDescentLineSearch(
        5.0, 100, convergenceMinimumValue=1e-4, convergenceWindowSize=5
    )

    R.SetInterpolator(sitk.sitkLinear)

    R.SetInitialTransformAsBSpline(tx, inPlace=True, scaleFactors=[1, 2, 5])
    R.SetShrinkFactorsPerLevel([4, 2, 1])
    R.SetSmoothingSigmasPerLevel([4, 2, 1])

    R.AddCommand(sitk.sitkIterationEvent, lambda: command_iteration(R, tx))
    R.AddCommand(
        sitk.sitkMultiResolutionIterationEvent, lambda: command_multi_iteration(R)
    )

    outTx = R.Execute(fixed, moving)

    print("-------")
    print(tx)
    print(outTx)
    print(f"Optimizer stop condition: {R.GetOptimizerStopConditionDescription()}")
    print(f" Iteration: {R.GetOptimizerIteration()}")
    print(f" Metric value: {R.GetMetricValue()}")

    # sitk.WriteTransform(outTx, out_path)
    if "SITK_NOSHOW" not in os.environ:
        resampler = sitk.ResampleImageFilter()
        resampler.SetReferenceImage(fixed)
        resampler.SetInterpolator(sitk.sitkLinear)
        resampler.SetDefaultPixelValue(100)
        resampler.SetTransform(outTx)

        out = resampler.Execute(moving)
        sitk.WriteImage(out, out_path)

上面的代码大部分来自官网的例子，本文进行了简单的修改。

out_path = './data/after_reg'
registration(fix_img, mov_img, out_path)

官方关于配准的链接