first break том 24, декабрь 2006
техническая статья
Resolution and uncertainty in spectral decomposition
Разрешение и неопределенность спектрального разложения
Matt Hall*
За последние пять лет спектральное разложение стало одним из основных способов интерпретации сейсмических данных. Этот метод, обычно используемый для качественной оценки сейсмического геоморфологического анализа (напр., Marfurt and Kirlin, 2001), все чаще используется количественно при расчете толщины слоев (описано в Partyka et al., 1999). Он также применяется для прямого обнаружения углеводородов (Castagna et al., 2003). В широком смысле мы используем принципы спектрального разложения всякий раз, когда мы исследуем спектр или выделяем отдельный импульс.
На данный момент для разложения обычно используется «кратковременное» преобразование Фурье, но в последнее время к списку возможностей присоединились и другие методы: S-преобразование (Stockwell et al., 1996), преобразование элементарных волн, и matching pursuit (напр., Castagna and Sun, 2006). В литературе часто упоминается «разрешение» этих различных способов. Неожиданностью может явиться то, что принцип неопределенности Гейзенберга иногда используется для подтверждения превосходства одного метода перед другим. В этой статье мы пытаемся ответить на вопрос, что конкретно имеется в виду под «разрешением» в этом контексте, каковы его ограничения, и какое отношение к этому имеет квантовая теория.
Естественное свойство
Принцип неопределенности Гейзенберга следует из классического неравенства Коши-Шварца и является одним из основных принципов квантовой теории. Сам Вернер объяснил его так:
не ограничивает предел наших знаний, он описывает реальное положение вещей: электрон не может одновременно обладать точным местоположением и моментом. Это тревожное наблюдение является основополагающим в так называемой копенгагенской интерпретации квантовой теории, которая так расстроила Эйнштейна (который оказался неправ насчет нее).
Габор (1946) был первым, кто обнаружил, что принцип неопределенности имеет отношение к теории информации и передаче сигнала. Благодаря корпускулярно-волновому дуализму, сигналы ведут себя точно так же, как и квантовые системы. В результате, точное время и частота сигнала не могут быть известны одновременно: сигнал невозможно начертить в виде точки на временно-частотной плоскости. Эта неопределенность является свойством сигналов, а не ограничением, накладываемым математикой.
Принцип неопределенности Гейзенберга описывается через стандартное отклонение позиции ox, стандартное отклонение момента a, и постоянную Планка h:
Другими словами, произведение неопределенностей позиции и момента не может быть менее некой ненулевой постоянной. Для сигналов, нам не нужна постоянная Планка для определения соотношения квантовых величин, но форма уравнения остается такой же. Если обозначить стандартные отклонения времени и частоты как, соответственно, ot и a,, тогда мы можем переписать принцип неопределенности Габора в виде:
'В момент времени, когда позиция определена, то есть, когда фотон рассеян у электрона, электрон испытывает прерывное изменение момента. Это изменение обратно пропорционально длине волны света; при наименьшей длине волны расположение частицы определяется с наибольшей точностью. Таким образом, во время, когда расположение электрона известно, его момент может быть определен лишь до величины, при которой происходит прерывное изменение. Соответственно, чем точнее определена позиция электрона, тем менее точно может быть определен его момент. Верно и обратное.'
Heisenberg (1927), стр. 174-5.
Произведение стандартных отклонений времени в миллисекундах, и частоты в герцах должно равняться как минимум 80 (мс·Гц, или миллициклов). Мы позже вернемся к этой неопределенности. Если мы можем количественно описать неопределенность с такой точностью, в чем же проблема? Отчасти она заключается в том, что понятие, обычно описываемое словом «разрешение» состоит как минимум из двух частей: частота дискретизации и локализация. Эти тесно связанные понятия существуют как во временной, так и в частотной области.
Самым важным аспектом принципа неопределенности является то, что хотя он изначально был выражен через наблюдения и измерения, он не является следствием каких-либо ограничений измерительного оборудования или математических приемов, с помощью которых мы описываем результаты. Принцип неопределенности * hallmt@conocophillips.com
Дискретизация
Живя в цифровом веке, мы все знакомы с понятием дискретизации. Оно аналогично разрешению в цифровой фотографии, которое обычно измеряется в точках на дюйм (Рис.1), и описывает, насколько детально наше знание чего-либо.
© 2006 EAGE
43
техническая статья
first break том 24, декабрь 2006
Как при построении цифровых изображений дискретизации подвергается пространство, а цифровых сигналов – время, так цифровые спектры требуют дискретизации частот. Интуитивно можно предположить, что дискретизация частот Af зависит от дискретизации времени At. На самом деле, временная дискретизация определяет максимальную полосу пропускания (частота Найквиста равняется половине частоты дискретизации), но не имеет влияния на интервал между замерами в частотной области. Единственное, что контролирует частоту дискретизации при преобразовании Фурье – длина окна анализа (Рис.2). Конкретнее, интервал между замерами частоты Af равняется 1/T, где T – длина окна в секундах. Окно в 100 мс означает, что мы получаем информацию о частоте через каждые 1/0.100 = 10 Гц; если же сократить окно до 40 мс, интервалы становятся значительно больше: 1/0.040 = 25 Гц.
В чем же практическое значение, скажем, скверной дискретизации частоты? Во-первых, она может ограничить нашу возможность аккуратного определения частоты сигнала. Истинная частота любого сигнала находится, скорее всего, между двумя замерами (Рис.3). Различия лежат в пределах плюс или минус одной четвертой интервала, или 1/4T, и таким образом могут быть очень большими при коротком временном окне. Например, в Рис.3 погрешность составляет около ±4 Гц.
Помимо ухудшения точности оценки частоты одного сигнала, низкая частота дискретизации также затрудняет возможность определения присутствия в сигнале двух
или более различных частот. Продолжая аналогию с временной областью, или цифровыми изображениями, возможность такого разделения – то, что мы обычно называем «разрешением». Рис.4 демонстрирует, что для того, чтобы можно было аккуратно различить частоты двух сигналов, они должны быть хотя бы вдвое выше, нежели частота замеров. Интуиция подсказывает, что локальные максимумы должны быть разделены минимумом. Если они слишком близко, их сложно выделить в спектре (хотя точные результаты зависят от функции окна и реальных частот). 64-миллисекундное окно Ханна, показанное в Рис.4, вызывает дискретизацию в 15.6 Гц; мы не можем использовать преобразование Фурье для обнаружения двух сигналов, различающихся менее чем на 32 Гц. Подобные сигналы выглядят как один широкий пик в спектре.
А как же методы, не использующие преобразование Фурье? Другие алгоритмы спектрального разложения приобретают все большую популярность, особенно Sпреобразование и различные преобразования элементарных волн (Castagna and Sun, 2006). К несчастью, жаргон и математика, стоящая за этими методами, часто мешают понять их суть. Я рекомендую всем интересующимся книгу Hubbard (1998); она очень хорошо написана и должна быть доступна любому ученому.
Рис.1. Разрешение, локализация в цифровой фотографии. (a) Изображение с максимальным разрешением: 350 точек на дюйм (dpi), 512 x 512 пикселей. Различаются чрезвычайно мелкие детали. (b) То же изображение после «размывки по Гауссу» 25 x 25 пикселей. (c) То же изображение, где разрешение уменьшено фильтром 10 x 10 пикселей до 51 x 51, или 35 dpi. (d) Изображение после размывки (b) и уменьшения разрешения. Дискретизация сглаженного изображения не имеет такого эффекта, как с оригинальным изображением.
44
Рис.2. Дискретизация и локализация. (a) Сигнал в 128 мс, 32 Гц, прошедший через окно 256 мс. (b) Его спектр, рассчитанный с помощью быстрого преобразования Фурье (FFT); короткое окно означает, что спектр снят с большим интервалом замеров, т.к. Af = 1/0.256 = 3.9 Гц. (c) Тот же самый сигнал (128 мс, 32 Гц), окно 1024 мс. (d) Его спектр; длинное окно улучшило частоту дискретизации (Af = 1/1.024 = 0.977 Гц), но неопределенность остается такой же, судя по широкому центральному пику в спектре. (e) Сигнал в 768 мс, 32 Гц, окно 1024 мс. (f) Его спектр, с такой же частотой дискретизации, что и (d), но с гораздо лучшей локализацией. Адаптировано из Stockwell (1999).
© 2006 EAGE
first break том 24, декабрь 2006
Рис.3. Измерение одиночной частоты. (a) Три сигнала в 64 мс; частота дискретизации – 1000 Гц; интервал замеров 1/0.064 = 15.6 Гц. Сигналы – синусоиды с частотой 43 Гц (синий), 47 Гц (красный), и 51 Гц (зеленый). Были использованы окна Ханна. (b) Спектры чрезвычайно похожи, и, скорее всего, не могут быть различены, особенно в присутствии шума.
Важным аспектом этих алгоритмов с точки зрения понимания разрешения является то, что они – лишь модификации анализа Фурье, отличаясь от него лишь анализирующей функцией. В анализе Фурье и Sпреобразовании, это функция Гаусса; в анализе элементарных волн это набор волновых цугов. Стоить отметить, что один из этих наборов (часто используемые волны Morlet) использует функцию Гаусса в качестве основной формы. Отличие этих методов от анализа Фурье заключается не в математике, а в способе, которым анализирующая функция – окно – применяется к сигналу. Иногда утверждается, что метод элементарных волн не использует окна: это не совсем так.
Дискретизация в частотной области не очень хорошо исследована для альтернативных способов. Тем не менее, S-преобразование, преобразования элементарных волн и matching pursuit объединяет зависимость разрешения от частоты (интервал между замерами зависит от частоты сигнала). Это проще всего объяснить на примере Sпреобразования (Stockwell et al., 1996) – по сути, преобразования Фурье, для которого длина окна изменяется в соответствии с частотой сигнала, становясь короче для высоких частот и длиннее – для низких. Это означает, что высокие частоты получают худшую дискретизацию, а низкие частоты – лучшую (Рис.5). Преобразования элементарных волн используют более сложный (и трудоемкий) подход, но основной принцип остается тем же. Эти схемы часто называются «адаптивными» или, более описательно, «много-разрешительными».
В типичной интерпретации сигнала с использованием метода элементарных волн, мы получаем сигнал на частоте Найквиста fw другой на частоте fN/2, затем fN/4, fN/8, и так далее. (Это на самом деле не так просто: элементарные волны не обладают точными частотами). Таким образом, для сейсмического сигнала, снятого при 500 Гц (через каждые 2 мс), частота Найквиста – 250 Гц, а другие частоты – около 250, 125, 63, 32, 16, 8, 4, 2, и 1 Гц. Высочайшая частота при этом будет около 0.8 x Найквист = 200 Гц, и наименьшие частоты – в районе 10 Гц. Это означает, как минимум, пять полезных коэффициентов (и также объясняет, почему этот метод так полезен для сжатия данных: сигнал может быть воспроизведен из относительно малого количества данных).
техническая статья
Локализация
Менее конкретная, нежели частота дискретизации, локализация все-таки знакома геофизикам. Широко говоря, она аналогична фокусу в цифровой фотографии: измерение степени «размазанности» информации. Стоит отметить, что это отличается от разрешения и частоты дискретизации: данные могут быть размазаны, несмотря на короткий интервал между замерами (Рис. 1b и 2b).
Локализацию во временной области легко понять интуитивно: она зависит от широты полосы пропускания, эффектов настройки, и стабильности данных. Она ограничивает аккуратность и постоянность интерпретации пика отражаемости сейсмических данных.
Локализация в частотной области зависит от утечки. При проведении преобразования Фурье любого конечного сигнала его спектр демонстрирует феномен, известный, как «утечка». Так, хорошо известно, что преобразование Фурье, проведенное на пачке импульсов (квадратном окне) – гармоническая функция бесконечной длины (Рис.6). Бесконечные колебания этой функции, и тот факт, что ее пик не острый, а скругленный и широкий, и составляют утечку. Ее вызывает выбор части данных для анализа, так что она не является артефактом собственно анализа Фурье.
Как уже было сказано, локализация в частотной области характеризуется боковыми спектральными колебаниями и шириной пика. Рассмотрим каждое из этих свойств по отдельности.
Изменение окна уменьшает колебания. Спектр квадратной функции, равно как и любого сигнала с использованием квадратного окна, обладает бесконечными боковыми колебаниями; тем не менее, они сильно сокращаются при использовании конической функции (Рис.6). Этот процесс называется аподизация, что буквально означает «лишение ног». Наиболее часто используемые аподизующие функции – Гаусса, Ханна, и Хемминга, но существуют и многие другие. Как и в других аспектах обработки сигнала, за уменьшение колебаний приходится платить: спектральный пик уменьшается и расширяется. Расширение и уменьшение амплитуды зависит от конкретной выбранной функции.
Рис.4. Анализ двух сигналов. (a) Сигнал в 64 мс; частота дискретизации 1000 Гц, интервал между частотами 1/0.064 = 15.6 Гц. Сигнал состоит из синусоид 32 Гц и 63 Гц; использовано окно Ханна. (b) В спектре две частоты неразличимы. (c) Сигнал 32 Гц и 65 Гц, также с окном Ханна. (d) Сигналы различимы, так как разница между ними больше, чем два частотных интервала.
© 2006 EAGE
45
техническая статья
first break том 24, декабрь 2006
Рис.6. Leakage. (a) Относительно длинный сигнал и (b) его спектр. (c) Тот же сигнал с наложением квадратной функции и (d) его спектр, демонстрирующий боковые колебания и увеличенную ширину пика. (e) Тот же сигнал с наложением функции Ханна, и (f) его спектр, без колебаний, но с еще более широким и низким пиком.
Рис.5. Логоны Габора объясняют концепцию локализации во времени и частоте с помощью изображения идеализированной частотно-временной плоскости. (a) В так называемой системе Дирака логоны превращаются в бесконечные вертикальные линии. (b) С другой стороны, в системе Фурье спектр не привязан ко времени. (c) Преобразование Фурье с коротким временным окном вызывает дурную локализацию в частотной области. (d) Длинное окно предоставляет лучшую частотную локализацию за счет ухудшенной локализации во временной области. (e) Преобразования элементарных волн и S-преобразование дают локализацию, зависящую от частоты: низкие частоты хорошо локализованы в частотной, но не во временной области; для высоких частот верно обратное.
Ширина спектрального пика – второй ключевой компонент локализации. Она изменяется в зависимости от выбранного окна, но в первую очередь зависит от самого сигнала. Она приблизительно равняется 2/S, где S – длина ненулевой части сигнала, а не длина окна. В случае с сейсмическими отражениями, мы имеем дело с пучком перекрывающихся сигналов, представляющих отражения с подповерхностного слоя. Длина S каждого индивидуального сигнала – из тех, что мы пытаемся описать – это длина сейсмического цуга. Чему же она равняется? В теории – бесконечности; на самом деле в какой-то момент уровень энергии падает ниже уровня шума. Опыт показывает, что цуги обычно 64-128 мс длиной. Так что для типичных сейсмических данных, ширина спектрального пика равняется примерно 16-31 Гц, с большей неопределенностью для более коротких волн.
Конечно, если окно короче сейсмического цуга, пик становится еще шире, т.к. максимум не может быть шире двух интервалов дискретизации. Таким образом, не стоит делать окно уже, чем сейсмический цуг.
Познавательным является изучение эффекта увеличения длины окна на спектр короткого сигнала, и его сравнение с окном, наложенным на длинный сигнал (Рис.2). Мы находим, что короткие сигналы обладают неважной локализацией в частотной области. Размывание частоты не зависит от выбора окна: это естественное свойство.
Измерение неопределенности
Соотношение между разрешением во временной и частотной областях проиллюстрировано в Рис.5 с помощью «квадратов Гейзенберга». Габор назвал эти квадраты «логоны», что является более удачным названием, так как они описывают принцип неопределенности Габора. Эти идеализированные схемы демонстрируют, как частотно-временная плоскость может быть «вымощена» прямоугольниками размером at x a,, где at – стандартное отклонение времени и a, - стандартное отклонение частоты. Как мы видели раньше, их произведение должно быть как минимум \IAn. На практике, нижняя граница обычно не достигается, хотя некоторые преобразования подходят ближе к ней, нежели другие. Интересно было бы увидеть публикации об областях логонов Габора, которыми их алгоритмы «мостят» частотно-временную плоскость.
Рассмотрим конкретный пример неопределенности Габора. Если мы знаем временную локализацию до ±5 мс, тогда локализация в частотной области ±0.005 x An= ±16 Hz. Узнать частоту сигнала с большей точностью, при этом не теряя в разрешении временной области,
46
© 2006 EAGE
first break том 24, декабрь 2006
невозможно. Подобные неопределенности важны в интерпретации сейсмических данных, особенно при использовании количественных спектральных методов.
Эффекты изменения дискретизации также возможно измерить. Представьте себе, что вы пытаетесь использовать спектральное разложение для количественного анализа тонких слоев в сейсмических данных. Как было упомянуто ранее, длина окна будет определять интервал между замерами частот. Так, если вы выберете окно в 50 мс, тогда у вас будет информация об ответе при 20 Гц, 40 Гц, 60 Гц, и так далее. Любая другая информация (то есть, другие срезы в так называемом «кубе настройки»), будет неизвестна с окном этой длины и, таким образом, интерполируется. Окно в 20 мс, которое предлагает улучшенное разрешение во времени, предоставляет спектральную информацию о 50 Гц, 100 Гц, 150 Гц, и т.д. В попытках описать эти тонкие слои, вы можете остаться всего лишь с тремя или четырьмя точками внутри линии пропускания. Соответственно, попытки описать толщину слоя не пойдут далее определений «тонкий», «толстый» и «очень толстый». Скорее всего, вы могли различить «толстый» и «очень толстый» слои, и не прибегая к методу спектрального разложения.
Заключение
Спектральный анализ сигналов управляется принципом неопределенности Габора, выведенного из неравенства Коши-Шварца и аналогичного принципу неопределенности Гейзенберга. Он не ограничивает предел возможностей измерения, а, скорее, описывает фундаментальный аспект реальности: сигналы не имеют произвольно точной локализации времени и частоты. Неважно, как рассчитывается спектр частот – если вам нужна информация о времени, за нее придется платить потерей информации о частоте. Конкретнее, произведение неопределенностей во временной и частотной областях будет как минимум 1/4n.
Для хорошей локализации во временной области, используйте короткое окно T для метода Фурье (для метода элементарных волн, используйте короткие цуги с высокими частотами). Таким образом, вы увеличите интервалы дискретизации частот и ухудшите разрешение в частотной области. Из-за уменьшения частоты дискретизации станет сложнее точно определить частоту сигнала (±1/2T), и разделить похожие частоты (ближе, чем 2/T).
Для хорошей локализации в частотной области, используйте длинное окно (или анализируйте с длинными цугами). При этом вы снизите интервал дискретизации частот и получите лучшее разрешение в частотной области за счет уменьшения локализации во временной.
Вне зависимости от того, используете ли вы окна Фурье или элементарные волны, короткие сигналы обладают широкими спектральными пиками и, таким образом, неважной локализацией в частотной области. Лишь долговременные сигналы дают узкие спектральные пики и, соответственно, хорошую частотную локализацию. В этом суть принципа неопределенности Габора.
техническая статья
Ссылки
Castagna, J., Sun, S., and Siegfried, R. [2003] Instantaneous spectral analysis: Detection of low-frequency shadows associated with hydrocarbons. The Leading Edge, 22, 127-129. Castagna, J. and Sun, S. [2006] Comparison of spectral decomposition methods. First Break, 24 (3), 75-79. Gabor, D. [1946] Theory of communication. Journal of the Institute of Electrical Engineering, 93, 429-457. Heisenberg, W. [1927] Über den anschaulichen Inhalt der quantentheoretischen Kinematik und Mechanik, Zeitschrift für Physik, 43, 172-198. English translation: Quantum Theory and Measurement, J. Wheeler and H. Zurek [1983]. Princeton University Press, Princeton. Hubbard, B. [1998] The world according to wavelets: the story of a mathematical technique in the making. A.K. Peters, Natick. Marfurt, K. and Kirlin, R. [2001] Narrowband spectral analysis and thin-bed tuning. Geophysics, 66, 1274-1283. Partyka, G., Gridley, J., and Lopez, J. [1999] Interpretational applications of spectral decomposition in reservoir characterization. The Leading Edge, 18, 353-360. Stockwell, R., Mansinha, L., and Lowe, R. [1996] Localization of the complex spectrum: the S transform. Institute of Electrical and Electronic Engineers Transactions on Signal Processing, 44,
998-1001.
Stockwell, R. [1999]. S-transform analysis of gravity wave activity from a small scale network of airglow imagers. PhD thesis, University of Western Ontario, London, Canada.
Благодарности
Мы выражаем свою признательность Tooney Fink, Sami
Elkurdi и Peter Rowbotham за их вклад, и анонимному
рецензенту за несколько важных улучшений.
© 2006 EAGE
47