rm(list=ls())
library(lattice)
setwd("C:\Users\Public\Documents\PC-Progress\Hydrus-1D 4.xx\Examples\Direct\ATNKCatotelmonlyWTD5cm") # Set home directory
getwd() #Check where R is looking
num <- 5000 #Number of realisations
threshold <- -100 #Threshold tension
output <- array(0, dim=c(num,3)) #set output directory
for(i in 1:2)
{INFO <-file.info (paste0("obs_theta_",sprintf("%04.f",i),".out"))
print (INFO$size)}
tension <- read.table(paste0("C:\\Users\\Public\\Documents\\PC-Progress\\Hydrus-1D 4.xx\\Examples\\Direct\\ATNKCatotelmonlyWTD5cm//obs_tension_",sprintf("%04.f",i),".out")) #load table of increasing file numbers
output[i,1] <- min(tension[,2]) #minimum pressure of surface node
output[i,2] <- tension[which.min(tension[,2]),1] #time of minimum tension
x <- dim(tension) #Number of rows of data in output file
value <- 0 #Starting tension
if (min(tension[,2])<threshold){ #if threshold tension is exceed at some point during the simulations
for (k in 1:x[1]) { #Go through each tension measurement in time order to find when the threshold is exceeded
if(tension[k,2]<threshold) { #If pressure is less than the threshold value
output[i,3] <- tension[k,1] # output time when threshold is crossed
break
}
}
}
else{output[i,3] <- 50*24} #output zero if pressure is not exceeded on that run.
rm("tension") # close current tension file ready to open the enxt one
} else{
output[i,2] <- "=NA()"
output[i,3] <- "=NA()"
}
}
colnames(output) <- c("MinPressure","TimeOfMinPressure","TimeOfThresholdPressure") #Define column headings of the output file
write.table(output,file = "Output.txt",sep = " ", quote = FALSE, append = FALSE, na = "NA") #output data
There are 5000 text files, I had written this code 5 years ago and then did not used it for long time. Each of the txt file containes time and theta and tension values i.e. time in the first column and tension in the first, second, third and fourth column. I want the R to out the first column of the theta and tension and time at which the tension reaches -100 cm. I have inserted that info you suggested and worked but not the way I wanted.