CC BY
35
Радиофизика
Вестник Нижегородского университета им. Н.И. Лобачевского, 2014, № 4 (1), с. 66-72
УДК 621.391.037.372
ПРИМЕНЕНИЕ АЛГОРИТМА ЭКСТРАПОЛЯЦИИ ПРИ АНАЛИЗЕ НАЧАЛЬНЫХ УЧАСТКОВ РАДИОИНТЕРФЕРОГРАММ, ПОЛУЧЕННЫХ В ГАЗОДИНАМИЧЕСКИХ ЭКСПЕРИМЕНТАХ
© 2014 г. С.Ю. Лупов, О.А. Морозов
Нижегородский госуниверситет им. Н.И. Лобачевского [email protected]
Пкступила в редакцию 02.06.2014
Представлен алгоритм экстраполяции начального участка дискретного сигнала, основанный на оценке оптимальных параметров синусоидальной и комплексной экспоненциальной моделей. Приведены результаты вычисления спектрально-временного распределения энергии радиоинтерферограммы, полученной в опыте по метанию стальной пластины продуктами взрыва, с применением и без применения алгоритма экстраполяции.
Ключевые слква: спектрально-временное распределение, экстраполяция, газодинамический процесс, оптимизация.
Введение
К числу важнейших характеристик, измеряемых при исследовании быстропротекающих газодинамических процессов, относятся скорости и перемещения границ раздела сред, в которых развивается процесс. Радиоинтерферомет-рический метод измерений, основы которого представлены в [1], получил в настоящее время широкое распространение. Достоинствами метода являются непрерывность измерений и отсутствие возмущений, вносимых измерительными средствами в объект исследования [2-7].
Методы цифровой обработки отсчетов ин-терферограммы, основанные на понятии «мгновенная частота аналитического сигнала» [2; 4], позволяют вычислять мгновенную скорость движущейся границы. Однако, как показал опыт проведения измерений, простые алгоритмы обработки далеко не всегда дают удовлетворительные результаты. Основными причинами этого являются нарушение квадратуры канальных сигналов, возникающее при неидеальном согласовании антенны с зондируемой средой, действие шума приемника и присутствие в спектре входного сигнала компонент, возникающих из-за многократных отражений радиоволн. Надежные результаты измерений при анализе сигналов в сложных условиях регистрации могут быть получены при использовании различных методов вычисления спектрально-временного распределения энергии (СВРЭ), в основе которых лежат параметрические методы спектрального оценивания [8].
Газодинамические эксперименты, проводимые по оценке метательной способности взрыв-
чатых веществ, интересны тем, что скорость пластины, метаемой продуктами взрыва заряда из тротила, за несколько микросекунд достигает нескольких километров в секунду. Изменение скорости передней поверхности пластины происходит не плавно, а скачками [9]. Длительность этого участка на интерферограмме составляет менее 10 периодов характерного времени изменения сигнала, и именно он представляет наибольший интерес для исследователей. Первый скачок происходит в начале процесса, далее скорость поверхности пластины до следующего скачка меняется незначительно, и время между скачками соответствует на интер-ферограмме примерно половине периода.
Вычисление СВРЭ позволяет оценить скорость пластины, однако оказывается не всегда возможно оценить скорость поверхности пластины от начала процесса до следующего скачка. Начальный участок СВРЭ получается размытым вследствие того, что используемое в методах вычисления СВРЭ скользящее окно в этот момент практически полностью перекрывается с началом процесса. Увеличить точность оценки скорости пластины на начальном участке можно, если уменьшить амплитуду скачка частоты в начале процесса.
Учитывая, что частота интерферограммы после начального скачка меняется незначительно в течение времени, соответствующего примерно половине периода, предложено использовать этот участок для экстраполяции дополнительных начальных отсчетов, убрав тем самым начальный скачок частоты. Данный подход позволяет улучшить оценку СВРЭ этого участка.
1. Описание алгоритма экстраполяции
Для экстраполяции начального отрезка ин-терферограммы, содержащего N отсчетов действительного сигнала s(n), будем аппроксимировать его линейной комбинацией М синусоид
М
~(п) = 2 А С08(ю*п + Фк), п = 0,1,...,N -1, (1)
к=1
или (при использовании сигнала в виде отсчетов квадратурных компонент) комплексных экспонент
(т =1,2, ... ,М) и приравняем их к нулю. В результате получится система 2М линейных уравнений:
' М М N-1
2 Хк ' СС1к ¥к ' ™1к = 2 ^(п)' С08®1П
1V1
S (n)=2 Ak
j (ш4И+ф4)
, n = 0,1,...,N-1, (2)
(не обязательно ортогональных на данном отрезке) с различными частотами гак, амплитудами Ак и фазами фк, так, чтобы ошибка, вычисляемая по формуле:
2 Xk ' "M — 2Y • cSMk s(n)'
k=1 k=1 n=0
M M N-1
k ' SC1k — 2 Y ' '
k=1 k=1 n=0
(5)
2 Xk • sc1k— 2 Y • ss1k=2s(n)' sin®1n
2 Xk' scMk— 2 Y' ssMk=2s(n)'si
sinra Mn
где
s =
j\ —1
21s(n)—5(n)2,
(3)
cc
N -1
mk =2(C0S атП ' C0S ®kn) =
n=0
была минимальной.
Вычислив оптимальные параметры экспоненциальной модели, можно экстраполировать недостающие начальные отсчеты, подставив в = < модель (1) или (2) значения индекса времени n = -1, -2, ...
Следует заметить, что такая постановка задачи в точности соответствует статистической формулировке задачи нахождения максимально правдоподобной оценки неизвестных параметров модели сигнала в аддитивной смеси его с белым гауссовым шумом [10, 11].
Минимизировать ошибку аппроксимации (3), зависящую от 3M параметров rak, Ak и ф^ аналитически невозможно, поэтому для поиска оптимальных параметров модели необходимо ис- = пользовать алгоритмы оптимизации. При этом можно использовать тот факт, что при фиксированных частотах rak значения амплитуд Ak и фаз ф^: вычисляются аналитически, что сокращает количество поисковых параметров в 3 раза.
Вычисление значений амплитуд и фаз в модели действительного сигнала. Для действительного сигнала значения частот rak задаются в пределах от 0 до л.
Подставим выражение (1) в (3). Проделав несложные тригонометрические преобразования, получим
N-1 í M Л2
s=2ls(n)—2(Xkcosrakn-Yksin®kn) , (4) =
n=0 V k=1 )
где Xk = Ak cos и Yk = Ak sin .
Для нахождения оптимальных параметров Xk и Yk вычислим производные dsjdXm , ds/dYm
Л + sin((( - 0.5)'(cam-ra)) + 2 4' sin(05m-®k )2) + sin(((-0.5)'(ram » ra *ra
+ 4' sin (ram +»k )/2) ' m * N +1 — sin ra m - sin((2N - ) 2 4' sin ram
N> ra m =rak = 0
ra m =rak * 0
cs
N —1
mk =2(C0s Ü^mn ' sin rakn) =
n=0
cos1
(ram — rak))) — cos((N — 0.5)' (ram — rak)) + 4' sin(ram — Щ V2) чч
+ cos(ram +rak V^- cos(N — 0.5Hram +rak )/ ' 4' sin((ram +rak)2)
ra m
cos ra m — cos(2N — l)'ram)
4 ' sin(ra m )
ram =rak * 0
0 ra m =rak = 0
N—1
Scmk = 2 ^in Ü^mn ' cos ran = cSkm
n=0
SSmk =2(sin ramn ' sin rakn) =
sin ((( — 0.5)' ( — )) __
4' sin(ra m —rak )2)
sin((( — 0.5)'(ю„ )) ra *
--77-\T"\—, ram *(
4' sin((co„ +rat ))) m
N +1 sin ram — sin((2 N — 1)-rom )
4 ' sin(ram )
ram =rak * 0
0 ram =rak = 0,
k=1
k=1
n=0
cosraMn
k=1
k=1
k=1
n=0
n=0
n=0
2
из решения которой получаем Хк и Ук, такие, чтобы ошибка аппроксимации е была минимальной для заданных частот гак. Оптимальные значения амплитуд Ак и фаз фк будут получены из выражений
п=0 (
Ак
,1 (кп+Фк)
Л
"(О-Е Ак
-I(кп+Фк)
Л
Сделав замену \ = АкеЛк, получим
=Е
I=0 (
?(и)-Е к
к=1
М
* («)-Е К*е -
Ъ1Я21 + Ъ2Я22 + - + ЪмЯ2М = Е 5(и)е-^
Ъ1 ЯМ 1 + Ъ2ЯМ2 + - + КМ§ММ = Е ^ ,
где
Ятк Е \
-гат )п _
1 - е
1(к-®т
-т , га, Ф га
I (к-»т ) ' к т
1-е N, ®к =®т ,
-Е
Яе ^(и)-Е(Хк С08 гакп - Ук 81п гакп)
(
+
V
(8)
1т*(п)-Е(Хк §1пгакп + Ук 008га кп)
V к=1 У
где Хк = Яе Ък = Ак 008 фк и Ук = 1т Ък = Ак 81п фк.
Тогда система линейных уравнений (7) примет вид:
' М М
Е Хк • с1к+Еук ■ 51к =
Ак =а/ Хк2 + Ук2, Фк = Лго1Е Х-.
X.
Вычисление сначений амплитуд и фас в мк-дели ккмплекснкгк сигнала. Подставим выражение (2) в (3), получим
(6)
Для нахождения оптимальных параметров Ък вычислим производные ёе/ёЪт (т=1,2, ...,М) и приравняем их к нулю, получится система Р линейных уравнений:
Ъ1 Яп + Ъ2Я12 + ... + Кёш =Е ^У
п=0
N-1
к=1
N-1
= Е ((е ¿"(п) • 008 га1п + 1т ¿(п) • 81п га 1п)
ММ
Е Х к • СМк + Е Ук • 5Мк =
к=1 к=1
N-1
= Е((е ¿"(п)- 008 гаМп + 1т ¿-(п)- 81п гаМп)
п=0
ММ
ЕХк • ¿1к -ЕУк • С1к =
к=1 к=1
N-1
= Е((е ¿(п)- 81п га1п - 1т ¿(п)- 008 га1п)
(9)
где
(7)
Е Хк • *Мк -ЕУк
к=1 к=1
N-1
— Е ((е¿(п) • 81п гаМп - 1т ¿(п) • 008 гаМп),
п—0
N-1
Стк =Е 008 (га т -гак )-1) —
I=0
1 + 81п(((-0.5)-(гат-гак», га ^
2 2 - 81п(га т -гак)2) ' т '
N, гат =га.
тк
N-1
¿тк =Е 81п (05т -гак )' ) =
1=0
008(гат -гак )))-008(( - 0.5)(т -гак ^
2 - 81п((га т -гак )2)
0, га т =гак .
из решения которой получаем значения Ък, такие, чтобы ошибка аппроксимации е была минимальной для заданных частот гак.
Для численной реализации алгоритма выражение (6) удобнее представить в тригонометрическом виде. Формула (6) примет вид:
М У
Оптимальные значения амплитуд Ак и фаз фк будут получены из выражений
Ак = Ы = д/ Хк + Ук2,
1тй, У,
Фк = Лгс1§-- = Лгс1§
ЯеЪ,
Хк
Блкк-схема алгкритма. На рис. 1 представлена блок-схема работы алгоритма. С помощью генератора случайных чисел (блок 1) задаются начальные значения всех цифровых частот гак в диапазоне от 0 до п (для действительного сигнала) или от 0 до 2п (для комплексного сигна-
к=1
п=0
5
X
е
к=1
п=0
X
X
4'
к=1
п=0
п=0
га т Фгак
п=0
+
п=0
ла). В блоке 2 производится поиск оптимальных параметров юк, при которых ошибка аппроксимации, вычисляемая в целевой функции (блок 3), имеет минимальное значение. В блоке 4 реализована визуализация исследуемого и аппроксимированного сигналов (в виде графиков), найденных значений частот, амплитуд и фаз (в виде таблицы) и т.д. Аппроксимированный сигнал отображается вместе с экстраполированными отсчетами. В процессе поиска значения цифровых частот могут перейти за границы заданного диапазона, поэтому полученные данные при отображении приводятся в диапазон от 0 до п (для действительного сигнала) или от 0 до 2п (для комплексного сигнала).
Весь этот процесс происходит многократно, сохраняя значения поисковых параметров, при которых целевая функция имеет наименьшее значение. Заканчивать процесс следует, если значение погрешности е перестало уменьшаться.
Особенности реализации целевой функции. При вычислении целевой функции предвари-
тельно сравниваются друг с другом входные параметры. При совпадении двух или более параметров системы линейных уравнений (5) и (7) становятся вырожденными и требуют применения специальных методов решения. Чтобы избежать этого, предлагается находить параметры с близкими (с точностью 10-12) значениями, оставляя только один и сокращая количество уравнений в системе.
Дальнейшее изменение целевой функции связано с особенностями исследуемых сигналов. Как правило, в сигналах присутствует постоянная составляющая. Чтобы ее учесть, в целевой функции добавляется параметр юм+1=0, который не является поисковым. Это позволяет сократить количество поисковых параметров ещё на один.
Для поиска оптимальных параметров в блоке 2 (рис. 1) можно использовать подпрограмму E04CCF [12] библиотеки NAGLIB или виртуальный прибор LabVIEW Downhill Simplex nD, модифицированный для работы со сложными функциями [13]. В обоих реализован симплексный [14] метод оптимизации.
Выбкр ккличества частктных ккмпкнент в отсчетов короткого сигнала и основан на оценке
мкдели. Количество частотных компонент в модели выбирается исходя из априорной информации о сигнале. Для экстраполяции начальных отчетов интерферометрических данных, полученных в экспериментах по метанию пластины продуктами взрыва заряда из тротила, количество частотных компонент следует выбирать не более 3. При увеличении их амплитуда экстраполированных отсчётов сигнала может резко возрасти по сравнению с известными отсчётами. Это проявляется из-за неполного соответствия модели и сигнала.
2. Результаты применения алгоритма экстраполяции
Осциллограмма квадратурных составляющих сигнала, полученного в опыте по оценке метательной способности взрывчатых веществ, представлена на рис. 2. На рис. 3 представлено СВРЭ, вычисленное модифицированным преобразованием Вигнера-Вилля [8] с размером окна, равным 300 отсчетам, и количеством поисковых параметров, равным 5. По оси абсцисс отложено время в отсчетах (1 отсчет = = 0.016 мкс), по оси ординат отложены безразмерные цифровые частоты от 0 до 0.025, которые соответствуют скорости пластины от 0 до 10 км/с. На рисунке хорошо выделяется основная частотная компонента, соответствующая скорости движения передней поверхности пластины. Качественно можно видеть, что начальный участок распределения получился размытым.
На рис. 4 представлено СВРЭ интерферограм-мы, вычисленное с использованием дополнительных экстраполированных начальных отсчетов. На начальном участке отчетливо стала видна ступенька в районе цифровой частоты 0.009.
Выводы
Описанный в данной статье алгоритм предназначен для экстраполяции дополнительных
оптимальных параметров синусоидальной (2) и экспоненциальной (3) моделей сигнала, заданного на коротком отрезке.
Применение алгоритма дает возможность экстраполировать дополнительные начальные отсчеты интерферограммы, полученной в опыте по метанию стальной пластины продуктами взрыва заряда из тротила, что позволяет улучшить оценку спектрально-временного распределения энергии этого участка интерферограммы.
Спискк литературы
1. Горелик Г.С. Колебания и волны. 2-е изд. М.: Физматлит, 1950. 572 с.
2. Поршнев С.В. Радиолокационные методы измерений экспериментальной баллистики. Екатеринбург: УрО РАН, 1999. 211 0.
3. Баталов С.В., Филин В.П., Шапошников В.В. Радиоволновой метод исследования физических явлений и химических превращений в гетерогенных ВВ под действием УВ // ФГВ. 1991. Т. 27. № 6. С. 107-109.
4. Канаков В.А., Лупов С. Ю., Орехов Ю.И., Родионов А.В. Методы извлечения информации о перемещении границ раздела в газодинамических экспериментах с использованием радиоинтерферометров миллиметрового диапазона длин волн // Известия вузов. Радиофизика. 2008. Т. Ь1. № 3. С. 234-246.
5. Взятышев В.Ф., Михайлов А.Л., Орехов Ю.И., Родионов А.В. О возможности повышения качества радиоинтерферометрии при диагностике газодинамических процессов специально сформированными зондирующими волновыми образованиями // Экстремальные состояния вещества. Детонация. Ударные волны: Труды Междунар. конф. «IX Харитоновские тематические научные чтения». Саров, 2007. С. 643-647.
6. Богданов Е.Н., Бельский В.М., Родионов А.В. О влиянии ударосжатого слоя воздуха перед метаемой пластиной на измерения её скорости радиоин-терферометрическим методом // Экстремальные состояния вещества. Детонация. Ударные волны: Труды Междунар. конф. «IX Харитоновские тематические научные чтения». Саров, 2007. С. 680-685.
7. Родионов А.В., Канаков В.А., Лупов С.Ю. Методы обработки результатов радиоинтерферометри-
ческих измерений параметров газодинамических процессов // Экстремальные состояния вещества. Детонация. Ударные волны: Тр. Междунар. конф. «7 Харитоновские тематические научные чтения». Са-ров, 2005. С. 680-685.
8. Лупов С.Ю., Кривошеев В.И. Модификация преобразования Вигнера - Виля для анализа интер-ферометрических данных газодинамических процессов // Вестник Нижегородского университета им. Н.И. Лобачевского. 2011. № 5(3). С. 95-103.
9. Грязнов Е.Ф., Колпаков В.И., Уткин А.В. Особенности математического моделирования процессов метания цилиндрических оболочек и пластин, нагруженных взрывом // Труды VII Харитоновских чтений. Саров, 2005. C. 560-565.
10. Advance Topics in Signal Processing / Ed. by J. Lim, A. Oppenheim. Prentice Hall, Englewood Cliffs, NJ, 1998. 515 p.
11. Kay S.M. Fundamentals of Statistical Signal Processing: Estimation Theory (v.1). Prentice Hall, Upper Sadle River, NJ, 1998. 595 p.
12. NAG Fortran Library Routine Document. E04CCF/E04CC A. http ://www.nag. co .uk/numeric/fl/ manual/pdf/E04/e04ccf.pdf (дата обращения: 20.05.2014).
13. Лупов С.Ю., Фрадкина Е.П. Модификация функции LabVIEW «Downhill Simplex nD» и тестирование её на примере аппроксимации тестовых сигналов суммой синусоид // Труды Международной научно-практической конференции «Образовательные, научные и инженерные приложения в среде LabVIEW и технологии National Instruments». М.: Изд-во РУДН, 2006. C. 278-281.
14. Nelder J.A., Mead R. A Simplex Method for Function Minimization // Computer Journal. 1965. V. 7. P. 308-313.
APPLICATION OF AN EXTRAPOLATION ALGORITHM IN ANALYSIS OF RADIO INTERFEROGRAM INITIAL PARTS OBTAINED IN GASDYNAMIC EXPERIMENTS
S. Yu. Lupov, O.A.Morozov
An extrapolation algorithm of the discrete signal initial part is presented, based on an estimation of optimal parameters of sinusoidal and complex exponential models. Calculation results with and without application of the extrapolation algorithm are given for the spectral-time distribution of the radio interferogram energy obtained in the experiments on throwing a metal plate by TNT explosion products.
Keywords: spectral-time distribution, extrapolation, gasdynamic, optimization.
References
1. Gorelik G.S. Kolebaniya i volny. izd. M.: Fiz-matlit, 1950. 572 s.
2. Porshnev S.V. Radiolokacionnye metody izmeren-ij ehksperimental'noj ballistiki. Ekaterinburg: UrO RAN, 1999. 211 s.
3. Batalov S.V., Filin V.P., Shaposhnikov V.V. Ra-diovolnovoj metod issledovaniya fizicheskih yavlenij i himicheskih prevrashchenij v geterogennyh VV pod dej stviem UV // FGV. 1991. T. 27. № 6. S. 107-109.
4. Kanakov V.A., Lupov S.Yu., Orekhov Yu.I., Ro-dionov A.V. Metody izvlecheniya informacii o peremeshchenii granic razdela v gazodinamicheskih ehksperimentah s ispol'zovaniem radiointerferometrov millimetrovogo diapazona dlin voln // Izvestiya vuzov. Radiofizika. 2008. T. LI. № 3. S. 234-246.
5. Vzyatyshev V.F., Mihajlov A.L., Orekhov Yu.I., Rodionov A.V. O vozmozhnosti povysheniya kachestva radiointerferometrii pri diagnostike gazodinamicheskih processov special'no sformirovannymi zondiruyushchimi volnovymi obrazovaniyami // Ehkstremal'nye sostoyani-ya veshchestva. Detonaciya. Udarnye volny: Trudy Mezhdunar. konf. «IX Haritonovskie tematicheskie nauchnye chteniya». Sarov, 2007. S. 643-647.
6. Bogdanov E.N., Bel'skij V.M., Rodionov A.V. O vliyanii udaroszhatogo sloya vozduha pered metaemoj plastinoj na izmereniya eyo skorosti radiointerfero-metricheskim metodom // Ehkstremal'nye sostoyaniya veshchestva. Detonaciya. Udarnye volny: Trudy
Mezhdunar. konf. «IX Haritonovskie tematicheskie nauchnye chteniya». Sarov, 2007. S. 680-685.
7. Rodionov A.V., Kanakov V.A., Lupov S.Yu. Metody obrabotki rezul'tatov radiointerferometricheskih izmerenij parametrov gazodinamicheskih processov // Ehkstremal'nye sostoyaniya veshchestva. Detonaciya. Udarnye volny: Tr. Mezhdunar. konf. «7 Haritonovskie tematicheskie nauchnye chteniya». Sarov, 2005. S. 680-685.
8. Lupov S.Yu., Krivosheev V.I. Modifikaciya preobrazovaniya Vignera - Vilya dlya analiza interfero-metricheskih dannyh gazodinamicheskih processov // Vestnik Nizhegorodskogo universiteta im. N.I. Loba-chevskogo. 2011. № 5(3). S. 95-103.
9. Gryaznov E.F., Kolpakov V.I., Utkin A.V. Oso-bennosti matematicheskogo modelirovaniya processov metaniya cilindricheskih obolochek i plastin, nagruzhen-nyh vzryvom // Trudy VII Haritonovskih chtenij. Sarov, 2005. S. 560-565.
10. Advance Topics in Signal Processing / Ed. by J. Lim, A. Oppenheim. Prentice Hall, Englewood Cliffs, NJ, 1998. 515 p.
11. Kay S.M. Fundamentals of Statistical Signal Processing: Estimation Theory (v.1). Prentice Hall, Upper Sadle River, NJ, 1998. 595 p.
12. NAG Fortran Library Routine Document. E04CCF/E04CC A. http ://www. nag. co .uk/numeric/fl/ manual/pdf/E04/e04ccf.pdf (data obrashcheniya: 20.05.2014).
13. Lupov S.Yu., Fradkina E.P. Modifikaciya funkcii LabVIEW «Downhill Simplex nD» i testirovanie eyo na primere approksimacii testovyh signalov summoj sinusoid // Trudy Mezhdunarodnoj nauchno-prakticheskoj
konferencii «Obrazovatel'nye, nauchnye i inzhenernye 14. Nelder J.A., Mead R. A Simplex Method for prilozheniya v srede LabVIEW i tekhnologii National Function Minimization // Computer Journal. 1965. V. 7. Instruments». M.: Izd-vo RUDN, 2006. S. 278-281. P. 308-313.
