SimpleITK 图像配准
在网上找的资源,效果不佳,等清楚了函数和原理再细改,调试效果。
1 # -*- coding : UTF-8 -*- 2 # @file : regist.py 3 # @Time : 2021-11-12 17:00 4 # @Author : wmz 5 6 import SimpleITK as sitk 7 8 # Utility method that either downloads data from the MIDAS repository or 9 # if already downloaded returns the file name for reading from disk (cached data). 10 # %run update_path_to_download_script 11 # from downloaddata import fetch_data as fdata 12 13 # Always write output to a separate directory, we don't want to pollute the source directory. 14 import os 15 OUTPUT_DIR = 'Output' 16 17 import matplotlib.pyplot as plt 18 # % matplotlib 19 # inline 20 21 from ipywidgets import interact, fixed 22 from IPython.display import clear_output 23 24 25 # Callback invoked by the interact IPython method for scrolling through the image stacks of 26 # the two images (moving and fixed). 27 def display_images(fixed_image_z, moving_image_z, fixed_npa, moving_npa): 28 # Create a figure with two subplots and the specified size. 29 plt.subplots(1, 2, figsize=(10, 8)) 30 31 # Draw the fixed image in the first subplot. 32 plt.subplot(1, 2, 1) 33 plt.imshow(fixed_npa[fixed_image_z, :, :], cmap=plt.cm.Greys_r); 34 plt.title('fixed image') 35 plt.axis('off') 36 37 # Draw the moving image in the second subplot. 38 plt.subplot(1, 2, 2) 39 plt.imshow(moving_npa[moving_image_z, :, :], cmap=plt.cm.Greys_r); 40 plt.title('moving image') 41 plt.axis('off') 42 43 plt.show() 44 45 46 # Callback invoked by the IPython interact method for scrolling and modifying the alpha blending 47 # of an image stack of two images that occupy the same physical space. 48 def display_images_with_alpha(image_z, alpha, fixed, moving): 49 img = (1.0 - alpha) * fixed[:, :, image_z] + alpha * moving[:, :, image_z] 50 plt.imshow(sitk.GetArrayViewFromImage(img), cmap=plt.cm.Greys_r); 51 plt.axis('off') 52 plt.show() 53 54 55 # Callback invoked when the StartEvent happens, sets up our new data. 56 def start_plot(): 57 global metric_values, multires_iterations 58 59 metric_values = [] 60 multires_iterations = [] 61 62 63 # Callback invoked when the EndEvent happens, do cleanup of data and figure. 64 def end_plot(): 65 global metric_values, multires_iterations 66 67 del metric_values 68 del multires_iterations 69 # Close figure, we don't want to get a duplicate of the plot latter on. 70 plt.close() 71 72 73 # Callback invoked when the IterationEvent happens, update our data and display new figure. 74 def plot_values(registration_method): 75 global metric_values, multires_iterations 76 77 metric_values.append(registration_method.GetMetricValue()) 78 # Clear the output area (wait=True, to reduce flickering), and plot current data 79 clear_output(wait=True) 80 # Plot the similarity metric values 81 plt.plot(metric_values, 'r') 82 plt.plot(multires_iterations, [metric_values[index] for index in multires_iterations], 'b*') 83 plt.xlabel('Iteration Number', fontsize=12) 84 plt.ylabel('Metric Value', fontsize=12) 85 plt.show() 86 87 88 # Callback invoked when the sitkMultiResolutionIterationEvent happens, update the index into the 89 # metric_values list. 90 def update_multires_iterations(): 91 global metric_values, multires_iterations 92 multires_iterations.append(len(metric_values)) 93 94 95 fixed_image = sitk.ReadImage(r"E:\Data\right_knee\01-YCQ-cl-L.nrrd", sitk.sitkFloat32) 96 moving_image = sitk.ReadImage(r"E:\Data\right_knee\02-LAXI-r-cl lv.nrrd", sitk.sitkFloat32) 97 98 interact(display_images, fixed_image_z=(0,fixed_image.GetSize()[2]-1), moving_image_z=(0,moving_image.GetSize()[2]-1), fixed_npa = fixed(sitk.GetArrayViewFromImage(fixed_image)), moving_npa=fixed(sitk.GetArrayViewFromImage(moving_image))); 99 100 101 initial_transform = sitk.CenteredTransformInitializer(fixed_image, 102 moving_image, 103 sitk.Euler3DTransform(), 104 sitk.CenteredTransformInitializerFilter.GEOMETRY) 105 106 moving_resampled = sitk.Resample(moving_image, fixed_image, initial_transform, sitk.sitkLinear, 0.0, moving_image.GetPixelID()) 107 108 interact(display_images_with_alpha, image_z=(0,fixed_image.GetSize()[2]), alpha=(0.0,1.0,0.05), fixed = fixed(fixed_image), moving=fixed(moving_resampled)); 109 110 registration_method = sitk.ImageRegistrationMethod() 111 112 # Similarity metric settings. 113 registration_method.SetMetricAsMattesMutualInformation(numberOfHistogramBins=50) 114 registration_method.SetMetricSamplingStrategy(registration_method.RANDOM) 115 registration_method.SetMetricSamplingPercentage(0.01) 116 117 registration_method.SetInterpolator(sitk.sitkLinear) 118 119 # Optimizer settings. 120 registration_method.SetOptimizerAsGradientDescent(learningRate=1.0, numberOfIterations=100, convergenceMinimumValue=1e-6, convergenceWindowSize=10) 121 registration_method.SetOptimizerScalesFromPhysicalShift() 122 123 # Setup for the multi-resolution framework. 124 registration_method.SetShrinkFactorsPerLevel(shrinkFactors = [4,2,1]) 125 registration_method.SetSmoothingSigmasPerLevel(smoothingSigmas=[2,1,0]) 126 registration_method.SmoothingSigmasAreSpecifiedInPhysicalUnitsOn() 127 128 # Don't optimize in-place, we would possibly like to run this cell multiple times. 129 registration_method.SetInitialTransform(initial_transform, inPlace=False) 130 131 # Connect all of the observers so that we can perform plotting during registration. 132 registration_method.AddCommand(sitk.sitkStartEvent, start_plot) 133 registration_method.AddCommand(sitk.sitkEndEvent, end_plot) 134 registration_method.AddCommand(sitk.sitkMultiResolutionIterationEvent, update_multires_iterations) 135 registration_method.AddCommand(sitk.sitkIterationEvent, lambda: plot_values(registration_method)) 136 137 final_transform = registration_method.Execute(sitk.Cast(fixed_image, sitk.sitkFloat32), 138 sitk.Cast(moving_image, sitk.sitkFloat32)) 139 140 print('Final metric value: {0}'.format(registration_method.GetMetricValue())) 141 print('Optimizer\'s stopping condition, {0}'.format(registration_method.GetOptimizerStopConditionDescription())) 142 143 moving_resampled = sitk.Resample(moving_image, fixed_image, final_transform, sitk.sitkLinear, 0.0, moving_image.GetPixelID()) 144 145 interact(display_images_with_alpha, image_z=(0,fixed_image.GetSize()[2]), alpha=(0.0,1.0,0.05), fixed = fixed(fixed_image), moving=fixed(moving_resampled)); 146 147 sitk.WriteImage(moving_resampled, os.path.join(OUTPUT_DIR, 'RIRE_training_001_mr_T1_resampled.mha')) 148 sitk.WriteTransform(final_transform, os.path.join(OUTPUT_DIR, 'RIRE_training_001_CT_2_mr_T1.tfm'))
标签:配准,image,sitk,moving,SimpleITK,plt,图像,fixed,method From: https://www.cnblogs.com/ybqjymy/p/17550348.html