ПРИМЕНЕНИЕ СПОСОБА ЭТАЛОННЫХ ПРИМЕРОВ ПРИ РЕШЕНИИ ОБРАТНОЙ ЗАДАЧИ СПЕКТРОСКОПИИ МЕТОДОМ РЕГУЛЯРИЗАЦИИ
44 В. С. Сизиков, А. В. Кривых
УДК 535.33+519.642.3+517.983.54+519.682
В. С. СИЗИКОВ, А. В. КРИВЫХ
ПРИМЕНЕНИЕ СПОСОБА ЭТАЛОННЫХ ПРИМЕРОВ ПРИ РЕШЕНИИ ОБРАТНОЙ ЗАДАЧИ СПЕКТРОСКОПИИ
МЕТОДОМ РЕГУЛЯРИЗАЦИИ
Рассматривается обратная задача спектроскопии — восстановление истинного спектра по измеренному спектру и аппаратной функции спектрометра путем решения интегрального уравнения Фредгольма I рода. Поскольку решается некорректная задача, то для получения устойчивого численного решения используется метод регуляризации Тихонова, при этом параметр регуляризации выбирается на основе решения эталонных примеров. Эффективное решение данной задачи позволит повысить разрешающую способность спектрометра. Приведены численные примеры.
Ключевые слова: обратная задача спектроскопии, интегральное уравнение, метод регуляризации Тихонова, способ эталонных примеров, разрешающая способность спектрометра, MatLab.
Введение. Для качественного и количественного исследования веществ широко используется спектральный анализ [1—6], основанный на изучении спектров излучения, поглощения, отражения, комбинационного рассеяния света и люминесценции. Под спектром u(ν) будем подразумевать зависимость интенсивности излучения u от частоты ν. Областями применения спектрального анализа являются физика, астрофизика, томография, механика, металлургия, химия и т.д.
Различают следующие виды спектров [3, 5, 6]: — непрерывные или сплошные — спектры, ширина спектральных линий которых заметно превышает естественную ширину [4, 6]; спектры с медленно меняющейся интенсивностью [1]; спектры, интенсивность которых, в терминах математики, изменяется непрерывно (однако могут присутствовать и „пустóты“ в спектре); — дискретные — спектры, состоящие из отдельных почти монохроматических спектральных линий [3]; — полосатые — спектры, состоящие из ряда полос, каждая из которых состоит из набора близких дискретных линий [3—6], и др. Для разложения излучения в спектр и его регистрации используются спектральные приборы (см. обзор В. А. Никитина: [5, с. 703—707]). Спектральный анализ можно условно разделить на два типа: — широкополосный анализ — изучение спектра в широкой полосе частот ν ∈[νmin , νmax ] , например изучение спектра звезды во всем видимом диапазоне; — узкополосный анализ — изучение спектра в узкой полосе частот, например изучение сверхтонкой структуры отдельной линии [1, 3]. Постановка задачи. Измеренный спектрометром (например, интерферометром Фабри — Перо [5, 6]) спектр u(ν) обычно отличается от истинного спектра z(ν) . Это проявляется, во-первых, в большей сглаженности спектра u(ν) по сравнению с z(ν) (не разрешены близкие линии, сглажена тонкая структура спектральной линии, что является результатом воздействия аппаратной функции спектрального прибора [2—4, 6]) и, во-вторых, в зашумленности спектра u(ν) (слабые линии „тонут“ в шуме, что является результатом случайных погрешностей измерений [1—3]).
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 45
Приведем определение аппаратной функции (АФ) [3]: аппаратной функцией K (ν, ν′)
спектрометра называется его реакция (в виде измеренной интенсивности) на дискретную линию единичной интенсивности и частоты ν′ при настройке спектрометра на частоту ν. Другие названия АФ [1, 3, 4, 6]: спектральная чувствительность, инструментальный контур, функция щели и т.д.
Аппаратная функция — это, вообще говоря, функция двух переменных, если ее форма (ширина и т.д.) изменяется с изменением частоты настройки (регистрации) ν ∈[νmin , νmax ] . Обычно с увеличением частоты ν АФ становится ýже, что характерно для широкополосной спектрометрии. Если же АФ практически не изменяется при изменении ν, то АФ является разностной (частотно-инвариантной) функцией: K (ν, ν′) = K (ν − ν′) , что имеет место, напри-
мер, при изучении тонкой структуры отдельной линии. На рис. 1 приведен непрерывный широкополосный измеренный спектр u(ν) , сглаженный
аппаратной функцией спектрометра K (ν, ν′) [3, 5], зашумленный измеренный спектр u(ν) и АФ
спектрометра, причем поскольку K (ν, ν′) — функция двух переменных, то приведено одно из ее
„сечений“ K (2, ν′) ; на рис. 1 обозначения ν, ν′ , u, u и K даны в условных единицах.
u, u , K, у.е.
6
5 4
3
2
1 0
1234567
8 9 10 ν, ν′, у.е.
Рис. 1
Следует отметить, что подобный вид имеет непрерывный узкополосный спектр [3], напри-
мер сверхтонкая структура отдельной линии, обусловленная магнитными или электрическими
полями (эффект Зеемана или Штарка), а также тепловым уширением (эффект Доплера) [5].
Как видно из рис. 1, в измеренном спектре u(ν) (тем более, в зашумленном спектре
u(ν) ) могут быть не разрешены близкие линии и не выявлены слабые, причем этот эффект
тем сильнее, чем шире АФ K (ν, ν′) (а также чем выше уровень шумов), другими словами, чем
меньше разрешающая способность спектрометра [1, 5]. Математическая формулировка обратной задачи спектроскопии. Рассмотрим ис-
тинный непрерывный спектр z(ν′) , ν′∈[a, b] , характерный для веществ с повышенной плот-
ностью (жидкий металл).
Измеренная интенсивность u(ν) при настройке спектрометра на частоту ν равна инте-
гралу (сумме) по всем интенсивностям z(ν′) с весовой функцией K :
b
u(ν) = ∫ z(ν′) K (ν, ν′) dν′ ,
(1)
a
откуда, варьируя значение ν (т.е. выполняя сканирование по спектру) и учитывая зашумлен-
ность спектра u(ν) , получаем
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
46 В. С. Сизиков, А. В. Кривых
b
∫ K (ν, ν′) z(ν′) dν′ = u(ν), c ≤ ν ≤ d ,
(2)
a
где [c, d] — пределы изменения ν (обычно более широкие, чем [a,b]).
В соотношении (2) известны (измерены или заданы) u(ν) , K (ν, ν′) , a , b , c , d , а z(ν′)
является искомой функцией.
Соотношение (2) есть интегральное уравнение Фредгольма I рода относительно z(ν′) .
Если K (ν, ν′) = K (ν − ν′) , то
∞
∫ K (ν − ν′) z(ν′) dν′ = u(ν), 0 ≤ ν < ∞ .
(3)
0
Соотношение (3) есть интегральное уравнение Фредгольма I рода типа свертки на полу-
оси. Задача решения уравнений (2) и (3) является некорректной [2, 3, 4, 7]: если решать урав-
нение (2), например, методом квадратур, а уравнение (3) — методом преобразования Фурье,
то в качестве решения получим так называемую „пилу“ [3] — крайне неустойчивое решение.
Поэтому для устойчивого решения этих уравнений необходимо применение устойчивых ме-
тодов, например метода регуляризации Тихонова [3, 4, 7].
При решении задачи обработки спектра в широкой полосе частот следует учитывать изме-
нение формы АФ K (ν, ν′) с изменением частоты ν . При решении задачи обработки спектра в
узкой полосе частот следует использовать уравнение Фредгольма I рода с разностным ядром:
b
∫ K (ν − ν′) z(ν′) dν′ = u(ν), c ≤ ν ≤ d .
(4)
a
Задача решения уравнений (2)—(4) связана с задачей редукции к идеальному спек-
тральному прибору [1—5, 8] — одним из вариантов редукционной проблемы Рэлея [9]. Ус-
пешное решение задачи редукции позволит путем численной обработки результатов повы-
сить разрешающую способность спектрального прибора.
Метод регуляризации Тихонова. Применительно к интегральному уравнению Фред-
гольма I рода (типа соотношения (2)) [3, 7]
b
A y ≡ ∫ K (x, s) y(s) ds = f (x), c ≤ x ≤ d ,
(5)
a
метод регуляризации Тихонова сводится к решению интегрального уравнения
b
α yα (t) + ∫ R(t, s) yα (s) ds = F (t), a ≤ t ≤ b ,
(6)
a
где α > 0 — параметр регуляризации, а новое ядро и новая правая часть равны
dd
R(t, s) = R(s,t) = ∫ K (x,t) K (x, s) dx, F (t) = ∫ K (x,t) f (x) dx .
(7)
cc
Способ эталонных примеров для выбора параметра регуляризации. Существует ряд
способов выбора параметра регуляризации α: способ невязки, обобщенный принцип невязки,
метод перекрестной значимости, локальный регуляризирующий алгоритм, способ подбора и
др. [3, 4, 7, 10—13].
В данной статье предлагается использовать способ эталонных примеров для выбора парамет-
ра α [7, 14, 15] с учетом дополнительной (априорной) информации о решении (искомом спектре).
Рассмотрим операторное уравнение I рода: A y = f (ср. с уравнением (5)). Пусть вместо
точных значений f и A известны f и A , такие что || f − f ||≤ δ , || A − A ||≤ ξ . При использова-
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 47
нии метода регуляризации Тихонова решается уравнение α yα + A*A yα = A* f (ср. с уравнением (6)). Пусть ∆ yα = yα − y — погрешность регуляризованного решения yα , а y — точное решение (нормальное псевдорешение [7, 12, 13]). В работах [7, 14, 15] получена следующая
оценка относительной погрешности регуляризованного решения по норме:
||
∆ ||
yα y ||
||
≤
ε(α)
,
(8)
где
ε(α)
=
|| A 2
|| η α
+
pα p α +1
.
(9)
Здесь η = δотн + ξотн , причем δотн = δ || f || и ξотн = ξ || A || — относительные погрешности
исходных данных; p =|| A+ ||2 , A+ — псевдообратный оператор: A+ f = y [12].
Функция ε(α) является огибающей для истинной относительной погрешности
σ(α)
=
||
∆ ||
yα y ||
||
.
(10)
В работах [7, 14, 15] показано, что функция ε(α) имеет (единственный) минимум при
( )p || A || η 2 < 27 16 ≈ 1, 69 . Согласно соотношениям (8), (9) оценка относительной погрешно-
сти || ∆ yα || || y || регуляризованного решения yα зависит от A и η. Поэтому если решается
несколько задач (примеров) с одинаковыми A и η, то для них оценки погрешности
|| ∆ yα || || y || будут одинаковыми. Отсюда следует, что при решении некоторого исходного примера P с неизвестным (ис-
комым) решением yP можно использовать результаты решения другого, модельного (эта-
лонного), примера Q с известным (заданным) точным решением yQ , но с такими же A и η,
что и в примере P. При решении примера Q можно рассчитать функцию
σ(α)Q = || ∆ yαQ ||
|| yQ ||
и по ней найти
αопт Q
—
значение α, при
котором
σ(α) = min .
α
Это
значение αопт Q может быть использовано при решении исходного примера P.
Отметим, что для повышения эффективности предлагаемого способа модельные примеры Q должны содержать дополнительную, априорную, информацию об исходном примере P.
Применительно к спектроскопии следует на основе анализа функции fP произвести предва-
рительную оценку количества спектральных линий (максимумов) в искомом спектре yP , соотношений их интенсивностей и значений их частот. Использование такой информации в
модельных примерах Q позволит более точно оценить параметр α. Примеры. В пакете прикладных программ MatLab были решены два численных примера.
Оба примера описываются интегральным уравнением типа (2), которое решается методом квадратур с регуляризацией Тихонова [3, с. 195] с помощью разработанной m-функции Tikh.m.
Первый пример — это исходный пример P (оригинал), для которого известны лишь зашумленный измеренный спектр u(ν) на сетке узлов ν = c, c + h,…, d , причем c = 0,5 ,
d = 10, 5 , h = const = 0, 05 (в условных единицах), и аппаратная функция
K (ν, ν′) =
q
π
⋅
exp
⎡ ⎢− ⎣⎢
1
q(ν − ν′)2 − 0, 2(ν 10)2
⎤ ⎥ ⎦⎥
,
(11)
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
48 В. С. Сизиков, А. В. Кривых
где q = 15 (ширина K (ν, ν′) убывает с ростом частоты ν). При этом истинный спектр z(ν)
неизвестен. Напомним (см. рис. 1), что измеренный спектр u(ν) имеет сложную структуру, а
именно: содержит шесть флуктуаций, две из которых (при ν ≈ 2 и ν ≈ 6,8) скорее всего состоят каждая из двух линий, но они не разрешены из-за того, что АФ имеет некоторую ширину и, тем самым, ограничивает разрешающую способность спектрометра. Из этого следует, что на самом деле в спектре имеется восемь или более спектральных линий. Поэтому в качестве второго (модельного) примера Q был рассмотрен близкий к оригиналу P пример, истинный спектр z(ν) которого состоит из 10 спектральных линий, имеющих вид гауссианы:
z(ν) = 3, 0 exp{−[(ν − 1, 7) / 0, 25]2} + 4,5exp{−[(ν − 2,3) / 0,35]2} + +2,5exp{−[(ν − 3, 4) / 0, 20]2} + 0, 4 exp{−[(ν − 3,9) / 0,10]2} + +8, 2 exp{−[(ν − 4,7) / 0,30]2} + 0,5exp{−[(ν − 5,9) / 0, 25]2} + +9, 2 exp{−[(ν − 6,5) / 0,10]2} + 8,5exp{−[(ν − 6,9) / 0,10]2} + +0, 4 exp{−[(ν − 7, 4) / 0, 25]2} + 2, 0 exp{−[(ν − 8, 7) / 0,50]2}.
(12)
Измеренный спектр u(ν) для примера Q был рассчитан численно по формуле
b
u(ν) = ∫ K (ν, ν′) z(ν′) dν′, c ≤ ν ≤ d ,
(13)
a
при a = 1, b = 10 (в условных единицах). К значениям u были добавлены случайные нормаль-
ные погрешности с СКО = 0,02 (что соответствует δотн ≈ 1 % ) и СКО = 0,002 ( δотн ≈ 0,1 % ). АФ спектрометра для примера Q была взята в виде (11), причем (чтобы промоделировать не-
точность знания АФ) параметр q был принят равным 15 (что соответствует ξотн = 0 ), а также q = 15 ⋅ (1− 0, 001) при ξотн ≈ 0,1 % и q = 15 ⋅ (1− 0, 03) при ξотн ≈ 3 % .
На рис. 2 для примера Q представлены истинный (точный) спектр z(ν) — кривая 1, из-
меренный (экспериментальный) спектр u(ν) — кривая 2, зашумленный спектр u(ν) — отме-
чен точками 2, и АФ K (2, ν′) — кривая 3.
u, u , z, K, у.е.
91 8 71 62
51 4
32
2 13
2
0
1 2 3 4 5 6 7 8 9 10 ν, ν′, у.е. Рис. 2
Далее, модельный пример Q был решен методом квадратур с регуляризацией Тихонова
с помощью m-функции Tikh.m для ряда значений параметра регуляризации α и построена за-
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 49
висимость относительной погрешности регуляризованного решения zα (ν) по отношению к точному решению z(ν) (см. формулу (10)):
σ(α)
=
||
zα (ν) − z(ν) || z(ν) ||
||
.
(14)
На рис. 3 представлены зависимости σ(α) для трех вариантов погрешностей, отобра-
женных соответственно следующими кривыми: 1 — при δотн ≈ 1 % , ξотн = 0 ; 2 — при δотн ≈ ξотн ≈ 0,1 % ; 3 — при δотн ≈ 1 % , ξотн ≈ 3 % , поскольку в примере-оригинале P обычно погрешности известны неточно. На рис. 3 представлена также огибающая ε(α) , при построе-
нии которой было положено: η = 10−2 , || A || = || u ||L2 || z ||L2 = 0,823 , а значение p было подобрано таким образом, чтобы кривая ε(α) касалась кривых 1—3.
σ, у.е.
0,9 0,8 0,7 0,6 0,5 0,4 0,3 0,2 0,1 0
–9
2 –8 –7
ε(α)
3 1
–6 –5 –4 –3 –2 –1
0
1 lg α
Рис. 3
Как видно из рис. 3, оптимальное значение параметра регуляризации α, например, для варианта 1 по минимуму кривой σ(α) равно приблизительно 10−3,7 . Однако при та-
ком значении α в спектре zα (ν) хорошо восстанавливаются сильные линии, а слабые линии „отягощены“ ложными флуктуациями и лучше восстанавливаются при большем значении α.
На рис. 4 представлены истинный спектр z(ν) — кривая 1, измеренный спектр
u(ν) — кривая 2 и восстановленный спектр zα (ν) — отмечен точками 3, при δотн ≈ 1 % ,
ξотн = 0 и α = 10−3 ; кривой 4 отмечена АФ K (2, ν′) . Кроме того, огибающая ε(α) (см. рис. 3)
имеет минимум также при 10−3 . Анализ графиков показывает, что спектр в модельном примере Q восстановился вполне удовлетворительно: разрешены близкие линии и восстановлены слабые линии, лишь на краю области (при ν > 9,5 ) наблюдаются слабые лож-
ные флуктуации в спектре zα (ν) — краевой эффект Гиббса.
При значении α = 10−3 , выбранном для модельного примера Q, восстановлен спектр также в примере-оригинале P (см. рис. 5, где обозначения соответствуют принятым для рис. 4). Как и в примере Q, в примере P разрешены близкие линии и восстановлены слабые линии, но также проявился краевой эффект Гиббса (при ν > 9,5 ), причем не только на краю, но и в
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
50 В. С. Сизиков, А. В. Кривых
середине спектра: при ν ≈ 5,8 слабая линия показывает ложное расщепление на две линии, однако это слабые эффекты (на уровне погрешностей метода).
z, zα, u, K, у.е.
9
8
7
6 53 41 32
2 14
3
1 2
1 3
2
0
12 z, zα, u, K, у.е.
34
5 6 7 8 9 ν, ν′, у.е.
Рис. 4
9 8 7 6 53 4 21 3 24 1
0
3 21
1 3
2
1 2 3 4 5 6 7 8 9 ν, ν′, у.е. Рис. 5
Практическое использование изложенной методики позволит повысить разрешающую способность спектрометра.
Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант № 09-08-00034.
СПИСОК ЛИТЕРАТУРЫ
1. Раутиан С. Г. Реальные спектральные приборы // Успехи физических наук. 1958. Т. 66, вып. 3. С. 475—517.
2. Кочиков И. В., Курамшина Г. М., Пентин Ю. А., Ягола А. Г. Обратные задачи колебательной спектроскопии. М.: Изд-во МГУ, 1993. 204 с.
3. Сизиков В. С. Математические методы обработки результатов измерений. СПб: Политехника, 2001. 240 с.
4. Старков В. Н. Конструктивные методы вычислительной физики в задачах интерпретации. Киев: Наукова думка, 2002. 264 с.
5. Физический энциклопедический словарь / Гл. ред. А. М. Прохоров. М.: Сов. энциклопедия, 1984. 944 с.
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 51
6. Ландсберг Г. С. Оптика: Учеб. пособие. М.: Физматлит, 2006. 848 с.
7. Верлань А. Ф., Сизиков В. С. Интегральные уравнения: методы, алгоритмы, программы. Киев: Наукова думка, 1986. 544 с.
8. Краулиня Э. К., Лиепа С. Я., Пикалов В. В., Скудра А. Я. К проблеме исследования атомной сенсибилизированной флуоресценции по контурам спектральных линий // Некорректные обратные задачи атомной физики / Под ред. Н. Г. Преображенского. Новосибирск: ИТПМ, 1976. С. 61—72.
9. Глазов М. В., Болохова Т. А. Решение редукционной проблемы Рэлея с использованием различных модификаций метода регуляризации // Оптика и спектроскопия. 1989. Т. 67, вып. 3. С. 533—537.
10. Воскобойников Ю. Е., Преображенский Н. Г., Седельников А. И. Математическая обработка эксперимента в молекулярной газодинамике. Новосибирск: Наука, 1984. 240 с.
11. Воскобойников Ю. Е., Мухина И. Н. Локальный регуляризирующий алгоритм восстановления контрастных сигналов и изображений // Автометрия. 2000. № 3. С. 45—53.
12. Морозов В. А. Регулярные методы решения некорректно поставленных задач. М.: Наука, 1987. 240 с.
13. Тихонов А. Н., Гончарский А. В., Степанов В. В., Ягола А. Г. Численные методы решения некорректных задач. М.: Наука, 1990. 232 с.
14. Сизиков В. С. Обобщенный метод редукции измерений // Электронное моделирование. 1991. Т. 13, № 4. С. 7—14.
15. Верлань А. Ф., Сизиков В. С., Мосенцова Л. В. Метод вычислительных экспериментов для решения интегральных уравнений в обратной задаче спектроскопии // Там же. 2011. Т. 33, № 2. С. 3—12.
Сведения об авторах
Валерий Сергеевич Сизиков
— д-р техн. наук, профессор; Санкт-Петербургский государственный
университет информационных технологий, механики и оптики, ка-
федра измерительных технологий и компьютерной томографии;
E-mail: sizikov2000@mail.ru
Александр Владимирович Кривых — магистр; Санкт-Петербургский государственный университет инфор-
мационных технологий, механики и оптики, кафедра измерительных
технологий и компьютерной томографии;
E-mail: krivykh1987@mail.ru
Рекомендована кафедрой измерительных технологий и компьютерной томографии
Поступила в редакцию 24.01.11 г.
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
УДК 535.33+519.642.3+517.983.54+519.682
В. С. СИЗИКОВ, А. В. КРИВЫХ
ПРИМЕНЕНИЕ СПОСОБА ЭТАЛОННЫХ ПРИМЕРОВ ПРИ РЕШЕНИИ ОБРАТНОЙ ЗАДАЧИ СПЕКТРОСКОПИИ
МЕТОДОМ РЕГУЛЯРИЗАЦИИ
Рассматривается обратная задача спектроскопии — восстановление истинного спектра по измеренному спектру и аппаратной функции спектрометра путем решения интегрального уравнения Фредгольма I рода. Поскольку решается некорректная задача, то для получения устойчивого численного решения используется метод регуляризации Тихонова, при этом параметр регуляризации выбирается на основе решения эталонных примеров. Эффективное решение данной задачи позволит повысить разрешающую способность спектрометра. Приведены численные примеры.
Ключевые слова: обратная задача спектроскопии, интегральное уравнение, метод регуляризации Тихонова, способ эталонных примеров, разрешающая способность спектрометра, MatLab.
Введение. Для качественного и количественного исследования веществ широко используется спектральный анализ [1—6], основанный на изучении спектров излучения, поглощения, отражения, комбинационного рассеяния света и люминесценции. Под спектром u(ν) будем подразумевать зависимость интенсивности излучения u от частоты ν. Областями применения спектрального анализа являются физика, астрофизика, томография, механика, металлургия, химия и т.д.
Различают следующие виды спектров [3, 5, 6]: — непрерывные или сплошные — спектры, ширина спектральных линий которых заметно превышает естественную ширину [4, 6]; спектры с медленно меняющейся интенсивностью [1]; спектры, интенсивность которых, в терминах математики, изменяется непрерывно (однако могут присутствовать и „пустóты“ в спектре); — дискретные — спектры, состоящие из отдельных почти монохроматических спектральных линий [3]; — полосатые — спектры, состоящие из ряда полос, каждая из которых состоит из набора близких дискретных линий [3—6], и др. Для разложения излучения в спектр и его регистрации используются спектральные приборы (см. обзор В. А. Никитина: [5, с. 703—707]). Спектральный анализ можно условно разделить на два типа: — широкополосный анализ — изучение спектра в широкой полосе частот ν ∈[νmin , νmax ] , например изучение спектра звезды во всем видимом диапазоне; — узкополосный анализ — изучение спектра в узкой полосе частот, например изучение сверхтонкой структуры отдельной линии [1, 3]. Постановка задачи. Измеренный спектрометром (например, интерферометром Фабри — Перо [5, 6]) спектр u(ν) обычно отличается от истинного спектра z(ν) . Это проявляется, во-первых, в большей сглаженности спектра u(ν) по сравнению с z(ν) (не разрешены близкие линии, сглажена тонкая структура спектральной линии, что является результатом воздействия аппаратной функции спектрального прибора [2—4, 6]) и, во-вторых, в зашумленности спектра u(ν) (слабые линии „тонут“ в шуме, что является результатом случайных погрешностей измерений [1—3]).
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 45
Приведем определение аппаратной функции (АФ) [3]: аппаратной функцией K (ν, ν′)
спектрометра называется его реакция (в виде измеренной интенсивности) на дискретную линию единичной интенсивности и частоты ν′ при настройке спектрометра на частоту ν. Другие названия АФ [1, 3, 4, 6]: спектральная чувствительность, инструментальный контур, функция щели и т.д.
Аппаратная функция — это, вообще говоря, функция двух переменных, если ее форма (ширина и т.д.) изменяется с изменением частоты настройки (регистрации) ν ∈[νmin , νmax ] . Обычно с увеличением частоты ν АФ становится ýже, что характерно для широкополосной спектрометрии. Если же АФ практически не изменяется при изменении ν, то АФ является разностной (частотно-инвариантной) функцией: K (ν, ν′) = K (ν − ν′) , что имеет место, напри-
мер, при изучении тонкой структуры отдельной линии. На рис. 1 приведен непрерывный широкополосный измеренный спектр u(ν) , сглаженный
аппаратной функцией спектрометра K (ν, ν′) [3, 5], зашумленный измеренный спектр u(ν) и АФ
спектрометра, причем поскольку K (ν, ν′) — функция двух переменных, то приведено одно из ее
„сечений“ K (2, ν′) ; на рис. 1 обозначения ν, ν′ , u, u и K даны в условных единицах.
u, u , K, у.е.
6
5 4
3
2
1 0
1234567
8 9 10 ν, ν′, у.е.
Рис. 1
Следует отметить, что подобный вид имеет непрерывный узкополосный спектр [3], напри-
мер сверхтонкая структура отдельной линии, обусловленная магнитными или электрическими
полями (эффект Зеемана или Штарка), а также тепловым уширением (эффект Доплера) [5].
Как видно из рис. 1, в измеренном спектре u(ν) (тем более, в зашумленном спектре
u(ν) ) могут быть не разрешены близкие линии и не выявлены слабые, причем этот эффект
тем сильнее, чем шире АФ K (ν, ν′) (а также чем выше уровень шумов), другими словами, чем
меньше разрешающая способность спектрометра [1, 5]. Математическая формулировка обратной задачи спектроскопии. Рассмотрим ис-
тинный непрерывный спектр z(ν′) , ν′∈[a, b] , характерный для веществ с повышенной плот-
ностью (жидкий металл).
Измеренная интенсивность u(ν) при настройке спектрометра на частоту ν равна инте-
гралу (сумме) по всем интенсивностям z(ν′) с весовой функцией K :
b
u(ν) = ∫ z(ν′) K (ν, ν′) dν′ ,
(1)
a
откуда, варьируя значение ν (т.е. выполняя сканирование по спектру) и учитывая зашумлен-
ность спектра u(ν) , получаем
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
46 В. С. Сизиков, А. В. Кривых
b
∫ K (ν, ν′) z(ν′) dν′ = u(ν), c ≤ ν ≤ d ,
(2)
a
где [c, d] — пределы изменения ν (обычно более широкие, чем [a,b]).
В соотношении (2) известны (измерены или заданы) u(ν) , K (ν, ν′) , a , b , c , d , а z(ν′)
является искомой функцией.
Соотношение (2) есть интегральное уравнение Фредгольма I рода относительно z(ν′) .
Если K (ν, ν′) = K (ν − ν′) , то
∞
∫ K (ν − ν′) z(ν′) dν′ = u(ν), 0 ≤ ν < ∞ .
(3)
0
Соотношение (3) есть интегральное уравнение Фредгольма I рода типа свертки на полу-
оси. Задача решения уравнений (2) и (3) является некорректной [2, 3, 4, 7]: если решать урав-
нение (2), например, методом квадратур, а уравнение (3) — методом преобразования Фурье,
то в качестве решения получим так называемую „пилу“ [3] — крайне неустойчивое решение.
Поэтому для устойчивого решения этих уравнений необходимо применение устойчивых ме-
тодов, например метода регуляризации Тихонова [3, 4, 7].
При решении задачи обработки спектра в широкой полосе частот следует учитывать изме-
нение формы АФ K (ν, ν′) с изменением частоты ν . При решении задачи обработки спектра в
узкой полосе частот следует использовать уравнение Фредгольма I рода с разностным ядром:
b
∫ K (ν − ν′) z(ν′) dν′ = u(ν), c ≤ ν ≤ d .
(4)
a
Задача решения уравнений (2)—(4) связана с задачей редукции к идеальному спек-
тральному прибору [1—5, 8] — одним из вариантов редукционной проблемы Рэлея [9]. Ус-
пешное решение задачи редукции позволит путем численной обработки результатов повы-
сить разрешающую способность спектрального прибора.
Метод регуляризации Тихонова. Применительно к интегральному уравнению Фред-
гольма I рода (типа соотношения (2)) [3, 7]
b
A y ≡ ∫ K (x, s) y(s) ds = f (x), c ≤ x ≤ d ,
(5)
a
метод регуляризации Тихонова сводится к решению интегрального уравнения
b
α yα (t) + ∫ R(t, s) yα (s) ds = F (t), a ≤ t ≤ b ,
(6)
a
где α > 0 — параметр регуляризации, а новое ядро и новая правая часть равны
dd
R(t, s) = R(s,t) = ∫ K (x,t) K (x, s) dx, F (t) = ∫ K (x,t) f (x) dx .
(7)
cc
Способ эталонных примеров для выбора параметра регуляризации. Существует ряд
способов выбора параметра регуляризации α: способ невязки, обобщенный принцип невязки,
метод перекрестной значимости, локальный регуляризирующий алгоритм, способ подбора и
др. [3, 4, 7, 10—13].
В данной статье предлагается использовать способ эталонных примеров для выбора парамет-
ра α [7, 14, 15] с учетом дополнительной (априорной) информации о решении (искомом спектре).
Рассмотрим операторное уравнение I рода: A y = f (ср. с уравнением (5)). Пусть вместо
точных значений f и A известны f и A , такие что || f − f ||≤ δ , || A − A ||≤ ξ . При использова-
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 47
нии метода регуляризации Тихонова решается уравнение α yα + A*A yα = A* f (ср. с уравнением (6)). Пусть ∆ yα = yα − y — погрешность регуляризованного решения yα , а y — точное решение (нормальное псевдорешение [7, 12, 13]). В работах [7, 14, 15] получена следующая
оценка относительной погрешности регуляризованного решения по норме:
||
∆ ||
yα y ||
||
≤
ε(α)
,
(8)
где
ε(α)
=
|| A 2
|| η α
+
pα p α +1
.
(9)
Здесь η = δотн + ξотн , причем δотн = δ || f || и ξотн = ξ || A || — относительные погрешности
исходных данных; p =|| A+ ||2 , A+ — псевдообратный оператор: A+ f = y [12].
Функция ε(α) является огибающей для истинной относительной погрешности
σ(α)
=
||
∆ ||
yα y ||
||
.
(10)
В работах [7, 14, 15] показано, что функция ε(α) имеет (единственный) минимум при
( )p || A || η 2 < 27 16 ≈ 1, 69 . Согласно соотношениям (8), (9) оценка относительной погрешно-
сти || ∆ yα || || y || регуляризованного решения yα зависит от A и η. Поэтому если решается
несколько задач (примеров) с одинаковыми A и η, то для них оценки погрешности
|| ∆ yα || || y || будут одинаковыми. Отсюда следует, что при решении некоторого исходного примера P с неизвестным (ис-
комым) решением yP можно использовать результаты решения другого, модельного (эта-
лонного), примера Q с известным (заданным) точным решением yQ , но с такими же A и η,
что и в примере P. При решении примера Q можно рассчитать функцию
σ(α)Q = || ∆ yαQ ||
|| yQ ||
и по ней найти
αопт Q
—
значение α, при
котором
σ(α) = min .
α
Это
значение αопт Q может быть использовано при решении исходного примера P.
Отметим, что для повышения эффективности предлагаемого способа модельные примеры Q должны содержать дополнительную, априорную, информацию об исходном примере P.
Применительно к спектроскопии следует на основе анализа функции fP произвести предва-
рительную оценку количества спектральных линий (максимумов) в искомом спектре yP , соотношений их интенсивностей и значений их частот. Использование такой информации в
модельных примерах Q позволит более точно оценить параметр α. Примеры. В пакете прикладных программ MatLab были решены два численных примера.
Оба примера описываются интегральным уравнением типа (2), которое решается методом квадратур с регуляризацией Тихонова [3, с. 195] с помощью разработанной m-функции Tikh.m.
Первый пример — это исходный пример P (оригинал), для которого известны лишь зашумленный измеренный спектр u(ν) на сетке узлов ν = c, c + h,…, d , причем c = 0,5 ,
d = 10, 5 , h = const = 0, 05 (в условных единицах), и аппаратная функция
K (ν, ν′) =
q
π
⋅
exp
⎡ ⎢− ⎣⎢
1
q(ν − ν′)2 − 0, 2(ν 10)2
⎤ ⎥ ⎦⎥
,
(11)
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
48 В. С. Сизиков, А. В. Кривых
где q = 15 (ширина K (ν, ν′) убывает с ростом частоты ν). При этом истинный спектр z(ν)
неизвестен. Напомним (см. рис. 1), что измеренный спектр u(ν) имеет сложную структуру, а
именно: содержит шесть флуктуаций, две из которых (при ν ≈ 2 и ν ≈ 6,8) скорее всего состоят каждая из двух линий, но они не разрешены из-за того, что АФ имеет некоторую ширину и, тем самым, ограничивает разрешающую способность спектрометра. Из этого следует, что на самом деле в спектре имеется восемь или более спектральных линий. Поэтому в качестве второго (модельного) примера Q был рассмотрен близкий к оригиналу P пример, истинный спектр z(ν) которого состоит из 10 спектральных линий, имеющих вид гауссианы:
z(ν) = 3, 0 exp{−[(ν − 1, 7) / 0, 25]2} + 4,5exp{−[(ν − 2,3) / 0,35]2} + +2,5exp{−[(ν − 3, 4) / 0, 20]2} + 0, 4 exp{−[(ν − 3,9) / 0,10]2} + +8, 2 exp{−[(ν − 4,7) / 0,30]2} + 0,5exp{−[(ν − 5,9) / 0, 25]2} + +9, 2 exp{−[(ν − 6,5) / 0,10]2} + 8,5exp{−[(ν − 6,9) / 0,10]2} + +0, 4 exp{−[(ν − 7, 4) / 0, 25]2} + 2, 0 exp{−[(ν − 8, 7) / 0,50]2}.
(12)
Измеренный спектр u(ν) для примера Q был рассчитан численно по формуле
b
u(ν) = ∫ K (ν, ν′) z(ν′) dν′, c ≤ ν ≤ d ,
(13)
a
при a = 1, b = 10 (в условных единицах). К значениям u были добавлены случайные нормаль-
ные погрешности с СКО = 0,02 (что соответствует δотн ≈ 1 % ) и СКО = 0,002 ( δотн ≈ 0,1 % ). АФ спектрометра для примера Q была взята в виде (11), причем (чтобы промоделировать не-
точность знания АФ) параметр q был принят равным 15 (что соответствует ξотн = 0 ), а также q = 15 ⋅ (1− 0, 001) при ξотн ≈ 0,1 % и q = 15 ⋅ (1− 0, 03) при ξотн ≈ 3 % .
На рис. 2 для примера Q представлены истинный (точный) спектр z(ν) — кривая 1, из-
меренный (экспериментальный) спектр u(ν) — кривая 2, зашумленный спектр u(ν) — отме-
чен точками 2, и АФ K (2, ν′) — кривая 3.
u, u , z, K, у.е.
91 8 71 62
51 4
32
2 13
2
0
1 2 3 4 5 6 7 8 9 10 ν, ν′, у.е. Рис. 2
Далее, модельный пример Q был решен методом квадратур с регуляризацией Тихонова
с помощью m-функции Tikh.m для ряда значений параметра регуляризации α и построена за-
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 49
висимость относительной погрешности регуляризованного решения zα (ν) по отношению к точному решению z(ν) (см. формулу (10)):
σ(α)
=
||
zα (ν) − z(ν) || z(ν) ||
||
.
(14)
На рис. 3 представлены зависимости σ(α) для трех вариантов погрешностей, отобра-
женных соответственно следующими кривыми: 1 — при δотн ≈ 1 % , ξотн = 0 ; 2 — при δотн ≈ ξотн ≈ 0,1 % ; 3 — при δотн ≈ 1 % , ξотн ≈ 3 % , поскольку в примере-оригинале P обычно погрешности известны неточно. На рис. 3 представлена также огибающая ε(α) , при построе-
нии которой было положено: η = 10−2 , || A || = || u ||L2 || z ||L2 = 0,823 , а значение p было подобрано таким образом, чтобы кривая ε(α) касалась кривых 1—3.
σ, у.е.
0,9 0,8 0,7 0,6 0,5 0,4 0,3 0,2 0,1 0
–9
2 –8 –7
ε(α)
3 1
–6 –5 –4 –3 –2 –1
0
1 lg α
Рис. 3
Как видно из рис. 3, оптимальное значение параметра регуляризации α, например, для варианта 1 по минимуму кривой σ(α) равно приблизительно 10−3,7 . Однако при та-
ком значении α в спектре zα (ν) хорошо восстанавливаются сильные линии, а слабые линии „отягощены“ ложными флуктуациями и лучше восстанавливаются при большем значении α.
На рис. 4 представлены истинный спектр z(ν) — кривая 1, измеренный спектр
u(ν) — кривая 2 и восстановленный спектр zα (ν) — отмечен точками 3, при δотн ≈ 1 % ,
ξотн = 0 и α = 10−3 ; кривой 4 отмечена АФ K (2, ν′) . Кроме того, огибающая ε(α) (см. рис. 3)
имеет минимум также при 10−3 . Анализ графиков показывает, что спектр в модельном примере Q восстановился вполне удовлетворительно: разрешены близкие линии и восстановлены слабые линии, лишь на краю области (при ν > 9,5 ) наблюдаются слабые лож-
ные флуктуации в спектре zα (ν) — краевой эффект Гиббса.
При значении α = 10−3 , выбранном для модельного примера Q, восстановлен спектр также в примере-оригинале P (см. рис. 5, где обозначения соответствуют принятым для рис. 4). Как и в примере Q, в примере P разрешены близкие линии и восстановлены слабые линии, но также проявился краевой эффект Гиббса (при ν > 9,5 ), причем не только на краю, но и в
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
50 В. С. Сизиков, А. В. Кривых
середине спектра: при ν ≈ 5,8 слабая линия показывает ложное расщепление на две линии, однако это слабые эффекты (на уровне погрешностей метода).
z, zα, u, K, у.е.
9
8
7
6 53 41 32
2 14
3
1 2
1 3
2
0
12 z, zα, u, K, у.е.
34
5 6 7 8 9 ν, ν′, у.е.
Рис. 4
9 8 7 6 53 4 21 3 24 1
0
3 21
1 3
2
1 2 3 4 5 6 7 8 9 ν, ν′, у.е. Рис. 5
Практическое использование изложенной методики позволит повысить разрешающую способность спектрометра.
Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант № 09-08-00034.
СПИСОК ЛИТЕРАТУРЫ
1. Раутиан С. Г. Реальные спектральные приборы // Успехи физических наук. 1958. Т. 66, вып. 3. С. 475—517.
2. Кочиков И. В., Курамшина Г. М., Пентин Ю. А., Ягола А. Г. Обратные задачи колебательной спектроскопии. М.: Изд-во МГУ, 1993. 204 с.
3. Сизиков В. С. Математические методы обработки результатов измерений. СПб: Политехника, 2001. 240 с.
4. Старков В. Н. Конструктивные методы вычислительной физики в задачах интерпретации. Киев: Наукова думка, 2002. 264 с.
5. Физический энциклопедический словарь / Гл. ред. А. М. Прохоров. М.: Сов. энциклопедия, 1984. 944 с.
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9
Применение способа эталонных примеров при решении обратной задачи спектроскопии 51
6. Ландсберг Г. С. Оптика: Учеб. пособие. М.: Физматлит, 2006. 848 с.
7. Верлань А. Ф., Сизиков В. С. Интегральные уравнения: методы, алгоритмы, программы. Киев: Наукова думка, 1986. 544 с.
8. Краулиня Э. К., Лиепа С. Я., Пикалов В. В., Скудра А. Я. К проблеме исследования атомной сенсибилизированной флуоресценции по контурам спектральных линий // Некорректные обратные задачи атомной физики / Под ред. Н. Г. Преображенского. Новосибирск: ИТПМ, 1976. С. 61—72.
9. Глазов М. В., Болохова Т. А. Решение редукционной проблемы Рэлея с использованием различных модификаций метода регуляризации // Оптика и спектроскопия. 1989. Т. 67, вып. 3. С. 533—537.
10. Воскобойников Ю. Е., Преображенский Н. Г., Седельников А. И. Математическая обработка эксперимента в молекулярной газодинамике. Новосибирск: Наука, 1984. 240 с.
11. Воскобойников Ю. Е., Мухина И. Н. Локальный регуляризирующий алгоритм восстановления контрастных сигналов и изображений // Автометрия. 2000. № 3. С. 45—53.
12. Морозов В. А. Регулярные методы решения некорректно поставленных задач. М.: Наука, 1987. 240 с.
13. Тихонов А. Н., Гончарский А. В., Степанов В. В., Ягола А. Г. Численные методы решения некорректных задач. М.: Наука, 1990. 232 с.
14. Сизиков В. С. Обобщенный метод редукции измерений // Электронное моделирование. 1991. Т. 13, № 4. С. 7—14.
15. Верлань А. Ф., Сизиков В. С., Мосенцова Л. В. Метод вычислительных экспериментов для решения интегральных уравнений в обратной задаче спектроскопии // Там же. 2011. Т. 33, № 2. С. 3—12.
Сведения об авторах
Валерий Сергеевич Сизиков
— д-р техн. наук, профессор; Санкт-Петербургский государственный
университет информационных технологий, механики и оптики, ка-
федра измерительных технологий и компьютерной томографии;
E-mail: sizikov2000@mail.ru
Александр Владимирович Кривых — магистр; Санкт-Петербургский государственный университет инфор-
мационных технологий, механики и оптики, кафедра измерительных
технологий и компьютерной томографии;
E-mail: krivykh1987@mail.ru
Рекомендована кафедрой измерительных технологий и компьютерной томографии
Поступила в редакцию 24.01.11 г.
ИЗВ. ВУЗОВ. ПРИБОРОСТРОЕНИЕ. 2011. Т. 54, № 9