R/areacalculator.R

#' Works out the area shared between shapefiles, used to convert between GADM version
#' @param shp First shapefile.
#' @param shp2.8 Second shapefile.
#' @param adm The administrative level used.
#' @examples
#' areacalculator(AGO_v1,AGO_v2.8)
#' @keywords internal

areacalculator<-function(shp1,shp2,adm="adm2"){
  
  admval<-if(adm=="adm1") "ID_1" else "ID_2"
  
  #Set up shape 1
  shp1att<-shp1
  shp1att<-as(shp1att, 'SpatialPolygons')
  shp1_ID2<-SpatialPolygonsDataFrame(shp1att,data.frame("ID_2"=shp1[,admval]), match.ID=F)
  #Set up shape 2
  shp2_ID2<-SpatialPolygonsDataFrame(as(shp2, 'SpatialPolygons'), data.frame(field=shp2[,admval]), match.ID=F)
  #Make sure projections are the same
  projection(shp2_ID2)<-projection(shp1_ID2)
  #Intersecting points
  shp1shp2int<-raster::intersect(shp1_ID2, shp2_ID2)
  #Extract areas from polygon objects then attach as attribute
  areas<-data.frame(area=sapply(shp1shp2int@polygons, FUN=function(x) {slot(x, 'area')}))
  row.names(areas)<-sapply(shp1shp2int@polygons, FUN=function(x) {slot(x, 'ID')})
  #Combine attributes info and areas
  attArea<-spCbind(shp1shp2int, areas)
  #Now to apply this stuff
  df1<-data.frame(attArea)
  df1<-df1[order(df1[,1]),]
  #Fine areas
  shp1_areas<-data.frame(id2=shp1_ID2[,admval],area=sapply(shp1_ID2@polygons, FUN=function(x) {slot(x, 'area')}))
  shp1_areas<-aggregate(shp1_areas,by=list(shp1_areas[,admval]),FUN=sum)[,c(1,3)]
  names(shp1_areas)<-c("id2","area")
  
  df1$prop<-sapply(1:dim(df1)[1], function(x){
    oldid2<-df1[,1][x]
    all_area<-shp1_areas[which(shp1_areas$id2 %in% oldid2),"area"]
    df1[x,"area"]/all_area
  })
  
  df1
  
}
arranhamlet/popvac_package documentation built on May 10, 2019, 1:48 p.m.