\documentclass{cmi}
\usepackage{upgreek}
\usepackage{multirow}
\usepackage{graphicx}

% Данные поля заполняются редакцией журнала
%\makeatletter
%\newcommand{\@volume}{}
%\newcommand{\@issue}{}
%\newcommand{\@num}{}
%\newcommand{\@pages}{}

\begin{document}

\classify{УДК 519.6}

\renewcommand*{\thefootnote}{\fnsymbol{footnote}} %to make a number in end of title delete this line

\title {РЕКОНСТРУКЦИЯ ИЗОБРАЖЕНИЯ ПО МЕТОДУ ОБРАТНОГО ПРОЕЦИРО-ВАНИЯ  С ИСПОЛЬЗОВАНИЕМ ВЕЙВЛЕТ-ФИЛЬТРАЦИИ ПРОЕКЦИОННЫХ ДАННЫХ В РЕНТГЕНОВСКОЙ КОМПЬЮТЕРНОЙ ТОМОГРАФИИ}
\footnote{Если статья рекомендована к публикации программным
комитетом научной конференции, это указывается в сноске к названию статьи. Благодарность за финансовую поддержку в подготовке статьи необходимо поместить после заключения.}%


\authors{Е.Н. Симонов$^{1}$, К.М. Виноградов$^{2}$}

\organizations{$^{1}$Южно-Уральский государственный университет\\
	(454080 Челябинск, пр.~им.~В.И.~Ленина, д.~76),\\ 
	$^{2}$Южно-Уральский государственный университет\\
	(454080 Челябинск, пр.~им.~В.И.~Ленина, д.~76)}

\newcommand{\emailsinsert}{\href{mailto:e.n.simonov@yandex.ru}{e.n.simonov@yandex.ru}, \href{mailto:vinogradovkm@susu.ru}{vinogradovkm@susu.ru}}

% Дата поступления статьи в редакцию
\newcommand{\rdate}{ДД.ММ.ГГГГ}
\emails{\emailsinsert}
\redactiondate{\rdate}

\renewcommand*{\thefootnote}{\arabic{footnote}}
\setcounter{footnote}{0}

% оформление аннотации
\begin{abstract}%
В статье представлен метод уменьшения ошибки реконструкции изображения для рентгеновской компьютер-ной томографии путем применения вейвлет-фильтрации зашумленных проекционных данных. Вейвлет преобразование и основанное на нем вейвлет-фильтрация одномерных сигналов дает возможность определять конкретное место соот-ветствия частотной и временной (в данном случае пространственной по координате детекторов) области. Это позволяет однозначно определять переход  из частотной области в пространственную и обратно. Для фильтрации проекционных данных используется вейвлет преобразование, которое дает возможность через коэффициенты, определяющие мас-штабирующие функции и функции вейвлетов определять в частотной и пространственной области место шума в за-шумленном сигнале и осуществлять выделение не зашумленного сигнала путем назначения порогов фильтрации на вышеуказанные коэффициенты. Для усиления фильтрующих свойств вейвлет преобразования предложено разбивать проекционные данные на интервалы, для каждого из которых определяются свои коэффициенты.     Вейвлет-фильтрация проводится с использованием  вейвлетов Добеши. Результаты исследований были подтверждены  матема-тическим моделированием  зашумленных проекционных данных, их вейвлет-фильтрации и реконструкции по ним те-стового томографического  изображения.  Математическая модель тестового объекта исследования  и разработанный авторами программный реконструктор томографического изображения позволили осуществлять моделирование пря-мой (получение проекционных данных  по тестовому объекту), обратной (получение тестового томографического изображения по проекционным данным объекта) задач томографии  и  осуществлять сравнительный анализ качества реконструкции изображения с «идеальными» и зашумленными проекционными данными.

\keywords{рентгеновская компьютерная томография, проекционные данные, вейвлеты}
\end{abstract}

%Образец цитирования статьи
\forcitation{Симонов~Е.Н., Виноградов~К.М. Реконструкция изображения по методу обратного проецирования с использованием вейвлет-фильтрации проекци-онных данных в рентгеновской компьютерной томографии}

%Верхние колонтитулы
\markboth{Симонов~Е.Н., Виноградов~К.М. Реконструкция изображения по методу обратного проецирования с использованием вейвлет-фильтрации проекци-онных данных в рентгеновской компьютерной томографии\ldots}
% Заголовки разделов формируются при помощи команд  \section{}, \subsection{}, \subsubsection{}

\section{Введение}
\label{sec-intro}
Повышение точности реконструкции изображений в рентгеновской компьютерной томографии, со времени ее открытия и получения Нобелевской премии,  является одной из важных задач в достижении пространственного и плотностного разрешения томогра-фических систем. Точность реконструкции зависит от факторов:
\begin{enumerate}
\item Устойчивости решения томографического уравнения через обратное преобразова-ние Радона (обратная задача томографии), т.к. эта задача относится к классу не-корректно поставленных  задач и требует регуляризации в решении.

\item Погрешности рентгено-оптического тракта томографа при измерении сигналов де-текторной системой, вызванной различными причинами: погрешностью самой электронной системы измерения, квантовыми шумами детектора, неоднородно-стью поля излучения источника, изменением геометрии тракта при сканировании объекта.
\end{enumerate}


В данной статье рассматривается задача уменьшения влияния погрешности второго фактора, который непосредственно определяет  точность формирования проекционных данных, участвующих в реконструкции томографического изображения. Влияние первого фактора в статье не рассматривается. Считается, что алгоритм реконструкции является корректным и регуляризированным.

Структура данной работы такова. Раздел ~\ref{sec-first} содержит обзор работ, посвященных фильтрации томографических проекционных данных.  В разделе ~\ref{sec-second} дано  описание разра-ботанного программного реконструктора, как инструмента  исследования влияния шумов на качество томографического изображения. Раздел ~\ref{sec-third} содержит предлагаемый метод фильтрации проекционных данных. Вычислительные эксперименты и их результаты пред-ставлены в разделе ~\ref{sec-fourth}. Последний раздел~\ref{sec-conclusion} содержит заключение и выводы исследования.

\section{Обзор работ}
\label{sec-first}
Вопросам фильтрации в рентгеновской компьютерной томографии посвящено до-статочно много работ. Однако основная их часть посвящена вопросам фильтрации самого алгоритма реконструкции, чтобы он был устойчивым к малым флюктуациям проекцион-ных данных, т.е. фильтрация используется для регуляризации обратного преобразования Радона при решении основного уравнения рентгеновской томографии [\cite{Luitt01}, \cite{LuisNetter02}, \cite{BarrettJE03}, \cite{ArseninVYa04}, \cite{ShiH22}], или вопросам  фильтрации уже полученного томографического изображения [\cite{PikalovVV05}, \cite{VoskoboinikovYuE06}, \cite{VoskoboinikovYuE07}, \cite{VoskoboinikovYuE08}, \cite{Simonov09}, \cite{Simonov10}, \cite{Simonov11}, \cite{Bessonov12}]. 

В работах [\cite{Simonov09}, \cite{Simonov10}, \cite{Simonov11}, \cite{Bessonov12}, \cite{Bessonov13}, \cite{Klonov14}, \cite{StaroverovNE15}, \cite{ObodovskiyAV16}] было показано, что  погрешности детекторной системы измерения проекционных данных такие, как квантовый шум детекторов, статистические скачки изме-рения единичный детекторов или группы детекторов (единичных проекционных данных или группы), статистические  пространственные колебания поля источника излучения или колебания геометрических параметров рентгенооптического тракта (горизонтальное или вертикальное смещение детекторов) приводят к различным шумам и артефактам на томо-графическом изображении. Уменьшить их влияние на качество изображения с применени-ем пространственных и частотных методов фильтрации изображений практически не представляется возможным [\cite{VoskoboinikovYuE06}, \cite{VoskoboinikovYuE07}, \cite{VoskoboinikovYuE08},\cite{Simonov09}]. Поэтому с таким влиянием указанных погрешностей на качество изображения необходимо бороться на самом начальном этапе томографического процесса, на этапе формирования проекционных данных.

В работах [\cite{Simonov10}, \cite{Simonov11}]  предлагалось  осуществлять фильтрацию проекционных данных  с применением преобразования Фурье (фурье-фильтрация). Однако результаты исследования 
неудовлетворительные. Причинами  неточной фильтрации проекционных данных могут  быть: ---  преобразование Фурье, базисной функцией которого являются гармонические колебания, может давать значительные погрешности  при скачкообразных измерительных данных с детекторов (эффект Гиббса);

Фурье преобразование имеет и другой фундаментальный недостаток, как  невозмож-ность указать соответствие локальной частотной области сигнала его локальной временной (пространственной) области. Это отражается  на диагностике (определении) локаль-ных областей зашумленных детекторов (измерений), определении частотного спектра  их фильтрации. 

Работ по фильтрации проекционных данных в компьютерной томографии доста-точно мало [\cite{Simonov09}, \cite{Simonov11}, \cite{Bessonov12}, \cite{Bessonov13}, \cite{Klonov14}, \cite{StaroverovNE15}, \cite{ObodovskiyAV16}]. Они имеют специальный характер, направлены на уменьшение влияния конкретного фактора, приводящего к специфическому артефакту на изображении, как правило, физическими методами.  

Главной причиной проблемы фильтрации проекционных данных в томографии яв-ляется сложный многомерный  алгоритм реконструкции томографического изображения [\cite{Simonov09}]. Проследить путь от шумящего по какой-либо причине конкретного детектора (группы детекторов) до какого-то конкретного пикселя (группы пикселей) томографического изоб-ражения не представляется возможным в силу  многократного преобразования проекцион-ных данных (сотни тысяч операций) с различными их весами в силу физики поглощения излучения,  их усиления, сверточного    преобразования (усечения высоких частот). На погрешность конкретного пикселя изображения влияют все проекции (шумящие и не шу-мящие), полученные по всем детекторам и ракурсам, и определение  влияния шума про-екционных данных (или какой-то их части) на точность реконструкции изображения явля-ется сложной задачей. Это является главным отличием от «классического» измерительно-го канала с помехой [\cite{Simonov10}].

Общего математического подхода к фильтрации  проекционных данных, который бы давал возможность выделять требуемую структуру изменения сигнала, определять ло-кальные области зашумленных детекторов и применять к ним локальные фильтры, в настоящее время не выработано. Такими возможностями,  в определенной мере, обладает фильтрация на основе вейвлет преобразования [\cite{DaubechiesI17}]. Фильтрация одномерных сигналов на основе вейвлет преобразования пока не нашла должного применения в практических при-ложениях, как фурье-фильтрация. Основная причина, видимо, достаточно сложный математический аппарат вейвлет преобразования и практическое отсутствие разработанных прикладных методик использования его уникальных свойств.

\section{Описание инструмента исследования – программного реконструктора томографического изображения}
\label{sec-second}

Для исследования влияния погрешности второго фактора на томографическое изображение авторами статьи был разработан программный реконструктор томографиче-ских рентгеновских изображений, позволяющий решать прямую задачу томографии, т.е. определять массив проекционных данных для различных ракурсов облучения и отсчетов линейки детекторов  по модели объекта исследования, и обратную задачу, т.е. по задан-ному массиву проекционных данных проводить томографическую реконструкцию изоб-ражения модели объекта. 

Таким образом, программный реконструктор, как инструмент исследования, дает возможность по модели объекта исследования получать «идеальные» проекционные дан-ные, накладывать на них общие и локальные  (по ряду ракурсов или отсчетов) шумы, проводить реконструкцию томографического изображения на зашумленных проекционных данных, оценивать качество изображения (ошибку реконструкции) относительно «идеаль-ного» изображения.

При разработке программного реконструктора  использовались теоретические ос-новы рентгеновской компьютерной томографии [\cite{Luitt01}, \cite{Simonov09}]. С помощью рентгеновской трубки (источника) получают пучки рентгеновских лучей интенсивности $I_0$, которые проходят через объект исследования и регистрируются детекторами.  На рис. \ref{fig1} показана параллель-ная схема сканирования рентгеновского пучка излучения.

\begin{figure}[h]%
	\centering%
	\includegraphics[scale=0.8]{pic/fig-1.eps}%
	\caption{Параллельная схема сканирования рентгеновского пучка излучения}%
	\label{fig1}%
\end{figure}%


Закон Бугера-Ламберта-Бера определяет ослабление монохроматического пучка излучения при распространении его в объекте исследования:

\begin{equation}
	\label{eq1}
	I(l) = I_0exp(-\int_{-\infty}^{+\infty} \ \mu_l dl)
\end{equation}

где $I_0$~--- интенсивность входящего пучка в объект исследования,
$I(l)$~--- интенсивность выходного пучка из объекта исследования,
$dl$~--- элементарное расстояние, которое проходит рентгеновский луч через вещество по прямой l,
$\mu_l$~---  линейный коэффициент ослабления рентгеновского излучения на элементар-ном расстоянии $dl$.

Обозначив через $q(l, \theta)$ проекционные данные для положения $l$ детектора и угла поворота $\theta$ относительно началь-ного положения, получим из (\ref{eq1}):

\begin{equation}
	\label{eq2}
	q(l, \theta) = - ln(\frac{I(l,\theta)}{I_0})=(\int_{-\infty}^{+\infty} \ \mu_l dl)_l,_\theta
\end{equation}

Запись (\ref{eq2}) есть основное уравнение рентгеновской компьютерной томографии и  означает, что интеграл определяется для каждого значения $l$ и $\theta$.

В полученных проекционных данных присутствует погрешность, потому что опре-деленный уровень погрешности всегда имеет место в любом детекторной электронной системе, которая измеряет интенсивность рентгеновского пучка на входе $I_0$ и выходе  $I(l)=I(l, \theta)$ объекта исследования. Учитывая, что множество вышеуказанных факторов, вли-яющих на уровень погрешности, являются независимыми и ни один из них не преобладает над остальными, распределение погрешности целесообразно выбрать нормальным, то есть в виде гауссовского шума.  

Согласно работе [\cite{Simonov09}], зашумленные проекционные данные можно представить в виде суммы проекционных данных без шума и аддитивного гауссовского шума 

\begin{equation}
	\label{eq3}
	`q = q+\eta
\end{equation}

где	$q`$~--- зашумленные проекционные данные;
$q$~--- проекционные данные без шума;
$\eta$~--- гауссовский шум.
Реконструкция томографического изображения в параллельных лучах для инте-грального аналитического алгоритма обратного проецирования с фильтрацией сверткой (Filtered Back Projection FBP) [\cite{Luitt01}] представляет собой определение  в (\ref{eq2}) $\mu_l$  в координатах $(x,y)$. Такое решение представляется, как обратное преобразование Радона уравнения (\ref{eq2})

\begin{equation}
	\label{eq4}
	\mu(x,y) = \int_{0}^{x} \int_{-\infty}^{+\infty}g(l`-l)q`(l,\theta)dl d \theta
\end{equation}

где  $g(l`-l)$~--- ядро свертки (сворачивающая функция) для параллельной геометрии, где
$q(l)=\frac{1}{4\pi}\int_{-\infty}^{+\infty}|w|W(w)e^{iws} dw$~, где $w$~--- пространственная частота, $W(w)$~--- окно фильтрации,
$q`(l,\theta)$---  исходные зашумленные проекционные данные,
$l`$~--- расстояние между точкой $Q (x,y)$ и прямой $ОЕ$ (рис. \ref{fig1}).
В работе [\cite{Simonov09}] автором данной статьи для множества различных окон фильтрации $W(w)$ (прямоугольного окна, синусного, параболического и др.) определены сворачивающие функции $g(l)$ (ядра свертки). Так, например, для прямоугольного окна сворачивающая функция будет иметь вид:

\begin{equation}
	\label{eq5}
	g(l)=\frac{1}{2\pi^2}(\frac{w_{np}}{l}sin(w_{np}l-\frac{2}{l^2}sin^2(\frac{w_{np} l}{2})),
\end{equation}

где пространственная частота $w_{np}=\frac{1}{2\Delta l}$.

Этот вид сворачивающей функции используется в программном реконструкторе для по-лучения томографического изображения.
Схема получения томографического изображения функции (\ref{eq4}) (обратная задача томографии) по заданным проекционным данным (\ref{eq2}) показана на рис. \ref{fig2}.

\begin{figure}[h]%
	\centering%
	\includegraphics[scale=0.25]{pic/fig-2.eps}%
	\caption{Схема получения томографического изображения}%
	\label{fig2}%
\end{figure}%

Программный код реконструкции томографического изображения по заданным проекционным данным на языке С++ отражен в [\cite{Simonov21}].
Схема получения проекционных данных (\ref{eq2}) (прямая задача томографии) показана на рис. \ref{fig3}.

\begin{figure}[h]%
	\centering%
	\includegraphics[scale=0.33]{pic/fig-3.eps}%
	\caption{Схема получения проекционных данных}%
	\label{fig3}%
\end{figure}%

Программный код моделирования проекционных данных для параллельной геометрии сканирования на языке С++ отражен в [\cite{Simonov21}].

\section{Предлагаемый метод фильтрации}

Для фильтрации зашумленных проекционных данных использовалось их   вейвлет  преоб-разование. Вейвлеты --– это обобщенное название семейств математических функций определенной формы, которые локальны во времени (в пространстве) по амплитуде и  частоте и в которых все функции получаются из одной базовой (порождающей) функции посредством ее сдвигов и растяжений по оси времени (пространства) [\cite{DaubechiesI17}]. Дискретное вейвлет преобразование любой функции  $f(x)$ можно представить, как

\begin{equation}
	\label{eq6}
	P_jf(x)=\sum_{n}a_{j,n} \phi_{j,n}(x)=\sum_{n}a_{j+1,n}(x)+\sum_{n}d_{j+1,n}\Psi_{j+1,n}(x)=
	\sum_{n}a_{j+2,n}\Psi_{j+2,n}(x)+\sum_{n}d_{j+2,n}\Psi_{j+2,n}(x)+\sum_{n}d_{j+1,n}\Psi_{j+1,n}
\end{equation}

причем  $\lim\limits_{j\to-\infty} P_jf(x)=f(x)$

Такое преобразование можно продолжить далее (рис. \ref{fig4}).

\begin{figure}[h]%
	\centering%
	\includegraphics[scale=0.5]{pic/fig-4.eps}%
	\caption{Схема вейвлет преобразования}%
	\label{fig4}%
\end{figure}%

Функция $\varphi(x)$ называется  масштабирующей функцией  и определяется, как

\begin{equation}
	\label{eq7}
	\varphi(x)=\sum_{n}h_n \varphi_{-1,n}(x)=\sqrt{2}\sum_{n}h_n\varphi(2x-n),
\end{equation}
где $h_n$~--- коэффициенты фильтра для применяемого вейвлета, $n={0,1,2…}.$

Функция $\varphi(x)$ называется функцией вейвлета и определяется, как:

\begin{equation}
	\label{eq8}
	\psi(x)=\sqrt{2}\sum_{n}g_n\varphi(2x-n)
\end{equation}

где $g_n=(-1)^2h_{1-n}$.

{$a_j$},{$d_j$} в (\ref{eq6}) --- коэффициенты дискретного вейвлет преобразования,  где $j$~–-- уровень разложения сигнала на низкочастотные и высокочастотные составляющие.
Соответственно, 

\begin{equation}
	\label{eq9}
	\phi_{j,n}=2^{-j/2}\phi(2^{-j} x-n)=\sum_{n}h_n\phi_{j-1,1+2n}(x)
\end{equation}

\begin{equation}
	\label{eq10}
	\psi_{j,n}=2^{-j/2}\psi(2^{-j} x-n)=\sum_{n}g_n\psi_{j-1,1+2n}(x)
\end{equation}

Дискретное вейвлет преобразование проводит перевод последовательности  ${a_{j,n}}$ в последовательности  ${a_{j+1,n}, d_{j+1,n}}$. Множества  ${h_n}$, ${g_n}$   позволяют по формулам (\ref{eq8}, \ref{eq9}, \ref{eq10}) численно находить коэффициенты дискретного вейвлет преобразования. 
Применяя преобразование (\ref{eq6}) для нашей задачи, получим, что проекционные данные $q`(l,\theta)$ , можно представить, как

\begin{equation}
	\label{eq11}
	q`(l,\theta)=\sum_{n=1}^{\infty}a_{j,n}\phi_{j,n}(l,\theta)+\sum_{j=1}^{\infty}\sum_{n=1}^{\infty}d_{j,n}\Psi_{j,n}(l,\theta),
\end{equation}

где: $a_{j,n}, d_{j,n}$~--- соответственно, аппроксимирующие и детализирующие коэффициенты вейвлет преобра-зования;
$\phi_{j,n}$~--- масштабирующая функция;
$\psi_{j,n}$~--- функция вейвлета.

Таким образом, анализируя (\ref{eq6}) и рис. \ref{fig3}, аппроксимирующую   часть функции мож-но многократно разложить на еще более сглаженную (аппроксимирующую) часть и де-тализирующую часть.

Считаем, что шум содержится в высокочастотных компонентах, т.е. находится в коэффициентах вейвлет преобразования, отвечающих за малые масштабы $(d)$. Фильтрация сводится к обрезанию “высоких частот” всего массива проекционных данных или его ло-кальной части, т.е. приравниванию к нулю коэффициентов $d$ ниже заданного значения (порога). В результате ожидается, что структура обрабатываемых данных, которая лежит в основе проекционного сигнала, сохранится и проявится после очищения от общего или локального  шума [\cite{DonohoD18}, \cite{BirgeL19}, \cite{ChangS20}]. 

Учитывая  фундаментальное свойство вейвлет преобразования, как соответствие временной (пространственной) и частотной областей  сигнала в каждой точке, имеется возможность локального анализа сигнала и его шума.  Чтобы усилить локальные свойства вейвлет преобразования, авторами работы было предложено разбиение проекционного сигнала на интервалы с высоким изменением  амплитуды, где амплитуда шума мало из-меняет общую тенденцию сигнала, и фильтрация шума возможна без искажения сигнала, и интервалы с малым изменением амплитуды сигнала, где амплитуда шума изменяет общую тенденцию сигнала, и фильтрация шума может быть  не возможна без искажения сигнала. Для каждого интервала проводится вейвлет преобразование и определяются свои аппроксимирующие а и детализирующие $d$ коэффициенты.

Критерием разбиения проекционного сигнала на интервалы является отношение

\begin{equation}
	\label{eq12}
	K=\frac{{A_N}^{max}}{\Delta A},
\end{equation}
	
где  ${A_N}^{max}$ - максимальная амплитуда шума в интервале,
${\Delta A}$ – изменение амплитуды  проекционного сигнала без шума  в интервале.

Для интервала, где коэффициент $0<K\leq1$, фильтрация шума возможна без искажения сигна-ла более жесткими фильтрами и порогами, дающие меньший эффект сглаживания сигнала. Для интервала, где $K>1$, фильтрация шума требует применение мягких фильтров и поро-гов.  

Таким образом, алгоритм очистка проекционных данных с применением вейвлет преобразования сводится к следующим шагам:
\begin{enumerate}
	\item Проведение вейвлет преобразования для каждого интервала проекционных данных.
	\item Задание порога для коэффициентов вейвлет преобразования $d$ в соответствии с их уровнем разложения $j$ для всего сигнала или его интервала.
	\item Проведение фильтрации для интервалов проекционных данных  с  заданным порогом.
	\item Восстановление проекционных данных по измененным коэффициентам $d$. 
\end{enumerate}

Существуют различные виды задания порога для коэффициентов вейвлет преобразования  $d$ [\cite{DonohoD18}]. Пусть $d_j,{d`}_j$  – вейвлет-коэффициенты некоторого уровня  $j$ до фильтрации и после, соответственно, $t$ – значение порога, тогда:

мягкая вейвлет-фильтрация определяется следующим образом

\begin{equation}
	\label{eq13}
	{d`}_j=sign(d_j) max(0,|d_j|-t);
\end{equation}

жесткая вейвлет-фильтрация

\begin{equation}
	\label{eq14}
	{d`}_j=
	 \begin{cases}
	 	(0,|d_j|\leq{t});\\
	 	d_j,|d_j|>t
	 \end{cases}
\end{equation}
аффинная вейвлет-фильтрация

\begin{equation}
	\label{eq15}
	{d`}_j=
	\begin{cases}
		0,|d_j|<t/2\\
		2d_j+t, -t\leq{d_j}\leq{-t/2}\\
		2d_j-t, t\leq{d_j}\leq{t}\\	
		d_j,|d_j|>t
	\end{cases}
\end{equation}

Значение порога $t$ играет важную роль. Слишком малое значение порога не может удалить шумовые составляющие, слишком большое значение может удалить полезные составля-ющие проекционных данных. Существуют следующие способы  поиска порога [\cite{DonohoD18}, \cite{BirgeL19}, \cite{ChangS20}]:

\begin{enumerate}
	\item Универсальный порог
	\begin{equation}
		\label{eq16}
		t=\sigma\sqrt{2log_2(n)},
	\end{equation}
	
	где $n$~--- число отсчетов сигнала, 
	$\sigma$~--- уровень шума. 
	\item Минимаксный порог 
	\begin{equation}
		\label{eq17}
		t=\sigma t_n,
	\end{equation}
	где $t_n$~--- выбирается, исходя из минимаксного правила, т.е. $inf(sup(E(f,f`)) , где E(f,f`)$~--- оценка среднеквадратичной ошибки между искомой и найденными функциями.
	\item Порог, устанавливаемый на основании правила $Birge˗Massart$,
	\begin{equation}
		\label{eq18}
		t=|c_k|,
	\end{equation}
	где $k=min(-sum({c_j}^2)+2{\sigma}^2 n`(\alpha+log_2(n/n`))$, 
	$c_j$-   вейвлет-коэффициенты $d_j$, отсортированные в порядке убывания своих абсолютных значений, $n`$--- число отсчетов сигнала, $1\leq{n`}\leq{n}$.
	\item Порог, устанавливаемый на основании правила $Donoho˗Johnstone$
	\begin{equation}
		\label{eq19}
		t=\frac{\sigma}{{\sigma}_x},
	\end{equation}
\end{enumerate}
	где ${\sigma}_x$~--- локальная ошибка исследуемой функции, $\sigma$ определяется, как 
	\begin{equation}
		\label{eq20}
		\sigma=\frac{med(|W_1(f`)|)}{0,6745},
	\end{equation}
	где  $W_1({f`})$~--- нижний уровень детализирующих вейвлет-коэффициентов.

Реконструируя функцию  $\mu(x,y)$ по формуле (\ref{eq4}) по зашумленным проекциям ${q`}(l,\theta)$  и фильтрованным $q(l,\theta)$ , можно оценить качество томографического изображения и сделать сравнительный анализ.

Таким образом, последовательность реконструкции томографического изображе-ния с использованием вейвлет-фильтрации проекционных данных будет следующей (рис. \ref{fig5}).

\begin{figure}[h]%
	\centering
	\includegraphics[scale=0.5]{pic/fig-5.eps}
	\caption{Схема реконструкции томографического изображения с использованием вейвлет-фильтрации проекционных данных}
	\label{fig5}
\end{figure}

Поясним  принципиальные особенности вейвлет-фильтрации проекционных данных (шаг 4 схемы рис. \ref{fig5}). Проекционные данные $q (l,\theta)$  в параллельных рентгеновских лучах (рис. \ref{fig1}) представляются в виде матрицы $N×M$ (рис. \ref{fig6),
где  $N$~--- количество детекторов, $l= \deltal×n$, $n$~--- порядковый индекс (номер) детектора, $n={1,2,3… N}$,
$M$~--- количество ракурсов облучения, 
$\theta=\delta\theta×m$, $m$~--- порядковый индекс (номер) ракурса, $m={1,2,3… M}$.

\begin{figure}[h]%
	\centering
	\includegraphics[scale=0.5]{pic/fig-6.eps}
	\caption{Матрица проекционных данных}
	\label{fig6}
\end{figure}

Из матрицы проекционных данных объекта исследования (рис. \ref(fig6)) для каждого ракурса $m={1,2,3… M}$ выделяют одномерный проекционный сигнал $q(l)$  и определяют интервалы для вейвлет анализа в соответствии с критерием (\ref{eq10}).

На рис. \ref(fig7) показан типовой одномерный проекционный сигнал $q (l)$ с шумом и ин-тервалы разбиения сигнала в соответствии с критерием (\ref{eq10}).  В интервалах $1-7,9$ значения $К<1$, в интервале $8$ $К>1$. Это значит, что для интервалов $1-7,9$ целесообразно применять жесткую или аффинную  вейвлет-фильтрацию (\ref{eq9}), (\ref{eq10}),  т.к. шум мало оказывает влияние на общую тенденцию изменения сигнала.  Порог для этих интервалов возможен по (\ref{eq14}), (\ref{eq15}), т.к. здесь не требуется адаптации к изменению  сигнала, ибо это изменение явно вы-ражено. Для интервала 8 целесообразно применять мягкую вейвлет-фильтрацию (\ref{eq11}), а порог возможен по правилу $Birge и Massart$ (\ref{eq16}) или правилу $Donoho и  Johnstone$ (\ref{eq17}), т.к. сигнал сильно зашумлен, изменения его малозаметны на фоне шума, и здесь требуются некоторые элементы адаптации к изменению сигнала при выборе типа вейвлет-фильтрации и порога, что и заложены в (\ref{eq11})), (\ref{eq16}) и (\ref{eq17}). В противном случае,  жесткая фильтрация  и «жесткие» пороги  могут удалить полезные составляющие сигнала.

\begin{figure}[h]%
	\centering
	\includegraphics[scale=2.1]{pic/fig-7.eps}
	\caption{Типовой одномерный проекционный сигнал $q (l)$ с шумом}
	\label{fig7}
\end{figure}

\section{Вычислительные эксперименты }

Для проведения вычислительного эксперимента использовался разработанный авторами статьи  программный реконструктор томографического изображения.  Реконструктор позволяет получать матрицы изображений  размером  $512×512$ пикселей.  Объем проекционных данных может моделироваться для $M=600$ ракурсов облучения (угол $/theta$), количества единичных отсчетов ($\delta l$) в линейке детекторов $N=512$, т.е. получать матрицу проекционных данных  $M×N (512×600)$.  

Для количественной оценки шума на томографическом изображении использова-лась среднеквадратическая ошибка:

\begin{equation}
	\label{eq21}
	E=(\frac{\sum_{i}{\epsilon}\sum_{j}{\epsilon}(A_{i,j}-B_{i,j})^2}{\sum_{i}{\epsilon}\sum_{j}{\epsilon}(A_{i,j}})^{1/2}
\end{equation}

где $A, B$~--- значение матрицы исходного и реконструируемого изображений,
$i, j$~--- номер соответствующего элемента  (пикселя) в матрице изображения,
$\epsilon$ - количество элементов  (пикселей) в матрице.

Анализ влияния вейвлет-фильтрации проекционных данных на томографическое изобра-жение был проведен для следующих вариантов:

\begin{enumerate}
	\item Вейвлет-фильтрация без разбиения одномерного проекционного сигнала на интер-валы с заданием мягкой (Мвф) и жесткой (Жвф) фильтрации, с заданием порогов по Birge˗  Massart (BM) и по Donoho˗Johnstone (DJ).
	\item Вейвлет-фильтрация с разбиением одномерного проекционного сигнала на интер-валы в соответствии с критерием (\ref{eq10}с заданием типов фильтрации и порогов, что в п.1.
\end{enumerate}

Для фильтрации в обоих вариантах использовался вейвлет Добеши типа 8 [\cite{DaubechiesI17}], уровень разложения проекционных данных взят $j=3$.

Для анализа была разработана тестовая математическая модель объекта с высоко контрастными прямоугольными и низко контрастными круговыми вставками  (рис.\ref(fig8)). На модели объекта  проверялась  точность работы программного реконструктора.  Из рис. \ref(fig8) видно, что точность реконструкции высокая: отсутствуют на томографическом изображении шумы и артефакты.  

\begin{figure}[!ht]
	\begin{minipage}[t]{1.0\textwidth}
		\centering
		\includegraphics[width=0.35\linewidth]{pic/fig-8-1.eps}
		\includegraphics[width=0.35\linewidth]{pic/fig-8-2.eps}
		\\ а. изображение модели объекта исследования  б. реконструированное томографическое 
		изображение модели	
	\end{minipage}
	\caption{Изображение модели объекта исследования  и реконструированное томографическое изображение модели  }
	\label{fig8}
\end{figure}
		
Последовательность операций по оценке влияния вейвлет-фильтрации зашумленных про-екционных данных (одномерных проекционных сигналов $q (l)$ для каждого ракурса) на качество томографического изображения была следующей:

\begin{enumerate}
	\item  По исходному тестовому изображению  (рис. \ref{fig8}) моделировались идеальные проекционные данные в соответствии с (\ref{eq2}).
	\item  На идеальные проекции аддитивно накладывался гауссовский шум с   с получением зашумленных проекционных данных в соответствии с (\ref{eq3}).
	\item  Проводилась реконструкция изображения с зашумленными данными в соответствии с (\ref{eq4}	\item  Проекционные данные разбивались на интервалы в соответствии с критерием (\ref{eq10}).
	\item  Для  каждого интервала определялись детализирующие коэффициенты  , тип их фильтрации  в соответствии с (\ref{eq11}) - (\ref{eq14}) и пороги в соответствии с (\ref{eq14}) --– (\ref{eq18})).
	\item  Проводилась реконструкция изображения в соответствии с (\ref{eq4}) для отфильтрованных проекционных данных.
	\item  Проводился сравнительный анализ с оценкой (\ref{eq19}томографических изображений, ре-конструируемых без применения вейвлет-фильтрации (п.3), с вейвлет-фильтрацией без разбиения проекционного сигнала на интервалы и с вейвлет-фильтрацией с разбиением проекционного сигнала на интервалы.
\end{enumerate}

Результаты реконструкции и оценка качества томографического изображения – среднеквадратическая ошибка $Е$, определенная по (\ref{eq19}), приведены на рис. \ref{fig9}.

\begin{figure}[!ht]
	\begin{minipage}[t]{1.0\textwidth}
		\centering
		\includegraphics[width=0.3\linewidth]{pic/fig-9-1.eps}
		\includegraphics[width=0.3\linewidth]{pic/fig-9-2.eps}
		\includegraphics[width=0.3\linewidth]{pic/fig-9-3.eps}
		\\ а. Без вейвлет-фильтрации ($Е=20,86\%$)  	 
		\\ б. С использованием вейвлет-фильтрации без разбиения проекционного сигнала для мягкой фильтрации	с заданием порогов по Birge˗Massart ($Е=7,42\%$)     
		\\ с. С использованием вейвлет-фильтрации без разбиения проекционного сигнала для мягкой фильтрации	по по Donoho˗Johnstone ($Е=10,25\%$)  
	\end{minipage}
	
		\begin{minipage}[t]{1.0\textwidth}
		\centering
		\includegraphics[width=0.3\linewidth]{pic/fig-9-4.eps}
		\includegraphics[width=0.3\linewidth]{pic/fig-9-5.eps}
		\includegraphics[width=0.3\linewidth]{pic/fig-9-6.eps}
		\\ а. Для жесткой фильтрации с заданием  порогов проекционного сигнала по Birge˗Massart ($Е=7,86\%$)   
		\\ б. Для жесткой фильтрации с заданием  порогов проекционного сигнала по Donoho˗Johnstone ($Е=10,01\%$)     
		\\ с. С разбиением проекционного сигнала  на интервалы по критерию (\ref{eq10}($Е=6,43\%$)  
	\end{minipage}
	
	\caption{Реконструированные изображения}
	\label{fig9}
\end{figure}

Результаты оценки качества изображений без применения, с применением вейвлет-фильтрации  и с применением фурье-фильтрации приведены в таблице \ref{tab1}.


\begin{table}[]
	%\resizebox{\textwidth}
	\centering
	\caption{Результаты оценки качества изображений}\label{tab1}}{
		\begin{tabular}{|llllllll|}
			\hline
			\multicolumn{8}{|l|}{Погрешность реконструкции}                                                                                                                                                                                                                                                                                                                        \\ \hline
			\multicolumn{4}{|l|}{Вейвлет-фильтрация проекционных данных без разбиения проекционного сигнала на интервалы}    & \multicolumn{2}{l|}{Вейвлет-фильтрация проекционных данных с разбиением проекционного сигнала на интервалы} & \multicolumn{1}{l|}{\multirow{5}{*}{Без  фильтрации проекционного сигнала}} & \multirow{5}{*}{Фурье-фильтрация проекционного сигнала} \\ \cline{1-6}
			\multicolumn{4}{|l|}{Тип вейвлет-фильтрации}                                                                     & \multicolumn{2}{l|}{Мвф}                                                                                    & \multicolumn{1}{l|}{}                                                       &                                                         \\ \cline{1-6}
			\multicolumn{2}{|l|}{Мвф}                              & \multicolumn{2}{l|}{Жвф}                                & \multicolumn{1}{l|}{\multirow{3}{*}{DJ}}             & \multicolumn{1}{l|}{\multirow{3}{*}{BM}}             & \multicolumn{1}{l|}{}                                                       &                                                         \\ \cline{1-4}
			\multicolumn{2}{|l|}{Тип порога}                       & \multicolumn{2}{l|}{Тип порога}                         & \multicolumn{1}{l|}{}                                & \multicolumn{1}{l|}{}                                & \multicolumn{1}{l|}{}                                                       &                                                         \\ \cline{1-4}
			\multicolumn{1}{|l|}{DJ}   & \multicolumn{1}{l|}{BM}   & \multicolumn{1}{l|}{DJ}    & \multicolumn{1}{l|}{BM}    & \multicolumn{1}{l|}{}                                & \multicolumn{1}{l|}{}                                & \multicolumn{1}{l|}{}                                                       &                                                         \\ \hline
			\multicolumn{1}{|l|}{7.42} & \multicolumn{1}{l|}{7.86} & \multicolumn{1}{l|}{10.25} & \multicolumn{1}{l|}{10.01} & \multicolumn{1}{l|}{6.43}                            & \multicolumn{1}{l|}{6.48}                            & \multicolumn{1}{l|}{20.86}                                                  & 20.4                                                    \\ \hline
		\end{tabular}%
	}
\end{table}

При анализе исходного тестового томографического изображения и изображений, полученных при реконструкции с зашумленными проекционными данными, видно, что применение вейвлет-фильтрации проекционных данных  может снизить ошибку рекон-струкции (с 20.86 \%  до 6.43 \%). Наилучшие результаты получены с вейвлет-фильтрацией проекционных данных с разбиением проекционного сигнала на интервалы. Фурье-фильтрация проекционных сигналов дает, в лучшем случае, качество томографического изображения на уровне без фильтрации [\cite{Simonov10}, \cite{Simonov11}]. Это объясняется  большой степенью сглаживания сигнала, а также фундаментальными недостатками Фурье преобразования.
Методология практического применения вейвлет-фильтрации проекционных дан-ных для конкретных рентгеновских томографических систем может быть следующей:

\begin{enumerate}
	\item На тестовом калибровочном физическом фантоме плотностного разрешения [\cite{Simonov09}] определяют:
	\subitem идеальные не зашумленные проекционные данные по математической модели калибровочного физического фантома.
	\subitem иеальные зашумленные проекционные данные калибровочного физического фантома на реальной  томографической системе.
	\item Зная изменения амплитуды не зашумленных проекционных данных (а).), реальные зашумленные проекционные данные (б).) разбивают на интервалы.
	\item В каждом интервале  проводят вейвлет-фильтрацию. С отфильтрованными проек-ционными данными проводят реконструкцию томографического изображения  те-стового калибровочного физического фантома на реконструкторе реальной  томо-графической системы.
	\item Проводят анализ реального томографического изображения:  
	\subitem если реальное томографическое изображение  удовлетворяет требуемым харак-теристикам  (разрешение, нелинейность, шум) [\cite{Simonov09}], то параметры вейвлет-фильтрации проекционных данных, как тип фильтрации, тип порога, уровень раз-ложения, задают в программное обеспечение обработки проекционных данных ре-альной рентгеновской томографической системы,
	\subitem если реальное томографическое изображение  не удовлетворяет требуемым ха-рактеристикам, то проводится изменение параметров вейвлет-фильтрации и пункты 3,4 повторяют. Как правило, требуется 1-2  итерации.
\end{enumerate}



\section*{Заключение}
\label{sec-conclusion}

Повышение точности реконструкции изображений в рентгеновской компьютерной 
томографии является одной из важных задач в достижении пространственного и плот-ностного разрешения реальных томографических систем. Одним из главных факторов влияния на достижение этих физических характеристик являются погрешности рентгено-оптического тракта томографа при измерении сигналов детекторной системой, вызванной различными причинами: погрешностью самой электронной системы измерения, кванто-выми шумами детектора, неоднородностью поля излучения источника, изменением гео-метрии тракта при сканировании объекта. Все эти причины влияют на погрешность полу-чения проекционных данных, по которым проводится реконструкция томографического изображения с применением  сложного многомерного   алгоритма. 

На погрешность конкретного пикселя изображения влияют все проекции (шумящие и не шумящие), полученные по всем детекторам и ракурсам, и на этапе реконструкции томографического изображения фильтрацию шумящих детекторов практически осуществить невозможно. Существующие пространственные и частотные методы фильтрации изображений, в том числе и на основе  Фурье преобразования, не дают желаемого резуль-тата  при фильтрации томографических изображений и проекционных данных.

В работе представлен метод фильтрации сигналов на этапе получения проекцион-ных данных по сигналам с детекторов с применением вейвлет преобразования. Учитывая фундаментальные свойства вейвлет-преобразования, проекционный  сигнал возможно  представить  в пространственной и частотной областях и определить соответствие этим областям для  конкретной точки сигнала. Это соответствие определяется через апрокси-мирующие и детализирующие коэффициенты пространственной области, которые ответственны, соответственно, за низкие и высокие частоты сигнала  в частотной области. Фильтруя детализирующие коэффициенты, ответственные за высокие частоты сигнала,  посредством применения определенного вида фильтрации и типа порога, мы тем самым регулируем  их влияние на исходный сигнал.

Исследования показали, что если проекционный сигнал  разбить на интервалы с большими и малыми изменениями сигнала и провести отдельно в каждом интервале филь-трацию детализирующих коэффициентов, то адаптивные и фильтрующие свойства вейвлет преобразования значительно усиливаются. В этом заключается суть представленного ме-тода фильтрации проекционных сигналов. Представленный метод фильтрации проекцион-ных данных с применением вейвлет анализа показал, что ошибка реконструкции томогра-фического изображения может быть снижена в 2-3 раза. 

Метод фильтрации проекционных данных с применением вейвлет анализа может найти практическое применение при проектировании томографических систем медицин-ского, промышленного и военного назначения. Требуются проведение дальнейших иссле-дований по  применения вейвлет-фильтрации проекционных данных, как двумерного массива сигналов, а также  для уменьшения специфических артефактов на томографическом изображении в  виде расплывчатых колец, полос, затенений.


\vspace{1em}
{\it Исследование выполнено за счет гранта Российского научного фонда №23-21-10051 https://rscf.ru/project/23-21-10051/.}


\begin{biblio}	
\bibitem{Luitt01}	
Луитт~Р.М. . Алгоритмы реконструкции  с использованием интегральных преобразований // ТИИЭР. 1983. Т.~71. №~3. С.~125-148.

\bibitem{LuisNetter02}
Луис~А.К. , Неттерер~Ф. Математические проблемы реконструктивной вычислительной томографии // ТИИЭР. 1983. Т~71. №~3. С.~111-125.

\bibitem{BarrettJE03}
Barrett~J.E., Keant~N. Artifacts in CT: Recognition and Avoidance // Radio Graphics. 2004. Vol.~24. P.~1679-1691. DOI:~\href{http://dx.doi.org/10.1148/rg.246045065}{10.1148/rg.246045065}

\bibitem{ArseninVYa04}	
Арсенин~В. Я., Криксин~Ю. А., Тимонов~А. А. Метод локальной регуляризации линейных операторных уравнений I рода и его приложения // Вычислительная математика и математическая  физика. 1988. №~28, С.~793–808.


\bibitem{PikalovVV05}
Пикалов~В.В., Непомнящий~А.В. Итерационный алгоритм с вэйвлет-фильтрацией в задаче двумерной томографии // Вычислительные методы и программирование. 2003. №.~4. С.~244–253

\bibitem{VoskoboinikovYuE06}
Воскобойников~Ю.Е., Колкер~А.Б. Комбинированные алгоритмы фильтрации зашумленных сигналов и изображений //Автометрия. 2002. №~4. С.~51-60.

\bibitem{VoskoboinikovYuE07}
Воскобойников~Ю.Е., Бронников~А.В. Адаптивный алгоритм фильтрации изображений и преобразование изображений в векторный формат // Автометрия. 1990. №~1. С~124-132.

\bibitem{VoskoboinikovYuE08}
Воскобойников~Ю.Е., Белявцев~В.Г. Алгоритмы фильтрации изображений с адаптацией размеров апертуры // Автометрия. 1998. №~3. С.~81-89.

\bibitem{Simonov09}
Симонов~Е.Н., Физика визуализации изображений в рентгеновской компьютерной томографии. Челябинск: Изд-во ЮУрГУ, 2014. 479 с.

\bibitem{Simonov10}
Ласьков~В.В., Симонов~Е.Н. Методы фильтрации изображений в рентгеновской
компьютерной томографии // Вестник Южно-Уральского государственного университета, Серия «Компьютерные технологии, управление, радиоэлектроника». 2014.  №~3. С.~29-33.


\bibitem{Simonov11}
Laskov~V.V., Simonov~E.N. Reduction of ring artifacts in  computer tomography //
Biomedical Engineering. 2016. Т.~49. №~5. С.~274-277. DOI:~\href{http://dx.doi.org/10.1007/s10527-016-9547-9}{10.1007/s10527-016-9547-9}.


\bibitem{Bessonov12}
Бессонов~В. Б., Клонов~В. В., Ларионов~И. А., Староверов~Н. Е. Разработка метода коррекции металлических артефактов при томографических исследованиях  // Физические основы приборостроения.  2020. №~4(38). С.~54-59. DOI:~\href{http://dx.doi.org/10.25210/jfop-2004-054059}{10.25210/jfop-2004-054059}


\bibitem{Bessonov13}
Бессонов~В.Б.,  Потрахов~Н.Н., Ободовский~А.В. Рентгеновская томография // Фотоника. 2019. №~7. С.~688-693.DOI:~\href{http://dx.doi.org/10.22184/1992-7296.FRos.2019.13.7.688.692}{10.22184/1992-7296.FRos.2019.13.7.688.692}

\bibitem{Klonov14}
Klonov~V.V., Larionov~I.A., Bessonov~V.B., Baksheev~I.K. Development of x-ray dose sensor // AIP Conference Proceedings : 7, Saint Petersburg, 26–27 ноября 2020 го-да. Saint Petersburg, 2021.  P.~020013. DOI:~\href{https://doi.org/10.1063/5.0053146}{10.1063/5.0053146}

\bibitem{StaroverovNE15}
Staroverov~N.E., Gryaznov~A.Y., Bessonov~V.B. Research of the possibility of using neural networks to identify areas of interest in tomographic data  // AIP Conference Pro-ceedings : 6, Saint Petersburg, 28–29 ноября 2019 года. Saint Petersburg, 2020.  P.~020027. DOI:~\href{https://doi.org/10.1063/5.0013424}{10.1063/5.0013424}
 
\bibitem{ObodovskiyAV16}
Obodovskiy~A.V., Bessonov~V.B., Larionov~I.A. Features of the practical application of microfocus x-ray tomograph in biomedical engineering  // AIP Conference Proceedings, St. Petersburg, American Ins. of Physics Inc., 2019.  P.~020049. DOI:\href{https://doi.org/10.1063/1.5121974}{10.1063/1.5121974}
 
\bibitem{DaubechiesI17}
Daubechies~I. The wavelet transform,  time-frequency localization and signal analyses// JEEE Trans. Inform. Theory. 1990. no.~36. P.~961-1005.
 
\bibitem{DonohoD18}            
Donoho~D. Nonlinear solution of linear inverse problems by wavelet-vaguelette decompositions // Journal of Applied and Computational Harmonic Analysis. 1995. Vol.~2. no.~2. P.~101-126.
 
\bibitem{BirgeL19} 
Birge~L., Massart~P. From model selection to adaptive estimation // in D. Pollard (ed), Festchrift for L. Le Cam, Springer. 1997. Р.~55-88. 

\bibitem{ChangS20}
Chang~S., Yu~B., Vetterli~M. Spatially adaptive wavelet thresholding with context model-ing for image denoising // IEEE Transactions on Image Processing. 2000. Vol.~9. no.~9. Р.~1522-1531. DOI:\href{https://doi.org/10.1109/83.862630}{10.1109/83.862630}

\bibitem{Simonov21}
Симонов~Е.Н. Реконструктор томографического изображения. Свидетельство о Государственной регистрации программы для ЭВМ №2011612631. Зарегистри-ровано 31.04.2011г.

\bibitem{ShiH22}	
Shi~H., Luo~S.,~Yang Z., Wu~G. A Novel Iterative CT Reconstruction Approach Based on FBP Algorithm // PLOS One. 2015.~10(9):~e0138498. DOI: \href{https://doi.org/10.1371/journal.pone.0138498}{10.1371/journal.pone.0138498}
 
 
\end{biblio}


% Сведения об авторах (Ф.И.О.(полностью), степень, звание, кафедра (отдел), вуз, e-mail)
Симонов Евгений Николаевич, д.т.н., профессор, кафедра техники, технологии и строительства, Южно-Уральский государственный университет (национальный исследовательский университет) (Челябинск, Российская Федерация)

Виноградов Константин Михайлович, к.т.н, доцент, кафедра техники, технологии и строительства, Южно-Уральский государственный университет (национальный исследовательский университет) (Челябинск, Российская Федерация)


%\newpage

{\phantom{eee} \hrule \vskip 4 pt}
\classify{} % Код MSC не нужен
%\classify{MSC 00000} % Коды MSC см. http://www.ams.org/msc/

\title{IMAGE RECONSTRUCTION BY THE METHOD OF REVERSE PROJECTION USING WAVELET FILTERING OF PROJECTION DATA IN X-RAY COMPUTED TOMOGRAPHY}

\authorsen{E.N.~Simonov$^{1}$, V.K.~Vinogradov $^{2}$}

\organizations{$^{1}$South Ural State University (76 V.I. Lenin Ave., Chelyabinsk, 454080),\\$^{2}$South Ural State University (76 V.I. Lenin Ave., Chelyabinsk, 454080)}

\emails{\emailsinsert}
\redactiondateen{\rdate}

\begin{abstract}%
The article presents a method for reducing the error of image reconstruction for X-ray com-puter tomography by using wavelet filtering of noisy projection data. The wavelet transfor-mation and the wavelet filtering of one-dimensional signals based on it makes it possible to determine a specific place of correspondence between the frequency and time (in this case, the spatial coordinate of the detectors) region. This makes it possible to uniquely determine the transition from the frequency domain to the spatial domain and vice versa. To filter the projec-tion data, the wavelet transform is used, which makes it possible, through coefficients defin-ing scaling functions and wavelet functions, to determine in the frequency and spatial domain the place of noise in a noisy signal and to isolate a non-noisy signal by assigning filtering thresholds to the above coefficients. To enhance the filtering properties of the wavelet trans-form, it is proposed to divide the projection data into intervals, for each of which its coeffi-cients are determined. Wavelet filtering is carried out using Daubeshi wavelets. The research results were confirmed by mathematical modeling of noisy projection data, their wavelet fil-tering and reconstruction of the test tomographic image based on them. The mathematical model of the test object of the study and the software reconstructor of the tomographic image developed by the authors made it possible to simulate direct (obtaining projection data on the test object), reverse (obtaining a test tomographic image from the projection data of the ob-ject) tomography tasks and to carry out a comparative analysis of the quality of image recon-struction with "ideal" and noisy projection data.

\keywordsen{X-ray computed tomography, projection data, wavelets}
\end{abstract}

\forcitationen{Simonov E.N., Vinogradov K.M. Image reconstruction by the method of reverse projection using wavelet filtering of projection data in X-ray computed tomography }

\openaccess

\begin{biblio_lat}
\bibitem{Luitt01} 
Luitt~R.M. Reconstruction algorithms using integral transformations // TIIER. 1983. Vol.~71. No.~3. pp.~125-148.

\bibitem{LuisNetter02}
Louis~A.K., Netterer~F. Mathematical problems of reconstructive computational tomography // TIIER. 1983. Vol.~71. No.~3. pp.~111-125. 

\bibitem{BarrettJE03}
Barrett~J.E., Keant~N. Artifacts in CT: Recognition and Avoidance // Radio Graphics. 2004. Vol.~24. P.~ 1679-1691. DOI:~\href{http://dx.doi.org/10.1148/rg.246045065}{10.1148/rg.246045065}.

\bibitem{ArseninVYa04}
Arsenin~V.Ya., Kriksin~Yu.A., Timonov~A.A. Method of local regularization of linear operator equations of the first kind and its applications // Computational mathematics and mathematical physics. 1988. No.~28, pp.~793-808.

\bibitem{PikalovVV05}
Pikalov~V.V., Nepomnyashchy~A.V. Iterative algorithm with wavelet filtering in a two-dimensional tomography problem // Computational methods and programming. 2003. No.~4. pp.~244-253. 

\bibitem{VoskoboinikovYuE06}
Voskoboynikov~Yu.E., Kolker~A.B. Combined algorithms for filtering noisy signals and im-ages //Autometry. 2002. No.~4. pp.~51-60.

\bibitem{VoskoboinikovYuE07}
Voskoboinikov~Yu.E., Bronnikov~A.V. Adaptive image filtering algorithm and image con-version to vector format // Autometry. 1990. No.~1. pp.~124-132. 

\bibitem{VoskoboinikovYuE08}
Voskoboynikov~Yu.E., Belyavtsev~V.G. Image filtering algorithms with aperture size adap-tation // Autometry. 1998. No.~3. pp.~81-89.

\bibitem{Simonov09}
Simonov~E.N., Physics of image visualization in X-ray computed tomography. Chelyabinsk: SUSU Publishing House, 2014. 479~p.

\bibitem{Simonov10}
Laskov~V.V., Simonov~E.N. Methods of image filtering in X-ray computed tomography // Bulletin of the South Ural State University, Series "Computer technologies, control, radio electronics". 2014. No.~3. pp.~29-33. 

\bibitem{Simonov11}	
Laskov~V.V., Simonov~E.N. Reduction of ring artifacts in computer tomography // Biomedical Engineering. 2016. Vol.~49. No.~5. pp.~274-277. DOI:~\href{http://dx.doi.org/10.1007/s10527-016-9547-9}{10.1007/s10527-016-9547-9}.

\bibitem{Bessonov12}
Bessonov V.B., Klonov V.V., Larionov I.A., Staroverov N.E. Development of a method for correcting metal artifacts in tomographic studies // Physical fundamentals of instrumenta-tion. 2020. No.~4(38). pp.~54-59. DOI:~\href{http://dx.doi.org/10.25210/jfop-2004-054059}{10.25210/jfop-2004-054059} 

\bibitem{Bessonov13}
Bessonov~V.B., Potrakhov N.N., Obodovsky~A.V. X-ray tomography // Photonics. 2019. No.~7. pp.~688-693. DOI:~\href{http://dx.doi.org/10.22184/1992-7296.FRos.2019.13.7.688.692}{10.22184/1992-7296.FRos.2019.13.7.688.692}


\bibitem{Klonov14}
Klonov~V.V., Larionov~I.A., Bessonov~V.B., Baksheev~I.K. Development of x-ray dose sensor // AIP Conference Proceedings : 7, Saint Petersburg, November 26-27, 2020. Saint Petersburg, 2021. P.~020013.DOI:~\href{https://doi.org/10.1063/5.0053146}{10.1063/5.0053146} 

\bibitem{StaroverovNE15}
Staroverov~N.E., Gryaznov~A.Y., Bessonov V.B. Research of the possibility of using neural networks to identify areas of interest in tomographic data // AIP Conference Proceedings : 6, Saint Petersburg, November 28-29, 2019. Saint Petersburg, 2020. P~ 020027.DOI:~\href{https://doi.org/10.1063/5.0013424}{10.1063/5.0013424}

\bibitem{ObodovskiyAV16}
Obodovskiy~A.V., Bessonov~V.B., Larionov~I.A. Features of the practical application of microfocus x-ray tomograph in biomedical engineering  // AIP Conference Proceedings, St. Petersburg, American Ins. of Physics Inc., 2019.  P.~020049. DOI:\href{https://doi.org/10.1063/1.5121974}{10.1063/1.5121974}

\bibitem{DaubechiesI17}
Daubechies~I. The wavelet transform,  time-frequency localization and signal analyses// JEEE Trans. Inform. Theory. 1990. no.~36. P.~961-1005.

\bibitem{DonohoD18}            
Donoho~D. Nonlinear solution of linear inverse problems by wavelet-vaguelette decompositions // Journal of Applied and Computational Harmonic Analysis. 1995. Vol.~2. no.~2. P.~101-126.

\bibitem{BirgeL19} 
Birge~L., Massart~P. From model selection to adaptive estimation // in D. Pollard (ed), Festchrift for L. Le Cam, Springer. 1997. Р.~55-88. 

\bibitem{ChangS20}
Chang~S., Yu~B., Vetterli~M. Spatially adaptive wavelet thresholding with context model-ing for image denoising // IEEE Transactions on Image Processing. 2000. Vol.~9. no.~9. Р.~1522-1531. DOI:\href{https://doi.org/10.1109/83.862630}{10.1109/83.862630}

\bibitem{Simonov21}
Симонов~Е.Н. Реконструктор томографического изображения. Свидетельство о Государственной регистрации программы для ЭВМ №2011612631. Зарегистри-ровано 31.04.2011г.

\bibitem{ShiH22}	
Shi~H., Luo~S.,~Yang Z., Wu~G. A Novel Iterative CT Reconstruction Approach Based on FBP Algorithm // PLOS One. 2015.~10(9):~e0138498. DOI: \href{https://doi.org/10.1371/journal.pone.0138498}{10.1371/journal.pone.0138498}
\end{biblio_lat}

\end{document}
