О РЕШЕНИИ ЗАДАЧИ СИНТЕЗА МАГНИТНОГО ПОЛЯ В МР-ТОМОГРАФЕ МЕТОДОМ РЕГУЛЯРИЗАЦИИ С ОГРАНИЧЕНИЯМИ И МЕТОДОМ АППРОКСИМАЦИИ δ-ФУНКЦИЕЙ

О РЕШЕНИИ ЗАДАЧИ СИНТЕЗА МАГНИТНОГО ПОЛЯ

В МР-ТОМОГРАФЕ МЕТОДОМ РЕГУЛЯРИЗАЦИИ С ОГРАНИЧЕНИЯМИ И МЕТОДОМ АППРОКСИМАЦИИ

5-ФУНКЦИЕЙ С.В. Иванов, Н.Г. Рущенко, В.С. Сизиков, Д.Ю. Соколов, Е.В. Хомутникова

Среди многочисленных задач магнитно-резонансной томографии (МР-томографии) [1-7], [8, с. 33-62] рассмотрим задачу синтеза магнитного поля в МР-томографе, уже рассмотренную в ряде работ (1, 2, 4, 7-12 и др.).

Формулировка задачи

В МР-томографии исключительно важной является задача создания высокооднородного магнитного поля внутри МР-томографа. Относительная неоднородность напряженности Н основного (статического) поля в пределах рабочего

объема томографа должна быть АН / Н ~10-5 -10-6 [4-7]. Это позволит с высоким разрешением решать основную задачу МР-томографии - задачу реконструкции МР-изображений [5], [6, р. 270-274], [8, с. 45-51].

Как правило, задача создания высокооднородного магнитного поля решается путем введения соленоидальных компенсирующих (корректирующих) катушек различного порядка [1, с. 274], [7]. При этом параметры катушек определяются из условия компенсации членов ряда Тейлора, в который разложено магнитное поле. Однако это довольно сложная техническая (а также математическая) задача.

Рассмотрим другой подход к решению задачи создания высокооднородного поля в МР-томографе. Этот подход основан на расчете такого закона распределения тока вдоль обмотки катушки томографа, что он обеспечит высокую однородность магнитного поля внутри катушки. Та и другая задача может ставиться как прямая, или задача анализа (расчет магнитного поля по распределению тока), и обратная, или задача синтеза (расчет распределения тока и конфигурации катушек по заданному магнитному полю) [10, с. 17]. Обратная задача является значительно более сложной, чем прямая. При этом рассматривается локальный метод синтеза (к нему относится метод компенсирующих катушек [1,7]) и интегральный метод синтеза на оси катушки [2, 9, 11, 12] и во всем объеме [7].

Задача синтеза

Рассмотрим задачу синтеза магнитного поля внутри катушки МР-томографа [4], [8, с. 55-61], [9], [10, с. 17], [11, 12], причем сначала рассмотрим ее простейший вариант - синтез на оси цилиндрической катушки с бесконечно тонкой обмоткой.

Одним из первых задачу интегрального синтеза магнитного поля сформулировал с учетом ее некорректности и решил с применением метода регуляризации Тихонова К. Адамиак [9]. Эту задачу рассмотрел также Л. Луганский [2] и др. [11,12]. В данной работе дается дальнейшее развитие этой задачи.

Итак, рассмотрим следующую обратную задачу МР-томографии: определить распределение плотности тока J(а) вдоль бесконечно тонкой обмотки цилиндрической катушки по заданному магнитному полю (напряженности) Н(г) на ее оси. Это - задача интегрального синтеза магнитного поля на оси катушки МР-томографа [4]. Здесь а -

координата вдоль обмотки катушки, z - координата вдоль ее оси. Рассмотрим случай H(z ) = const.

На рис. 1 представлена цилиндрическая катушка с бесконечно тонкой обмоткой, где a е [-/, l] - расстояние вдоль обмотки катушки, z е [-/, l] - расстояние от центра катушки вдоль ее оси, l - полудлина катушки, H(z) = H = const, z е [-/, l] - заданная напряженность магнитного поля на оси катушки, J(a), a е [-l, l] - искомое распределение тока вдоль обмотки катушки, R - радиус катушки.

J

/_ a

н ( z) • / / ^R

-l 0 l z

Рис. 1. Цилиндрическая катушка с бесконечно тонкой обмоткой

Согласно [1, с. 37], обмотка в случае J (a) ^ const называется соленоидом с переменной плотностью тока или с неоднородным распределением тока. Реализовать неоднородное распределение тока можно, например, путем разбиения соленоида на несколько отдельных секций [1, с. 37]. Другие способы (создание изолированных витков с индивидуальной подводкой тока через единый источник или создание соленоида с единым, намотанным на цилиндр проводом, но имеющим на каждом витке свое сопротивление) предложены в работах [11, 12].

Вывод интегрального уравнения задачи синтеза

Как известно [1, с. 23], [7, 9], [10, с. 20], [13, с. 330-333], [14, с. 327], напряженность магнитного поля на оси тонкого витка

3 К2

Н(г) = к ,_ _=, -х< z <х,

V

R2 + ( z - a )2

где к - коэффициент пропорциональности, 3 - ток в витке, К - радиус витка, z -координата точки, в которой вычисляется Н, а - z-координата центра витка (ось z направлена перпендикулярно плоскости витка).

R^ н

0 а 1 z z

Рис. 2. Напряженность на оси тонкого витка

Для упрощения дальнейших записей будем полагать к = 1. Тогда

H ( z ) =

JR2

[2 + (z - a )2

— x < z < X.

Далее полагаем, что имеется бесконечное множество витков, намотанных на цилиндр радиуса Я и полудлины I и имеющих плотность тока J = J(а). В этом случае, интегрируя (1) по виткам, получим следующее выражение для суммарной напряженности магнитного поля в точке с координатой г на оси цилиндра с бесконечно тонкой обмоткой (см. рис. 1).

H (z)=/

J (a) R da

V[R2 + (z — a)2

— X < z < X .

Вычисление поля H (z) по заданному току J (a), согласно (2), есть прямая задача. Сделаем важное замечание: формула (2) справедлива как для точек внутри цилиндра ( | z | < l), так и вне его ( | z | > l). Выполним краткий анализ прямой задачи. Если J(a) = J = const , то интегрирование (2) дает:

H ( z ) =

z+l

z — l

Jr2 + (z +1 )2 д/R2 + (z — l )2

J , — x < z < x .

Частные случаи (3) - поле в центре катушки h (0) = 2 J + (R /1 )2

, поле на краю

катушки H(l) = J+ (R/2l)2 . Если R = l, то H(0) = 2 JV2, H(l) = 2 JV5, т.е.

напряженность в центре катушки в д/5/2 «1.58 раз больше, чем на ее краях. Если

2 3

z —^ x, то [13, с. 333] H(z) — 2 JR l/z , т.е. вне катушки поле H(z) падает как

— 3

~ z . Такая асимптотика справедлива и для тонкого витка (см. (1)), а также будет иметь место и для цилиндрической катушки при J(a) ^ const .

Если l >> R (бесконечно длинная катушка), то H(0) = 2 J, H(l) = J .

Если, например, l = 10R, то H(0) = 1.9901 J, H(R) = 1.9898 J, т.е.

| H (R ) — H (0)|/H (0) = 1.5-10 —4. Видим, что в случае длинной катушки в принципе

можно получить неоднородность поля порядка 10—4 в пределах | z | < R, однако такую катушку технически трудно изготовить.

Вышеприведенный краткий анализ прямой задачи показывает, что при J(a) = const поле H(z) падает от центра катушки к ее краям, а вне катушки

_ 3

Перейдем к рассмотрению обратной задачи. Запишем (2) в виде:

R

-l^2 + ( z — a )2

-J ( a ) da = H ( z ),

— l < z < l.

Соотношение (4) есть интегральное уравнение Фредгольма I рода, где H(z)-задаваемая правая часть (напряженность магнитного поля на оси катушки), например, H(z) = H = const, а J(a) - искомая функция (распределение тока вдоль бесконечно тонкой обмотки цилиндрической катушки МР-томографа). Задача решения уравнения (4) является некорректной [8, с. 178].

Используем безразмерные переменные: s = ^R, х = z/R, ^ = l/R . Тогда (4) можно записать в виде

l

s0

J K ( x, s ) J ( s ) ds = H ( x ),

- s0 ^ x ^ s0,

s0

где

K (x, s) =

VR+(x - s)2 ]3

- ядро интегрального уравнения.

Решение уравнения (5) и техническая реализация решения позволят в принципе создать на оси катушки МР-томографа поляризующее поле с заданным законом изменения его напряженности, в частности, H(x) = H = const .

Решение уравнения методом регуляризации с ограничениями

Рассмотрим вопрос о решении уравнения (5).

Вначале выполним качественный анализ. Заметим, что запись (5)-(6) справедлива не только для x е [-s0,s0], но и для x £ [-s0,s0], т.е. для x е (-да, да).

Пусть H(x) = const при x е [-s0, s0]. Тогда при x £ [-s0, s0] поле H(x) ^ const и

будет монотонно падать с ростом | x |, переходя при | x да в асимптотику

H(x) ~| x |-3 (см. рис. 3).

—о 0 ¿0 х

Рис. 3. Поведение напряженности поля Н (х) внутри и вне катушки

Видим, что функция Н(х) имеет разрывные производные Н'(х), Н"(х), ... при х = ± . Отсюда следует, что искомая функция 3(¿) уравнения (5) должна быть сингулярной при ^ = -¿0 + 0 и ^ = ¿0 -0, а при | 51>¿0 ток 3(5) = 0, как это следует из физической постановки задачи. Кроме того, 3(-¿) = 3(¿) , если Н(-х) = Н(х) .

Таким образом, функция 3(¿) должна быть сингулярной и симметричной, т.е. иметь качественный вид, изображенный на рис. 4.

Теперь перейдем к количественному алгоритму решения.

Уравнение (5) в принципе может быть решено аналитически. Однако аналитического решения получить не удается. Поэтому рассмотрим вопрос о численном решении уравнения (5). Учитывая некорректность уравнения (5), воспользуемся методом регуляризации Тихонова [8, 15].

В работе [9] также использовался метод регуляризации Тихонова. При этом, как показало численное решение, при малых значениях параметра регуляризации а регуляризованное решение 3а (¿) имеет большие знакопеременные флуктуации. В то же время параметр регуляризации а должен быть действительно очень мал (порядка

1

10 5 -10 10). Это обусловлено тем, что погрешность задания Н(х) отсутствует (есть только погрешность численного алгоритма). Кроме того, регуляризованное решение Аа (з) должно иметь очень большой диапазон значений (из-за сингулярности

точного решения). При этом решение Аа (з) должно быть монотонно возрастающим от центра катушки (з = 0) к ее краям (з = ± з0 ).

А

—0 0 з0 з

Рис. 4. Принципиальный вид функции А

Чтобы регуляризованное решение Аа (з) было таковым, а именно,

неотрицательным и монотонно возрастающим от центра катушки к ее краям, нужно воспользоваться методом регуляризации Тихонова с ограничениями на решение [15, с. 118].

Ограничения в данном случае состоят в том, что решение Аа (з) ищется на множестве неотрицательных монотонно невозрастающих функций [15, с. 118]. Однако такому условию удовлетворяет лишь левая половина решения Аа(з), т.е. при з е (-з0,0]. Поэтому модифицируем уравнение (5), учитывая симметрию решения: А(-з) = А(з) . Для этого запишем (5) в виде:

Н (х)=Й

0

= 1

А (з) йз

(х - з)2]; А (з) йз

з0 I

А (з) йз

Л

1 + (х - з)"

и

I

1 + (х - з) • А (-з) йз

1 + (х + з) '

Поскольку А ( з ) = А (-з ), т.е. распределение тока симметрично, то Н ( х ) = Н ( - х ), т. е. магнитное поле Н(х) также симметрично. Справедливо и обратное: если Н ( х ) = Н ( - х ), то А ( з ) = А (-з ).

В результате уравнение (5) может быть записано в виде

| Я( х, з ) А ( з ) йз = Н ( х ),

_з0

где новое ядро равно Я( х, з) =

1

- з0 < х < 0,

1

Я

[1 + (х + з)2 ]3 д/[1 + (х - з)2 ]3 Итак, А (з) будем искать в левом полупространстве путем решения уравнения (7), т.е. при з е (-з0,0], после чего в правом полупространстве А(з) = А(-з), з е [0,з0). При

этом дополнительно будем учитывать, что функция 3(5), 5 е [—,0], является неотрицательной монотонно невозрастающей (типа изображенной на рис. 5).

J

-0

0 s

Рис. 5. Ветвь функции J

Для получения численного решения Ja (s) методом регуляризации Тихонова на

множестве неотрицательных монотонно невозрастающих функций можно воспользоваться программой PTIPR на Фортране [15, с. 118, 174]. В этой программе выполняется минимизация сглаживающего функционала методом проекции сопряженных градиентов на множество неотрицательных монотонно невозрастающих функций при каждом значении параметра регуляризации а > 0 .

В [11, 12] приведены результаты численного решения уравнения (7) при H (х) = const методом регуляризации Тихонова с ограничениями на решение.

Регуляризованное решение Ja (s) при s ^ -s0 + 0 имеет крутое нарастание (тем

больше, чем меньше а ), а в остальной области практически Ja (s) = const. Технически

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

Однако, как видно из работ [11, 12], напряженность магнитного поля

не получается высокооднородной ни при каких значениях а.

К более однородному полю Н(х), как следует из [11, 12], приводит решение задачи синтеза путем приближения 3а (5) 5 -функцией. Рассмотрим этот вопрос подробнее.

Положим, что искомое распределение тока 3(5) представляется суммой двух 5 -функций, точнее, полу-5 -функций:

0

Решение путем приближения к 5 -функции

J(s) = 51 /2 (s + s0) + 51/2 (s - s0) (см. рис. 6). Тогда в соответствии с (5) и (6),

2LV[1 +

1

+

1

А А

-1-

-з0 0 з0 з

Рис. 6. Представление функции ^в) в виде суммы полу- 5 -функций

На рис. 7. представлена Н(х) для ряда значений з0, вычисленная согласно (10). Видим, что при з0 = 0.5 поле Н(х) имеет для такой простой катушки степень однородности | АН | / Н ~ 10-4 для | х |< 0.4з0.

-1.1

1.1

Рис. 7. Напряженность Н(х) вдоль оси Данный результат можно улучшить, если представить распределение тока А (з) в

виде

А(з) = А0 • [51/2 (з + з0) + 51/2 (з - з0)] + А

Подбирая з0 и отношение А / А0 можно с помощью (11) получить однородность

В работах [11, 12] А(з) представлен суммой "усеченных" 5 -функций:

|А0, з = ±з0,

А (з) =

где А0 и А - конечные константы, причем А0 >> А. В этом случае, подбирая А0 / А, также можно добиться | АН |/Н ~ 10 4 -И0 5 для | х |< 0.4з0.

X

Можно сделать следующий вывод. Если распределение тока J(s) вдоль обмотки катушки МР-томографа сделать однородным, за исключением крайних витков (т. е. типа (9), (11), или (13)), то это позволит повысить однородность магнитного поля H (х) на оси катушки. Данный способ (приближение к 5 -функции), хотя и является более простым, чем метод регуляризации Тихонова с ограничениями на решение, но приводит к более точному результату (к более однородному полю), т. е. дает более адекватное математическое описание задачи и, как следствие, более точное решение.

Литература

1. Монтгомери Д.Б. Получение сильных магнитных полей с помощью соленоидов. М.: Мир, 1971.

2. Луганский Л.Б. Расчет соленоида с заданной формой магнитного поля // Журн. техн. физики. 1985. Т. 55. Вып. 7. С. 1263-1271.

3. Тихонов А.Н., Арсенин В.Я., Тимонов А.А. Математические задачи компьютерной томографии. М.: Наука, 1987.

4. Лухвич А.А., Чурило В.Р. Источники поляризующего магнитного поля и его градиентов для я.м.р.-томографии (обзор) // ПТЭ. 1987. № 5. С. 172-173.

5. Лич М. Получение ЯМР-изображений с пространственной локализацией // Физика визуализации изображений в медицине. / Под ред. С. Уэбба. М.: Мир, 1991. Т. 2. С. 105-231.

6. Cho Z.H., Jones J.P., Singh M. Foundations of medical imaging. NY: Wiley, 1993.

7. Галайдин П.А., Замятин А.И., Иванов В.А. Расчет и проектирование электромагнитных систем магниторезонансных томографов. / Уч. пособие. СПб: ИТМО, 1998.

8. Сизиков В.С. Математические методы обработки результатов измерений. СПб: Политехника, 2001 (Электронный вариант: Сизиков В.С. Устойчивые методы обработки результатов измерений. Эл. учебник. http://de.ifmo.ru/library/sizikov.pdf) .

9. Adamiak K. Method of the magnetic field synthesis on the axis of cylinder solenoid // Appl. Phys. 1978. Vol. 16. P. 417-423.

10. Афанасьев Ю.В., Студенцов Н.В., Хорев В.Н. и др. Средства измерений параметров магнитного поля. Л.: Энергия, 1979.

11. Сизиков В.С., Марусина М.Я., Иванов С.В., Колобухова Т.Б., Николаев Д.Б., Соколов Д.Ю., Хомутникова Е.В. Прямая и обратная задачи синтеза магнитного поля в ЯМР-томографе. // Научно-техн. вестник СПбГИТМО(ТУ). 2001. Вып. 3. С. 209-214.

12. Сизиков В.С., Ахмадулин Р.И., Николаев Д.Б. Синтез магнитного поля вдоль оси катушки ЯМР-томографа. // Изв. вузов. Приборостроение. 2002. Т. 45, № 1. С. 5257.

13. Фриш С.Э., Тиморева А.В. Курс общей физики. Т. 2. МЛ: ГИТТЛ, 1952.

14. Дружкин Л.А. Задачи теории поля. М.: МИРГЭ, 1964.

15. 15.Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. М.: Наука, 1990.