Создание карты высот на основании растра с цветной палитрой

[Zeloni]

Добрый день. Помогите, пожалуйста, решить задачу в QGis или GlobalMapper.

У меня есть цветной растр RGB обозначающий высоты местности (пример):



Растр содержит 113 цветов, каждый из которых соответствует определенной высоте местности.
В зависимости от проекта диапазон высот может меняться. К примеру я делаю проект местности с высотой от 0 до 150 метров. (синий цвет - 150 метров, красный 0 метров).

Я знаю значения RGB каждого цвета и соответствующая ему высоту например:
MaxElevation="15.93" red="0" blue="222" green="125" MinElevation="14.6" или
0 222 125 15,93 14,6 (без разницы какой формат, смогу сделать любой)

Задача получить DEM модель или GeoTiff с картой высоты чтобы ее можно было использовать в любом другом софте.

Для себя пока выработал следующий алгоритм:
1. Делаю из RGB с помощью специального алгоритма grayscale.
2. Привязываю карту к местности с помощью Georeferencer.
3. Создаю контуры произвольных высот с шагом например 10 с помощью Raster - Extraction - Contour.
4. В Attribute Table полученного слоя преобразовываю значения высот в поле ELEV в нужный диапазон (например я получил значения поля ELEV от 30 до 1150, преобразовал их в диапазон значений от 0 до 150).
5. Далее с помощью Raster - Analysis - Grid (Interpolation) делаю нужный мне файл.

У меня вопрос, а нет ли способа попроще.
Спасибо.

[Максим Дубинин]

http://gis-lab.info/qa/gdal-hillshade-colormap.html

[Zeloni]

Максим, спасибо за ответ.

В указанном Вами примере (я его ранее читал):
В качестве единственного источника выступает ЦМР, растровая матрица, где каждому пикселю соответствует значение его высоты над уровнем моря. Это может быть SRTM, ASTER GDEM и другие источники.

Я умею доставать эти данные (strm, aster и т.д. из других источников), а мне нужно решить обратную задачу - получить растровую матрицу, где каждому пикселю соответствует значение его высоты над уровнем моря из моего источника, в котором нет информации какой высоте над уровнем моря соответствует каждый пиксель. Мне нужно как то преобразовать мой источник, чтобы эта информация появилась.

Я не могу понять каков формат этой матрицы. Какой цвет соответствует определенной высоте над уровнем моря. и каким образом можно ее сгенерировать используя исходные данные что у меня есть.

Мучаюсь уже не первый день. Допускаю что все может быть достаточно просто и не хватает какихто базовых знаний.

[gamm]

Растр содержит 113 цветов, каждый из которых соответствует определенной высоте местности.

не очень понятно, в чем проблема: загрузите растр и таблицу перекодировки в Python/R, напишите 113 if(...) ... для перекодировки (или используйте словарь в Питоне/match() в R, или сгенерируте программу с этими if() из таблицы перекодировки), и перекодируйте цвета в высоты. Потом сохраните максимальные/средние/минимальные высоты в растр. Все это будет работать, если у вас действительно растр с цветами, а не пережатый jpeg, в последнем случае ничего умного сделать не получится, разве что искать приблизительно похожие цвета.

[Zeloni]

Спасибо. Посоветуйте, пожалуйста, документацию по Python и R буду изучать.
Просто для меня "перекодируйте цвета в высоты" - это темный лес

Если не сложно можно для примера код одного
if(...) ..

Спасибо.

[lam]

Инструмент для Арки: ArcToolbox-Spatial Analyst Tools-Reclass-Reclassify. Растр должен иметь цветовую модель "8-bit Palette"

[gamm]

выложите небольшой фрагмент растра, попробую написать программку

[Zeloni]

Добрый день.СПАСИБО!!
Прилагаю файлы файлы.
и таблицу соответствия высот цветам.

[gamm]

сам рельеф не фонтан (явно TIN), на местах, где были надписи - дырки, некоторых цветов в таблице не хватает, эти данные выводятся в файл colors_not_found.txt. Нужно поставвить R (скачать с https://cran.r-project.org/), потом пакеты:
install.packages("rgdal")
install.packages("raster")

запустить R, установить рабочую папку (где файлы - в пути не должно быть русских букв, пробелов, и пр.), потом поправить имена и скопипастить программу в консоль R

И если понравится, поддержите gis-lab.

prg.rar

program

(325.68 КБ) 422 скачивания

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

library(rgdal)

# --- files
elev.fn<-"Elev r g b.txt"
log.fn<-"colors_not_found.txt"
in.fn <-"SNAG-00149.png"
out.fn<-"SNAG-00149-Elev.tif"

# --- load recode table and create color keys
p.recode<-read.table(elev.fn,header=FALSE)
names(p.recode)<-c("Elev","r","g","b")
p.key<-sprintf("%03d_%03d_%03d",p.recode$r,p.recode$g,p.recode$b)

# --- load image
work.region.tif <- new("GDALReadOnlyDataset", in.fn)
p.b<-getRasterTable(work.region.tif)        # by rows, from the upper-left corner
p.bbox<-work.region.tif[,,1]@bbox           # bounding box
p.csize<-work.region.tif[,,1]@grid@cellsize # cell (pixel) size
p.b.prj<-getProjectionRef(work.region.tif)  
p.b.dim<-dim(work.region.tif)               # c(ny,nx,nband)

# --- get grid
tmp<-work.region.tif[,,1:4]

# --- close image
GDAL.close(work.region.tif)

# --- create color keys
img.key<-sprintf("%03d_%03d_%03d",tmp$band1,tmp$band2,tmp$band3)

# --- get positions in recode table
pos<-match(img.key,p.key)
if(any(ind<-is.na(pos))) {
  write.table(tmp[ind,],log.fn,row.names=FALSE,sep="\t")
}

# --- recode positions to elevation
tmp@data$Elev<-NA
tmp@data$Elev[!ind]<-p.recode$Elev[pos[!ind]]

# --- save elevation as TIFF
p.out<-data.frame(x=p.b[,1],y=p.b[,2],Elev=tmp@data[,"Elev"])
coordinates(p.out) = ~x+y
gridded(p.out) = TRUE
p.out.grd=as(p.out,"SpatialGridDataFrame")
p.out.grd@proj4string@projargs<-as.character(p.b.prj)
p.out.tif<-create2GDAL(p.out.grd,drivername="GTiff",type="Float32")
saveDataset(p.out.tif,out.fn)
GDAL.close(p.out.tif)

[Zeloni]

Спасибо огромное. Все работает.

И если понравится, поддержите gis-lab.

С удовольствием