УДК 519
АВТОМАТИЗАЦИЯ ВЫЧИСЛЕНИЯ К-ИНДЕКСА С ИСПОЛЬЗОВАНИЕМ ПАКЕТНОГО ВЕЙВЛЕТ-ПРЕ ОБРАЗОВАНИЯ
© 2008 г. О.В. Мандрикова, С.Э. Смирнов
Институт космофизических исследований и распространения радиоволн ДВО РАН, 684034, Камчатская обл., Елизовский р-он, с. Паратунка, ул. Мирная, д. 7, 1кк@кгг.катсЬа1к:а. ги
Institute of Cosmophysical Research and Radio Wave Propagation FEB RAS, 684034, Kamchatskiy obl., v. Paratunka, Mirnaya St., 7, [email protected]
Предлагается автоматический алгоритм вычисления К-индекса с использованием пакетного вейвлет-преобразования. Обработка временных рядов производилась на основе алгоритма быстрого вейвлет-преобразования. Приводится описание способа определения К-индекса и результатов обработки магнитных данных Камчатской области.
Ключевые слова: вейвлет-преобразование, цифровая обработка сигналов, индекс геомагнитной активности.
In the paper an automatic algorithm of calculation ofK-index applying wavelet transform package is suggested. Time series processing was performed on the basis of the algorithm of fast wavelet-transform. The method of estimation ofK-index and the results ofprocessing of magnetic data for Kamchatka region are described.
Keywords: wavelet-transformation, digital processing of signals, index of geomagnetic activity.
Постановка задачи
Индекс геомагнитной активности К характеризует меру интенсивности геомагнитных возмущений в месте регистрации и используется при решении широкого ряда научных и прикладных задач. Индекс был введен Дж. Бартельсом в 1938 г. и представляет собой значения от 0 до 9 для каждого трехчасового интервала мирового времени. К-индекс геомагнитной активности вычисляется из магнитограмм по специальной методике, которая
предполагает высокий уровень профессиональной подготовки персонала магнитной обсерватории. Основные требования данной методики предъявляются к выделению невозмущенной вариации геомагнитного поля, которая называется Sq-кривой (Solar quiet-вариация), обусловленной токовой системой, возникающей на освещенной стороне Земли на высоте слоя E ионосферы. На основе разницы между наибольшим и наименьшим отклонениями в течение трехчасового интервала реальной
магнитограммы от значения Бд-кривой определяют амплитуду возмущения в нТл для данной обсерватории. Эта амплитуда переводится по квазилогарифмической шкале в К-индекс. Шкала для каждой обсерватории определяется в соответствии с геомагнитной широтой индивидуально.
Автоматизация процесса определения индекса геомагнитной активности позволит более эффективно решать задачи обработки и исследования геофизических данных, обеспечит автоматический анализ регистрируемых сигналов. В настоящее время эта задача актуальна, т.к. есть нерешенные методологические вопросы определения К-индекса, особенно для нулевых дней [1]. Создание автоматических методов определения К-индекса является серьезной проблемой, трудности решения которой связаны с реализацией процедуры определения Бд-кривой, максимально приближенной к ручному способу ее выделения. Существующие методы автоматизации данной процедуры не удовлетворяют основному требованию, поскольку не включают в себя средства адаптации, тем самым не учитывают возможность изменчивости Sq-кривой с течением времени. В работе [2] предложен способ определения «не-К-вариации» на основе метода естественных ортогональных составляющих, который позволяет идентифицировать данную кривую на любой из дней текущей десятидневки. Данный метод хотя и является удовлетворительным, но не обеспечивает требуемую точность определения К-индекса. Недостаток этого метода заключается в отсутствии во-первых, критерия для выборки дней десятидневки, во-вторых, адаптации алгоритма к периодам влияния токовой системы DP0, постоянно существующей в высоких широтах и генерируемой процессами квазивязкого обтекания магнитосферы солнечным ветром [2]. В случае попадания возмущенных дней в эту десятку Sq-вариация подвергается искажениям, что влияет на результат определения К-индекса, особенно в дневные часы.
В данной работе предлагается способ автоматического определения Бд-вариации на основе разложения составляющих напряженности магнитного поля в вейвлет-пакеты. В настоящее время вейвлет-преобразование широко применяется для решения многих задач в радиотехнике, электронике, сейсмо-акустике, биологии, ядерной физике и других областях науки и техники [3, 4]. В данном случае эта математическая конструкция - наиболее удобный инструмент для реализации процедуры определения Sq-вариации ввиду следующих ее свойств. Разложение сигнала по вейвлет-базису позволяет осуществить разнесение в пространстве разномасштабных его составляющих, что значительно облегчает процедуру выделения информативных компонент сигнала. Одни из выделенных компонент является регулярными составляющими процесса, другие содержат в себе информацию о величине отклонений от среднего его уровня в локальные моменты времени. На основе анализа и обработки этих компонент может быть реализована процедура идентификации спокойного (либо неспокойного) дня. В свою очередь выбор базисной функции (из большого числа
возможных) и определение наилучшего пути по вейвлет-дереву с учетом определенной функции стоимости обеспечивают получение алгоритма, адаптируемого к сигналу по мере его изменения. Преимуществом данного способа адаптации служит также тот факт, что при этом не требуется процедуры обучения и знания статистических свойств сигнала. Таким образом, может быть получен алгоритм, во-первых, позволяющий выполнить идентификацию возмущенной составляющей напряженности магнитного поля, во-вторых, реализующий определение Sq-вариации на конкретный момент времени.
Описание алгоритма
В процессе работы были рассмотрены узлы дерева вейвлет-пакета для составляющих напряженности магнитного поля. При проведении анализа использовались разные базисные функции класса Добеши [3]. Разложение осуществлялось до 6-го уровня полного вейвлет-дерева включительно. Путем обработки статистических данных было обнаружено, что для выполнения процедуры идентификации спокойного дня информативными являются ветки дерева вейвлет-пакета, показанные на рис. 1. Эти ветки обозначим AAAD, AADD, AAADD. В качестве базисной функции была выбрана функция Добеши 3-го порядка. На основе анализа спокойных дней для каждой из выделенных компонент были введены в рассмотрение так называемые пороговые значения (ПЗ), равные суммарным значениям вейвлет-коэф-фициентов данной ветки, не превышение которых служило характеристикой того, что компонента является невозмущенной.
Способ идентификации спокойного дня заключался в следующем:
1. На основе обработки ближайшего к текущей дате спокойного дня производился расчет суммы ПЗ.
2. Для составляющей магнитного поля текущего дня производился расчет суммы значений вейвлет-коэффи-циентов компонент разложения АААБ, ААББ, АААББ, обозначаемой далее Б, и выполнялось ее сравнение с суммой ПЗ (п. 1).
3. Если найденная сумма Б была меньше суммы ПЗ, то текущий день определялся как спокойный.
На основе описанной методики был предложен следующий автоматический алгоритм определения К-индекса:
Рис. 1
1. В память загружается сигнал составляющей магнитного поля текущего дня, и на основе описанной выше схемы осуществляется проверка, является ли он спокойным.
2. Если текущий день был определен как спокойный, то он замещает первый день в тройке спокойных дней, и с учетом изменений рассчитывается Sq-вариация.
3. Для расчета Sq-вариации используются значения трех спокойных дней, ближайших по дате к текущему дню (на рис. 2 обозначены соответственно X, У, 2): их соответствующие значения складываются и делятся на 3 (вектор 8 на рис. 2).
4. На основе полученной Sq-вариации далее определяется К-индекс по методике, соответствующей ручному способу (на основе разницы между наибольшим и наименьшим отклонениями в течение трехчасового интервала реальной магнитограммы от значения Sq-кривой определяется амплитуда К-возмущения, которая потом переводится по соответствующей геомагнитной широте шкале в К-индекс).
Результаты экспериментов
На основе предложенной методики были обработаны магнитные данные, полученные на обс. «Пара-тунка» (с. Паратунка, Камчатская область, ИКИР ДВО РАН). В таблице для сравнения приведены результаты определения К-индекса соответственно ручным и автоматическим методами за январь 2005 г.
Анализ результатов расчета данного месяца показывает, что отличия составляют малый процент (45 значений из 248) и в своем большинстве не превышают 1 балла. Данный результат позволяет считать предлагаемую методику вполне удовлетворительной по точности и допустимой для проведения автоматических расчетов по определению индекса геомагнитной активности К.
Результаты расчета К-индекса двумя различными методиками за январь 2005 г.
X
Y
Z
S
Число Ручной расчет Автоматический расчет
1 2 3 3 3 3 3 4 4 2 3 3 2 3 3 3 4
2 3 3 5 4 3 5 3 5 3 3 4 4 3 5 1 5
3 3 3 3 4 4 4 3 2 3 3 2 4 4 4 3 1
4 3 3 3 3 4 3 3 3 3 3 5 3 4 3 3 3
5 3 2 4 4 4 4 2 2 3 2 4 4 4 3 2 2
6 1 1 0 1 1 1 2 1 1 2 0 1 0 1 2 1
7 0 1 0 3 5 5 4 7 0 1 0 3 5 5 4 6
8 5 5 3 4 3 2 2 2 5 2 3 4 3 2 1 2
9 1 1 0 3 1 1 2 0 1 1 0 3 1 0 2 0
10 1 1 1 1 2 3 3 2 1 1 1 1 2 3 3 2
11 1 3 4 2 2 3 2 3 1 3 4 2 2 2 2 3
12 3 3 3 5 5 4 4 3 3 3 3 5 5 4 4 3
13 2 3 2 3 3 3 3 2 1 3 2 3 5 3 3 2
14 2 2 2 1 1 1 4 4 2 2 2 1 3 2 4 4
15 3 3 4 3 3 3 2 3 3 3 4 3 2 3 2 2
16 2 1 2 3 3 4 2 3 2 1 2 3 3 4 2 2
17 4 3 4 5 5 6 5 6 4 3 4 5 5 6 6 6
18 6 5 7 4 5 4 4 4 6 5 6 4 5 4 4 3
19 5 4 5 5 4 4 2 3 5 4 5 4 4 4 2 3
20 2 2 2 2 3 4 3 2 2 2 1 2 3 4 3 2
21 2 2 2 3 2 6 7 6 2 2 2 2 2 6 7 5
22 5 4 4 4 4 4 4 3 5 4 4 4 3 3 3 3
23 2 3 2 3 4 4 2 2 2 3 2 3 4 2 1 2
24 2 1 1 2 3 3 2 2 3 1 1 2 3 3 2 2
25 2 1 1 3 3 2 2 1 2 1 1 3 3 2 2 1
26 0 1 0 1 2 0 1 2 0 1 0 1 2 1 1 2
27 2 2 1 1 1 1 2 2 2 2 1 1 1 1 2 2
28 2 0 1 1 1 2 3 3 2 0 1 3 2 2 3 2
29 2 2 2 4 3 3 4 3 2 2 1 4 3 3 4 2
30 2 3 2 5 3 3 2 2 2 3 2 5 2 3 2 1
31 2 2 3 5 5 4 3 2 2 1 3 4 5 4 2 2
ч. s
K'IVIII Iii
MAU
tat Mw К11Ж май НИ Я*
1гшу7 L/^v-Im. fr\ гшгя 21МИ дрш
0 720 i« о 72t 1«
tM 14* 11 II MS nmtmwoglMi
|шйм i uW^ 2117241 21JWL Zj/t.OS
1 0 720 IM 0 72Ö 1«
Рис. 2
Из вышеизложенного следует, что на основе новой математической конструкции -вейвлет-преобразовании предложен способ автоматического определения К-индекса по регистрируемым данным магнитного поля Земли. В основе алгоритма лежит конструкция разложения сигнала в вейвлет-пакеты.
На примере обработки магнитных данных Камчатского региона выполнена апробация метода и получены положительные результаты.
Литература
1. Будько Н., Зайцев А., Карпачев А., Козлов А., Филиппов Б. Космическая среда вокруг нас. Троицк, 2006.
2. Головков В.П., Папиташвили В.О., Папи-ташвили Н.Е. Автоматизированное вычисление К-индексов с использованием метода естественных ортогональных составляющих // Геомагнетизм и аэрономия. 1989. Т. 29. № 4. С. 667-670.
3. Добеши И. Десять лекций по вейвлетам: Пер. с англ. Ижевск, 2001.
4. Малла С. Вейвлеты в обработке сигналов: Пер. с англ. М., 2005.
Поступила в редакцию
24 марта 2008 г.