9.4. Lesson: 補足の実習¶

9.4.1. 問題の状態¶

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

• それは、東向き斜面に生えています。

• It grows on slopes with a gradient between 15% and 60%.
• It grows in areas that have a total annual rainfall of > 1200 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 6000m2 in area.

As a volunteer for Cape Nature, you have agreed to search for the plant on the closest suitable piece of land to your house. Use your GIS skills to determine where you should go to look.

9.4.2. ソリューションの概要¶

In order to solve this problem, you will have to use the available data (available in exercise_data/more_analysis) to find the candidate area that is closest to your house. If you don’t live in Cape Town (where this problem is based) you can choose any house in the Cape Town region. The solution will involve:

• analysing the DEM to find the east facing slopes and the slopes with the correct gradients;
• analysing the rainfall raster to find the areas with the correct amount of rainfall;
• analysing the Zoning vector layer to find areas that are away from human settlement and are of the correct size.

9.4.3. マップの設定¶

• Click on the “CRS status” button in the extreme lower right corner of the screen. Under the CRS tab of the screen that appears, you will see the box Coordinate reference systems of the world.
• In this box, navigate to Projected Coordinate Systems ‣ Universal Transverse Mercator (UTM).
• Select the entry WGS 84 / UTM zone 33S (with the EPSG code 32733).

• OK をクリックします。マップは UTM33S のCRSです。

• Save the map by clicking on the Save Project As toolbar button, or use the File ‣ Save Project As... menu item.
• Save the map in a directory called Rasterprac that you should create somewhere on your computer. You will save whatever layers you create in this directory as well.

9.4.4. マップへの地図の読み込み¶

In order to process the data, you will need to load the necessary layers (street names, zones, rainfall, DEM) into the map canvas.

9.4.5. ベクタレイヤのシンボロジの変更¶

• In the Layers list, right-click on the Streets layer.
• Select Properties from the menu that appears.

• 表示されているダイアログで スタイル タブに切り替えます。

• Click on the button labeled Change, with a square showing the current color of the Streets layer.

• 表示されるダイアログで新しい色を選択します。

• Click OK.
• Click OK again in the Layer Properties dialog. This will change the color of the Streets layer.
• Follow a similar process for the Zoning layer and choose an appropriate color for it.

9.4.6. ラスタレイヤのシンボロジの変更¶

ラスタレイヤシンボロジはどこか違います。

• Open the Properties dialog for the Rainfall raster.
• Switch to the Style tab. You’ll notice that this style dialog is very different from the version used for vector layers.
• Ensure that the button Use standard deviation is selected.
• Change the value in the associated box to 2.00 (it should be set to 0.00 by default).
• Under the heading Contrast enhancement, change the value of the Current dropdown list to Stretch to MinMax.
• Click OK. The “Rainfall” raster, if visible, should change colors, allowing you to see different brightness values for each pixel.
• Repeat this process for the DEM, but set the standard deviations used for stretching to 4.00.

9.4.7. レイヤ順の変更¶

• In the Layers list, click and drag layers up and down to change the order they appear in on the map.
• Newer versions of QGIS may have a Control rendering order checkbox beneath the Layers list. Ensure that it is checked.

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 aren’t going to be used anyway.

9.4.8. 正しい地区の検索¶

• Load the vector layer admin_boundaries/Western_Cape_UTM33S.shp into your map.

• Districts にリネームします。

• Right-click on the Districts layer in the Layers list.
• In the menu that appears, select the Query... menu item. The Query Builder dialog appears.

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

• Bellville,
• Cape,
• Goodwood,
• Kuils River,
• Mitchells Plain,

• Simons Town, そして

• Wynberg.
• In the Fields list, double-click on the NAME_2 field. It appears in the SQL where clause text field below.
• Click the = button; an = sign is added to the SQL query.
• 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).
• Double-click the value Bellville in the Values list. As before, this will be added to the SQL query.

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

• Click the OR button and it will be added to the SQL query.

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

  "NAME_2" = 'Cape'
  

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

• 最終的なクエリはこうであるべきです

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

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

9.4.9. ラスタのクリップ¶

Now that you have an area of interest, you can clip the rasters to this area.

• Ensure that the only layers that are visible are the DEM, Rainfall and Districts layers.
• Districts must be on top so that they are visible.
• Open the clipping dialog by selecting the menu item Raster ‣ Extraction ‣ Clipper.
• In the Input file (raster) dropdown list, select the DEM layer.
• Specify an output location in the Output file text field by clicking the Select... button.

• Rasterprac ディレクトリに移動します。

• ファイル名を入力します。

• Save the file. Leave the No data value checkbox unchecked.
• Use the Extent clipping mode by ensuring the correct radio button is selected.
• Click and drag an area in the canvas, so that the area which includes the districts is selected.
• Check the Load into canvas when finished box.
• Click OK.
• After the clipping operation is completed, DO NOT CLOSE the Clipper dialog. (Doing so would cause you to lose the clipping area that you have already defined.)
• Select the Rainfall raster in the Input file (raster) dropdown list and choose a different output file name.
• Do not change any other options. Do not alter the existing clipping area which you drew previously. Leave everything the same and click OK.
• After the second clipping operation has completed, you may close the Clipper dialog.

• マップを保存します。

9.4.10. マップをクリーンアップします。¶

• Remove the original Rainfall and DEM layers from the Layers list:

• これらのレイヤ上で右クリックし、 削除 を選択します。

  • This will not remove the data from your storage device, it will merely take it out of your map.

• Streets レイヤのラベリングを無効にします。

  • ラベリング ボタンをクリックします。

  • Uncheck the Label this layer with box.
  • Click OK.

• 再度 Streets をすべて表示します。

  • Right-click on the layer in the Layers list.

  • クエリ を選択します。

• In the Query dialog that appears, click the Clear button, then click OK.
• Wait while the data is loaded. All the streets will now be visible.
• Change the raster symbology as before (see Changing the symbology of raster layers).

• マップを保存します。

• You can now hide the vector layers by unchecking the box next to them in the Layers list. This will make the map render faster and will save you some time.

In order to create the hillshade, you will need to use a plugin that was written for this purpose.

9.4.11. Activating the Raster Terrain Analysis plugin¶

This plugin is included by default in QGIS 1.8. However, it may not be immediately visible. To check if it is accessible on your system:

• Click on the menu item Plugins –> Manage Plugins....
• Ensure that the box next to Raster Terrain Analysis plugin is selected.
• Click OK.

You will now have access to this plugin via the Raster ‣ Terrain analysis menu item.

Remember that plugins may sometimes depend on certain Python modules being installed on your system. Should a plugin refuse to work while complaining of missing dependencies, please ask your tutor or lecturer for assistance.

9.4.12. 陰影図の作成¶

• In the Layers list, ensure that the DEM is the active layer (i.e., it is highlighted by having been clicked on).
• Click on the Raster ‣ Terrain analysis ‣ Hillshade menu item to open the Hillshade dialog.
• Specify an appropriate location for the output layer and call it hillshade.

• 結果をプロジェクトに追加する をチェックします。

• Click OK.

• 処理が完了するのを待ちます。

The new hillshade layer has appeared in your Layers list.

• Right-click on the hillshade layer in your Layers list and bring up the Properties dialog.
• Click on the Transparency tab and set the transparency slider to 80%.

• ダイアログで OK をクリックします。

• Note the effect when the transparent hillshade is superimposed over the clipped DEM.

9.4.13. 傾斜¶

• Click on the menu item Raster ‣ Terrain analysis.
• Select the Slope analysis type, with the clipped DEM as the input layer.
• Specify an appropriate file name and location for output purposes.

• 結果をプロジェクトに追加する をチェックします。

• Click OK.

The slope image has been calculated and added to the map. However, as usual it is just a gray rectangle. To properly see what’s going on, change the symbology as follows.

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

• スタイル タブをクリックします。

• Where it says Grayscale (in the Color map dropdown menu), change it to Pseudocolor.
• Ensure that the Use standard deviation radio button is selected.

9.4.14. 傾斜方位¶

• Use the same approach as for calculating the slope, but select Aspect in the initial dialog box.

定期的にマップを保存することを忘れないで下さい。

9.4.15. ラスタの再階級¶

• Click the menu item Raster ‣ Raster calculator.
• Specify your Rasterprac directory as the location for the output layer.
• Ensure that the Add result to project box is selected.

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

The slope needs to be between 15 and 60 degrees. Everything less than 15 or greater than 60 must therefore be excluded.

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

  ((slope@1 < 15) OR (slope@1 > 60)) = 0

• Set the Output layer field to an appropriate location and file name.

• Click OK.

Now find the correct aspect (east-facing: between 45 and 135 degrees) using the same approach.

• 次の式を組み立てます:

  ((aspect@1 < 45) OR (aspect@1 > 135)) = 0

• Find the correct rainfall (greater than 1200mm) the same way. Build the following expression:

  (rainfall@1 < 1200) = 0

Having reclassified all the rasters, you will now see them displayed as gray rectangles in your map (assuming that they have been added to the map correctly). To properly display raster data with only two classes (1 and 0, meaning true or false), you will need to change their symbology.

9.4.16. 再階級されたレイヤのスタイルの設定¶

• Open the Style tab in the layer’s Properties dialog as usual.
• Under the heading Load min / max values from band, select the Actual (slower) radio button.

• 読み込み ボタンをクリックします:

The Custom min / max values fields should now populate with 0 and 1, respectively. (If they do not, then there was a mistake with your reclassification of the data, and you will need to go over that part again.)

• Under the heading Contrast enhancement, set the Current dropdown list to Stretch To MinMax.
• Click OK.
• Do this for all three reclassified rasters, and remember to save your work!

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

9.4.17. 農村地域の検索¶

• Hide all layers in the Layers list.

• Unhide the Zoning vector layer.

• Right-click on it and bring up the Query dialog.

• 次のクエリをビルドします:

  "Gen_Zoning" = 'Rural'
  

  See the earlier instructions for building the Streets query if you get stuck.

• それを実行したら :guilabel:`クエリ`ダイアログを閉じます。

You should see a collection of polygons from the Zoning layer. You will need to save these to a new layer file.

• On the right-click menu for Zoning, select Save as....
• Save your layer under the Zoning directory.

• 出力ファイル名を rural.shp とします。

• Click OK.

• レイヤをマップに追加します。

• Click the menu item Vector ‣ Geoprocessing Tools ‣ Dissolve.
• Select the rural layer as your input vector layer, while leaving the Use only selected features box unchecked.
• Under Dissolve field, select — Dissolve all —.
• Save your layer under the Zoning directory.
• Click OK. A dialog will appear asking whether you want to add the new layer to the TOC (“Table of Contents”, referring to the Layers list).

• はい をクリックします。

• ディゾルブ ダイアログを閉じます。

• Remove the rural and Zoning layers.

• マップを保存します。

Now you need to exclude the areas that are within 250m from the edge of the rural areas. Do this by creating a negative buffer, as explained below.

9.4.18. Creating a negative buffer¶

• Click the menu item Vector ‣ Geoprocessing Tools ‣ Buffer(s).
• In the dialog that appears, select the rural_dissolve layer as your input vector layer (Use only selected features should not be checked).
• Select the Buffer distance button and enter the value -250 into the associated field; the negative value means that the buffer must be an internal buffer.

• バッファ結果をディゾルブする ボックスにチェックします。

• Set the output file to the same directory as the other rural vector files.

• 出力ファイル名を :kbd:`rural_buffer.shp`とします。

• Click Save.

• OK をクリックして処理が完了するのを待ちます。

• ダイアログに現れる:guilabel:Yes を選択します。

• Close the Buffer dialog.
• Remove the rural_dissolve layer.

• マップを保存します。

In order to incorporate the rural zones into the same analysis with the three existing rasters, it will need to be rasterized as well. But in order for the rasters to be compatible for analysis, they will need to be the same size. Therefore, before you can rasterize, you’ll need to clip the vector to the same area as the three rasters. A vector can only be clipped by another vector, so you will first need to create a bounding box polygon the same size as the rasters.

9.4.19. バウンディングボックスベクタの作成¶

• Click on the menu item Layer –> New –> New Shapefile Layer....
• Under the Type heading, select the Polygon button.
• Click Specify CRS and set the coordinate reference system WGS 84 / UTM zone 33S : EPSG:32733.

• OKをクリックします。

• Click OK on the New Vector Layer dialog as well.
• Save the vector in the Zoning directory.
• Name the output file bbox.shp.
• Hide all layers except the new bbox layer and one of the reclassified rasters.
• Ensure that the bbox layer is highlighted in the Layers list.
• Navigate to the View > Toolbars menu item and ensure that Digitizing is selected. You should then see a toolbar icon with a pencil or koki on it. This is the Toggle editing button.
• Click the Toggle editing button to enter edit mode. This allows you to edit a vector layer.
• Click the Add feature button, which should be nearby the Toggle editing button. It may be hidden behind a double arrow button; if so, click the double arrows to show the Digitizing toolbar’s hidden buttons.
• With the Add feature tool activated, left-click on the corners of the raster. You may need to zoom in with the mouse wheel to ensure that it is accurate. To pan across the map in this mode, click and drag in the map with the middle mouse button or mouse wheel.
• For the fourth and final point, right-click to finalize the shape.

• シェープIDの任意の番号を入力します。

• Click OK.

• 編集を保存 ボタンをクリックします。

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

• マップを保存します。

Now that you have a bounding box, you can use it to clip the rural buffer layer.

9.4.20. ベクタレイヤのクリップ¶

• Ensure that only the bbox and rural_buffer layers are visible, with the latter on top.
• Click the menu item Vector > Geoprocessing Tools > Clip.
• In the dialog that appears, set the input vector layer to rural_buffer and the clip layer to bbox, with both Use only selected features boxes unchecked.
• Put the output file under the Zoning directory.

• 出力ファイル名を file :kbd:`rural_clipped`とします。

• Click OK.
• When prompted to add the layer to the TOC, click Yes.

• ダイアログを閉じます。

• Compare the three vectors and see the results for yourself.
• Remove the bbox and rural_buffer layers, then save your map.

これでラスタライズする準備ができました。

9.4.21. ベクタレイヤのラスタ化¶

You’ll need to specify a pixel size for a new raster that you create, so first you’ll need to know the size of one of your existing rasters.

• Open the Properties dialog of any of the three existing rasters.

• :guilabel:`メタデータ`タブに切り替えます。

• Make a note of the X and Y values under the heading Dimensions in the Metadata table.

• プロパティ ダイアログを閉じます。

• Click on the Raster ‣ Conversion ‣ Rasterize menu item. You may receive a warning about a dataset being unsupported. Click it away and ignore it.

• 入力レイヤとして:guilabel:`rural_clipped`を選択します。

• Set an output file location inside the Zoning directory.

• 出力ファイル名を :kbd:`rural_raster.tif``とします。

• Check the New size box and enter the X and Y values you made a note of earlier.

• キャンバスにロードします ボタンをチェックします。

• Click the pencil icon next to the text field which shows the command that will be run. At the end of the existing text, add a space and then the text -burn 1. This tells the Rasterize function to “burn” the existing vector into the new raster and give the areas covered by the vector the new value of 1 (as opposed to the rest of the image, which will automatically be 0).
• Click OK.
• The new raster should show up in your map once it has been computed.
• The new raster will look like a grey rectangle – you may change the display style as you did for the reclassified rasters.

• マップを保存します。

Now that you have all four criteria each in a separate raster, 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, but if a pixel has the value of 0 in any of the four rasters, then it will be 0 in the result. In this way, the result will contain only the overlapping areas.

9.4.22. ラスタの結合¶

• Click the Raster ‣ Raster calculator menu item.

• Build the following expression (with the appropriate names for your layers, depending on what you called them):

  [Rural raster] * [Reclassified aspect] * [Reclassified slope] *
  [Reclassified rainfall]

• 出力場所を :kbd:`Rasterprac`ディレクトリに設定します。

• Name the output raster cross_product.tif.

• 結果をプロジェクトに追加する ボックスがチェックされていることを確認します。

• OKをクリックします。

• Change the symbology of the new raster in the same way as you set the style for the other reclassified rasters. The new raster now properly displays the areas where all the criteria are satisfied.

To get the final result, you need to select the areas that are greater than 6000m^2. However, computing these areas accurately is only possible for a vector layer, so you will need to vectorize the raster.

9.4.23. ラスタのベクタライズ¶

• Click on the menu item Raster ‣ Conversion ‣ Polygonize.
• Select the cross_product raster.

• 出力場所を:kbd:`Rasterprac`に設定します。

• :kbd:`candidate_areas.shp`と改名します。

• 終了時にキャンバスにロードします ボタンがチェックされていることを確認します。

• OKをクリックします。

• 処理が完了したらダイアログを閉じます。

All areas of the raster have been vectorized, so you need to select only the areas that have a value of 1.

• Open the Query dialog for the new vector.

• クエリのビルド:

  "DN" = 1
  

• Click OK.

• Create a new vector file from the results by saving the candidate_areas vector after the query is complete (and only the areas with a value of 1 are visible). Use the Save as... function in the layer’s right-click menu for this.

• :kbd:`Rasterprac`ディレクトリにファイルを保存します。

• :guilabel:`candidate_areas_only.shp`と名前を付けます。

• マップを保存します。

9.4.24. 各ポリゴンの面積の計算¶

• 新規ベクタレイヤの右クリックメニューを表示します。

• 属性テーブルを開く を選択します。

• Click the Toggle editing mode button along the bottom of the table, or press Ctrl+E.

• Click the Open field calculator button along the bottom of the table, or press Ctrl+I.

• Under the New field heading in the dialog that appears, enter the field name area. The output field type should be an integer, and the field width should be 10.

• :guilabel:`フィールド演算式`で、タイプ:

  $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.

• Click OK.

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

  $id

  This ensures that each polygon has a unique ID for identification purposes.

• Click Toggle editing mode again, and save your edits if prompted to do so.

9.4.25. 与えられたサイズの面積を選択¶

今このエリアは知られています：

• Build a query (as usual) to select only the polygons larger than 6000m^2. The query is:

  "area" > 6000
  

• Save the selection as a new vector layer called solution.shp.

You now have your solution areas, from which you will pick the one nearest to your house.

9.4.26. あなたの家のデジタイズ¶

• Create a new vector layer as before, but this time, select the Type value as being a Point.

• それが正しいCRSであることを確認します！

• 新規レイヤを:kbd:`house.shp`とします。

• 新規レイヤ作成の完了

• （新規レイヤが選択された状態で）編集モードにします。

• Click the point where your house or other current place of residence is, using the streets as a guide. You might have to open other layers to help you find your house. If you don’t live anywhere nearby, just click somewhere among the streets where a house could conceivably be.

• シェープIDの任意の番号を入力します。

• Click OK.

• 編集内容を保存し、編集モードを終了します。

• マップを保存します。

You will need to find the centroids (“centers of mass”) for the solution area polygons in order to decide which is closest to your house.

9.4.27. ポリゴンの中心点の算出¶

• Click on the Vector ‣ Geometry Tools ‣ Polygon centroids menu item.
• Specify the input layer as solution.shp.
• Provide the output location as Rasterprac.
• Call the destination file solution_centroids.shp.
• Click OK and add the result to the TOC (Layers list), then close the dialog.
• Drag the new layer to the top of the layer order so that you can see it.

9.4.28. あなたの家に最も近い重心の計算¶

• Click on the menu item Vector –> Analysis Tools –> Distance matrix.
• The input layer should be your house, and the target layer solution_centroids. Both of these should use the id field as their unique ID field.
• The output matrix type should be linear.

• 適切な出力場所と名前を設定します。

• Click OK.
• Open the file in a text editor (or import it into a spreadsheet). Note which target ID is associated with the shortest Distance. There may be more than one at the same distance.
• Build a query in QGIS to select only the solution areas closest to your house (selecting it using the id field).

これは、研究の質問への最終的な答えです。

For your submission, include 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 polygon of the closest solution area(s), as well as your house. Follow all the best practices for cartography in creating your output map.