УДК 622.831.32:681.5.08
Константинов Александр Владимирович
младший научный сотрудник, аспирант, Институт горного дела ДВО РАН, 680000, г. Хабаровск, ул. Тургенева, 51 e-mail: [email protected]
Гладырь Андрей Владимирович
старший научный сотрудник, Институт горного дела ДВО РАН e-mail: [email protected]
Ломов Михаил Андреевич
младший научный сотрудник, Институт горного дела ДВО РАН
РАЗРАБОТКА АЛГОРИТМА АВТОМАТИЧЕСКОЙ ИДЕНТИФИКАЦИИ СЕЙСМОАКУСТИЧЕСКИХ СИГНАЛОВ СРЕДСТВАМИ ЛОКАЛЬНОГО МОНИТОРИНГА
DOI: 10.25635/2313-1586.2019.02.043
Konstantinov Alexander V.
Senior Researcher Institute of Mining, FEB RAS, Khabarovsk, 51 Turgenev Str., e-mail: alex-sdt@yandex. ru
Gladyr Andrey V.
Senior Researcher, Institute of Mining, FEB RAS, e-mail: [email protected]
Lomov Mikhail A.
Junior Researcher, Institute of Mining, FEB RAS
DEVELOPMENT OF ALGORITHM FOR AUTOMATIC IDENTIFICATION OF SEISMIC ACOUSTIC SIGNALS BY MEANS OF LOCAL MONITORING
Аннотация:
Рассматривается вопрос создания эффективного алгоритма идентификации сейсмоакусти-ческих импульсов. Разработанный алгоритм станет частью программных средств многоканального прибора локального контроля ударо-опасности в массиве горных пород. Проблема выделения полезных сигналов из потока регистрируемых данных о геомеханическом состоянии решается путем определения моментов времени начала и окончания импульсов. Начальный момент времени определяется с использованием быстрой и медленной экспоненциальных скользящих средних, служащих для сглаживания вариаций данных и характеризующих отношение сигнал - шум. Метод скользящих средних применяется совместно с адаптивным порогом, позволяющим повысить точность определения начального момента, подстраиваясь под шумовой фон. Конечный момент времени полезного сигнала детектируется при помощи кумулятивной огибающей, рассчитываемой с момента времени обнаружения начала импульса и прекращается при достижении порогового значения. Разработанный математический алгоритм позволяет с достаточной точностью идентифицировать сейсмоакустические сигналы и является нетребовательным к аппаратным ресурсам прибора, что позволяет использовать его в режиме реального времени.
Abstract:
The paper deals with the issue of creating an effective algorithm to identify seismoacoustic impulses. The developed algorithm will become a part of the software tools of the multichannel device of local impact control in the rock mass. The problem of separating the useful signals from the stream of registered data on the geomechanical state is solved by determining the moment of start and end of impulses. The initial moment of time is determined using the fast and slow exponential moving averages, which are used to smooth the data variations and characterize the signal-to-noise ratio. The method of moving averages is used in conjunction with the adaptive threshold, which allows to improve the accuracy of the initial moment detection, adjusting to the noise background. The end point of time of the useful signal is detected by the cumulative envelope calculatedfrom the moment of detection of the pulse start and terminated when the threshold value is reached. The developed mathematical algorithm allows to identify seismic acoustic signals with sufficient accuracy and does not require the hardware resources of the device, which allows it to be used in real time.
2 S X
re £
a
<u EC О
U it
Ключевые слова: горный массив, удароопасное состояние, обработка цифровых сигналов, локальный контроль, акустическая эмиссия, геомеханический мониторинг
Key words: rock mass, rock-bump hazardous state, digital signal processing, local control, acoustic emission, geomechanical monitoring
2
Введение
Для обеспечения безопасного ведения горных работ процесс разработки должен сопровождаться непрерывным мониторингом состояния горного массива [1 - 3]. Одним из наиболее перспективных направлений по осуществлению такого контроля в подземных выработках являются геомеханические системы мониторинга, которые в совокупности с высокой оперативностью и низкой трудоемкостью допускают применение в различных горно-геологических условиях и обеспечивают достаточную точность снимаемых показаний.
Среди широкого спектра применяемых средств выделим класс приборов локального контроля, решающий задачи оперативного краткосрочного прогноза удароопасно-сти рабочих и нерабочих выработок. Средства локального контроля удароопасности, относящиеся к геомеханическим системам мониторинга, имеют набор характерных особенностей:
- охват локального участка (в пределах 30 - 50 м) горного массива с возможностью произвольного выбора контролируемого участка на территории шахтного поля в зоне доступа оператора;
- мобильность в процессе проведения измерений и гибкость по настройке параметров конкретного устройства;
- легкая замена чувствительных элементов, совместное использование их разных типов.
Ярким представителем класса средств контроля удароопасности является прибор локального мониторинга PrognozL, разработанный в Институте горного дела ДВО РАН и сочетающий в себе современные технические средства и гибкое программное обеспечение [4]. Данный прибор позволяет выполнять краткосрочный прогноз удароопасности локального участка горного массива, а также формировать заключение о безопасности дальнейшего производства горных работ в контролируемом объеме. Практическое использование и лабораторные испытания прибора локального контроля показали достаточную точность и адекватность оценки категории удароопасности, соответствующей начальному моменту разрушения для крепких вмещающих пород, а гибкость настройки параметров расчета позволяет проводить корректировку аналитических алгоритмов с учетом особенностей практически любого горного предприятия.
Однако с учетом достоинств прибора PrognozL в процессе многолетнего использования была сформулирована потребность в модернизации, включающей конструктивные, программные и аппаратные улучшения. Задача, рассматриваемая в данной статье, решается на программном уровне. В работе [5] предлагается использовать полнофункциональную операционную систему с высокой степенью интеграции и широкими возможностями по встраиванию Windows 10 IoTEnterprise в отличие от применяемой в настоящее время управляющей программы, написанной на языке низкого уровня Си [6]. Предлагаемая операционная система позволяет задействовать возможности продвинутой графической среды, расширенные средства объектно-ориентированных языков с высоким уровнем абстракции и разработанные сторонними разработчиками модули и библиотеки, включая пакеты математического анализа.
Постановка задачи
Основная цель модернизации прибора локального контроля состоит в повышении качества и надежности прогнозирования удароопасности горного массива с возможно-
2 стью локации источников сейсмоакустических событий и определении параметров очага
ге их возникновения [7 - 9]. Данная задача решается путем расширения количества каналов
о. регистрации геомеханического состояния, расширения динамического и частотного диа-
сс пазонов и разработки соответствующего математического аппарата, позволяющего ана-
о и
лизировать сейсмоакустические данные в реальном времени и определять параметры от-
дельных импульсов, групп импульсов, очаговых зон. Создание качественно новой и эффективной методики оценки удароопасности средствами локального контроля позволит обеспечить безопасную обстановку для работы шахтного персонала.
Упрощенная структурная схема обработки сейсмоакустических данных разрабатываемым прибором представлена на рис. 1.
Д1 ПУ1
Аппаратная
Программная часть
Д2
ПУ2
Математический модуль
ПУ8
АЦП
Цифрово
й фильтр
Буфер данных
Карта памяти
____________I
2 X
ге *
а
ф
ЕС
о и
Рис. 1 - Упрощенная структурная схема обработки сейсмоакустических данных
Обработка сейсмоакустической информации состоит из следующих этапов:
1. Сигнал регистрируется специализированными датчиками (Д1-Д8).
2. Информация о параметрах волновых процессов в массиве горных пород, регистрируемая датчиками, поступает на предварительные усилители (ПУ1-ПУ8), где решается задача согласования уровней датчика и АЦП.
3. АЦП переводит аналоговую информацию в цифровую форму путем дискретизации и квантования.
4. Оцифрованные данные поступают в модуль цифровой фильтрации, в котором при помощи механизма полосового фильтра производится подавление нежелательных низких и высоких частот.
5. В буфере данных накапливаются оцифрованные данные, которые по мере необходимости сохраняются в виде многоканального аудиофайла на карте памяти прибора и последовательно передаются в модуль математического анализа для дальнейшей обработки.
6. В математическом модуле выполняется идентификация импульсов акустической эмиссии и расчет их параметров; фильтрация шумов, определение принадлежности импульсов к источнику, выделение акустически активных зон, определение удароопас-ности участков горного массива.
В настоящей работе ставится задача разработки алгоритмов эффективной идентификации сейсмоакустических сигналов с расчетом их параметров в условиях шумовой обстановки природного и техногенного характера.
Для обработки цифровых данных воспользуемся пакетом прикладных программ для решения задач технических вычислений МАТЬАБ и одноименным языком программирования, применяемым в этом пакете.
Параметрическое описание сейсмоакустического сигнала
На первоначальном этапе необходимо обозначить характеристики импульса. Для этого рассмотрим рис. 2, представляющий участок сейсмоакустических данных с ярко выраженным фронтом.
Рис. 2 - Волновая форма сейсмоакустических импульсов: 11 и ¿2 - соответственно, моменты начала и окончания сигнала; ta - момент времени с максимальной амплитудой, A - амплитуда сигнала, Ao - амплитуда шума
2 X
ге *
а
ф
ЕС
о и
В качестве параметрического описания сигнала выделим следующие характеристики:
¿1, ¿2 -момент начала и окончания сигнала; ¿2 - ¿1 - длительность сигнала; ¿а - момент времени с максимальной амплитудой; ¿а - ¿1 - длительность фронта; А - амплитуда сигнала; Ао - средняя амплитуда шума;
/ 3 Л(^количество осцилляций, равное числу пересечения сигналом нулевого ^ 1
значения амплитуды; • площадь сигнала, вычисляемая как .
$ТЛ/ЬТЛ детектор с адаптивным порогом
Чтобы успешно идентифицировать сейсмоакустический импульс, требуется определить момент вступления волн полезного сигнала. Указанная задача осложняется тем, что среда проводимых измерений является зашумленной со средней амплитудой Ао, и поэтому необходимо эффективно отфильтровывать помехи, образуемые в результате наложения шумов природного и техногенного характера.
С целью нахождения начала сейсмоакустического сигнала рассмотрим участок данных, предшествующих его появлению (рис. 3).
На рис. 3 можно отчетливо выделить начало формирования полезного сигнала. Для автоматического определения воспользуемся методом БТА/ЬТА (8Ьог1:Т1шеАуега§е1:оЬоп§Т1шеАуега§еХ описанном в работах Фрейбергера и его последователей [10]. Метод заключается в расчете средних значений амплитуд в коротком и длинном временных окнах:
п-1
БМА(г) =^^-1)1
П ¿—I
¿ = 0
2 s
X
ге £
а
ф
ЕС
о и
где 5MЛ(t) - значение скользящей средней в отчет времени t,
п - число отсчетов в коротком (длинном) временном окне, называемое «сглаживающий интервал»;
у^ — I) - значение амплитуды сигнала в отчет t — /.
Данный метод является нетребовательным к вычислительным ресурсам и позволяет его использовать в системах реального времени [11]. Однако, чтобы еще более снизить требования к оборудованию и повысить производительность вычислений, предлагается обрабатывать поступающий от АЦП сигнал частями, используя экспоненциальные скользящие средние, рассчитываемые по формуле
ЕМА(г) = «• у(0 + (1—к) • ЕМА(г — 1), где ЕМА(€) - значение экспоненциальной скользящей средней в отчет времени ^
а - коэффициент, характеризующий скорость уменьшения значения амплитуды, называемый также сглаживающей константой, 0 <<х < 1; у(I) - значение амплитуды сигнала в текущем отчете ^
ЕМА(1 — 1) - значение экспоненциальной скользящей средней в предыдущем отчете времени t — 1. .
Рис. 3 - Волновая форма начального участка сейсмоакустического сигнала с предысторией
Использование экспоненциальных скользящих средних позволяет снизить сложность алгоритма с 0(N) для простых скользящих средних, где N- число отчетов анализируемой записи, до 0(11) для экспоненциальных, т. е. алгоритм будет работать за константное время.
Таким образом, путем выбора значений сглаживающих констант можно определить значение быстрой (STA) и медленной (LTA) экспоненциальных средних. Константы для определения STA, LTA были подобраны эмпирически путем обработки множества сей-смоакустических данных и составляют для частоты дискретизации 100 кГц 0,001 и 0,2, соответственно.
Параметр позволяет оценить отношение сигнал - шум в конкретный момент времени и считать началом полезного сигнала момент времени, когда выполняется условие r(t{) > c > 1 , где c - пороговое значение, определяемое эмпирическим путем и зависящее от условий проведения измерения; r(t1) = STA(t1)/LTA(t1) - отношение STA/LTA
2 X
ге *
а
ф
ЕС
о и
в момент времени являющегося началом сейсмоакустического сигнала. На рис. 4 представлен начальный участок сигналограммы импульса, графики функций STA/LTA детектора и адаптивного порога, а также выделена область отчетов полезного сигнала.
Рис. 4 - Волновая форма начального участка сейсмоакустического сигнала с отношением сигнал/шум и порогом
Согласно работам Фрейбергера значение пороговой линии задается как константа. С целью повышения точности проводимых расчетов требуется понижать данное значение, что может приводить к увеличению частоты ложных срабатываний.
Чтобы решить указанную проблему, предлагается использовать адаптивное пороговое значение [12,13], определяемое по формуле
с(0 = ЕМА06(г) + ЕМА02(г) • 0.05 + сп, где ЕМАОвЮ, ЕМА02(€) - значения экспоненциальных скользящих средних, сглаживающие коэффициенты которых определяются в зависимости от константного значения предшествующих секунд и частоты дискретизации сигнала;
сп - смещение порога детектора, предотвращающее ложные срабатывания и определяемое опытным путем.
Описанный STA/LTA детектор делает возможным автоматически детектировать момент начала формирования импульса акустической эмиссии в режиме реального времени либо по уже сформированным данным, например, wav-файлу.
Определение продолжительности акустического импульса
Следующим этапом алгоритма идентификации полезного сигнала выступает определение момента окончания импульса. Для достижения данной цели предлагается воспользоваться понятием кумулятивной огибающей отношения STA/LTA [11], рассчитываемой по формуле:
сз(г) = • Г(I - (•)],
1=1
где Ъ < 1 - константа, гарантирующая убывание функции сб^) на участках, где не наблюдается полезный сигнал.
2 s
X
ге £
а
ф
ЕС
о и
Тогда моментом окончания полезного сигнала ^ можно считать отсчет времени, когда выполняются следующие условия: t > t0 и <0. С этого момента происходит прекращение расчета функции и производится поиск начала следующего полезного сигнала. Алгоритм поиска окончания импульса выполняется за амортизационное время и его сложность О (а(п)).
На рис. 5 приведен пример расчета функции кумулятивной огибающей.
Рис. 5 - Волновая форма сейсмоакустического сигнала с функцией кумулятивной огибающей
Стадии идентификации сейсмоакустического импульса
Подведем итоги проделанной работы и перечислим основные стадии алгоритма идентификации полезного сигнала:
1. Инициализация константных параметров: сглаживающих констант, характеризующих степень влияния предшествующих отчетов, и величины смещения порога детектора;
2. Цикл перебора сэмплов и аккумулирование импульсов с параметрическим описанием акустической эмиссии в отдельном массиве:
• на каждой итерации цикла производится расчет экспоненциальных STA, LTA и от абсолютных значений амплитуды текущего сэмпла;
• проверка флага, имеющего отрицательное значение в случае, если начало импульса еще не определено:
- если начало импульса не определено выполняется расчет адаптивного порога и сравниваются отношение STA/LTA и значение порога. В случае выполнения условия r(t) > c(t) начало импульса считается найденным, флаг переводится в положительное состояние;
- если условие проверки обнаружения полезного сигнала не выполняется, алгоритм переходит к следующей итерации;
• если флаг находится в положительном состоянии, то продолжается расчет параметров импульса и кумулятивной огибающей cs(t). Отчет с удовлетворяющимся условием cs(t) < 0 считается моментом окончания импульса акустической эмиссии, и флаг переводится в неактивное состояние.
3. Параметрическое описание обнаруженных импульсов сохраняется в файл с расширением .х1бх.
Схематическое представление алгоритма идентификации импульсов акустической эмиссии приведено на рис. 6.
Рис. 6 - Блок-схема алгоритма идентификации сейсмоакустических импульсов
Заключение
Разработанный алгоритм идентификации сейсмоакустических импульсов позволяет эффективно с достаточной точностью выявлять моменты начала и окончания полезного сигнала; является нетребовательным к аппаратным ресурсам прибора и может быть масштабирован к любому количеству входных данных. В дальнейших
работах планируется его доработка и усовершенствование, в частности, более точное определение момента начала сигнала путем вариации константных значений скользящих экспоненциальных средних, а также расчет дополнительных характеристик импульсов, например, площадь под огибающей сигнала (MARSE), представляющей собой аналог энергии.
Идентификация полезных сигналов и определение их параметрических характеристик послужат важным аспектом усовершенствованной методики контроля ударо-опасности в массиве горных пород и станут неотъемлемой частью программного обеспечения разрабатываемого прибора.
Литература
1. Концепция комплексного геодинамического мониторинга на подземных горных работах / А.Н. Шабаров, С.В. Цирель, К.В. Морозов, И.Ю. Рассказов // Физико-технические проблемы разработки полезных ископаемых. - 2017. - № 9. - С. 59 - 64.
2. Чебан А.Ю. Способ доработки глубокого карьера с применением фрезерных машин / А.Ю. Чебан // Маркшейдерия и недропользование. - 2017. - № 4. - С. 23 - 29.
3. Чебан А.Ю. Совершенствование техники и технологий безвзрывной разработки горных пород: моногр. / А. Ю. Чебан. - Хабаровск: ИГД ДВО РАН, 2017. -260 с.
4. Результаты применения геоакустического метода локального контроля ударо-опасности на рудниках Дальнего Востока / А.А. Терёшкин, И.Ю. Рассказов, П.А. Аникин, Д.С. Мигунов // Горный информационно-аналитический бюллетень. -2017. - № 24. - С. 338 - 347.
5. Konstantinov A., Development of Multi-Channel Portable Impact Control Device for Local Assessment of The State of The Edge Parts of The Rock Massif / VII International Scientific Conference "Problems of Complex Development of Georesources" (Khabarovsk). - 2018. - P. 02024.
6. Свидетельство о государственной регистрации программы для ЭВМ RUS 2017663503 (2017) / И.Ю. Рассказов, А.В. Гладырь, Д.С. Мигунов, П.А. Аникин, А.А. Терёшкин.
7. Использование данных сейсмоакустических наблюдений для определения характера развития очага разрушения породного массива / И.Ю. Рассказов, С.В. Ци-рель, А.О. Розанов, А.А. Терешкин, А.В. Гладырь // Физико-технические проблемы разработки полезных ископаемых. - 2017. - № 2. - С. 29 - 37.
8. RozanovA.O., Petrov D.N., RozenbaumA.M., TereshkinA.A., IlinovM.D.Acoustic Emission Precursor Criteria of Rock Damage / Eurock 2018. Geomechanics and Geody-namics of Rock Masses (set of 2 volumes). - 2018. - P. 669 - 672.
9. Rasskazov I.Y., Tereshkin A.A., Gladyr A.V., TsirelS. V., RozanovA.O. Application of Acoustic Measurement Data to Characterize Initiation and Development of Disintegration Focus in A Rock Mass // Journal of Mining Science. - 2017. - V. 53. - N. 2. -P. 224 - 231.
10. Joswig M. Pattern Recognition for Earthquake Detection // Bulletin of Seismologi-cal Society of America. 1990 - V. 80. -N 1. - P. 170 - 186.
11. Баранов С.В. Автоматическое определение длительности сейсмического события в режиме реального времени / С.В. Баранов // Региональный вестник молодых учёных. Геофизика. - 2004. - № 3.
s 12. Красовский А.А. Цифровая обработка в ZETLAB при идентификации сейсми-
« ческого сигнала / А.А. Красовский // Цифровая обработка сигналов. - 2010. - № 3. -
С. 70 - 76.
13. Сергиенко А.Б. Алгоритмы адаптивной фильтрации: особенности реализации 8 в MATLAB // ExponentaPro. Математика в приложениях. - 2003. - № 1. - С. 18 - 28.