3. Merging fluvial clusters and creating lakes and seas

This is trickiest part of the procedure and it has to be done with care. But no worries, the concept is very clear. I tried many options before accepting that this was the easiest and more effective way. However, if you find an even easiest one, I will be very happy to learn! This is how I proceed:

1. I define a fluvial cluster for every end-point. As you can see in the fluvial network above, it is composed by unconnected sub-networks (the big one that includes every river that ends at the ocean, and the two that end at each end-point). Identify them using standard graph functions and libraries. I am using the next trick: I build the transpose of the adjacency matrix of the fluvial network. It is a NxN sparse matrix that I will call A. I then create a vector v with N components in such a way that all components are zero excepting the one at the position of the end-point (remember, his ID) which is 1. With the dot product A.v I got a new vector which components are zero excepting those at the first neighbors of the end-point (according to the reverse fluvial network). Thus, if I iterate with v = v + A.v I reach a limit where v doesn't change anymore and it's 1 for those nodes that belong to the cluster. You see what I'm doing here: I'm going upstream until reaching every point belonging to the river basin.

2. I include as first cluster that of the ocean. I just put together all these nodes with all the basins that end directly at the ocean. It can be done using the same trick as before, just starting in v with those points at the sea. So, at this point, we have a list of clusters, starting with the one of the ocean, with a length of number of end-points plus one. In our example, we have two end-points, that is, three clusters. The image below (1) show the three clusters: the big thing is the one of the ocean, and the other two those in each caldera.