How can I get weighted Flow Accumulation to deal with negative values? [duplicate]

Question

here's one I haven't seen anywhere, to my surprise. I am trying to do something a little different with Flow Accumulation, which I would think has been attempted before, but it's stumping me for now.

I want to weight flow accumulation by a loading map of pollutant z to get its accumulation as it's transported downstream. This is simple enough, except I'd like some areas to be net sinks -- that is, where a fixed amount of the pollutant is removed rather than accumulated. So, say our flow path contains cells with values 1 1 -1 1. -1 is a net sink, and so our accumulation here, with flow going left to right, would be 1 2 1 2. This too seems to be doable. I can make one weighting raster of 1 1 0 1 and another of 0 0 1 0, run the flow accumulation with both, then subtract the 2nd from the 1st. OK.

However, imagine this flowpath scenario -- 1 1 -5 1. Here I have a powerful "sink" at -5, which is in a position to remove everything that has flowed to it. Accumulation at that cell would be 2, and with 5 removed, this leaves us with -3 using the above method. But of course, having a quantity of -3 is nonsensical, and instead it should simply drop to 0, as everything is removed. At the end of this little 4-cell flowpath, accumulation should thus be 1.

Is there any way to weight or otherwise get Flow Accumulation to intelligently do the opposite of accumulating (de ... cumulating?) in certain specific areas, while also respecting a lower bound of 0? I'm hoping the method, if it exists, will be implementable in Python as well.

UPDATE:
If anyone is interested in this problem in the future, I may have figured out how to do it (though am still open to correction, as I'm not 100% sure). By default, the Flow Accumulation tool treats negative weightings as 0.

The steps are as follows:
1). Run Flow Accumulation with your weighting grid that includes negatives.
2). Apply the following conditional statement to the result:

Con(("FlowWeight" < 0) & ("FlowAccum" > 0),Con(-"FlowWeight" < "FlowAccum",
-"FlowWeight","FlowAccum"),0)

This assigns a non-zero value to all cells with a negative weight and a non-zero accumulation, such that if the weight represents uptake of some but not all of what has accumulated there, the cell gets the value of that amount, as a positive number, and if uptake is more than what has accumulated, the value is simply the accumulated amount.

3). Run Flow Accumulation weighted by the grid resulting from the above. This is like the "loss accumulation".

4). Subtract the result of step 3 from the result of step 1. Force any negative values to 0.

Answer 1

I do not work with ArcGIS, but I use SAGA mainly for hydrologic modelling. What I can suggest here is: - fill the drainage depressions, and then subtract the filled DEM from the original DEM, getting a raster with the depression depth; - compute flow accumulation, which will flow over the filled depressions; - compute an accumulation function, depending on the depression depth, to get a raster of accumulated polluant; - use the accumulation polluant raster in the modelling process (in SAGA there is an option for both material and weight in such a modelling).

Answer 2

I dealt with a similar problem here. Like you speculated, you have to split your weight raster into positive and negative values, flow accumulate them separately, and recombine them. Here's a python code snippet from my answer that will do it for you:

def angleAccum(flowdir,angle):
    pos_angle = Con(angle,angle,0,"Value > 0")
    neg_angle = Abs(Con(angle,angle,0,"Value < 0"))
    pos_accum = FlowAccumulation(flowdir, pos_angle,"FLOAT")
    neg_accum = FlowAccumulation(flowdir, neg_angle,"FLOAT")
    return Minus(pos_accum,neg_accum)

It was written with averaging upslope aspect in mind (hence theangle argument), but it will work on any raster with negative values.

To use it, open your python window and enter:

from arcpy.sa import *

then copy and paste the function into the window and hit enter. Finally, call the angleAccum function passing in your flow direction raster and your weight raster like so:

weightAccum = angleAccum(inFlowDir,inWeightRaster)
weightAccum.save("C:/Path/to/outputraster")

Another option that should theoretically work is to transform all of your numbers to be positive numbers, then subtract the unweighted accumulation from your weighted accumulation. However, in practice this didn't seem to work quite right (See comments in linked question).