Способ подготовки высокодетальных гравиметрических карт
METHOD FOR PREPARATION OF HIGHLY DETAILED GRAVIMETRIC MAPS

Изобретение относится к области гравиметрии и может быть использовано при оперативной подготовке карт в интересах разведки полезных ископаемых, подготовки строительства больших гидросооружений, подготовки карт для обеспечения навигации по гравитационному полю Земли (ГПЗ) и т.д.
Известен способ подготовки высокодетальных гравиметрических карт, который основан на измерении параметров гравитационного поля Земли в узлах сетки размерами от первых десятков метров до пяти километров [1]. Измерения выполняются с помощью абсолютных и относительных гравиметров (зарубежных Scintrex CG-5 Autograv, Micro-g LaCoste A-10, отечественных - ГНУ-КВ, ГАБЛ), а также гравитационных градиентометров, например Е-60. В результате создаются карты масштабов от 1:5 000 до 1:500 000.
Способ обладает следующими недостатками:
- низкая оперативность, поскольку для обеспечения создания карт, например масштаба 1:200 000 (расстояние между узлами сетки 1 км), на площади 100×100 км необходимо произвести не менее 10000 измерений, для чего потребуется не менее 1000 приборо-смен, продолжительностью от 4 часов до 8 часов каждая;
- высокая стоимость полевых работ, поскольку для создания карт масштаба 1:200 000 и крупнее требуется от 50 млн.руб. и более в ценах 2021 года;
- способ неприменим для создания карт на недоступных территориях.
Известен способ подготовки гравиметрических карт, основанный на использовании моделей ГПЗ, по которым в узловых точках вычисляются значения параметров гравитационного поля Земли. В частности, наиболее совершенная модель ГПЗ EGM2008 [2] имеет размер сетки узловых точек около 20×20 км.
Недостатком такого способа является низкая точность и детальность создаваемых карт, что не удовлетворяет требованиям потребителей, нуждающихся в высокоточных картах детальностью 1 км и менее (например, прикладная геодезия, поиск и разведка полезных ископаемых и др.).
Известен способ-прототип, основанный на интерполяции значений параметров ГПЗ, известных для карт с низкой детальностью [3]. Например, по известной гравиметрической карте, подготовленной по модели ГПЗ с детальностью значений ускорения свободного падения (УСП) и гравитационных градиентов (ГТ) 20 км (масштаб 1:2 000 000), путем интерполяции вычисляют значения УСП и ГТ в промежуточных точках, отстоящих на расстоянии 1 км и менее. На основе этого можно создать карты обладающие большей детальностью, например, более мелкого масштаба 1:200 000.
Сущность способа-прототипа поясняется фиг. 1 на примере карт УСП. Способ включает следующие операции:
N
N+1
1. Определяют по модели ГПЗ значения УСП или гравитационного градиента в узловых точках сетки размерами 20 км (значения g и S в соответствии с фиг. 1).
2. Интерполируют вычисленные значения УСП на промежуточные точки карты, отстоящие друг от друга, например, на 1 км по формуле:
где  - интерполированное значение УСП,
N
N+1
g и g - известные из модели значения УСП или гравитационного градиента
в крайних точках профиля интерполяции,
L - длина профиля интерполяции,
- интервал интерполяции.
3. Результаты интерполяции сводят в сетку карты нового масштаба с размерами ячейки сетки 1 км (масштаб 1:200 000).
Недостатком способа-прототипа является низкая точность, поскольку значительные погрешности в интерполяцию вносят неопределенности в знании рельефа и плотности пород на дистанции между исходными узловыми точками. В зависимости от типа местности (равнинная, холмистая или горная) и аномальности гравитационного поля (слабо- и сильно аномальная) погрешность способа-прототипа может достигать значений от 3 до 9 мГал и более, что является практически невостребованной величиной для всех типов потребителей.
Техническим результатом изобретения является повышение точности создания гравиметрических карт.
Сущность предлагаемого способа подготовки высокодетальных карт УСП и гравитационных градиентов состоит в следующем:
1. при построении карт ускорения свободного падения на основе известных его значений в узловых точках модели дополнительно используют данные о координатах узловых точек пространственной сетки заданной детальности, данные о модели рельефа местности, данные о модели плотности подстилающих пород;
2. по всем этим данным для всех точек пространственной сетки заданной детальности вычисляют аномальные поправки в значения ускорения свободного падения нормального поля, вызванные влиянием геодезической высоты точек над референц-эллипсоидом, влиянием промежуточного слоя пород между поверхностью референц-эллипсоида и физической поверхностью Земли;
3. для узловых точек модели определяют значения аномальной составляющей, вызванной влиянием неизвестной неоднородностью плотности пород под поверхностью референц-эллипсоида и интерполируют эти значения на точки высокодетальной пространственной сетки;
4. по результатам интерполяции уточняют вычисленные значения ускорения свободного падения для всех точек высокодетальной пространственной сетки;
5. для точек высокодетальной пространственной сетки, высота которых превышает заданную, по данным об их координатах, данным о модели рельефа местности и модели плотности пород вычисляют значения гравитационных градиентов;
6. по результатам вычисления значений ускорения свободного падения и гравитационных градиентов для узловых точек высокодетальной пространственной сетки подготавливают карту параметров ГПЗ с заданной детальностью.
Таким образом, заявленный технический результат изобретения, состоящий в повышении точности создания гравиметрических карт, а также возможности подготовки карт на недоступные территории, достигнут:
- при подготовке карт УСП на доступные территории выигрыш в точности по отношению к прототипу составляет от 3,3 до 0,7 мГал (до 5 раз);
- при подготовке карт УСП на недоступные территории выигрыш в точности по отношению к прототипу составляет от 3,3 до 1,9 мГал (до 1,5 раз);
- при подготовке карт гравитационных градиентов выигрыш в точности по отношению к прототипу составляет от 19 до 6 Этвеш (до 3 раз).
Для подготовки карт УСП и гравитационных градиентов для недоступных территорий необходима дополнительная информация о моделях (картах) УСП, о рельефе (модель SRTM) и о плотности пород (модель TopoDensT), которая имеется на открытых ресурсах https://www.usgs.gov и https://www.unb.ca.
Изобретение поясняется графическими изображениями и таблицами:
Фиг. 1 - Схема способа-прототипа;
Фиг. 2 - Схема предлагаемого способа;
Фиг. 3 - Схема расположения экспериментальных профилей в Московской области;
Фиг. 4 - Результаты оценки погрешности предлагаемого способа создания карт УСП;
Фиг. 5 - Разность  и  значений УСП для варианта труднодоступной территории;
Фиг. 6 - Схема эксперимента для проверки величины аномальной составляющей  на различных высотах;
Фиг. 7 - Результаты экспериментальной оценки предельной высоты влияния составляющей ;
xy
Δ
Фиг. 8 - Схема района аэро-градиентометрической съемки: а) -карта рельефа, б) - составляющих W, W;
Фиг. 9 - Результаты оценки погрешности способа-прототипа;
Фиг. 10 - Результаты оценки погрешности вычислений. Как следует из отчета аэро-градиентометрической съемки, зоны высоких значений  и  связаны исключительно с высокой помеховой обстановкой во время съемки (турбулентность): а) карта разности . СКО разности 6 Этвеш, среднее значение  Этвеш; б) карта разности. . СКО разности 6 Этвеш, среднее значение  Этвеш.
1. Обоснование точности предлагаемого способа
1.1. Обоснование операций способа повышения точности создания карт ускорения свободного падения (УСП)
Сущность предлагаемого способа создания карт УСП поясняется фиг. 2.
Ускорение свободного падения g в точке на поверхности Земли представим в виде [4]:
0
где γ - нормальное гравитационное поле Земли,  - аномальная поправка за счет геодезической высоты точки,  - аномальная поправка за счет влияния промежуточного слоя пород, заключенных между поверхностью земного эллипсоида и физической поверхностью Земли,  - аномальная составляющая УСП, вызванная влиянием неоднородности плотности пород под поверхностью эллипсоида.
0
Величина γ определяется по известной формуле:
e
p
где g- значение УСП на экваторе,  g - значение УСП на полюсе, ϕ - широта пункта наблюдений,  α - сжатие общеземного эллипсоида.
Величина  вычисляется по формуле:
где h - геодезическая высота точки.
Аномальная поправка за счет влияния промежуточного слоя  в общем виде определяется по формуле:
i
i
i
где G - гравитационная постоянная, ρ - плотность элементарной ячейки пород промежуточного слоя, ν - объем элементарной ячейки пород промежуточного слоя, r - расстояние до элементарной ячейки пород промежуточного слоя.
Δρ
uзм
Аномальная поправка Δg неизвестна, однако в i-й точке с известным значениям g (взятым из карты, либо измеренным) эту величину можно вычислить на основе формулы (2):
где сумма  вычисляется по вышеперечисленным формулам (3), (4), (5).
выч
Слагаемое  вычисляется по формуле (5) по известной информации о высотах рельефа местности и значениях плотности пород вокруг точки g. Эти данные находятся в свободном доступе в виде цифровых моделей рельефа, например модель Shuttle radar topographic mission (SRTM) [5] и цифровых моделей плотности пород, например модель TopoDensT [6].
Информация о высотах рельефа можно представить в виде набора прямоугольных призм. Поэтому составляющая  в формуле (5) вычисляется следующим образом:
i
1
2
1
2
1
2
где G - гравитационная постоянная, ρ - плотность пород, ξ и ξ - координаты углов призмы по оси X, η и η - координаты углов призмы по оси Y, ζ и ζ - координаты углов призмы по оси Z. , где х, у и z - координаты точки вычислений.
N
N+1
Для двух соседних точек g и g по формулам (3), (4), (7) и (6) определяются значения  и , на основе которых проводится интерполирование в промежуточную i-ю точку по формуле:
где  - интерполированное значение аномальной составляющей УСП,  и  - вычисленные по формуле (6) аномальные составляющее УСП в крайних точках профиля интерполяции, L - длина профиля интерполяции  - интервал интерполяции.
В результате, значение УСП в промежуточной i-й точке  определяется по формуле:
1.2. Оценка точности способа повышения точности создания карт УСП
Способ, по сравнению с прототипом, обладает более высокой точностью поскольку учитывает изменения рельефа и плотности пород между узловыми точками исходной карты.
N
N+1
Оценка повышения точности предлагаемого способа в части создания карт УСП на доступных территориях приведем на примере трех экспериментальных профилей в Московской области (фиг. 3). Для каждого профиля расстояние L между крайними точками g и g составило ≈5 км, количество промежуточных точек равно пяти, расстояние  между промежуточными точками 800-900 м,. Итого для каждого измерительного профиля интерполяции получен набор из семи точек с измеренными координатами, каждой из которых соответствует значение . Измерения выполнялись при помощи гравиметра CG-5 Autograv, координатная привязка вьшолнялась с помощью геодезического приемника Javad Sigma.
Погрешность измерений  составила 0,005 мГал.
Сравнение способа прототипа и предлагаемого способа выполнено по разности  и , соответственно. Результаты оценки этих разностей, которые представляются как погрешности интерполяции, представлены в таблице (столбцы 6 и 7, фиг. 4.)
Как видно из таблицы, погрешность интерполяции предлагаемого способа (столбец 7) составила от 0,1 до 0,7 мГал, что в более чем в 5 раз лучше, чем при использовании способа прототипа, для которого погрешность интерполяции составляет от 0,1 до 33 мГал (столбец 6).
Способ применим и для создания карт УСП на недоступные территории. Проверка предлагаемого способа для создания карт УСП на труднодоступные территории выполнена на примере трех экспериментальных профилей в Московской области. Для определения значений  в промежуточных точках профилей, высоты этих точек определялись по цифровой модели рельефа.
Результаты оценки погрешности интерполяции  представлены в таблице, представленной на фиг. 5.
Из таблицы (фиг. 5, столбец 2) следует, что погрешность созданных карт УСП не превышает 1,9 мГал, тогда как в способе-прототипе для этого же участка погрешность интерполяции достигает 3,3 мГал (фиг. 4, столбец 6).
1.3. Обоснование операций способа повышения точности создания карт гравитационных градиентов
αβ
Гравитационный градиент W в любой точке поверхности Земли можно представить в виде суммы [4]:
где α,β принимают значения x,y,z;  - составляющие гравитационных градиентов нормального гравитационного поля, создаваемого земным эллипсоидом,  - аномальные составляющие гравитационных градиентов, вызванные влиянием пород промежуточного слоя, заключенного между поверхностью земного эллипсоида и физической поверхностью Земли,  - аномальные составляющие гравитационных градиентов, вызванные влиянием плотностных неоднородностей под поверхностью эллипсоида.
Δ
xy
Применительно к составляющим гравитационного градиента уравнение (10) на примере W и W представим в виде:
Составляющие гравитационных градиентов нормального гравитационного поля Земли находятся по известным формулам:
0
где γ - нормальное гравитационное поле Земли, вычисляемое по формуле (12); М и N - радиус кривизны меридиана и первого вертикала соответственно.
Эти составляющие вычисляются по координатам наземной точки с погрешностью тысячные доли Этвеш.
Слагаемые , ,  в формуле (11) в общем виде вычисляются по формуле:
i
i
i
где G - гравитационная постоянная, ρ - плотность элементарной ячейки пород промежуточного слоя, ν - объем элементарной ячейки пород промежуточного слоя, r - расстояние до элементарной ячейки пород промежуточного слоя.
Учитывая форму представления информации о рельефе и плотности пород, составляющая  в формуле (11) вычисляется в виде определенного интеграла при известных размерах призм:
1
2
1
2
1
2
где G - гравитационная постоянная, ρ - плотность пород, ξ и ξ - координаты углов призмы по оси X, η и η - координаты углов призмы по оси Y, ζ и η - координаты углов призмы по оси Z. , где х, у и z координаты точки вычислений.
Опубликованных сведений о величине аномальных составляющих  нет. Выполненная экспериментальным путем оценка величины аномальной составляющей  (фиг. 6, 7) показала, что при высоте 250 м и выше над эллипсоидом (геодезическая высота) влиянием плотностных неоднородностей под поверхностью эллипсоида (т.е. членом ) можно пренебречь.
При соблюдении вышеуказанных условий, формула (11) упрощается:
Δ
xy
Таким образом, для подготовки карт гравитационных градиентов составляющих W и W необходимо использовать формулы (12), (13) и (14) - (18) при геодезической высоте в районе подготовки карты 250 м и выше.
1.4. Оценка точности способа повышения точности создания карт гравитационных градиентов
Проверка применимости предлагаемого способа создания карт гравитационных градиентов на недоступные территории выполнена на примере опубликованных результатов воздушных геофизических исследований для территории штата Миссури, США [7], представленных на фиг. 8. Для указанной территории доступны по публикациям координаты и высоты измерительных профилей, а также измеренные значения градиентов  и .
Оценка погрешности интерполяции способа-прототипа на примере указанной территории выполнена по разности  и  Результаты оценки погрешности интерполяции составили 19 и 13 Е для  и  соответственно (фиг. 9 столбцы 6 и 9).
Оценка погрешности предлагаемого способа выполнена на примере указанной территории выполнена по разности  и  для всех измерительных профилей. Значения  и  не превысили 6 Этвеш (фиг. 10).
Таким образом, выигрыш в точности предлагаемого способа подготовки карт гравитационных градиентов составляет более 3 раз (6 Этвеш против 19 Этвеш в прототипе).
7.5. Приведенные обоснования можно представить в виде последовательности действий:
В части создания карт УСП
1. По известным моделям (картам) ГПЗ вычисляют сетку значений УСП на заданную территорию поверхности Земли с детальностью, которую обеспечивает выбранная модель (карта).
2. По координатам точек, модели рельефа местности и модели плотности пород по формулам (3), (4), (7) и (6) вычисляют аномальную поправку  в узловых точках сетки.
3. В соответствии с заданной детальностью, например, 1 км выбирают промежуточные точки сетки, в которые интерполируют аномальную поправку  по формуле (8).
4. В промежуточных точках сетки по координатам точек, модели рельефа местности и модели плотности пород по формулам (3), (4), (7) и (9) вычисляют значения УСП  и создают карту с детальностью 1 км.
В части создания карт гравитационных градиентов
1. Выполняют оценку величин геодезической высоты на заданной территории. Убеждаются, что значения высот больше 250 м.
Δ
xy
2. В соответствии с заданной детальностью выбирают промежуточные точки сетки, для которых по координатам точек, модели рельефа местности и модели плотности пород по формулам (12), (14-17) и (11) вычисляют значения гравитационных градиентов W и W.
Δ
xy
3. Создают карту гравитационных градиентов W и W с детальностью 1 км.
Источники информации
1. Инструкция по гравиразведке, М.: Недра, 1980 г.
2. N.K. Pavlis, S.A. Holmes, S.C. Kenyon, J.K. Factor // The development and evaluation of the Earth Gravitational Model 2008 (EGM2008). Journal of Geophysical Research, vol. 117, B04406, DOI: 10.1029/2011JB008916, 2012.
3. A. Lindau // A map based alternative for the determination of the local gravity, PTB, Braunschweig, Germany.
4. Огородова Л.В., Шимбирев Б.П., Юзефович А.П. Гравиметрия. М.:, «Недра», 1978, 325 с.
5. https://www.usgs.gov/centers/eros/science/usgs-eros-archive-digital-elevation-shuttle-radar-topography-mission-srtm-1-arc?qt-science_center_obiects=0#qt-science_center_objects
6. https://www.unb.ca/fredericton/engineering/depts/gge/resources.html
7. GEOPHYSICAL SURVEY REPORT Helifalcon™ AIRBORNE GRAVITY GRADIOMETER AND MAGNETIC STINGER SURVEY SULLIVAN NORTH, MISSOURI PROJECT 14012 USGS. 2014.