Skip to content
GitLab
Commit 6964072e authored by Georges Kunstler's avatar Georges Kunstler
Browse files

moved formating script for France in ./R/

parent 735f4964
master BRIEFCOM ctfs package test.dplyr
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
Showing
with 7 additions and 687 deletions
+7 -687
.gitignore
+ 4
− 1
View file @ 6964072e
output
output_cluster
./data/raw/
data/raw/
data/species.list
figs
*.xls
*.xlsx
......@@ -10,3 +11,5 @@ data/process
*.pdf
*.doc
*.html
*.Rdata
*.rds
CLIMATE.FRANCE.R deleted 100644 → 0
+ 0
− 242
View file @ 735f4964
################################################
################################################
################################################
################################################
### LOAD and PROCEED CLIMATIC DATA
################################################
################################################
################################################
################################################
source("./R/FUN.climate.R")
## Data extracted from C. Piedallu data base with
## Radition with cloud cover for each plots
## Piedallu, C, and J Gegout. 2008. “Efficient Assessment of Topographic Solar Radiation to Improve Plant Distribution Models.” Agricultural and Forest Meteorology 148 (11) (October): 1696–1706. doi:10.1016/j.agrformet.2008.06.001.
## Temperature and precipitation with temperature corrected by elevation based on
## Piedallu, C., J. C. Gégout, V. Perez, F. Lebourgeois, and R. Field. 2012. “Soil Water Balance Performs Better Than Climatic Water Variables in Tree Species Distribution Modelling.” Global Ecology and Biogeography. http://onlinelibrary.wiley.com/doi/10.1111/geb.12012/full.
######
## LOAD PIDEALLU DATA
data.CLIM <- read.csv("./data/raw/DataFrance/climate_piedallu/placettesGK_avec_2011.csv",
header=T,sep='\t',stringsAsFactors=FALSE,dec=",",na.strings="")
data.CLIM$rumkg_500 <- as.numeric(gsub(",",".",data.CLIM$rumkg_500))
#########################################
#########################################
#### COMPUTE CLIMATIC VARIABLES !!!
### NEED TO COMPUTE DDG5
### NEED TO COMPUTE SOIL MAX WATER CONTENT
### NEED TO COMPUTE PET
### NEED TO COMPUTE WATER BUDGET
## remove NA lines
data.CLIM.na <- is.na(data.CLIM[["tmoy6190_1_cor"]])
data.CLIM2 <- data.CLIM[!data.CLIM.na,]
### apply function sgdd
sgdd.vec <- apply(data.CLIM2[,82:93],MARGIN=1,FUN=fun.sgdd)
MeanT.vec <- apply(data.CLIM2[,82:93],MARGIN=1,FUN=mean)
plot(MeanT.vec,sgdd.vec)
dim(data.CLIM)
#################################################
### Max Soil Water Content
### compute based on
## Piedallu, C., J. C. Gégout, A. Bruand, and I. Seynave. 2011. “Mapping Soil Water Holding Capacity over Large Areas to Predict Potential Production of Forest Stands.” Geoderma 160 (3): 355–366.
### read data texture as in Piedallu et al. 2011
data.texture <- read.table("./data/raw/DataFrance/climate_piedallu/texture.txt",header=T)
### load ecological data
load(file="./data/process/ecologie_tot.Rdata")
head(ecologie_tot)
codeprof <- c(2.5,10,20,30,40,50,60,70,80,92.5) # code prof en cm
codecaillou <- c(2.5,10,20,30,40,50,60,70,80,90,97.5) ## code percentage en %
#compute swhc
swhc <- fun.swhc(affroc=ecologie_tot$affroc,cailloux=ecologie_tot$cailloux,text2=ecologie_tot$text2,text1=ecologie_tot$text1,prof2=ecologie_tot$prof2,prof1=ecologie_tot$prof1,codeprof,codecaillou,data.texture)
swhc[is.na(ecologie_tot$prof2)] <- NA
### add to ecological data based
ecologie_tot$SWHC <- swhc
##################################
##################################
#### COMPUTE PET WITH TURC FORMULA
### unit to convert radneb61_ from J/cm2/month in MJ/m2/day /100/30 then to convert in KJ/m2/day as in Nick formula *1000
## RAD and Temperature of the 12 months
radneb61_1to12 <- 66:77
tmoy6190_1_cor.1to12 <- 82:93
### apply in parallel
library(doParallel)
registerDoParallel(cores=6) ## affect automaticaly half of the core detected to the foreach
getDoParWorkers() ## here 8 core so 4 core if want to use more registerDoParallel(cores=6)
PET.matrix <-
foreach(i=1:length(data.CLIM$idp), .combine=rbind) %dopar%
{
fun.PET(i,rad=data.CLIM[,radneb61_1to12],temp=data.CLIM[,tmoy6190_1_cor.1to12])
}
PET.matrix <- PET.matrix[!data.CLIM.na,]
PET.matrix[PET.matrix<0] <- 0 ## affect zero if negative PET
## plot to check PET computed by me vs PET of Christian
## par(mfrow=c(3,4))
## for (i in 1:12) {plot(PET.matrix[,i],data.CLIM2[,i+48]) ; lines(0:100,0:100)}
colnames(PET.matrix) <- paste("PET.cor.",1:12,sep="")
##############
### MERGE CLIMATE and PET
data.CLIM2 <- cbind(data.CLIM2,PET.matrix,sgdd=sgdd.vec)
### MERGE WITH ECOLOGICAL DATA
ecologie.clim <- merge(ecologie_tot,data.CLIM2,by="idp",all.x=T)
dim(ecologie_tot)
dim(ecologie.clim)
names(ecologie.clim)
ecologie.clim2 <- ecologie.clim[!is.na(ecologie.clim$PET.cor.1) &!is.na(ecologie.clim$SWHC),]
###########################################################################################
###### COMPUTE WATER BUDGET
## ### test function
## fun.WaterBudget(i=16,prcp.m=(ecologie.clim2[,paste("prec6190_",1:12,sep="")]),
## PET.m=(ecologie.clim2[,paste("PET.cor.",1:12,sep="")]),
## Ta.m=(ecologie.clim2[,c("tmoy6190_1_cor",paste("tmoy6190_1_cor.",1:11,sep=""))]),
## SWHC.v=(ecologie.clim2[ ,"SWHC"]),n=2)
### apply function in parallel
library(doParallel)
registerDoParallel(cores=6) ## affect automaticaly half of the core detected to the foreach
getDoParWorkers() ## here 8 core so 4 core if want to use more registerDoParallel(cores=6)
WBWS.matrix <-
foreach(i=1:length(ecologie.clim2$idp), .combine=rbind) %dopar%
{
fun.WaterBudget(i,prcp.m=(ecologie.clim2[,paste("prec6190_",1:12,sep="")]),
PET.m=(ecologie.clim2[,paste("PET.cor.",1:12,sep="")]),
Ta.m=(ecologie.clim2[,c("tmoy6190_1_cor",paste("tmoy6190_1_cor.",1:11,sep=""))]),
SWHC.v=(ecologie.clim2[ ,"SWHC"]),n=2)
}
### MERGE all and saved
WBWS.matrix2 <- cbind("idp"=ecologie.clim2$idp,WBWS.matrix)
ecologie.clim.data <- merge(ecologie.clim,WBWS.matrix2,by="idp",all.x=T)
#####
# change tplant not good format
load("./data/process/placette_tot.Rdata")
ecologie.clim.data$tplant <- placette_tot$tplant
### Compute mean T and annual sum of precip
ecologie.clim.data$MAT <- apply(ecologie.clim.data[,114:125],MARGIN=1,FUN=mean)
ecologie.clim.data$SAP <- apply(ecologie.clim.data[,69:80],MARGIN=1,FUN=sum)
saveRDS(ecologie.clim.data,file="./data/process/ecologie.clim.data.rds")
#########################################################
#########################################################
###### FIGURES OF CLIMATIC DATA
#############################
ecologie.clim.data <- readRDS("./data/process/ecologie.clim.data.rds")
names(ecologie.clim.data)
###check climatic data
pdf("./figs/sgdd.tmin.map.pdf")
par(mfrow=c(2,2))
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(heat.colors(10))[cut(ecologie.clim.data$sgdd,quantile(ecologie.clim.data$sgdd,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA,main="SGDD")
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(heat.colors(10))[cut(ecologie.clim.data$tmin6190_min_cor,quantile(ecologie.clim.data$tmin6190_min_cor,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA,main="Tmin")
plot( ecologie.clim.data$tmin6190_min_cor,ecologie.clim.data$sgdd,cex=0.2,xlab="Tmin",ylab="SGDD")
dev.off()
pdf("./figs/Water.map.pdf")
par(mfrow=c(2,2),mar=c(1,1,1,1))
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(colorRampPalette(c("blue", "red"))( 10 ))[cut(ecologie.clim.data$WB.y,quantile(ecologie.clim.data$WB.y,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA)
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(colorRampPalette(c("blue", "red"))( 10 ))[cut(ecologie.clim.data$WB.s,quantile(ecologie.clim.data$WB.s,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA)
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(colorRampPalette(c("blue", "red"))( 10 ))[cut(ecologie.clim.data$WS.y,quantile(ecologie.clim.data$WS.y,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA)
plot(ecologie.clim.data$xl93,ecologie.clim.data$yl93,
col=rev(colorRampPalette(c("blue", "red"))( 10 ))[cut(ecologie.clim.data$WS.s,quantile(ecologie.clim.data$WS.s,probs=(0:10)/10,na.rm=T),labels=F)],
cex=0.2,pty="s",xlab=NA,ylab=NA)
dev.off()
pdf("./figs/Water.var,cor.pdf")
par(mfrow=c(2,2))
plot(ecologie.clim.data[['WB.y']],ecologie.clim.data[['WS.y']],xlab="WB year",ylab="ratio AET/D year",cex=0.2)
plot(ecologie.clim.data[['WB.s']],ecologie.clim.data[['WS.s']],xlab="WB growing season",ylab="ratio AET/D growing season",cex=0.2)
plot(ecologie.clim.data[['WB.s']],ecologie.clim.data[['WB.y']],xlab="WB growing season",ylab="WB growing year",cex=0.2)
plot(ecologie.clim.data[['WS.s']],ecologie.clim.data[['WS.y']],xlab="ratio AET/D growing season",ylab="ratio AET/D year",cex=0.2)
dev.off()
pdf("./figs/sgdd.water.cor.pdf")
plot(ecologie.clim.data[['sgdd']],ecologie.clim.data[['WS.y']],xlab="SGDD",ylab="ratio AET/D year",cex=0.2)
dev.off()
### DO PCA OF CLIMATIC DATA
data.prc <- prcomp(ecologie.clim.data[apply(is.na(ecologie.clim.data[,c('sgdd','tmin6190_min_cor','WB.s','WS.s')]),MARGIN=1,FUN=sum)<1,c('sgdd','tmin6190_min_cor','WB.s','WS.s')],na.action=na.omit,scale=T)
summary(data.prc)
scores <-data.prc$x
pdf("./figs/climate.pca.pdf")
biplot(data.prc,pch=1,xlabs=rep('.',length(data.prc$x[,1])),ylabs=c('SGDD','Tmin','WB','WS'))
text(150,180,labels=paste("Var PC1 ", ((summary(data.prc)))$importance[2,1],sep=""))
text(150,165,labels=paste("Var PC2 ", ((summary(data.prc)))$importance[2,2],sep=""))
dev.off()
##### FIGS on elevation correction of temperature done by Christian
### map of diff between mean altitude of the 1x1km cell and actual elevation of the plots
pdf("./figs/map.error.alti.pdf")
plot(data.CLIM$xl2,data.CLIM$yl2,col="grey",cex=0.3,xlab="X",ylab="Y")
data.temp <- data.CLIM[data.CLIM$diff..altIFN...mnt50 < (-100) | (data.CLIM$diff..altIFN...mnt50>100),]
points(data.temp$xl2,data.temp$yl2,col=c("red","coral3","coral","grey","coral","coral3","red")[cut(data.temp$diff..altIFN...mnt50,breaks=c(-500,-300,-200,-100,100,200,300,500),labels=F)],cex=0.5,pch=1)
dev.off()
##### PLOT HIST OF ALL VARIABLES
pdf('./figs/climatvar.hist.pdf')
par(mfrow=c(3,4))
for (i in 37:48) hist(data.CLIM[,i],xlab='Monthly Precip (mm)',main=names(data.CLIM)[i])
par(mfrow=c(3,4))
for (i in 49:60) hist(data.CLIM[,i],xlab='PET (mm)',main=names(data.CLIM)[i])
par(mfrow=c(3,4))
for (i in 66:77) hist(data.CLIM[,i],xlab='Radiation (?)',main=names(data.CLIM)[i])
par(mfrow=c(3,4))
for (i in 82:93) hist(data.CLIM[,i],xlab='temperature (C)',main=names(data.CLIM)[i])
dev.off()
pdf('./figs/xyplot.Tunco.Tcor.pdf')
par(mfrow=c(3,4))
for (i in 1:12) plot(data.CLIM[,(23:34)[i]],data.CLIM[,(82:94)[i]],xlab='temperature (C) uncorrected',ylab='temperature (C) corrected',
main=names(data.CLIM)[(23:34)[i]])
dev.off()
READ.DATA.NFI.FRANCE.R deleted 100644 → 0
+ 0
− 443
View file @ 735f4964
This diff is collapsed.
merge.data.FRANCE.R
+ 2
− 0
View file @ 6964072e
......@@ -3,6 +3,8 @@
#############################################
### MERGE FRENCH DATA
source("./R/format.function.R")
## source("./R/READ.DATA.NFI.FRANCE.R") need to be run if no already doen (long ...)
## source("./R/CLIMATE.FRANCE.R")
###########################################################################
######################### READ DATA
......
species.list.R
+ 1
− 1
View file @ 6964072e
......@@ -5,7 +5,7 @@
#################################
#################################
#### READ SPECIES LIST I HAVE FOR NOW FOR EACH SITES NON TROPICAL
#### READ SPECIES LIST I HAVE NOW FOR EACH SITES NON TROPICAL
#### (ASSUMING THAT TROPICAL SITES COMES WITH TRAITS DATA
#############################3
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment