Crear un shapefile para huella de trama ' áreas de datos válidos de s con GDAL

Question

Me gustaría mostrar las fronteras de las áreas de un determinado ASCII DTM archivo para el que hay datos. Este problema ha sido discutido aquí anteriormente, pero la solución recomendada, los bordes de la Imagen con el plugin de QGIS, ahora se considera obsoleto. Además, tengo que ser capaz de scripts, por lo que el uso de la GDAL herramientas desde la línea de comandos es probablemente un plan mejor para mí.

Yo también tengo una restricción adicional: tengo que ser capaz de servir a este archivo en GeoServer, y en última instancia, la muestra mediante Prospecto, por lo que un shapefile o GeoJSON podrían ser buenas opciones para el formato de archivo de la salida final.

Qué comando debo usar?

Answer 1

Tengo un PASTO/R flujo de trabajo a continuación que implica GDAL. Es mucho más complicado y desordenado de lo que yo recuerdo, pero funciona para lo que sea. Blog aquí con los comentarios adicionales, llanura instrucciones para un entorno Windows con OSGeo4W64, R 3.3.1 y la HIERBA 7.0.5 a continuación.

1: Obtener todos sus rásteres en una carpeta. Asegúrese de que:

• los archivos se nombran en una secuencia lógica que no empezar con un número.
• los archivos están todos en la misma proyección
• los archivos deben ser espacial alineados; es decir, los píxeles deben estar todos en el mismo general el sistema de la red.

2: Vaya a Sistema > Configuración Avanzada del Sistema > Variables de Entorno y establezca las siguientes (ajustar según sea necesario, si usted instalado en otro lugar, o tener una copia anterior de la HIERBA):

• GDAL_DATA: C:\OSGeo4W64\share\gdal
• GISBASE: C:\OSGeo4W64\apps\grass\grass-7.0.5
• GISRC: C:\Users\INICIO de sesión AQUÍ\AppData\Roaming\GRASS7\rc
• GRASS_PYTHON: C:\OSGeo4W64\bin\pythonw.exe
• PROJ_LIB: C:\OSGeo4W64\share\proj
• PYTHONPATH: C:\OSGeo4W64\apps\grass\grass-7.0.5\etc\python

• RUTA de acceso: agregar las siguientes rutas de acceso, en el mismo orden a continuación

  • C:\OSGeo4W64\bin
  • C:\OSGeo4W64\apps\grass\grass-7.0.5\bin
  • C:\OSGeo4W64\apps\grass\grass-7.0.5\lib
  • C:\OSGeo4W64\apps\qgis
  • C:\OSGeo4W64\apps\Python27

• Quitar GRASS_ADDON_BASE.

3: Configurar la ubicación de la HIERBA

Abrir el CÉSPED y establecer una nueva ubicación con el mismo CRS, ya que su objetivo rásteres. No se preocupe acerca de la configuración de la región, por ahora, y no hay necesidad de crear cualquier extra mapsets más allá de la 'PERMANENTE'. Abra la PERMANENTE mapset en su nueva ubicación y, a continuación, salga de la HIERBA. Esto actualizará el archivo que GISRC puntos.

R código:

La carga de las bibliotecas y de la lista de los archivos de destino

library(sp)  
library(raster)  
library(rgdal)  
library(XML)  
library(rgrass7)  

setwd('C:/footprint_test')  
rasters <- list.files(getwd(), pattern = '\\.tif$', full.names = T)  

Importación de datos

grassnames_list <- list()  
for (item in rasters) {  
  grass_name <- substr(basename(item), 1, nchar(basename(item)) - 4)  
  execGRASS('r.in.gdal',   
            flags = 'overwrite',  
            parameters = list(input  = item,  
                              output = grass_name))  
  grassnames_list[[item]] <- grass_name  
}  

extender la región para cubrir los datos importados

execGRASS('g.region', parameters= list(raster=paste(grassnames_list, collapse=",")))

escribir un archivo de instrucciones para r.mapcalc

exp_list <- list()
fpname_list <- list()
for (item in rasters) {
  grassname <- substr(basename(item), 1, nchar(basename(item)) -4)
  fp_name   <- paste0(grassname, '_foot')
  exp       <- paste0(fp_name, ' = if(isnull(', grassname, '), null() , 1)')

  exp_list[[item]]    <- exp
  fpname_list[[item]] <- fp_name
}  

write(c(paste(exp_list, collapse = '\n'), 'end'),
      file = paste0(getwd(), '/rmc_instructions.txt'),
      ncolumns = 1)

utilice el archivo de instrucciones a los datos de presencia/ausencia de la huella de rásteres

rmc   <- 'C:/OSGeo4W64/apps/grass/grass-7.0.5/bin/r.mapcalc.exe'
calcs <- 'C:/footprint_test/rmc_instructions.txt'
system2(rmc, paste0('file=', calcs, ' --overwrite'))

vectorise las huellas

vname_list <- list()
for (fp in fpname_list) {
  v_name <- paste0(fp, '_v')
  execGRASS('r.to.vect', 
            flags = 'overwrite',
            parameters = list(input  = fp,
                              output = v_name,
                              type   = 'area'))
 vname_list[[fp]] <- v_name  
}

escribir las huellas de archivo JSON (también es compatible)

outname_list <- list()
for (vn in vname_list) {
  out_name <- paste0(getwd(), '/', vn, '.shp')
  execGRASS('v.out.ogr', 
            flags = 'overwrite', 
            parameters = list(input  = vn,
                              output = out_name,
                              format = 'ESRI_Shapefile'))
  outname_list[[vn]] <- out_name
}

Si desea multipart polígonos

dissname_list <- list()
for (vn in vname_list) {
  diss_name <- paste0(vn, '_sp')
  execGRASS('v.dissolve', 
            flags = 'overwrite',
            parameters = list(input  = vn,
                              column = 'value',
                              output = diss_name))
  dissname_list[[vn]] <- diss_name
}

outname_list2 <- list()
for (dn in dissname_list) {
  out_name2 <- paste0(getwd(), '/', dn, '.shp')
  execGRASS('v.out.ogr', 
            flags = c('m', 'overwrite'), 
            parameters = list(input  = dn,
                              output = out_name2,
                              format = 'ESRI_Shapefile'))
  outname_list2[[dn]] <- out_name2
}

Déjame saber cómo te va con eso.