-
Notifications
You must be signed in to change notification settings - Fork 0
/
21LCMAP_CleanUp.R
88 lines (72 loc) · 3.55 KB
/
21LCMAP_CleanUp.R
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
# this is to mosaic LCMAP data and clip them with HR data
library(raster)
library(rgdal)
setwd("C:/Users/wenjing.xu/Google Drive/RESEARCH/Elk/data_backup")
# # same crs as LCMAP data
target.crs <- CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=23 +lon_0=-96 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0")
# read in shapefiles
AOI <- readOGR("General_AOI.shp")
AOI <- spTransform(AOI, target.crs)
elk_HR_25ind <- readOGR("C:/Users/wenjing.xu/Google Drive/RESEARCH/Elk/Analysis/GYE-Elk-LULC/AllHerd_25ind_w_Agency_clean.shp")
elk_HR_25ind <- spTransform(elk_HR_25ind, target.crs)
###############################################################
############ cleaning LCMAP data ##############################
###############################################################
### mosaic data together #####################################
# mosaicLCMAP <- function(year) {
# dir.i <- paste0("D:/LargeDataSet/LCMAP_1985-2017_GYE/", year, "/")
# filenames <- list.files(dir.i, pattern = "*\\LCPRI.tif$", recursive = TRUE) #only read the primary land cover files for the year
# files <- lapply(paste0(dir.i, filenames), raster)
# files$fun <- mean
# rast.mosaic <- do.call(mosaic,files)
# writeRaster(rast.mosaic, paste0(getwd(),"/LCMAP/LCMAP_LCPRI_", year), format = "GTiff")
# }
# for (i in 1985:2017) {
# mosaicLCMAP(i)
# }
############# clip data into for each HR #######################
### conducting in ArcGIS is a bit easier #######################
setwd("C:\\Users\\wenjing.xu\\Google Drive\\RESEARCH\\Elk\\Analysis\\GYE-Elk-LULC")
LCMAP.dir <- paste0(getwd(), "/Data/LCPRI_Mosaic_GYE_1985-2017/")
filenames <- list.files(LCMAP.dir, pattern = "*\\.tif$")
LCPRI <- lapply(paste0(LCMAP.dir, filenames), raster)
#shorten the list to help with the memory limitation
#LCPRI.1 <- LCPRI[27]
#LCPRI.2 <- LCPRI[17:32]
LCPRI_HR.1 <- lapply(LCPRI.1, crop, y=elk_HR_25ind)
#LCPRI_HR.2 <- lapply(LCPRI.2, crop, y=elk_HR_25ind)
# LCPRI_HR <- lapply(LCPRI, mask, y=elk_HR_25ind)
# save the rasters
writeRaster.LCMAP <- function(x) {
name <- names(x)
name <- substr(name, 7,10)
writeRaster(x, filename = paste0(getwd(), "/Data/LCPRI_HR_25ind_1985-2017/LCPRI_HR_25ind_", name), format ="GTiff")
}
# lapply(LCPRI_HR, writeRaster.LCMAP) # cannot allocate RAM, try for loop
for (i in 1:length(LCPRI_HR.1)) {
LCPRI_HR.i <- LCPRI_HR.1[[i]]
writeRaster.LCMAP(LCPRI_HR.i)
}
########################################################################
############ start here for recalculation ##############################
########################################################################
############ calculate LCLU in each herd winter HR #####################
LCMAP.dir <- paste0(getwd(), "/Data/LCPRI_HR_25ind_1985-2017/")
LCPRI_HR <- list.files(LCMAP.dir, pattern = "*\\.tif$")
countpixel <- function(x, ...){
v <- c(sum(x==1), sum(x==2), sum(x==3), sum(x==4), sum(x==5), sum(x==6), sum(x==7), sum(x==8))
return(v)
}
LCPRI.df <- data.frame()
for (i in 1:length(LCPRI_HR)) {
LCPRI_HR.i <- raster(paste0(LCMAP.dir,LCPRI_HR[[i]]))
cal <- as.data.frame(raster::extract(LCPRI_HR.i, elk_HR_25ind, fun = countpixel))
year <- substr(names(LCPRI_HR.i), 17, 20)
cal$year <- year
cal <- cbind(cal, herd = elk_HR_25ind$id, agency = elk_HR_25ind$AGENCY)
LCPRI.df <- rbind(LCPRI.df, cal)
}
LCPRI.df <- cbind(LCPRI.df[,10:11], LCPRI.df[,9], LCPRI.df[,1:8])
names(LCPRI.df) <- c( "herd", "agency", "year", "developed", "cropland", "grass_shrub", "tree_cover", "water", "wetland", "ice_snow", "barrern")
# LCPRI.df.25ind <- LCPRI.df
# write.csv(LCPRI.df, "LCPRI_allHerds_25ind.csv")