The Grand Canyon DEM colored by topographic gradient (based on SRTM-3 downloaded from

What is TopoToolbox?

TopoToolbox provides a set of MATLAB functions that support the analysis of relief and flow pathways in digital elevation models. The major aim of TopoToolbox is to offer helpful analytical GIS utilities in a non-GIS environment in order to support the simultaneous application of GIS-specific and other quantitative methods.

Who should use TopoToolbox?

TopoToolbox is not a click & go software. You should be familiar with MATLAB and its programming syntax. If you feel well versed in Matlab, you will like TopoToolbox for its versatility, functionality and ability to be customized.

How do I install TopoToolbox?

The folder with TopoToolbox as well as its subfolders (e.g. utilities, GIStools) must be on the MATLAB search path. Unzip the TopoToolbox file and add the folder and its subfolders to the search path with the function addpath (see readme file) or the Set Path dialog box.

Where can I find a documentation?

“What makes code bad?” “No comment!” Your job as a programmer is only half done if there is no documentation. We have done our job and there is a documentation that you’ll find by typing

doc topotoolbox

in the command window. This will open TopoToolbox help pages in the MATLAB’s documentation browser. This help features a Getting Started with TopoToolbox section, a user guide, and a function reference.

In addition, you will find a help section in each function including syntax, input and output arguments, and an example. For example, to access the help of the function gradient8, just type

help gradient8

in the command line.

This blog is intended as an additional resource of more comprehensive examples. An overview on the topics covered is here.

In addition, you might want to take a look at two papers about TopoToolbox published in Environmental Modelling and Software and Earth Surface Dynamics.

Where can I report bugs or suggestions for enhancement?

You found a bug? Sorry for the inconvenience. Please tell me more about it here. The same holds for enhancement requests. I cannot guarantee that I’ll have the time to deliver, but I’ll do my best to keep TopoToolbox up-to-date.

How can I stay informed about updates?

There is no mailing list. All news about TopoToolbox will be posted in this blog. In addition, you may want to follow the software on GitHub.


109 thoughts on “TopoToolbox

    Jean-Marie AUGUSTIN said:
    May 3, 2015 at 12:41 pm

    I have just discovered you toolbox which seems very interesting. I am the developer of a software suite written in Matlab language (SonarScope :
    I would be very interested to propose DEM processings using TopoToolbox in SonarScope.
    In addition to SonarScope I developped a 3DViewer based on Nasa World Wind that could interest you. If you are interested in getting copies of these software suites please tell me,

    Jean-Marie Augustin, Ifremer, France

    Martine Collaud Coen said:
    August 25, 2016 at 11:28 am

    I do use TopoToolbox since several months and have a Problem with the coor2ind function. There is almost always a shift between the IX linear index calculated by coor2ind leading to an error between 0 and 50 on IX. If I plot on the topography the x,y coordinate (latitude and longitude), the shift can be visualised. Up to now I always plot them and caculate manually the shift. I have now to make calculation for a whole grid. Do you have a solution ?

    wschwanghart responded:
    August 25, 2016 at 1:56 pm

    Hi Martine, I do not fully understand the error. coord2ind snaps x and y coordinates to the center of the nearest cell, and returns this cell’s linear index. Thus, there will usually always a shift between the x and y coordinates that you supply to coord2ind, and those coordinates that you obtain from ind2coord.

    Without knowing exactly your application, would the function GRIDobj/interp do the job? This function interpolates values in a grid to your coordinates.

    Best, Wolfgang

    Martine Collaud Coen said:
    August 25, 2016 at 2:28 pm

    Hi Wolfgang,

    I do not know if I’m doing something wrong. But coord2ind does not give me the nearest cell, but sometimes some 10 to 30 cells before or after the real cell. So that I implemented my Points with latitude, longitude and with the number of cells that have to be shifted to have the real Point. For example, using
    IX=coord2ind(X,Y, lon_st, lat_st),
    hold on;

    then I Count on the map the difference d and if I plot


    then xx2 and yy2 are on the right place.

    d can be up to 45 sometimes. There is gridcells where this shift is small, then it increases and becomes suddently negatives and then decrease again. I asked my self is this can be due to a wrong rounding of the coordinates.
    I hope that you understand my Problem and many thanks for your answer
    Best regards

    wschwanghart responded:
    August 26, 2016 at 7:34 am

    Hi Martine,
    quite difficult to judge from distance. But following possible errors come to my mind:
    – Do the GRIDobj and lon_st, lat_st share the same coordinates?
    – isequal(DEM.size,size(DEM.Z))
    If that doesn’t help, please write me an email (

    Best regards, Wolfgang

    Francisco Peña said:
    November 4, 2016 at 5:05 pm

    Hi Wolfgang,
    I’m very interested in using your software as part of my PhD research.
    I’ve recently downloaded Topotoolbox and I have experience some issues using basic commands during the 1st tutorial such as crop, resample or imageschs. For instance, after writing the command imageschs(DEM,min(gradient8(DEM),1)) I received the following error:
    ??? Undefined function or method ‘narginchk’ for input arguments of type ‘double’.

    Error in ==> GRIDobj.min at 5
    narginchk(1, 2)

    The paths have been correctly assigned so basically I don’t understand what I’m doing wrong.
    Could you please provide some advice in order to solve the problem?
    Thank you in advance,

    floodriskgis said:
    November 4, 2016 at 5:26 pm

    Dear Wolfgang,
    I’m very interested in using your software for my PhD research project.
    I recently downloaded Topotoolbox and unfortunately I’ve experienced some issues during the 1st tutorial.
    This happens when using the command >> imageschs(DEM,min(gradient8(DEM),1))
    ??? Undefined function or method ‘narginchk’ for input arguments of type ‘double’.

    Error in ==> GRIDobj.min at 5
    narginchk(1, 2)
    The paths were added correctly and I’ve followed the user guide step-by-step.
    Could you please provide some advice to solve this issue?
    Thanks in advance,

    Francisco Peña

    wschwanghart responded:
    November 4, 2016 at 6:36 pm

    Dear Francisco,

    I guess that these issues arise due to an outdated version of matlab. narginchk was introduced in 2011b ( It seems your version is older than that.

    Best, Wolfgang

    ahmed said:
    November 26, 2016 at 12:01 pm

    Hi There, I have just start to work on TopoToolbox that i worked on TecDEM and stream profiler. Kindly I would like to ask if there someone in the UK working on TopoToolbox to contact him and see how to work. Another question is there a window appears to to choose the function(similler to TecDEM) or i have to write a command for each process in the MATLAB window.

    Many thanks

    VBuford said:
    December 4, 2016 at 3:43 am

    Hi there,

    I am currently using topotools to run ksn analyses. My tif is from Peru, so both negative lat and long. When I run the ksn help matlab code, and I get to L86 (the STATS section), I get this error:
    Error using histc
    Edge vector must be monotonically non-decreasing.

    Error in slopearea (line 141)
    [~,ix] = histc(a,edges);
    Do you have any idea why this is ocurring or how to fix this?

      wschwanghart responded:
      December 4, 2016 at 9:30 pm

      Hi, I had this issue before. It happened because a user worked with a geographic coordinate system (latitude, longitude). Please work with projected coordinate systems (e.g. UTM) with the horizontal and vertical coordinates having the same unit, e.g. [m]. You can use the function

      DEM = reproject2utm(DEM);

      Hope this works out. Please let me know if it doesn’t.
      Best regards, Wolfgang

        Ramendra Sahoo said:
        March 28, 2017 at 4:51 am

        Dear Wolfgang, I have recently started using topotoolbox for my research, which mainly focuses on quantitative attributes of drainage network. I was going through the userguide_XX.m files for acquaintance with the available functions. One of the example in that is the derivation of flow direction (both single and multiple) from a dem(baranja_hill.mat). Upon executing the line 28, I got an error message. Following are the line and the associated error message.

        load baranja_hill
        demf = fillsinks(dem);
        Undefined function ‘fillsinks’ for input arguments of type ‘double’.

        Later I checked and saw that ‘fillsinks’ use a GRIDobj as an input. Also, in case of the ‘flowdir’ function, I tried both the ‘GRIDobj’ type argument and ‘double’ type argument. For both of them, it returns the same type of error.

        Kindly look into this. I am using MATLAB R2015a, and I have licenses for both image processing and mapping toolbox. Please help me with the proble.


    wschwanghart responded:
    March 28, 2017 at 7:37 am

    Dear Ramendra, many thanks for your comment. You are right. The baranja_hill.mat contained the matrices X,Y and dem as they were used by TopoToolbox before we introduced GRIDobj. I have now updated the file which now contains the GRIDobj. At the same time, I realized that userguide_2_flats.m is completely outdated and incompatible with version 2 of TopoToolbox. At has been a while since I have looked at these files. Thanks for these hints!
    Best regards, Wolfgang

      Ramendra Sahoo said:
      March 28, 2017 at 8:25 am

      Dear Wolfgang, thank you for addressing my concern. I have downloaded the latest work-in-progress copy of topotoolbox only. But as per your suggestion, I again downloaded it and checked for the banjara_hills file. The same problem persists and there is no GRIDobj type file. And can you also kindly check the ‘flowdir’ tool, as it is taking neither mat file nor GRIDobj file as input. All the examples concerning the flowdir function has shown a mat file as input. But it is throwing the same error.
      ‘Undefined function ‘fillsinks’ for input arguments of type ‘double’.’


    Ramendra Sahoo said:
    March 28, 2017 at 9:36 am

    Dear Wolfgang, thank you very much. All the confusions seems to have been sorted out for now.

    Thanks and Regards

    Ramendra Sahoo said:
    April 4, 2017 at 10:09 am

    Dear Wolfgang, can you kindly suggest any way to write a 2-d matrix into a geotiff file? I have been trying to do that using the geotiifwrite. But the function doesn’t perform consistently and throws some errors now and then.


      wschwanghart responded:
      April 8, 2017 at 8:06 pm

      Dear Ramendra, without seeing the error message, it is difficult to help out here. Unless you are writing a GRIDobj to a geotiff, please check the help of geotiffwrite ( ). It might also help to take a look at the code of GRIDobj2geotiff.
      Best regards, Wolfgang

    Caterina Samela said:
    May 4, 2017 at 1:54 pm

    Dear Wolfgang, thanks for these super-useful functions! Can you kindly suggest me the proper one to extract the main channel (e.g. the longest one) in a pre-identified river network?

      wschwanghart responded:
      May 4, 2017 at 2:09 pm

      Hi Caterina, thanks for your kind comment! trunk is the function you need.
      Best regard, Wolfgang

    Sotiris Karalis said:
    July 9, 2017 at 10:32 am

    Wolfgang Hi!
    I am going through the TTLEM_usersguide_1_intro, to see what is in there, but when i give

    p = ttlemset

    i get an error saying
    Undefined function ‘addParameter’ for input arguments of type ‘inputParser’.

    Error in ttlemset (line 160)
    addParameter(p,’TimeSpan’,2000000,@(x) isscalar(x) && x>0);

    Error in TTLEM_usersguide_1_intro (line 138)
    p = ttlemset %#ok show parameters
    What is it that i am doing wrong?
    Thank you in advance!

    Sotiris Karalis

      wschwanghart responded:
      July 10, 2017 at 11:19 am

      Hi Sotiris, replace all calls to addParameter by addParamValue. addParamValue was replaced a couple of years ago with addParameter. Hope this works.
      Best regards, Wolfgang

        Sotiris Karalis said:
        July 11, 2017 at 12:06 pm

        Thank you Wolfgang! I followed your instructions and run all 3 simulations in TTLEM, which are very impressive i must add! Fyi i don’t see the box depicting the status of the simulation, but this is not a problem. Is it because i have MATLAB 2012b?

        Now i have a new question regarding TTbox usersguide1 and 2. In particular, when i give

        i get:
        Attempt to reference field of non-structure array.

        Error in STREAMobj/plotdz (line 102)
        colororderindex = mod(ax.ColorOrderIndex, size(ax.ColorOrder,1));
        What is the problem now?

        Warm regars.

    wschwanghart responded:
    July 11, 2017 at 12:15 pm

    Hi Sotiris, sorry for those problems. I guess that you downloaded topotoolbox from here: . Yet, this repository is not as frequently updated as my own repository, where this bug should be removed:

    Cheers, W

      Sotiris Karalis said:
      July 12, 2017 at 7:37 am

      Hi Wolfgang, even though i downloaded from the other site, the same problem remains..
      How can i resolve it?
      Could you please give instructions? This is a very usefull function..

      Warm Regards

        wschwanghart responded:
        July 12, 2017 at 8:35 am

        Hi Sotiris, ah, that seems also a problem related to an older version of MATLAB. MATLAB changed its graphics engine with version 2014b (?). With this change, axes and lines have become objects rather than handles with meaningless numbers. Let me fix this in a minute.

        wschwanghart responded:
        July 12, 2017 at 8:49 am

        Sotiris, plotdz should be able to run on versions before 2014b now. Note, however, that I do not have an old version of MATLAB to test the changes. In addition, other functions may have similar backward compatibility issues that I won’t be able to resolve.
        Cheers, Wolfgang

    Sotiris Karalis said:
    July 12, 2017 at 6:31 pm

    Thanks Wolfgang. Probably it’s time to update..
    I ll get back to you when i ll start a real world application with my data!
    Good job!

    Saleh Alqahtani said:
    August 17, 2017 at 2:13 pm

    Thanks for the great program.
    I have tried to use the topotoolbox functions on a DSM tif file from a drainage basin I’m working on. I am trying to calculate the ksn for the basin streams but the result map shows overlapping ksn values even at very local scale. I used the STREAMobj2mapstruct function to have the ksn values generalised for 10km long segments of the stream but the map stays the same.
    I am not expert in matlab so I appreciate that the mistake might be a trivial one to avoid.

      wschwanghart responded:
      August 18, 2017 at 6:51 am

      Saleh, thanks for your comment. Without seeing the code, I am unfortunately unable to resolve this issue. However, please check the new function STREAMobj/aggregate ( Alternatively, you can also use STREAMobj/smooth with a large value of K (e.g. K=100) to smooth your ksn values along the river network.
      Best, Wolfgang

        Saleh Alqahtani said:
        August 21, 2017 at 10:13 am

        HI Wolfgang. Thanks for your reply and sorry for my late reply. I have tried it again and it worked on DEM data I got from another source but for the same region.

        Saleh Alqahtani said:
        August 24, 2017 at 10:37 am

        Hi dear Wolfgang,
        May I ask if there is a way to extract the ks and theta values for just one stream instead of all of the streams in one drainage catchment? And the same also about plotting the elevation vs distance for just one stream?
        Best, Saleh

    wschwanghart responded:
    August 24, 2017 at 11:03 am

    Hi Saleh, sure, this is possible. You just need to adjust your stream network so that it contains only one stream. You can do this interactively using the modify function:
    Sn = modify(S,’interactive’,’reach’);
    Alternatively, construct a one-stream STREAMobj as follows:
    S = STREAMobj(FD,’channelheads’,IX);
    where IX can be obtained from the channelhead coordinates x, y and the function IX = coord2ind(DEM,x,y).
    Best, Wolfgang

      Saleh Alqahtani said:
      August 25, 2017 at 10:15 am

      Brilliant.. the Sn worked very well.. Many thanks Wolfgang

    Xiaoli LIU said:
    March 6, 2018 at 3:07 am

    Dear Wolfgang,
    Thanks for your great program very much. I have download the latest copy of Topotoolbox from last week. Unfortunately I’ve experienced some issues when I go through ‘topoapp’.
    1)About the button ‘Import outlets from ascii file(x,y)’

    My outlets is a ascii file, it contains a pair of project coordinate value of the outlet. When I input a .txt file, the following error will occur:

    Undefined function ‘max’ for input arguments of type ‘struct’.
    Error in importoutlets (line 34)
    if max(outlets(:,1))<=180 && max(outlets(:,2))<=90 % probably lat,lon

    Error while evaluating PushTool ClickedCallback

    2) About the button ‘Pick channel reach’

    After I pick two points along the streams on the main windows, the following error will occur:

    Undefined function ‘snap2stream’ for input arguments of type ‘GRIDobj’.
    Error in pickreach (line 38)
    [ix1,~] = snap2stream(W,ix1);

    Error while evaluating PushTool ClickedCallback

    I am using MATLAB R2014b, the paths were added correctly and I’ve followed the user guide step-by-step.
    I would deeply appreciate if you can provide me some advices to solve the above issues.
    Thank you for your time and consideration. I look forward to hearing from you at your earliest convenience.

    Xiaoli LIU


      wschwanghart responded:
      March 6, 2018 at 11:57 am

      Dear Xiaoli Liu, thanks for your message and details about those errors. I’ll forward them to Dirk, who developed topoapp. The second error is obviously my fault, because I removed the function GRIDobj/snap2stream.m function. I’ll add it again to warrant that the tool runs properly.
      Cheers, Wolfgang

        Xiaoli LIU said:
        March 7, 2018 at 12:50 am

        Hi Wolfgang.
        Thanks for your reply and Dirk’ contribution. I’ll do it again and look forward to a new copy will bring new experiences. I’ll also feedback results timely.
        Best, Xiaoli

    Vivien Mai Yung Sen said:
    May 28, 2018 at 9:04 am

    Hi Wolfgang,
    I begin to use Topotoolbox for my research work. I was wondering if there was a way to convert shapefiles to streamobject? I calculated streams before on Qgis and I would like to use them on Topotoolbox.
    Thank you very much for your work


      wschwanghart responded:
      May 28, 2018 at 9:15 am

      Hi Vivien,
      unfortunately such a tool does not exist. And this is for a reason: The stream network (STREAMobj) should/must be a subgraph of the flow network (FLOWobj) because otherwise some algorithms will throw errors or return nonsense.
      Why did you derive the stream network in QGIS? Mapping channelheads? If so, then consider deriving the FLOWobj in TopoToolbox and use the mapped channelheads as input into the STREAMobj function

      S = STREAMobj(FD,'channelhead',ix);


      ix = coord2ind(DEM,x,y)

      where x and y are the coordinates of channelheads.
      Does this help?
      Cheers, Wolfgang

        Vivien Mai Yung Sen said:
        May 30, 2018 at 12:41 pm

        Hi Wolfgang,
        Thank you for your answer. I derived the stream network in Qgis because I didn’t know until two weeks Topotoolbox!
        Ever since I have derived the streams in Topotoolbox but the problem is for some low slope areas. Sometimes, the calculated channels don’t match with real ones. So I drew some channels from satellite images and dug the DEM using a raster calculator in order to force the flow. It worked well on Qgis but not in Topotoolbox. The flow ignores the dug pixels and the calculated streams are the same than without digging.
        Cheers, Vivien

    Federico Carballo said:
    May 30, 2018 at 11:51 am

    Hello, Wolfgang,
    Here in Argentina, we’re trying to use the topotoolbox too. We worked with it in a Matlab R2012a with a version call “master”, there was no problem with it.
    We upgrade the Matlab version to a R2016a and the topotoolbox to 2.2 one, and we get the message error:
    About the button ‘Pick channel reach’

    After I pick two points along the streams on the main windows, the following error will occur:

    Undefined function ‘snap2stream’ for input arguments of type ‘GRIDobj’.
    Error in pickreach (line 38)
    [ix1,~] = snap2stream(W,ix1);

    Error while evaluating PushTool ClickedCallback

    Reading the blog, you explain it could be related with the absent of the GRIDobj/snap2stream.m function. So we replace the entire folder of the @GridObj with the one of the master version, and reinstalled the Topotoolbox 2.2. This problem was repeated again. Then, we tried to change the version of the Topotoolbox to the master one and the issue apear again. Finally, we change the Matlab version to the R2015a and the problem continued.

    Following the issue, in the snap2stream.m, line 69, parse(p,varargin{:}); an error message appears, it tries to use the varagin as a variable:
    Undefined variable varargin.
    Could be this related with the issue?
    Well, We hope that this problem has a solution.
    We find the topotoolbox very usefull.



    wschwanghart responded:
    May 30, 2018 at 1:05 pm

    Hi Fede, thanks for your comment. I assume that you refer to errors that occur when using topoapp, right? And you are right. Dirk and I have still not addressed some of the bugs that were introduced by the new graphics system in R2014b and we apologize for the inconveniences that this causes.

    Moreover, I recently removed the GRIDobj/snap2stream function because of redundancy (there is a better function STREAMobj/snap2stream). However, I then realized that the function is used by topoapp and added it back to topotoolbox. With the latest version 2.3pre Version available on my github account (, the error message should not appear anymore.

    However, I can only encourage not to use topoapp anymore but to use the command line. This is much more stable and versatile.

    Cheers, Wolfgang

    DSG said:
    June 11, 2018 at 9:24 pm


    I am new to using TopoToolbox 2 and have been using STREAMobj with some success. However, I am trying to extract the stream networks separately after using STREAMobj/split for processing in other programs. Split works nicely, but I am unable to separate the streams with their x and y locations


      wschwanghart responded:
      June 12, 2018 at 7:00 am

      HI DSG, split may be the wrong function here. You can use STREAMobj2mapstruct to create a mapping structure that you can export as shapefile using the function shapewrite.

    Helbert Garcia said:
    August 20, 2018 at 3:01 pm

    Dear Wolfgang,

    I have been working with the ksn tutorial and it all worked out well. At the moment I want to export the mapping structure using shapewrite I am getting this error message:

    Error using shapewrite>writeSHP (line 94)
    Unable to open ksn.shp for writing.

    Error in shapewrite (line 80)
    [shapeType, boundingBox, index] = writeSHP(S,basename);

    I am a very new user of MATLAB and TopoToolbox so this could be probably some beginner mistake. Thanks for your help.

    Kind regards,


      wschwanghart responded:
      August 20, 2018 at 5:14 pm

      Hi Helbert, I suspect that the file is locked by another program, e.g. ArcGIS.
      Best regards, W.

    Emily Schottenfels said:
    August 22, 2018 at 8:33 pm

    Hello, I use this toolbox for a number of functions (thank you!), one of which is to fill erroneous pits/sinks in a bathymetry DEM before performing basin metrics/analyses. The fillsinks (& depth threshold) method is useful, however I would like to select specific sinks to be filled/left unfilled. I have figured out a method to create a sinks raster in ArcMap (find and select the deepest cell in the sinks) and use that as an input in MATLAB, however there are problems with this and some sinks are still filled.

    example code:
    S = GRIDobj(‘sinks.tif’); %tiff from ArcMap of the sinks to fill (select deepest cells in a sink)
    SINKS = GRIDobj(DEM, ‘logical’);
    DEMfs = fillsinks(DEM, SINKS);

    I was wondering if there is a method to select sinks/create a sink raster all in MATLAB to avoid switching between programs and the issues that arise from that.

    Please let me know if you need clarification or have additional questions.

      wschwanghart responded:
      August 22, 2018 at 8:57 pm

      Hi Emily, now that’s an excellent question that may be interesting to many that work with bathymetry or in areas with endorheic basins. However, I don’t have a straight-forward answer, except that TopoToolbox currently doesn’t have the tools to do this. There used to be such tool (preprocessapp in v 2.1 (, but I removed it as it was rather buggy. Yet, you may still want to try it.

      Yet, the MATLAB image processing toolbox has everything you need. imregionalmin, imhmin or imextendedmin could provide the functionality to build the available tools. bwselect may help in selecting the sinks with mouse clicks.

      Hope this helps.
      Best, Wolfgang

    Mark Jwaideh said:
    December 1, 2018 at 4:11 pm

    Hello Wolfgang,

    I am currently trying to use topotoolbox to do a global weighted accumulated flow analysis. Trying to analyse contaminant loads in rivers. However, when processing the code on a global raster 4320 x 2160 – it seems further downstream of the river no data is given. Stream flow accumulation get cut off further up stream without flowing to the mouth. Any help would be much appreciated.



      wschwanghart responded:
      December 1, 2018 at 4:58 pm

      Hi Mark, there can be several reasons for this:
      – NaNs in the DEM.
      – Topographic depressions and no preprocessing. Use FLOWobj(DEM,’preprocess’,’carve’) if you are using an outdated version of TopoToolbox.
      Best, Wolfgang

    Clément Desormeaux said:
    March 5, 2019 at 4:05 pm

    Hello Wolfgang

    I am currently trying to use topotoolbox to develop chi maps and KSN maps from DEMs (30m and 10m) but I have some problems on DEMs. A solution to my problem could be to select an outlet and recalculate the corresponding watershed as well as drainage networks, mn ratio, chi map etc … Is there a function to give an outlet (from its coordinates) and not the whole DEM to topotoobox?

    Thanks for your help,
    Kind regards,


      wschwanghart responded:
      March 5, 2019 at 4:41 pm

      Hi Clément, not sure how that could/should work. Why not loading the entire DEM?

        Clement Desormeaux said:
        March 5, 2019 at 6:09 pm

        I use DEMs cut from the watersheds (shp format) with a different initial resolution. This creates problems during processing (especially in the calculation of flows direction and stream networks) and despite the many tools provided by topotoolbox (smooth operator, minerea etc…) problems persist.

          wschwanghart responded:
          March 6, 2019 at 8:04 am

          Hi Clement, if you need a solution to the problems of your analysis, you need to be much more specific. Articulating the problem right is the first step to its solution. 😉 Cheers, Wolfgang

    Clément Desormeaux said:
    March 8, 2019 at 9:18 am

    Hi Wolfgang,
    ok sorry, I’ll try to be clearer.
    the thing we try to do really is to calculate the m/n ratio for all of the basins. But with the scripts we found, it seems to be possible only on the biggest basin in the DEM. What we would like to do, is to give a list of outlets, and then make the calculation of the m/n on each corresponding basin.
    Is it possible with topotoolbox?

    Cheers, Clément

      wschwanghart responded:
      March 11, 2019 at 8:16 am

      Dear Clément,
      I often use the function klargestconncomps to apply some of the methods to only the largest drainage basin. The function chiplot that I’d recommend to find the mn ratio, for example, can be applied to one drainage basin only. If you want to find the mn ratio for several drainage basins, I would proceed as follows:

      DEM = GRIDobj('srtm_bigtujunga30m_utm11.tif');
      FD     = FLOWobj(DEM);
      A       = flowacc(FD);
      S      = STREAMobj(FD,'minarea',1000);
      CS    = STREAMobj2cell(S);
      mn   = zeros(size(CS));
      MN   = GRIDobj(DEM);
      for r = 1:numel(CS)
           C = chiplot(CS{r},DEM,A,'plot',false);
           mn(r) =;
           D = drainagebasins(FD,CS{r}) > 0;
           MN.Z(D.Z) = mn(r);

      If you run that script, you will notice that there are numerous high and negative values of mn. Preprocessing the DEM (e.g. using crs) and removing edge effects in the DEM (e.g. using removeedgeeffects) will likely improve the results. Moreover, difficulties will arise if rivers have large knickpoints, i.e. they are transient.
      Good luck, Wolfgang

    Mark Jwaideh said:
    March 11, 2019 at 11:43 am

    Hi Wolfgang,

    I’m wondering if topotoolbox is able to conduct a problem I am facing within its functions. I am trying to calculate the product of all downstream weighted fractions of flow. Therefore the upstream cell should be the product of all downstream cell weightings.

    Flow Accumulation does not bring the desired output as it is the addition of upstream flows. Is there a function that works backwards and calculates the products of fractions of flows?

      wschwanghart responded:
      March 11, 2019 at 12:01 pm

      Hi Mark,
      I had a similar request recently. I came up with the function below that sums up values from downstream.

      function G = downslopesum(FD,G)
      ix = FD.ix;
      ixc = FD.ixc;
      g = G.Z;
      for r = numel(ix):-1:1
          g(ix(r)) = g(ixc(r)) + g(ix(r));
      G.Z = g;

      I am wondering what addition/multiplication from downstream could be good for? I might add this as a function if I can think of some use cases.
      Cheers, Wolfgang

        mark Jwaideh said:
        March 11, 2019 at 12:10 pm

        Hi Wolfgang,

        Many thanks – I will give it a go.

        I am using this to determine spatially explicit fate factors within watersheds – useful when looking into contaminant fate and transport.

        mark Jwaideh said:
        March 11, 2019 at 2:11 pm

        Apologies, may I confirm that within the function downslopesum(FD,G) is G meant to be the flowacc value or the weights raster?

          wschwanghart responded:
          March 11, 2019 at 2:47 pm

          Hi Mark, no, in this case G contains just values that you may want to accumulate in upstream direction. These could also be a GRIDobj of ones. In this case, the resulting grid contains for each cell the number of downstream cells.

    mark Jwaideh said:
    March 14, 2019 at 11:42 am

    Apologies for yet another question. Just wondering if there is a function to convert flow directions obtained from ERSI or the HYDROSheds database to FLOWobj format. I can see there is a function for the reverse.


    Konstantin Hohmann said:
    March 25, 2019 at 2:19 pm

    Hi Wolfgang,

    i want to extract river terraces/relict landscapes from my 10×10 DEM.
    During field work (and DEM slopeanalysis) i mapped several flat areas along the maintrunk of my Streamnetwork. I am looking for a way to track the elevation of these mapped polygons with respect to the baselevel of the mainriver to extract different errosional levels.
    This sounds like a job for SWATHobj.?
    Note: Some of these polygons stretch along the mainriver and so the mean elevation is not a good metric. Is a running mean with a moving window of lets say 50-100 m a good approach?
    Another question: Is it possible to map those relict landscapes via SWATH profiles with slope as parameter so that in a elevation-distance-plot i can qualitatively interpret different levels?

    Can you think of a solution? Thank you very much!

    BR, Konstantin

      wschwanghart responded:
      March 27, 2019 at 3:52 pm

      Hi Konstantin,
      SWATHobj is probably the best approach. Alternatively, you can use maplateral:
      Did you differentiate the terrace levels in the field? If yes, I’d clip the DEM to the different terraces (e.g. polygon2GRIDobj), map the clipped DEM to the trunk river (e.g. maplateral), and then use the function STREAMobj/smooth to filter the data. The smoothing factor K will depend on the amount of scatter.
      Cheers, Wolfgang

    Kai said:
    May 15, 2019 at 6:21 am

    Hi Wolfgang,

    Is there any way to plot NaNs as transparent (no color) when using imageschs.m to plot up GRIDobj. It would be very helpful for me, if it could be done.



      wschwanghart responded:
      June 3, 2019 at 12:43 pm

      By, default this should be done.

      DEM = GRIDobj('srtm_bigtujunga30m_utm11.tif');
      FD  = FLOWobj(DEM);
      S   = STREAMobj(FD,'minarea',1000);
      S   = klargestconncomps(S);
      D   = drainagebasins(FD,S);
      DEM = clip(DEM,D>0);
      A   = flowacc(FD);
      A   = clip(A,D>0);

      Best regards, Wolfgang

    Patricio Castillo said:
    June 1, 2019 at 3:51 pm

    hi, Topotoolbox works in the student licenced Matlab Version??


    thanks for you time 🙂

      wschwanghart responded:
      June 3, 2019 at 8:35 am

      Yes, it should. But remember to install following Toolboxes that are required: Image Processing Toolbox and Mapping Toolbox. Parallel Toolbox and Optimization Toolbox is also good to have.

    Lingyun said:
    September 11, 2019 at 12:30 pm

    Hello Wolfgang,
    I’m currently learning how to use topotoolbox, but this time I ran into a bit of trouble. When I used GRIDobj, the system indicated that there was no projection , and when I followed its instructions, it said, “undefined function ‘reproject2utm’ corresponding to input parameters of type ‘char'”.Could you help me to solve this problem?Thank you.
    Kind regard,

      wschwanghart responded:
      September 11, 2019 at 1:32 pm

      Hi Lingyun, is it possible that you wrote DEM = reproject2utm(‘DEM’,30) instead of DEM = reproject2utm(DEM,30)?

        Lingyun said:
        September 12, 2019 at 11:27 am

        I wrote DEM = reproject2utm(DEM,30) and it now shows undefined function variables

          wschwanghart responded:
          September 12, 2019 at 11:52 am

          Please check again whether TopoToolbox is on your MATLAB path, e.g.

          which reproject2utm

          should return a path.

    Lingyun said:
    September 15, 2019 at 12:39 pm

    Hi Wolfgang.

    Thanks for your reply and help . I will continue to learn this nice toolbox.


    Liu Li said:
    September 21, 2019 at 11:21 am

    Hello Wolfgang,

    Recently, I wanted to use the TopoToolbox to extract the specified two rivers from the DEM data, and then do the river profile analysis. Although there have trunk stream and maximum watershed extraction methods in the tutorial, there are many unnecessary tributaries in this way. Could you help me to solve this problem?


      wschwanghart responded:
      September 21, 2019 at 12:28 pm

      Hi Liu, I have no idea why there are unneccessary tributaries. The function trunk extracts the main channel in a stream network. If there is more than one drainage basin, then there will be many trunk streams. Please also check the STREAMobj methods klargestconncomps, extractconncomps, and modify.
      Best, Wolfgang

        Liu Li said:
        September 25, 2019 at 11:44 am

        Hi, my question is this: I have two rivers here, they have a river capture phenomenon. When doing river profile analysis, I need to put two rivers on the same screen for consideration, rather than just extract the river profile of the largest trunk stream, so how should I do to achieve this step?

          wschwanghart responded:
          September 26, 2019 at 4:11 pm

          Well, without seeing the data, it’s a bit challenging to help. Did you check the function flowpathapp to extract the streams?

    Audrey said:
    December 2, 2019 at 7:57 am

    Hi Wolfgang,
    I want to extract a swath profile from a DEM using the function SWATHobj. I am getting this error after drawing the line on the DEM:

    DEM = GRIDobj(‘srtm_bigtujunga30m_utm11.tif’);
    Undefined function ‘getdistance’ for input arguments of type ‘double’.

    Error in SWATHobj (line 195)
    d0 = getdistance(x0,y0);

    I am working with MATLAB R2017b, I tried other functions to analyse river profiles and they are working pretty well.

      wschwanghart responded:
      December 2, 2019 at 8:29 am

      Hi Audrey, I guess that some folders are not on your MATLAB search path. getdistance is the utilities folder.
      Best, Wolfgang

    Liu said:
    December 9, 2019 at 3:12 am

    Hi Wolfgang ,
    I would want to know if there is any way to extract only the tributaries of interest and number them. It would be very helpful for me, if it could be done.

      wschwanghart responded:
      December 9, 2019 at 7:13 am

      Hi Liu, yes, check the function STREAMobj/modify and STREAMobj2cell.
      Best, Wolfgang

        Liu said:
        December 9, 2019 at 2:32 pm

        I used modify to extract the tributary, but I encountered another problem: how to change the flowacc to the flow of tributary Sm when calculating Chi value?

        DEM = GRIDobj(‘srtm_bigtujunga30m_utm11.tif’);
        FD = FLOWobj(DEM,’preprocess’,’carve’);
        S = STREAMobj(FD,’minarea’,1000);
        Sm = modify(S,’interactive’,’outletselect’);
        A = flowacc(FD);
        c = chitransform(Sm,A,’mn’,0.45);

          wschwanghart responded:
          December 9, 2019 at 3:44 pm

          Why would flow accumulations change? chitransform extracts the drainage areas from the flow accumulation grid at each node of the stream network Sm. Hence, I don’t see a problem with using the code as is.

    William Medwedeff said:
    January 23, 2020 at 4:27 pm

    Hi Wolfgang,

    I have been working with TTLEM and I ran into a problem setting the hillslope diffusion scheme to nonlinear, with a critical slope threshold (i.e. p.diffScheme = ‘imp_nonlin_sc’).

    I get an error in the function nonlinDiff_Sc: “Unrecognized function or variable ‘mexGetsMats_BC_v1′” on line 57.

    It looks like the function mexGetsMats_BC_v1.c is written in C, so it is not recognized by Matlab. Is there an easy way to execute this function from Matlab? I see that MEX files are used to implement Matlab functions and variables in C++. I’m guessing nonlinDiff_Sc.m is just missing a few lines to build a MEX file (?) but this is really pushing the limits of my programming knowledge 🙂

    Sorry if you have already answered this question – I did some searching but couldn’t find it below. FYI I am on OSX, Matlab 2019b, and I have the most recent version of topotoolbox from your github.

    Have any insight?

    PS: Just want to say that I am always amazed at how clean/fast/useful your software is, and I owe you many thanks as I use it all the time for my research!! (I am a PhD student at University of Michigan). Your work really inspires me!

      wschwanghart responded:
      January 24, 2020 at 7:37 am

      Hi William, thanks for your comment!
      My guess is that mexGetMats_BC_v1.c needs to be compiled. We distribute windows binaries of our mex-functions (there are some in the @FLOWobj/private folder and in ttlem), but not for mac. This has not been reported so far, because there are also m-files that will be used if MATLAB cannot find the mex functions. However, mexGetMats_BC_v1 apparently only exists as windows binary (*.mexw64). Hence, you would need to compile it yourself. Check the mcc function which should do the job.
      Cheers, Wolfgang

        William Medwedeff said:
        January 29, 2020 at 4:26 pm

        Thanks for the quick reply. Now all fixed!

    Federica Ferrarini said:
    March 5, 2020 at 6:19 am

    Hi Wolfgang,

    I’m using Topotoolbox for fluvial network analysis and I would like to exploit it also to plot terrace surfaces.

    I tried at the beginning using both Topotoolbox and GIS (ArcMap): I extracted terrace (point) elevations for the different clipped surfaces (different terrace orders),
    I computed a near distance between points and the river mouth and I plotted, to follow, the points along the river. That was not a good attempt because the river has meanders and the along-distance of the river and of the terraces do not step correctly together.

    I luckily run into a reply that you gave to another guy on Topotoolbox WordPress where you suggested to use ‘maplateral’ starting from clipped dem.

    I tried using these code lines (only for one terrace order):

    DEM = GRIDobj (”….path\wideA_Lidar_prj.tif’);
    FD= FLOWobj(DEM,’preprocess’,’fill’, ‘mex’, true);
    A= flowacc(FD);
    S=STREAMobj (FD, A >800000);
    river_lid = klargestconncomps(S,1);
    river_lid_TR = trunk (river_lid);
    river_lid_TR40=modify(river_lid_TR,’distance’,[0 40000]);

    DEMACT = GRIDobj (‘….path\ACT_left_side_UTM33.tif’); %to load the clipped terraces

    MSACT = GRIDobj2polygon(DEMACT);
    PACT = polygon2GRIDobj(DEMACT,MSACT,’ID’);

    PACTr = resample(PACT,DEM) % I had ‘validate aligments’ troubles so I tried to resample PACTr
    VACT=validatealignment(river_lid_TR40,PACTr) % to verify the alignmment

    gACT = maplateral(river_lid_TR40,PACTr,500,{@mean}); % try mapping within a buffer of 500m
    plotdz (river_lid_TR40,DEM)
    ylim([0 400])
    hold on
    plotdz(river_lid_TR40, gACT);

    So, admitting that the code is right, what I obtained was a sort of constant value (mean) plotted onto the longitudinal profile.
    Since ‘maplateral’ request an aggregate function, what type of function could be adequate to plot simply the elevation of points in the gridobj PACTr?
    Could you give, please, some suggestion?

    Apologies for the long message.

    Thank you in advance for the help that you would grant me!

      wschwanghart responded:
      March 5, 2020 at 9:48 am

      Hi Federica, sorry for my late reply. Your message landed in the spam folder.

      As far as I can see, your code looks quite alright. One thing I am unsure about is whether the grid PACTr actually contains the clipped elevations, or just the ID of the gridded polygons. Thus, maplateral should actually be called like
      gACT = maplateral(river_lid_TR40,DEMACT,500,{@mean});
      The aggregation function takes the values from all pixels that are nearest to a stream pixel, and maps them to that particular pixel. This means, the function takes a vector of values and returns a scalar. @mean does that, but you may also choose @median or other metrics.

      Hope this helps. Best regards, Wolfgang

        ferrarinifederica said:
        March 11, 2020 at 5:05 am

        Thank you very much.
        Now fixed!

    Mehdi said:
    March 6, 2020 at 10:59 am

    Hi Wolfgang.
    I’m trying to load my aster GDEMV3 with GRIDobj but i get this error:
    GRIDobj cannot read the TIF-file because it does not have a tfw-file.

      wschwanghart responded:
      March 6, 2020 at 11:12 am

      Do you have the Mapping Toolbox? Check the command


      . If you don’t have the mapping toolbox, then TopoToolbox cannot read geotiffs. As a solution, GRIDobj then searches for a .tfw file, which contains the coordinates of the grid. You can write world-files using ArcGIS or QGIS.

        Mehdi said:
        March 6, 2020 at 2:06 pm

        Thank you for your answer.
        I have it. It just happens when i use aster GDEMV3. other DEMs load correctly.

    lucianopav96 said:
    April 9, 2020 at 10:14 pm

    Hi Wolfgang!

    I don’t know if this is the rigth section for this question.

    I have a DEM.
    I’m trying to copy from a drainagebasins GRIDobj matrix (ex. DB. Z) the DEM. Z values. The problem is that all the values are integers because the DB. Z matrix is type uint32 of default. How can i convert the matrix from uint32 to double?

    Best Reguards,

      wschwanghart responded:
      April 10, 2020 at 7:02 pm

      Hi Luciano,
      that’s easy.

      DB.Z = double(DB.Z);

      Hope this works. Best, Wolfgang

        lucianopav96 said:
        April 11, 2020 at 10:05 am

        In this way, i don’t know why doesn’t work for my problem, because DB.Z wasn’t a vector. But anyway I solved with a loop for.

        Thank you for your time,
        Best, Luciano

    Seb said:
    April 29, 2020 at 7:00 am

    Hi Wolfgang,
    I was wondering how I can run fillsinks (and flowacc etc) for an area that is not rectangular. I dont have a DEM in certain areas and therefor only want to run my code for a polygon-shaped area. If the surrounding areas of the polygon in the GridObj are NaN or 0 this leads to erroneous sinks at the “boundary”. How could I handle this?

    I tried using the clip function (not sure if its meant for this):
    Poly_shp = shaperead(‘….shp’);
    P = polygon2GRIDobj(DEM,Poly_shp);
    DEMc = clip(DEM,P.Z);
    But I receive the following error (also when using DEMc = clip(DEM,Poly_shp):
    Undefined function ‘clip’ for input arguments of type ‘GRIDobj’.


      wschwanghart responded:
      April 29, 2020 at 9:01 am

      Hi Sebastian, it’s strange that there seems to be no clip function. Do you have the latest version of TopoToolbox (
      Setting values to zero outside the mask is not a good idea, because a flow network will be computed for these areas. However, setting these values to nan should be just fine. That’s actually what the function clip does. So, the problem seems to be the missing function clip, and my guess is that you do not have the latest version of TT.
      Hope this helps,

        Seb said:
        April 29, 2020 at 12:33 pm

        Thanks Wolfgang, seems like addpath didn’t work quite right. Clip and everything else works fine now. Thanks for the great tools! 🙂

          wschwanghart responded:
          April 29, 2020 at 2:11 pm

          Super. Glad to hear everything works well. Cheers, Wolfgang

    Ahmed said:
    October 29, 2020 at 1:40 pm

    Hi W. Schwanghart, First of all I would like to thank you for colossal effort to provide such a useful Toolbox,
    It was recommended for me by a friend who has been working with it.
    I’m just a beginner in Matlab, but I’m following courses to get used to,
    I’m writing this comment to get some help because I have some issues in installing the toolbok, I looked up for some tuto, but I did’t find. After having downloaded the folders Topotoolbox-master ; Topotoolbox 2-3, and Tototoolbox.mltbx, I managed to put the folders in the Matlab search path, but I still don’t know where to put Topotoolbox.mltbx.
    Could you please show me how to make these things right,

    Oscar Camargo said:
    May 28, 2021 at 5:46 pm

    Dear Wolfgang

    I hope you are well, we are using topotoolbox and calculating a heavy flow accumulation using the precipitation, both the DEM, and the precipitation we created at the same resolution, projection, etc and did resampling. However we get the following validatealignment error:

    Error using FLOWobj/validatealignment (line 42) FLOWobj and GRIDobj do not align each other. Make sure that both instances have the same spatial reference. Both variables are deemed to have the same reference if their properties ‘size’ and ‘refmat’ are both equal.

    We checked and indeed there are differences in the values but they are very small, is there any recommendation you can give us to solve this problem.

    Thanks in advance, all the best


      wschwanghart responded:
      May 28, 2021 at 6:15 pm

      Hi Oscar, yes, that’s a common problem. I commonly then resample the grid so that it perfectly aligns the DEM and FLOWobj.

      RAIN = resample(RAIN,DEM);

      Hope this works for you.
      Best, Wolfgang

    Vishal Vats said:
    July 3, 2021 at 4:49 pm

    Hello Dr. Schwanghart,
    Congratulation for this amazing software. I have one querey that can we use topotoolbox to extract the value of elevation along with lat lon from DEM in an excel file?

      wschwanghart responded:
      July 6, 2021 at 8:26 am

      Dear Vishal,

      if you have some points as x,y-coordinates, then extract z values with the function

      z = interp(DEM,x,y);

      Latitude and longitude are retrieved by transforming x,y-coordinates to a geographic coordinate system.

      [lat,lon] = minvtran(DEM.georef.mstruct, x,y);

      Best regards, Wolfgang

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

This site uses Akismet to reduce spam. Learn how your comment data is processed.