Вестник КРАУНЦ. Физ.-мат. науки. 2019. Т. 29. № 4. C. 87-97. ISSN 2079-6641
DOI: 10.26117/2079-6641-2019-29-4-87-97
УДК 516
МЕТОД ОБНАРУЖЕНИЯ ПОМЕХ В ГЕОМАГНИТНЫХ ДАННЫХ
С. Ю. Папшева, О. В. Мандрикова, С. Ю. Хомутов
Институт космофизических исследований и распространения радиоволн ДВО РАН, 684034, Камчатский край, с. Паратунка, ул. Мирная, 7 E-mail: 75sp19@mail.ru, oksanam1@mail.ru, khomutov@ikir.ru
Предложен вычислительный метод обнаружения помех в геомагнитных данных, основанный на вейвлет-преобразовании и пороговых функциях. Эффективность метода показана на примере обработки результатов измерений с помощью вариационного феррозон-дового магнитометра FGE-DTU (обсерватория «Паратунка», Камчатский край, ИКИР ДВО РАН). Рассмотрены некоторые виды помех от естественных источников (при землетрясениях на Камчатке) и техногенных помех, связанных с работой ионозонда. Детально изучена частотно-временная структура помех в вариациях Z- и D- составляющих геомагнитного поля (частота измерений 2 Гц). Для рассмотренного вида помех определены информативные масштабные уровни вейвлет-преобразования и оценены параметры алгоритма реализации метода.
Ключевые слова: помехи в геомагнитных данных, вейвлет-преобразование, анализ данных
(с) Папшева С. Ю., Мандрикова О. В., Хомутов С. Ю., 2019
Введение
При регистрации вариаций геомагнитного поля возникают помехи различного характера, которые оказывают негативное влияние на качество получаемых данных и вносят существенные погрешности в результаты их обработки и научного анализа [1]. Задача обнаружения помех возникает во многих областях науки и техники и является одной из основных проблем в различных технических системах -радиотехнических, гидроакустических, навигационных спутниковых и др. Методы, направленные на решение данной задачи, используют разные подходы, основанные, как на применении классических методов анализа данных (например, факторный анализ, метод главных компонент [2]), так и современных математических аппаратов (нейронных сетей, вейвлет-преобразования, теории нечетких множеств) и их комбинаций (например, [3]- [5] и др.). В данной работе используется подход, основанный на применении методов вейвлет-преобразования. Вейвлет-преобразование
широко используется в задачах обработки данных в различных прикладных областях, в т.ч. составляет основу современных методов обнаружения сигналов на фоне помех ( [2]- [3], [6]- [9] и др.). В работе предложено использовать конструкцию непрерывного вейвлет-преобразования [11] совместно с пороговыми функциями. Данный подход был впервые предложен в [9] и показана его эффективность в задаче первичной обработки магнитных данных. Статья является продолжением исследований в этом направлении. В работе рассмотрены виды помех, связанных с работой ионозонда, и от естественных источников - землетрясений на Камчатке. Выполнен анализ частотно-временной структуры помех в вариациях 2- и Э-составляющих геомагнитного поля. В исследовании использовались геомагнитные данные (вариации составляющих Н, Э, 2) Геофизической обсерватории «Паратунка» (сертифицирована как обсерватория INTERMAGNET), полученные с помощью вариационного ферро-зондового магнитометра FGE-DTU (частота измерений 2 Гц).
Описание метода
Следуя работе [8], для отображения данных в вейвлет-пространство использовалось непрерывное вейвлет- преобразование
(Щ/)(Ь,а) = |а|-1/2 | /(г^ ^, (1)
—^
где ^ - вейвлет, / е Ь2(Я), а,Ь е Я'а = 0, а- масштаб, Ь- время. В качестве базисного вейвлета использовался ортогональный вейвлет Добеши порядка 3. За меру интенсивности помехи на масштабе а принята величина
еь1 = !№/)(Ь, а)|, (2)
где г = Ь/- момент возникновения помехи. Интенсивность помехи в момент времени г = Ь определялась как
Еь = £ еЬ,а, (3)
а
Операция обнаружения помехи основана на применении пороговой функции:
Рт(Еь) = ( 1 если ЕЬ > Т (4)
4 ' [0, если ЕЬ < Т,
где Т - порог. Значение Рт (Еь) = 1 свидетельствует о наличии помехи в момент времени г = Ь.
Выбор порога Т. Пороговая функция в (4) определяет правила выбора решения о наличии/отсутствии помехи в сигнале в момент времени г = Ь. В силу случайной природы возникновения помех использование любого правила связано с возможностью ошибочных решений. Определение порога в работе выполнялось путем оценки апостериорного риска.
Порог разбивает пространство значений X анализируемой функции на две непересекающиеся области Хо и Х1 . Тогда правило выбора решения устанавливает соответствие между решениями о наличии/отсутствии помехи и областями. При использовании определенного правила выбора решения для заданного состояния Н^
(характеризует наличие/отсутствие помехи) средняя величина потерь может быть определена как
1
I) (х) = £п ар{х е X/к)}, 1=0
где Пц - функция потерь, Р{х е Х[/к)} - условная вероятность попадания выборки в область Х1, если в действительности имеет место состояние к), ) = I, ), I- индексы состояний (знак "/"означает условную вероятность). Усредняя условную функцию риска по всем состояниям к), ) = 0,1, получаем средний риск [9]:
1 = £ р!),
1=0
где р)- априорная вероятность состояния к).
Наилучшим правилом будет такое, для которого средний риск будет наименьшим (байесовский риск, [9]).
Поскольку мы не знаем априорное распределение состояний Р), для выбора наилучшего правила будем использовать апостериорный риск [9]. Апостериорные вероятности Р{к)/х}, I = 0,1 представляют наиболее полную характеристику состояний к)
при располагаемых априорных данных. Для простой функции потерь П = | ^ , | = ^ апостериорный риск II (х) равен
I(х)= £Р{к)/х е XI}.
В этом случае критерием качества выбора решения является критерий наименьшей частоты ошибок. Порог определяется наилучшим правилом выбора решения, обеспечивающим наименьшее значение апостериорного риска I(х). Критерием оценки состояния к) являлось решение эксперта о наличии/отсутствии помехи. Состояния несут информацию о наличии/отсутствии помехи и их оценка (на основе принятого правила выбора решения) позволяет фиксировать моменты возникновения помех в сигнале.
Определение информативных масштабов а). Операция обнаружения помех (соотношение (4)) требует определения информативных масштабов а), вносящих основной вклад в интенсивность помехи Е^. Определение информативных масштабов а) выполнялось на основе оценки величин:
• средней интенсивности фона для Н-, Э-, 2- компонент геомагнитного поля
аеЬ,а /г\
Е/ап = -, (5)
7 п
где еьа = /)(Ь,а)|, п- количество временных отсчетов 9-часовой выборки фона;
• отношения сигнал/шум в моменты возникновения помех ? = I, I = 1,Ь, где Ь- количество помех
к = ЕЕ-; (6)
Е/ап
• максимальных Е[тах и минимальных Етп значений интенсивностей помех.
Алгоритм обнаружения помехи в момент времени г = г/ включал следующие шаги:
1) отображение данных в вейвлет-пространство (операция (1)) до масштаба а/, где а/ - наибольший информативный масштаб;
2) оценка интенсивности функции Е{, в момент времени г = г/ (операция (3));
3) применение пороговой функции Рт (операция (4)).
Результаты экспериментов
В работе изучены вариации геомагнитного поля за период с 01.01.2016 г. по 31.12.2016 г. Анализируемый период содержал 56 суточных интервалов с помехами, которые вызваны землетрясениями на Камчатке, из них 49 произошли в спокойные (индекс геомагнитной активности К не превышал значения 1) и слабовозмущенные (К-индекс имел значения 2 и 3) дни. Идентификация помех в магнитных записях и определение вероятного их источника выполнялись вручную магнитологами обсерватории «Паратунка» ИКИР ДВО РАН.
На рис. 1 в качестве примера показаны нормированные вариации Н-, Э-, 2- составляющих геомагнитного поля за 01.12.2016 г.
Рис. 1. Нормированные вариации составляющих Н, Э, 2 геомагнитного поля за 01.12.2016 г.
Нормировка данных выполнялась на основе операции
xi ximin
Xl =
Ximax ximin
(7)
где Xi - ¡-ое значение исходной выборки, ximin и Ximax- минимальное и максимальное значения исходной выборки.
Землетрясение, произошедшее в 05:16:24.8иТ (широта 52.2408Ы, долгота 158.1810Е, 96 км к северо-востоку от Петропавловска-Камчатского, глубина 107 км, Кб=Ю.8, МЬ=3.6, (отмечено на рис. 1 овалом), проявилось во всех компонентах - Н, Б и 2, но наиболее выражено в вариации склонения Б. Магнитометр РОЕ-БТи является чувствительным к землетрясениям, поскольку блок его феррозондовых датчиков закреплён на карданном подвесе и испытывает механические колебания при прохождении сейсмической волны через место установки прибора. Поэтому наблюдаемые в магнитной записи осцилляции размахом от долей до десятков нТл не связаны с вариациями магнитного поля и являются помехой.
Кроме того, в Б- и 2- составляющих (рис. 1) хорошо заметны импульсные помехи, возникающие при сеансах вертикального зондирования ионосферы ионозондом "Парус" каждые 15 мин. и достигающие 1 нТл.
На рис. 2, 3 и 4 показаны результаты расчета интенсивностей (операция 3) техногенных и естественных возмущений в Н-, Б- и 2- компонентах в 9-часовом интервале.
Результаты оценок (операции (5), (6)) показали, что наиболее информативными масштабами являются масштабы а = 1...5. Интенсивность помех на данных масштабах значительно возрастает, что позволяет их обнаружить на основе операции (4). На более крупных масштабах возрастают амплитуды вариаций естественных возмущений (фона).
Для оценки порогов Т^, Г® и Г/ использовались все суточные вариации с землетрясениями, произошедшими в спокойные (К-индекс = 1) и слабовозмущенные (К-индекс = 3) периоды. Результаты расчета средней интенсивности фона Е/оп, максимальной Е^птах и минимальной Еопт1п интенсивности помех от ионозонда, а также максимальной Е:еттах и минимальной Е:етт1п интенсивности помех от землетрясений в периоды спокойного и слабовозмущенного геомагнитного поля представлены в Табл. 1, 2. Анализ результатов в Табл. 1,2 подтверждает эффективность предлагаемого метода и показывает наибольшую информативность Б-компоненты.
Таблица 1
Интенсивности помех и фона в периоды спокойного геомагнитного поля (К- индекс = 1) в D-, Z-компонентах для различных масштабов а (указаны
в скобках)
K = 1 H (1 - 5) H (1 -10) H (1 - 15) D (1 - 5) D (1 -10) D (1 - 15) Z (1 - 5) Z (1 -10) Z (1 - 15)
Eionmax 0.0579 0.1733 0.2345 0.235 0.4402 0.7007 0.7887 1.5484 2.4203
Elonmin 0.0099 ' 0.0247 0.0339 0.0113 0.0239 0.0339 0.0502 0.883 0.1234
Ezemmax 0.4767 0.5519 0.634 1.913 2.5188 3.6875 0.1599 0.5227 0.7396
Ezemmin 0,0045 0,0215 0,0377 0,0024 0,0148 0,0276 0,0144 0,0269 0,0544
Efon 0.0373 0,1331 0.1852 0.0447 0.1727 0.2415 0.1049 0.4214 0.5896
На рис. 5 представлены результаты оценки максимальных и минимальных значений интенсивностей для помех от ионозонда, землетрясений и для естественного фона в Н-, Б- и 2- компонентах для различных диапазонов масштабов вейвлет-разложения.
Н- компонента
10«ЦхО,5 с)
Рис. 2. Интенсивности техногенных и естественных возмущений в Н- составляющей, рассчитанные для различных диапазонов масштабов a
Таблица 2
Интенсивности помех и фона в периоды слабовозмущенного геомагнитного поля (К- индекс = 3) в Н-, Э-, 2-компонентах для различных масштабов a
(указаны в скобках)
K = 3 H (1 - 5) H (1 -10) H (1 - 15) D (1 - 5) D (1 -10) (1 - 15) Z (1 - 5) Z (1 -10) Z (1 - 15)
Eionmax 0.0363 0.0873 0.1097 0.0953 0.18 0.2721 0.3836 0.7209 1.1254
Eionmin 0.0086 0.0324 0.0484 0.0134 0.0265 0.0311 0.0405 0.0698 0.0874
Ezemmax 0.48 0.1102 0.1529 0.4871 0.9702 1.4762 0.0847 0.2847 0.378
Ezemmin 0.0059 0.0128 0.0285 0.0115 0.0287 0.0456 0.0174 0.0593 0.11084
Efon 0.0173 0.0626 0.0874 0.0143 0.0547 0.0764 0.0525 0.2152 0.2999
Результаты показывают широкий разброс интенсивностей помех, и в некоторые моменты времени интенсивность фона может превышать интенсивность помехи.
Рис. 3. Интенсивности техногенных и естественных возмущений в Б- компоненте, рассчитанные для различных диапазонов масштабов а
Рис. 4. Интенсивности техногенных и естественных возмущений в 2- компоненте, рассчитанные для различных диапазонов масштабов а
Н- компонента Б- компонента Z- компонента
0.7 - 4 3
ИЗ 1-Ю 1-15 1-5 МО 1-15 1-5 1-10 1-15
Рис. 5. Интенсивности помех от ионозонда (Егоп), землетрясений (Е:ет) и естественного фона (Еуоп) в Н-, Б- и 2-компонентах для масштабов вейвлет-разложения а=1^5, а=1^10, а=1^15
На рис. 6 показаны результаты обнаружения помех предлагаемым методом. Для масштабов разложения а = 1,..., 5 использовались значения порогов Т^ =0.042, Гр =0.07 и Т^ =0.15 , определенные путем оценки апостериорного риска. В результате обработки идентифицированно соответственно 61363, 16098 и 9055 измерений, содержащих помехи.
clear norm Н
1 t
Г 1 1
О 2 4 6 8 10 12 14 16 18
х104 t(x 0,5 С)
clear norm D
_I_u_i_i_i_i_i_i_
0 2 4 6 8 10 12 14 16 18
x io* t(x 0,5 c)
clear norm Z
0 2 4 6 8 10 12 14 16 18
хЮ" I (X 0,5 С)
Рис. 6. Результаты идентификации помех в период 01.12.2016 (красным цветом отображены данные, в которых удалены помехи, синим - удаленные фрагменты записей с помехами)
Заключение
Предложенный метод выделения помех в геомагнитных данных показал свою эффективность и может применяться для первичной обработки данных магнитометров. Эффективность способа показана на примере анализа результатов измерений с помощью вариационного феррозондового магнитометра РОЕ-ОТи (обсерватория «Па-ратунка» ИКИР ДВО РАН, Камчатский край). Рассмотренные помехи от естественных источников (землетрясения на Камчатке) и техногенные, связанные с работой ионозонда, наиболее выражены в вариациях О- и 2-составляющих магнитного поля. Оценки показали широкий разброс интенсивностей помех, и в некоторые моменты времени интенсивность фона может превышать интенсивность помехи. Наиболее информативными масштабными уровнями вейвлет-преобразования являются =Ы-5.
Авторы планируют продолжить исследование в указанном направлении и на основе комплексной обработки различных составляющих магнитного поля оптимизировать работу метода.
Список литературы/References
[1] Khomutov S. Y. et al., "Noise in raw data from magnetic observatories", Geosci. Instrum. Method. Data Syst., 6:2 (2017), 329-343.
[2] Lewicki M.S., "A review of methods for spike sorting: the detection and classification of neural action potentials", Network: Computation in Neural Systems, 9:4 (1998), R53-R78.
[3] Hulata E. et al., "Detection and sorting of neural spikes using wavelet packets", Physical Review Letters, 85:21 (2000), 4637-4640.
[4] Pavlov A.N. et al., "Sorting of neural spikes: When wavelet based methods outperform principal component analysis", Natural Computing, 6 (2007), 269-281.
[5] Soloviev A., Agayan S., Bogoutdinov S., "Estimation of geomagnetic activity using measure of anomalousness", Annals of Geophysics, 59:6 (2016).
[6] Letelier J., Weber P., "Spike sorting based on discrete wavelet transform coefficients", Journal of neuroscience methods, 101 (2000), 93-106.
[7] Quian Quiroqa R., Nadasdy Z., Ben-Shaul Y., "Unsupervised spike detection and sorting with wavelets and superparamagnetic clustering", Neural Computation, 16 (2004), 16611687.
[8] Короновский A.A. и др., Вейвлеты в нейродинамике и нейрофизиологии, Физматлит, М., 2013, 272 с. [Koronovskiy A.A. i dr., Veyvlety v neyrodinamike i neyrofiziologii, Fizmatlit, M., 2013, 272 pp., (in Russian)].
[9] Жижикина Е. А., Мандрикова О. В., Хомутов С. Ю., "Алгоритм выделения техногенных помех в геомагнитных данных", Вестник КамчатГТУ, 35 (2016), 21-26. [Zhizhikina Ye. A., Mandrikova O. V., Khomutov S.YU., "Algoritm vydeleniya tekhnogennykh pomekh v geomagnitnykh dannykh", Vestnik KamchatGTU, 35 (2016), 21-26, (in Russian)].
[10] Левин Б. Р., Теоретические основы статистической радиотехники. В трех книгах. Книга вторая. Изд. 2-е, перераб. и дополнен., Сов. Радио, М., 1975, 392 с. [Levin B. R., Teoreticheskiye osnovy statisticheskoy radiotekhniki. V trekh knigakh. Kniga vtoraya. Izd. 2-ye, pererab. i dopolnen., Sov. Radio, M., 1975, 392 pp., (in Russian)].
[11] Добеши И., Десять лекций по вейвлетам, НИЦ «Регулярная и хаотическая динамика», Ижевск, 2001, 464 с. [Dobeshi I., Desyat' lektsiy po veyvletam, NITS «Regulyarnaya i khaoticheskaya dinamika», Izhevsk, 2001, 464 pp., (in Russian)].
Список литературы (ГОСТ)
[1] Khomutov S. Y. et al. Noise in raw data from magnetic observatories // Geosci. Instrum. Method. Data Syst. 2017. vol. 6. no. 2. pp. 329-343.
[2] Lewicki M. S. A review of methods for spike sorting: the detection and classification of neural action potentials // Network: Computation in Neural Systems. 1998. vol. 9. no. 4. R53-R78.
[3] Hulata E. et al. Detection and sorting of neural spikes using wavelet packets // Physical Review Letters. 2000. vol. 85. no. 21. pp. 4637-4640.
[4] Pavlov A. N. et al. Sorting of neural spikes: When wavelet based methods outperform principal component analysis // Natural Computing. 2007. no. 6. pp. 269-281.
[5] Soloviev A., Agayan S., Bogoutdinov S. Estimation of geomagnetic activity using measure of anomalousness // Annals of Geophysics. 2016. vol. 59. no. 6.
[6] Letelier J., Weber P. Spike sorting based on discrete wavelet transform coefficients // Journal of neuroscience methods. 2000. vol. 101. pp. 93-106.
[7] Quian Quiroqa R., Nadasdy Z., Ben-Shaul Y. Unsupervised spike detection and sorting with wavelets and superparamagnetic clustering // Neural Computation. 2004. no. 16. pp. 1661-1687.
[8] Короновский A.A. и др. Вейвлеты в нейродинамике и нейрофизиологии М.: Физматлит. 2013, 272 c.
[9] Жижикина Е. А., Мандрикова О. В., Хомутов С. Ю. Алгоритм выделения техногенных помех в геомагнитных данных // Вестник КамчатГТУ. 2016. №35. С. 21-26.
[10] Левин Б. Р. Теоретические основы статистической радиотехники. В трех книгах. Книга вторая. Изд. 2-е, перераб. и дополнен. М.: Сов. Радио, 1975. 392 c.
[11] Добеши И. Десять лекций по вейвлетам. Ижевск: НИЦ «Регулярная и хаотическая динамика», 2001. 464 c.
Для цитирования: Папшева С. Ю., Мандрикова О. В., Хомутов С. Ю. Метод обнаружения помех в геомагнитных данных // Вестник КРАУНЦ. Физ.-мат. науки. 2019. Т. 29.
№4. C. 87-97. DOI: 10.26117/2079-6641-2019-29-4-87-97
For citation: Papsheva S. Y., Mandrikova O. V., Khomutov S.Y. Method of noise detection
in magnetic data, Vestnik KRAUNC. Fiz.-mat. nauki. 2019, 29: 4, 87-97. DOI: 10.26117/20796641-2019-29-4-87-97
Поступила в редакцию / Original article submitted: 26.11.2019
Vestnik KRAUNC. Fiz.-Mat. Nauki. 2019. vol. 29. no.4. pp. 87-97.
DOI: 10.26117/2079-6641-2019-29-4-87-97
MSC 86-05
METHOD OF NOISE DETECTION IN MAGNETIC
DATA
S. Y. Papsheva, O. V. Mandrikova, S. Y. Khomutov
Institute of Cosmophysical Research and Radio Wave Propagation FEB RAS,
684034, Paratunka, Mirnaya st., 7, Russia
E-mail: 75sp19@mail.ru, oksanam1@mail.ru, khomutov@ikir.ru
The method of detection of noise in magnetic data based on wavelet transformation and threshold functions is considered. Efficiency is shown by the results of analysis of fluxgate magnetometer FGE-DTU measurements at Observatory Paratunka, Kamchatka, IKIR FEB RAS. The noise from natural sources such as earthquakes in Kamchatka region and from artificial sources such as the vertical sounding of ionosphere by the ionosonde near Observatory is considered. Detailed time-frequency structure of noise in 2 Hz records of Z and D components is investigated. To automation the method for considered examples of noise the informative scale levels of wavelet-transformation are determined and parameters of threshold functions appreciated.
Key words: noise in magnetic data, wavelet transform, data analysis
© Papsheva S.Y., Mandrikova O.V., Khomutov S.Y., 2019