riverflow=scan("dischargesynt.txt",skip=1)
#upload of observed river flow
rainfall=scan("rain-evaposynt.txt",skip=1)
#upload of observed rainfall
rainfall=array(rainfall,dim=c(2,8760))
rainfall=rainfall[1,]
#eliminating the column of evapotranspiration

#Definition of vectors and input data
riverflowLR=rep(0,8760)
#tentative value for k
k=300*3600
areab=1214000000
tstep=3600
ndeltat=8760
w=0

#Start of the computational cycle
for(i in 1:ndeltat)
{
  w=w+rainfall[i]/1000/tstep*tstep
  riverflowLR[i]=w/k
  w=w-riverflowLR[i]*tstep
  if(w<0) w=0
}
riverflowLR=riverflowLR*areab
#rescaling to the catchment area

plot(riverflowLR,type="l")
lines(riverflow,col="red",lwd=2)
#Plot of the output