8.4. Lesson: Exercício Suplementar

En esta lección, serás guiado a través de un análisis SIG completo en QGIS.

Nota

Lesson developed by Linfiniti Consulting (South Africa) and Siddique Motala (Cape Peninsula University of Technology)

8.4.1. Planteamiento del Problema

You are tasked with finding areas in and around the Cape Peninsula that are suitable habitats for a rare fynbos plant species. The extent of your area of investigation covers Cape Town and the Cape Peninsula between Melkbosstrand in the north and Strand in the south. Botanists have provided you with the following preferences exhibited by the species in question:

  • It grows on east facing slopes

  • It grows on slopes with a gradient between 15% and 60%

  • It grows in areas that have a total annual rainfall of > 1000 mm

  • It will only be found at least 250 m away from any human settlement

  • The area of vegetation in which it occurs should be at least 6000 ㎡ in area

As a student at the University, you have agreed to search for the plant in four different suitable areas of land. You want those four suitable areas to be the ones that are closest to the University of Cape Town where you live. Use your GIS skills to determine where you should go to look.

8.4.2. Esquema de la Solución

The data for this exercise can be found in the exercise_data/more_analysis folder.

You are going to find the four suitable areas that are closest to the University of Cape Town.

The solution will involve:

  1. Analysing a DEM raster layer to find the east facing slopes and the slopes with the correct gradients

  2. Analysing a rainfall raster layer to find the areas with the correct amount of rainfall

  3. Analysing a zoning vector layer to find areas that are away from human settlement and are of the correct size

8.4.3. Follow Along: Setting up the Map

  1. Click on the projectionEnabled Current CRS button in the lower right corner of the screen. Under the CRS tab of the dialog that appears, use the “Filter” tool to search for “33S”. Select the entry WGS 84 / UTM zone 33S (with EPSG code 32733).

  2. Clique OK.

  3. Save the project file by clicking on the fileSave Save Project toolbar button, or use the File ▶ Save As… menu item.

    Save it in a new directory called Rasterprac, that you should create somewhere on your computer. You will save whatever layers you create in this directory as well. Save the project as your_name_fynbos.qgs.

8.4.4. Carregar Dados para o Mapa

Para procesar los datos, necesitarás cargar las capas necesarias (nombres de calles, zonas, precipitación, MDE) en la vista del mapa.

8.4.4.1. Para vetores…

  1. Click on the dataSourceManager Open Data Source Manager button in the Data Source Manager Toolbar, and enable the addOgrLayer Vector tab in the dialog that appears, or use the Layer ▶ Add Layer ▶ addOgrLayer Add Vector Layer… menu item

  2. Ensure that radioButtonOn File is selected

  3. Click on the … button to browse for vector dataset(s)

  4. In the dialog that appears, open the exercise_data/more_analysis/Streets directory

  5. Select the file Street_Names_UTM33S.shp

  6. Haz clic en Abrir.

    The dialog closes and shows the original dialog, with the file path specified in the text field next to Vector dataset(s). This allows you to ensure that the correct file is selected. It is also possible to enter the file path in this field manually, should you wish to do so.

  7. Click Add. The vector layer will be loaded into your map. Its color is automatically assigned. You will change it later.

  8. Rename the layer to Streets

    1. Right-click on it in the Layers panel (by default, the pane along the left-hand side of the screen)

    2. Click Rename in the dialog that appears and rename it, pressing the Enter key when done

  9. Repita o processo de adição de vetores, mas desta vez selecione o arquivo Generalised_Zoning_Dissolve_UTM33S.shp no diretório Zoning.

  10. Renomeie para Zoneamento.

8.4.4.2. Para rasters…

  1. Click on the dataSourceManager Open Data Source Manager button and enable the addRasterLayer Raster tab in the dialog that appears, or use the Layer ▶ Add Layer ▶ addRasterLayer Add Raster Layer… menu item

  2. Ensure that radioButtonOn File is selected

  3. Navigate to the appropriate file, select it, and click Open

  4. Do this for each of the following two raster files, DEM/SRTM.tif and rainfall/reprojected/rainfall.tif

  5. Rename the SRTM raster to DEM and the rainfall raster to Rainfall (with an initial capital)

8.4.5. Cambio de orden de capas

Click and drag layers up and down in the Layers panel to change the order they appear in on the map so that you can see as many of the layers as possible.

Now that all the data is loaded and properly visible, the analysis can begin. It is best if the clipping operation is done first. This is so that no processing power is wasted on computing values in areas that are not going to be used anyway.

8.4.6. Encuentra los Distritos Correctos

  1. Carregue a camada vetorial admin_boundaries/Western_Cape_UTM33S.shp no seu mapa.

  2. Rename it to Districts.

  3. Right-click on the Districts layer in the Layers panel.

  4. In the menu that appears, select the Filter… menu item. The Query Builder dialog appears.

    You will now build a query to select only the following districts:

    • Bellville

    • Cabo

    • Goodwood

    • Kuils River

    • Mitchells Plain

    • Simon Town

    • Wynberg

    1. In the Fields list, double-click on the NAME_2 field to make it appear in the SQL where clause text field below

    2. Click the = button; an = sign is appended to the SQL query

    3. Click the All button below the (currently empty) Values list.

      After a short delay, this will populate the Values list with the values of the selected field (NAME_2).

    4. Double-click the value Bellville in the Values list to append it to the SQL query.

      In order to select more than one district, you’ll need to use the OR boolean operator.

    5. Click the OR button and it will be appended to the SQL query

    6. Using a process similar to the above, add the following to the existing SQL query:

      "NAME_2" = 'Cape'

    7. Add another OR operator, then work your way through the list of districts above in a similar fashion.

      The final query should be:

      "NAME_2" = 'Bellville' OR "NAME_2" = 'Cape' OR
"NAME_2" = 'Goodwood' OR "NAME_2" = 'Kuils River' OR
"NAME_2" = 'Mitchells Plain' OR "NAME_2" = 'Simon Town' OR
"NAME_2" = 'Wynberg'

      Nota

      By using the IN operator, the query would look like this:

      "NAME_2" in ('Bellville', 'Cape', 'Goodwood', 'Kuils River',
             'Mitchells Plain', 'Simon Town', 'Wynberg')

    8. Clique em OK.

      The districts shown in your map are now limited to those in the list above.

8.4.7. Recorta los Ráster

Ahora que tienes un área de interés, puedes recortar los ráster a esa área.

  1. Ensure that the only layers that are visible are the DEM, Rainfall and Districts layers

  2. Districts must be on top to be visible

  3. Open the clipping dialog by selecting the menu item Raster ▶ Extraction ▶ Clip Raster by Mask Layer…

  4. In the Input layer dropdown list, select the DEM layer

  5. In the Mask layer dropdown list, select the Districts layer

  6. Scroll down and specify an output location in the Clipped (mask) text field by clicking the … button and choosing Save to File…

    1. Navigate to the Rasterprac directory

    2. Enter a file name - DEM_clipped.tif

    3. Save

  7. Make sure that checkbox Open output file after running algorithm is checked

  8. Click Run

    After the clipping operation has completed, leave the Clip Raster by Mask Layer dialog open, to be able to reuse the clipping area

  9. Select the Rainfall raster layer in the Input layer dropdown list and save your output as Rainfall_clipped.tif

  10. Do not change any other options. Leave everything the same and click Run.

  11. After the second clipping operation has completed, you may close the Clip Raster by Mask Layer dialog

  12. Save the map

8.4.7.1. Align the rasters

For our analysis we need the rasters to have the same CRS and they have to be aligned.

First we change the resolution of our rainfall data to 30 meters (pixel size):

  1. Right-click on the Rainfall_clipped layer and select Export▶ Save As… in the context menu.

  2. Set the Horizontal and Vertical resolution to 30 (meters).

  3. Save the file as Rainfall30.tif in rainfall/reprojected (File name)`

Then we align the DEM:

  1. Right-click on the DEM_clipped layer and select Export▶ Save As… in the context menu

  2. For CRS, choose WGS 84 / UTM zone 33S (EPSG code 32733)

  3. Set the Horizontal and Vertical resolution to 30 (meters)

  4. Under Extent, click on Calculate from Layer and choose rainfall30

  5. Save the file as DEM30.tif in DEM/reprojected (File name)

Para ver correctamente qué está pasando, se necesita cambiar la simbología para las capas.

8.4.8. Cambio de simbología de capas vectoriales

  1. In the Layers panel, right-click on the Streets layer

  2. Select Properties from the menu that appears

  3. Switch to the Symbology tab in the dialog that appears

  4. Click on the Color dropdown

  5. Select a new color in the dialog that appears

  6. Clique OK.

  7. Vuelve a hacer clic en Aceptar en el cuadro de diálogo Propiedades de la capa. Esto cambiará el color de la capa Streets.

  8. Follow a similar process for the Zoning layer and choose an appropriate color for it

8.4.9. Cambio de simbología de capas ráster

La simbología de capas ráster es algo diferente.

  1. Open the Properties dialog for the Rainfall30 raster layer

  2. Alterne para a guia Simbologia. Você vai notarque esse diálogo é muito diferente da versão usada para as camadas de vetor.

  3. Expand Min/Max Value Settings

  4. Ensure that the button Mean +/- standard deviation is selected

  5. Make sure that the value in the associated box is 2.00

  6. For Contrast enhancement, make sure it says Stretch to MinMax

  7. For Color gradient, change it to White to Black

  8. Clique OK.

    The Rainfall30 raster, if visible, should change colors, allowing you to see different brightness values for each pixel

  9. Repeat this process for the DEM30 layer, but set the standard deviations used for stretching to 4.00

8.4.10. Limpia el mapa

  1. Remove the original Rainfall and DEM layers, as well as Rainfall_clipped and DEM_clipped from the Layers panel:

    • Haz clic derecho en esas capas y selecciona Eliminar.

      Nota

      Esto no borrará los datos de tu dispositivo de almacenamiento, solamente lo quitará de tu mapa.

  2. Save the map

  3. Agora você pode ocultar as camadas vetoriais desmarcando a caixa ao lado delas no painel Camadas. Isso fará o mapa renderizar mais rápido e economizará algum tempo.

8.4.11. Crear el sombreado del relieve

Para criar o sombreamento, você precisará usar um algoritmo que foi escrito para esse fim.

  1. In the Layers panel, ensure that DEM30 is the active layer (i.e., it is highlighted by having been clicked on)

  2. Click on the Raster ▶ Analysis ▶ Hillshade… menu item to open the Hillshade dialog

  3. Scroll down to Hillshade and save the output in your Rasterprac directory as hillshade.tif

  4. Make sure that checkbox Open output file after running algorithm is checked

  5. Click Run

  6. Aguarde pelo final do processamento.

The new hillshade layer has appeared in the Layers panel.

  1. Right-click on the hillshade layer in the Layers panel and bring up the Properties dialog

  2. Click on the Transparency tab and set the Global Opacity slider to 20%

  3. Clique OK.

  4. Note the effect when the transparent hillshade is superimposed over the clipped DEM. You may have to change the order of your layers, or click off the rainfall30 layer in order to see the effect.

8.4.12. Inclinação

  1. Click on the Raster ▶ Analysis ▶ Slope… menu item to open the Slope algorithm dialog

  2. Select DEM30 as Input layer

  3. Check checkbox Slope expressed as percent instead of degrees. Slope can be expressed in different units (percent or degrees). Our criteria suggest that the plant of interest grows on slopes with a gradient between 15% and 60%. So we need to make sure our slope data is expressed as a percent.

  4. Specify an appropriate file name and location for your output.

  5. Make sure that checkbox Open output file after running algorithm is checked

  6. Click Run

The slope image has been calculated and added to the map. As usual, it is rendered in grayscale. Change the symbology to a more colorful one:

  1. Open the layer Properties dialog (as usual, via the right-click menu of the layer)

  2. Click on the Symbology tab

  3. Where it says Singleband gray (in the Render type dropdown menu), change it to Singleband pseudocolor

  4. Choose Mean +/- standard deviation x for Min / Max Value Settings with a value of 2.0

  5. Select a suitable Color ramp

  6. Click Run

8.4.13. Try Yourself Aspect

Use the same approach as for calculating the slope, choosing Aspect… in the Raster ▶ Analysis menu.

Remember to save the project periodically.

8.4.14. Reclassificar rasters

  1. Choose Raster ▶ Raster calculator…

  2. Specify your Rasterprac directory as the location for the Output layer (click on the … button), and save it as slope15_60.tif

  3. Verifique se a caixa :guilabel:`Abrir arquivo de saída após executar o algoritmo’ está selecionada.

    In the Raster bands list on the left, you will see all the raster layers in your Layers panel. If your Slope layer is called slope, it will be listed as slope@1. Indicating band 1 of the slope raster.

  4. The slope needs to be between 15 and 60 degrees.

    Using the list items and buttons in the interface, build the following expression:

    (slope@1 > 15) AND (slope@1 < 60)

  5. Ajusta el campo Capa de salida a un nombre y localización adecuados.

  6. Haz clic en Run.

Agora encontre o aspecto correto (voltado para o leste: entre 45 e 135 graus) usando a mesma abordagem.

  1. Build the following expression:

    (aspect@1 > 45) AND (aspect@1 < 135)

You will know it worked when all of the east-facing slopes are white in the resulting raster (it’s almost as if they are being lit by the morning sunlight).

Find the correct rainfall (greater than 1000 mm) the same way. Use the following expression:

rainfall30@1 > 1000

Now that you have all three criteria each in separate rasters, you need to combine them to see which areas satisfy all the criteria. To do so, the rasters will be multiplied with each other. When this happens, all overlapping pixels with a value of 1 will retain the value of 1 (i.e. the location meets the criteria), but if a pixel in any of the three rasters has the value of 0 (i.e. the location does not meet the criteria), then it will be 0 in the result. In this way, the result will contain only the overlapping areas that meet all of the appropriate criteria.

8.4.15. Combinación de rásters

  1. Open the Raster Calculator (Raster ▶ Raster Calculator…)

  2. Build the following expression (with the appropriate names for your layers):

    [aspect45_135] * [slope15_60] * [rainfall_1000]

  3. Set the output location to the Rasterprac directory

  4. Name the output raster aspect_slope_rainfall.tif

  5. Ensure that checkbox Open output file after running algorithm is checked

  6. Click Run

The new raster now properly displays the areas where all three criteria are satisfied.

Save the project.

The next criterion that needs to be satisfied is that the area must be 250 m away from urban areas. We will satisfy this requirement by ensuring that the areas we compute are inside rural areas, and are 250 m or more from the edge of the area. Hence, we need to find all rural areas first.

8.4.16. Encontrar áreas rurales

  1. Hide all layers in the Layers panel

  2. Unhide the Zoning vector layer

  3. Right-click on it and bring up the Attribute Table dialog. Note the many different ways that the land is zoned here. We want to isolate the rural areas. Close the Attribute table.

  4. Right-click on the Zoning layer and select Filter… to bring up the Query Builder dialog

  5. Build the following query:

    "Gen_Zoning" = 'Rural'

    See the earlier instructions if you get stuck.

  6. Click OK to close the Query Builder dialog. The query should return one feature.

You should see the rural polygons from the Zoning layer. You will need to save these.

  1. In the right-click menu for Zoning, select Export ▶ Save Features As….

  2. Save your layer under the Rasterprac directory

  3. Name the output file rural.shp

  4. Clique OK.

  5. Save the project

Agora você precisa excluir as áreas que estão dentro de 250m da borda das áreas rurais. Faça isso criando um buffer negativo, conforme explicado abaixo.

8.4.17. Crear un buffer negativo

  1. Click the menu item Vector ▶ Geoprocessing Tools ▶ Buffer…

  2. In the dialog that appears, select the rural layer as your input vector layer (Selected features only should not be checked)

  3. Set Distance to -250. The negative value means that the buffer will be an internal buffer. Make sure that the units are meters in the dropdown menu.

  4. Check checkbox Dissolve result

  5. In Buffered, place the output file in the Rasterprac directory, and name it rural_buffer.shp

  6. Click Save

  7. Click Run and wait for the processing to complete

  8. Feche a caixa de diálogo Buffer.

    Make sure that your buffer worked correctly by noting how the rural_buffer layer is different from the rural layer. You may need to change the drawing order in order to observe the difference.

  9. Remove the rural layer

  10. Save the project

Now you need to combine your rural_buffer vector layer with the aspect_slope_rainfall raster. To combine them, we will need to change the data format of one of the layers. In this case, you will vectorize the raster, since vector layers are more convenient when we want to calculate areas.

8.4.18. Vectorizar el ráster

  1. Click on the menu item Raster ▶ Conversion ▶ Polygonize (Raster to Vector)…

  2. Select the aspect_slope_rainfall raster as Input layer

  3. Set Name of the field to create to suitable (the default field name is DN - Digital number data)

  4. Save the output. Under Vectorized, select Save file as. Set the location to Rasterprac and name the file aspect_slope_rainfall_all.shp.

  5. Ensure that checkbox Open output file after running algorithm is checked

  6. Click Run

  7. Close the dialog when processing is complete

All areas of the raster have been vectorized, so you need to select only the areas that have a value of 1 in the suitable field. (Digital Number.

  1. Open the Query Builder dialog (right-click - Filter…) for the new vector layer

  2. Build this query:

    "suitable" = 1

  3. Clique OK.

  4. After you are sure the query is complete (and only the areas that meet all three criteria, i.e. with a value of 1 are visible), create a new vector file from the results, using the Export –> Save Features As… in the layer’s right-click menu

  5. Save the file in the Rasterprac directory

  6. Name the file aspect_slope_rainfall_1.shp

  7. Remove the aspect_slope_rainfall_all layer from your map

  8. Save your project

When we use an algorithm to vectorize a raster, sometimes the algorithm yields what is called “Invalid geometries”, i.e. there are empty polygons, or polygons with mistakes in them, that will be difficult to analyze in the future. So, we need to use the “Fix Geometry” tool.

8.4.19. Fixing geometry

  1. In the Processing Toolbox, search for “Fix geometries”, and Execute… it

  2. For the Input layer, select aspect_slope_rainfall_1

  3. Under Fixed geometries, select Save file as, and save the output to Rasterprac and name the file fixed_aspect_slope_rainfall.shp.

  4. Ensure that checkbox Open output file after running algorithm is checked

  5. Click Run

  6. Close the dialog when processing is complete

Now that you have vectorized the raster, and fixed the resulting geometry, you can combine the aspect, slope, and rainfall criteria with the distance from human settlement criteria by finding the intersection of the fixed_aspect_slope_rainfall layer and the rural_buffer layer.

8.4.20. Determining the Intersection of vectors

  1. Click the menu item Vector ▶ Geoprocessing Tools ▶ Intersection…

  2. In the dialog that appears, select the rural_buffer layer as Input layer

  3. For the Overlay layer, select the fixed_aspect_slope_rainfall layer

  4. In Intersection, place the output file in the Rasterprac directory

  5. Name the output file rural_aspect_slope_rainfall.shp

  6. Click Save

  7. Click Run and wait for the processing to complete

  8. Close the Intersection dialog.

    Make sure that your intersection worked correctly by noting that only the overlapping areas remain.

  9. Save the project

The next criteria on the list is that the area must be greater than 6000 ㎡. You will now calculate the polygon areas in order to identify the areas that are the appropriate size for this project.

8.4.21. Cálgulo del área para cada polígono

  1. Open the new vector layer’s right-click menu

  2. Select Open attribute table

  3. Click the toggleEditing Toggle editing button in the top left corner of the table, or press Ctrl+e

  4. Click the calculateField Open field calculator button in the toolbar along the top of the table, or press Ctrl+i

  5. In the dialog that appears, make sure that checkbox Create new field is checked, and set the Output field name to area The output field type should be a decimal number (real). Set Precision to 1 (one decimal).

  6. In the Expression area, type:

    $area

    This means that the field calculator will calculate the area of each polygon in the vector layer and will then populate a new integer column (called area) with the computed value.

  7. Clique OK.

  8. Do the same thing for another new field called id. In Field calculator expression, type:

    $id

    Eso asegura que cada polígono tiene una ID única para su identificación.

  9. Click toggleEditing Toggle editing again, and save your edits if prompted to do so

8.4.22. Selección de áreas para un tamaño dado

Ahora que las áreas son conocidas:

  1. Build a query (as usual) to select only the polygons that are larger than 6000 ㎡. The query is:

    "area" > 6000

  2. Save the selection in the Rasterprac directory as a new vector layer called suitable_areas.shp.

You now have the suitable areas that meet all of the habitat criteria for the rare fynbos plant, from which you will pick the four areas that are nearest to the University of Cape Town.

8.4.23. Digitize the University of Cape Town

  1. Create a new vector layer in the Rasterprac directory as before, but this time, use Point as Geometry type and name it university.shp

  2. Ensure that it is in the correct CRS (Project CRS:EPSG:32733 - WGS 84 / UTM zone 33S)

  3. Finish creating the new layer (click OK)

  4. Hide all layers except the new university layer and the Streets layer.

  5. Add a background map (OSM):

    1. Go to the Browser panel and navigate to XYZ Tiles ▶ OpenStreetMap

    2. Drag and drop the OpenStreetMap entry to the bottom of the Layers panel

    Using your internet browser, look up the location of the University of Cape Town. Given Cape Town’s unique topography, the university is in a very recognizable location. Before you return to QGIS, take note of where the university is located, and what is nearby.

  6. Ensure that the Streets layer clicked on, and that the university layer is highlighted in the Layers panel

  7. Navigate to the View ▶ Toolbars menu item and ensure that Digitizing is selected. You should then see a toolbar icon with a pencil on it (toggleEditing Toggle editing). This is the Toggle Editing button.

  8. Click the Toggle editing button to enter edit mode. This allows you to edit a vector layer

  9. Click the capturePoint Add Point Feature button, which should be nearby the toggleEditing Toggle editing button

  10. With the Add feature tool activated, left-click on your best estimate of the location of the University of Cape Town

  11. Supply an arbitrary integer when asked for the id

  12. Clique OK.

  13. Click the saveEdits Save Layer Edits button

  14. Click the Toggle editing button to stop your editing session

  15. Save the project

8.4.24. Find the locations that are closest to the University of Cape Town

  1. Go to the Processing Toolbox, locate the Join Attributes by Nearest algoritm (Vector general ▶ Join Attributes by Nearest) and execute it

  2. Input layer should be university, and Input layer 2 suitable_areas

  3. Set an appropriate output location and name (Joined layer)

  4. Set the Maximum nearest neighbors to 4

  5. Ensure that checkbox Open output file after running algorithm is checked

  6. Leave the rest of the parameters with their default values

  7. Click Run

The resulting point layer will contain four features - they will all have the location of the university and its attributes, and in addition, the attributes of the nearby suitable areas (including the id), and the distance to that location.

  1. Open the attribute table of the result of the join

  2. Note the id of the four nearest suitable areas, and then close the attribute table

  3. Open the attribute table of the suitable_areas layer

  4. Build a query to select the four suitable areas closest to the university (selecting them using the id field)

Esta es la respuesta final a la pregunta investigada.

For your submission, create a fully labeled layout that includes the semi-transparent hillshade layer over an appealing raster of your choice (such as the DEM or the slope raster, for example). Also include the university and the suitable_areas layer, with the four suitable areas that are closest to the university highlighted. Follow all the best practices for cartography in creating your output map.