first break том 25, март 2007 особая тема Визуализация и интерпретация Автоматическое выделение элементов геологического строения из 3D сейсмических данных с помощью спектрального разложения многоканальных объемных данных. Automated delineation of geological elements from 3D seismic data through analysis of multichannel, volumetric spectral decomposition data Джонатан Гендерсон (Jonathan Henderson),1 Стивен Ж. Пурвс (Stephen J Purves)1 и Крис Лепард (Chris Leppard2) описывают некоторые новые методы объемного спектрального разложения, созданные для оптимизации обработки геологической информации, полученной из сейсмических данных. Спектральный состав регистрируемых сейсмических данных зависит от акустических свойств пород вдоль пути распространения луча. Вариации в спектре возникают вследствие вариаций в интерференции, природы границ между различными слоями (грубозернистые/тонкозернистые пачки слоев) и изменениями в заполнителе порового пространства, и особенно присутствие газа. Хотя в большинстве случаев при недостатке информации все же удается определить причину изменений в спектре сигнала, значительно больше информации можно получить из анализа пространственных изменений в спектре на нескольких частотах. Это обстоятельство привело к существенному интересу к методам спектрального разложения как к инструменту интерпретации сейсмических данных в задачах стратиграфического анализа (Partyka et al., 1999, Marfurt and Kirlin 2001, Johann et al., 2003) и при поиске углеводородов (Castagna et al., 2003). В целом, приложение спектрального анализа для интерпретации 3D сейсмических данных включает выделение данных, ассоциируемых с 2D плоскостями или горизонтами или построение 3D поверхностей для одной частоты. Для реализации полного потенциала методов спектрального разложения нужно использовать метод 3D пространственного разложения на нескольких частотах, результаты которого применяются для выделения трехмерных элементов в геологическом строении (GeoBodies). В этой статье описан набор новых методов 3D визуализации и выделения геологических тел на основе спектрального разложения в процессе интерпретации, разработанный в результате совместной исследовательской работы ffA и Norsk Hydro. Среди этих методов 2D и 3D спектральное разложение, генерация объемных спектральных атрибутов, 2D и 3D RGB визуализация в нескольких частотных диапазонах, 3D визуализация с помощью построения оптических RGB плотностей и выявление геологических тел (GeoBodies). Множество современных мощных инструментов для интерпретации в геофизике считаются слишком сложными для использования и поэтому неоправданно редко используются. Одна из целей этой работы, таким образом, состоит во внедрении инструментов объемного спектрального анализа в простой и интуитивный процесс обработки, чтобы они могли использоваться в составе стандартного набора методов интерпретации, применяемых в компании Hydro (Норвегия). Применение этого комплекса методов для интерпретации данных для нескольких районов, отчетливо демонстрируют богатство геологической информации, получаемой с помощью методов объемного спектрального разложения. 1 ffA, Newcastle upon Tyne, UK. 2 Hydro, Bergen, Norway. Спектральное разложение Спектральным разложением называется процесс анализа того, каким образом данный сигнал может быть представлен в виде суммы простых, хорошо известных базисных функций. Наиболее известным методом спектрального разложения является Фурье-анализ, в котором сигнал представляется в виде взвешенной суммы синусоидальных и косинусоидальных колебаний на различных частотах. Применяя преобразование Фурье, можно определить амплитуду и фазу косинусоидальных колебаний на всех частотах. Для преобразования Фурье предполагается, что спектральный состав сигнала постоянен 30 Hz 50 Hz во времени. Эти ограничения привели к разработке множества методов разложения нестационарных сигналов (таких как сейсмические данные), в числе которых: быстрое преобразование Фурье (Nawab and Quatieri, 1988); Wigner-Ville разложение (WVD) (Claasen and Mecklenbrauker, 1980); вейвлет-преобразование (Daubechies 1988; Mallat, 1989); пакетное вейвлет-преобразование (Coifman and Wicker-hauser, 1992) и разложение типа Matching Pursuit (MPD) (Mallat and Zhang, 1993, Chakraborty and Okaya, 1995). Временные/частотные представления, получаемые при применении этих методов неодинаковы и их использование зависит от решаемой задачи. Рис. 1 Вертикальный разрез амплитуд на 30 Гц и 50 Гц, рассчитанный с использованием постоянных полосовых фильтров (a) и постоянных Q-фильтров(b). Отметьте четкое выделение отклика на частоте 50 Гц с помощью постоянных Q-фильтров. © 2007 EAGE 87 особая тема Визуализация и интерпретация first break том 25, март 2007 При создании графа обработки методом спектрального разложения, который описывается ниже, основной целью было выбрать и сравнить отклики от сигнала возбуждения для ряда различных частот, а не построить сплошной спектр для входного 3D массива данных. Эффективный вычислительный метод для достижения этой цели состоит в оценке амплитудной и частотной компоненты отклика на данной частоте при помощи свертки входных данных с вейвлетом с известными характеристиками. Из множества известных вейвлетов был выбран вейвлет Габора (Gabor wavelet). Удобное свойство вейвлета Габора – это локальность его временных/частотных характеристик. Принцип неопределенности (Williams et al, 1991) выглядит следующим образом: (1) то есть, чем меньше ширина окна (∆t), тем хуже разрешающая способность по частоте и наоборот. Путем подбора соответствующих параметров вейвлета, неравенство в выражении 1 может быть минимизировано для всех частот с использованием вейвлета Габора. Помимо этого, вейвлета Габора является близкой аппроксимацией вейвлета Рикера (Ricker wavelet), часто используемого при обработке сейсмических сигналов. Выражение для вейвлета Габора: (2) Где: a - масштаб, T - временной сдвиг, f – частота модуляции и qp фаза атома (базисного вейвлета). Масштаб и частота модуляции определяют центральную частоту, разрешение по частоте и пространственную локализацию атома. Оптимальное соотношение между пространственной локализацией и разрешением по частоте достигается с помощью сохранения постоянной частотной модуляции (f) и изменением масштаба (a). Мы назвали это постоянным Q-методом, поскольку отношение между частотной полосой и центральной частотой остается постоянным. Постоянный Q-метод дает результат эквивалентный результату вейвлет-преобразования с вейвлета Габора. Если масштаб (a) остается постоянным, а частота модуляции (f) меняется, тогда полоса фильтра постоянна и не зависит от центральной частоты. Этот метод с постоянной полосой фильтра эквивалентен быстрому преобразованию Фурье с Гауссовым окном. Различия в пространственной локализации для постоянного Qметода и метода с постоянной полосой представлены на рис. 1. Вследствие преимуществ в локализованности частотно-временных параметров, постоянный Q-метод используется с целью анализа спектральных вариаций, связанных с отдельными сейсмическими волнами. Неудобство постоянного Q-метода заключается в том, что пространственная локализация зависит от частоты. Метод с постоянной полосой применяется, когда целью спектрального разложения является различение геологических элементов на основании их объемных свойств, например выделение больших соляных тел или газовых резервуаров. В подобных случаях это могут быть свойства всего частотного спектра, наиболее удобного для анализа, при этом могут быть полезны такие свойства метода с постоянной полосой, как дискретизация в частотной области и постоянство времени. Общность между кратковременным преобразованием Фурье и вейвлет преобразованием с вейвлетом Габора-Морле (Gabor-Morlet) предоставляют большую гибкость и контролируемость спектрального разложения, что позволяет создать простой механизм работы с 3D частотными вариациями с использованием многоканальной частотной информации для выделения геологических тел. Спектральные атрибуты Ценность методов спектрального разложения в интерпретации данных сейсморазведки дает возможность сравнивать 3D вариации в отклике на сейсмический сигнал на разных частотах. Хотя это 4D задача, она может быть преобразована в 3D задачу, путем создания множества спектральных атрибутов, дающих информацию о форме и спектре колебания в каждой точке сейсмического разреза (рис. 2). Метод спектрального разложения с постоянной полосой создает базу для вычисления спектральных атрибутов, позволяя создавать широкий диапазон описаний частотных спектров для всего объема, в том числе: ■ Основные параметры мощности Полная спектральная мощность Средняя мощность в полосе ■ Квартильная мощность Это измерение мощности для каждого квартиля в полосе частот. Каждый квартиль имя одинаковую длительность по частоте, как показано для положительной части спектра. 88 © 2007 EAGE first break том 25, март 2007 особая тема Визуализация и интерпретация ■ Отношение квартилей Отношение мощности первого и четвертого квартиля является индикатором симметрии распределения. ■ Частота для максимальной мощности ■ Статистические оценки Средняя частота Дисперсия (Полоса) Асимметрия ■ Q-фактор 'Q' или фактор качества – это мера достоверности (сфокусированности) сигнала. Она основа на соотношении между полосой и пиковой мощностью распределения. Частотные атрибуты применяются уже некоторое время в интерпретации данных сейсморазведки и могут быть использованы для подчеркивания особенностей в сейсмических данных, трудно выделяемых другими способами. Поскольку они вычисляются для полного частотного спектра, они менее подвержены влиянию шума, по сравнению с атрибутами, вычисленными для узкой полосы частот. Спектральные атрибуты дополняют методы многоканальной 3D визуализации и сегментации, которые являются основными целями данной работы. RGB визуализация Наиболее важное свойство методов спектрального разложения заключается в возможности получать и сравнивать отклик для разных полос частот. Для этого в разных областях, а также в спектральном разложении сейсмических данных (Hall and Trouillot, 2004), широко применяется метод RGB (Red-Green-Blue, красный-зеленый-синий) смешивания 2D изображений. Методы объемной RGB визуализации и построения оптических RGB плотностей для 3D сейсмических (c) Рис. 2 Временной срез (a) по кубу мгновенной амплитуды в сравнении со спектральными атрибутами, (b) 4th квартиль мощности, (c) спектральная асимметрия, (d) средняя частота. данных до настоящего времени не были опубликованы. Рис. 3 – цветовое пространство RGB. Метод RGB смешивания объединяет информацию из трех входных изображений для образования “полноцветного” выходного изображения. Для RGB изображения используется 3D цветовое пространство, при этом по осям отложены относительные цветов: красного, зеленого и синего (рис. 3). RGB изображение образуется в соответствии с позицией в 3D цветовом пространстве, при этом координаты по осям определяются значениями для трех входных изображений. (3) Здесь Cout(x,y,z) цвет для точки (x,y,z) в выходном изображении и IR(x,y,z), IG(x,y,z) и IB(x,y,z) значения для элемента объема в точке (x,y,z) во входных изображениях, которые присваиваются распределениям Red, Green и Blue соответственно. Качество смешанного RGB изображения зависит от доступного цветового разрешения. Многие системы 3D визуализации используют 8 битные (256 уровней) цветовые таблицы, такие дисплеи не могут адекватно представлять изображения в методе смешивания RGB. В большинстве случаев смешанные RGB изображения хорошего качества могут быть получены с 18 битным цветовым разрешением (262144 цветов), увеличение цветового разрешения до 21 бита (2097152 цветов) в некоторых ситуациях может быть полезным, но крайне небольшие улучшения заметны при дальнейшем увеличении цветового разрешения. В том случае, когда скорость построения изображений существенна, может быть полезно понизить цветовое разрешение до 15 бит (32768 цветов). Для создания 18 битного RGB изображения 6 бит (64 уровня) используются для каждого из каналов. Для оптимизации использования динамического диапазона были разработаны методы автоматического отображения входных данных в 3D цветовое пространство. Они основаны как на использовании анализа данных во всем 3D кубе так и для данных, соответствующих заданному горизонту; кроме этого доступны пользовательские функции трансформации. © 2007 EAGE 89 особая тема Визуализация и интерпретация first break том 25, март 2007 Построение оптических RGB плотностей и выделение геологических тел Объемные изображения, построенные методом смешения RGB, являются инструментом для простого и быстрого доступа к результатам спектрального разложения (рис. 4). Однако, как только интересующий объект идентифицирован при помощи информации от спектрального разложения, целесообразно использовать эту информацию для выделения и оконтуривания соответствующего 3D геологического тела (GeoBody). Общий подход для определения геологического тела – использовать пороговые соотношения и анализ взаимосвязи, при котором все элементы объема в пределах заданного диапазона соединяются между собой с присвоением общей метки. Преимущество этого подхода в том, что прямой визуальный эффект может быть получен при помощи соединением кривой оптической плотности с заданными параметрами элементов объема интересующего объекта. Это процесс гораздо сложнее применять к RGB изображениям, поскольку кривая оптической плотности должна быть построена в 3D цветовом пространстве, а не в 1D. Разработаны два способа построения кривой оптической плотности в 3D пространстве. Первый использует определение интересующей 3D области в RGB цветовом пространстве путем отображения 1D кривых плотностей индивидуально для каждого из каналов Red, Green и Blue (рис. 5). Преимущество данного подхода в том, что пользователь знает, какие значения он выбирает для каждого из трех входных объемов и он хорошо работает, когда все интересующие элементы объема непрерывны в 3D RGB пространстве. Этот метод также может использоваться для ступенчатых RGB кривых оптической плотности. Недостаток этого метода в том, что трудно построить кривую оптической плотности для выделения интересующих объектов в случае, когда элементы объема объекта исследования широко распределены или образуют множество кластеров в 3D цветовом пространстве. Вторая методика была разработана для устранения этого недостатка. Прямое наблюдение RGB изображения является наиболее простым методом для определения цветов, описывающих объект исследования, таким образом, представляется логичным простой выбор цветов в интерактивном режиме из самого RGB изображения. Во втором методе, разработанном для построения кривой оптической плотности, пользователь должен выделить интересующую область (Region of Interest, ROI), выделяя на 2D плоскости полигон, содержащий цвета, для которых установлен атрибут непрозрачности (рис. 6a). Затем, одно значения оптической плотности устанавливается для всех цветов в этой области и дополнительное значение устанавливается для остальных цветов. Далее, кривая оптической плотности используется для отображения всего RGB объема данных. Эта технология позволяет пользователю видеть только выбранные им цвета, при этом не нужно точно выделять границы объекта исследования. Она также позволяет подобрать Рис. 5 Выбор цветов с использованием индивидуальных кривых оптической плотности для каналов red, green и blue. Рис. 4 Объемное RGB изображение для различных сейсмических отражений (разный цвета) в разных сдвиговых блоках. (d) Рис. 6 Настройки цветовой оптической плотности для RGB изображения объекта исследования. (a) объект исследования в 2D сечении, (b) цвета, выбранные с настройками непрозрачности, (c) цвета, выбранные с настройками диффузной прозрачности, (d) Объем данных, построенный с выбором диффузно-прозрачных цветов. (RGB настройки объекта исследования могут быть выбраны для любого изображения, полученного методом смешения RGB.) 90 © 2007 EAGE first break том 25, март 2007 особая тема Визуализация и интерпретация наиболее подходящую кривую оптической плотности для объектов в соответствии с их цветами, распределенными в 3D цветовом пространстве. Если RGB изображение содержит большое количество цветов, может быть сложно гарантировать, что выделенный объект исследования содержит все мелкие особенности, выражаемые в переходах цвета. Вариации цвета в пределах изучаемого объекта могут быть учтены с помощью расширения границы объекта таким образом, чтобы они полностью окружали объект исследования. Дополнительно можно установить параметр диффузности цвета. Параметр диффузности определяет расстояние в 3D цветовом, при этом вес цвета, находящиеся в пределах этого расстояния от границы объекта, выбранной пользователем, попадают в область выделения (рис. 6b, 6c). Это дает большую степень контроля при выделении объектов. Эффект при изменении вертикальной протяженности области исследования или параметра диффузности проявляется немедленно при построении 3D плотностного изображения, что дает пользователю визуальную обратную связь (рис. 6d). Как только выбраны цветовые плотностные параметры, можно применить процесс сегментации для выявления элементов объема, образующих интересующие объекты. После этого к результату процесса сегментации изображения применяется анализ взаимосвязи элементов объема и алгоритм маркировки для выделения геологических тел (рис. 7). Рис. 7 Сравнение традиционных методов и метода RGB для выделения тел (a) применение метода мгновенных амплитуд, (b) тела выделены на основании дифференциации в пределах диапазона амплитуд в западной части изучаемой области, (c) метод RGB смешения для трех разных частот для одного массива данных, (d) тела выделены с помощью оптической плотностной RGB технологии. Процесс RGB спектрального разложения и его результаты На основании этапов обработки и визуализации, описанных выше, был построен простой в использовании процесс обработки (рис. 8), позволяющий эффективно применять спектральное разложение, метод 3D RGB смешивания и сегментацию для больших объемов данных. Этот процесс охватывает предварительную оценку 2D параметров, 3D обработку и анализ и выделение геологических тел. Фаза предварительной оценки 2D параметров (фаза I) процесса спектрального разложения (рис. 9) проводится с использованием программы-“помощника” спектрального разложения, в которой пользователь проводит пошаговую оптимизацию информации, полученной в результате процесса спектрального разложения и RGB смешивания. Фаза I процесса спектрального разложения включает следующие шаги: ■ Выбор данных для теста. Результаты теста получаются для временного или глубинного сечения объема данных. ■ Выбор параметров спектрального разложения (минимальная и максимальная частоты, количество частотных интервалов, используемый фильтр: постоянный полосовой или постоянный Q-фильтр). ■ Два режима предварительного просмотра результатов. Режим для одной частоты/атрибута и RGB режим, при котором отображаются отображения для трех часто или атрибутов одновременно. ■ Выбор частотных компонентов или атрибутов для расчета по всему объему данных. Рис .8 Схема процесса спектрального разложения. © 2007 EAGE 91 особая тема Визуализация и интерпретация first break том 25, март 2007 Цель фазы I состоит в изучении частотных интервалов и атрибутов и выборе тех из них, что наиболее ясно предоставляют информацию в результате спектрального разложения. Это позволяет произвести эффективные вычисления для всего объема данных без избыточных расчетов. Фаза II процесса спектральной декомпозиции служит для поддержки объемного анализа RGB смешанных данных (обратите внимание, что объемный RGB массив может быть получен из любых трех массивов сейсмических данных, не обязательно полученных в фазе I процесса спектрального разложения). 3D RGB изображение в главном окне программы сопровождается дополнительным окном 2D просмотра для контроля параметров цветовой плотности изучаемых объектов и процесса сегментации. Фаза II процесса спектрального разложения включает следующие шаги: ■ 3D объемное или послойное изображение RGB смешанных данных ■ Регулировка масштабных коэффициентов метода RGB смешивания. ■ Построения RGB с учетов выбранной модели цветовой плотности ■ Сегментация RGB изображений в соответствии с кривой цветовой плотности, определенной пользователем ■ Выделение геологических тел при помощи анализа взаимосвязи элементов объема сегментированного RGB изображения Сравнение результатов для осадочных систем Послойные сейсмические атрибуты, такие как амплитуда огибающей (коэффициент отражения), часто используются для определения литологии, а шаблоны атрибутов, пространственно распределенные в сейсмическом разрезе, часто интерпретируются как признаки осадочных систем. Сравнения карт, построенных с помощью этих традиционных методов и карт, построенных с помощью метода RGB спектрального разложения, показывают важность метода RGB спектрального разложения в получении геологической информации прямо из сейсмических данных. Послойные карты, показанные на рис. 10 демонстрируют сравнение сейсмических данных на шельфе (для полного набора удалений, с миграцией по времени). Входные сейсмические данные были очищены от шума по методике ffA для устранения влияния шума при построении последующих изображений. Карта горизонта для амплитуды огибающей (рис. 10a) четко выделяет меандрирующую канальную систему, протягивающуюся из верхнего правого угла в нижний левый угол карты. Однако трудно выделить детали канальной системы, такие как геометрию отдельных каналов и их внутреннюю архитектуру и морфологию, в силу отсутствия контраста амплитуд между заполнителем канала, кромкой канала и засклоновыми осадками. По сравнению с этой картой, карта горизонта, полученная по данным RGB спектрального разложения (рис. 10b) дает гораздо большую Рис . 9 Фаза I процесса частотного разложения. 92 © 2007 EAGE first break том 25, март 2007 особая тема Визуализация и интерпретация инструменты RGB визуализации эффективно решают проблемы получения информации для 3D структурного анализа. Инструменты метода цветовой RGB плотности и сегментации предоставляют механизм построения 3D объектов, что позволяет использовать информацию, полученную при объемном спектральном разложении на последующих стадиях в процессе моделирования и интерпретации. Рис. 10 Сравнение послойных карт построенных по (a) коэффициентам отражения и (b) RGB спектральному разложению. четкость в изображении геометрии, морфологии и внутренней архитектуры канальной системы. Применяя возможности RGB сегментации, могут быть выделены интересующие элементы осадочной системы и сегментированы в виде 3D геологических тел для дальнейшего использования для вычислений для всего объема в качестве входной геологической модели. Выводы Метод объемного спектрального разложения может использоваться для получения информации из 3D сейсмических данных, недоступной для получения другими методами. Данная работа имеет практические результаты применения метода объемного спектрального разложения в процессе 3D интерпретации. Связь между 2D и 3D методами спектрального разложения и новые 3D Литература Castagna, J.P., Sun, S.J., and Siegfried, R.W. [2003] Instantaneous spectral analysis: Detection of low-frequency shadows associated with hydrocarbons. The Leading Edge, 22, 120-127. Chakraborty, A. and Okaya, D. [1995] Frequency-time decomposition of seismic data using wavelet-based methods. Geophysics, 60, 1906-1916. Classen T. and Macklenbrauker W. [1980] The Wigner distribution- A tool for time-frequency signal analysis (3 parts). Philips J. Res, 35, 217-250, 276-300, 372-389. Coifman, R. R. and Wickerhauser, M. V. [1992] Entropy-based algorithms for best basis selection. IEEE Trans. Inf. Th., 38,713-719. Daubechies, I., [1988] Orthonormal bases of compactly supported wavelets. Comm. in Pure Applied Math, 41, 909-996. Hall. M. and Trouillot, E. [2004] Predicting stratigraphy with spectral decomposition. Canadian Society of Exploration Geophysics, Expanded Abstracts. Johann, P., Ragagnin, G., and Spinola, M. [2003] Spectral Decomposition Reveals Hidden Features in Amplitude Maps from a Deep Water Reservoir in the Campos Basin. 73rd Meeting, Society of Exploration Geophysics, Expanded Abstracts. Mallat, S. and Zhang, Z. [1993] Matching Pursuit with Time-Frequency Dictionaries. IEEE Trans Sig. Proc., 41, 3397-3415. Mallat, S. [1989] A theory of multi resolution signal decomposition, the wavelet representation: IEEE Trans Pat. Anal, and Mach. Intell., 11, 674693. Marfurt, K.J. and Kirlin, R.L. [2001] Narrow-band spectral analysis and thin-bed tuning. Geophysics, 66, 1274-1283. Nawab, S. and Quatieri, T. [1988] Short-time Fourier transform. In Lim, J., and Oppenheim, A., (Eds) Advanced topics in signal processing. Prentice Hall Signal Processing Series, 289-337. Partyka, G.A., Gridlay, J.M., and Lopez, J [1999] Interpretational Applications of Spectral Decomposition in Reservoir Characterization. The Leading Edge, 18, 3, 353-360. Williams W. J., Brown M. L., and Hero A. O. [1991] Uncertainty, information, and time-frequency distributions: Advanced signal processing algorithms, architectures and implementations II, SPIE 1566, 144-156. © 2007 EAGE 93