Hi! I have R script (found in an open source) scripted for raster files. I tried to use it in R studio and think the code works but did not work some part of the output raster as I opened and saw the output raster in ArcGIS. I have not used R studio before and therefore do not pose knowledge about it. All I know that the script would solve my problem. What I need is to calculate sum of cell of a raster file along flow path from each cell. Is there anyone kind enough to help me out to get my desire output? I am attaching a a jpg of my desire solution of the problem and coding below. I have two input files (r1 & rdir) needed for running this code.However could not attach the two files here because GTiff file seems to be not allowed to attach here.If anyone need I can provide them.
library(raster)
library(sp)
library(rgdal)
Charge your raster
r1<-raster( ) #include the complete path of your raster file inside the paenthesis
rdir<-raster( ) # Add flow direction raster
r1.matrix<-as.matrix(r1)
mdir<-as.matrix(rdir)
r2.matrix<-r1.matrix
for (i in 1:ncol(r1.matrix)){
for (j in 1:nrow(r1.matrix)){
x<-i
y<-j
valor<-0
while ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if (mdir[y,x]==1){
valor<-valor+r1.matrix[y,x]
x<-x+1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==16){break}
}
next
}
if (mdir[y,x]==2){
valor<-valor+r1.matrix[y,x]
x<-x+1
y<-y+1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==32){break}
}
next
}
if (mdir[y,x]==4){
valor<-valor+r1.matrix[y,x]
y<-y+1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==64){break}
}
next
}
if (mdir[y,x]==8){
valor<-valor+r1.matrix[y,x]
x<-x-1
y<-y+1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==128){break}
}
next
}
if (mdir[y,x]==16){
valor<-valor+r1.matrix[y,x]
x<-x-1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==1){break}
}
next
}
if (mdir[y,x]==32){
valor<-valor+r1.matrix[y,x]
x<-x-1
y<-y-1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==2){break}
}
next
}
if (mdir[y,x]==64){
valor<-valor+r1.matrix[y,x]
y<-y-1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==4){break}
}
next
}
if (mdir[y,x]==128){
valor<-valor+r1.matrix[y,x]
x<-x+1
y<-y-1
if ((x>=1) & (x<=ncol(r1.matrix)) & (y>=1) & (y<=nrow(r1.matrix))){
if (is.na(mdir[y,x])==TRUE){
break
}
if(mdir[y,x]==8){break}
}
next
}
}
r2.matrix[j,i]<-valor
}
}
r2<-r1
r2<-r2.matrix
writeRaster(r2, filename="r2.tif", format="GTiff", overwrite=TRUE) #change the path with the name of your raster including the .asc (if you use the ascii format)