new branch forked

.Rbuildignore

@@ -0,0 +1,2 @@
+^.*\.Rproj$
+^\.Rproj\.user$

.travis.yml

@@ -0,0 +1,5 @@
+# R for travis: see documentation at https://docs.travis-ci.com/user/languages/r
+
+language: R
+sudo: false
+cache: packages

DESCRIPTION

@@ -0,0 +1,28 @@
+Package: svfAHN2point
+Type: Package
+Title: Calculate the Sky view factor from las files
+Version: 1.0
+Date: 2017-06-15
+Authors@R:person("Marieke", "Dirksen", email = "marieke.dirksen@knmi.nl",
+                  role = c("aut", "cre"))
+Maintainer: Who to complain to <marieke.dirksen@knmi.nl>
+Description: More about what it does (maybe more than one line)
+License: What license is it under?
+Imports:  data.table, 
+          doParallel,
+          horizon,
+          leaflet,
+          magrittr,
+          methods,
+          raster,
+          RCurl,
+          rgdal,
+          rgeos,
+          rLiDAR,
+          R.utils,
+          shiny,
+          sp,
+          stringr,
+          stringi
+RoxygenNote: 6.0.1
+

NAMESPACE

@@ -0,0 +1,7 @@
+exportPattern("^[[:alpha:]]+")
+importFrom("methods", "as")
+importFrom("stats", "df")
+importFrom("utils", "write.table")
+importFrom("stringr","str_pad")
+importFrom("raster","crop","extent","rasterize","res","crs","mask")
+importFrom("rgeos","gBuffer")

R/SVF.R

@@ -0,0 +1,69 @@
+#' Sky view factor
+#' 
+#' \code{SVF} Calculates the sky view factor from las point data and save raster files. 
+#' @description calculates the SVF
+#' 
+#' 
+#' @param pointX x-coordinates
+#' @param pointY y-coordinates
+#' @param maxView viewing distance of the sky view factor
+#' @param proj projection
+#' 
+#' @seealso This function calls upon the following functions \code{\link{loadTile}}, \code{\link{makeSpatialDF}},
+#' \code{\link{mergeNeighborTiles}}, \code{\link{makeRaster}},and from the horizon package \code{\link{svf}}.
+
+SVF<-function(pointX, pointY, maxView, proj){
+
+
+  #tileNumberXCoord<-str_pad(as.integer(floor(pointX/1000)*1000), 6, pad = "0")
+  #tileNumberYCoord<-str_pad(as.integer(floor(pointY/1000)*1000), 6, pad = "0")
+
+  tileNumberXCoord<-pointX
+  tileNumberYCoord<-pointY
+
+  mainTile<-loadTile(lazFolder, as.integer(tileNumberXCoord), as.integer(tileNumberYCoord))
+  mainTile<-makeSpatialDF(mainTile,projection = pro)
+  extensionMainTile<-extent(mainTile)
+
+  neighbors<-mergeNeighborTiles(lazFolder, as.integer(tileNumberXCoord), as.integer(tileNumberYCoord), extensionMainTile, maxView, pro)
+
+  rasterizedNeighbors<-lapply(neighbors, makeRaster, Xres, Yres, pro)
+  mergedNeighbors<-do.call(merge, c(rasterizedNeighbors, tolerance =10))
+  rm(neighbors)
+  rm(rasterizedNeighbors)
+  rasterizedMainTile<-makeRaster(mainTile,Xres,Yres,pro)
+  rm(mainTile)
+  #rasterOptions(tolerance = 0.1)
+  fullRaster<-merge(rasterizedMainTile, mergedNeighbors, tolerance = 10)
+  rm(mergedNeighbors)
+  gc()
+
+  r.svf<-svf(fullRaster, nAngles=16, maxDist= maxView, ll=F)
+  out<-crop(r.svf,extent(rasterizedMainTile))
+
+  #to fix the small inconsistencies in the extention (sometimes small approx are introduced in the extent)
+  out<-fix_extent(rasterizedMainTile,out)
+
+  r.b<-brick(rasterizedMainTile,out[[1]])
+  names(r.b)<-c("Z","SVF")
+  r.df<-as.data.frame(r.b,xy=TRUE)
+
+  ##############################################
+  ###############SAVE FUNCTION##################
+  ##############################################
+
+  #Writing the raster file with tile names
+  writeRaster(r.b,filename=paste0(output_dir,
+                                  str_pad(as.integer(floor(tileNumberXCoord/1000)*1000), 6, pad = "0"), "_",
+                                  str_pad(as.integer(floor(tileNumberYCoord/1000)*1000), 6, pad = "0"),".grd"),
+              format="raster",
+              overwrite=TRUE)
+
+  #Writing a table which will be appended if exists
+  write.table(r.df,file="NLSVF.txt",sep=",",row.names = FALSE, append = TRUE, col.names = !file.exists("NLSVF.txt"))
+  ##############################################
+  ##############################################
+  ##############################################
+  rm(r.svf,r.b,r.df)
+  gc()
+}

R/SVFGMStest.R

@@ -0,0 +1,64 @@
+#' Test function
+#' @title Test function for a point dataset of the GMS
+#' @description test function for the GMS dataset to calculate the SVF
+#' @param pointX x-coordinates of the point in RD coordinates
+#' @param pointY y-coordinates of the point in RD coordinates
+#' @param proj projection
+#' @param maxView maximum viewing distance of the svf calculation
+#' @export
+
+SVFGMSTest<-function(pointX, pointY, maxView, proj){
+
+
+  #tileNumberXCoord<-str_pad(as.integer(floor(pointX/1000)*1000), 6, pad = "0")
+  #tileNumberYCoord<-str_pad(as.integer(floor(pointY/1000)*1000), 6, pad = "0")
+
+  tileNumberXCoord<-pointX
+  tileNumberYCoord<-pointY
+
+  mainTile<-loadTileGMSFromLas(lazFolder, as.integer(tileNumberXCoord), as.integer(tileNumberYCoord))
+  mainTile<-makeSpatialDF(mainTile,projection = pro)
+  extensionMainTile<-extent(mainTile)
+
+  neighbors<-mergeNeighborTilesGSMTest(lazFolder, as.integer(tileNumberXCoord), as.integer(tileNumberYCoord), extensionMainTile, maxView, pro)
+
+  rasterizedNeighbors<-lapply(neighbors, makeRaster, Xres, Yres, pro)
+  mergedNeighbors<-do.call(merge, c(rasterizedNeighbors, tolerance =10))
+  rm(neighbors)
+  rm(rasterizedNeighbors)
+  rasterizedMainTile<-makeRaster(mainTile,Xres,Yres,pro)
+  rm(mainTile)
+  #rasterOptions(tolerance = 0.1)
+  fullRaster<-merge(rasterizedMainTile, mergedNeighbors, tolerance = 10)
+  rm(mergedNeighbors)
+  gc()
+
+  r.svf<-svf(fullRaster, nAngles=16, maxDist= maxView, ll=F)
+  out<-crop(r.svf,extent(rasterizedMainTile))
+
+  #to fix the small inconsistencies in the extention (sometimes small approx are introduced in the extent)
+  out<-fix_extent(rasterizedMainTile,out)
+
+  r.b<-brick(rasterizedMainTile,out[[1]])
+  names(r.b)<-c("Z","SVF")
+  r.df<-as.data.frame(r.b,xy=TRUE)
+
+  ##############################################
+  ###############SAVE FUNCTION##################
+  ##############################################
+
+  #Writing the raster file with tile names
+  writeRaster(r.b,filename=paste0(output_dir,
+                                  str_pad(as.integer(floor(tileNumberXCoord/1000)*1000), 6, pad = "0"), "_",
+                                  str_pad(as.integer(floor(tileNumberYCoord/1000)*1000), 6, pad = "0"),"-GMS.grd"),
+              format="raster",
+              overwrite=TRUE)
+
+  #Writing a table which will be appended if exists
+  write.table(r.df,file="NLSVFGMS.txt",sep=",",row.names = FALSE, append = TRUE, col.names = !file.exists("NLSVFGMS.txt"))
+  ##############################################
+  ##############################################
+  ##############################################
+  rm(r.svf,r.b,r.df)
+  gc()
+}

R/SVFWholeNL.R

@@ -0,0 +1,79 @@
+#' Calculate the SVF for the whole Netherlands
+#' @title SVF Netherlands
+#' @param filepath path to the LAS files
+#' @param maxView maximum viewing distance for the sky view factor calculations
+#' @export
+
+SVFWholeNL<-function(filepath, maxView){
+  file<-basename(filepath)
+  splits<-unlist(strsplit(file, c("\\.")))
+  splits<-unlist(strsplit(splits[[1]], "_"))
+
+  tileNumberXCoord<-as.integer(splits[[2]])
+  tileNumberYCoord<-as.integer(splits[[3]])
+
+  mainTile<-loadTileWholeNL(filepath)
+  mainTile<-makeSpatialDF(mainTile,projection = pro)
+  extensionMainTile<-extent(mainTile)
+
+  if(xmax(extensionMainTile)!=xmin(extensionMainTile) & ymax(extensionMainTile)!=ymin(extensionMainTile))
+     {
+
+  neighbors<-mergeNeighborTiles(lazFolder, tileNumberXCoord, tileNumberYCoord, extensionMainTile, maxView, pro)
+
+  #neighbors<-lapply(neighbors,checkCoordinates)
+
+  rasterizedNeighbors<-lapply(neighbors, makeRaster, Xres, Yres, pro)
+  if(length(rasterizedNeighbors)==1){
+    mergedNeighbors<-rasterizedNeighbors[[1]]
+  }
+  else if (length(rasterizedNeighbors)==0){
+    return(-1)  
+  }
+  else{
+    mergedNeighbors<-do.call(merge, c(rasterizedNeighbors, tolerance =10))  
+  }
+  rm(neighbors)
+  rm(rasterizedNeighbors)
+  rasterizedMainTile<-makeRaster(mainTile,Xres,Yres,pro)
+  rm(mainTile)
+  #rasterOptions(tolerance = 0.1)
+  fullRaster<-merge(rasterizedMainTile, mergedNeighbors, tolerance = 10)
+  rm(mergedNeighbors)
+  gc()
+
+  r.svf<-svf(fullRaster, nAngles=16, maxDist= maxView, ll=F)
+  out<-crop(r.svf,extent(rasterizedMainTile))
+
+  #to fix the small inconsistencies in the extention (sometimes small approx are introduced in the extent)
+  out<-fix_extent(rasterizedMainTile,out)
+
+  r.b<-brick(rasterizedMainTile,out[[1]])
+  names(r.b)<-c("Z","SVF")
+  r.df<-as.data.frame(r.b,xy=TRUE)
+
+  ##############################################
+  ###############SAVE FUNCTION##################
+  ##############################################
+
+  #Writing the raster file with tile names
+  writeRaster(r.b,filename=paste0(output_dir,
+                                  str_pad(as.integer(floor(tileNumberXCoord/1000)*1000), 6, pad = "0"), "_",
+                                  str_pad(as.integer(floor(tileNumberYCoord/1000)*1000), 6, pad = "0"),".grd"),
+              format="raster",
+              overwrite=TRUE)
+
+  #Writing a table which will be appended if exists
+  write.table(r.df,file="NLSVF.txt",sep=",",row.names = FALSE, append = TRUE, col.names = !file.exists("NLSVF.txt"))
+
+  }
+  else {
+    print(paste0("tile ", filepath, " with odd extention: ", as.character(extensionMainTile)))
+  }
+  ##############################################
+  ##############################################
+  ##############################################
+  rm(r.svf,r.b,r.df,splits,file)
+  gc()
+}
+

R/calculate_path.R

@@ -0,0 +1,38 @@
+#' Calculate paths from a raster
+#' 
+#' @title Calculate Path
+#' \code{calculate_path} calculates the path for given distance, angles and coordinates from a raster file.
+#' 
+#' @description calculates the path for different angles
+#' @seealso This function is called in \code{\link{lines_from_las}}
+#' @param theta angle or angles in degrees
+#' @param x x-coordinates in meters (RDcoordinates)
+#' @param y y-coordinates in meters (RDcoordinates)
+#' @param distance distance for which the path should be calculated
+#' @param rbuffer raster which is masked in \code{\link{lines_from_las}}
+#' 
+#' 
+#' @export
+
+calculate_path<-function(theta,x=x,y=y,distance,rbuffer){
+  dx<-cos(theta)*distance
+  dy<-sin(theta)*distance
+
+  x1<-x+dx
+  y1<-y+dy
+
+  xy1<-data.frame(x1,y1)
+
+  if(length(xy1$x1)!=length(distance)){
+    warning('Your vectors have different lengths, returning NULL')
+    return(NULL)
+  }
+  #XY<-rbind(c(x,y),xy1) #Bug!
+  XY<-cbind(xy1,distance)
+  XY<-data.frame(XY)
+
+  values<-raster::extract(rbuffer,XY[1:2],buffer=3,fun=mean,na.rm=TRUE)
+  XY<-cbind(XY,values)
+  names(XY)<-c("x","y","distance","height")
+  return(XY)
+}