writeOutSamples.r 3.13 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
#' Sample Collection for Habitat Types
#'
#'Writes out a set of samples (SpatialPointsDataFrame) into ESRI shapefiles or a GeoJSON file for a selected habitat type. Each point represents a valid sample location that identifies the selected habitat type.
#'
#' @param inPath file path (character) for results of habitat type sampling and probability mapping (same as outPath from function multi_Class_Sampling)
#' @param step step number (numeric)
#' @param className name (character) of habitat type for which samples should be selected
#' @param output_format format (character) of output; whether shp (default) or geojson
#'
#' @return ESRI shapefiles/GeoJSON with name: SamplePoints_step_classname.shp/SamplePoints_step_classname.geojson
#' 1) Point Shape/GeoJSON represents pixel that belong to selected habitat type and can be used as reference for further model building.
#'    ESRI shapefiles have the same CRS as the input raster. GeoJSON files are in the standard CRS of GeoJSON (EPSG:4326).
#'
#'
#' @export

###write out selected samples
writeOutSamples <- function (inPath, step, className, output_format = c("shp", "geojson")) {

  paste(inPath, "step_", step, "_", className, ".tif", sep = "")
  run1 <- get(load(paste(inPath, "Run", step, sep = "")))
  load(paste(inPath, "threshold_step_", step, sep = ""))
  dummy_sample <-
    raster::raster(paste(inPath, "step_", step, "_", className, ".tif", sep =
                           ""))

  length_threshold <- length(threshold)
  thres <- threshold[length_threshold]
  dummy_sample[dummy_sample < thres] <- NA
  dummy_sample[dummy_sample >= thres] <- 1

  collect <- list()
  j <- 0

  ###extract only class samples
  for (i in 1:length(run1@ref_samples)) {
    if (length(dim(run1@ref_samples[[i]])) != 0)
    {
      if (is.na(run1@switch[i]) == F) {
        j = j + 1
        collect[[j]] <-
          run1@ref_samples[[i]][which(run1@ref_samples[[i]]@data == 1), ]
      } else
      {
        j = j + 1
        collect[[j]] <-
          run1@ref_samples[[i]][which(run1@ref_samples[[i]]@data == 2), ]
      }
    }
  }
  
  result <- do.call(rbind, collect)

  res <- raster::extract(dummy_sample, result)
  if (length(which(is.na(res))) > 0) {
    res <- result[-which(is.na(res)), ]
  }
  
  output_format <- match.arg(output_format)
  
  if (output_format == "geojson") {
    crs_dummy <- sp::proj4string(dummy_sample)
    raster::crs(res) <- crs_dummy
    res <- sp::spTransform(res, sp::CRS("+init=epsg:4326"))

    rgdal::writeOGR(
      obj = res,
      layer = paste("SamplePoints_step_", step, "_", className, sep = ""),
      dsn = paste(inPath, "SamplePoints_step_", step, "_", className, ".geojson", sep = ""),
      driver = "GeoJSON",
      check_exists = TRUE,
      overwrite_layer = TRUE
    ) 
   } else {
      crs_dummy <- sp::proj4string(dummy_sample)
      raster::crs(res) <- crs_dummy
      
      rgdal::writeOGR(
        obj = res,
        layer = paste("SamplePoints_step_", step, "_", className, sep = ""),
        dsn = paste(inPath, "SamplePoints_step_", step, "_", className, ".shp", sep = ""),
        driver = "ESRI Shapefile",
        check_exists = TRUE,
        overwrite_layer = TRUE
      )
     }
}