Физика Земли, 2020, № 2, стр. 127-147
Различение типов микросейсмических источников по данным малоапертурных сейсмических групп
И. А. Санина 1, О. Ю. Ризниченко 1, А. Ф. Кушнир 2, А. В. Варыпаев 2, *, С. И. Сергеев 3, С. Г. Волосов 1
1 Институт динамики геосфер РАН
г. Москва, Россия
2 Институт теории прогноза землетрясений и математической геофизики РАН
г. Москва, Россия
3 Российские сети вещания и оповещения
г. Москва, Россия
* E-mail: avalex89@gmail.com
Поступила в редакцию 25.04.2019
После доработки 17.06.2019
Принята к публикации 24.06.2019
Аннотация
В настоящей работе рассматривается задача регистрации и обработки микросейсмических данных с дальнейшей их интерпретацией и принятием решений о возможных угрозах, связанных с нарушениями целостности земной среды в результате горных выработок. Указанные решения принимаются на основе обнаружения микросейсмических источников и распознавания типов их очагов в результате обработки многоканальных сейсмограмм малоапертурных сейсмических групп. Процедуры обработки сейсмограмм содержат алгоритмы обнаружения сигналов микросейсмических источников в зашумленных сейсмограммах датчиков групп и процедуры оценивания различных параметров этих источников. В работе предложен новый алгоритм обнаружения, основанный на оценке множественной когерентности сейсмограмм датчиков группы, который позволяет выявлять волновое излучение микросейсмических источников со сложными механизмами очагов на фоне диффузных сейсмических помех. В отличие от широко применяемого STA/LTA группового детектора, ориентированного на обнаружение плоских сейсмических волн, предлагаемый детектор предназначен для обнаружения волн от локальных сейсмических событий в сложно устроенной приповерхностной земной среде. С помощью численных экспериментов с модельными и реальными данными сейсмических групп, содержащими сигналы источников со сложными механизмами очагов, показано, что распознавание типов микросейсмических источников можно производить, используя различные алгоритмы анализа данных сейсмических групп: как новые алгоритмы оценивания координат микросейсмических источников, так и традиционные алгоритмы пространственно-временного спектрального анализа.
1. ВВЕДЕНИЕ
Проблема обнаружения сигналов от слабых сейсмических источников и оценки их параметров стала популярной в связи с необходимостью микросейсмического мониторинга опасных промышленных или природных объектов, таких как горные работы [Maochen, 2005; Malovichko, Lynch, 2006], вулканические зоны [Droznin et al., 2015], площадки добычи углеводородов и гидротермы [Eisner, 2010; Cros et al., 2017]. Такой мониторинг включает регистрацию и обработку сейсмических данных с дальнейшей их интерпретацией и принятием решений о возможных угрозах, связанных с нарушениями целостности земной среды.
Предметом настоящей работы было создание методики автоматизированного микросейсмического мониторинга земной среды в окрестностях шахт для разработки залежей полезных ископаемых, таких как железная руда или руды других металлов, которые содержатся в кристаллических породах. При разработке таких ископаемых порода дробится направленными взрывами, которые нарушают стабильность напряженного состояния земной среды, ослабленной шахтными выработками. Следствием таких взрывов является наведенная микросейсмическая активность в среде, которую можно интерпретировать как поток афтершоков, сопровождающих каждый взрыв.
Регистрация сейсмических волн, вызванных взрывами и наведенными микросейсмическими источниками, производится группой сейсмических станций, расположенных на земной поверхности и представляющих собой “сейсмическую антенну”, обычно называемую малоапертурной сейсмической группой. Главные особенности малоапертурной группы следующие:
а) сейсмометры (датчики) группы регистрируют сейсмические волны синхронно во времени;
б) максимальное расстояние между датчиками достаточно мало, так что земную среду под группой можно считать горизонтально однородной.
Поэтому совокупность сейсмограмм, зарегистрированных датчиками группы, могут математически интерпретироваться как многомерный стационарный процесс, и полезные сигналы от сейсмических событий, зарегистрированные различными датчиками, являются когерентными (точный смысл этого термина будет ясен из дальнейшего).
Целью работы было выделить из потока микросейсмических событий те события, которые наводятся сильными шахтными взрывами и свидетельствуют о наведенной микросейсмической активности среды. Наведенные события являются своего рода атершоками взрывов, и свидетельствуют о вызванной взрывами перестройке механических напряжений в среде, которая может привести даже к разрушению участков среды, ослабленных проходческими работами в шахтах.
С точки зрения физики очаги взрывных и наведенных микросейсмических источников различаются типами их тензоров сейсмического момента [Аки, Ричардс, 1983], которые определяют пространственные диаграммы излучения этими источниками сейсмических волн. Взрывной источник (всестороннее расширение) имеет диагональный тензор сейсмического момента и излучает P-волну с сферически однородным фронтом. Большинство природных микросейсмических событий имеют сложный тензор сейсмического момента и излучают и P-, и S-волны с существенно неоднородными волновыми фронтами. Примеры подобных диаграмм излучения приведены, в частности, в публикации [Кушнир, Варыпаев и др., 2014].
Возможности использования сейсмических групп для измерения характеристик природных и искусственных сейсмических источников изучаются c 70-х гг. прошлого века [Capon, 1969; 1973; Davies et al., 1971; Posmentier, Herrmann, 1971] и продолжаются до настоящего времени. Этому вопросу посвящена обширная литература. Были предложены различные методы направленного приема сейсмических волн, излучаемых источниками, т.е. формирования пространственной диаграммы направленности сейсмической группы, которая позволяет подавлять природные сейсмические помехи или мешающие волны от “посторонних” событий.
Процедура формирования диаграммы направленности группы, которая позволяет наилучшим образом регистрировать временные колебания плоской сейсмической волны, приходящей на группу с заданного азимута и с заданным углом выхода, и подавлять другие волны, называется процедурой формирования луча (beamforming) сейсмической группы. Процедура формирования луча позволяет измерять мощность плоской волны, приходящей с заданного направления. Измерение с помощью “веера” лучей мощности плоских волн, приходящих на группу с разных направлений, называют пространственным спектрально-временным анализом (F-K-analysis).
Ниже мы будем использовать термины: “процедура фокусировки” или “процедура стекинга (stacking) сейсмограмм группы” для обозначения таких алгоритмов обработки многоканальной сейсмограммы малоапертурной сейсмической группы, которые позволяют измерять мощности сейсмических волн, излучаемые точечным источником, расположенным в заданной точке упругой среды. Частный случай такой процедуры был назван в работе [Kiselevitch et al., 1991] сейсмической эмиссионной томографией (СЭТ) и этот термин закрепился в литературе по микросейсмическому мониторингу. Отметим, что все процедуры фокусировки требуют задания скоростной модели упругой среды (в простейшем случае – указания, что среда предполагается однородной).
В настоящей работе мы использовали оба метода измерения мощности сейсмических волн, приходящих на группу: фазовый вариант пространственного спектрального анализа [Варыпаев и др., 2018] и фазовый вариант сейсмической эмиссионной томографии [Kushnir, Varypaev, 2017].
Автоматическая обработка сейсмограмм групп предполагает последовательное решение двух задач. Первая – это обнаружение слабых сейсмических сигналов, наблюдаемых на фоне сейсмических помех. Вторая задача – оценка параметров сейсмического источника, таких как его пространственные координаты и тензор сейсмического момента или (в случае плоских сейсмических волн) как азимут и кажущаяся медленность плоской волны. Очевидно, что эти две задачи тесно связаны, поскольку оценка параметров источника может быть произведена только в течение интервала времени, в котором наблюдаются его сигналы.
Как уже обсуждалось в работе [Кушнир, Варыпаев и др., 2014], при разработке систем микросейсмического мониторинга инженеры-геофизики сталкиваются с такими проблемами как низкое отношение сигнал/шум (SNR) в регистрируемых сигналах и сложные механизмы очагов сейсмических источников. В этом случае ряд широко используемых процедур фокусировки, например, сейсмическая эмиссионная томография [Duncan et al., 2010], становятся непригодными для оценки параметров сейсмического источника с необходимой точностью [Anikeev et al., 2014; Кушнир, Варыпаев и др., 2014]. Более того, низкое значение SNR вообще не позволяет обнаруживать сигналы источника со сложными диаграммами излучения в зашумленных сейсмограммах группы с использованием стандартных алгоритмов формирования лучей группы и STA/LTA детектора.
В данной работе мы разработали новый алгоритм многоканального обнаружения сигналов сейсмических источников, основанный на оценивании меры когерентности многомерной сейсмограммы группы. Этот алгоритм является эффективным даже при коротких интервалах существования полезных сигналов, содержащихся в сейсмограммах группы. Концепция этого алгоритма основана на оценке сингулярных значений матрицы спектральной плотности мощности (МСПМ) сейсмограммы группы.
Также использовались фазовые алгоритмы оценивания положения сейсмического источника (фазовые процедуры фокусировки) и фазовые алгоритмы F-K-анализа. Фазовые алгоритмы применялись потому, что они обладают свойством робастности, то есть устойчивости к изменениям статистических характеристик полезных сигналов и маскирующих их помех [Kushnir, Varypaev, 2017]. Указанные методы обнаружения сигналов сейсмических источников и оценки параметров микросейсмических источников со сложным тензором сейсмического момента позволяют создать компьютерную систему обработки данных малоапертурных групп, которая может использоваться для микросейсмического мониторинга с целью как распознавания типов микросейсмических источников, так и оценки параметров этих источников.
Предложенная методика автоматической обработки данных сейсмических групп была опробована на реальных сейсмических записях поверхностной малоапертурной группы, установленной в районе Курской магнитной аномалии в непосредственной близости от шахт по добыче железной руды. Основной целью этой экспериментальной работы было исследование микросейсмической активности среды, которая наводится подземными шахтными взрывами, производимыми при добыче железной руды.
С помощью описываемых ниже алгоритмов были обработаны многоканальные сейсмограммы группы, зарегистрированные в течение 20 ч. Сначала производилось детектирование сигналов от микроземлетрясений и выделялись короткие временные интервалы сейсмограмм группы, имеющие высокие значения меры когерентности. Затем эти интервалы записей обрабатывались с целью оценивания параметров микросейсмических источников, таких как координаты источников и направления прихода на группу P-волн, излучаемых этими источниками [Zhang et al., 2008; Варыпаев и др., 2017].
За 20 ч регистрации данных группы были обнаружены сильные локальные подземные взрывы и слабые микроземлетрясения. У некоторых землетрясений был выявлен сложный механизм их очагов. Для оценивания указанных параметров использовались робастные фазовые алгоритмы. С использованием многоканального фильтра Кейпона [Capon, 1996] было также продемонстрировано сходство временных форм сигналов от микроземлетрясений в диапазоне частот от 15 до 30 Гц. А именно, многоканальная фильтрация сейсмограмм группы показала, что длительность наблюдаемых сигналов микросейсмических источников очень коротка (менее 0.2 с), а их форма схожа с хорошо известной формой Риккер-вейвлета [Bording, 1996].
2. МЕТОДЫ ОБНАРУЖЕНИЯ СИГНАЛОВ МИКРОСЕЙСМИЧЕСКИХ ИСТОЧНИКОВ И ОЦЕНИВАНИЯ ИХ ПАРАМЕТРОВ ПРИ СЛОЖНЫХ МЕХАНИЗМАХ ОЧАГОВ
Математическая модель наблюдаемых данных сейсмической группы, используемая в настоящей работе, определяется следующими предположениями:
А1. Пространственная область очага микросейсмического события предполагается настолько малой, что возможно использовать простую теорию точечного сейсмического источника. Физически это означает, что размеры очага микросейсмического источника значительно меньше, чем апертура сейсмической группы и расстояний от источника до датчиков группы.
А2. Земная среда между сейсмическим источником и датчиками группы рассматривается как горизонтально-однородная упругая среда, поскольку предполагается, что глубина микросейсмических очагов достаточно мала. Это позволяет использовать лучевую теорию распространения сейсмических волн и рассматривать результат распространения сейсмической волны от источника к l-му датчику группы как результат изменений амплитуд и фаз частотных составляющих временной функции волны.
А3. Пространственная диаграмма излучения микросейсмического источника полностью определяется его тензором сейсмического момента [Аки, Ричардс, 1982].
А4. Для целей работы достаточно анализировать многоканальную сейсмограмму группы на достаточно малых интервалах наблюдений, содержащих только колебания P-волны микросейсмического источника, т.е. на тех интервалах, которые следуют за первым вступлением последовательности волн, генерируемых источником.
А5. Сейсмические помехи, искажающие многомерный сигнал сейсмического источника, регистрируемый сейсмической группой, могут описываться как многомерный стационарный случайный процесс.
Предположения А.1–А.5. приводят к следующей математической модели оператора преобразования временной функции ${{s}_{0}}\left( t \right)$ точечного источника в сейсмограмму ${{s}_{l}}\left( t \right)$ l-го датчика группы:
(1)
$\begin{gathered} {{s}_{l}}\left( {t\,{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right) = {{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right){{\gamma }_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right){{a}_{l}}\left( {\mathbf{r}} \right){{s}_{0}}\left( {t - {{\tau }_{l}}\left( {\mathbf{r}} \right)} \right), \\ l \in \overline {1,m} , \\ \end{gathered} $Сейсмограммы ${{s}_{l}}\left( {t\,{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right)$ маскируются сейсмическими помехами ${{\xi }_{k}}\left( t \right),$ аддитивно воздействующими на датчики группы, в результате чего наблюдаемые сейсмограммы на выходах датчиков группы имеют вид:
(2)
$\begin{gathered} {{x}_{l}}\left( {t\,{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right) = {{d}_{l}}\left( {{\mathbf{\theta }},{\mathbf{r}}} \right)\gamma \left( {{\mathbf{\theta }},{\mathbf{r}}} \right){{a}_{l}}\left( {\mathbf{r}} \right) \times \\ \times \,\,{{s}_{0}}\left( {t - {{\tau }_{l}}\left( {\mathbf{r}} \right)} \right) + {{\xi }_{l}}\left( t \right),\,\,\,\,l \in \overline {1,m} , \\ \end{gathered} $Непрерывные сейсмограммы датчиков (2) квантуются по времени с частотой дискретизации ${{f}_{s}}$, т.е. каждым датчиком регистрируются отсчеты:
Обнаружение микросейсмического источника, т.е. обнаружение его сигналов ${{s}_{l}}\left( {t{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right)$ в сейсмограммах группы ${{x}_{l}}\left( {t\,{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right),$ и оценивание значений векторных параметров ${\mathbf{r}},\theta $ производится по дискретным выборкам отсчетов многоканальной сейсмограммы группы ${\mathbf{X}}_{{\tau }}^{{m,n}}$ = $\left\{ {{{x}_{{l,k}}},l \in \overline {1,m} ,k \in \overline {\tau + 1,\tau + n} } \right\},$ где $\tau \in {{Z}^{ + }}$ – начальные моменты выборок, n – число отсчетов (размер выборки).
Для построения алгоритмов обнаружения многоканального сигнала источника и оценивания его параметров мы используем модель наблюдений в частотной области диcкретного конечного преобразования Фурье, которая эквивалентна модели (2) при достаточно больших n:
(3)
$\begin{gathered} {{{\dot {X}}}_{l}}\left( {{{f}_{j}}\,;\,\,{\mathbf{\theta }},{\mathbf{r}}} \right) = {{d}_{l}}\left( {{\mathbf{\theta }},{\mathbf{r}}} \right){{a}_{l}}\left( {\mathbf{r}} \right){{S}_{0}}\left( {{{f}_{j}}} \right) \times \\ \times \,\,\exp \left( {i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right) + \gamma \left( {{\mathbf{\theta }},{\mathbf{r}}} \right)\pi } \right) + {{{\dot {\zeta }}}_{l}}\left( {{{f}_{j}}} \right), \\ \end{gathered} $Обнаружение сигнала источника, маскируемого случайной помехой, трактуется как проверка по выборкам ${\mathbf{X}}_{\vartheta }^{{m,n}}$ многоканальной сейсмограммы группы в каждом скользящем временном окне:
(4)
${\mathbf{w}}\left( \vartheta \right) = \left( {{{t}_{{\vartheta + 1}}},...,{{t}_{{\vartheta + n}}}} \right),\,\,\,\,\vartheta \in {{Z}^{ + }}$Широко используется алгоритм обнаружения многомерного сигнала сейсмического источника в сейсмограмме группы, который состоит из следующих процедур:
1. Процедуры $b\left( \rho \right)$ формирования луча группы (beamforming procedure), где $\rho = \left( {\delta ,\beta } \right);$ $\delta $ – обратный азимут (back azimuth) луча; $\beta $ – угол входа (angle of incidence) луча, т.е. угол между вертикалью к земной поверхности и нормалью к фронту волны, приходящей на каждый датчик группы. С помощью этой процедуры рассчитываются одномерные временные ряды $\eta _{k}^{{\rho }},$ $k \in {{Z}^{ + }},$ $\rho \in \overline {1,L} ,$ энергия которых, в основном, определяется волновым излучением, приходящим на группу с заданного направления $\delta ,$ $\beta .$
2. Процессы $\eta _{k}^{{\rho }}$ сканируется временным окном (4), и для каждого окна ${{{\mathbf{w}}}_{{\tau }}}$ рассчитываются временные ряды $p_{\tau }^{\rho }$, $\tau \in {{{\mathbf{Z}}}^{ + }}$, $\rho \in \overline {1,L} $ значений известной статистики STA/LTA.
3. Процедуры сравнения с порогом ${{K}_{\alpha }}$, где $\alpha $ – заданная вероятность обнаружения источника, величин ${{\lambda }_{{\tau }}} = \mathop {\max }\limits_{\rho } p_{{\tau }}^{{\rho }}$ – максимальных значений временных рядов $p_{{\tau }}^{{\rho }}$ в момент $\tau .$
4. Если оказывается, что ${{\lambda }_{{\tau }}} > {{K}_{{\alpha }}}$ для заданного числа N последовательных моментов $\tau ,$ и $\mathop {\max }\limits_{\rho } p_{{\tau }}^{{\rho }}$ достигается в каждый из этих моментов для одного и того же луча $b\left( {{{\rho }_{0}}} \right),$ то выносится решение об обнаружении достаточно сильного волнового излучения, приходящего с направления $\left( {{{\delta }_{0}},{{\beta }_{0}}} \right).$
Этот алгоритм, традиционный в сейсмологии малоапертурных групп, хорошо работает для обнаружения телесейсмических и региональных событий, волны от которых приходят на группу с фронтами, несильно отличающихся от плоских. Кроме того, этот алгоритм требует задания достаточно точной скоростной модели земной среды.
Подобные алгоритмы обнаружения сигналов телесейсмических и региональных событий до сих пор используются для обработки сигналов малоапертурных сейсмических групп, расположенных в разных районах земного шара и работающих в рамках Международной системы контроля за выполнением Международного договора о всеобъемлющем запрещении ядерных испытаний, заключенном на Генеральной Ассамблее ООН 24.09.1996 г.
Аналогичные алгоритмы обнаружения можно реализовать, используя вместо процедуры формирования лучей группы различные процедуры стекинга сейсмограмм группы для измерения энергии, излучаемой из различных точек среды под группой. Эти процедуры являются основными компонентами алгоритмов оценивания координат микросейсмических источников. Ряд таких процедур описан в работе [Кушнир, Варыпаев, 2014]. Простейшей процедурой стекинга является статистика широко используемой оценки координат сейсмического источника по методу сейсмической эмиссионной томографии (СЭТ). Процедура стекинга по методу СЭТ в настоящее время широко применяется в системах микросейсмического мониторинга процессов гидроразрыва пластов земной коры при добыче сланцевого газа и сланцевой нефти.
Описанные традиционные методы обнаружения, основанные на процедурах формирования лучей группы или стекинга сейсмограмм группы, при обнаружении локальных микроскопических событий имеют существенный недостаток. Для этого есть две причины:
1. Процедура формирования лучей и процедуры стекинга используют скоростные модели земной коры между сейсмической группой и сейсмическим источником. Когда поверхностная среда под группой горизонтально неоднородна, создание достаточно точных скоростных моделей и/или их использование в компьютерных программах затруднительно.
2. И процедура формирования лучей, и традиционные процедуры стекинга, такие как СЭТ, предполагают, что знаки вступлений сейсмических волн от источника одинаковы на всех датчиках группы. Иными словами, что сигналы источника ${{s}_{l}}\left( {t\,{\text{;}}\,{\mathbf{r}},{\mathbf{\theta }}} \right),$ описываемые выражением (1), имеют одинаковые значения коэффициентов ${{\gamma }_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right)$ для всех номеров датчиков $l \in \overline {1,m} .$ Однако при малых расстояниях от источника до датчиков группы и при механизме источника, отличном от механизма взрыва, это предположение далеко не всегда выполняется.
Поэтому в ряде задач микросейсмического мониторинга представляется целесообразным применять другие методы обнаружения источников, не требующие стекинга сейсмограмм группы.
В настоящей работе мы использовали детектор, который в рамках модели наблюдений (1), (2) позволяет сделать вероятности ошибок обнаружения микросейсмических событий менее зависимыми от расстояния между источником и малоапертурной группой и от скоростных характеристик земной среды.
Новый алгоритм обнаружения основан на оценке меры когерентности многоканальных сигналов, генерируемых в датчиках группы микросейсмическим источникам. Эта мера вводится из следующих соображений. Если (без потери общности) рассматривать временную функцию источника как “отрезок” стационарного случайного процесса, то можно показать [Кушнир, 2012], что многоканальный сигнал ${\mathbf{s}}\left( {{\mathbf{r}},\theta } \right)$ = $\left( {{{s}_{1}}\left( {{\mathbf{r}},\theta } \right),...,{{s}_{m}}\left( {{\mathbf{r}},\theta } \right)} \right),$ возбуждаемый источником в датчиках группы и определяемый формулой (1), имеет вырожденную на всех частотах матричную спектральную плотность мощности (МСПМ) с рангом, равным единице. То есть, отношение максимального собственного числа этой матрицы к сумме всех остальных ее собственных чисел бесконечно велико. В то же время, МСПМ естественных случайных помех, воздействующих на датчики группы, практически всегда является невырожденной, т.е. имеет полный ранг. Поэтому многоканальная аддитивная смесь полезных сигналов с помехами (2) имеет МСПМ, для которой (при не слишком малом отношении сигнал–помеха) максимальное собственное число при всех частотах значительно превосходит сумму остальных собственных чисел. Причем, это свойство МСПМ многоканальной сейсмограммы группы не зависит от задержек временной функции сейсмической волны, излучаемой источником, которые приобретаются при распространении волны от очага до различных датчиков группы. Т.е. это свойство МСПМ не зависит от скоростной модели среды и от расстояний между источником до датчиков группы.
Предлагаемый алгоритм обнаружения многоканального сигнала микросейсмического источника на фоне помех основан на сравнении с порогом усредненного по всем частотам отношения квадрата главного собственного числа оценки МСПМ к сумме квадратов остальных собственных чисел.
Оценку МСПМ спектральных наблюдений группы (2) в скользящем временном окне (4) можно вычислять, например, по методу Блэкмана–Тьюки [Марпл, 1990; Дженкинс, Ватс, 1971]:
(5)
$\begin{gathered} {{{{\mathbf{\hat {F}}}}}_{x}}({{f}_{j}}) = \left[ {{{{\hat {F}}}_{{l,q}}}\left( {{{f}_{j}}} \right),\,\,\,l,q \in \overline {1,m} \,} \right], \\ {{{\hat {F}}}_{{l,q}}}\left( {{{f}_{j}}} \right) = \sum\limits_{k = 1}^n {w\left( {{{f}_{j}} - {{f}_{k}}} \right)\left( {\dot {X}_{l}^{*}\left( {{{f}_{k}}} \right){{{\dot {X}}}_{q}}\left( {{{f}_{k}}} \right)} \right)} , \\ {{f}_{j}} = j\frac{{{{f}_{{smp}}}}}{{2n}},\,\,\,\,j = \overline {1,n} , \\ \end{gathered} $Статистика обнаружения алгоритма имеет следующий вид:
(6)
$\overline {{{\Lambda }_{x}}} = \sum\limits_{j = {{j}_{1}}}^{{{j}_{2}}} {{{\Lambda }_{x}}\left( {{{f}_{j}}} \right)} ,\,\,\,\,{{j}_{1}} < {{j}_{2}},$(6a)
${{\Lambda }_{x}}\left( {{{f}_{j}}} \right) = \frac{{\lambda _{1}^{2}\left( {{{f}_{j}}} \right)}}{{\sum\limits_{q = 1}^{m - 1} {\lambda _{{q + 1}}^{2}\left( {{{f}_{j}}} \right)} }},$Статистика (6а) является мерой близости ранга МСПМ к единице [Jolllife, 1986].
Алгоритм обнаружения сигнала источника по наблюдениям в скользящем окне (4) может быть выражен следующим соотношением:
(7)
$\begin{gathered} {\text{Сигнал}}\,{\text{источника}}\,{\text{присутствует,}}\,{\text{если}}~\,\overline {{{\Lambda }_{x}}} \, \geqslant \,{{K}_{{\alpha }}}, \\ {\text{Сигнал}}\,{\text{источника}}\,{\text{присутствует,}}\,{\text{если}}\,\overline {{{\Lambda }_{x}}} \, < \,{{K}_{{\alpha }}}, \\ \end{gathered} $Порог обнаружения ${{K}_{{\alpha }}}$ выбирается так, чтобы обеспечить заданную вероятность $\alpha $ ложного обнаружения при некотором фиксированном среднем значении отношения сигнал/шум (СОСШ) в датчиках группы. Точный выбор ${{K}_{{\alpha }}}$ возможен лишь при известной функции распределения статистики (6) при гипотезе, что сигнал источника отсутствует. Теоретическое определение этого распределения весьма затруднительно, поэтому в настоящей работе мы оценивали необходимое значения порога методом Монте-Карло с использованием массива записей “чистого” шума, воздействовавшего на группу, когда микросейсмические источники (взрывы, землетрясения) отсутствовали. В экспериментальной части эта процедура описана более подробно.
В общем случае, ранг МСПМ ${{{\mathbf{F}}}_{x}}(f)$ многоканальной сейсмограммы группы, порожденной k < m статистически независимыми точечными источниками, равен k; при числе таких источников, большем числа датчиков группы, матричная функция ${{{\mathbf{F}}}_{x}}(f)$ становится невырожденной [Кушнир, 2012]. В настоящей работе основным предположением является наличие в каждом скользящем окне (4) только одного источника. Обнаружение в сейсмограмме группы сигналов, генерируемых несколькими одновременно действующими микросейсмическими источниками, является “ложной тревогой”, но вероятность такого события достаточно мала.
Заметим, что мера (6а) когерентности многомерного временного ряда не является единственной. В принципе, любую функцию, характеризующую ранг матричной функции ${{{\mathbf{F}}}_{x}}(f),$ можно использовать как меру когерентности сейсмограмм группы (2). Подобные меры когерентности представлены, например, в работах [Samson, Olson, 1982; Любушин, 2014], посвященных исследованию свойств сейсмических помех.
На втором этапе обработки данных малоапертурной группы фрагмент многоканальной сейсмограммы (2), соответствующий некоторому положению $\tilde {\vartheta }$ скользящего временного окна ${\mathbf{w}}\left( \vartheta \right),$ в котором был обнаружен сигнал сейсмического источника, анализируется с целью оценивания параметров ${\mathbf{r}},{\mathbf{\theta }}$ сейсмического очага. Большинство алгоритмов оценивания параметров удобнее всего представлять в частотной форме как функции спектральной выборки (3). Так, например, наиболее используемый в настоящее время алгоритм оценивания координат микросейсмических источников по методу СЭТ имеет следующую спектральную форму:
(8)
${{{\mathbf{\hat {r}}}}_{{set}}}\left( {{\mathbf{\dot {X}}}_{{m,n}}^{{\tilde {\vartheta }}}} \right)\, = \,\mathop {\arg \max }\limits_{{\mathbf{r}} \in Q} {{\sum\limits_{j = 1}^n {\left| {\sum\limits_{l = 1}^m {\dot {X}_{{l,j}}^{{\tilde {\vartheta }}}\exp \left( { - i\, \times \,2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right)} \right)} } \right|} }^{2}},$Алгоритм (8) содержит процедуру стекинга сейсмограмм группы в скользящем временном окне ${\mathbf{w}}\left( \vartheta \right){\text{:}}$
(9)
$\begin{gathered} \sum\limits_{l = 1}^m {\dot {X}_{{l,j}}^{\vartheta }\exp \left( { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right)} \right)} \mathop \Leftrightarrow \limits^{{\text{DFFT}}} \hat {s}_{k}^{\vartheta }\left( {\mathbf{r}} \right), \\ j \in \overline {1,n} ,\,\,\,\,k \in \overline {1,n} , \\ \end{gathered} $Заметим, что функция времени $\hat {s}_{k}^{{\tilde {\vartheta }}}\left( {{{{{\mathbf{\hat {r}}}}}_{{set}}}} \right)$ является оценкой временной функции источника ${{s}_{{0,k,\,\,k \in \overline {1,n} }}}.$
Алгоритм (8) оценивания вектора ${\mathbf{r}}$ и алгоритм стекинга (9) соответствуют следующим предположениям о сейсмическом источнике, среде и помехах:
Б1. Сейсмический источник имеет характер взрыва – однородного по пространству растяжения (или сжатия) среды, т.е. тензор его сейсмического момента – скалярная матрица. Это предположение соответствует модели сейсмограмм группы (2), в которой коэффициенты ${{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right),$ ${{\gamma }_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right)$ постоянны, т.е. не зависят от $l,$ ${\mathbf{r}},$ ${\mathbf{\theta }}.$
Б2. Ослабление амплитуды сейсмической волны средой при распространении от источника до датчиков группы одинаково для всех датчиков.
Б3. Сейсмические помехи, воздействующие на датчики группы, представляют собой белые Гауссовские шумы, которые статистически независимы для разных датчиков и имеют мощность, одинаковую для всех датчиков.
Эти предположения, во многих случаях естественные при анализе электромагнитных или акустических волн, никогда не выполняются в сейсмологии. Особенно критично первое предположение. Если оно неверно, то, как показывают и теоретический анализ, и эксперименты, статистическое качество оценки координат источника (8) и оценки (9) временной функции источника становится неприемлемым для задач микросейсмического мониторинга [Кушнир, Варыпаев и др., 2014; Eisner et al., 2014].
В работе [Кушнир, Варыпаев и др., 2014] предложена модификация СЭТ-оценки (8), которая позволяет обеспечить ее работоспособность при оценивании координат источников со сложными механизмами очагов, т.е. в случае, когда сейсмограммы группы описываются моделью наблюдений (2) с функциями ${{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right),$ ${{\gamma }_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right),$ ${{a}_{l}}\left( {\mathbf{r}} \right),$ соответствующими широкому классу физических моделей среды и тензора сейсмического момента миросейсмических источников:
(10)
$\begin{gathered} {{{{\mathbf{\hat {r}}}}}_{{mset}}}\left( {{{{{\mathbf{\dot {X}}}}}^{{m,n}}}} \right) = \mathop {\arg \max }\limits_{{\mathbf{r}} \in Q;\,\,{\mathbf{\theta }} \in \Theta } \sum\limits_{j = {{j}_{1}}}^{{{j}_{2}}} {\left| {\sum\limits_{l = 1}^m {{{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right){{a}_{l}}\left( {\mathbf{r}} \right){{{\dot {X}}}_{{l,j}}} \times } } \right.} \\ {{\left. {{{{_{{_{{_{{}}}}}}}}^{{^{{^{{}}}}}}} \times \,\,\exp \left( { - 2\pi i{{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right) - \pi {{{\bar {\gamma }}}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right)} \right)} \right|}^{2}}, \\ \end{gathered} $В выражении (10) для простоты опущен индекс $\tilde {\vartheta }$ временного окна.
Алгоритм (10) одновременно определяет и координаты источника r , и параметры ${\mathbf{\theta }}$ его тензора сейсмического момента, что с точки зрения математической статистики обеспечивает наименьшие ошибки их оценивания (при условия выполнения предположений Б2, Б3). Ясно, что алгоритм (10) вычислительно значительно сложнее, чем алгоритм СЭТ, но он позволяет получить полную информацию о характеристиках микросейсмического очага.
Используя метод, впервые примененный в работе [Zhang et al., 2008], нетрудно построить для алгоритмов (9), (10) их фазовые варианты, которые сохраняют высокую точность оценивания координат микросейсмических источников при отклонении характеристик среды и помех от строгих ограничений Б2, Б3 [Kushnir, Varupaev, 2017]:
(11)
$\begin{gathered} {{{{\mathbf{\hat {r}}}}}_{{pset}}}\left( {{{{{\mathbf{\dot {X}}}}}_{{m,n}}}} \right) = \\ = \mathop {\arg \max }\limits_{{\mathbf{r}} \in Q} {{\sum\limits_{j = {{j}_{1}}}^{{{j}_{2}}} {\left| {\sum\limits_{l = 1}^m {\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}}\exp \left( { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right)} \right)} } \right|} }^{2}}, \\ \end{gathered} $(12)
$\begin{gathered} {{{{\mathbf{\hat {r}}}}}_{{pmset}}}\left( {{\mathbf{\dot {X}}}_{{m,n}}^{{}}} \right) = \sum\limits_{l = 1}^m {\mathop {\arg \max }\limits_{{\mathbf{r}} \in Q;\,{\mathbf{\theta }} \in \Theta } } \sum\limits_{j = 1}^n {\left| {\sum\limits_{l = 1}^m {{{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right){{a}_{l}}{{{\left( {\mathbf{r}} \right)}}^{{^{{^{{^{{^{{}}}}}}}}}}}} } \right.} \\ \times \,\,{{\left. {\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}}\exp \left( { - 2\pi i{{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right) - \pi {{{\bar {\gamma }}}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right)} \right)} \right|}^{2}}, \\ \end{gathered} $Уточним, что фазовые алгоритмы (11) и (12) сохраняют качество оценивания параметров источника при изменении статистических характеристик помех ${{\xi }_{l}}\left( t \right),$ воздействующих на датчики группы (формула (2)) при сохранении их статистической независимости. Только это следует из модельных экспериментов, приведенных в работе [Kushnir, Varypaev, 2017]. Изучение других свойств фазовых алгоритмов микросейсмического мониторинга требует дальнейших теоретических и экспериментальных исследований.
В алгоритмах (11) и (12) обратное преобразование ДКПФ каждой из сумм по l (рассматриваемых как функции от $j \in \overline {1,n} $) определяет варианты процедуры стекинга сейсмограмм группы, подобные процедуре (10):
(13)
$\sum\limits_{l = 1}^m {\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}}\exp \left( { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right)} \right)} \mathop \Leftrightarrow \limits^{{\text{DFFT}}} {{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{s} }_{k}}\left( {\mathbf{r}} \right),\,\,k \in \overline {1,n} ,$(14)
$\begin{gathered} \sum\limits_{l = 1}^m {{{d}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right){{a}_{l}}\left( {\mathbf{r}} \right)\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}} \times } \\ \times \,\,\exp \left( { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right) - \pi {{{\bar {\gamma }}}_{l}}\left( {{\mathbf{r}},{\mathbf{\theta }}} \right)} \right)\mathop \Leftrightarrow \limits^{{\text{DFFT}}} {{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{s} }}_{{0,k}}},\,\,\,\,k \in \overline {1,n} . \\ \end{gathered} $Функции ${{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{s} }_{k}}\left( {{{{{\mathbf{\hat {r}}}}}_{{pset}}}} \right)$ и ${{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{s} }_{k}}\left( {{{{{\mathbf{\hat {r}}}}}_{{mpset}}}} \right),$ $k \in \overline {1,n} ,$ конечно, не являются, в строгом смысле, оценками временной функции источника ${{s}_{{0,k}}},$ $k \in \overline {1,n} ,$ но они сохраняют некоторые ее характерные черты.
Основная цель данной работы заключалась в создании методики обработки данных сейсмических групп для различения типов очагов микросейсмических источников по критерию: взрывной или сдвиговый тип. Для этого мы исследовали возможность использования вычислительно простой процедуры робастного к помехам стекинга (13), и вычисляли следующую функцию координат ${\mathbf{r}} \in Q{\text{:}}$
(15)
${{K}_{{pset}}}\left( {\mathbf{r}} \right) = \,{{\sum\limits_{j = {{j}_{1}}}^{{{j}_{2}}} {\left| {\sum\limits_{l = 1}^m {\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}}\exp \left( { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {\mathbf{r}} \right)} \right)} } \right|} }^{2}},\,\,{\mathbf{r}} \in Q.$Результаты моделирования и обработки реальных данных, зарегистрированных малоапертурной сейсмической группой в районе добычи железной руды, также показали, что простые алгоритмы пространственно-временного (F-K)-анализа [Capon, 1973; Кушнир, 2012], широко используемые для обработки данных сейсмических групп тоже во многих случаях также обеспечивают возможность распознавания типов микросейсмических очагов. В настоящей работе мы использовали робастную фазовую версию алгоритма F-K-анализа, устойчивую к изменениям статистических характеристик сейсмических помех ${{\xi }_{l}}\left( t \right)$ в выражении (2):
(16)
$W\left( {\mathbf{p}} \right) = \sum\limits_{j = {{j}_{1}}}^{j = {{j}_{2}}} {{{{\left| {\sum\limits_{l = 1}^m {\frac{{{{{\dot {X}}}_{{l,j}}}}}{{\left| {{{{\dot {X}}}_{{l,j}}}} \right|}}\,} {\text{exp}}\left( { - i \times 2\pi {{f}_{j}}{\mathbf{u}}_{l}^{T}{\mathbf{p}}} \right)} \right|}}^{{\text{2}}}}} ,\,\,\,\,{\mathbf{p}} \in P,$Карты F-K-анализа (16), также как и диаграммы интенсивности (15), имеют различный графический вид в зависимости от того, является ли источник взрывным или его механизм имеет сильные сдвиговые компоненты.
Процедуры оценивания координат источника и параметров тензора сейсмического момента, подобные (10) и (12), открывают возможность точного оценивания временной функции источника $s\left( {k{{\Delta }_{{smp}}}} \right),$ $k \in \overline {1,n} $ путем модификации известного фильтра Кейпона [Capon, 1969]. Этот алгоритм представляет собой многоканальный фильтр с m входами и одним выходом, который минимизирует для всех $k \in \overline {1,n} $ среднеквадратическое отклонение оценки ${{\hat {s}}_{k}}$ от истинного значения сигнала $s\left( {k{{\Delta }_{{smp}}}} \right)$ при условии несмещенности оценки сигнала. В частотной области ДКПФ эти условия могут быть записаны следующим образом:
(17)
$\begin{gathered} {{E}_{s}}\left\{ {{{{\left| {\hat {s}\left( {{{f}_{j}}} \right) - s\left( {{{f}_{j}}} \right)} \right|}}^{{\text{2}}}}} \right\} = \\ = {{E}_{s}}\left\{ {{{{\left| {{\mathbf{\dot {\varphi }}}_{С}^{*}\left( {{{f}_{j}}} \right){\mathbf{\dot {X}}}\left( {{{f}_{j}}} \right) - s\left( {{{f}_{j}}} \right)} \right|}}^{2}}} \right\} = \mathop {{\text{min}}}\limits_{{\mathbf{\dot {\varphi }}}\left( f \right)} , \\ \end{gathered} $(18)
${{E}_{s}}\left\{ {{\mathbf{\dot {\varphi }}}{\text{*}}\left( {{{f}_{j}}} \right){\mathbf{\dot {X}}}\left( {{{f}_{j}}} \right)} \right\} = s\left( f \right),$В рассматриваемом случае сложного механизма сейсмического источника фильтр Кейпона описывается следующим соотношением:
(19)
$\begin{gathered} {{{\hat {s}}}_{{0,j}}} = {\mathbf{\dot {\varphi }}}_{{c,j}}^{*}{{{{\mathbf{\dot {X}}}}}_{j}} = \frac{{{\mathbf{\dot {h}}}_{j}^{ * }\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right){\mathbf{\dot {F}}}_{{{\text{ш}},j}}^{{ - 1}}{{{{\mathbf{\dot {X}}}}}_{j}}}}{{{\mathbf{\dot {h}}}_{j}^{ * }\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right){\mathbf{\dot {F}}}_{{{\text{ш}},j}}^{{ - 1}}{{{{\mathbf{\dot {h}}}}}_{j}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right)}}, \\ j \in \overline {1,n} , \\ \end{gathered} $(20)
$\begin{gathered} {{{{\mathbf{\dot {h}}}}}_{j}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right) = \left( {{{d}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right){{a}_{l}}\left( {{{{\mathbf{r}}}_{0}}} \right) \times } \right. \\ \left. { \times \,\,\exp \left\{ { - i \times 2\pi {{f}_{j}}{{\tau }_{l}}\left( {{{{\mathbf{r}}}_{0}}} \right) - \pi {{{\bar {\gamma }}}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right)} \right\},\,\,\,\,l \in \overline {1,m} } \right), \\ \end{gathered} $В случае источника взрывного типа в выражении (20) следует положить ${{d}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right) = 1,$ ${{\bar {\gamma }}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right) = 0,$ а в случае плоских сейсмических волн от телесейсмических источников – еще и ${{a}_{l}}\left( {{{{\mathbf{r}}}_{0}}} \right) = 1.$ Заметим, что последний случай как раз и рассматривался в классической работе [Capon, 1969].
Подстановка в выражения (19), (20) вместо истинных значений параметров ${{{\mathbf{r}}}_{0}},$ ${{{\mathbf{\theta }}}_{0}}$ источника их оценок с помощью описных выше алгоритмов, позволяет получить представление о длительности и частотном составе микросейсмических источников. Это особенно просто сделать в случае источников взрывного типа. Если с помощью диаграммы интенсивности, подобной (15) или даже (16), установлено, что источник имеет взрывной механизм, то оценивание его координат ${{{\mathbf{r}}}_{0}}$ – гораздо более простая вычислительная процедура. В результате можно определить временную функцию с помощью ДКПФ–преобразования выражения (19), в котором ${{d}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right) = 1,$ ${{\bar {\gamma }}_{l}}\left( {{{{\mathbf{r}}}_{0}},{{{\mathbf{\theta }}}_{0}}} \right) = 0,$ а вместо ${{\tau }_{l}}\left( {{{{\mathbf{r}}}_{0}}} \right)$ – величины ${{\tau }_{l}}\left( {{\mathbf{\hat {r}}}} \right),$ где ${\mathbf{\hat {r}}}$ – оценка координат источника, полученная по тем же сейсмограммам группы.
3. ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ И ОБРАБОТКА РЕАЛЬНЫХ ДАННЫХ
3.1. Моделирование и обработка синтетических сейсмограмм
В этом разделе приведены результаты моделирования сейсмограмм поверхностной малоапертурной группы из 5-ти сейсмических датчиков, развернутой в окрестности шахты для добычи железной руды в районе Курской магнитной аномалии, и результаты обработки этих сейсмограмм с помощью робастной диаграммы локации (15) и робастной диаграммы F-K-анализа (16). Целью этого модельного исследования было продемонстрировать, что диаграммы стекинга и диаграммы F-K-анализа позволяют различать по сейсмограммам малоапертурной группы микросейсмические источники с взрывными и сложными сдвиговыми механизмами очагов.
Моделирование сейсмограмм группы осуществлялось в рамках ограничений А.1–А.5 Раздела 2 с использованием различных вариантов тензора сейсмического момента источника. Параметры слоистой горизонально-однородной скоростной модели среды представлены в табл. 1, расположение датчиков группы показано на рис. 1. Черным кругом на этом рисунке отмечена проекция источника на дневную поверхность. Глубина источника была выбрана равной 0.35 км. Использовавшиеся при моделировании характеристики среды, конфигурация группы и положение относительно нее контролируемых микросейсмических источников, типичны для мониторинга подземных горных работ.
Таблица 1.
Параметры слоистой горизонтально-однородной упругой среды
| Толщины слоев, км | Скорости продольных волн, км/с | Скорости поперечных, волн, км/с | Плотности, gm/cc |
|---|---|---|---|
| 0.15 | 3.3 | 1.34 | 2.11 |
| 0.18 | 3.5 | 1.44 | 2.6 |
| (Полупространство) 3.2 | 4.0 | 1.73 | 2.8 |
Рис. 1.
Конфигурация поверхностной сейсмической группы: черные треугольники – вертикальные сейсмометры группы; черная точка – проекция эпицентра сейсмического события на дневную поверхность.
Моделирование сигналов группы, порожденных источником, выполнялось на основе метода построения обобщенных сейсмических лучей, реализованного в свободном программном обеспечении (подпрограммы gprep96 и gpulse96 [Herrmann, 2013]). Моделировались сейсмограммы датчиков группы как проекции колебаний продольных сейсмических волн на нормаль к поверхности в точках расположения вертикальных сейсмометров группы. Сейсмический точечный источник описывался при моделировании тремя тензорами сейсмического момента: а) однородного расширения; б) чистого сдвига; в) сложного типа с расширениями и сдвигами:
Временная функция источника моделировалась в виде Рикер-вейвлета:
Рис. 2.
Фокальные механизмы сейсмических источников с вышеописанными тензорами ${{{\mathbf{M}}}_{1}},$ ${{{\mathbf{M}}}_{2}},$ ${{{\mathbf{M}}}_{3}}{\text{:}}$ интенсивности P-волн в дальней зоне и их стереографические проекции на нижнюю полусферу.
Модельные сейсмические сигналы, порожденные в разных датчиках группы источниками с тензорами ${{{\mathbf{M}}}_{1}}$, ${{{\mathbf{M}}}_{2}}$ и ${{{\mathbf{M}}}_{3}}$ и временной функцией $s\left( t \right),$ показаны на рис. 3.
Рис. 3.
Модельные сейсмограммы датчиков группы длительностью 0.5 с для тензоров сейсмического момента источника ${{{\mathbf{M}}}_{1}},$ ${{{\mathbf{M}}}_{2}}$ и ${{{\mathbf{M}}}_{3}}.$
Сейсмические сигналы рис. 3. были обработаны с помощью алгоритмов (15) и (16) во временном окне (0.08–0.2 с), охватывающем только сейсмограммы P-волн на всех датчиках группы. При данных параметрах моделирования это было возможно, поскольку достаточно малы относительные временные сдвиги P-волн на сейсмограммах всех датчиков, определяемые скоростной моделью среды (табл. 1) и взаимным расположением сейсмического источника и датчиков группы.
На рис. 4 представлены: в первом ряду – диаграммы пространственной интенсивности сейсмических волн (15) как функции горизонтальных координат источника, во втором ряду – диаграммы F-K-анализа (16) как функции горизонтальных медленностей.
Рис. 4.
Первый ряд – диаграммы локации для различных фокальных механизмов. Второй ряд – диаграммы F-K-анализа для тех же источников в зависимости от кажущихся медленностей волн: (а) – механизм равномерного расширения (взрыва) с тензором ${{{\mathbf{M}}}_{1}};$ (б), (в) – сдвиговые механизмы с тензорами ${{{\mathbf{M}}}_{2}},$ ${{{\mathbf{M}}}_{3}}.$
В случае механизма источника типа равномерного расширения (взрыва) диаграмма излучения источника свидетельствует, что очаг сейсмической эмиссии локализован в некоторой точке пространства, а диаграмма F-K-анализа – что сейсмическая волна приходит с определенного направления. Поскольку положение источника при моделировании было выбрано довольно близким по отношению к сейсмической группе, кажущаяся скорость сейсмической волны должна быть довольно высокой, что и показывает F-K-диаграмма рис. 4а, где эта скорость равна 17.24 км/с.
Как показывают диаграммы локации рис. 4б, 4в, в случае сложных механизмов источника невозможно достаточно точно оценить координаты источника, так же как и оценить вектор кажущихся медленностей P-волны от источника по диаграммам F-K-анализа. Как обсуждалось выше, это связано с различными полярностями (знаками вступления) P-волн от источника, наблюдаемыми на различных датчиках сейсмической группы.
Полученные результаты хорошо согласуются с результатами работ [Кушнир, Варыпаев и др., 2014; Eisner et al., 2014], несмотря на то, что в указанных работах моделирование производилось для сейсмических групп со значительно большей апертурой и большим числом датчиков. Однако эти результаты показывают, что симметричное расположение экстремумов в случае диаграмм излучения сложных источников, отмечавшееся в указанных работах, необязательно проявляется при малых расстояниях от источника до датчиков группы и/или при любых конфигурациях группы.
Таким образом, из результатов моделирования следует вывод, что и при малой апертуре сейсмической группы, малом числом датчиков в ней и близком расположении источника к датчикам группы диаграммы кажущегося излучения и F-K-анализа позволяют надежно отличать источники с простым взрывным механизмом от источников со сложными сдвиговыми механизмами. Компьютерный эксперимент также подтверждает вывод, ранее сделанный в указанных выше работах, что наличие сложных механизмов сейсмических источников не позволяет использовать стандартный алгоритм сейсмической эмиссионной томографии (СЭТ) для локализации источников. Можно также утверждать, что статистку СЭТ нецелесообразно использовать для обнаружения микросейсмических источников со сложными механизмами. Это связано с тем, что сложные механизмы источников значительно уменьшают максимальное значение статистики СЭТ и тем самым уменьшают вероятность обнаружения подобных источников в условиях сильных сейсмических помех.
3.2. Регистрация и обработка реальных сейсмограмм, зарегистрированных датчиками малоапертурой группы
Сейсмические измерения производились с 28 по 29.10.2017 г. с помощью малоапертурной группы сейсмических датчиков, развернутой в районе Курской магнитной аномалии в непосредственной близости от шахт для добычи железной руды. В качестве датчиков группы использовались сейсмометры Reftek-130 с частотой дискретизации 500 Гц, синхронизацию которых обеспечивали приемники спутниковой системы GPS. На рис. 1 показано расположение датчиков сейсмической группы с апертурой около 300 м. Краткая информация о параметрах датчиков сейсмической группы представлена в табл. 2, где отражены только вертикальные датчики группы.
Таблица 2.
Параметры сейсмических датчиков
| Станция | Широта | Долгота | Высота над уровнем моря, км | Тип сейсмометра |
|---|---|---|---|---|
| “Север” | 51:18:28.14 | 37:33:51.13 | 0.172 | SM-3KV (вертикальный) |
| “Центр” | 51:18:25.30 | 37:33:50.47 | 0.17 | SM-3KV (вертикальный) |
| “Запад” | 51:18:23.74 | 37:33:46.39 | 0.172 | SM-3KV (вертикальный) |
| “Восток” | 51:18:23.69 | 37:33:54.60 | 0.171 | SM-3KV (вертикальный) |
| “Бочка” | 51:18:22.57 | 37:33:38.37 | 0.179 | SPV-3K |
Цифровая регистрация сейсмических данных группы производилась непрерывно в течение 20 часов. Основной целью этих наблюдений являлось исследование микросейсмической активности земной среды после сильных подземных взрывов, которые производятся в шахтах для добычи железной руды. То есть, необходимо было зарегистрировать сейсмические волны от микросейсмических источников, маскируемых техногенными сейсмическими помехами, которые достаточно интенсивны в этом индустриальном регионе, для того чтобы исследовать характеристики микроземлетрясений, возбуждаемых сильными взрывами в шахтах, которые приводят к изменению напряженного состояния земной среды.
3.3. Определение параметров детектора сигналов микросейсмических источников с помощью численного моделирования
Как указывалось в Разделе 2, для обнаружения многоканальных сигналов, возбуждаемых в датчиках группы слабыми микросейсмическими источниками и маскируемых сейсмическими помехами, целесообразно использовать алгоритм детектирования (6). Такой детектор необходим в рассматриваемом случае, поскольку отсутствовала какая-либо информация о скоростных характеристиках среды распространения сейсмических волн от микросейсмических источников до датчиков группы.
В соответствии с методикой обнаружения сейсмических сигналов, основанной на общей теории проверки статистических гипотез, необходимо установить пороговое значение статистики (6) – функции от наблюдений, обеспечивающее вероятность ложной тревоги достаточно близкой к нулю, и оценить, какова для такого порога вероятность правильного обнаружения при заданном отношении сигнал–шум в наблюдаемых сейсмограммах. Для этой цели мы провели моделирование по методу Монте-Карло, которое заключалось в обработке большого количества смесей синтетического сигнала микросейсмического источника с последовательными фрагментами 120‑секундной записи реальных “чистых” помех. Подобное моделирование использовалось ранее в работах [Kushnir et al., 2014; Kushnir, Varypaev, 2017] для сравнения качества оценивания координат микросейсмических источников с помощью различных алгоритмов при разных соотношениях сигнал/помеха и разных статистических свойствах сейсмических помех. Процедуры создания синтетических смесей сигналов с помехами и определения вероятностей ошибок по методу Монте-Карло, используемые при этом моделировании, подробно описаны в указанных работах. Ниже укажем только формулу, которая использовалась для вычисления усредненного по каналам группы отношения сигнал/шум (ASNR) в синтетических смесях сигналов с помехами:
На рис. 5 показаны две синтетические многоканальные сейсмограммы группы, содержащие каждая по 51 смесей с различными значениями ASNR в полосе частот (10 -30) Гц: а) – 4; б) – 0.5. Число смесей в каждой многоканальной сейсмограмме равно 51. Расстояние во времени между смесями равно 2 с. Отметим, что на сейсмограммах рис. 5б, где ASNR смесей равно 0.5, трудно визуально распознать временные интервалы, в которых присутствуют сигналы источника.
Рис. 5.
Синтетические смеси коротких сигналов источника (с длительностью 0.4 с) и реальных помех, воздействующих на датчики группы. Число смесей – 51. Расстояние во времени между смесями – 2 с. (а) – ASNR = 4; (б) – ASNR = 0.5. Значения ASNR заданы в полосе частот (10–30) Гц.
Основная идея моделирования по методу Монте-Карло с целью определения вероятностных характеристик детектора состоит в том, чтобы вычислять статистику (6) и выполнять ее сравнение с некоторым порогом достаточно много раз, используя статистически независимые смеси сигналов с помехами, имеющими заданное ASNR. Вероятность правильного обнаружения далее оценивается как число превышений порога в процедуре (6), деленное на общее количество операций. Естественно, количество таких операций равно количеству имеющихся смесей.
Та же стратегия применяется для оценивания вероятности ложной тревоги, но вместо смесей сигналов с помехами берутся интервалы “чистых” помех. Вероятность ложного срабатывания детектора оценивается как количество превышений порога, деленное на общее количество обработанных временных интервалов “чистых” помех.
В данной работе расчеты статистики (6) выполнялись 51 раз для синтетических смесей модельного сигнала с реальными помехами при ASNR = 0.5. и 51 раз – для интервалов “чистых помех” в промежутках между интервалами смесей. Для наглядности также расчеты значений статистики (6) производились в скользящем временном окне с шагом 0.1 с и длительностью 0.4 с. При расчетах учитывался только интервал частот (10–30) Гц, поскольку ожидаемые спектры сигналов микросейсмических источников преимущественно находятся в этом диапазоне.
Результаты расчетов статистики (6) представлены на рис. 6.
Рис. 6.
Значения статистики детектирования (6) в скользящем временном окне с длительностью 0.4 с и с шагом 0.1 с. Горизонтальными линиями на рисунках отмечены значения порога обнаружения в масштабах вертикальных осей рисунков: (а) – Смеси синтетических сигналов источника с реальными помехами. Число смесей равно 51, расстояния во времени между смесями равно 2 с, ASNR = 0.5; (б) запись реальных “чистых” помех.
Из рис. 6б следует, что для выбранного порогового значения статистики (отмеченного горизонтальной чертой) вероятность ложного срабатывания очень близка к нулю, поскольку на временном интервале записи “чистых помех”, равном 120 с, нет пересечений значений статистики (6) с порогом. Из рис. 4а (2 ряд) следует, что вероятность правильного обнаружения близка 0.94, поскольку в 48 интервалах из 51-го интервала, где сигнал источника присутствовал, значение статистики (6) превышало порог обнаружения.
Отметим, что при моделировании использовалась статистика обнаружения (6), “реагировавшая” только на степень взаимной когерентности сейсмограмм датчиков группы, содержащих сигналы источника, которая отражалось в близости матричной спектральной плотности этих сейсмограмм к сингулярной матричной функции. Указанное свойство статистики обнаружения не зависит от положения источника относительно датчиков сейсмической группы. Для эффективной работы детектора существенно лишь то, чтобы расстояние от сейсмического источника до датчиков группы позволяло в одном и том же достаточно длительном временном окне “захватывать” на всех датчиках группы сигналы только от P-волны источника.
Итогом описанного численного моделирования является определение таких параметров детектора, как:
а) пороговое значение детекторной статистики, обеспечивающее необходимую вероятность ложного обнаружения;
б) вероятность пропуска сигналов источников, получаемая при этом пороге и заданном значении ASNR.
Ниже в Разделе 3.4 приводятся результаты анализа микросейсмичности в районе шахтных выработок, в котором использовался детектор (6) с порогом обнаружения ${{K}_{{\alpha }}}$ = 7, который, согласно проведенному модельному исследованию, обеспечивает вероятность ложной тревоги, меньшую 0.1, и вероятность правильного обнаружения, превышающую 0.9 при ASNR > 0.5.
3.4. Обнаружение микросейсмических источников и анализ типов их очагов по записям сейсмической группы в районе горнорудных работ на Курской магнитной аномалии
Приведенный в Разделе 2 и исследованный в Разделе 3.3 алгоритм обнаружения (6) был использован для анализа микросейсмичности в районе горнорудных работ на основе обработки записей сейсмической группы. Было обработано двадцать часов записей, в которых были выделены интервалы, на которых многоканальные записи группы имели сильную когерентность. Поэтому записи группы на этих интервалах могли соответствовать локализованным микросейсмическим событиям в окрестности группы. На двух интервалах были зафиксированы сигналы, вызванные промышленными взрывами в шахтах. Моменты этих взрывов были известны. Сейсмограммы взрывов и моменты времени их обнаружения детектором представлены на рис. 7.
Расстояние от гипоцентров обоих взрывов до группы составляло около 3 км, и направления от центра группы на эпицентры этих взрывов были известны. Мы провели F-K-анализ многоканальных записей группы с помощью робастного алгоритма (16), чтобы оценить кажущиеся медленности P-волн от взрывов. Это позволило проверить качество наших данных, поскольку вектор медленностей сейсмических волн позволяет определить направление на очаг сейсмического события. На F-K-диаграммах рис. 8 приведены результаты оценок обратных азимутов и кажущихся медленностей P-волн от обоих взрывов, которые хорошо согласуются с теоретическими значениями.
Отметим, что из-за некоторых инструментальных проблем, которые в течение короткого времени привели к сбою в регистрации данных группы, мы не смогли обработать сейсмограммы третьего взрыва, эпицентр которого был расположен вблизи от сейсмической группы.
С помощью F-K-анализа данных с использованием робастного алгоритма (15) сейсмические события на всех интервалах высокой когерентности записей группы были разделены на 3 группы: А, Б, В. На рис. 9 показаны диаграммы F-K-анализа, соответствующие каждой из групп: столбцы рис. 9 содержат по две диаграммы из каждой группы. Первый столбец соответствует группе Б, в которую мы отнесли события с высоким значением кажущейся скорости волны (>10 км/с). Такая кажущаяся скорость означает, что относительные задержки P-волн, вызванных сейсмическими событиями, довольно малы (в пересчете на плоский волновой фронт). Следовательно, источники этих событий расположены недалеко от сейсмической группы – в пределах 3-х ее апертур. С большой вероятностью эти события имеют механизм очага типа всестороннего расширения.
Рис. 9.
F-K-диаграммы временных интервалов записей сейсмической группы, обладающих высокой когерентностью. Диаграммы в столбцах (а), (б) и (в) соответствуют группам событий А, Б, В.
События из группы В невозможно идентифицировать, поскольку на их F-K-диаграммах слишком много областей с приблизительно одинаковой энергией волн (рис. 9в).
Существует простой критерий, который позволяет различать по F-K-диаграммам, порождаются ли наблюдаемые сейсмические волны некогерентными шумовыми источниками или они порождаются, в основном, одним точечным источником (хотя и искаженым некогерентными шумовыми волнами). Таким критерием является отношение величин двух самых высоких максимумов F-K-диаграммы. Для событий из группы В это отношение очень близко к 1, в отличие от событий группы А и группы Б, для которых это отношение много больше единицы. Типичные для событий из групп А и Б изображены на рис. 9а и 9б. Отметим, что F-K-диаграммы рис. 9б типичны для событий, которые могут быть идентифицированы как точечные источники со сложными сдвиговыми тензорами сейсмического момента, но для подтверждения этого факта требуется дополнительная информация. Подобные F-K-диаграммы мы получали при анализе синтетических сейсмограмм от точечных источников со сложными сдвиговыми тензорами сейсмических моментов (рис. 4), но теоретически не исключена возможность, что они могут порождаться и другими сейсмическими источниками.
Симметричное поведение главных максимумов F-K-диаграмм рис. 9б для точечных источников теоретически объясняется тем, что колебания сейсмограмм различных датчиков группы при таких источниках являются когерентными, но имеют разные полярности (разные знаки вступлений). Это и приводит к усложнению F-K-диаграммы по сравнению с простым точечным источником взрывного типа.
Чтобы проверить, соответствуют ли колебания наблюдаемых сейсмограмм из групп Б этому теоретическому выводу, мы многократно повторили следующую процедуру: изменяли полярности сейсмограмм датчиков группы с +1 на –1 в разных сочетаниях, и рассчитывали F-K-диаграммы этих модифицированных сейсмограмм. Эта процедура повторялась для всевозможных сочетаний –1 и +1 в множествах из 5 таких чисел. Если некоторая “оптимальная” F-K-диаграмма, у которой глобальный максимум имеет самое высокое значение, обладает только одним выраженным максимумом, то это означает, что анализируемые сейсмограммы датчиков группы действительно порождены точечным источником.
На рис. 10 показаны F-K-диаграммы, полученные в результате “оптимального” распределения полярностей P-волны в исходных сейсмограммах датчиков группы от двух источников из группы Б (для которых F-K-диаграммы представлены на рис. 9б). Важно отметить, что максимальные значения диаграмм рис. 10 больше, чем те же значения на диаграммах рис. 9б. Это означает, что мы определили “правильные координаты” анализируемых сейсмических источников из группы Б в пространстве векторов медленности.
Рис. 10.
Диаграммы F-K-анализа сейсмограмм источников из группы Б при “оптимальном” распределении полярностей P-волны на датчиках группы.
“Оптимальная” F-K-диаграмма в множестве F-K-диаграмм, рассчитанных для различных комбинаций полярностей наблюденных сейсмограмм группы, позволяет также определить истинное распределение знаков вступления P-волны от анализируемого источника на различные датчики группы. Действительно, зная последовательность распределений чисел +1 и –1 в множестве из 5 таких чисел, которая привела к “оптимальной” F-K-диаграмме, можно получить распределение этих чисел в исходных сейсмограммах. Заметим, что это распределение знаков вступлений P-волны на датчики группы практически невозможно определить непосредственно по наблюденным сейсмограммам из-за малых отношений сигнал–помеха в этих сейсмограммах. Таким образом, используя достаточно простые манипуляции с сейсмограммами группы от анализируемого микросейсмического источника, то есть, изменяя полярности этих сейсмограмм (умножая некоторые из них на –1), можно выявить “истинное” распределение знаков вступления слабой P‑волны на различные датчики группы. А это распределение физически согласуется с локальным напряжениями и сжатиями земной среды, определяемыми тензором сейсмического момента источника.
Отсюда следует вывод, что события из группы Б (представленная F-K-диаграммами рис. 9б), так же как события из группы А (рис. 9а), могут описываться как события с точечными источниками, но имеют сложные тензоры сейсмического момента, и, следовательно, излучают одновременно и продольные, и поперечные сейсмические волны во взаимно перпендикулярных направлениях (направления также определяются этими тензорами).
Отметим, что на основании анализа характеристик P-волн можно доказать наличие сложного механизма источника не во всех случаях. Это зависит от направления вектора сдвиговых напряжений в источнике и, следовательно, от диаграммы волнового излучения источника. Если направление этого вектора почти горизонтальное и источник находится вблизи группы, то когерентные сейсмограммы датчиков будут иметь одинаковые полярности, и на основе анализа этих сейсмограмм невозможно определить, имеет ли источник сдвиговый или взрывной механизм. Отсюда также следует, что невозможно однозначно утверждать, что события группы А имеют механизм очага взрывного типа, поскольку ориентация диаграммы направленности и сдвигового источника может быть такой, что на всех сейсмических датчиках группы будет наблюдаться одинаковая полярность когерентных колебаний P-волны. Также трудно категорично утверждать, что все выявленные слабые события из групп А и Б находятся недалеко от апертуры группы только на основании, оценка кажущейся скорости сейсмической волны превышает 10 км/с.
Важно отметить, что сходство F-K-диаграмм на рис. 9а, 9б и рис. 4 можно интерпретировать как физическую идентичность синтетических сейсмограмм группы, вызванных сейсмическими источниками со сложными механизмами очагов и наблюдаемых сейсмограмм от реальных микросейсмических событий. То есть, все атрибуты описанного в Разделе 3.1 моделирования, а именно: относительное расположение сейсмических источников и датчиков сейсмической группы, скоростная модель упругой среды, продолжительность времени импульса источника, сложность механизмов очагов сейсмических источников в определенной степени близки к атрибутам реального сейсмического мониторинга.
Заметим также, что большинство сейсмических событий, выявленных с помощью алгоритма детектирования, было обнаружено в течение двух с половиной часов после шахтных промышленных взрывов. Это соответствует распространенным представлениям о процессах наведенной микросейсмической активности, вызванных сильными воздействиями на земные породы.
Аналогичные результаты обнаружения сейсмических событий были получены также в ходе проведения измерительного эксперимента 20.10.2018 г. Это был в точности повторный эксперимент с использованием 4 датчиков группы рис. 1, за исключением датчика “Бочка”. Ниже, на рис. 11 приведен пример пары F-K-диаграмм, соответствующих очагам со сложным механизмом. На рис. 12 показаны F-K-диаграммы, рассчитанные для тех же временных интервалов сейсмограмм, что и на рис. 11, а также полученные в результате “оптимального” распределения полярностей P-волны в исходных сейсмограммах датчиков группы.
Рис. 11.
F-K-диаграммы временных интервалов записей сейсмической группы, обладающих высокой когерентностью.
Рис. 12.
Диаграммы F-K-анализа сейсмограмм источников из группы Б при “оптимальном” распределении полярностей P-волны на датчиках группы, соответствующие временным интервалам рис. 11.
Очень важно отметить, что вид F-K-диаграмм рис. 11 соответствует общим представлениям о симметрии главных экстремумов относительно истинного значения вектора кажущейся медленности, поэтому события, сигналы которых были обнаружены на данных временных промежутках относятся к группе Б рис. 9. Однако ввиду того, что в данном эксперименте было использовано всего 4 датчика, то F-K-диаграммы рис. 11 отличны от F-K-диаграмм рис. 9б, а также от рис. 4. Повторный эксперимент позволил также убедиться в устойчивости полученных ранее результатов по детектированию и распознаванию временных интервалов, содержащих многомерный сигнал сейсмического события.
3.5. Оценивание временной функции P-волны источника с помощью фильтра Кейпона
Выше мы оценили направления прихода (векторы кажущихся медленностей) P-волн от двух событий из группы А и группы Б (рис. 9а, 9б). Это позволило нам произвести еще одну проверку, что эти события действительно были вызваны точечными микросейсмическими источниками. С этой целью мы применили алгоритм оценивания временной формы сейсмической волны на основе сейсмограмм события, зарегистрированными датчиками группы. Алгоритм оценивания описан в Разделе 2 как фильтр Кейпона и определяется выражением (9). Чтобы исключить влияние низкочастотного сейсмического шума, мы применили фильтра Кейпона в диапазоне частот от 15 до 30 Гц. Были обработаны события из группы Б, F-K-диаграммы которых изображены на рис. 9а, 9б. Оценки временных функций сейсмических волн, зарегистрированных для этих событий датчиками группы, показаны на рис. 13. Как обсуждалось ранее, в случае слабого сейсмического события отсутствует информация о полярностях P-волны, воздействующей на разные датчики группы, поэтому временная функция волны может быть восстановлена с точностью до знака. Колебания на графиках рис. 13, ограниченные красными линиями, и представляют оценки временных форм P‑волн от двух анализируемых микросейсмических источников. Эти колебания имеют длительность менее 0.2 с и имеют похожую (с точностью до полярности) форму. Эти симметричные функции времени похожи на функцию Рикер-вейвлета [Bording, 1996], которая широко используется при моделировании синтетических сейсмограмм от микросейсмических источников. Отметим, что оценки временных функций сейсмических источников, представленные на рис. 13, близки к типичной форме многократно регистрировавшихся сейсмических колебаний от слабых источников в полосе частот (15–30) Гц.
4. ЗАКЛЮЧЕНИЕ
Основным результатом данной работы является разработка методики многоканальной обработки данных сейсмических групп при микросейсмическом мониторинге и апробирование этой методики как на реальных, так и на модельных данных малоапертурных сейсмических групп. Результаты численного моделирования, проделанного перед обработкой реальных данных, оказались весьма сходными с результатами обработки сейсмограмм группы, зарегистрированными при мониторинге микросейсмической активности, вызванной шахтными промышленными взрывами.
Разработанная методика микросейсмического мониторинга позволяет распознавать и анализировать когерентные сейсмические сигналы от микросейсмических источников, скрытые естественным сейсмическим шумом при отношении сигнал/помеха, существенно меньшем единицы, т.е. когда энергия P-волн от микросейсмического источника намного меньше энергии сейсмических помех. Особенности разработанной методики микросейсмического мониторинга состоят в следующем:
1) используется алгоритм обнаружения сейсмических сигналов, основанный только на анализе когерентности сейсмограмм датчиков группы, что позволяет выделять сейсмические сигналы микросейсмических источников в условиях сильных сейсмических помех искусственного происхождения, сейсмограммы которых, как правило, не являются когерентными;
2) используются фазовые робастные алгоритмы анализа сигналов от микросейсмических источников для определения координат и характеристик их очагов, качество которых слабо зависит от статистических характеристик промышленных сейсмических помех.
Результаты обработки реальных данных показывают, что обнаруженные детектором многоканальные сигналы группы, скорее всего, порождены микросейсмическими событиями, т.е. не имеют техногенного происхождения по следующим причинам:
– согласно принципу работы детектора эти сигналы имеют высокую когерентность;
– пространственный спектральный анализ этих сигналов показывает, что зарегистрированные группой сейсмические волны приходят на группу с определенного направления;
– спектры этих сигналов лежат, в основном, в полосе частот (10–30) Гц;
– временные функции сейсмических волн, порождающих эти сигналы, имеют длительность порядка 0.2 с;
– в ряде случаев полярности этих когерентных сигналов разные и соответствуют сдвиговым механизмам точечного сейсмического источника, что не типично для сейсмических помех.
Авторы настоящей работы полагают, что полученные результаты могут способствовать усовершенствованию существующей технологии микросейсмического мониторинга с помощью малоапертурных сейсмических групп, что позволит улучшить надежность обнаружения слабых событий.
Список литературы
Аки К., Ричардс П. Количественная сейсмология. Т. 1, 2. М.: МИР. 1983.
Кушнир А.Ф. Статистические и вычислительные методы сейсмического мониторинга. М.: КРАСАНД. 2012. 464 с. ISBN: 978-5-396-00450-4.
Кушнир А.Ф, Варыпаев A.В., Рожков М.В., Епифанский А.Г., Дрикер И. Определение параметров очагов микросейсмических событий по данным поверхностных сейсмических групп при сильных коррелированных помехах и сложных механизмах источников излучения // Физика Земли. 2014. № 3. С. 28–50. УДК: 550.34.016.
Варыпаев А.В., Санина И.А., Чулков А.Б., Кушнир А.Ф. Применение робастных фазовых алгоритмов для выявления сейсмической эмиссии в районе проведения взрывных работ в шахтах // Сейсмические приборы. 2018. Т. 54. № 2. С. 5–18. https://doi.org/10.21455/si2018.2-2
Anikiev D., Valenta J., Staněk F., Eisner L. Joint location and source mechanism inversion of microseismic events: benchmarking on seismicity induced by hydraulic fracturing // Geophysical J. International. 2014. V. 198. Is. 1, July. P. 249–258.
Bendat S., Piersol G. Random Data: Analysis and Measurement Procedures, 4th Edition. Wiley Press. 2010.
Bording P. SeismicWave Propagation – Modeling and Inversion. Society of Exploration Geophysics, Tulsa, Oklahoma. 1996.
Capon J. High-resolution frequency-wavenumber spectrum analysis. Proceedings of the IEEE 57. 1969. P. 1408–1418.
Cros E., Roux P., Vandemeulebrouck J., Kedar S. Locating hydrothermal acoustic sources at Old Faithful Geyser using Matched Field Processing // Geophys. J. Int. 2017. V. 187(1). P. 385–393.
Droznin D., Shapiro N., Droznina S., Senyukov S., Chebrov V., Gordeev E. Detecting and locating volcanic tremors on the Klyuchevskoy group of volcanoes (Kamchatka) based on correlations of continuous seismic records // Geophys. J. Int. 2015. V. 203. P. 1001–1010.
Duncan P.M., Lakings J.D., Flores R.A. Method for passive seismic emission tomography. US Patent #7,663,970. 2010.
Eisner L., Williams-Stroud S., Hill A., Duncan P., Thornton M. Beyond the dots in the box: Microseismicity-constrained fracture models for reservoir simulation. The Leading Edge 29. 2010. P. 326–333.
Herrmann R.B. Computer programs in seismology: an evolving tool for instruction and research // Seismol. Res. Lett. 2013. V. 84. P. 1081–1088.
Hannan E.J. Multiple Time Series Analysis. John Willey and Sons, Inc., N.Y. 1970.
Jollife I. Principal component analyses. Springer-Verlag. 1986.
Kiselevitch V.L., Nikolaev A.V., Troitskiy P.A., Shubik B.M. Emission tomography: main ideas, results, and prospects. Proceedings of the 61st annual international meeting, SEG, expanded abstracts, 1602. 1991.
Krim Hamid, Viberg Mats. Sensor Array Signal Processing: Two Decades Later. 1995.
Kushnir A.F. Varypaev A.V. Robustness of statistical algorithms for location of microseismic sources based on surface array data. Comp. Geoscience. 2017. https://doi.org/10.1007/s10596-017-9623-6
Kushnir A., Varypaev A., Dricker I., Rozhkov M., Rozhkov N. Passive surface microseismic monitoring as a statistical problem: location of weak microseismic signals in the presence of strongly correlated noise // Geophysical Prospecting. 2014. V. 62(4). P. 819–833.
Lyubushin A. Global seismic noise synchronization and seismic danger. Second European conference on earthquake engineering and seismology, Istanbul, Aug. 2014. P. 25–29.
Malovichko D.A., Lynch R.A. Micro-seismic monitoring of open-pit slopes // Mining Echo. 2006. V. 24(2). P. 21–30.
Maochen Ge. Efficient mine microseismic monitoring // Internat. J. Coal Geol. 2005. V. 64. P. 44–56.
Marple S. Digital Spectral Analysis. Englewood Cliffs. N.J.: Prentice-Hall. 1985.
Peter Stoica, Moses R. Spectral Analysis of Signals. N.J.: Prentice Hall. 2005.
Samson J.C., Olson J.V. Data – adaptive polarization filter for multichannel geophysical data // Geophysics. 1982. V. 46. № 10. P. 1423–1431.
Singh Hema, Jha, Rakesh Mohan. Trends in Adaptive Array Processing. 2012.
Stanek Frantisek, Valenta J., Anikiev Denis, Eisner L. Semblance for Microseismic Event Detection. 2014.
Zhang C., Florêncio D., Ba D.E., Zhang Z. Maximum likelihood sound source localization and beam forming for directional microphone arrays in distributed meetings // IEEE Trans. Multimedia. 2008. V. 10(3). P. 538–548.
Дополнительные материалы отсутствуют.
Пн.-пт.: 10.00-19.00