Nothing
rsmlToDART <- function(rsml.path, final.date, connect){
# Create the plant object
rsml <- xmlToList(xmlParse(rsml.path))
plants <- rsml$scene
resolution<-rsml$metadata$resolution
if (is.null(resolution)==TRUE){resolution<-1}
unitlength1<-rsml$metadata$unit
# Find unittime and unitdiameter in 'property-definition'
unittime1<-NULL
unitdiameter1<-NULL
#Check unittime
property<-rsml$metadata$'property-definition'
if (is.character(final.date)==TRUE) {
if (is.null(property)==FALSE){
for (i in 1:length(property)){if (TRUE %in% (final.date %in% property[[i]]$label)){unittime1<-as.character(property[[i]]$unit)}}}
if (is.null(unittime1)==TRUE){message(paste("No time unit found in ", sub(basename(rsml.path), pattern=".rsml", replacement=""), " metadata (property-definition)", sep=""))}}
#Check unitdiameter
if (is.null(property)==FALSE){
for (i in 1:length(property)){if (property[[i]]$label=="diameter"|property[[i]]$label=="rootDiameter"){unitdiameter1<-as.character(property[[i]]$unit)}}}
if (is.null(unitdiameter1)==TRUE){message(paste("No diameter unit found in ", sub(basename(rsml.path), pattern=".rsml", replacement=""), " metadata (property-definition)", sep=""))}
# Create LIE and RAC files for each root system
n <- 0 # Initialise number of LIE/RAC files (as 1 RSML file can contain more than 1 root system)
lie.all<-list()
rac.all<-list()
tps.all<-list()
for (r0 in plants){# For each plant
#Calculate total number of nodes in scene (! add additional nodes if connect=TRUE)
nodes<-length(grep(x=names(rapply(r0, length, how="unlist")), pattern="point"))
#Calculate total number of roots in scene
roots<-length(grep(x=names(rapply(r0, length, how="unlist")), pattern="root..attrs"))
if (connect==TRUE){nodes<-nodes+roots}
n<-n+1 #Add one unit for each root system
r <- 0 # Initialise the number of roots consisting a root system
timeserie<-FALSE #Does the rsml file contain time series data? (Ex: the root system age for each node)
#For each plant root system, create LIE, RAC and TPS files
lie<-matrix(nrow=nodes, ncol=13)
rac<-matrix(nrow=roots, ncol=6)
lie.lines<-0
for (r1 in r0){ # For each first-order root
if (class(r1)=="list"){
ns<-r1$geometry$polyline
# Extract plant age when nodes appeared for LIE and TPS files
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r1)){
age<-r1$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
timeserie<-TRUE
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r1)){
age<-r1$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
length1<-length(ns)
r<-r+1 #Add one unit for each root
if (r==1){#For the first first-order root only
currentMother<-r
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[1:length1,1:6]<-c(c(1:length1),rep(1, length1),rep(0, length1), rep(0, length1), 0:(length1-1), 2:(length1+1))}
if (timeserie==TRUE) {lie[1:length1,1:6]<-c(c(1:length1),time,rep(0, length1), rep(0, length1), 0:(length1-1), 2:(length1+1))}
lie[1:length1,11]<-r
if (is.null(diameter)==TRUE){} else {lie[1:length1, 12]<-diameter}
lie[1:length1,13]<-1
#c(X,Y,Z)
lie[1:length1,7]<-sapply(ns, xnodes)
lie[1:length1,8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[1:length1,9]<-sapply(ns, znodes)} else {lie[1:length1,9]<-0}
#c(dist)
lie[1:length1, 10]<-c(0, cumsum(sqrt((diff(lie[1:length1, 7]))^2+(diff(lie[1:length1, 8]))^2+(diff(lie[1:length1, 9]))^2)))
#For the first point of the first order root
if (timeserie==FALSE) {lie[1,c(2:3)]<-c(0,1)} else {lie[1,3]<-1}
start1<-as.numeric(lie[1, 1])
#For the last point of the first order root
lie[length1,c(4,6)]<-c(1,0)
stop1<-as.numeric(lie[length1, 1])
# Fill RAC file for the first order root
#c(Root, Mother, Ord, DBase, DApp, Length)
cumulDist<-sum(sqrt((diff(lie[1:length1, 7]))^2+(diff(lie[1:length1, 8]))^2+(diff(lie[1:length1, 9]))^2))
rac[r,1:6]<-c(r,-1,1,0,0,cumulDist)
lie.lines<-lie.lines+length1}
else {#For the other first-order roots only
currentMother<-r
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length1),1:6]<-c(c((lie.lines+1):(lie.lines+length1)),rep(1, length1),rep(0, length1), rep(0, length1), c(0,(lie.lines+1):(lie.lines+length1-1)), (lie.lines+2):(lie.lines+length1+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length1),1:6]<-c(c((lie.lines+1):(lie.lines+length1)),time,rep(0, length1), rep(0, length1), c(0,(lie.lines+1):(lie.lines+length1-1)), (lie.lines+2):(lie.lines+length1+1))}
lie[(lie.lines+1):(lie.lines+length1),11]<-r
if (is.null(diameter)==TRUE){} else {lie[(lie.lines+1):(lie.lines+length1),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length1),13]<-1
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length1),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length1),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length1),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length1),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length1), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length1), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length1), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length1), 9]))^2)))
#For the first point of the first order root
if (timeserie==FALSE) {lie[lie.lines+1,c(2:3)]<-c(0,1)} else {lie[lie.lines+1,3]<-1}
start1<-as.numeric(lie[lie.lines+1, 1])
#For the last point of the first order root
lie[(lie.lines+length1),c(4,6)]<-c(1,0)
stop1<-as.numeric(lie[lie.lines+length1, 1])
# Fill RAC file for the first order root
#c(Root, Mother, Ord, DBase, DApp, Length)
cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length1), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length1), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length1), 9]))^2))
rac[r,1:6]<-c(r,-1,1,0,0,cumulDist)
lie.lines<-lie.lines+length1}
#----------------------------------------------------------------------------------------------------
# If there are lateral roots
if ("root" %in% names(r1)){
for (r2 in r1){# For each 2-order root
if ("geometry" %in% names(r2)){
r<-r+1
ns <- r2$geometry$polyline
length2<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r2)){
age<-r2$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r2)){
age<-r2$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length2),1:6]<-c((lie.lines+1):(lie.lines+length2), rep(1, length2), rep(0, length2), rep(0, length2), lie.lines:(lie.lines+length2-1), (lie.lines+2):(lie.lines+length2+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length2),1:6]<-c((lie.lines+1):(lie.lines+length2), time, rep(0, length2), rep(0, length2), lie.lines:(lie.lines+length2-1), (lie.lines+2):(lie.lines+length2+1))}
lie[(lie.lines+1):(lie.lines+length2),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length2),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length2),13]<-2 #Order 2
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length2),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length2),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length2),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length2),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length2), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length2), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start1:stop1, 7]==lie[lie.lines+1, 7] & lie[start1:stop1, 8]==lie[lie.lines+1, 8] & lie[start1:stop1, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start1+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start1+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start1:(stop1-1), 7]-lie[lie.lines+1, 7])^2+(lie[start1:(stop1-1), 8]-lie[lie.lines+1, 8])^2+(lie[start1:(stop1-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start1:stop1, 7])*(lie[start1:(stop1-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start1:stop1, 8])*(lie[start1:(stop1-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start1:stop1, 9])*(lie[start1:(stop1-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start1:stop1, 7])^2+diff(lie[start1:stop1, 8])^2+diff(lie[start1:stop1, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start1+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start1+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length2), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length2), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start1:stop1, 7])*t+lie[start1:(stop1-1), 7]
yn<-diff(lie[start1:stop1, 8])*t+lie[start1:(stop1-1), 8]
zn<-diff(lie[start1:stop1, 9])*t+lie[start1:(stop1-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
dist1<-min(dist1)
index<-as.numeric(match(min(dist1), dist1)[1])}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start1+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start1+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length2), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length2), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length2), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length2), 10]
dbase<-lie[start1+index-1, 10]+distance3D(x1=lie[start1+index-1, 7], x2=xn[index], y1=lie[start1+index-1, 8], y2=yn[index], z1=lie[start1+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length2+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start1+index-1, 11], NA, 1)
lie[1:(lie.lines+length2+1),]<-lie[order(lie[1:(lie.lines+length2+1), 11], lie[1:(lie.lines+length2+1), 10]),]
length1<-length1+1
lie.lines<-lie.lines+1
stop1<-stop1+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length2),5]<-match(lie[1:(lie.lines+length2),5], lie[1:(lie.lines+length2),1])
lie[1:(lie.lines+length2),6]<-match(lie[1:(lie.lines+length2),6], lie[1:(lie.lines+length2),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length2),1]<-match(lie[1:(lie.lines+length2),1], lie[1:(lie.lines+length2),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start2<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length2,c(4,6)]<-c(1, 0)
stop2<-as.numeric(lie[lie.lines+length2, 1])
# Fill RAC file for the 2-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length2), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length2), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length2), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother, 2, dbase, 1, cumulDist)
currentMother2<-rac[r,1]
lie.lines<-lie.lines+length2
#----------------------------------------------------------------------------------------------------
if ("root" %in% names(r2)){
for (r3 in r2){# For each 3-order root
if ("geometry" %in% names(r3)){
r<-r+1
ns <- r3$geometry$polyline
length3<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r3)){
age<-r3$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r3)){
age<-r3$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length3),1:6]<-c((lie.lines+1):(lie.lines+length3), rep(1, length3), rep(0, length3), rep(0, length3), lie.lines:(lie.lines+length3-1), (lie.lines+2):(lie.lines+length3+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length3),1:6]<-c((lie.lines+1):(lie.lines+length3), time, rep(0, length3), rep(0, length3), lie.lines:(lie.lines+length3-1), (lie.lines+2):(lie.lines+length3+1))}
lie[(lie.lines+1):(lie.lines+length3),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length3),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length3),13]<-3
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length3),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length3),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length3),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length3),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length3), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length3), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start2:stop2, 7]==lie[lie.lines+1, 7] & lie[start2:stop2, 8]==lie[lie.lines+1, 8] & lie[start2:stop2, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start2+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start2+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start2:(stop2-1), 7]-lie[lie.lines+1, 7])^2+(lie[start2:(stop2-1), 8]-lie[lie.lines+1, 8])^2+(lie[start2:(stop2-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start2:stop2, 7])*(lie[start2:(stop2-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start2:stop2, 8])*(lie[start2:(stop2-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start2:stop2, 9])*(lie[start2:(stop2-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start2:stop2, 7])^2+diff(lie[start2:stop2, 8])^2+diff(lie[start2:stop2, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start2+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start2+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length3), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length3), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start2:stop2, 7])*t+lie[start2:(stop2-1), 7]
yn<-diff(lie[start2:stop2, 8])*t+lie[start2:(stop2-1), 8]
zn<-diff(lie[start2:stop2, 9])*t+lie[start2:(stop2-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
index<-as.numeric(match(min(dist1), dist1)[1])
dist1<-min(dist1)}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start2+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start2+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length3), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length3), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length3), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length3), 10]
dbase<-lie[start2+index-1, 10]+distance3D(x1=lie[start2+index-1, 7], x2=xn[index], y1=lie[start2+index-1, 8], y2=yn[index], z1=lie[start2+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length3+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start2+index-1, 11], NA, 2)
lie[1:(lie.lines+length3+1),]<-lie[order(lie[1:(lie.lines+length3+1), 11], lie[1:(lie.lines+length3+1), 10]),]
length2<-length2+1
lie.lines<-lie.lines+1
stop2<-stop2+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length3),5]<-match(lie[1:(lie.lines+length3),5], lie[1:(lie.lines+length3),1])
lie[1:(lie.lines+length3),6]<-match(lie[1:(lie.lines+length3),6], lie[1:(lie.lines+length3),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length3),1]<-match(lie[1:(lie.lines+length3),1], lie[1:(lie.lines+length3),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start3<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length3,c(4,6)]<-c(1, 0)
stop3<-as.numeric(lie[lie.lines+length3, 1])
# Fill RAC file for the 3-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length3), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length3), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length3), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother2, 3, dbase, 1, cumulDist)
currentMother3<-rac[r,1]
lie.lines<-lie.lines+length3
#----------------------------------------------------------------------------------------------------
if ("root" %in% names(r3)){
for (r4 in r3){# For each 4-order root
if ("geometry" %in% names(r4)){
r<-r+1
ns <- r4$geometry$polyline
length4<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r4)){
age<-r4$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r4)){
age<-r4$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length4),1:6]<-c((lie.lines+1):(lie.lines+length4), rep(1, length4), rep(0, length4), rep(0, length4), lie.lines:(lie.lines+length4-1), (lie.lines+2):(lie.lines+length4+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length4),1:6]<-c((lie.lines+1):(lie.lines+length4), time, rep(0, length4), rep(0, length4), lie.lines:(lie.lines+length4-1), (lie.lines+2):(lie.lines+length4+1))}
lie[(lie.lines+1):(lie.lines+length4),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length4),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length4),13]<-4
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length4),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length4),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length4),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length4),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length4), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length4), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start3:stop3, 7]==lie[lie.lines+1, 7] & lie[start3:stop3, 8]==lie[lie.lines+1, 8] & lie[start3:stop3, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start3+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start3+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start3:(stop3-1), 7]-lie[lie.lines+1, 7])^2+(lie[start3:(stop3-1), 8]-lie[lie.lines+1, 8])^2+(lie[start3:(stop3-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start3:stop3, 7])*(lie[start3:(stop3-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start3:stop3, 8])*(lie[start3:(stop3-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start3:stop3, 9])*(lie[start3:(stop3-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start3:stop3, 7])^2+diff(lie[start3:stop3, 8])^2+diff(lie[start3:stop3, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start3+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start3+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length4), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length4), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start3:stop3, 7])*t+lie[start3:(stop3-1), 7]
yn<-diff(lie[start3:stop3, 8])*t+lie[start3:(stop3-1), 8]
zn<-diff(lie[start3:stop3, 9])*t+lie[start3:(stop3-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
index<-as.numeric(match(min(dist1), dist1)[1])
dist1<-min(dist1)}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start3+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start3+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length4), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length4), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length4), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length4), 10]
dbase<-lie[start3+index-1, 10]+distance3D(x1=lie[start3+index-1, 7], x2=xn[index], y1=lie[start3+index-1, 8], y2=yn[index], z1=lie[start3+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length4+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start3+index-1, 11], NA, 3)
lie[1:(lie.lines+length4+1),]<-lie[order(lie[1:(lie.lines+length4+1), 11], lie[1:(lie.lines+length4+1), 10]),]
length3<-length3+1
lie.lines<-lie.lines+1
stop3<-stop3+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length4),5]<-match(lie[1:(lie.lines+length4),5], lie[1:(lie.lines+length4),1])
lie[1:(lie.lines+length4),6]<-match(lie[1:(lie.lines+length4),6], lie[1:(lie.lines+length4),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length4),1]<-match(lie[1:(lie.lines+length4),1], lie[1:(lie.lines+length4),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start4<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length4,c(4,6)]<-c(1, 0)
stop4<-as.numeric(lie[lie.lines+length4, 1])
# Fill RAC file for the 4-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length4), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length4), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length4), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother3, 4, dbase, 1, cumulDist)
currentMother4<-rac[r,1]
lie.lines<-lie.lines+length4
#----------------------------------------------------------------------------------------------------
if ("root" %in% names(r4)){
for (r5 in r4){# For each 5-order root
if ("geometry" %in% names(r5)){
r<-r+1
ns <- r5$geometry$polyline
length5<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r5)){
age<-r5$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r5)){
age<-r5$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length5),1:6]<-c((lie.lines+1):(lie.lines+length5), rep(1, length5), rep(0, length5), rep(0, length5), lie.lines:(lie.lines+length5-1), (lie.lines+2):(lie.lines+length5+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length5),1:6]<-c((lie.lines+1):(lie.lines+length5), time, rep(0, length5), rep(0, length5), lie.lines:(lie.lines+length5-1), (lie.lines+2):(lie.lines+length5+1))}
lie[(lie.lines+1):(lie.lines+length5),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length5),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length5),13]<-5
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length5),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length5),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length5),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length5),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length5), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length5), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start4:stop4, 7]==lie[lie.lines+1, 7] & lie[start4:stop4, 8]==lie[lie.lines+1, 8] & lie[start4:stop4, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start4+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start4+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start4:(stop4-1), 7]-lie[lie.lines+1, 7])^2+(lie[start4:(stop4-1), 8]-lie[lie.lines+1, 8])^2+(lie[start4:(stop4-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start4:stop4, 7])*(lie[start4:(stop4-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start4:stop4, 8])*(lie[start4:(stop4-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start4:stop4, 9])*(lie[start4:(stop4-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start4:stop4, 7])^2+diff(lie[start4:stop4, 8])^2+diff(lie[start4:stop4, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start4+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start4+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length5), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length5), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start4:stop4, 7])*t+lie[start4:(stop4-1), 7]
yn<-diff(lie[start4:stop4, 8])*t+lie[start4:(stop4-1), 8]
zn<-diff(lie[start4:stop4, 9])*t+lie[start4:(stop4-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
index<-as.numeric(match(min(dist1), dist1)[1])
dist1<-min(dist1)}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start4+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start4+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length5), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length5), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length5), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length5), 10]
dbase<-lie[start4+index-1, 10]+distance3D(x1=lie[start4+index-1, 7], x2=xn[index], y1=lie[start4+index-1, 8], y2=yn[index], z1=lie[start4+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length5+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start4+index-1, 11], NA, 4)
lie[1:(lie.lines+length5+1),]<-lie[order(lie[1:(lie.lines+length5+1), 11], lie[1:(lie.lines+length5+1), 10]),]
length4<-length4+1
lie.lines<-lie.lines+1
stop4<-stop4+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length5),5]<-match(lie[1:(lie.lines+length5),5], lie[1:(lie.lines+length5),1])
lie[1:(lie.lines+length5),6]<-match(lie[1:(lie.lines+length5),6], lie[1:(lie.lines+length5),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length5),1]<-match(lie[1:(lie.lines+length5),1], lie[1:(lie.lines+length5),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start5<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length5,c(4,6)]<-c(1, 0)
stop5<-as.numeric(lie[lie.lines+length5, 1])
# Fill RAC file for the 5-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length5), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length5), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length5), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother4, 5, dbase, 1, cumulDist)
currentMother5<-rac[r,1]
lie.lines<-lie.lines+length5
#------------------------------------------------------------------------------------
if ("root" %in% names(r5)){
for (r6 in r5){# For each 6-order root
if ("geometry" %in% names(r6)){
r<-r+1
ns <- r6$geometry$polyline
length6<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r6)){
age<-r6$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r6)){
age<-r6$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length6),1:6]<-c((lie.lines+1):(lie.lines+length6), rep(1, length6), rep(0, length6), rep(0, length6), lie.lines:(lie.lines+length6-1), (lie.lines+2):(lie.lines+length6+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length6),1:6]<-c((lie.lines+1):(lie.lines+length6), time, rep(0, length6), rep(0, length6), lie.lines:(lie.lines+length6-1), (lie.lines+2):(lie.lines+length6+1))}
lie[(lie.lines+1):(lie.lines+length6),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length6),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length6),13]<-6
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length6),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length6),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length6),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length6),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length6), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length6), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start5:stop5, 7]==lie[lie.lines+1, 7] & lie[start5:stop5, 8]==lie[lie.lines+1, 8] & lie[start5:stop5, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start5+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start5+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start5:(stop5-1), 7]-lie[lie.lines+1, 7])^2+(lie[start5:(stop5-1), 8]-lie[lie.lines+1, 8])^2+(lie[start5:(stop5-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start5:stop5, 7])*(lie[start5:(stop5-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start5:stop5, 8])*(lie[start5:(stop5-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start5:stop5, 9])*(lie[start5:(stop5-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start5:stop5, 7])^2+diff(lie[start5:stop5, 8])^2+diff(lie[start5:stop5, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start5+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start5+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length6), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length6), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start5:stop5, 7])*t+lie[start5:(stop5-1), 7]
yn<-diff(lie[start5:stop5, 8])*t+lie[start5:(stop5-1), 8]
zn<-diff(lie[start5:stop5, 9])*t+lie[start5:(stop5-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
index<-as.numeric(match(min(dist1), dist1)[1])
dist1<-min(dist1)}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start5+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start5+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length6), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length6), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length6), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length6), 10]
dbase<-lie[start5+index-1, 10]+distance3D(x1=lie[start5+index-1, 7], x2=xn[index], y1=lie[start5+index-1, 8], y2=yn[index], z1=lie[start5+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length6+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start5+index-1, 11], NA, 5)
lie[1:(lie.lines+length6+1),]<-lie[order(lie[1:(lie.lines+length6+1), 11], lie[1:(lie.lines+length6+1), 10]),]
length5<-length5+1
lie.lines<-lie.lines+1
stop5<-stop5+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length6),5]<-match(lie[1:(lie.lines+length6),5], lie[1:(lie.lines+length6),1])
lie[1:(lie.lines+length6),6]<-match(lie[1:(lie.lines+length6),6], lie[1:(lie.lines+length6),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length6),1]<-match(lie[1:(lie.lines+length6),1], lie[1:(lie.lines+length6),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start6<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length6,c(4,6)]<-c(1, 0)
stop6<-as.numeric(lie[lie.lines+length6, 1])
# Fill RAC file for the 6-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length6), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length6), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length6), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother5, 6, dbase, 1, cumulDist)
lie.lines<-lie.lines+length6
#----------------------------------------------
if ("root" %in% names(r6)){
for (r7 in r6){# For each 7-order root
if ("geometry" %in% names(r7)){
r<-r+1
ns <- r7$geometry$polyline
length7<-length(ns)
age=NULL
diameter=NULL
if (is.character(final.date)==TRUE & "functions" %in% names(r7)){
age<-r7$functions
for (i in 1:length(age)){
if (TRUE %in% (final.date %in% age[[i]]$.attrs)) {
time<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)
if (length(time)>1 & time[1]<time[2]){time[1]<-time[2]}}}}
if ("functions" %in% names(r7)){
age<-r7$functions
for (i in 1:length(age)){if ("diameter" %in% age[[i]]$.attrs|"rootDiameter" %in% age[[i]]$.attrs) {diameter<-sapply(age[[i]][1:(length(age[[i]])-1)], xnodes)}}}
#c(Num, Date, Bran, Apic, Prec, Suiv)
if (timeserie==FALSE) {lie[(lie.lines+1):(lie.lines+length7),1:6]<-c((lie.lines+1):(lie.lines+length7), rep(1, length7), rep(0, length7), rep(0, length7), lie.lines:(lie.lines+length7-1), (lie.lines+2):(lie.lines+length7+1))}
if (timeserie==TRUE) {lie[(lie.lines+1):(lie.lines+length7),1:6]<-c((lie.lines+1):(lie.lines+length7), time, rep(0, length7), rep(0, length7), lie.lines:(lie.lines+length7-1), (lie.lines+2):(lie.lines+length7+1))}
lie[(lie.lines+1):(lie.lines+length7),11]<-r
if(is.null(diameter)==TRUE) {} else {lie[(lie.lines+1):(lie.lines+length7),12]<-diameter}
lie[(lie.lines+1):(lie.lines+length7),13]<-7 #order
#c(X,Y,Z)
lie[(lie.lines+1):(lie.lines+length7),7]<-sapply(ns, xnodes)
lie[(lie.lines+1):(lie.lines+length7),8]<-sapply(ns, ynodes)
if (length(ns[[1]])==3) {lie[(lie.lines+1):(lie.lines+length7),9]<-sapply(ns, znodes)} else {lie[(lie.lines+1):(lie.lines+length7),9]<-0}
#c(dist)
lie[(lie.lines+1):(lie.lines+length7), 10]<-c(0, cumsum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length7), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 9]))^2)))
# Search the closest point on the mother root (calculate DBase)
parentnode<-which(lie[start6:stop6, 7]==lie[lie.lines+1, 7] & lie[start6:stop6, 8]==lie[lie.lines+1, 8] & lie[start6:stop6, 9]==lie[lie.lines+1, 9])
if (length(parentnode)>0){
dbase<-lie[start6+parentnode[1]-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start6+parentnode[1]-1, 1]
dist1<-0}}
else { #If no physical connection between parent and daughter root. Interpolate new point or find closest point.
alldist<-sqrt((lie[start6:(stop6-1), 7]-lie[lie.lines+1, 7])^2+(lie[start6:(stop6-1), 8]-lie[lie.lines+1, 8])^2+(lie[start6:(stop6-1), 9]-lie[lie.lines+1, 9])^2)
scalx<-diff(lie[start6:stop6, 7])*(lie[start6:(stop6-1), 7]-lie[lie.lines+1, 7])
scaly<-diff(lie[start6:stop6, 8])*(lie[start6:(stop6-1), 8]-lie[lie.lines+1, 8])
scalz<-diff(lie[start6:stop6, 9])*(lie[start6:(stop6-1), 9]-lie[lie.lines+1, 9])
d2<-diff(lie[start6:stop6, 7])^2+diff(lie[start6:stop6, 8])^2+diff(lie[start6:stop6, 9])^2
t<-(-(scalx+scaly+scalz)/d2)
if (length(which(t>=0 & t<=1))==0){
index<-which(alldist==min(alldist))[1]
dbase<-lie[start6+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start6+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length7), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length7), 10]}}
else {
t[t<0]<-NA
t[t>1]<-NA
xn<-diff(lie[start6:stop6, 7])*t+lie[start6:(stop6-1), 7]
yn<-diff(lie[start6:stop6, 8])*t+lie[start6:(stop6-1), 8]
zn<-diff(lie[start6:stop6, 9])*t+lie[start6:(stop6-1), 9]
dist1<-sqrt((xn-lie[lie.lines+1, 7])^2+(yn-lie[lie.lines+1, 8])^2+(zn-lie[lie.lines+1, 9])^2)
if (sum(is.na(dist1)==T)>0) {
index<-as.numeric(match(min(dist1, na.rm=T), dist1)[1])
dist1<-min(dist1, na.rm=T)}
else {
index<-as.numeric(match(min(dist1), dist1)[1])
dist1<-min(dist1)}
if (dist1>min(alldist)){
index<-which(alldist==min(alldist))[1]
dist1<-min(alldist)
dbase<-lie[start6+index-1, 10]
if (connect==TRUE){
lie[lie.lines+1,5]<-lie[start6+index-1, 1] #Update Prec
dist1<-min(alldist)
lie[(lie.lines+1):(lie.lines+length7), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length7), 10]}}
else{
lie[(lie.lines+1):(lie.lines+length7), 10]<-dist1+lie[(lie.lines+1):(lie.lines+length7), 10]
dbase<-lie[start6+index-1, 10]+distance3D(x1=lie[start6+index-1, 7], x2=xn[index], y1=lie[start6+index-1, 8], y2=yn[index], z1=lie[start6+index-1, 9], z2=zn[index])
if (connect==TRUE){
lie[lie.lines+length7+1,1:13]<-c(NA, NA, 0, 0, NA, NA, xn[index], yn[index], zn[index], dbase, lie[start6+index-1, 11], NA, 6)
lie[1:(lie.lines+length7+1),]<-lie[order(lie[1:(lie.lines+length7+1), 11], lie[1:(lie.lines+length7+1), 10]),]
length6<-length6+1
lie.lines<-lie.lines+1
stop6<-stop6+1
pos<-match(NA, lie[1:(lie.lines),1])
lie[1:(lie.lines+length7),5]<-match(lie[1:(lie.lines+length7),5], lie[1:(lie.lines+length7),1])
lie[1:(lie.lines+length7),6]<-match(lie[1:(lie.lines+length7),6], lie[1:(lie.lines+length7),1])
lie[pos, 5]<-pos-1
lie[pos, 6]<-pos+1
lie[pos, 2]<-lie[lie[pos, 6], 2]
lie[pos+1, 5]<-pos
lie[pos-1, 6]<-pos
lie[which(is.na(lie[1:(lie.lines),5])==TRUE), 5]<-0
lie[which(is.na(lie[1:(lie.lines),6])==TRUE), 6]<-0
lie[1:(lie.lines+length7),1]<-match(lie[1:(lie.lines+length7),1], lie[1:(lie.lines+length7),1])
lie[lie.lines+1,5]<-pos
#Calculate diameter for interpolated point (linear function)
if (is.null(diameter)==FALSE){
prec<-lie[pos, 5]
suiv<-lie[pos, 6]
slope<-(lie[suiv, 12]-lie[prec, 12])/distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[suiv, 7], y2=lie[suiv, 8], z2=lie[suiv, 9])
intercept<-lie[prec, 12]
lie[pos, 12]<-intercept+slope*distance3D(x1=lie[prec, 7], y1=lie[prec, 8], z1=lie[prec, 9], x2=lie[pos, 7], y2=lie[pos, 8], z2=lie[pos, 9])}}}}}
if (connect==FALSE) {lie[lie.lines+1,5]<-lie.lines+1}
lie[lie.lines+1,3]<-1
start7<-as.numeric(lie[lie.lines+1, 1])
# Change Suiv and Apic values for the last point of a lateral root
lie[lie.lines+length7,c(4,6)]<-c(1, 0)
stop7<-as.numeric(lie[lie.lines+length7, 1])
# Fill RAC file for the 7-order root
if (connect==TRUE){cumulDist<-dist1 + sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length7), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 9]))^2))}
else {cumulDist<-sum(sqrt((diff(lie[(lie.lines+1):(lie.lines+length7), 7]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 8]))^2+(diff(lie[(lie.lines+1):(lie.lines+length7), 9]))^2))}
rac[r, 1:6]<-c(max(rac[,1], na.rm=TRUE)+1, currentMother5, 7, dbase, 1, cumulDist)
lie.lines<-lie.lines+length7
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
}
#Remove lines with NA values in lie
index<-which(is.na(lie[,1])==TRUE)
if (length(index)>0) {lie<-lie[-index,]}
#Remove lines with NA values in rac
index<-which(is.na(rac[,1])==TRUE)
if (length(index)>0) {
rac<-rac[-index,]
message(paste("Roots with a branching order greater than 7 have been skipped in ", sub(basename(rsml.path), pattern=".rsml", replacement=""), sep=""))}
#Convert matrices to data frames
lie<-as.data.frame(lie)
rac<-as.data.frame(rac)
colnames(lie)<-c("Num","Date", "Bran", "Apic", "Prec", "Suiv", "X", "Y", "Z", "dist", "root", "diameter", "ord")
colnames(rac)<-c("Root", "Mother", "Ord", "DBase", "DApp", "Lengths1")
lie$Bran[lie$Bran==0]<-"false"
lie$Bran[lie$Bran==1]<-"true"
lie$Apic[lie$Apic==0]<-"false"
lie$Apic[lie$Apic==1]<-"true"
if (sum(lie$Z)==0){lie<-lie[,-9]}
#Round dates to integers
lie[,2]<-ceiling(lie[,2])
#Create TPS file
dates<-as.numeric(unique(lie[,2]))
dates<-sort(dates)
if (is.null(unittime1)==TRUE){unittime1<-"unittime"}
if (timeserie==FALSE) {
if (is.null(final.date)==TRUE) {tps<-data.frame(Num=1, Date=1)}
if (is.null(final.date)==FALSE & is.numeric(final.date)==TRUE) {tps<-data.frame(Num=1, Date=final.date)}
if (is.null(final.date)==FALSE & is.character(final.date)==TRUE) {tps<-data.frame(Num=1, Date=1)}}
if (timeserie==TRUE) {tps<-data.frame(Num=c(1:length(dates)), Date=dates)}
#Replace Date by Num in LIE file
date.lie<-lie$Date
lie$Date<-match(lie$Date, tps$Date)
lie$Date[which(is.na(lie$Date)==TRUE)]<-0
#Make RAC file for time series
rac$Root<-rac$Root-1
rac$Mother[rac$Mother!=-1]<-rac$Mother[rac$Mother!=-1]-1
if (timeserie==TRUE){
rac<-rac[,-6] #Delete Lengths column because must be replaced by a matrix
cols<-nrow(tps)
rows<-nrow(rac)
length1<-matrix(0, nrow=rows, ncol=cols)
maxlength<-aggregate(lie$dist, by=list(lie$root, lie$Date), max)
length1[as.matrix(maxlength[,1:2])]<-maxlength$x
length1<-t(apply(length1, 1, function(x){for (i in 2:length(x)){if (x[i]<x[i-1]){x[i]<-x[i-1]}};return(x)}))
colnames(length1)<-paste(rep("Lengths", cols), c(1:cols), sep="")
rac<-cbind(rac, as.data.frame(length1))}
#Export RAC, LIE, TPS
lie.all[[n]]<-lie
rac.all[[n]]<-rac
tps.all[[n]]<-tps
}
result<-list(resolution=resolution, length=unitlength1, diameter=unitdiameter1, time=unittime1, lie=lie.all, rac=rac.all, tps=tps.all)
return(result)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.