# Creating the Van Wirdum function

#You first need to instal the package "RCurl".
#  	For R go to Packages>Install package(s)...
#	For Rstudio go to Tools>Install packages and type RCurl
#	Or run once the next code
#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()