Paso 1: Cargar imágenes

Este ejemplo utiliza dos imágenes tridimensionales de la cabeza del mismo paciente. En los problemas de registro, consideramos que una imagen es la imagen fija y la otra imagen es la imagen en movimiento. El objetivo del registro es alinear la imagen en movimiento con la imagen fija. En este ejemplo, la imagen fija es una imagen RMN ponderada T1. La imagen en movimiento que queremos registrar es una imagen de TC. Los datos se almacenan en el formato de archivo utilizado por el proyecto retrospectivo de evaluación de registro de imágenes (RIRE). Se utiliza para leer los archivos binarios que contienen datos de imagen.multibandread Utilice la función para obtener los metadatos asociados a cada imagen.helperReadHeaderRIRE Puede utilizar el siguiente enlace para encontrar más información sobre el formato de archivo RIRE:Formato de datos RIRE

fixedHeader  = helperReadHeaderRIRE('rirePatient007MRT1.header'); movingHeader = helperReadHeaderRIRE('rirePatient007CT.header');  fixedVolume  = multibandread('rirePatient007MRT1.bin',...                             [fixedHeader.Rows, fixedHeader.Columns, fixedHeader.Slices],...                             'int16=>single', 0, 'bsq', 'ieee-be' );                          movingVolume = multibandread('rirePatient007CT.bin',...                             [movingHeader.Rows, movingHeader.Columns, movingHeader.Slices],...                             'int16=>single', 0, 'bsq', 'ieee-be' );

La función es una función auxiliar que se proporciona para ayudar a juzgar la calidad de los resultados de registro 3D.helperVolumeRegistration Puede girar interactivamente la vista y ambos ejes permanecerán sincronizados.

helperVolumeRegistration(fixedVolume,movingVolume);

También puede utilizar para mirar planos individuales de los volúmenes fijos y en movimiento para obtener una idea de la alineación general de los volúmenes.imshowpair En la imagen superpuesta de , las áreas grises corresponden a áreas que tienen intensidades similares, mientras que las áreas magenta y verde muestran lugares donde una imagen es más brillante que la otra.imshowpair Se utiliza para observar el registro incorrecto de los volúmenes de imagen a lo largo de un segmento axial tomado a través del centro de cada volumen.imshowpair

centerFixed = size(fixedVolume)/2; centerMoving = size(movingVolume)/2; figure imshowpair(movingVolume(:,:,centerMoving(3)), fixedVolume(:,:,centerFixed(3))); title('Unregistered Axial Slice')

Paso 2: Configurar el Registro Inicial

La función facilita la elección del optimizador y la configuración de métricas correctos para usar con .imregconfigimregister Estas dos imágenes proceden de dos modalidades diferentes, la RMN y la TC, por lo que la opción "multimodal" es adecuada.

[optimizer,metric] = imregconfig('multimodal');

El algoritmo utilizado por convergerá para obtener mejores resultados más rápidamente cuando se especifique información de referencia espacial sobre la resolución y/o ubicación de las imágenes de entrada.imregister En este caso, la resolución de los conjuntos de datos CT y MRI se define en los metadatos de la imagen. Utilice estos metadatos para construir objetos de referencia espacial que pasaremos como argumentos de entrada para el registro.imref3d

Rfixed  = imref3d(size(fixedVolume),fixedHeader.PixelSize(2),fixedHeader.PixelSize(1),fixedHeader.SliceThickness); Rmoving = imref3d(size(movingVolume),movingHeader.PixelSize(2),movingHeader.PixelSize(1),movingHeader.SliceThickness);

Las propiedades de los objetos de referencia espacial definen dónde se encuentran los volúmenes de imagen asociados en el sistema de coordenadas universales y cuál es la extensión de píxel en cada dimensión. La propiedad XWorldLimits de Rmoving define la posición del volumen móvil en la dimensión X. La propiedad PixelExtentInWorld define el tamaño de cada píxel en unidades del mundo en la dimensión X (a lo largo de las columnas). El volumen en movimiento se extiende de 0.3269 a 334.97 en el sistema de coordenadas X mundial y cada píxel tiene una extensión de 0.6536mm. Las unidades están en milímetros porque la información de encabezado utilizada para construir la referencia espacial estaba en milímetros. La propiedad ImageExtentInWorldX determina la extensión completa de las unidades mundiales del volumen de imagen en movimiento en unidades del mundo.

Rmoving.XWorldLimits

ans = 1×2

    0.3268  334.9674

Rmoving.PixelExtentInWorldX

ans = 0.6536 

Rmoving.ImageExtentInWorldX

ans = 334.6406 

La desalineación entre los dos volúmenes incluye la conversión, el escalado y la rotación. Utilice una transformación de similitud para registrar las imágenes.

Comience utilizando para obtener un volumen de imagen de salida registrado que puede ver y observar directamente para acceder a la calidad de los resultados de registro.imregister

Especifique una configuración no predeterminada para la propiedad InitialRadius del optimizador para lograr una mejor convergencia en los resultados de registro.

optimizer.InitialRadius = 0.004; movingRegisteredVolume = imregister(movingVolume,Rmoving, fixedVolume,Rfixed, 'rigid', optimizer, metric);

Utilice de nuevo y repita el proceso de examinar la alineación de un segmento axial tomado a través del centro de los volúmenes registrados para tener una idea de lo exitoso que es el registro.imshowpair

figure imshowpair(movingRegisteredVolume(:,:,centerFixed(3)), fixedVolume(:,:,centerFixed(3))); title('Axial Slice of Registered Volume')

Desde el sector axial anterior, parece que el registro se realizó correctamente. Utilícelo de nuevo para ver el volumen registrado para seguir juzgando el éxito del registro.helperVolumeRegistration

helperVolumeRegistration(fixedVolume,movingRegisteredVolume);

Paso 3: Obtener transformación geométrica 3D que alinea el movimiento con fijo.

La función se puede utilizar cuando está interesado en la estimación de transformación geométrica que se utiliza para formar la imagen de salida registrada. utiliza el mismo algoritmo que y toma los mismos argumentos de entrada que .imregtformimregisterimregtformimregisterimregister Puesto que la inspección visual del volumen resultante indica que el registro se realizó correctamente, puede llamar con los mismos argumentos de entrada para obtener la transformación geométrica asociada a este resultado de registro.imregisterimregtform

geomtform = imregtform(movingVolume,Rmoving, fixedVolume,Rfixed, 'rigid', optimizer, metric)

geomtform =    affine3d with properties:      Dimensionality: 3                  T: [4x4 double]  

El resultado de imregtform es un objeto de transformación geométrica. Este objeto incluye una propiedad, T, que define la matriz de transformación afín 3D.

geomtform.T

ans = 4×4

    0.9704   -0.0143   -0.2410         0
    0.0228    0.9992    0.0324         0
    0.2404   -0.0369    0.9700         0
  -15.8648  -17.5692   29.1830    1.0000

El método transformPointsForward de la transformación geométrica se puede utilizar para determinar dónde se puede determinar dónde se puede utilizar un punto [u,v,w] en los mapas de imágenes en movimiento como resultado del registro. Dado que las entradas referenciadas espacialmente se especificaron para imregtform, la transformación geométrica asigna puntos en el sistema de coordenadas universales de moverse a fijo. El método transformPointsForward se utiliza a continuación para determinar la ubicación transformada del centro de la imagen en movimiento en el sistema de coordenadas universales.

centerXWorld = mean(Rmoving.XWorldLimits); centerYWorld = mean(Rmoving.YWorldLimits); centerZWorld = mean(Rmoving.ZWorldLimits); [xWorld,yWorld,zWorld] = transformPointsForward(geomtform,centerXWorld,centerYWorld,centerZWorld);

Puede utilizar el método worldToSubscript de Rfixed para determinar el elemento del volumen fijo que se alinea con el centro del volumen en movimiento.

[r,c,p] = worldToSubscript(Rfixed,xWorld,yWorld,zWorld)

r = 116 

c = 132 

p = 13 

Paso 4: Aplique la estimación de transformación geométrica al volumen de la imagen en movimiento.

La función se puede utilizar para aplicar la estimación de transformación geométrica de imregtform a un volumen 3D.imwarp El nombre/valor 'OutputView' se utiliza para definir un argumento de referencia espacial que determina los límites mundiales y la resolución de la imagen remuestreada de salida. Puede producir los mismos resultados dados mediante el objeto de referencia espacial asociado a la imagen fija.imregister Esto crea un volumen de salida en el que los límites del mundo y la resolución de la imagen fija y en movimiento son los mismos. Una vez que los límites mundiales y la resolución de ambos volúmenes son los mismos, hay correspondencia de píxeles a píxeles entre cada muestra de los volúmenes en movimiento y fijos.

movingRegisteredVolume = imwarp(movingVolume,Rmoving,geomtform,'bicubic','OutputView',Rfixed);

Utilice imshowpair de nuevo para ver un sector axial a través del centro del volumen registrado producido por .imwarp

figure  imshowpair(movingRegisteredVolume(:,:,centerFixed(3)), fixedVolume(:,:,centerFixed(3))); title('Axial Slice of Registered Volume')