prot1=read.table("outputcharmm/calc_bindingEnergy_alongMD_aspp_prot_1_anion+ini_8.out")

#**********************************************
#Separate les energies pr3-mb et pr3-eau

mb.prot1=prot1[,1][seq(1,1999,2)]
eau.prot1=prot1[,1][seq(2,2000,2)]
mb.prot1.Et=prot1[,2][seq(1,1999,2)]
eau.prot1.Et=prot1[,2][seq(2,2000,2)]

#***********************************************
#Calcul la moyenne et l'erreur standard de l'energie effective et de l'energie totale

stderr <- function(x) sqrt(var(x)/length(x))

mean(mb.prot1)
stderr(mb.prot1)

mean(mb.prot1.Et)
stderr(mb.prot1.Et)

#*******************************************************
#Calcul l'energie de liason et sa moyenne sur les 500 dernieres ps

WWimm.prot1=mb.prot1-eau.prot1
mean(WWimm.prot1[500:1000])
stderr(WWimm.prot1[500:1000])