実習です！！ 実習です！！今回は講義で説明した地すべりを例として，サイトとラスターレイヤー間の簡単な統計と演算を行ってみましょう．
なお，地域は今までと同じですが，前回を済ませていない方も実習可能なように地形を含めた新しいデータセットを準備しました．実習の準備として下記をクリックして，データをダウンロードして下さい．新しいMAPSETで作業を行ってください．
Step０．データのダウンロードと解凍
vunivdata.tar.gz （←をクリックし，ダウンロード後ファイル名を確認）ここでは，上記のファイルが
/home/grass のディレクトリにダウンロードされたとして説明をします．各自で/home/grass のディレクトリ名を読み替えて以下の作業を行なってください．データを解凍します．
tar -zxvf vunivdata.tar.gz ↓
(↓は，リターンを意味します)
Step１．GRASSの起動
以下のLOCATIONとMAPSET名でGRASSを立ち上げます．なお，
/home/grass は，各自で先ほどのディレクトリ名に変更してください．(↓は，リターンを意味します)>
grass4.2 ↓LOCATION=
vunivdataMAPSET=
vuniv99DATABASE=
/home/grass
Step２．地形データと地すべり地データの確認
grass4.2.1>
d.mon x1 ↓(モニターを開きます)
grass4.2.1>
d.rast dem ↓（地形データ（dem）を表示します）
demの表示例（163KB）
grass4.2.1>
d.sites landslide ↓（地すべり地を示すサイトデータを表示します）
地すべり地の表示例（54KB）
Step３．地すべりデータの構築
ASCII形式で書かれたベクトル型の地すべり地データをGRASSのバイナリー形式のベクトルファイルに変換します．ASCII形式のベクトルファイルは，現在のMAPSETの下の
dig_ascii ディレクトリにある必要があり，ここでは jisuberi というファイル名でデータを入れておきました．この変換は，v.in.asciiコマンドにより行ないます．grass4.2.1>
v.in.ascii input=jisuberi output=landslide ↓（jisuberiというASCII形式のベクトルファイルを入力し，landslideというベクトルレイヤーを作成する）
次に，このlandslideというレイヤーのサポートファイルの作成とラベル化を行ないます．
grass4.2.1>
v.support map=landslide option=build ↓grass4.2.1>
v.alabel map=landslide value=1 ↓作成したベクトルレイヤー（landslide）を表示してみます．
grass4.2.1>
d.erase ↓grass4.2.1>
d.vect landslide ↓
Step４．ベクトル→ラスター変換
Step３で作成した地すべりのベクトルデータをラスターデータに変換します．
grass4.2.1>
v.to.rast input=landslide output=landslide ↓（ベクトルデータと同じ名前でラスターデータに変換しています）
変換後のラスター化された地すべりデータを表示してみましょう．
grass4.2.1>
d.erase ↓grass4.2.1>
d.rast landslide ↓
地すべり地域の表示例（38KB）
Step５．傾斜量図の作成
地形の傾斜の大きさ（傾斜量図；slope）を地形データ（dem）から作成します．
grass4.2.1>
r.slope.aspect elevation=dem slope=slope ↓（標高データ（elevation）としてdemを指定し，計算結果の傾斜量図（slope）としてslopeというファイル名をつけます）
地形の傾斜量図を表示してみましょう．
grass4.2.1>
d.rast slope ↓
Step６．起伏量図の作成
近傍点の計算を行なうr.neighbors コマンドを用いて，地形データから各点の近傍の最大値と最小値を求めます．
grass4.2.1>
r.neighbors input=dem output=max size=5 method=maximum ↓grass4.2.1>
r.neighbors input=dem output=min size=5 method=minimum ↓（各点の周囲５×５の近傍領域の最大値（method=maximun）と最小値（method=minimun）をそれぞれ求め，maxとminというラスターを作成します（５×５近傍の接峰面と接谷面図に相当します））．求めた結果からr.mapcalcコマンドを用いて起伏量に変えます．
grass4.2.1>
r.mapcalc relief=max-min ↓（最大値（max）−最小値（min）の計算を行ない，結果をreliefというファイルに出力します．）
地形の起伏量図を表示してみましょう．
grass4.2.1>
d.rast relief ↓
Step７．傾斜量と起伏量の再分類
傾斜量と起伏量のデータを再分類するためのテーブルをあらかじめ作成し，r.reclassコマンドを用いて再分類します．２つの再分類テーブルは
viやxeditなどを用いて作成してください．・傾斜量の再分類
傾斜を５°単位で区切ったファイルを傾斜量のファイルから再分類して作成します．そのために傾斜量を再分類するための再分類テーブルとして，slope5degというファイル名で下記の内容のファイルを作成します．
0 = 1
1 thru 5 = 2
6 thru 10 = 3
11 thru 15 = 4
16 thru 20 = 5
21 thru 25 = 6
26 thru 30 = 7
31 thru 35 = 8
36 thru 40 = 9
41 thru 45 = 10
46 thru 48 = 11
（傾斜は，48°が最大のため，最後は48にしています）
grass4.2.1>
r.reclass input=slope output=slopereclass < slope5deg ↓（傾斜量のファイル（slope）を入力し，再分類結果をslopereclassというファイル名で出力します．< はUnixのリダイレクションで再分類テーブルをslope5degというファイルから入力することを示します）
傾斜量の再分類結果表示例（56KB）
・起伏量の再分類
起伏量を約10ｍ単位で区切ったファイルを起伏量のファイルから再分類して作成します．そのために起伏量を再分類するための再分類テーブルとして，relief10というファイル名で下記の内容のファイルを作成します．
0 = 1
1 thru 10 = 2
11 thru 20 = 3
21 thru 30 = 4
31 thru 40 = 5
41 thru 50 = 6
51 thru 60 = 7
61 thru 70 = 8
71 thru 80 = 9
81 thru 90 = 10
91 thru 113 = 11
（起伏を約10ｍごとに分けています）
grass4.2.1>
r.reclass input=relief output=reliefreclass < relief10 ↓（起伏量のファイル（relief）を入力し，再分類結果をreliefreclassというファイル名で出力します．< はUnixのリダイレクションで再分類テーブルをrelief10というファイルから入力することを示します）
起伏量の再分類結果表示例（52KB）
Step８．統計量の計算
サイト形式で表した地すべり地（landslide）とラスター形式の再分類された傾斜量（slopereclass）・起伏量（reliefreclass）との関係をs.menuコマンドを用いて調べます．
grass4.2.1>
s.menu ↓Please select one of the following
1 Read an existing site list
2 Mask current site list
3 Save the current site list in your mapset
4 Check site list for duplicate sites
5 Edit site list using a unix editor
6 Convert site list to raster file (0/1)
7 Convert site list to raster file (frequency of occurence)
8 Run reports on the current site list
stop Leave the SITES program
choice>
1 ↓ （サイトファイルを選択するために１を選びます）
Read an existing site list
Enter the name of an existing site list file
Enter 'list' for a list of existing site list files
Hit RETURN to cancel request
>
landslide ↓ （サイトファイルとして地すべり地（landslide）を指定します）
Please select one of the following
1 Read an existing site list
2 Mask current site list
3 Save the current site list in your mapset
4 Check site list for duplicate sites
5 Edit site list using a unix editor
6 Convert site list to raster file (0/1)
7 Convert site list to raster file (frequency of occurence)
8 Run reports on the current site list
stop Leave the SITES program
choice>
8 ↓ （指定したサイトファイルに関するレポートの作成作業を選択します）
Run reports on the current site list
You may select one or more map layers to be analyzed
select first map layer
Enter 'list' for a list of existing raster files
Enter 'list -f' for a list with titles
Hit RETURN to cancel request
>
reliefreclass ↓ （ラスターファイルとして再分類した起伏量のファイルを指定します）<reliefreclass>
selected layers: reliefreclass
select another map layer
Enter 'list' for a list of existing raster files
Enter 'list -f' for a list with titles
Hit RETURN to cancel request
>
slopereclass ↓ （ラスターファイルとして再分類した傾斜量のファイルを指定します）<slopereclass>
selected layers: reliefreclass slopereclass
select another map layer
Enter 'list' for a list of existing raster files
Enter 'list -f' for a list with titles
Hit RETURN to cancel request
>
↓ （２つを指定し終えたので，リターンします）
SITE REPORTS
Please set the site quad mask size.
This is the number of cells away from the site to be included in the analysis
eg, enter
0 if you just want to include the site itself
1 if you want to include the 8 cells surrounding the site
2 if you want to include both the above 8 cells and the 14 cells
surrounding those 8, etc.
enter <stop> if you want to stop
site quad mask size>
0 ↓ （レポートを作成する際に，サイトデータの１地点に対するものか，周囲のどれくらいを取るかを指示します．０はサイトのみを考えます）
SITES REPORT MENU
Please select one of the following
1 Site characteristics report
2 Site occurrence report
3 Convert data to S input format
4 Produce machine readable data file
stop return to SITES MAIN MENU
choice>
2 ↓ （再分類した両ラスターに関するサイトのレポートを見るために，２を選択）
SITE REPORTS
report format phase
reliefreclass
slopereclass
Report ready
would you like to see the report? (y/n)
y ↓ （yを入力します）（以下のように計算結果が表示されます）
表−１ 再分類した起伏量と地すべり地点との関係
表−２ 再分類した傾斜量と地すべり地点との関係
would you like the report sent to the printer? (y/n)
n ↓（プリンタへの出力はしない）
would you like to save the report in a file? (y/n)
y ↓ （ファイルへの保存は行なう）
SITE REPORTS
enter the name of a file to hold report
(this file will be relative to your home directory)
or hit RETURN to cancel
file>
landslide.rep ↓ （レポートの出力ファイル名を指定します）
SITES REPORT MENU
Please select one of the following
1 Site characteristics report
2 Site occurrence report
3 Convert data to S input format
4 Produce machine readable data file
stop return to SITES MAIN MENU
choice>
stop ↓ （s.menuの終了を指示します）choice>
stop ↓ （再度，s.menuの終了を指示します）
Step９．統計量の計算
表−１，表−２に示した再分類した起伏量・傾斜量と地すべり地点の関係を考えます．両表中の黄色で塗ったexpected sitesというのは，もし，その地域に地すべりが，一様に分布すると考えた場合，そのそれぞれのカテゴリーには，どれくらいの地すべり地があることが期待されるかを示したものです．計算は
expected sites＝地すべり地総数×（カテゴリーの面積（セル数）／（本地域の全面積（セル数））
で，例えば，起伏量のカテゴリー１の面積は848（セル数）で，全面積は45269（セル数），地すべり総数は123個から，expected sites＝2.3となります．actual sitesは実際にそのカテゴリーに含まれる地すべり地の数です．
ここで，そのカテゴリーに含まれる地すべり地が多いか，少ないかを決める判断基準として，
（actual sites／expected sites）
を考えてみましょう．expected sitesというのは，説明したように，もし地すべり地が一様にあるとしたときの期待される数ですから，それを基準に多いか少ないかを判断します．もし，この値が１より大きければ，一様にあるより多いといえます．起伏量・傾斜量とも，全てのカテゴリーについて，この計算を行ないます（電卓）．さらに，この計算結果を整数化して重ね合せることができるようにします．これらの結果は次のようになります．
Step10．２つの情報の結合
この傾斜量・起伏量と地すべり地の関係を結合する方法は２通りあります．１つは，それぞれ再分類した傾斜量(slopereclass)・起伏量(reliefreclass)を計算結果の整数化したものにもとづき，再度再分類（r.reclass）を行ない，その２つの結果をr.mapcalcで単純に足し算する方法です．もう１つは，GRASSのr.weightコマンドを用いて結合する方法です．
どちらの結果とも同じものが得られます．この結果を適当に色づけして表示したものを下に示します．
分類結果の図（59KB）
この図では，凡例にあるように０〜２までを白〜灰色で示し，３以上を黄色〜赤で示しています．これは２以下が一様に地すべり地があることよりも地すべりが少ないと判断し，３以上は高く，値が高いほど（赤に近づくほど）高いと考えたからです．
なお，このラスターファイルのカラーテーブルは以下の通りです．カラーテーブルを直接作成する方法は，Mapsetディレクトリの下に
colrというカラーテーブルのディレクトリがあります（この場合は，/home/grass/vunivdata/vuniv99/colr）．この中に作成したラスターファイル名と同じ名前でカラーテーブルをvi等で直接作ります．% 1 22
0:255:255:255
1:200:200:200
2:128:128:128
3:255:255:0 22:255:0:0
最後に，ここで求めた傾斜と起伏を結合した情報と地すべり地との関係をs.menuコマンドで調べて見ましょう．
