adjustments to get the sample points saved as geojson in epsg 4326; work in progress

R-package/R/WriteOutSamples.r

@@ -12,48 +12,104 @@
 #'
 #' @export
 
-###write out selected samples for step = 6 using
+###write out selected samples
 write_Out_Samples <- function (inPath, step, className) {
+  
+  print("function started")
+  
   paste(inPath, "step_", step, "_", className, ".tif", sep = "")
   run1 <- get(load(paste(inPath, "Run", step, sep = "")))
   load(paste(inPath, "threshold_step_", step, sep = ""))
-  dummy <-
+  dummy_sample <-
     raster::raster(paste(inPath, "step_", step, "_", className, ".tif", sep =
                            ""))
-  thres <- threshold[6]
-  dummy[dummy < thres] <- NA
-  dummy[dummy >= thres] <- 1
+  
+  print("problem zone 0")
+  
+  length_threshold <- length(threshold)
+  thres <- threshold[length_threshold]
+  dummy_sample[dummy_sample < thres] <- NA
+  dummy_sample[dummy_sample >= thres] <- 1
+  
+  print("problem zone 1")
 
   collect <- list()
   j <- 0
 
   ###extract only class samples
   for (i in 1:length(run1@ref_samples)) {
+    print("problem zone 1.1")
     if (length(dim(run1@ref_samples[[i]])) != 0)
     {
+      print("problem zone 1.2")
       if (is.na(run1@switch[i]) == F) {
         j = j + 1
         collect[[j]] <-
           run1@ref_samples[[i]][which(run1@ref_samples[[i]]@data == 1), ]
       } else
       {
+        print("problem zone 1.3")
         j = j + 1
         collect[[j]] <-
           run1@ref_samples[[i]][which(run1@ref_samples[[i]]@data == 2), ]
       }
     }
   }
+  
+  print("problem zone 2")
+  
   result <- do.call(rbind, collect)
 
-  res <- raster::extract(dummy, result)
+  res <- raster::extract(dummy_sample, result)
   if (length(which(is.na(res))) > 0) {
     res <- result[-which(is.na(res)), ]
   }
+  
+  print("problem zone 3")
+  ### save the points as .shp with crs
+  # crs_dummy <- proj4string(dummy_sample)
+  # proj4string(res) <- crs_dummy
+  # print(str(res))
+  # print(proj4string(res))
+  # 
+  # rgdal::writeOGR(
+  #   res,
+  #   layer = "result",
+  #   dsn = paste(inPath, "RefHaSa_", className, "_", step, ".shp", sep = ""),
+  #   driver = "ESRI Shapefile",
+  #   check_exists = TRUE,
+  #   overwrite_layer = TRUE
+  # )
 
+  #### this approach is for a shape file in a Geopackage (.dpkg) container
+  # print(proj4string(res))
+  # 
+  # crs_dummy <- proj4string(dummy_sample)
+  # proj4string(res) <- crs_dummy
+  # 
+  # print(proj4string(res))
+  # 
+  # rgdal::writeOGR(
+  #   res,
+  #   layer = paste("RefHaSa_step_", step, "_", className, sep = ""),
+  #   dsn = paste(inPath, "RefHaSa_step_", step, "_", className, ".gpkg", sep = ""),
+  #   driver = "ESRI Shapefile",
+  #   check_exists = TRUE,
+  #   overwrite_layer = TRUE
+  # )
+  
+  crs_dummy <- proj4string(dummy_sample)
+  crs(res) <- crs_dummy
+  res <- spTransform(res, CRS("+proj=longlat +datum=WGS84 +init=epsg:4326"))
+  
   rgdal::writeOGR(
     res,
     layer = "result",
-    dsn = paste(inPath, "RefHaSa_", className, "_", step, ".shp", sep = ""),
-    driver = "ESRI Shapefile"
+    dsn = paste(inPath, "RefHaSa_step_", step, "_", className, ".geojson", sep = ""),
+    driver = "GeoJSON",
+    check_exists = TRUE,
+    overwrite_layer = TRUE
   )
+
+  print("finished")
 }