СПОСОБ МОНИТОРИНГА РЕГИОНА СЕТЬЮ СЕЙСМОСТАНЦИЙ
METHOD OF MONITORING REGION WITH NETWORK OF SEISMIC STATIONS

Изобретение относится к области сейсмологии и может найти применение в национальных системах наблюдения и обработки данных геофизических измерений для прогнозирования землетрясений.
2
4
Достоверный прогноз землетрясений возможен при установлении их первопричины. Одной из первопричин, вызывающих землетрясения, является накачка земной коры дополнительной энергией восходящего потока газов (Н, Не, СН). Факт эманации газов из земной коры накануне сейсмического удара [см., например, Патенты RU №2204852, 2003 г.; №2275659, 2006 г.; №2302020, 2007 г.]. Накачка земной коры дополнительной энергией приводит к раскачке очага землетрясений, сопровождаемой распространением от гипоцентра сверхнизких литосферных волн [см. Научное открытие №365. «Явление раскачки очага землетрясений перед сейсмическим ударом», РАЕН, 2008 г.].
Известен «Способ предсказания землетрясений». Патент RU №2170446, 2001 г., путем измерения сверхнизких литосферных волн раскачки очага землетрясений космическими средствами системы Navstar (GPS) - аналог.
i
i
i
В способе-аналоге размещают в сейсмическом районе приемные станции космической навигационной системы, разнесенные на протяженной измерительной базе, осуществляют непрерывное высокоточное измерение координат точек размещения приемных станций, регистрируют момент появления периодических отклонений Δx, Δy, Δz координат точек и отслеживают изменения этих отклонений во времени, вычисляют гипотетический центр очага как точку пересечения радиус-векторов в пространстве, длину и направляющие косинусов, которые определяют из соотношений
рассчитывают время удара, отсчитываемое от момента появления периодических отклонений координат:
и магнитуду
где Т - период отклонения координат, час, dek r - натуральный логарифм отношения амплитуд отклонения координат двух смежных периодов; d, I - коэффициенты регрессии.
Недостатками способа-аналога являются:
- наличие скрытого участка нечувствительности (мертвой зоны), ограниченной среднеквадратичной ошибкой средств GPS, снижающей интервал прогнозирования;
- неточность регрессионных зависимостей расчета времени удара и магнитуды ввиду зависимости периода (Т) отклонения координат от времени.
Ближайшим аналогом к заявленном способу является «Способ краткосрочного предсказания землетрясений». Патент RU №2181205, G01V, 9/00, 2002 г.
max
В способе ближайшего аналога осуществляют регистрацию сейсмического фона в виде непрерывной последовательности дискретных отсчетов амплитуды сигнала Δ(t) в двух разнесенных по координатам пунктах, рассчитывают спектр Фурье от последовательности выборок измерений с объемом отчетов в каждой выборке N≥2F/σ, вычисляют автокорреляционные функции В(τ) сигналов выборок и определяют интервал корреляции τ, регистрируют начало изменения параметра τ и при его непрерывном отслеживании фиксируют время запаздывания Δτ изменения фазы данного признака между двумя пунктами, рассчитывают направляющий косинус траверзы прихода сверхнизких волн очага
определяют гипотетический центр очага как точку пересечения на траверзе радиус-векторов пунктов с косинусом угла при вершине:
0
x
0
mах
1
2
вычисляют период Т параметра τ и по его величине прогнозируют магнитуду  (час) и время удара t=2,3Т, где F - максимальная частота спектра сейсмического фона, Гц; σ - среднеквадратическая ошибка вычисления спектра Фурье по дискретной выборке измерений; α - длина базы между двумя пунктами, м; B(0), В(0) -значения автокорреляционных функций в нуле для каждого пункта; V -скорость сейсмических волн в земной коре, м/с.
К недостаткам ближайшего аналога следует отнести:
- интервал автокорреляции сигнала τ определяется в первую очередь шириной спектра сигнала сейсмического фона, поэтому измерение параметра Δτ технически трудно реализуемо, поскольку оно практически не зависит от размера базы;
0
0
0
- неточность регрессионных зависимостей определения магнитуды и времени удара, в частности, известна из соотношений Гутенберга - Рихтера, что чем больше время Т, тем магнитуда ожидаемого удара должна быть больше, к тому же расчетный период Т не является постоянным, а является функцией времени T(t) раскачки очага землетрясения.
Задача, решаемая заявленным способом, состоит в расширении интервала времени упреждающего прогноза, повышении точности и достоверности определения прогнозируемых параметров сейсмического удара существующими наземными средствами путем выделения скрытой информации из соотношения скоростей продольных и поперечных волн на выбранных трассах их распространения.
p
s
p
s
p
s
р
s
p
s
i
Технический результат достигается тем, что способ мониторинга региона сетью сейсмостанций включает трассовые измерения соотношения скоростей распространения продольных V и поперечных V волн в земной коре (V/V) под различными азимутами трасс каждой из станций сети, использующих в качестве первичных датчиков скоростей датчики давлений, размещенные во взаимно ортогональных плоскостях измерений V и V, получение дискретных выборок измерений (V/V)для каждой из трасс в одни и те же интервалы времени [V/V (t)], расчет спектров Фурье полученных выборок с определением средневзвешенных частот выборок, отождествление момента изменения средневзвешенных частот с возникновением аномальной зоны подготавливаемого землетрясения на измеряемой трассе, пеленгацию аномальной зоны как точки пересечения радиус-векторов скоростей нескольких станций на измеряемых трассах, определение энергии сигнала расчетных спектров как
0
уст
y
и отслеживание динамики ее изменения в виде экспоненты E(t)=E * exp(-t/T), расчет размеров аномальной зоны R, км ≂ 1/E прогнозирование магнитуды ожидаемого сейсмического удара М≂ln[R,км] и времени удара, отсчитываемого от момента изменения средневзвешенной частоты спектра Фурье: t~4,7T, где  - постоянная времени затухания сигнала на измеряемой трассе;
1
2
Δt - интервал времени между расчетными определениями энергии сигнала E и E;
0
E - энергия сигнала до момента возникновения аномалии;
уст
0
-4,7
E≅E * е - установившееся значение энергии сигнала, прошедшего через аномалию накануне сейсмического удара.
Изобретение поясняется чертежами, где
фиг.1 - схема сети сейсмостанций Камчатского геофизического полигона, на котором отрабатывался заявленный способ и пеленгация аномалии как точки пересечения радиус - векторов сейсмостанций;
фиг.2 - реализации функций Vp / Vs (t) для отдельных трасс в одних и тех же интервалах измерений;
фиг.3 - спектры Фурье выборок измерений отслеживаемых функций;
фиг.4 - динамика изменения энергии сигналов на измеряемых трассах;
фиг.5 - функциональная схема устройства, реализующего способ;
фиг.6 - сводный график результатов мониторинга.
p
s
Техническая сущность изобретения заключается в следующем. Достоверное предсказание могут обеспечить те методы, которые основаны на измерении предвестников первопричины землетрясений. Одной из первопричин землетрясений является дегазация земной коры и ее накачка дополнительной энергией восходящего потока газов высокого давления. Насыщение земной коры газовой компонентой изменяет ее вязкоупругие характеристики. Дисперсия плотности породы приводит к дисперсии скоростей распространения литосферных волн, т.е. изменению соотношения между скоростями распространения продольных V и поперечных V волн.
Скорость распространения продольных волн в среде
поперечных
p
s
p
s
где Е - модуль упругости Юнга, G - модуль сдвига, ρ - плотность среды [см. Курс физики Г.А.Зисман, О.М.Тодес, т.I, стр.301-302, издательство Наука. М., 1964 г.]. В качестве функции чувствительности изменения вязкоупругих характеристик среды выбрана функция V/V, изменяющаяся во времени и в пространстве от азимута трассы измерений. Реализации функции V/V до и после состоявшихся сейсмических ударов на измеряемых трассах иллюстрируются графиками фиг.2.
Известно [см, например, Р.Дуда, П.Харт. «Распознавание образов и анализ сцен», перев. с англ. М.: Мир, 1976, стр.319-331], что максимум информации о происходящих процессах содержится в образе объекта. Образом объекта земной коры является спектр сигнала измеряемой функции.
Для получения спектра-образа используют математические процедуры спектрального анализа регистрируемых функций путем их разложения в тригонометрический ряд, обеспечивающий наивысшую точность отображения, [см, например, Н.С.Пискунов. Дифференциальное и интегральное исчисления для ВТУ3-ов, учебник, т. II, 5-е издание. М.: Наука, 1964 г., Ряды Фурье, стр.180-182, 218-221].
По определению, спектр Фурье F(jω) функции  рассчитывается из соотношения:
2
1
t-t - интервал выборки изменений. В соответствии с теоремой Котельникова, для адекватного представления выборки измерений частотным спектром F(jω) необходимо, чтобы интервал Δt дискретизации
отсчетов временной функции удовлетворял условию
max
(F - максимальная частота спектра сейсмического сигнала), а объем выборки составлял не менее 600 отсчетов.
При выполнении данного условия малейшие изменения спектра выбранной функции будут обнаружены.
Расчет спектра осуществляют алгоритмами быстрого Фурье-преобразования (БПФ) по стандартным программам, входящим в комплект специализированного программного обеспечения ПЭВМ типа MATHCAD, ER MAPPER, (см., например, специализированное программное обеспечение MATHCAD 6.0 PLUS, издание 2-е стереотипное. М.: Информационно-издательский дом «Филинъ», 1997 г., стр.441).
Расчетные спектры Фурье выборок от функций  (фиг.2) до и после состоявшихся сейсмических ударов иллюстрируются графиками фиг.3.
Возникновение на измеряемой трассе газонасыщенного объема земной коры зоны подготавливаемого землетрясения изменяет параметры регистрируемого сигнала: ширину спектра, пропорции в соотношении амплитуд спектральных линий, средневзвешенную частоту спектра, величину затухания. С началом переходного к сейсмическому удару процесса скорость изменения функции  увеличивается, что отражается в спектре как увеличение его средневзвешенной частоты.
Момент изменения средневзвешенной частоты спектра отождествляют с возникновением на измеряемой трассе аномальной зоны подготавливаемого землетрясения и началом сейсмического процесса.
После обнаружения аномалии на трассе измерений осуществляют ее пеленгацию. Известно (см. например, Корн Г., Корн Т. «Справочник по математике для научных работников и инженеров» перевод с англ. М.: Наука, 1971 г., раздел «Аналитическая геометрия», стр.73-74), что положение радиус вектора в пространстве полностью определяется его косинус направляющими.
дл
Из аналитической геометрии известно, что косинус направляющая вектора равна отношению его проекции на данную ось к длине вектора. Длина вектора полной скорости (V) находится, как корень квадратный из суммы его проекций:
x
y
где V,V - скорости продольных волн для данной сейсмостанции по координатам x, y.
Для каждой сейсмостанции измеряется направление
Координаты аномалии находят как точку пересечения радиус векторов скоростей от нескольких сейсмостанции, как это иллюстрируется фиг.1.
После обнаружения аномалии и ее пеленгации на измеряемых трассах осуществляют прогноз параметров ожидаемого сейсмического удара. Для чего отслеживают изменение энергии сигнала.
В соответствии с законом Рэлея, энергию спектра сигнала рассчитывают как:
2
1
1
2
Изменение затухания сигнала на измеряемой трассе во времени определяют как разницу |ΔE|=E(t) - E(t) где t, t - время выборок измерений, для которых рассчитывался спектр Фурье.
Известна (см. «Краткосрочный прогноз катастрофических землетрясений с помощью радиофизических наземно-космических методов». Доклады конференции, ОИФЗ им. О.Ю.Шмидта. М.: РАН, 1998 г., стр.9) примерная зависимость размеров зоны подготавливаемого землетрясения R[км] от магнитуды (М) ожидаемого сейсмического удара.
R[км]≈ехр(М)
Известно также, что сама физическая величина и скорость ее изменения связаны дифференциальным уравнением первой степени, общим решением которого является экспоненциальная зависимость (см. Н.С.Пискунов. «Дифференциальное и интегральное исчисление для ВТУ3-ов». М.: Наука, 1976 г., стр.458). Экспоненциальная зависимость обладает тем свойством, что по ее дискретным значениям можно восстановить всю функцию.
Экспоненциальная зависимость затухания энергии сигнала на измеряемых трассах при увеличении размеров (R) аномальной зоны иллюстрируется графиком фиг.4
По величине затухания сигнала судят о размерах аномальной зоны. Расчет размеров аномальной зоны [R,км] для оценки ожидаемой магнитуды (М) сейсмического удара приведен ниже в примере реализации способа.
Пример реализации способа
Реализация способа осуществлялась на Камчатском Геофизическом полигоне РАН. Схема полигона иллюстрируется фиг.1. Полигон включает 6 сейсмических станций: Русская (RUS), Петропавловск (PET), Шипунский (SPN), Мыс Козлова (MKZ), Крутоберегово (КВТ), Авача (AVH).
Сейсмические станции мониторинга фиксируют приход продольной и поперечной волн от каждого землетрясения, произошедшего в фокальной зоне расположенной в 150 км от побережья полуострова Камчатки.
Фокальная зона Камчатки представляет сочленение континентальной и океанической плит, последняя из которых осуществляет движение под континентальную плиту. К северу, в зоне мониторинга станций КВТ, имеет место сочленение этих плит с более северной, т.е. имеет место тройное сочленение, сейсмическая активность которого очень высока. Фокальная зона землетрясений - от севера Камчатки на юг, вплоть до Индонезии. В мониторинге среды станциями RUS, PET, КВТ, MKZ, SPN, AVH, использованы данные, полученные только от местных землетрясений.
Функциональная схема системы измерений, реализующей способ, иллюстрируется фиг.5.
Система измерений (фиг.5) содержит размещенные в Камчатском регионе шесть сейсмостанций (1…6). Каждая из сейсмостанций, для азимутов трасс измерений в фокальной зоне, содержит группу из двух трехкомпонентных сейсмоприемников, (7, 8) оси чувствительности одноименных компонент которых взаимно параллельны. Сейсмоприемники размещены на измерительной базе (9), кратной размерам средней длине волны измеряемого спектра. Базы располагают перпендикулярно измеряемым трассам. Сейсмоприемники устанавливаются заглубленными, соосно друг другу.
px
л
p
Через встроенные усилители 10, 11, 12 сигналы с выхода сейсмоприемников подают на соответствующие входы канального коммутатора 13. Канальный коммутатор реализует функцию отношений сумм сигналов и сигнал пеленга аномалии (V/gV). Результирующие сигналы от канального коммутатора 13 оцифровывают в АЦП 14, 15 записывают в буферное ЗУ 16 для последующего ввода данных в компьютер 17, в стандартном наборе элементов: процессор 18, оперативное ЗУ 19, винчестер 20, дисплей 21, принтер 22, клавиатура 23. Компьютеры 17 объединены в локальную вычислительную сеть и обеспечена синхронизация их работы во времени.
Все элементы устройства представляют собой существующие технические разработки. В устройстве использованы элементы виброизмерительной аппаратуры фирмы Bruel & Kjairi ENDEVCO, Дания, следующих моделей: встроенные усилители 10, 11, 12 модели 2626, 2628; канальный коммутатор, АЦП, буферное ЗУ (элементы 13, 14, 15, 16) - многофункциональный блок 3560-L.
2
Сейсмическая волна характеризуется избыточным давлением Р[н/м]=ρc|ν|, где ρ - плотность среды, с - скорость акустических волн,  - скорость колебательного движения, u - смещение частиц от положения равновесия в процессе колебательного движения. Поскольку давление во фронте сейсмической волны пропорционально скорости, для ее измерения на сейсмостанциях используют датчики давления. Известен класс устройств мембранного типа, позволяющих измерять инфразвуковые колебания земной коры (см., например, «Кондуктометрический датчик колебаний». Патент RU №2055352, G.01.N, 27/02, 1996 г.). Датчик имеет высокую чувствительность и линейность выходной характеристики в области инфрачастот 0,1…100 Гц, быстродействие 100 ms. Промышленная разработка кондуктометрического датчика представляет собой трехкомпонентный сейсмоприемник (см., например, «Способ определения гипоцентра нефтегазового месторождения». Патент RU №2150719, G01V 1/100, 2000 г.). Он может использоваться в любой пространственной ориентации, т.е. осуществлять измерения во взаимно ортогональных плоскостях.
Анализ сейсмических данных проводился для периода с 1.01.2001 по 31.12.2004 гг. Интервалы времени, в которые рассматривались изменения состояний среды, составляли ~ 16, 8, 4, 2 суток до состоявшегося землетрясения и соответственно 2, 4, 8, 16 суток после. Обработка сигналов на каждой станции проводилась по зарегистрированным сейсмическим волнам с учетам направления их прихода, т.е. выбирались различные азимуты (азимут может иметь значение от 0 до 360) и решение зависело уже от того, откуда пришла волна землетрясений. Для станции КВТ было выбрано два промежутка азимутов: первый от 0 до 130, второй от 50 до 110; по каждому интервалу азимутов было проведено независимое исследование.
Для станций MKZ - от 70 до 130 и от 160 до 220. SPN - 70-130 и 110-190. PET-90-150, 100-130 и 180-230. RUS-60-110, 90-200 и 180-230.
p
s
cp1
ср2
cp3
1
2
3
p
s
p
s
х
Реализация функций сигналов (V/V) для 2 станций (KВТ, MKZ) иллюстрируются графиком фиг.2. Выборки дискретных отсчетов сигналов для станций проводились в интервалах времени 16, 8, 4 суток до и после состоявшихся землетрясений. Спектры Фурье этих выборок иллюстрируются фиг.3. При насыщении земной коры газовой компонентой среда становится диспергирующей. Возникновение подобной аномалии на измеряемых трассах приводит к изменению амплитудных соотношений между гармониками спектра, как это иллюстрируется фиг.3. Средневзвешенное значение частоты расчетных спектров соответственно составили F (16 суток) ≂ 5 Гц, F (8 суток)≂Гц, F (4 суток) ≂ 8 Гц. Энергия сигналов выборок измерений, вычисленная по соотношению Рэлея, соответственно равнялась: E(16 суток)≂0,306 мВт, E(8 суток) ≂ 0,0407 мВт, E(4 суток) ≂ 0,0203 мВт. Поскольку измеряемый сигнал представляет отношение двух функций (V/V), то измеряемое значение энергии сигнала (Е) не зависит от абсолютных значений величин V, V, т.е. является устойчивым параметром. Постоянная времени затухания сигнала на измеряемой трассе, для экспоненциального процесса:
Расчетная величина установившегося значения энергии сигнала накануне сейсмического удара
min
max
Признаком-предвестником землетрясения является увеличение размеров зоны подготавливаемого землетрясения (аномальной зоны). В первом приближении можно считать, что энергия сигнала на трассе распространения убывает обратно пропорционально E~1/R радиусу R аномальной зоны, как это иллюстрируется фиг.1. Откуда размеры зоны подготавливаемого землетрясения
Ожидаемая магнитуда сейсмического удара ln[R, км]=М, ln[147 км] ≈ 5,1 балла.
го
Сводный график расчетных значений землетрясений по всем шести станциям сети в интервале времен наблюдения 2001 г. - 2003 г. (в наблюдаемых азимутах станций) иллюстрирует фиг.6. Сходимость прогнозируемых результатов (с точность до 2 знака) соответствуют апостериорным данным состоявшихся землетрясений. Эффективность способа характеризуется высокой достоверностью, масштабностью территорий, охваченной мониторингом и большим интервалом времени упреждающего прогноза.