Просмотр кода страницы Основы макроязыка EASI (PCI Geomatica)Часть 3. Географическая привязка и проекции

{{Статья|Опубликована|easi3}}
{{Аннотация|Третья часть серии статей о языке EASI, основные команды для географической привязки и работы проекциями, системами координат.}}

=== Поддерживаемые проекции, эллипсоиды и датумы ===

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

Проекции, поддерживаемые на эллипсоиде:

* ACEA — Albers Conical Equal-Area
* CASS — Cassini
* EC — Equidistant Conic
* KROV — Krovak
* LCC — Lambert Conformal Conic
* MER — Mercator
* MSC — Mod. Stereographic Conformal (Alaska)
* OM — Oblique Mercator
* PC — Polyconic
* PS — Polar Stereographic
* RSO — Rectified Skew Orthomorphic
* SOM — Space Oblique Mercator
* SPAF — State Plane (U.S. Feet)
* SPCS — State Plane (Metres)
* SPIF — State Plane (International Feet)
* TM — Transverse Mercator
* UPS — Universal Polar Stereographic
* UTM — Universal Transverse Mercato

Кроме того, географические координаты — LONG/LAT.

Проекции, поддерживаемые ТОЛЬКО на сфере:

* AE — Azimuthal Equidistant
* ER — Equirectangular
* GNO — Gnomonic
* GOOD — Goode's Homolosine
* GVNP — Gen. Vertical Near-Side Perspective
* LAEA — Lambert Azimuthal Equal-Area
* MC — Miller Cylindrical
* OG — Orthographic
* ROB — Robinson
* SG — Stereographic
* SIN — Sinusoidal
* VDG — Van der Grinten

Если в параметрах проекции указан другой эллипсоид, этот параметр игнорируется и выполняется пересчет со сферы.

Эллипсоиды обозначаются кодами, начинающимися с буквы «Е». Буква может быть опущена (тогда предполагается, что используется эллипсоид). Список эллипсоидов довольно большой. Вы можете сами найти код нужного Вам эллипсоида в справочной системе EASI. Приведем лишь наиболее часто используемые:

* Е000 (или просто 0) — Сlarke 1866 (этот эллипсоид будет подразумеваться в том случае, если никакого эллипсоида и датума не указано!)
* Е008 — GRS 1980
* Е012 — WGS 84
* Е015 — Krassovsky
* Е019 — Sphere of Radius 6370997 metres

Эллипсоид может быть задан полуосями («6378206.4 6356583.8»).

Датум задается числовым кодом, начинающимся с буквы «D». Информация о проекции содержит только что-то одно: или имя эллипсоида (тогда подразумевается геоцентрический датум) или имя датума (тогда имя эллипсоида однозначно определено кодом датума). Некоторые распространенные датумы:

* D800 – геоцентрический датум для сферы (=Е019)
* D000 – WGS 84 (=Е012)
* D203 – Pulkovo 1942 (для России)
* D-01 — NAD27 (USA, NADCON)
* D-02 — NAD83 (USA, NADCON)

Пользовательские датумы могут быть добавлены в таблицу датумов, хранящуюся в текстовом формате (<tt>$PCIHOME/etc/datum.txt</tt>). Обратите внимание, что бы номер нового датума не совпадал с уже имеющимися в таблице!

=== Вывод информации о проекции ===

Информация о привязке и проекции файла хранится в специальном сегменте данных, имеющем номер 1. Вывести на экран или в файл информацию о проекции можно с помощью команды <tt>GEOREP</tt>:

==== GEOREP – вывод информации о привязке ====

* <tt>file="</tt> (имя файла)
* <tt>report=</tt> (“term или имя файла для вывода)

<pre>
Georeference Units:                UTM    33   D000
Projection:                        Universal Transverse Mercator
Datum - Ellipsoid:                 WGS 1984 (Global Definition) - WGS 84
Grid units:                        METRE
	
Upper Left Corner          288804.750  E               108912.750  N
Upper Right Corner         533106.750  E               108912.750  N
Image Centre               410955.750  E               -99.750  N
Lower Left  Corner         288804.750  E               -109112.250  N
Lower Right Corner         533106.750  E               -109112.250  N

Pixel Size                 57.000  E                   57.000 N

Upper Left Corner          13d06'07.60" E  Lon         0d59'05.36" N  Lat
Upper Right Corner         15d17'51.23" E  Lon         0d59'07.27" N  Lat
Image Centre               14d11'59.31" E  Lon         0d00'03.25" S  Lat
Lower Left  Corner         13d06'07.60" E  Lon         0d59'11.86" S  Lat
Lower Right Corner         15d17'51.23" E  Lon         0d59'13.77" S  Lat
</pre>

Формат таблицы стандартный, что позволяет ее использовать при автоматизации процессов с помощью языка PERL (см. [http://gis-lab.info/qa/easi4.html статью 4]). Следует обратить внимание на замеченные проблемы с выводом таблицы. В ряде случаев (по неизвестной причине) координаты углов могут быть перепутаны (обычно правый и левый углы меняются местами).

=== Назначение параметров привязки и проекции растровой базе данных ===

==== Импорт из другого формата ====
Существуют несколько способов задания проекции для файла PCIDSK. Прежде всего, проекция может быть автоматически импортирована при импорте данных из другого формата командой FIMPORT. Следует обратить внимание на две замеченные проблемы, возникающие при импорте информации о проекции и привязке:

# При импорте из формата GeoTIF PCI подразумевает, что координаты верхнего левого угла являются координатами границ углового пикселя, в то время как ряд программ подразумевает, что это координаты центра углового пикселя. Ошибка легко может быть замечена при импорте файла GeoTIF одновременно в ERDAS Imagine и PCI, с последующим экспортом обратно в GeoTIF. При этом файл, экспортированный из PCI, будет смещен на пол-пикселя по обеим осям относительно оригинала и файла, экспортированного из ERDAS. Возможные пути решения проблемы: смещение импортированного в PCI файла или импорт в PCI через другие форматы (например, *.img).
# В ряде случаев параметры проекции могут быть прочитаны неправильно. Особенно часто проблема проявляется при импорте данных из южного полушария в проекции UTM. Для устранения ошибки, проекция может быть переписана в EASI или исправлена через другие приложения PCI (Focus, ImageWorks).

==== Создание нового файла с указанием проекции и привязки ====

Создать новый файл и назначить ему параметры привязки и проекцию можно с помощью команды <tt>CIMPRO2</tt>. Использование команды рассмотрено в [http://gis-lab.info/qa/easi1.html статье 1]. Здесь мы приведем пример использования команды для создания нового файла в проекции UTM:

* <tt>file="new</tt>(имя файла)
* <tt>dbnc=10,0,0,0</tt> (число каналов - 10 8-битных каналов)
* <tt>dblayout=”file</tt> (тип хранения каналов – по разным файлам)
* <tt>project="UTM</tt> (имя проекции)
* <tt>zone=33</tt> (зона проекции UTM)
* <tt>ellips="D000</tt> (датум WGS1984)
* <tt>llbound="n</tt> (координаты углов: если «y» – географические, если «n» - метрические)
* <tt>ulx="288804.750</tt>
* <tt>uly="108912.750</tt>
* <tt>lrx="533106.750</tt>
* <tt>lry="-109112.250</tt>
* <tt>bxpxsz="57.000 (размер пикселя)</tt>
* <tt>bypxsz="57.000</tt>
* <tt>extrainf=</tt> (дополнительная информация о проекции по умолчанию)

<pre>r cimpro2</pre>

==== Назначение проекции созданному файлу ====

Существует несколько команд задания проекции файлу PCIDSK. Мы рекомендуем по возможности всегда использовать команду <tt>SETRPO2</tt>. После ввода параметров и запуска команды, значения первого сегмента (сегмента привязки) файла PCIDSK будут изменены.

===== SETPRO2 – Ввод информации о привязке и проекции растрового файла =====

* <tt>file="</tt> (имя файла)
* <tt>project="</tt> (имя проекции)
* <tt>zone=</tt> (зона проекции UTM)
* <tt>ellips="</tt> (эллипсоид или датум)
* <tt>llbound="</tt> (координаты углов: если y – географические, если n - метрические)
* <tt>ulx="</tt> (координаты углов)
* <tt>uly="</tt>
* <tt>lrx="</tt>
* <tt>lry="</tt>
* <tt>extrainf=</tt> (дополнительная информация по проекции. Команда включает множество дополнительных параметров для разных проекций. Все эти параметры могут быть объявлены «по умолчанию»)

<pre>r setpro2</pre>

Существуют другие команды для ввода информации о проекции и привязке:

* <tt>SETPRO</tt> – команда аналогичная SETPO2, но с меньшими возможностями при вводе параметров проекции
* <tt>GEOSET</tt> – Простая команда, позволяющая назначить ограниченный набор проекций (UTM, LONG/LAT, State plane). Замечена следующая ошибка: при вводе координат углов значения по оси Y обычно округляются, что приводит к неправильному назначению информации о привязке.

==== Создание дополнительного сегмента привязки с информацией о проекции ====

Для использования ряда программ (например, <tt>DBREPRO</tt>) бывает необходимым создание дополнительного сегмента с указанием проекции, в которую Вы хотите спроецировать данные. Для этого используем команду <tt>PROSET</tt>.

===== PROSET – изменение значений или создание нового сегмента привязки с указанием проекции =====

* <tt>file="</tt> (имя файла)
* <tt>dbgeo=</tt> (номер сегмента. если номер=1, то будет изменена информация о проекции файла. Если номер не указан, будет создан новый сегмент)
* <tt>dbsn="</tt> (имя нового сегмента – обязательный параметр)
* <tt>mapunits="</tt> (Имя проекции. Если проекция зональная, так же и номер зоны. Например, “UTM 33”)
* <tt>ellips="</tt> (эллипсоид или датум)
* <tt>prolat="</tt> (стандартная параллель)
* <tt>prolong="</tt> (центральный меридиан)

<pre>r prosin</pre>

Команда включает множество дополнительных параметров для разных проекций. Для упрощения ввода только параметров нужной проекции, используются специальные процедуры, основанные на команде <tt>PROSET</tt>. Имя процедур для каждой проекции Вы можете найти в справочной системе EASI. Например, для определения синусоидальной проекции используется процедура <tt>PROSIN</tt>:

* <tt>file="</tt> (имя файла)
* <tt>dbgeo=</tt>
* <tt>dbsn="</tt> (имя нового сегмента – обязательный параметр)
* <tt>ellips="19</tt> (для синусоидальной проекции всегда сфера!)
* <tt>prolong="0.0</tt> (стандартные параметры)

<pre>r prosin</pre>

=== Изменение проекции ===

Изменение проекции (перепроецирование) растровых данных в EASI решено в виде серии команд:

# Создание нового сегмента с информацией о требуемой проекции (<tt>PROSET</tt>)
# Пересчет координат растрового файла и определение размеров выходного файла в требуемой проекции (<tt>DBREPRO</tt>)
# Создание файла с рассчитанными параметрами (<tt>CIM</tt>, <tt>CIMPRO2</tt>)
# Трансформация каналов растра из одной проекции в другую (<tt>REGPRO</tt>, <tt>MOSAIC</tt>).

==== DBREPRO – расчет параметров перепроецирования ====

Данная команда имеет дополнительные параметры, позволяющие вводить ограниченное число проекций для трансформации. Для трансформации в пользовательскую проекцию обычно бывает удобным сначала создать новый сегмент с параметрами проекции с помощью команды <tt>PROSET</tt>, а затем использовать в параметрах только имя нового сегмента:

* <tt>file="</tt> (имя файла)
* <tt>project=</tt> (оставляем пустым – этот параметр будет замещен информацией из нового сегмента)
* <tt>ellips=</tt> (оставляем пустым – этот параметр будет замещен информацией из нового сегмента)
* <tt>outgeo=</tt> (номер нового сегмента с параметрами проекции)
* <tt>report="</tt> (вывод параметров трансформации на экран или в файл)

<pre>r dbrepro</pre>

Дополнительно можно указать размер пикселя выходной базы данных (параметры <tt>bxpxsz</tt> и <tt>bypxsz</tt>).

После запуска команды <tt>DBREPRO</tt> в параметрах EASI сохраняются границы и параметры проекции, в которую Вы хотите перепроецировать Ваши данные. Вы можете воспользоваться командой <tt>CIMPRO2</tt> для создания новой базы данных, в которую будет произведена трансформация. Рассчитанные с помощью <tt>DBREPRO</tt> параметры нового файла будут автоматически использованы командой <tt>CIMPRO2</tt>. Обратите внимание, чтобы число каналов создаваемого файла соответствовало числу каналов, которые Вы хотите перепроецировать! Другой возможный вариант – вывод параметров в файл (опция report) и обработка файла программой на языке PERL. Этот путь может быть необходим, если Вы, допустим, хотите, что бы новый файл был «привязан» к стандартной растровой сетке.

После создания нового файла Вы можете произвести трансформацию исходного растра в полученный пустой файл с помощью команды <tt>REGPRO</tt>.

==== REGPRO – трансформация растрового изображения ====

* <tt>fili="</tt> (имя исходного файла)
* <tt>dbic=</tt> (список каналов для трансформации)
* <tt>ingeo=</tt> (номер сегмента привязки исходного файла – по умолчанию 1)
* <tt>filo="</tt> (имя выходного файла)
* <tt>dboc=</tt> (список каналов выходного файла)
* <tt>outgeo=</tt> (номер сегмента привязки выходного файла – по умолчанию 1)
* <tt>resample="</tt> (тип пересчета near/bilin/cubic)

<pre>r regpro</pre>

Дополнительные параметры программы смотрите в справочной системе EASI.

Вместо команды <tt>REGPRO</tt> для трансформации изображения может использоваться команды <tt>MOSAIC</tt>. Эта команда всегда будет использовать первые сегменты файлов для расчета параметров трансформации, и пересчет всегда будет производиться методом ближайшего соседа (см. [http://gis-lab.info/qa/easi1.html статью 1]).

=== Пример использования команд для перепроецирования снимка ===

Исходный снимок Landsat имеет проекцию UTM. Мы перепроецируем снимок в проекцию Albers Conical Equal Area с параметрами, настроенными для Европейской России.

# Создаем новый сегмент привязки с параметрами требуемой проекции:
#* <tt>file="img.pix</tt> (имя исходного снимка)
#* <tt>dbgeo= dbsn="albers</tt> (имя проекции)
#* <tt>ellips="D203</tt> (датум – Пулково 1942)
#* <tt>prolong="45</tt> (центральный меридиан)
#* <tt>prolat="0</tt> (широта начала отсчета координат)
#* <tt>stpar1="52</tt> (1 стандартная параллель)
#* <tt>stpar2="64</tt> (2 стандартная параллель)
#* <tt>efalse="8500000</tt> (смещение начала координат на запад)
#* <tt>nfalse=</tt><pre>r proacea</pre>
# Вычисляем параметры трансформации:
#* <tt>file="img</tt>
#* <tt>project=</tt>
#* <tt>ellips=</tt>
#* <tt>outgeo=2</tt> (новый сегмент с информацией о проекции)
#* <tt>report="</tt>
#* <tt>bxpxsz="28.5</tt> (требуемые размеры пикселя)
#* <tt>bypxsz="28.5</tt><pre>r dbrepro</pre>
# Создаем новый файл командой <tt>CIMPRO2</tt>. Все необходимые параметры команды будут прописаны автоматически:
#* <tt>file="img_alb</tt> (имя нового файла)
#* <tt>dbnc=3,0,0,0</tt> (число каналов нового файла)
#* <tt>dblayout="file</tt> (способ хранения каналов)<pre>r cimpro2</pre>
# Выполняем трансформацию:
#* <tt>fili="img</tt> (имя исходного файла)
#* <tt>dbic=1,2,3</tt> (список каналов для трансформации)
#* <tt>ingeo= filo="img_alb</tt> (имя выходного файла)
#* <tt>dboc=1,2,3</tt> (список каналов выходного файла)
#* <tt>outgeo=</tt>
#* <tt>resample="near</tt> (тип пересчета)<pre>r regpro</pre>

=== Калькулятор координат ===

<tt>PROCONV</tt> — простая команда EASI, позволяющая пересчитать координаты из одной проекции (или из географических координат) в другую проекцию (или географические координаты).

Ввод координат для пересчета обычно осуществляется из тестового файла. Пары координат в текстовом файле разделяются запятыми:

<pre>288804.750,108912.750</pre>

Если вводятся географические координаты, то они могут иметь вид десятичных градусов или градусы-минуты-секунды:

<pre>
-20.4560,35.4568 (десятичные градусы)
13 06 07.60 E ,0 59 05.36 N (градусы-минуты-секунды)
</pre>

Исходная и выходная проекции могут быть заданы по имеющимся сегментам привязки или как параметр команды (только для LONG/LAT, UTM и State Plane). Пример использования программы с заданием проекции из сегментов привязки:

* <tt>file="</tt> (имя файла, содержащего сегменты привязки)
* <tt>dbgeo=1,2</tt> (номера сегментов, описывающих входную и выходную проекции)
* <tt>tfile="</tt> (имя текстового файла с парами координат)
* <tt>ltyp="brief</tt> (вывод в короткой форме)
* <tt>dmsform="off</tt> (если выходные координаты географические, то они будут представлены в десятичных градусах)
* <tt>report="</tt> (имя файла для сохранения результата)

<pre>r proconv</pre>

Пример использования команды для пересчета координат с заданием проекции через параметры файла: <tt>file=</tt> (пустой)

* <tt>dbgeo=</tt> (пустой)
* <tt>iunits=lon e012</tt> (исходные координаты – географические)
* <tt>ounits="utm zone e012</tt> (выходные координаты – UTM)
* <tt>tfile="</tt> (имя текстового файла с парами координат)
* <tt>ltyp="brief</tt> (вывод в короткой форме)
* <tt>dmsform="off</tt> (если выходные координаты географические, то они будут представлены в десятичных градусах)
* <tt>report="</tt> (имя файла для сохранения результата)

<pre>r proconv</pre>

=== Краткий перечень команд привязки изображений по набору контрольных точек (GCP) ===

Расстановка контрольных точек (GCP) для привязки изображения может проводиться как в PCI, так и в других программа. Контрольные точки хранятся в специальном сегменте в файле PCIDSK. Существует возможность загрузки-выгрузки таблицы контрольных точек в текстовый файл. В PCI существуют модули автоматической расстановки контрольных точек (основанные на расчете корреляции между изображениями). Правда, они продаются отдельно от основного пакета программы. Помимо алгоритма с использованием целого исходного снимка для привязки, реализован подход для создания набора модельных участков (image chips) для автоматической привязки. Ниже приведен краткий перечень команд для работы с наборами контрольных точек и автоматической привязки. Трансформация снимка на основе набора контрольных точек производится командой REGPRO (см. выше) или REG. Описание команд смотрите в справочной системе EASI.

Работа с сегментом контрольных точек (GCP):

* <tt>GCPREAD</tt> – Прочитать набор контрольных точек из текстового файла
* <tt>GCPWRIT</tt> – Записать набор контрольных точек в текстовый файл
* <tt>GCPREP</tt> – Вывод описания сегмента контрольных точек
* <tt>GCPPRO</tt> – Проецирование сегмента контрольных точек в другую проекцию
* <tt>GCPREFN</tt> – Автоматическое удаление контрольных точек с большой ошибкой

Автоматическая расстановка контрольных точек и привязка снимка:

* <tt>AUTOGCP</tt> – Создание набора контрольных точек путем корреляционного анализа
* <tt>AUTOREG</tt> – Автоматическая ко-регистрация изображений
* <tt>CHIPEXT</tt> – Создание набора модельных участков на основе привязанного снимка
* <tt>CHIPMAT</tt> – Автоматический поиск контрольных точек на основе библиотеки модельных участков

=== Ссылки по теме ===

* [http://gis-lab.info/qa/easi1.html Основы макроязыка EASI (PCI Geomatica). Часть 1. Управление данными]
* [http://gis-lab.info/qa/easi2.html Основы макроязыка EASI (PCI Geomatica). Часть 2. Анализ данных]
* [http://gis-lab.info/qa/easi4.html Основы макроязыка EASI (PCI Geomatica). Часть 4. Примеры автоматизации процессов]