Kaynağa Gözat

Lisatud funktsioon 'piirkonnale_ruudustike_lisamine'.

Ardo Kubjas 5 yıl önce
ebeveyn
işleme
3f1d9c3ee9
2 değiştirilmiş dosya ile 135 ekleme ve 112 silme
  1. 2 112
      00_algandmed.R
  2. 133 0
      gpkg/01_piirkonnale_ruudustike_lisamine.R

+ 2 - 112
00_algandmed.R

@@ -29,6 +29,7 @@ objektid <- c("valga", "matsalu", "lahemaa")
 conn <- ruut::db_connect(conf = conf)
 ## Valitud objekti indeks
 i <- 1
+source("gpkg/01_piirkonnale_ruudustike_lisamine.R")
 
 for (i in 1:length(objektid)) {
   ## ---------------- 1. piirkonna piir ------------------
@@ -37,118 +38,7 @@ for (i in 1:length(objektid)) {
   pk <- pk_piir(obj = obj)
   sf::st_geometry(pk) %>% plot()
   gpkg_home <- "/data/gpkg/artiklid/artikkel_210127_valga_matsalu_lahemaa"
-  dsn <- sprintf("%s/%s.gpkg", gpkg_home, obj)
-  input_layer_name <- "piir"
-  input_layer <- sprintf("%s|layername=%s", dsn, input_layer_name)
-  output_layer_name <- "boundarybox_3301"
-  tmp_gpkg_file <- tempfile(fileext = ".gpkg")
-  # write to gpkg
-  class(pk)
-  sf::write_sf(pk,
-    dsn = dsn,
-    layer = "piir", driver = "gpkg", append = FALSE
-  )
-
-  ## ------------ 2. piirkonna 3301 projektsiooniga piirikast --------------
-  ## 2.1  Esmalt leiame milliset 500x500 meetri ruutudega on objekt kaetud.
-  # ruut::qgis_algorithm_search_by_word(str = "grid")
-  algorithm <- "native:creategrid"
-  # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
-  result <- qgisprocess::qgis_run_algorithm(
-    algorithm = algorithm,
-    CRS = "EPSG:3301",
-    # EXTENT = '25.454305528,26.259893095,59.454118579,59.714621582 [EPSG:4326]',
-    EXTENT = input_layer,
-    HOVERLAY = 0,
-    HSPACING = 500,
-    TYPE = 2,
-    VOVERLAY = 0,
-    VSPACING = 500,
-    OUTPUT = qgisprocess::qgis_tmp_file(ext = ".gpkg")
-    # .quiet = TRUE
-  )
-  result
-  pk_grid <- sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT"))
-  # sf::st_geometry(pk_grid) %>% plot()
-  ## 2.2 Leitud ruutudele leiame piirikasti.
-  algorithm <- "qgis:minimumboundinggeometry"
-  # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
-  result <- qgisprocess::qgis_run_algorithm(
-    algorithm = algorithm,
-    FIELD = "",
-    INPUT = pk_grid,
-    TYPE = 3,
-    OUTPUT = tmp_gpkg_file
-    # .quiet = TRUE
-  )
-  result
-  pk_grid_bb <- sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT"))
-  sf::st_geometry(pk_grid_bb) %>% plot(add = T)
-  system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
-
-  ## ------------- 3. piirkonna epk10t (5x5 km) ruudud --------------------
-  # 3.1 kogu ruutvõrgustik
-  ruudud <- c("epk10t_grid", "epk2t_grid", "epk02t_grid")
-  j <- 4
-  for (j in 1:length(ruudud)) {
-    ruut <- ruudud[j]
-    output_layer_name <- ruut
-    conf <- ruut::get_config()
-    pg <- ruut::construct_ogr2ogr_PG_connect_str(conf = conf)
-    input <- sprintf(
-      'postgres://dbname=\'%s\' host=%s port=%s user=\'%s\' sslmode=%s password=\'%s\' key=\'fid\' srid=4326 type=Polygon checkPrimaryKeyUnicity=\'1\' table=\"%s\".\"%s\" (geometry)',
-      conf$dbname, conf$host, conf$port, conf$user, conf$sslmode, conf$password,
-      "maaamet", ruut
-    )
-    # ruut::qgis_algorithm_search_by_word(str = "extract")
-    algorithm <- "native:extractbylocation"
-    # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
-    result <- qgisprocess::qgis_run_algorithm(
-      algorithm = algorithm,
-      INPUT = input,
-      INTERSECT = sprintf("%s|layername=%s", dsn, "boundarybox_3301"),
-      OUTPUT = tmp_gpkg_file,
-      PREDICATE = c(0, 1)
-      # .quiet = TRUE
-    )
-    result
-    assign(ruut, sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT")))
-    system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
-  }
-  sf::st_geometry(epk10t_grid) %>% plot(add = T)
-
-  # 3.2 ainult piirkonnaga seotud ning informatsiooni sisaldav ruutvõrgustik
-  ruudud <- c("epk10t", "epk2t")
-  j <- 4
-  for (j in 1:length(ruudud)) {
-    ruut <- ruudud[j]
-    output_layer_name <- ruut
-    conf <- ruut::get_config()
-    pg <- ruut::construct_ogr2ogr_PG_connect_str(conf = conf)
-    input <- sprintf(
-      'postgres://dbname=\'%s\' host=%s port=%s user=\'%s\' sslmode=%s password=\'%s\' key=\'fid\' srid=4326 type=Polygon checkPrimaryKeyUnicity=\'1\' table=\"%s\".\"%s\" (geometry)',
-      conf$dbname, conf$host, conf$port, conf$user, conf$sslmode, conf$password,
-      "maaamet", ruut
-    )
-    # ruut::qgis_algorithm_search_by_word(str = "extract")
-    algorithm <- "native:extractbylocation"
-    # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
-    result <- qgisprocess::qgis_run_algorithm(
-      algorithm = algorithm,
-      INPUT = input,
-      INTERSECT = input_layer,
-      OUTPUT = tmp_gpkg_file,
-      PREDICATE = c(0, 1)
-      # .quiet = TRUE
-    )
-    result
-    assign(ruut, sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT")))
-    system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
-  }
-  sf::st_geometry(epk10t) %>% plot(add = T)
-  ## ---------------------- vaata layer'id ----------------------
-  # Vaata layer'eid
-  sf::st_layers(dsn = dsn)
+  piirkonnale_ruudustike_lisamine(obj = obj, pk = pk, gpkg_home = gpkg_home)
 }
 
 

+ 133 - 0
gpkg/01_piirkonnale_ruudustike_lisamine.R

@@ -0,0 +1,133 @@
+#' Piirkonna ruudustike leidmine
+#'
+#' Etteantud piirkonna geomeetrilise piirjoone ('piir') järele leitakse selle piirkonna EPSG:3301 projektsioonile vastav piirikast ('boundarybox_3301'). Piirikasti järele leitakse selle kastiga kaetud 5x5 km ('epk10t_grid'), 1x1 km ('epk2t_grid') ja 100x100 meetrit ('epk02t_grid') ruudustikud. Samuti leitakse piirkonna piiri sees olevad ja statistilist ning muud infot sisaldavad 5x5 km ('epk10t') ja 1x1 km ('epk2t') sisaldavad ruudu. Andmed salvestatakse GPKG faili kihtidena. Faili nimeks on objekti nimi.
+#'
+#' @param obj str Objekti nimi.
+#' @param pk Piirkonna geomeetriline joon.
+#' @param gpkg_home path Salvestatavate GPKG faili asukoht.
+#'
+#'
+piirkonnale_ruudustike_lisamine <- function(obj = NULL, pk = NULL, gpkg_home = "/tmp") {
+  if (is.null(pk) || !sf::st_is_valid(pk)) {
+    cat("\nPalun kontrolli geomeetrilise kujundi õigsust.\n")
+    return()
+  }
+  if (is.null(obj)) {
+    cat("\nPuudu on objekti nimi.\n")
+    return()
+  }
+  cat(sprintf("\nAlgparameetrid: objekti nimi %s ja GPKG faili kataloog %s\n", obj, gpkg_home))
+  dsn <- sprintf("%s/%s.gpkg", gpkg_home, obj)
+  input_layer_name <- "piir"
+  input_layer <- sprintf("%s|layername=%s", dsn, input_layer_name)
+  output_layer_name <- "boundarybox_3301"
+  tmp_gpkg_file <- tempfile(fileext = ".gpkg")
+  # write to gpkg
+  class(pk)
+  sf::write_sf(pk,
+    dsn = dsn,
+    layer = "piir", driver = "gpkg", append = FALSE
+  )
+
+  ## ------------ 2. piirkonna 3301 projektsiooniga piirikast --------------
+  ## 2.1  Esmalt leiame milliset 500x500 meetri ruutudega on objekt kaetud.
+  # ruut::qgis_algorithm_search_by_word(str = "grid")
+  algorithm <- "native:creategrid"
+  # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
+  result <- qgisprocess::qgis_run_algorithm(
+    algorithm = algorithm,
+    CRS = "EPSG:3301",
+    # EXTENT = '25.454305528,26.259893095,59.454118579,59.714621582 [EPSG:4326]',
+    EXTENT = input_layer,
+    HOVERLAY = 0,
+    HSPACING = 500,
+    TYPE = 2,
+    VOVERLAY = 0,
+    VSPACING = 500,
+    OUTPUT = qgisprocess::qgis_tmp_file(ext = ".gpkg")
+    # .quiet = TRUE
+  )
+  result
+  pk_grid <- sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT"))
+  # sf::st_geometry(pk_grid) %>% plot()
+  ## 2.2 Leitud ruutudele leiame piirikasti.
+  algorithm <- "qgis:minimumboundinggeometry"
+  # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
+  result <- qgisprocess::qgis_run_algorithm(
+    algorithm = algorithm,
+    FIELD = "",
+    INPUT = pk_grid,
+    TYPE = 3,
+    OUTPUT = tmp_gpkg_file
+    # .quiet = TRUE
+  )
+  result
+  pk_grid_bb <- sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT"))
+  sf::st_geometry(pk_grid_bb) %>% plot(add = T)
+  system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
+
+  ## ------------- 3. piirkonna epk10t (5x5 km) ruudud --------------------
+  # 3.1 kogu ruutvõrgustik
+  ruudud <- c("epk10t_grid", "epk2t_grid", "epk02t_grid")
+  j <- 4
+  for (j in 1:length(ruudud)) {
+    ruut <- ruudud[j]
+    output_layer_name <- ruut
+    conf <- ruut::get_config()
+    pg <- ruut::construct_ogr2ogr_PG_connect_str(conf = conf)
+    input <- sprintf(
+      'postgres://dbname=\'%s\' host=%s port=%s user=\'%s\' sslmode=%s password=\'%s\' key=\'fid\' srid=4326 type=Polygon checkPrimaryKeyUnicity=\'1\' table=\"%s\".\"%s\" (geometry)',
+      conf$dbname, conf$host, conf$port, conf$user, conf$sslmode, conf$password,
+      "maaamet", ruut
+    )
+    # ruut::qgis_algorithm_search_by_word(str = "extract")
+    algorithm <- "native:extractbylocation"
+    # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
+    result <- qgisprocess::qgis_run_algorithm(
+      algorithm = algorithm,
+      INPUT = input,
+      INTERSECT = sprintf("%s|layername=%s", dsn, "boundarybox_3301"),
+      OUTPUT = tmp_gpkg_file,
+      PREDICATE = c(0, 1)
+      # .quiet = TRUE
+    )
+    result
+    assign(ruut, sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT")))
+    system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
+  }
+  sf::st_geometry(epk10t_grid) %>% plot(add = T)
+
+  # 3.2 ainult piirkonnaga seotud ning informatsiooni sisaldav ruutvõrgustik
+  ruudud <- c("epk10t", "epk2t")
+  j <- 4
+  for (j in 1:length(ruudud)) {
+    ruut <- ruudud[j]
+    output_layer_name <- ruut
+    conf <- ruut::get_config()
+    pg <- ruut::construct_ogr2ogr_PG_connect_str(conf = conf)
+    input <- sprintf(
+      'postgres://dbname=\'%s\' host=%s port=%s user=\'%s\' sslmode=%s password=\'%s\' key=\'fid\' srid=4326 type=Polygon checkPrimaryKeyUnicity=\'1\' table=\"%s\".\"%s\" (geometry)',
+      conf$dbname, conf$host, conf$port, conf$user, conf$sslmode, conf$password,
+      "maaamet", ruut
+    )
+    # ruut::qgis_algorithm_search_by_word(str = "extract")
+    algorithm <- "native:extractbylocation"
+    # cat(qgisprocess::qgis_show_help(algorithm = algorithm))
+    result <- qgisprocess::qgis_run_algorithm(
+      algorithm = algorithm,
+      INPUT = input,
+      INTERSECT = input_layer,
+      OUTPUT = tmp_gpkg_file,
+      PREDICATE = c(0, 1)
+      # .quiet = TRUE
+    )
+    result
+    assign(ruut, sf::read_sf(qgisprocess::qgis_output(result, "OUTPUT")))
+    system(sprintf("ogr2ogr -f GPKG -overwrite  %s %s -nln %s", dsn, tmp_gpkg_file, output_layer_name))
+  }
+  sf::st_geometry(epk10t) %>% plot(add = T)
+  ## ---------------------- vaata layer'id ----------------------
+  # Vaata layer'eid
+  sf::st_layers(dsn = dsn)
+}
+# piirkonnale_ruudustike_lisamine(obj = NULL, pk = NULL, gpkg_home = "/tmp")