GDAL SPOT-6 VRT

Ответить
juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

GDAL SPOT-6 VRT

Сообщение juffin_h » 04 фев 2015, 10:38

Есть снимок со спутника SPOT-6. Хотел его ортотрансформировать с помощью gdalwarp, но оказалось, что GDAL не “понимает” снимки SPOT-6.
Хочу попробовать RPC из метаданных снимка, преобразовать в формат VRT и подсунуть gdalwarp.
Такое возможно?
  Если да, то может ли кто поделится VRT файлом с RPC доменом для примера? В сети подобного VRT не нашел.
  Если нет, то направьте на путь истинный.
Последний раз редактировалось juffin_h 09 фев 2015, 10:38, всего редактировалось 3 раза.

alexandr cherepanov
Гуру
Сообщения: 534
Зарегистрирован: 30 ноя 2006, 13:31
Репутация: 116
Откуда: Moscow

Re: GDAL SPOT-6 VRT

Сообщение alexandr cherepanov » 04 фев 2015, 14:29

Код: Выделить всё

<VRTDataset rasterXSize="9720" rasterYSize="66383">
  <SRS>PROJCS["WGS 84 / UTM zone 37N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",39],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32637"]]</SRS>
  <GeoTransform> 4.2771600000000000e+005, 2.3999999999999999e+000, 0.0000000000000000e+000, 6.2590391999999993e+006, 0.0000000000000000e+000,-2.3999999999999999e+000</GeoTransform>
  <Metadata>
    <MDI key="AREA_OR_POINT">Area</MDI>
    <MDI key="TIFFTAG_XRESOLUTION">1</MDI>
    <MDI key="TIFFTAG_YRESOLUTION">1</MDI>
  </Metadata>
  <Metadata domain="RPC">
    <MDI key="HEIGHT_OFF">210</MDI>
    <MDI key="HEIGHT_SCALE">500</MDI>
    <MDI key="LAT_OFF">55.7571</MDI>
    <MDI key="LAT_SCALE">0.7430</MDI>
    <MDI key="LINE_DEN_COEFF">+1.000000E+00 -9.503236E-06 -6.386284E-04 +2.976162E-04 -1.904468E-07 -3.768348E-06 -3.588716E-06 +2.900618E-06 +2.049445E-07 +4.631482E-06 +4.712116E-07 +0.000000E+00 +1.548356E-08 +0.000000E+00 +5.036180E-08 +4.819471E-08 +4.843385E-08 +1.151909E-07 -4.660568E-06 +1.819432E-08</MDI>
    <MDI key="LINE_NUM_COEFF">-1.217299E-03 +1.969171E-03 -1.000925E+00 -1.019775E-04 -1.229830E-05 +5.464240E-06 +2.950766E-04 -1.953999E-04 +5.792297E-04 -1.083605E-07 -8.263748E-06 +0.000000E+00 -4.426955E-07 +1.433500E-08 -7.759899E-07 -1.435955E-07 -4.594503E-06 +2.071665E-08 -1.441813E-05 -1.001393E-08</MDI>
    <MDI key="LINE_OFF">33201</MDI>
    <MDI key="LINE_SCALE">34420</MDI>
    <MDI key="LONG_OFF">38.0486</MDI>
    <MDI key="LONG_SCALE">0.1888</MDI>
    <MDI key="SAMP_DEN_COEFF">+1.000000E+00 -2.492860E-03 +8.702597E-03 -4.031501E-04 +2.042522E-05 -2.319546E-05 +9.809980E-06 -5.156165E-05 +7.904697E-05 -9.248327E-05 -2.709752E-07 -1.701596E-08 +2.933039E-07 +3.913028E-08 -9.259787E-07 +1.708361E-06 -1.859015E-06 -1.031072E-08 +4.222594E-07 +1.142850E-07</MDI>
    <MDI key="SAMP_NUM_COEFF">+5.959196E-03 +9.446828E-01 +9.045097E-02 -2.201229E-02 -9.965023E-03 +8.334638E-04 -2.793291E-03 -2.345943E-03 +1.193017E-03 -1.723565E-05 +3.216177E-05 -4.933984E-05 -1.610468E-04 -8.584692E-05 +5.990033E-05 +8.395135E-06 -9.615944E-06 -2.591164E-05 +7.440910E-06 +1.979198E-06</MDI>
    <MDI key="SAMP_OFF">5204</MDI>
    <MDI key="SAMP_SCALE">5226</MDI>
  </Metadata>
  <VRTRasterBand dataType="UInt16" band="1">
    <Metadata />
    <ColorInterp>Red</ColorInterp>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">example.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="9720" RasterYSize="66383" DataType="UInt16" BlockXSize="9720" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
      <DstRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
    </SimpleSource>
  </VRTRasterBand>
  <VRTRasterBand dataType="UInt16" band="2">
    <Metadata />
    <ColorInterp>Green</ColorInterp>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">example.tif</SourceFilename>
      <SourceBand>2</SourceBand>
      <SourceProperties RasterXSize="9720" RasterYSize="66383" DataType="UInt16" BlockXSize="9720" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
      <DstRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
    </SimpleSource>
  </VRTRasterBand>
  <VRTRasterBand dataType="UInt16" band="3">
    <Metadata />
    <ColorInterp>Blue</ColorInterp>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">example.tif</SourceFilename>
      <SourceBand>3</SourceBand>
      <SourceProperties RasterXSize="9720" RasterYSize="66383" DataType="UInt16" BlockXSize="9720" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
      <DstRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
    </SimpleSource>
  </VRTRasterBand>
  <VRTRasterBand dataType="UInt16" band="4">
    <Metadata />
    <SimpleSource>
      <SourceFilename relativeToVRT="1">example.tif</SourceFilename>
      <SourceBand>4</SourceBand>
      <SourceProperties RasterXSize="9720" RasterYSize="66383" DataType="UInt16" BlockXSize="9720" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
      <DstRect xOff="0" yOff="0" xSize="9720" ySize="66383" />
    </SimpleSource>
  </VRTRasterBand>
</VRTDataset>

juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

Re: GDAL SPOT-6 VRT

Сообщение juffin_h » 04 фев 2015, 14:40

Спасибо большое! Буду пробовать.

juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

Re: GDAL SPOT-6 VRT

Сообщение juffin_h » 07 фев 2015, 13:48

С помощью gdalbuildvrt и файла метаданных SPOT-6 DIM*.XML создал файл vrt. В полученный vrt добавил RPC домен с данными, взятыми из файла RPC*.XML.
Далее делаю: gdalwarp –rpc in.vrt out.tif.
В результате координаты поученного растра врут. В одном углу на ~1000м, в другом на ~200м.
В чем моя принципиальная ошибка?

Аватара пользователя
Дмитрий Барышников
Гуру
Сообщения: 2572
Зарегистрирован: 17 ноя 2009, 19:17
Репутация: 261
Откуда: Москва

Re: GDAL SPOT-6 VRT

Сообщение Дмитрий Барышников » 07 фев 2015, 18:49

Вообще GDAL поддерживает rpc коэфициенты в формате DigitalGlobe. А SPOT там поддержки не было.
Для проверки что так и есть достаточно выполнить gdalinfo <путь до растра>
Вот у меня вывод для съемки с КА QuickBird вот такой и там есть секция про RPC. Если у вас не так, то ничего удивительного, что результат далек от ожидаемого.
$ gdalinfo 12JUN08182708-M2AS-052783824030_01_P001.TIL
Driver: TIL/EarthWatch .TIL
Files: 12JUN08182708-M2AS-052783824030_01_P001.TIL
./12JUN08182708-M2AS_R1C1-052783824030_01_P001.TIF
./12JUN08182708-M2AS_R1C2-052783824030_01_P001.TIF
./12JUN08182708-M2AS_R2C1-052783824030_01_P001.TIF
./12JUN08182708-M2AS_R2C2-052783824030_01_P001.TIF
12JUN08182708-M2AS-052783824030_01_P001.IMD
12JUN08182708-M2AS-052783824030_01_P001.RPB
Size is 2094, 2094
Coordinate System is:
PROJCS["WGS 84 / UTM zone 10N",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",-123],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AUTHORITY["EPSG","32610"]]
Origin = (542391.599999999976717,4185149.999999930150807)
Pixel Size = (2.400000000000000,-2.400000000000000)
RPC Metadata:
HEIGHT_OFF=43
HEIGHT_SCALE=501
LAT_OFF=37.7939
LAT_SCALE=0.0497
LINE_DEN_COEFF=+1.000000E+00 -5.924250E-06 +1.542447E-04 +3.469027E-04 +0.000000E+00 +1.836519E-06 -3.068663E-07 -7.924399E-07 +3.851062E-08 -9.499148E-07 -6.968846E-08 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +4.201451E-08 +0.000000E+00 +2.886710E-08 +0.000000E+00 +0.000000E+00
LINE_NUM_COEFF=+3.158127E-04 -1.007419E-02 -1.015209E+00 +8.584744E-03 +3.636214E-07 +4.502635E-05 +3.439088E-04 -9.018938E-04 -1.609656E-04 -3.146326E-06 +2.248350E-06 +0.000000E+00 +1.295983E-08 +2.167892E-08 +2.637776E-07 -4.052430E-08 +7.179718E-07 +7.509639E-07 -4.581548E-07 +0.000000E+00
LINE_OFF=855
LINE_SCALE=2263
LONG_OFF=-122.4301
LONG_SCALE=0.1020
SAMP_DEN_COEFF=+1.000000E+00 +3.939953E-04 +3.381574E-04 -5.806425E-04 +4.442241E-08 -7.145756E-08 +3.948624E-07 +1.140191E-08 +2.612522E-07 -3.726788E-07 +1.438509E-08 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 -1.659343E-08 -1.322187E-08 +0.000000E+00
SAMP_NUM_COEFF=-7.970981E-04 +1.008162E+00 -3.799454E-03 +9.610906E-03 -3.360840E-04 +5.825505E-04 -2.467101E-04 +3.994724E-04 -3.389467E-06 +6.916883E-06 -3.867707E-07 +1.494673E-07 -3.456001E-07 +2.320039E-07 -2.386330E-07 +0.000000E+00 +4.243663E-08 +5.620114E-07 +1.257254E-07 +0.000000E+00
SAMP_OFF=3250
SAMP_SCALE=3712

Corner Coordinates:
Upper Left ( 542391.600, 4185150.000) (122d31' 6.20"W, 37d48'45.94"N)
Lower Left ( 542391.600, 4180124.400) (122d31' 7.26"W, 37d46' 2.88"N)
Upper Right ( 547417.200, 4185150.000) (122d27'40.66"W, 37d48'45.05"N)
Lower Right ( 547417.200, 4180124.400) (122d27'41.85"W, 37d46' 1.99"N)
Center ( 544904.400, 4182637.200) (122d29'23.99"W, 37d47'23.98"N)
Band 1 Block=128x128 Type=Byte, ColorInterp=Undefined
Band 2 Block=128x128 Type=Byte, ColorInterp=Undefined
Band 3 Block=128x128 Type=Byte, ColorInterp=Undefined
Band 4 Block=128x128 Type=Byte, ColorInterp=Undefined
Еще использование rpc коэфицентов предполагет, что вы предоставляете матрицу рельефа. В вашей команде ( gdalwarp –rpc in.vrt out.tif) я этого не вижу.
Вот подробнее: http://gis-lab.info/qa/orbview3-ortho-gdal.html
Последний раз редактировалось Дмитрий Барышников 08 фев 2015, 03:18, всего редактировалось 1 раз.

juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

Re: GDAL SPOT-6 VRT

Сообщение juffin_h » 07 фев 2015, 22:31

Дмитрий Барышников писал(а):Вообще GDAL поддерживает rpc коэфициенты в формате DigitalGlobe. А SPOT там поддержки не было.
По этому я и спользовал VRT, чтобы GDAL увдел rpc коэфициенты. Иначе геопривязки нет вообще.
Дмитрий Барышников писал(а):Еще исползование rpc коэфицентов предполагет что вы предоставляете матрицу рельефа.
Матрцу я использовал. Но здесь не привел, поскольку при ошбке в 1000м она погоды не делает. Без матрицы gdalwarp по rpc делает георивязку без ортокоррекции.

В метаданных SPOT указано, что RPC в формате NITF RPC00B. Но GDAL, по моему, с этим вариантом и работает.

Аватара пользователя
Дмитрий Барышников
Гуру
Сообщения: 2572
Зарегистрирован: 17 ноя 2009, 19:17
Репутация: 261
Откуда: Москва

Re: GDAL SPOT-6 VRT

Сообщение Дмитрий Барышников » 08 фев 2015, 03:33

1. Я так и не понял у вас в выводе gdalinfo rpc секция есть или нет?
2. Вы уверены, что правильно создали vrt файл (в том же QGIS он ложится туда куда надо и правильно отображается)?

GDAL работает с форматом RPC и RPB. Но формат не причем, т.к. вы его уже изменили при записи в VRT.
По поводу привязки с орто и без - у вас в VRT написано:
<SRS>PROJCS["WGS 84 / UTM zone 37N",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",39],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AUTHORITY["EPSG","32637"]]</SRS>
<GeoTransform> 4.2771600000000000e+005, 2.3999999999999999e+000, 0.0000000000000000e+000, 6.2590391999999993e+006, 0.0000000000000000e+000,-2.3999999999999999e+000</GeoTransform>
Это и есть геопривязка. Тут и СК есть и положение верхней левой точки, и размер пиксела и поворот.
А rpc надо использовать либо с матрицей высот (-to "RPC_DEM..."), либо с какой-то базовой высотой (-to "RPC_HEIGHT..."). Иначе я не уверен, что вообще что-то считается.

juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

Re: GDAL SPOT-6 VRT

Сообщение juffin_h » 08 фев 2015, 15:04

Это не мой vrt файл. Это мне для примера выложили. В SPOT-6 (в том варианте, который у меня ) растры имеют систему координат сенсора. В файле DIM*.xml приведены метаданные. В том числе система координат, ссылки на растр и на файл RPC*.xml. Тип геопривязки указан – RPC. Выходит, что RPC и есть связь между СК сенсора и геодезической СК. Все метаданные в формате dimap. Gdalinfo по файлу tiff (драйвер geotiff) дает только привязку в системе координат сенсора. По файлу DIM*.xml (драйвер dimap) подхватывает еще геодезическую систему координат, но RPC не видит. По этому я их и добавил в vrt.
В понедельник приведу реальные файлы для ясности и то, что выдает gdalinfo.

Есть один нюанс. Файл RPC*.xml содержит два набора коэффициентов – прямые и обратные. Я использовал прямые.

Аватара пользователя
Дмитрий Барышников
Гуру
Сообщения: 2572
Зарегистрирован: 17 ноя 2009, 19:17
Репутация: 261
Откуда: Москва

Re: GDAL SPOT-6 VRT

Сообщение Дмитрий Барышников » 08 фев 2015, 16:38

На сколько я понимаю, RPC это математическая модель камеры. Т.е. вам не дали фокусное расстояние, размер элемента матрицы ПЗЦ, главную точку изображения и др. параметры, а представили это все в виде функции с набором параметров. Это помогает учесть геометрию сенсора при приведении снимка к ортографической проекции.
С собственно привязкой это связано опосредовано.
Я думаю, что фраза "Тип геопривязки указан – RPC" имеет аналогичное значение, что и у других поставщиков съемки - над снимком провели базовую операцию ортоисправления по какой-то средней высоте для региона съемки.
Обычно снимок имеет уже привязку, а при помощи rpc вы ее уточняете.
В принципе, по VRT видно, что привязка указана (я это уже отмечал) - главный вопрос насколько правильно.
Это проверяется через визуальную загрузку vrt растра в QGIS. Ну и все таки надо проконтролировать, что gdalinfo подхватывает RPC.
Если все ок, а орто все равно с такими гигантскими ошиками, тут либо действительно не те коэфициенты или gdal не поддерживает такой вариант (SPOT). Ну т.е. у спота несколько иная математика и надо пользоваться orfeo toolbox.

Аватара пользователя
Дмитрий Барышников
Гуру
Сообщения: 2572
Зарегистрирован: 17 ноя 2009, 19:17
Репутация: 261
Откуда: Москва

Re: GDAL SPOT-6 VRT

Сообщение Дмитрий Барышников » 08 фев 2015, 16:43

Я выполнил gdalinfo на вашем vrt подложив ему произвольный растр:
Секция RPC присутствует, так что проблема с наличием RPC отпадает.
Осталось проверить - правильность начальной привязки и собственное те ли коэфициенты.
$ gdalinfo example.vrt
Driver: VRT/Virtual Raster
Files: example.vrt
/home/bishop/example.tif
Size is 9720, 66383
Coordinate System is:
PROJCS["WGS 84 / UTM zone 37N",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",39],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AUTHORITY["EPSG","32637"]]
Origin = (427716.000000000000000,6259039.199999999254942)
Pixel Size = (2.400000000000000,-2.400000000000000)
Metadata:
AREA_OR_POINT=Area
TIFFTAG_XRESOLUTION=1
TIFFTAG_YRESOLUTION=1
RPC Metadata:
HEIGHT_OFF=210
HEIGHT_SCALE=500
LAT_OFF=55.7571
LAT_SCALE=0.7430
LINE_DEN_COEFF=+1.000000E+00 -9.503236E-06 -6.386284E-04 +2.976162E-04 -1.904468E-07 -3.768348E-06 -3.588716E-06 +2.900618E-06 +2.049445E-07 +4.631482E-06 +4.712116E-07 +0.000000E+00 +1.548356E-08 +0.000000E+00 +5.036180E-08 +4.819471E-08 +4.843385E-08 +1.151909E-07 -4.660568E-06 +1.819432E-08
LINE_NUM_COEFF=-1.217299E-03 +1.969171E-03 -1.000925E+00 -1.019775E-04 -1.229830E-05 +5.464240E-06 +2.950766E-04 -1.953999E-04 +5.792297E-04 -1.083605E-07 -8.263748E-06 +0.000000E+00 -4.426955E-07 +1.433500E-08 -7.759899E-07 -1.435955E-07 -4.594503E-06 +2.071665E-08 -1.441813E-05 -1.001393E-08
LINE_OFF=33201
LINE_SCALE=34420
LONG_OFF=38.0486
LONG_SCALE=0.1888
SAMP_DEN_COEFF=+1.000000E+00 -2.492860E-03 +8.702597E-03 -4.031501E-04 +2.042522E-05 -2.319546E-05 +9.809980E-06 -5.156165E-05 +7.904697E-05 -9.248327E-05 -2.709752E-07 -1.701596E-08 +2.933039E-07 +3.913028E-08 -9.259787E-07 +1.708361E-06 -1.859015E-06 -1.031072E-08 +4.222594E-07 +1.142850E-07
SAMP_NUM_COEFF=+5.959196E-03 +9.446828E-01 +9.045097E-02 -2.201229E-02 -9.965023E-03 +8.334638E-04 -2.793291E-03 -2.345943E-03 +1.193017E-03 -1.723565E-05 +3.216177E-05 -4.933984E-05 -1.610468E-04 -8.584692E-05 +5.990033E-05 +8.395135E-06 -9.615944E-06 -2.591164E-05 +7.440910E-06 +1.979198E-06
SAMP_OFF=5204
SAMP_SCALE=5226
Corner Coordinates:
Upper Left ( 427716.000, 6259039.200) ( 37d49'36.07"E, 56d28'13.00"N)
Lower Left ( 427716.000, 6099720.000) ( 37d52' 7.98"E, 55d 2'20.53"N)
Upper Right ( 451044.000, 6259039.200) ( 38d12'19.06"E, 56d28'23.80"N)
Lower Right ( 451044.000, 6099720.000) ( 38d14' 1.97"E, 55d 2'30.77"N)
Center ( 439380.000, 6179379.600) ( 38d 2' 2.68"E, 55d45'22.69"N)
Если и там все ОК, то gdal действительно с метаданными спота не работает.

alexandr cherepanov
Гуру
Сообщения: 534
Зарегистрирован: 30 ноя 2006, 13:31
Репутация: 116
Откуда: Moscow

Re: GDAL SPOT-6 VRT

Сообщение alexandr cherepanov » 08 фев 2015, 23:10

juffin_h писал(а): один нюанс. Файл RPC*.xml содержит два набора коэффициентов – прямые и обратные. Я использовал прямые.
И не угадали :) , нужны обратные насколько помню Inverse_Model.

juffin_h
Завсегдатай
Сообщения: 265
Зарегистрирован: 22 окт 2012, 08:35
Репутация: 49
Откуда: Нижний Новгород

Re: GDAL SPOT-6 VRT

Сообщение juffin_h » 09 фев 2015, 10:37

alexandr cherepanov
Да, в этом и есть моя принципиальная ошибка :) .

Сделал два файла aaa.vrt и bbb.vrt. C использованием Dirrect_Model и Inverse_Model коэффициентов соответственно. И попробовал пересчитать.
Предположительные координаты первой точки растра:
56.9260335419 51.1710700959

Результат работы gdaltransform:

Код: Выделить всё

>gdaltransform -rpc -to "RPC_HEIGHT=280.0" aaa.vrt 
input: 0.5 0.5
output: 56.9112644423733 51.1727931433332 0

Код: Выделить всё

>gdaltransform -rpc -to "RPC_HEIGHT=280.0" bbb.vrt 
input: 0.5 0.5
output: 56.9260330017306 51.1710635256391 0
Второй вариант явно лучше.

alexandr cherepanov
Гуру
Сообщения: 534
Зарегистрирован: 30 ноя 2006, 13:31
Репутация: 116
Откуда: Moscow

Re: GDAL SPOT-6 VRT

Сообщение alexandr cherepanov » 09 фев 2015, 10:54

Вот vrt по реальным RPC SPOT 6 только для PAN.

Код: Выделить всё

<VRTDataset rasterXSize="38609" rasterYSize="93093">
  <Metadata>
    <MDI key="TIFFTAG_XRESOLUTION">100</MDI>
    <MDI key="TIFFTAG_YRESOLUTION">100</MDI>
  </Metadata>
  <Metadata domain="RPC">
    <MDI key="HEIGHT_OFF">+500.0 meters</MDI>
    <MDI key="HEIGHT_SCALE">+500.0 meters</MDI>
    <MDI key="LAT_OFF">+38.14680766 degrees</MDI>
    <MDI key="LAT_SCALE">+0.65112231 degrees</MDI>
    <MDI key="LINE_DEN_COEFF">+1.0  +1.171970711705561E-8  +9.973772532013897E-9  +3.511684622474163E-10  +8.832303258929766E-8  +9.681464134541095E-11  +2.871180539102313E-9  -5.527444456947567E-9  +4.180193541283476E-8  -6.846687887830579E-12  -4.729135399547531E-12  -4.290152960551305E-10  +1.406997044694272E-8  +1.856199653627005E-13  -4.887411923344226E-10  -3.112499066440895E-9  -1.501083533861343E-13  +2.79327601197463E-11  +4.034054066320384E-10  +5.672286206760068E-15 </MDI>
    <MDI key="LINE_NUM_COEFF">+0.001276562702162548  +0.0001131283318517914  -1.000920103541099  -0.0003143472793749695  -0.0004704727976734494  -1.458719033516534E-7  -1.545949497336026E-5  -0.0007412028713334032  -0.000667087584202698  +1.993062791507356E-8  -2.968484313635449E-8  -3.476134708390151E-7  -2.281438390676952E-5  +4.890200607842853E-12  -4.135421386616048E-6  -2.180301820290683E-5  +1.431073958664704E-9  -1.416898803592059E-8  +6.826489454786002E-7  -1.193521682264104E-12 </MDI>
    <MDI key="LINE_OFF">+46546.5  pixels</MDI>
    <MDI key="LINE_SCALE">+46546.5  pixels</MDI>
    <MDI key="LONG_OFF">+71.95837831 degrees</MDI>
    <MDI key="LONG_SCALE">+0.34673781 degrees</MDI>
    <MDI key="SAMP_DEN_COEFF">+1.0  -4.477649065904744E-6  +1.265615480628109E-5  -4.889652298967428E-6  -1.547368875928256E-5  +4.121062570677158E-6  -2.119930069789608E-6  -2.345917347265104E-6  +9.812500637069498E-7  -8.273776625442081E-9  -6.28696796464711E-8  -8.147563669803938E-9  +2.055374127695115E-7  -3.758908482261207E-11  +4.449973618701111E-8  -7.537286287902737E-8  +2.864046099559989E-10  +2.702921830472764E-9  +5.289397551117258E-8  -1.040366854896172E-10 </MDI>
    <MDI key="SAMP_NUM_COEFF">-0.001882775812670119  +1.008367300680687  +0.003112176138881072  +0.0008088622253749416  -0.006682805177499917  +0.0007886064477391866  -0.0004465435510247784  -0.002248150685236973  -1.252716748716211E-5  +4.291810612294975E-7  -4.399974294921424E-6  -0.0001049794742026144  -0.0007831626929331334  +5.553394179977077E-7  +0.001265882577730264  +0.0001081137455484875  -3.241441475707899E-7  +7.138903940402302E-7  -7.79329484317618E-7  +3.080863041319949E-10 </MDI>
    <MDI key="SAMP_OFF">+19304.5  pixels</MDI>
    <MDI key="SAMP_SCALE">+19304.5  pixels</MDI>
  </Metadata>
  <VRTRasterBand dataType="Byte" band="1">
    <Metadata />
    <ColorInterp>Gray</ColorInterp>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">IMG_SPOT6_P.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="38609" RasterYSize="93093" DataType="Byte" BlockXSize="38609" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="38609" ySize="93093" />
      <DstRect xOff="0" yOff="0" xSize="38609" ySize="93093" />
    </SimpleSource>
  </VRTRasterBand>
</VRTDataset>

Соответственно, для орто нужны в команде gdalwarp ключи
-tr
-t_srs
-rpc -to "RPC_DEM=FILENAME" (высоты над эллипсоидом)
Ну и что то можно взять из примера приведенного выше с OrbView

gdalwarp -tr 2.0 2.0 -dstnodata 0 -srcnodata 0 -overwrite -r bilinear -t_srs epsg:32642 -wo "INIT_DEST=NO_DATA" -et 0.0 -rpc -to "RPC_DEM=c:\SPOT\dem.tif" c:\SPOT\IMG_SPOT6_P_vrt.vrt c:\SPOT\ortho_IMG_SPOT6_P_vrt.tif

Ответить

Вернуться в «GDAL/OGR»

Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей и 1 гость