# Creating the Van Wirdum function #------------------------------------------------------------ #You first need to instal the package "RCurl". #This function only works on R version 2.15.0 or higher. #This action only needs to be preformed once. # For R go to Packages>Install package(s)... # For Rstudio go to Tools>Install packages and type RCurl # Or run once the next code seperatly from the rest of the code (you still need to select a mirror). #install.packages("RCurl") #------------------------------------------------------------ require(RCurl) vanwirdum<-function(EC= 0,IR= 0,group = seq(from=1,to=1,length.out=length(EC)),pointcolor = FALSE,pointsize = 0.5,add2plot = FALSE) { #Test for add2plot if(add2plot != TRUE & add2plot != FALSE){stop("Error: Wrong value for add2plot. Should be TRUE or FALSE")} #Test for EC<>IR if(length(EC)!=length(IR)){stop("Error: EC and IR not same length of data")} #Test for EC if(sum(is.na(EC))>0){stop("Error: EC can not contain NA")} if(sum(EC<0 |EC>10000)>0){stop("Error: EC out of bounderies. Should be EC>= 0 and EC<= 10.000. EC should be measured at 25 degrees")} #Test for IR if(sum(is.na(IR))>0){stop("Error: IR can not contain NA")} if(sum(IR <0|IR>100)>0){stop("Error: IR out of bounderies. Should be IR>= 0 and IR<= 100. IR consists of Ca and CL")} #Test for group if(length(group)!=length(EC)){stop("Error: Not the same number of groups as values")} if(sum(class(group) != "factor" & class(group) != "numeric" & class(group) != "character")>0){stop("Error: Group should be factor,numeric or character")} if(length(unique(group))>44){stop("Error: Group can not have more than 44 dstinctions")} if(sum(is.na(group))>0){stop("Error: Group can not be NA")} #Test for pointcolor if(sum(class(pointcolor) != "numeric" & class(pointcolor) != "character" & pointcolor[1] != FALSE)>0){stop("Error: Pointcolor should be numeric, character or FALSE")} #Pointcolors nr.group<-match(group,unique(group)) col=c(2:7,"gold","darkorange2","purple1","pink2","darkolivegreen1","brown",10:30,3:7,45,100,150,200,250,300) if(pointcolor[1] == FALSE){pointcolor = col[nr.group]} #------------------------------------------------------------------------------------------------------------------------------------------------------------------# #Delete if you don't want to work by internet myCsv <- getURL("http://publicwiki.deltares.nl/download/attachments/90415361/coordinates_LAT_framework.csv?version=1&modificationDate=1353428316053") coordinates_LAT_framework.csv<-textConnection(myCsv) myCsv2 <- getURL("http://publicwiki.deltares.nl/download/attachments/90415361/reference.points.csv?version=1&modificationDate=1353428327557") reference.points.csv<-textConnection(myCsv2) #------------------------------------------------------------------------------------------------------------------------------------------------------------------# #Import LAT framework #setwd() LAT.framework<-read.csv(coordinates_LAT_framework.csv, header = TRUE, sep = ";") #Import reference points reference.points<-read.csv(reference.points.csv, header = TRUE, sep = ";") par(new = add2plot) plot(log10(EC),IR,type="n",xlab = "EC25 (mS/m)", ylab= "IR (%)",cex.main = 2, ylim=c(0,100),axes=F, xlim = log10(c(1,10000))) abline(h=c(20,40,60,80,100),v = c(0,10000),col="lightgrey") abline(h = c(0,100), v=log10(c(1,10,100,1000,10000)),col="lightgrey") lines(log10(LAT.framework$EC25),LAT.framework$IR, lty = "longdash") points(log10(reference.points$EC25),reference.points$IR, pch = 16) text(log10(reference.points$EC25) ,reference.points$IR, labels = reference.points$Name,adj = c(0,1.5), col = "blue") points(log10(EC),IR,col=pointcolor,cex=pointsize,pch=19) box() axis(1,at=log10(c(1,10,100,1000,10000)),label=c(1,10,100,1000,10000)) axis(2) } #------------------------------------------------------------------------------------------------------------------------------------------------------------------# #Delete if you don't want to work by internet myCsv3 <- getURL("http://publicwiki.deltares.nl/download/attachments/90415361/RondeHoep.csv?version=2&modificationDate=1353484644824") RondeHoep.csv<-textConnection(myCsv3) #------------------------------------------------------------------------------------------------------------------------------------------------------------------# #Import RHdata RHdata<-read.csv(RondeHoep.csv, header = TRUE) RHdata vanwirdum() ############################################################################################# #You just created the function vanwirdum() and plotted an empty # #vanwirdum diagram. If you want an example of the use of van wirdum, go to the website # # http://publicwiki.deltares.nl/display/VWD/Home and click on example. Paste this script # #underneath the script of the function. The function only exists # # as long as you don't clear your memory (When closing R you can choose to save the memory).# #############################################################################################