diff --git a/_freeze/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk/execute-results/html.json b/_freeze/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk/execute-results/html.json new file mode 100644 index 0000000..7203a76 --- /dev/null +++ b/_freeze/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk/execute-results/html.json @@ -0,0 +1,12 @@ +{ + "hash": "b7c0a59c0496131d3ba35ee7a1d01dc8", + "result": { + "engine": "jupyter", + "markdown": "---\ntitle: \"Modeling Movement in GRASS\"\nauthor: \"Michael Barton & Anna Petrasova\"\ndate: 2025-02-10\ndate-modified: today\nformat: \n html:\n embed-resources: true\n toc: true\n code-tools: true\n code-copy: true\n code-fold: false\npage-layout: article\ncategories: [intermediate, advanced, GUI, raster, cost surface, least cost path]\ndescription: Generating a cumulative cost surface and least cost path with *r.walk* and *r.path* to model movement by walking across a landscape.\nimage: images/thumbnail.webp\nengine: jupyter\nexecute: \n eval: false\njupyter: python3\n---\n\n\nGRASS has sophisticated tools to model movement across terrain, including\n[r.cost](https://grass.osgeo.org/grass-stable/manuals/r.cost.html),\n[r.walk](https://grass.osgeo.org/grass-stable/manuals/r.walk.html),\n[r.drain](https://grass.osgeo.org/grass-stable/manuals/r.drain.html), and\n[r.path](https://grass.osgeo.org/grass-stable/manuals/r.path.html).\nIn this tutorial, we will use *r.walk* and *r.path* to determine the best\nwalking path to reach a destination, such as a hospital.\n\n::: {.callout-note title=\"Dataset\"}\nThis tutorial uses data from Flagstaff, Arizona, USA,\nbut it can be completed with\na [North Carolina dataset](https://grass.osgeo.org/download/data/).\nWe will use the\n*elevation* DEM, the *hospitals* vector points file (any other vector point\nfile will serve), and the *landuse* raster map.\n:::\n\n::: {.callout-note title=\"Don't know how to get started?\"}\nIf you are not sure how to get started with GRASS using its garphical user interface\nor using Python, checkout the tutorials\n[Get started with GRASS & Python in Jupyter Notebooks](../get_started/fast_track_grass_and_python.qmd)\nand\n[Get started with GRASS GIS GUI](../get_started/fast_track.qmd)\n:::\n\n## Overview: cost surfaces\n\nA *cost surface* is a raster where each cell represents the \"cost\"\nor difficulty to move across landscape\n\nA *cumulative cost surface* shows the total accumulated cost of moving\nfrom a starting point to a location.\nCumulative cost surfaces are also used to find a *least cost path*\nbetween a location and the starting point.\n\n::::: grid\n::: g-col-6\n\nGRASS [r.walk](https://grass.osgeo.org/grass-stable/manuals/r.walk.html)\ngenerates a cumulative cost surface based on:\n\n- Naismith's rule for walking times based on slope (e.g., walking up\n the hill is slower) and\n- a cost surface called *friction*, where each cell's value\n represents additional walking time in seconds.\n:::\n\n::: g-col-6\n![Cumulative cost surface](images/GRASS_movement_0.webp)\n:::\n:::::\n\n## Overview of modeling movement with a cumulative cost surface and least cost path\n\n::::: grid\n::: g-col-6\n\n1. Start with a *digital elevation model (DEM)* to determine movement costs.\n2. Create a *friction map* with a value of 0 (or other values for additional\n movement costs).\n3. Select one or more *start point(s)*.\n4. Create a *cost surface*.\n5. A *least cost* path can then be calculated between any point on\n the cost surface and the start point.\n:::\n\n::: g-col-6\n![DEM](images/GRASS_movement_1.webp)\n:::\n:::::\n\n### Friction map\n\nMost cost surfaces are based on movement across terrain\n(e.g., travel up or down slopes or across flat areas).\nA *friction surface* in *r.walk* allows you to account for other costs to movement.\nThe friction surface cell value represents additional walking time in seconds\nto cross 1 meter distance (e.g., going through a forest has higher friction\nthan walking on a road).\n\nA friction surface can be made multiple ways. For our first cost surface,\nwe will focus only on terrain and will only need a friction surface\nwith the value of 0 in each cell.\nWe can easily make that with the GRASS *map calculator*, which\nlets you create a map, alter a map or combine maps using *map algebra*.\n\nThe map calculator can be used from a GUI or as a command\n[r.mapcalc](https://grass.osgeo.org/grass-stable/manuals/r.mapcalc.html)\nfrom shell or Python.\n\nMake sure the [computational region](https://grass.osgeo.org/grass-stable/manuals/g.region.html)\nmatches the DEM and then create the `friction0` map with 0 in every cell:\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n::::: grid\n::: g-col-6\n\n1. Add the *elevation* raster to the Layer Manager.\n\n2. Right click *elevation* in Layer Manager.\n\n3. Select **Set computational region from selected map**.\n\n4. Open the **Raster Map Calculator** from the top toolbar (\"abacus\" icon).\n\n5. Fill in the output map name as *friction0* and enter 0 in the expression field.\n\n6. Press Run.\n:::\n\n::: g-col-6\n![Raster map calculator](images/GRASS_movement_2.webp)\n:::\n:::::\n\n#### Command line\n\n\n```{bash}\n# Set the region to match the elevation map\ng.region raster=elevation\n\n# Create a friction map with a value of 0 in all cells in the region\nr.mapcalc \"friction0 = 0\"\n```\n\n\n#### Python\n\n::: {#a7ab0e15 .cell execution_count=1}\n``` {.python .cell-code}\n# Set the region to match the elevation map\ngs.run_command(\"g.region\", raster=\"elevation\")\n\n# Create a friction map with a value of 0 in all cells in the region\ngs.mapcalc(\"friction0 = 0\")\n```\n:::\n\n\n:::\n\n### Choose a start point\n\nNow that we have a DEM and a friction map, all we need is a start point:\nthe point from which to calculate movement costs.\n\nA start point can be one or more vector points, raster cells,\nor even pairs of coordinates.\nWe'll make a start point from the Flagstaff Medical Center,\nfound in the *hospitals* vector map.\n\nUse [v.extract]((https://grass.osgeo.org/grass-stable/manuals/v.extract.html))\nto select the Flagstaff Medical Center from the *hospitals* points file:\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n::::: grid\n::: g-col-6\n\nUse the **Attribute Table Manager** for *hospitals*:\n\n1. Display *hospitals* map by adding it to the Layer Manager from Data Catalog.\n2. Right click and open the Attribute Table Manager\n3. Select \"Flagstaff Medical Center\" record\n4. Right click and select *Extract selected features*\n5. Name the new map *FMC*\n:::\n\n::: g-col-6\n![Attribute Table Manager for hospitals vector points map](images/GRASS_movement_3.webp) \n\n![Table row context menu](images/GRASS_movement_4.webp){width=\"60%\"}\n:::\n:::::\n\n#### Command line\n\nUse the v.extract tool to create a vector point map named *FMC* from the *hospitals* map:\n\n\n```{bash}\nv.extract input=hospitals type=point where=\"FACILITY_N = 'FLAGSTAFF MEDICAL CENTER'\" output=FMC\n```\n\n\n#### Python\n\nUse the v.extract tool to create a vector point map named *FMC* from the *hospitals* map:\n\n::: {#184d4ce3 .cell execution_count=2}\n``` {.python .cell-code}\ngs.run_command(\"v.extract\", input=\"hospitals\", type=\"point\", where=\"FACILITY_N = 'FLAGSTAFF MEDICAL CENTER'\", output=\"FMC\")\n```\n:::\n\n\n:::\n\nDisplay FMC with a color and size to see it better\n\n![](images/GRASS_movement_5.webp)\n\n### Compute the cumulative cost surface\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\nUse the *r.walk* tool from the **Raster/Terrain Analysis** menu:\n\n::::: grid\n::: g-col-6\n\n1. Enter the input and output parameters: **elevation**=*elevation*, **friction**=*friction0*, and **output**=*FMC_cost_seconds*:\n\n![](images/GRASS_movement_6.webp)\n:::\n\n::: g-col-6\n2. Enter the name of a directions map **outdir**=*FMC_directions* for creating a least cost path:\n\n![](images/GRASS_movement_7.webp)\n\n:::\n:::::\n\n::::: grid\n::: g-col-6\n3. Enter the start point (**start_points**=*FMC*).\n\n![](images/GRASS_movement_8.webp){fig-align=\"left\" width=\"100%\"}\n:::\n::: g-col-6\n4. Optional: control the spatial extent of cost surface:\n\n![](images/GRASS_movement_9.webp)\n:::\n::::::\n\n::::: grid\n::: g-col-6\n5. Optional: adjust movement parameter settings:\n\n![](images/GRASS_movement_10.webp)\n:::\n::: g-col-6\n6. Recommended: select *knight's move* for more accurate cost and direction:\n\n![](images/GRASS_movement_11.webp)\n:::\n:::::\n\n::: {.callout-note title=\"Tip\"}\nClick the \"copy\" button to copy the GRASS command.\nYou can save it in a text file for later reuse or to document your work.\n:::\n\n#### Command line\n\nUse the r.walk command to generate the cumulative cost surface:\n\n\n```{bash}\nr.walk elevation=elevation friction=friction0 output=FMC_cost_seconds outdir=FMC_directions start_points=FMC -k\n```\n\n\n#### Python\n\nUse the r.walk command to generate the cumulative cost surface:\n\n::: {#ce7dc2b1 .cell execution_count=3}\n``` {.python .cell-code}\ngs.run_command(\"r.walk\",\n elevation=\"elevation\",\n friction=\"friction0\",\n output=\"FMC_cost_seconds\",\n outdir=\"FMC_directions\",\n start_points=\"FMC\",\n flags=\"k\")\n```\n:::\n\n\n:::\n\nEach raster cell value in the cost surface is the time in\nseconds to walk from FMC to that cell given the terrain.\n\n![Displaying the cost surface over a shaded relief map using the [d.shade](https://grass.osgeo.org/grass-stable/manuals/d.shade.html) tool and a relief map of elevation made with [r.relief](https://grass.osgeo.org/grass-stable/manuals/d.shade.html).](images/GRASS_movement_12.webp)\n\n### Movement across a cumulative cost surface\n\nYou can transform the walking time in seconds to hours by dividing the map by 3600 using the map calculator\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n::::: grid\n::: g-col-6\n\n1. Open the map calculator.\n2. Enter the name of the new map of walking time in hours: *FMC_cost_hours*.\n3. Enter `FMC_cost_seconds / 3600` into the expressions field\n4. Press Run.\n:::\n::: g-col-6\n![](images/GRASS_movement_13.webp)\n:::\n::::::\n\n::::: grid\n::: g-col-6\nWe can query or filter the cost surface to areas of equivalent walking time.\n\nFor example, to show the area within a 2 hour walk of FMC, filter the\ncumulative cost surface display using the *d.rast* tool in the *layer manager*\nto show 0-2 hours and set the opacity of the cost surface to 50% to see\nthe underlying terrain.\n:::\n::: g-col-6\n![](images/GRASS_movement_15.webp)\n:::\n:::::\n\n#### Command line\n\n\n```{shell}\nr.mapcalc \"FMC_cost_hours = FMC_cost_seconds / 3600\"\n```\n\n\n#### Python\n\n::: {#d719bb59 .cell execution_count=4}\n``` {.python .cell-code}\ngs.mapcalc(\"FMC_cost_hours = FMC_cost_seconds / 3600\")\n```\n:::\n\n\n:::\n\n![Shaded area represents terrain reached within a 2 hour walk from FMC](images/GRASS_movement_14.webp)\n\n### Least cost path\n\nWe can also plot a *least cost path* (LCP), which is the least costly\n(shortest time) route between any point on the cumulative cost surface and FMC\n\nImagine a stranded hiker NE of Flagstaff who has to walk to FMC.\nWhat path would take the least time? \n\n![Start and end point](images/GRASS_movement_16.webp)\n\nTo create a LCP in GRASS, we will use [r.path](https://grass.osgeo.org/grass-stable/manuals/r.path.html)\n(also under the Raster/Terrain analysis menu).\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n::::: grid\n::: g-col-6\n\n1. We need to use the directions map **input**=*FMC_directions*\nwe made when we created the cost surface.\n:::\n::: g-col-6\n![](images/GRASS_movement_17.webp)\n:::\n::::::\n\n::::: grid\n::: g-col-6\n2. The coordinates of the hiker **start_coordinates**=*477476,13914951*\n(we could also use a pre-defined vector point).\n:::\n::: g-col-6\n![](images/GRASS_movement_18.webp)\n:::\n:::::\n\n:::: grid\n::: g-col-6\n3. And specify the name of the output LCP **vector_path**=*LCP_cumulative*.\n:::\n::: g-col-6\n![](images/GRASS_movement_19.webp)\n:::\n:::::\n\n#### Command line\n\n\n```{shell}\nr.path input=FMC_directions format=auto vector_path=LCP_cumulative start_coordinates=477476,13914951\n```\n\n\n#### Python\n\n::: {#d39adc47 .cell execution_count=5}\n``` {.python .cell-code}\ngs.run_command(\"r.path\",\n input=\"FMC_directions\",\n format=\"auto\",\n vector_path=\"LCP_cumulative\",\n start_coordinates=[477476, 13914951])\n```\n:::\n\n\n:::\n\n![Resulting least cost path](images/GRASS_movement_20.webp)\n\n## Adding a friction map to movement costs\n\nA **friction map** can be used to incorporate other factors than\njust terrain into creating a cost surface and modeling movement.\nThe value of each cell in a friction surface raster is the amount\nof walking time, in seconds/meter, *in addition to* the time needed\nbecause of terrain.\n\nWe can reclassify the land use map (*landuse*) to create a friction surface\nof the amount of extra time it would take to walk through different kinds of land cover. \n\n![Reclassification of *landuse* to show major land cover types](images/GRASS_movement_21.webp)\n\n### Reclassification of landuse to a friction map\n\nA standard walking speed across flat terrain is about 5 km/hr = 0.72 sec/m.\nWe might then estimate that it would take an additional:\n\n- 3 sec/m to walk through dense pinyon-juniper woodland,\n- 1 sec/m to cross conifer forest,\n- 2 sec/m to wander through urban Flagstaff,\n- 5 sec/m to clamber over lava fields, and\n- 10 sec/m to try and cross water\n\nWe can create this friction map by reclassifying the *landuse* map using\n[r.reclass](https://grass.osgeo.org/grass-stable/manuals/r.reclass.html)\nto assign new friction values to the existing land cover categories.\n\n::: {.panel-tabset group=\"language\"}\n\n## GUI\n\nThe [r.reclass](https://grass.osgeo.org/grass-stable/manuals/r.reclass.html) tool\nis found under the Raster/Change category… menu.\n\n::::: grid\n::: g-col-6\n\n1. Enter *landuse* for the raster to be reclassified\n\n2. Enter *friction_reclassified* for the output reclassified map\n\n3. Enter the reclass rules directly in the text box or from a saved text file. Use and \\* symbol to represent everything not covered by the specific reclass rules\n\n > 11 90 95 = 10 water \n > 21 thru 24 = 2 urban \n > 31 = 5 lava \n > 41 thru 43 = 1 conifer forest \n > 52 = 3 pinyon juniper woodland \n > \\* = 0 no friction \n\n:::\n\n::: g-col-6\nThis creates a new friction map named *friction_landcover*.\n![](images/GRASS_movement_22.webp)\n:::\n:::::\n\n## Command line\n\n\n```{shell}\nr.reclass input=landuse output=friction_landcover rules=- << EOF\n11 90 95 = 10 water\n21 thru 24 = 2 urban\n31 = 5 lava\n41 thru 43 = 1 conifer forest\n52 = 3 pinyon juniper woodland\n* = 0 no friction\nEOF\n```\n\n\n## Python\n\n::: {#eefdb65d .cell execution_count=6}\n``` {.python .cell-code}\nrules = \"\"\"\\\n11 90 95 = 10 water\n21 thru 24 = 2 urban\n31 = 5 lava\n41 thru 43 = 1 conifer forest\n52 = 3 pinyon juniper woodland\n* = 0 no friction\n\"\"\"\ngs.write_command(\"r.reclass\",\n input=\"landuse\",\n output=\"friction_landcover\",\n rules=\"-\",\n stdin=rules)\n```\n:::\n\n\n:::\n\n![The new friction map](images/GRASS_movement_23.webp)\n\n### Modifying a cumulative cost surface with a friction map\n\nWe can now create a new cumulative cost surface using this new friction map\nand convert it from seconds to hours, as we did before:\n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n1. Follow the procedures in the [Compute the cumulative cost surface\nsection](#compute-the-cumulative-cost-surface) and substitute the new\n*friction_landcover* map for the *friction0* map used previously.\n\n2. Convert the cumulative cost surface to hours instead of seconds\nfollowing the procedures in the [Movement across a cumulative cost surface section](#movement-across-a-cumulative-cost-surface).\n\n#### Command line\n\n\n```{bash}\nr.walk elevation=elevation friction=friction_landcover output=FMC_vegcost_seconds outdir=FMC_vegcost_directions start_points=FMC -k\nr.mapcalc \"FMC_vegcost_hours = FMC_vegcost_seconds / 3600\"\n```\n\n\n#### Python\n\n::: {#1de5ca1b .cell execution_count=7}\n``` {.python .cell-code}\ngs.run_command(\"r.walk\",\n elevation=\"elevation\",\n friction=\"friction_landcover\",\n output=\"FMC_vegcost_seconds\",\n outdir=\"FMC_vegcost_directions\",\n start_points=\"FMC\",\n flags=\"k\")\ngs.mapcalc(\"FMC_vegcost_hours = FMC_vegcost_seconds / 3600\")\n```\n:::\n\n\n:::\n\n![Cumulative cost surface with landcover friction](images/GRASS_movement_24.webp)\n\n### Least cost paths and friction\n\nWe can create a new LCP from the stranded hiker to FMC over terrain\nwhere land cover also affects the cost of movement. \n\n::: {.panel-tabset group=\"language\"}\n\n#### GUI\n\n1. Follow the procedures in the [Least cost path section](#least-cost-path), substituting the new direction map made along with the cumulative cost surface with land cover friction map.\n\n2. Give this new LCP the name *LCP_vegcost* to distinguish from the previous LCP.\n\n#### Command line\n\n\n```{bash}\nr.path input=FMC_vegcost_directions format=auto vector_path=LCP_vegcost start_coordinates=477476,13914951\n```\n\n\n#### Python\n\n::: {#94266227 .cell execution_count=8}\n``` {.python .cell-code}\ngs.run_command(\"r.path\",\n input=\"FMC_vegcost_directions\",\n format=\"auto\",\n vector_path=\"LCP_vegcost\",\n start_coordinates=[477476, 13914951])\n```\n:::\n\n\n:::\n\nIt takes more time to reach FMC if the hiker has to navigate dense vegetation as well as terrain.\n\nIn the map below, terrain colored by land cover, the original LCP with terrain only is shown by the blue line. The LCP with a land cover friction map is shown by the heavier yellow line\n\n![Comparing LCP with terrain only and with land cover friction](images/GRASS_movement_25.webp)\n\n## Acknowledgement\n\nThe development of this tutorial was funded by the U.S.\n[National Science Foundation (NSF)](https://www.nsf.gov/),\naward [2303651](https://www.nsf.gov/awardsearch/showAward?AWD_ID=2303651).\n\n", + "supporting": [ + "least_cost_path_r_walk_files" + ], + "filters": [], + "includes": {} + } +} \ No newline at end of file diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_0.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_0.webp new file mode 100644 index 0000000..4df25ae Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_0.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_1.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_1.webp new file mode 100644 index 0000000..3ced818 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_1.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_10.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_10.webp new file mode 100644 index 0000000..42b0fff Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_10.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_11.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_11.webp new file mode 100644 index 0000000..1a36ae8 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_11.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_12.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_12.webp new file mode 100644 index 0000000..4df25ae Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_12.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_13.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_13.webp new file mode 100644 index 0000000..100c98d Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_13.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_14.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_14.webp new file mode 100644 index 0000000..a465640 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_14.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_15.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_15.webp new file mode 100644 index 0000000..758284d Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_15.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_16.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_16.webp new file mode 100644 index 0000000..241b2f5 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_16.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_17.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_17.webp new file mode 100644 index 0000000..1e5f860 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_17.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_18.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_18.webp new file mode 100644 index 0000000..a3c3d1b Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_18.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_19.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_19.webp new file mode 100644 index 0000000..11d127e Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_19.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_2.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_2.webp new file mode 100644 index 0000000..5b31aba Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_2.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_20.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_20.webp new file mode 100644 index 0000000..1d8afcf Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_20.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_21.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_21.webp new file mode 100644 index 0000000..32877ad Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_21.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_22.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_22.webp new file mode 100644 index 0000000..b75889f Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_22.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_23.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_23.webp new file mode 100644 index 0000000..f9d720a Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_23.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_24.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_24.webp new file mode 100644 index 0000000..3e832fb Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_24.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_25.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_25.webp new file mode 100644 index 0000000..8a09280 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_25.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_3.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_3.webp new file mode 100644 index 0000000..a75b1ff Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_3.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_4.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_4.webp new file mode 100644 index 0000000..61fef9a Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_4.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_5.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_5.webp new file mode 100644 index 0000000..98ee2ce Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_5.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_6.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_6.webp new file mode 100644 index 0000000..af8d51c Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_6.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_7.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_7.webp new file mode 100644 index 0000000..85a95f3 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_7.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_8.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_8.webp new file mode 100644 index 0000000..5c38723 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_8.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_9.webp b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_9.webp new file mode 100644 index 0000000..52416bc Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/GRASS_movement_9.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/images/thumbnail.webp b/content/tutorials/least_cost_path_r_walk/images/thumbnail.webp new file mode 100644 index 0000000..35c9ff1 Binary files /dev/null and b/content/tutorials/least_cost_path_r_walk/images/thumbnail.webp differ diff --git a/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk.qmd b/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk.qmd new file mode 100644 index 0000000..69f2d88 --- /dev/null +++ b/content/tutorials/least_cost_path_r_walk/least_cost_path_r_walk.qmd @@ -0,0 +1,598 @@ +--- +title: "Modeling Movement in GRASS" +author: "Michael Barton & Anna Petrasova" +date: 2025-02-10 +date-modified: today +format: + html: + embed-resources: true + toc: true + code-tools: true + code-copy: true + code-fold: false +page-layout: article +categories: [intermediate, advanced, GUI, raster, cost surface, least cost path] +description: Generating a cumulative cost surface and least cost path with *r.walk* and *r.path* to model movement by walking across a landscape. +image: images/thumbnail.webp +engine: jupyter +execute: + eval: false +jupyter: python3 +--- + +GRASS has sophisticated tools to model movement across terrain, including +[r.cost](https://grass.osgeo.org/grass-stable/manuals/r.cost.html), +[r.walk](https://grass.osgeo.org/grass-stable/manuals/r.walk.html), +[r.drain](https://grass.osgeo.org/grass-stable/manuals/r.drain.html), and +[r.path](https://grass.osgeo.org/grass-stable/manuals/r.path.html). +In this tutorial, we will use *r.walk* and *r.path* to determine the best +walking path to reach a destination, such as a hospital. + +::: {.callout-note title="Dataset"} +This tutorial uses data from Flagstaff, Arizona, USA, +but it can be completed with +a [North Carolina dataset](https://grass.osgeo.org/download/data/). +We will use the +*elevation* DEM, the *hospitals* vector points file (any other vector point +file will serve), and the *landuse* raster map. +::: + +::: {.callout-note title="Don't know how to get started?"} +If you are not sure how to get started with GRASS using its garphical user interface +or using Python, checkout the tutorials +[Get started with GRASS & Python in Jupyter Notebooks](../get_started/fast_track_grass_and_python.qmd) +and +[Get started with GRASS GIS GUI](../get_started/fast_track.qmd) +::: + +## Overview: cost surfaces + +A *cost surface* is a raster where each cell represents the "cost" +or difficulty to move across landscape + +A *cumulative cost surface* shows the total accumulated cost of moving +from a starting point to a location. +Cumulative cost surfaces are also used to find a *least cost path* +between a location and the starting point. + +::::: grid +::: g-col-6 + +GRASS [r.walk](https://grass.osgeo.org/grass-stable/manuals/r.walk.html) +generates a cumulative cost surface based on: + +- Naismith's rule for walking times based on slope (e.g., walking up + the hill is slower) and +- a cost surface called *friction*, where each cell's value + represents additional walking time in seconds. +::: + +::: g-col-6 +![Cumulative cost surface](images/GRASS_movement_0.webp) +::: +::::: + +## Overview of modeling movement with a cumulative cost surface and least cost path + +::::: grid +::: g-col-6 + +1. Start with a *digital elevation model (DEM)* to determine movement costs. +2. Create a *friction map* with a value of 0 (or other values for additional + movement costs). +3. Select one or more *start point(s)*. +4. Create a *cost surface*. +5. A *least cost* path can then be calculated between any point on + the cost surface and the start point. +::: + +::: g-col-6 +![DEM](images/GRASS_movement_1.webp) +::: +::::: + +### Friction map + +Most cost surfaces are based on movement across terrain +(e.g., travel up or down slopes or across flat areas). +A *friction surface* in *r.walk* allows you to account for other costs to movement. +The friction surface cell value represents additional walking time in seconds +to cross 1 meter distance (e.g., going through a forest has higher friction +than walking on a road). + +A friction surface can be made multiple ways. For our first cost surface, +we will focus only on terrain and will only need a friction surface +with the value of 0 in each cell. +We can easily make that with the GRASS *map calculator*, which +lets you create a map, alter a map or combine maps using *map algebra*. + +The map calculator can be used from a GUI or as a command +[r.mapcalc](https://grass.osgeo.org/grass-stable/manuals/r.mapcalc.html) +from shell or Python. + +Make sure the [computational region](https://grass.osgeo.org/grass-stable/manuals/g.region.html) +matches the DEM and then create the `friction0` map with 0 in every cell: + +::: {.panel-tabset group="language"} + +#### GUI + +::::: grid +::: g-col-6 + +1. Add the *elevation* raster to the Layer Manager. + +2. Right click *elevation* in Layer Manager. + +3. Select **Set computational region from selected map**. + +4. Open the **Raster Map Calculator** from the top toolbar ("abacus" icon). + +5. Fill in the output map name as *friction0* and enter 0 in the expression field. + +6. Press Run. +::: + +::: g-col-6 +![Raster map calculator](images/GRASS_movement_2.webp) +::: +::::: + +#### Command line + +```{bash} +# Set the region to match the elevation map +g.region raster=elevation + +# Create a friction map with a value of 0 in all cells in the region +r.mapcalc "friction0 = 0" +``` + +#### Python + +```{python} +# Set the region to match the elevation map +gs.run_command("g.region", raster="elevation") + +# Create a friction map with a value of 0 in all cells in the region +gs.mapcalc("friction0 = 0") +``` + +::: + +### Choose a start point + +Now that we have a DEM and a friction map, all we need is a start point: +the point from which to calculate movement costs. + +A start point can be one or more vector points, raster cells, +or even pairs of coordinates. +We'll make a start point from the Flagstaff Medical Center, +found in the *hospitals* vector map. + +Use [v.extract]((https://grass.osgeo.org/grass-stable/manuals/v.extract.html)) +to select the Flagstaff Medical Center from the *hospitals* points file: + +::: {.panel-tabset group="language"} + +#### GUI + +::::: grid +::: g-col-6 + +Use the **Attribute Table Manager** for *hospitals*: + +1. Display *hospitals* map by adding it to the Layer Manager from Data Catalog. +2. Right click and open the Attribute Table Manager +3. Select "Flagstaff Medical Center" record +4. Right click and select *Extract selected features* +5. Name the new map *FMC* +::: + +::: g-col-6 +![Attribute Table Manager for hospitals vector points map](images/GRASS_movement_3.webp) + +![Table row context menu](images/GRASS_movement_4.webp){width="60%"} +::: +::::: + +#### Command line + +Use the v.extract tool to create a vector point map named *FMC* from the *hospitals* map: + +```{bash} +v.extract input=hospitals type=point where="FACILITY_N = 'FLAGSTAFF MEDICAL CENTER'" output=FMC +``` + +#### Python + +Use the v.extract tool to create a vector point map named *FMC* from the *hospitals* map: + +```{python} +gs.run_command("v.extract", input="hospitals", type="point", where="FACILITY_N = 'FLAGSTAFF MEDICAL CENTER'", output="FMC") +``` + +::: + +Display FMC with a color and size to see it better + +![](images/GRASS_movement_5.webp) + +### Compute the cumulative cost surface + +::: {.panel-tabset group="language"} + +#### GUI + +Use the *r.walk* tool from the **Raster/Terrain Analysis** menu: + +::::: grid +::: g-col-6 + +1. Enter the input and output parameters: **elevation**=*elevation*, **friction**=*friction0*, and **output**=*FMC_cost_seconds*: + +![](images/GRASS_movement_6.webp) +::: + +::: g-col-6 +2. Enter the name of a directions map **outdir**=*FMC_directions* for creating a least cost path: + +![](images/GRASS_movement_7.webp) + +::: +::::: + +::::: grid +::: g-col-6 +3. Enter the start point (**start_points**=*FMC*). + +![](images/GRASS_movement_8.webp){fig-align="left" width="100%"} +::: +::: g-col-6 +4. Optional: control the spatial extent of cost surface: + +![](images/GRASS_movement_9.webp) +::: +:::::: + +::::: grid +::: g-col-6 +5. Optional: adjust movement parameter settings: + +![](images/GRASS_movement_10.webp) +::: +::: g-col-6 +6. Recommended: select *knight's move* for more accurate cost and direction: + +![](images/GRASS_movement_11.webp) +::: +::::: + +::: {.callout-note title="Tip"} +Click the "copy" button to copy the GRASS command. +You can save it in a text file for later reuse or to document your work. +::: + +#### Command line + +Use the r.walk command to generate the cumulative cost surface: + +```{bash} +r.walk elevation=elevation friction=friction0 output=FMC_cost_seconds outdir=FMC_directions start_points=FMC -k +``` + +#### Python + +Use the r.walk command to generate the cumulative cost surface: + +```{python} +gs.run_command("r.walk", + elevation="elevation", + friction="friction0", + output="FMC_cost_seconds", + outdir="FMC_directions", + start_points="FMC", + flags="k") +``` + +::: + +Each raster cell value in the cost surface is the time in +seconds to walk from FMC to that cell given the terrain. + +![Displaying the cost surface over a shaded relief map using the [d.shade](https://grass.osgeo.org/grass-stable/manuals/d.shade.html) tool and a relief map of elevation made with [r.relief](https://grass.osgeo.org/grass-stable/manuals/d.shade.html).](images/GRASS_movement_12.webp) + +### Movement across a cumulative cost surface + +You can transform the walking time in seconds to hours by dividing the map by 3600 using the map calculator + +::: {.panel-tabset group="language"} + +#### GUI + +::::: grid +::: g-col-6 + +1. Open the map calculator. +2. Enter the name of the new map of walking time in hours: *FMC_cost_hours*. +3. Enter `FMC_cost_seconds / 3600` into the expressions field +4. Press Run. +::: +::: g-col-6 +![](images/GRASS_movement_13.webp) +::: +:::::: + +::::: grid +::: g-col-6 +We can query or filter the cost surface to areas of equivalent walking time. + +For example, to show the area within a 2 hour walk of FMC, filter the +cumulative cost surface display using the *d.rast* tool in the *layer manager* +to show 0-2 hours and set the opacity of the cost surface to 50% to see +the underlying terrain. +::: +::: g-col-6 +![](images/GRASS_movement_15.webp) +::: +::::: + +#### Command line + +```{shell} +r.mapcalc "FMC_cost_hours = FMC_cost_seconds / 3600" +``` + +#### Python + +```{python} +gs.mapcalc("FMC_cost_hours = FMC_cost_seconds / 3600") +``` + +::: + +![Shaded area represents terrain reached within a 2 hour walk from FMC](images/GRASS_movement_14.webp) + +### Least cost path + +We can also plot a *least cost path* (LCP), which is the least costly +(shortest time) route between any point on the cumulative cost surface and FMC + +Imagine a stranded hiker NE of Flagstaff who has to walk to FMC. +What path would take the least time? + +![Start and end point](images/GRASS_movement_16.webp) + +To create a LCP in GRASS, we will use [r.path](https://grass.osgeo.org/grass-stable/manuals/r.path.html) +(also under the Raster/Terrain analysis menu). + +::: {.panel-tabset group="language"} + +#### GUI + +::::: grid +::: g-col-6 + +1. We need to use the directions map **input**=*FMC_directions* +we made when we created the cost surface. +::: +::: g-col-6 +![](images/GRASS_movement_17.webp) +::: +:::::: + +::::: grid +::: g-col-6 +2. The coordinates of the hiker **start_coordinates**=*477476,13914951* +(we could also use a pre-defined vector point). +::: +::: g-col-6 +![](images/GRASS_movement_18.webp) +::: +::::: + +:::: grid +::: g-col-6 +3. And specify the name of the output LCP **vector_path**=*LCP_cumulative*. +::: +::: g-col-6 +![](images/GRASS_movement_19.webp) +::: +::::: + +#### Command line + +```{shell} +r.path input=FMC_directions format=auto vector_path=LCP_cumulative start_coordinates=477476,13914951 +``` + +#### Python + +```{python} +gs.run_command("r.path", + input="FMC_directions", + format="auto", + vector_path="LCP_cumulative", + start_coordinates=[477476, 13914951]) +``` + +::: + +![Resulting least cost path](images/GRASS_movement_20.webp) + +## Adding a friction map to movement costs + +A **friction map** can be used to incorporate other factors than +just terrain into creating a cost surface and modeling movement. +The value of each cell in a friction surface raster is the amount +of walking time, in seconds/meter, *in addition to* the time needed +because of terrain. + +We can reclassify the land use map (*landuse*) to create a friction surface +of the amount of extra time it would take to walk through different kinds of land cover. + +![Reclassification of *landuse* to show major land cover types](images/GRASS_movement_21.webp) + +### Reclassification of landuse to a friction map + +A standard walking speed across flat terrain is about 5 km/hr = 0.72 sec/m. +We might then estimate that it would take an additional: + +- 3 sec/m to walk through dense pinyon-juniper woodland, +- 1 sec/m to cross conifer forest, +- 2 sec/m to wander through urban Flagstaff, +- 5 sec/m to clamber over lava fields, and +- 10 sec/m to try and cross water + +We can create this friction map by reclassifying the *landuse* map using +[r.reclass](https://grass.osgeo.org/grass-stable/manuals/r.reclass.html) +to assign new friction values to the existing land cover categories. + +::: {.panel-tabset group="language"} + +## GUI + +The [r.reclass](https://grass.osgeo.org/grass-stable/manuals/r.reclass.html) tool +is found under the Raster/Change category… menu. + +::::: grid +::: g-col-6 + +1. Enter *landuse* for the raster to be reclassified + +2. Enter *friction_reclassified* for the output reclassified map + +3. Enter the reclass rules directly in the text box or from a saved text file. Use and \* symbol to represent everything not covered by the specific reclass rules + + > 11 90 95 = 10 water + > 21 thru 24 = 2 urban + > 31 = 5 lava + > 41 thru 43 = 1 conifer forest + > 52 = 3 pinyon juniper woodland + > \* = 0 no friction + +::: + +::: g-col-6 +This creates a new friction map named *friction_landcover*. +![](images/GRASS_movement_22.webp) +::: +::::: + +## Command line + +```{shell} +r.reclass input=landuse output=friction_landcover rules=- << EOF +11 90 95 = 10 water +21 thru 24 = 2 urban +31 = 5 lava +41 thru 43 = 1 conifer forest +52 = 3 pinyon juniper woodland +* = 0 no friction +EOF +``` + +## Python + +```{python} +rules = """\ +11 90 95 = 10 water +21 thru 24 = 2 urban +31 = 5 lava +41 thru 43 = 1 conifer forest +52 = 3 pinyon juniper woodland +* = 0 no friction +""" +gs.write_command("r.reclass", + input="landuse", + output="friction_landcover", + rules="-", + stdin=rules) +``` + +::: + +![The new friction map](images/GRASS_movement_23.webp) + +### Modifying a cumulative cost surface with a friction map + +We can now create a new cumulative cost surface using this new friction map +and convert it from seconds to hours, as we did before: + +::: {.panel-tabset group="language"} + +#### GUI + +1. Follow the procedures in the [Compute the cumulative cost surface +section](#compute-the-cumulative-cost-surface) and substitute the new +*friction_landcover* map for the *friction0* map used previously. + +2. Convert the cumulative cost surface to hours instead of seconds +following the procedures in the [Movement across a cumulative cost surface section](#movement-across-a-cumulative-cost-surface). + +#### Command line + +```{bash} +r.walk elevation=elevation friction=friction_landcover output=FMC_vegcost_seconds outdir=FMC_vegcost_directions start_points=FMC -k +r.mapcalc "FMC_vegcost_hours = FMC_vegcost_seconds / 3600" +``` + +#### Python + +```{python} +gs.run_command("r.walk", + elevation="elevation", + friction="friction_landcover", + output="FMC_vegcost_seconds", + outdir="FMC_vegcost_directions", + start_points="FMC", + flags="k") +gs.mapcalc("FMC_vegcost_hours = FMC_vegcost_seconds / 3600") +``` + +::: + +![Cumulative cost surface with landcover friction](images/GRASS_movement_24.webp) + +### Least cost paths and friction + +We can create a new LCP from the stranded hiker to FMC over terrain +where land cover also affects the cost of movement. + +::: {.panel-tabset group="language"} + +#### GUI + +1. Follow the procedures in the [Least cost path section](#least-cost-path), substituting the new direction map made along with the cumulative cost surface with land cover friction map. + +2. Give this new LCP the name *LCP_vegcost* to distinguish from the previous LCP. + +#### Command line + +```{bash} +r.path input=FMC_vegcost_directions format=auto vector_path=LCP_vegcost start_coordinates=477476,13914951 +``` + +#### Python + +```{python} +gs.run_command("r.path", + input="FMC_vegcost_directions", + format="auto", + vector_path="LCP_vegcost", + start_coordinates=[477476, 13914951]) +``` + +::: + +It takes more time to reach FMC if the hiker has to navigate dense vegetation as well as terrain. + +In the map below, terrain colored by land cover, the original LCP with terrain only is shown by the blue line. The LCP with a land cover friction map is shown by the heavier yellow line + +![Comparing LCP with terrain only and with land cover friction](images/GRASS_movement_25.webp) + +## Acknowledgement + +The development of this tutorial was funded by the U.S. +[National Science Foundation (NSF)](https://www.nsf.gov/), +award [2303651](https://www.nsf.gov/awardsearch/showAward?AWD_ID=2303651).