Прямое и обратное преобразование Радона в среде MATLAB

Функция iradon выполняет обратное преобразование Радона, которое широко применяется в томографии. Это преобразование является инверсным к преобразованию Радона и выполняет реконструкцию изображений по их проекционным данным.

Функция прямых преобразований Радона radon для своего выполнения требует установки двух основных параметров – исходного изображения I и углов theta.

R = radon(I,theta);

Функция iradon реализует восстановление изображения I на основе проекционных данных.

IR = iradon(R,theta);

Функция iradon восстанавливает полутоновое изображение IR из матрицы проекций R, которые получены под соответствующими углами theta. Значения углов в векторе theta определяют направление осей. В общей сложности, функция iradon позволяет задать четыре параметра. Кроме описанных двух параметров, функция iradon задает также тип интерполяции, которая используется при восстановлении исходного изображения. Еще один параметр дает возможность определить фильтр предварительной обработки проекций, из которых будет восстанавливаться изображение. Если эти параметры не указываются, то используется соответствующее значение по умолчанию.

С помощью функции phantom сымитировали снимок головы человека. Выполнили прямое преобразование Радона (radon), выбрав различные углы. Восстановили изображение с помощью обратного преобразования Радона. Оценили качество восстановленного изображения.

close all; clear all; clc;

P=phantom(256); %Создание искусственного изображения среза головы

figure(1)

imagesc(P);

title('Head Phantom');

xlabel('x\prime (pixels)');

ylabel('y\prime (pixels)');

Рис. 6. Сгенерированный снимок головы.

%Прямые преобразования Радона с различными углами theta

theta1 = 0:10:170;

[R1,xp1] = radon(P,theta1);

theta2 = 0:5:175;

[R2,xp2] = radon(P,theta2);

theta3 = 0:2:178;

[R3,xp3] = radon(P,theta3);

theta4 = 0:0.1:179;

[R4,xp4] = radon(P,theta4);

figure(2)

subplot(2,2,1)

imagesc(theta1,xp1,R1)

title('Num angles = 18');

xlabel('Rotation Angle - \theta (degrees)');

ylabel('Sensor Position - x\prime (pixels)');

colormap(hot)

subplot(2,2,2)

imagesc(theta2,xp2,R2)

title('Num angles = 36');

xlabel('Rotation Angle - \theta (degrees)');

ylabel('Sensor Position - x\prime (pixels)');

colormap(hot)

subplot(2,2,3)

imagesc(theta3,xp3,R3)

title('Num angles = 90');

xlabel('Rotation Angle - \theta (degrees)');

ylabel('Sensor Position - x\prime (pixels)');

colormap(hot)

subplot(2,2,4)

imagesc(theta4,xp4,R4)

title('Num angles = 1791');

xlabel('Rotation Angle - \theta (degrees)');

ylabel('Sensor Position - x\prime (pixels)');

colormap(hot)

Рис. 7. Преобразования Радона с различными углами.

%Вычисление обратного преобразования Радона с различными углами theta

IR1=iradon(R1,theta1,'nearest','Ram-Lak',1,256);

IR2=iradon(R2,theta2,'nearest','Ram-Lak',1,256);

IR3=iradon(R3,theta3,'nearest','Ram-Lak',1,256);

IR4=iradon(R4,theta4,'nearest','Ram-Lak',1,256);

figure(3)

subplot(2,2,1)

imagesc(theta1,xp1,IR1)

title('Num angles = 18');

xlabel('x\prime (pixels)');

ylabel('y\prime (pixels)');

subplot(2,2,2)

imagesc(theta2,xp2,IR2)

title('Num angles = 36');

xlabel('x\prime (pixels)');

ylabel('y\prime (pixels)');

subplot(2,2,3)

imagesc(theta3,xp3,IR3)

title('Num angles = 90');

xlabel('x\prime (pixels)');

ylabel('y\prime (pixels)');

subplot(2,2,4)

imagesc(theta4,xp4,IR4)

title('Num angles = 1791');

xlabel('x\prime (pixels)');

ylabel('y\prime (pixels)');

Рис. 8. Обратные преобразования Радона с различными углами.

Таб. 1. Оценка качества восстановленных изображений.