Алгоритм 3D кинематической томографии на отраженных волнах для построения начальной скоростной модели для миграционного скоростного анализа

Гадыльшин Кирилл Геннадьевич; Неклюдов Дмитрий Александрович; Протасов Максим Игоревич; Исаков Николай Григорьевич; Правдухин Андрей Петрович

УДК 528.44

DOI: 1G.183G3/2618-981X-2G18-3-265-273

АЛГОРИТМ 3D КИНЕМАТИЧЕСКОЙ ТОМОГРАФИИ НА ОТРАЖЕННЫХ ВОЛНАХ ДЛЯ ПОСТРОЕНИЯ НАЧАЛЬНОЙ СКОРОСТНОЙ МОДЕЛИ ДЛЯ МИГРАЦИОННОГО СКОРОСТНОГО АНАЛИЗА

Кирилл Геннадьевич Гадыльшин

Институт нефтегазовой геологии и геофизики им. А. А. Трофимука СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Коптюга, 3, кандидат физико-математических наук, старший научный сотрудник, тел. (383)33G-27-96, e-mail: gadylshin@gmail.com

Дмитрий Александрович Неклюдов

Институт нефтегазовой геологии и геофизики им. А. А. Трофимука СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Коптюга, 3, кандидат физико-математических наук, старший научный сотрудник, тел. (383)330-27-9б, e-mail: dmitn@mail.ru

Максим Игоревич Протасов

Институт нефтегазовой геологии и геофизики им. А. А. Трофимука СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Коптюга, 3, кандидат физико-математических наук, старший научный сотрудник, тел. (383)33G-27-96, e-mail: protasovmi@ngs.ru

Николай Григорьевич Исаков

ООО «Ингеосервис», 625019, Россия, г. Тюмень, ул. Республики, 211, заместитель генерального директора по обработке материалов сейсморазведки, тел. (345)221-52-94, e-mail: isakov@ingeos.info

Андрей Петрович Правдухин

ООО «Ингеосервис», 625019, Россия, г. Тюмень, ул. Республики, 211, начальник отдела обработки сейсмических данных, тел. (345)221-52-94, e-mail: andrey.pravduhin@ingeos.info

В настоящей работе мы предлагаем алгоритм 3D кинематической томографии для построения начальных глубинно-скоростных моделей для миграционного скоростного анализа. Он основывается на использовании результатов «временной» обработки, а именно: стандартного скоростного анализа и PSTM. На реальном примере показывается его эффективность.

Ключевые слова: глубинная миграция, глубинно-скоростная модель, томография, отраженные волны.

ALGORYTHM OF 3D REFLECTION TOMOGRAPHY FOR INITIAL VELOCITY MODEL BUILDING FOR MVA

Kirill G. Gadylshin

Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, 3, Prospect Аkademik

Koptyug St., Novosibirsk, 63GG9G, Russia, Ph. D., Senior Researcher, phone: (383)33G-27-96, e-mail: gadylshin@gmail.com

Dmitry A. Neklyudov

Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, 3, Prospect Аkademik Koptyug St., Novosibirsk, 63GG9G, Russia, Ph. D., Senior Researcher, phone: (383)33G-27-96, e-mail: dmitn@mail.ru

Maxim I. Protasov

Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, 3, Prospect Akademik Koptyug St., Novosibirsk, 630090, Russia, Ph. D., Senior Researcher, phone: (383)330-27-96, e-mail: protasovmi@ngs.ru

Nikolay G. Isakov

Ingeoservis LTD, 211, Respubliki St., Tyumen, 625019, Russia, Chief Specialist for Seismic Data Processing, phone: (345)221-52-94, e-mail: isakov@ingeos.info

Andrey P. Pravduhin

Ingeoservis LTD, 211, Respubliki St., Tyumen, 625019, Russia, Head of Seismic Processing Department, phone: (345)221-52-94, e-mail: andrey.pravduhin@ingeos.info

In this paper, we propose an algorithm for 3D kinematic tomography for constructing initial-depth-velocity models for migration velocity analysis. It is based on using the results of "time" processing, namely standard speed analysis and PSTM. Real example demonstrates effectiveness of the proposed approach.

Key words: prestack depth migration, velocity model building, migration velocity analysis, reflected waves.

Введение

В настоящее время глубинная миграция до суммирования (англ: Pre-Stack Depth Migration, PSDM) стала необходимым инструментом в практике обработки сейсмических данных. Результаты глубинной миграции, по устоявшейся традиции, называются глубинными сейсмическими изображениями изучаемой геологической среды. Их ценность особенно отчетливо проявляется в районах со сложным строением, при наличии сильной латеральной изменчивости свойств геологического разреза. Степень достоверности «глубинных» изображений полностью зависит от адекватности глубинно-скоростной модели среды, в которую выполняется PSDM. Для построения глубинно-скоростных моделей, как правило, используются алгоритмы, основанные на использовании кинематических характеристик зарегистрированных волновых полей. Детальный обзор существующих методов построения глубинно-скоростных моделей выходит за рамки настоящей статьи, отметим лишь основные особенности применения алгоритмов кинематической томографии для данных сейсморазведки по методу отраженных волн (МОВ). Основной особенностью является тот факт, что получение входных данных - времен пробега отраженных волн - весьма затруднено (по сравнению с кинематической томографией, основанной на использовании проходящих, головных или рефрагированных волн, наблюдающихся в первых вступлениях сейсмических записей). Это связано с двумя факторами: 1) объем сейсмических данных, регистрируемых в 3D МОВ, очень велик, в настоящее время речь идет о десятках (и сотнях) терабайт; 2) качество самих данных -сейсмограммы МОВ, как правило, представляют собой сложную картину интерферирующих вступлений волн различных типов.

Сама идентификация отраженных волн даже на данных, прошедших предварительную обработку, весьма непроста: протяженные вступления отраженных волн трудно различимы «на глаз», зачастую маскируются сильным фоном когерентных и случайных помех. Все это практически полностью исключает возможность «ручного» пикирования времен вступлений отраженных волн.

Стандартным практическим методом построения глубинно-скоростных моделей в настоящее время стал «миграционный скоростной анализ» (англ. Migration Velocity Analysis, [1], MVA). Основная идея, заложенная в MVA, в упрощенном виде может быть описана следующим образом. При наличии некоторой, весьма приближенной, начальной глубинно-скоростной модели среды PSDM может быть проведена в нее независимо для набора сейсмограмм общего удаления (ОУ). Набор полученных изображений, соответствующих каждому выносу (и набору азимутов), упорядочиваются в так называемые сейсмограммы общей точки изображения (англ. Common Image Gathers, CIG). Каждая трасса CIG соответствует глубинному изображению в фиксированной точке площади при определенном выносе (и азимуте), т. е. CIG состоит из всех мигрированных трасс, полученных с разных удалений и изображающих одну и ту же «область» среды. Если PSDM была выполнена с «правильной» скоростной моделью, то мигрированные отражения на CIG будут выпрямлены, так как изображения одной и той же отражающей площадки должны находиться на одинаковой глубине, в независимости от того, с какого выноса строилось это изображение. Если же скоростная модель была неадекватна, мигрированные отражения на CIG будут иметь кривизну (или же, говоря другими словами, будут иметь остаточный годограф). Анализ остаточных годографов, проведенный по всему объему сейсмического изображения, позволяет скорректировать текущую глубинно-скоростную модель: кривизна каждого когерентного события на CIG, которое отождествляется с отражателем, пере-считывается в невязку времен пробега соответствующей отраженной волны от этой отражающей площадки (под невязками времен здесь подразумевается разница времен пробега отраженной волны в текущей и «правильной» скоростной модели). Полученные невязки идут на вход процедуре кинематической томографии, в результате работы которой получается уточнение скоростной модели, минимизирующее невязки времен пробега и, соответственно, уменьшающее кривизну остаточных годографов на CIG.

Недостаток MVA заключается в необходимости делать PSDM на каждом шаге уточнения скоростной модели. При обработке данных сейсморазведки 3D это может занимать очень много времени. Цель настоящей работы - предложить простой и относительно дешевый алгоритм 3D-томографии для того, чтобы сократить количество «дорогих» итераций MVA на первых этапах построения скоростной модели. Для этого мы предлагаем использовать данные, полученные при «временной» обработке.

Описание метода

PSDM всегда выполняется на заключительном этапе обработки сейсмических данных. На предшествующих этапах, как правило, всегда строятся «временные изображения» изучаемого района с помощью временной миграции до суммирования (англ. Pre-stack Time Migration, PSTM). На кубах изображений PSTM производится предварительная интерпретация, пикирование основных отражающих горизонтов. На этапах, предваряющих глубинную миграцию (мы будем называть этот этап этапом «временной» обработки) всегда строятся скорости суммирования ОГТ, V0 гт, и скорости временной миграции, VM j G (которые, как предполагается, аппроксимируют среднеквадратические скорости, VrM s). Имея эти данные, можно построить эффективный алгоритм кинематической томографии на отраженных волнах, с помощью которого существенно сократить число очень дорогих итераций MVA. Ниже мы опишем детально предлагаемый алгоритм. Вначале сформулируем, что же является входными данными для предлагаемого алгоритма.

1. Начальная глубинно-скоростная модель. Как правило, она строится из имеющихся скоростей временной миграции путем пересчета их в интервальные скорости по формуле Дикса. В большинстве случаев даже в районах с умеренными латеральными вариациями скоростей подобная модель весьма далека от приемлемой глубинно-скоростной модели, которая бы обеспечила построение достоверного глубинного изображения.

2. Набор снятых на суммарном кубе PSTM основных отражающих горизонтов Т/М(Х1М, YIM), j = 1,..N, задаваемых в виде таблицы (XIM YIM Т1М), где Хш, YIM - координаты точки на поверхности наблюдений в изображении PSTM, TIM - соответствующее время.

3. Куб скоростей ОГТ, Vo гт ((X cmp >Ycmp>To'>clz[) (возможно, зависящих от азимута, az). Xcmp, Ycmp - координаты средних точек на поверхности наблюдений, Т0 - время (двойное) по нормальному лучу.

Далее опишем сам алгоритм построения глубинно-скоростной модели по этапам. На первых этапах производится «численная демиграция через глубину» снятых «временных» отражающих поверхностей. Это происходит с использованием хорошо известной в сейсмике концепции «image ray» [3], т. е. «луча, строящего изображение». Этот термин используется для обозначения фиктивного луча, подходящего к поверхности наблюдения по нормали и соединяющего точку на поверхности наблюдения с отражающим горизонтом в глубине (далее для краткости мы будем обозначать его как «image»-луч).

Описанные ниже процедуры 1-5 производятся для каждой «временной» поверхности независимо.

1. Для каждой точки на поверхности Т/М(Х im, YIM) в начальную модель трассируется «image»-луч. Он выпускается с поверхности наблюдений из текущей точки ( , ) по нормали и продолжается в начальной глубинно-скоростной модели, пока не истечет время TIM/2. В итоге определяется положе-

ние соответствующей точки на отражающем горизонте в глубине (X, У, X). Продолжая таким образом «image»-лучи для всех точек текущей «временной» поверхности, строится ее отображение в глубине в текущей скоростной модели. Таким образом, происходит «глубинная» миграция «временной» поверхности с помощью «image»-лучей [3]. Отметим, что, трассируя «image»-лучи, мы не получаем никакой информации о наклонах отражающих площадок. Все, что можно сказать, - «image»-луч касается отражающей площадки по истечении времени Т1М /2.

2. Определение нормалей к мигрированной поверхности в глубине. На предыдущем этапе было получено положение всех точек отражателя в глубине, что можно рассматривать как задание набора значений некоторой функции 1] (X, У). Далее необходимо определить численно нормали к поверхности 1] (X, У) каждой точке. Для каждой точки в глубине в некоторой наперед заданной апертуре набираются ее соседи и строится локальная плоскость, наиболее близкая ко всем выбранным точкам в смысле наименьших квадратов. Таким образом, решается классическая задача трехмерной линейной регрессии. В результате по нерегулярному «облаку» мигрированных точек отражателя (X, У, X) строится набор локальных отражающих площадок, каждая из которых характеризуется своей нормалью. Отметим, что при миграции в неадекватную глубинно-скоростную модель изначально гладкие «временные» горизонты могут потерять гладкость в глубине. Мы предполагаем, что сохраняется локальная гладкость, что достигается отбрасыванием максимальных/минимальных значений

(X, У) при решении задачи локальной линейной регрессии и локальным сглаживанием полученных нормалей в некоторой апертуре.

3. Из каждой отражающей площадки к поверхности наблюдений трассируется нормальный луч. Таким образом численно строится соответствие мигри-рованной поверхности в глубине с кубом скоростей ОГТ, Уоп(ХСМР>УСМР,Т0). Этапы 1-3 фактически представляют собой «численную демиграцию через глубину» снятых «временных» отражающих поверхностей на кубе PSTM.

4. От каждой отражающей площадки с известной нормалью строится конус отраженных лучей. Для этого из глубины с некоторым шагом по углу раствора (угол отражения) и азимуту выпускается набор лучей. Каждые два сегмента, которые удовлетворяют в точке выхода (т. е. в точке отражения) закону Снеллиуса, образуют луч, соответствующий волне, отраженной от заданной площадки при заданном угле отражения и при заданном азимуте. Предположим, каждый из двух сегментов луча отраженной волны выходит на поверхность. Тогда отраженный луч на поверхности можно однозначно охарактеризовать следующими координатами: ХсМР,УсМР - координатами средней точки двух сегментов, абсолютным выносом И, т. е. расстоянием на поверхности между точками выхода и азимутом на поверхности а. При трассировке сегментов также вычисляется время пробега отраженной волны в текущей модели, которое можно записать как: Тс^с = ТсАЬс(ХСМР,¥СМР'к,а).

5. Далее попробуем определить «наблюденные» времена в точке выхода отраженных лучей. Для этого мы предлагаем использовать скорости суммирования ОГТ и в качестве «реальных» наблюденных времен взять их гиперболическую аппроксимацию:

Для каждого отраженного луча, вышедшего на поверхность, не составляет труда по координатам ХСМр ,Усмр найти соответствующее этой точке Т0 (для каждой мигрированной «временной» поверхности на шаге 3 мы уже определили «карту» Т0(Х, У); здесь же еще раз отметим, что описываемая процедура производится с каждой поверхностью независимо) и, следовательно, найти нужное значение У0 гт. И рассчитать время по формуле (1) в текущей точке на поверхности. Отметим, что подобная аппроксимация реальных времен отраженных волн нередко встречается на практике [2]. Это связано с тем, что скорости суммирования ОГТ представляют собой интегральный параметр, который определяется по значительной выборке трасс, т. е. он очень устойчив. Скорости ОГТ очень тщательно контролируются обработчиком. В принципе, это сугубо кинематический атрибут, который в общем случае может не иметь ничего общего с реальным распределением интервальных скоростей в среде, но всегда характеризует поведение годографа отраженных волн в смысле регрессии, т. е. характеризует гиперболу, наилучшим образом «ложащуюся» на реальный годограф отраженной волны. Таким образом, мы получаем входные данные для стандартной кинематической томографии на отраженных волнах, где минимизируются невязки времен пробега отраженных волн йТ = Тов$ - ТСА1^с, здесь времена вычисляются в точках выхода отраженных лучей с помощью лучевой трассировки для Тсдьс и по формуле (1) для Тов8.

6. Решение обратной кинематической задачи. Согласно стандартной схеме, принятой в кинематической томографии, составляется система линейных уравнений, связывающая наблюденные невязки отраженных волн йТ с искомой поправкой параметров модели Ла:

Здесь М обозначает томографическую матрицу, элементы которой - производные от времен пробега по параметрам модели. В нашем случае модель параметризуется значениями скоростей на трехмерной сетке, в каждой ячейке скорость предполагается постоянной. Матрица М заполняется на шаге 5, каждая поверхность соответствует блоку из строк в единой томографической матрице. Ее размерность в изотропном случае определяется как ,

где - число успешно протрассированных отраженных лучей, -

ТОВ8 (ХСМР, УСМР >к>а)- 4 Т0 (хсмр , Усмр ) +

2

(1)

МАа = dT.

(2)

количество ячеек сетки в области реконструкции. Отметим, что в нашем случае число строк на несколько порядков превосходит число столбцов.

Далее мы следуем схеме, предложенной в работе [1]. Вместо исходной системы (2) решается предобусловленная система

LMSRAa' = LdT, (3)

где L - диагональная матрица взвешивания строк, R - диагональная матрица взвешивания столбцов, S - сглаживатель, представляющий собой 3D-свертку по пространству с гладким ядром. Система (3) решается методом IRLS (англ. Iterative Reweighted Least-Squares) [1] в норме L1.5, что подразумевает использование итерационного метода решения СЛАУ, такого как LSQR. Решение томографической системы с помощью IRLS позволяет получать решения, гораздо более устойчивые к резким «скачкам» невязок (outliers), что является существенной проблемой в подобного рода задачах при решении их стандартным методом наименьших квадратов.

После решения линеаризованной обратной кинематической задачи (3) и получения поправки к начальной скоростной модели шаги 1-6 повторяются в скорректированную модель. Таким образом происходит достаточно быстрое итерационное уточнение 3D глубинно-скоростной модели без привлечения дорогой процедуры PSDM. После нескольких глобальных итераций уточнение скоростной модели необходимо прекратить, так как гиперболичная аппроксимация самих времен пробега отраженных волн может оказаться недостаточной для достоверного восстановления скоростной модели, особенно в глубокой ее части. Полученный результат можно использовать как начальное приближение для стандартного MVA, который основан на гораздо более слабом предположении о форме остаточного годографа на CIG. Отметим, что «томографическое ядро», т. е. часть программы, отвечающая за решение обратной кинематической задачи, может быть одинакова что для описанного выше подхода, что для стандартной CIG-томографии.

Рис. 1. Начальная глубинно-скоростная модель (скорости изменяются в диапазоне 1 500-4 400 м/с)

Численный эксперимент

Описанный выше подход использовался для построения глубинно-скоростной модели для участка одной из площадей Карского моря. Для эксперимента был выбран фрагмент с размерами ~19 х 5 км, общая площадь примерно равна 100 км . Глубина построения скоростной модели равна 6 км. Начальная скоростная модель, полученная конвертацией скоростей миграции в интервальные скорости по формуле Дикса, представлена на рис. 1. На кубе PSTM было снято 8 поверхностей. Результат PSDM в начальной скоростной модели представлен на рис. 2, A. Как можно видеть по сейсмограммам общей точки изображения, эта модель весьма далека от удовлетворительной. В нижней части изображений когерентные события не видны вовсе, что весьма затруднило бы применение стандартной CIG-томографии. PSDM, выполненная в уточненную модель, полученную предложенным выше способом, дает весьма удовлетворительный результат даже без применения CIG-томографии, что показывает достаточную эффективность алгоритма. В некоторых случаях он мог бы заменить стандартный MVA. Хотя, напомним, его основное предназначение - использоваться до MVA чтобы, обеспечить удовлетворительную начальную модель, которая могла бы быть быстро уточнена с минимальным количеством итераций MVA.

Рис. 2. СЮ в начальной модели (А). СЮ в модели, построенной с помощью

описанного подхода (Б)

Работа выполнена при частичной поддержке грантов РФФИ № 17-0500001 (М.П.) и 18-35-00253 (К.Г.)

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. A decade of tomography / V. Woodward, D. Nichols, O. Zdraveva et al. // Geophysics. -2008. - Vol. 5. - P. VE5-VE11.

2. Fagin S. Model based depth imaging: monograph. - Talsa, USA : SEG, 1999. - 120 s.

3. Hubral P., Krey T. Interval velocities from seismic reflection time measurements: monograph. - Talsa, USA : SEG, 1980. - 203 s.