ここでは、GRASSというGISソフトウェアを使った簡単な演習を説明しています。GRASSではリモートセンシング画像も表示、処理できます。このページは、GRASSのGUIでの操作を説明していきます。
演習内容: |
この演習1で最終的に得られる画像 ![]() |
GRASS起動後にWindow右下の「Projection values」をクリック
LOCATION: utm MAPSET: PERMANENT DATABASE: /home/student と入力したら、Escを押しながらEnter
Please specify the coordinate system for location <utm> A x,y B Latitude-Longitude C UTM D Other Projection RETURN to cancel >C(と入力してEnter)
Do you wish to specify a geodetic datum for this location? (y/n) [y] y Please specify datum name Enter 'list' for the list of available datums or 'custom' if you wish to enter custom parameters Hit RETURN to cancel request >list(と入力してEnter) (種々の測地系が表示される) Short Name Long Name / Description (中略) tokyo Tokyo (bessl ellipsoid) wgs72 WGS_1972 (wgs72 ellipsoid) wgs84 WGS_1984 (wgs84 ellipsoid) Please specify datum name Enter 'list' for the list of available datums or 'custom' if you wish to enter custom parameters Hit RETURN to cancel request >wgs84(と入力してEnter) Now select Datum Transformation Parameters Please think carefully about the area covered by your data and the accuracy you require before making your selection. Enter 'list' to see the list of available Parameter sets Enter the corresponding number, or <RETURN> to cancel request >list(と入力してEnter) Number Details --- 1 Used in whole wgs84 region (PROJ.4 Params towgs84=0.000,0.000,0.000) Default 3-Parameter Transformation (May not be optimum for older datums; use this only if no more appropriate options are available.) --- (中略) >1(と入力してEnter) Enter Zone: 53(と入力してEnter) Is this South Hemisphere (y/n) [n] n(と入力してEnter) DEFINE THE DEFAULT REGION ====== DEFAULT REGION ====== | NORTH EDGE: 4000000 | | | WEST EDGE | |EAST EDGE 380000 | |630000 | SOUTH EDGE: 3700000 | ============================ PROJECTION: 1 (UTM) ZONE: 53 GRID RESOLUTION East-West: 28.5 North-South: 28.5 と入力した後で、Escを押しながらEnter
→ GRASSのWindowが起動する
※ETM+ (Enhanced Thematic Mapper Plus) 画像のダウンロード:環境リモートセンシングの演習では既にダウンロード済みでローカルディスクにコピーしてあるが、今後は下記のWebサイトからダウンロードするといい。
Global Land Cover Facility (GLCF) Earth Science Data Interface, University of Maryland
UTM座標系のLandsat ETM+画像を使用して演習を行う。
次の3バンドのファイルに対し、下の画面の処理を繰り返す。
まずは、RGBの3層のレイヤーを指定し、redにb3、greenにb4、blueにb5を割り当てる。
次に、ディスプレイのウィンドウの上側にあるメニューバーの左端かその隣のアイコンをクリックして、画像を表示させる。
ディスプレイのウィンドウの上側にあるメニューバーの虫眼鏡のアイコンをクリックして、画像を拡大・縮小させる。
注意:上図のように画像が表示されないときは、「測地系の指定」での範囲の指定が間違っている可能性がある。その範囲を変更する手順を画像で示す手順を試すように。
上記の指定が終わったら、下図の左端のアイコンをクリックして新規にディスプレイを立ち上げて画像を表示させる。
植生の粗密と相関のあるNDVI (Normalized Difference Vegetation Index) と呼ばれる植生指数を表示させる。Landsat ETM+画像の場合、バンド3と4から計算できる。
ただし衛星画像に記録されているのはDN (Digital Number) と呼ばれる単なる数値であり、DNから放射輝度、反射率を経由してNDVIを計算するのが一般的である。また放射輝度あるいは反射率は、そのままでは大気上端での放射輝度や反射率を意味するので、大気補正を実施してからNDVIを計算することもある。ここでは大気補正をせず、DNから放射輝度を求め、NDVIを計算する手順を説明する。
/home/student/utm/PERMANENTの下にある「p110r036_7x20011015.met.txt」というファイルをテキストエディタで(Wordpadでも何でも)開いてみる。すると、下記の記述を見つけるだろう。
GROUP = MIN_MAX_RADIANCE (中略) LMAX_BAND3 = 152.900 LMIN_BAND3 = -5.000 LMAX_BAND4 = 157.400 LMIN_BAND4 = -5.100 LMAX_BAND5 = 31.060 LMIN_BAND5 = -1.000 (中略) END_GROUP = MIN_MAX_RADIANCE GROUP = MIN_MAX_PIXEL_VALUE (中略) QCALMAX_BAND3 = 255.0 QCALMIN_BAND3 = 1.0 QCALMAX_BAND4 = 255.0 QCALMIN_BAND4 = 1.0 QCALMAX_BAND5 = 255.0 QCALMIN_BAND5 = 1.0 (中略) END_GROUP = MIN_MAX_PIXEL_VALUE
DN値が0の時は放射輝度(L3やL4, L5)を0に設定する条件処理を行う
> r.mapcalc "L3=if(b3==0, 0, (b3-1)*(152.9+5.0)/(255.0-1.0)-5.0)" > r.mapcalc "L4=if(b4==0, 0, (b4-1)*(157.4+5.1)/(255.0-1.0)-5.1)" > r.mapcalc "L5=if(b5==0, 0, (b5-1)*(31.06+1.0)/(255.0-1.0)-1.0)"
上の式をGUIで処理させるには、下記の手順に従って行う。
手順1:「Map calculator」を選択
手順2:入力パラメータに例えば「b3」を、式の欄に「if(b3==0, 0, (b3-1)*(152.9+5.0)/(255.0-1.0)-5.0)」を入力する。
画像の端の領域が放射輝度が0になっているかを確認するには、画像を表示させた上で、
をクリックしてから画像上でクリックすれば、下図のように座標とRGBの値が表示される。
一面緑色の意味のないような画像が表示されるだろう。しかし、画像上でクリックすれば-1から1までの値が入っていることが分かる。この範囲に合わせたカラーテーブルを割り当てよう。
上図の「Name of input raster map」で指定した「ndvi」は生成した画像を、「Name of predefined rules file」の「ndvi」は予め用意されたカラーテーブルを指定することを意味している。
最終的に得られるNDVI画像を下記に示す。