Subject: [ts-mag:00097] OpenGIS Mail Magazine for T.S.U 【第 90号】
    Date: Fri, 14 Mar 2003 18:41:18 +0900
    From: Taichi FURUHASHI

OpenGIS Mail Magazine for Technical Support Users       2003/3/14
Keywords: OpenGIS,TNTmips,FocalMean,SML,IgnoreNull,Add All,Pyramid


オープンGIS
テクニカル・サポート専用メールマガジン

第90号「 ど う し て N u l l (ヌ ル) を 認 識 し な い の ?  」
#################################################################
                       株式会社 オープンGIS



接峰面(せっぽうめん)という言葉をご存知ですか?

地形を見る際、細かい谷に削られた場所を埋め戻すようにして、まだ侵食
されていない元の地形に近い形で復元する、地形解析手法の代表的なテク
ニックになります。


●


昔は、地形図を小さなメッシュに区切って、そのメッシュ内で最も標高の
高いところを抜き出し、その情報から大まかな等高線を手作業で引きなお
すことが主流でしたが、GIS を利用することで、ある程度、接峰面の自動
生成が可能になります。
 メルマガ第88号(2003/2/28 配信)にて「空間フィルタ」の話題をご紹介
しましたところ、ユーザーさまより「空間フィルタ」を利用して、接峰面
生成 SML を自作したとの連絡をいただきました。
 しかし同時に、9割がた完成したものの、どうしても海の部分に与えて
いたヌル(Null)値の処理がうまくいかないとのご相談も受けまして、今回
メルマガ上で、そのソースプログラムをご覧いただき、修正箇所などを、
皆さんにもご紹介したします。

TNTmips でのプログラミングの雰囲気と、注意する点などご覧いただけれ
ば幸いです。


●


まず、問題点と修正目標を先にご紹介しましょう。

・元の標高データ(DEM)は、このようなデータです。
陰影図から地形の起伏が細かいことがわかりますね。




・送っていただいた SML による出力結果はこのデータです。
空間フィルタによって、平滑化されているのですが、海岸線付近の平滑化
が影響して、陸地の面積が広がっていることがわかるかと思います。




・修正目標の出力結果はこのようなかたちです。平滑化はされつつも、
海の部分は元データと同じように出力されるというのが今回の目標です。




●


早速、送っていただいた SML を拝見させていただきましょう。
今回のヌル値処理以外につきましても、非常に参考になる点が多いスク
リプトですので、じっくりご覧ください。

尚、# の後はコメントとなりますので、処理には影響いたしません。





【送っていただいた SML】

clear()

## 必要なデータの設定
GetInputRaster(OrigDEM)
GetOutputRaster(DummySummit,NumLins(OrigDEM),NumCols(OrigDEM),RastType(OrigDEM))
GetOutputRaster(Summit,NumLins(OrigDEM),NumCols(OrigDEM),RastType(OrigDEM))

## 3つのデータとも Null を使います宣言
UseNull(OrigDEM)
UseNull(DummySummit)
UseNull(Summit)

## まずはダミーラスタと出力ラスタにコピー
for each OrigDEM begin
        DummySummit = OrigDEM
        Summit = OrigDEM
end

## 繰り返し、5x5 の窓枠平均を計算するためのループ
i = 1
while (i <= 20) begin       ### 繰り返し開始

   ## for each 文で出力ラスタを各ピクセルごとに空間フィルタをかける
   for each Summit begin
           ## 得られた値を四捨五入して整数値に
           ## ダミーラスタにその値を代入
           DummySummit = int(0.5 + FocalMean(Summit, 5, 5))
   end

   ## ダミーラスタとオリジナルラスタを比較
   ## もしオリジナルよりもダミーが高ければそのまま
   ## もしダミーが低ければオリジナルの値を代入
   ## この処理を繰り返すことで、小さな谷地形が埋められていきます
   for each DummySummit begin
           if (DummySummit < OrigDEM) then
                   DummySummit = OrigDEM
   end

   ## ダミーラスタの値を出力ラスタにコピー
   for each DummySummit begin
           Summit = DummySummit
   end

   ## 繰り返しの数をカウントしながらコンソールに表示
   printf("Iteration = %d \n" , i)
   i = i + 1

end      ### 繰り返し終了

   ## オリジナルラスタから海などの Null 値を検出
   for each OrigDEM begin
           ## もし Null だったら、出力ラスタも Null にする
           if (OrigDEM == NullValue(OrigDEM)) then
                   Summit = NullValue(OrigDEM)
           end

## 標高ラスタから出力ラスタへ
## 座標情報を司るジオリファレンスサブオブジェクトを
## まるごとコピー            
CopySubobjects(OrigDEM, Summit, "georef")

## 作業終了を教える音を鳴らす
beep()

## 使用したラスタをすべて閉じる
CloseRaster(OrigDEM)
CloseRaster(DummySummit)
CloseRaster(Summit)





さて、いかがでしょうか。

きちんと Null 値の場合の処理も行われているように見えるのですが、
実際に実行してみますと、なぜだかうまくいきません。

そこで、修正と詳細のコメント、更に僕の癖もあえて加えまして、プロ
グラムを変更してみました。といいましても、とても美しいとはいえな
い力技のプログラムですので、あえてこんな書き方をする場合もある程
度にお考えください。

では、お恥ずかしながら変更例をご覧ください。





【古橋の変更例】

clear(); ## コンソール画面を消去するおまじない
 ## ほとんどの SML はまずこのクリアから始まります。
 ## SML の場合セミコロンはあってもなくても良いのですが、
 ## C言語の癖を忘れない為にあえて付けています。

## 元データの設定
GetInputRaster(OrigDEM); ## 入力する標高ラスタを訊ねる関数です。
 lin = NumLins(OrigDEM); ## 入力した標高ラスタのライン数をチェック
 col = NumCols(OrigDEM); ## 入力した標高ラスタのカラム数をチェック

## 入力した標高ラスタのビット深度をチェック
rastertype$ = RastType(OrigDEM); 
 ## 文字列名は、頭文字が小文字、
 ## 一番最後に$を付けるというルールがあります。

## ダミーラスタの作成
CreateTempRaster(DummySummit, lin, col, rastertype$);
 ## 処理の途中だけ必要なラスタは、テンポラリラスタとして
 ## 利用したほうが、ハードディスクの節約にもなりますし、
 ## ラスタ名の入力の手間も省けます。

## 出力ラスタの設定。
GetOutputRaster(Summit, lin, col, rastertype$);
 ## なるべく1行を短く見やすくするために
 ## 引数はあらかじめ変数として準備しています。

## 元データを、ダミーラスタ、出力ラスタ両方へまるごとコピーします。
for each OrigDEM { 
 ## for each 文を使用するとラスタ内の各ピクセルごとに
 ## 処理をしてくれるのでかなり便利です。
   DummySummit = OrigDEM;
   Summit = OrigDEM;
 ## begin と end を使わずに、{}を使うのは、好みの問題ですが、
 ## 見かけ上シンプルになるため、僕はよく{}を使用します。
}

## カウント用の値を初期化
i = 1; 

## カウントしながら平滑化の計算を行い 20回実行すると処理から抜けます。
while (i <= 20){  ### start iteration
 ## ここから繰り返し平滑処理のスタートですね。
 ## for each 文と同じように、あえて {} に変えてみました。

   ## ここで 空間フィルタである FocalMean 関数を使用
   for each Summit {
           DummySummit = int(0.5 + FocalMean(Summit, 5, 5));
   }
    ## 引数の2つの 5 は窓領域(5x5)の大きさを意味します。

   ## 接峰面は、ある範囲の最高標高を使用しますので、
   ## 平均をとりながらも、もし元データより標高が低い場合は
   ## 値を入れ替えて、絶対に元データより低い値にはならないよう
   ## チェックしてあげます。
   for each DummySummit {
           if (DummySummit < OrigDEM){
                   DummySummit = OrigDEM;
           }
   }

   ## 元データとの高さチェックが終わると
   ## 出力ラスタへまるごとコピーしてあげます。
   for each DummySummit {
           Summit = DummySummit;
   }

   printf("Iteration = %d \n", i);
   i = i + 1;
    ## これで、平滑化の処理が1回完了です。
    ## カウントを増やして、2回目、3回目と処理をすすんでいきましょう。

}        ### end iteration
 ## ここまでが繰り返し処理です。
 ## タブをうまくあわせて{}を配置すると、
 ## 対応する括弧を見つけやすくなります。


##### ここからが今回の重要な修正箇所です

## まず、下データの Null 値設定されていた値を調べます。
orig_null = NullValue(OrigDEM);

## 元データのヌル値不使用宣言をします。
IgnoreNull(OrigDEM);
 ## こうしないと、次にでてきます if 文の判定で、
 ## 数値としてみてくれなくなります。

   for each OrigDEM {
           ### もし、元データのピクセル値が Null と設定されていた
           ### 値であったならば、同じ値を出力ラスタにコピーしなさい。
           if (OrigDEM == orig_null ){
                   Summit = orig_null;
           }
            ## これで、Null 値に使用されていた数値が
            ## 出力ラスタにコピーされました。
   }

## 最後に出力ラスタに対して、Null 値の設定をしてあげます。
SetNull(Summit, orig_null);

##### ここまでが今回の重要な修正箇所です。


## ジオリファレンスサブオブジェクトのコピー
CopySubobjects(OrigDEM, Summit, "georef");
 ## 元データから、ジオリファレンスサブオブジェクトを
 ## まるごとコピーしておくと、重ね合わせのとき便利です。

## 処理終了を教えるビープ音を鳴らす
beep();

## 作成したラスタをすべてクローズ。
CloseRaster(OrigDEM);
CloseRaster(DummySummit);
CloseRaster(Summit);

## 処理の終了宣言
Exit()
 ## とくになくてもいいのですが、癖としていつもつけています。





いかがでしょうか。

まず、ヌル(Null)という値は、処理の対象外として扱われますので、もし
Null値に設定した数値を利用して処理する場合には、一度 IgnoreNull()
関数にてヌル設定を外してあげる逃げ技が有効になります。

今回は処理が終わった後で、再度ヌルとして設定してみました。










今週の新プチ・テクニック
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
" 
" そんなに沢山クリックするのは面倒!!
" 
" Select Object ウィンドウでデータ(オブジェクト)を開くときに、
" データの名前をダブルクリックしたり、名前の先頭にあるマークを
" シングルクリックしますが、沢山のデータを何度も何度もクリック
" するのはとても面倒です。
" 
" そんなときは、下にあります黄色い十字の Add All ボタン。
" 
" RVC ファイル内にあるすべての対象データを選択してくれます。
" 
" さらにこのボタンは、RVC の外でも有効ですので、
" 同じフォルダ内にあるすべての RVC ファイル内データを選択すること
" だってできます。
" 
" 逆にすべてをキャンセルするには横にあります、Remove All ボタン
" 
" 是非一度お試しください。
" 慣れると病みつきになります :-)
" 
" 
""""  隔 週 で 新 / 過 去 テ ク ニ ッ ク を 交 互 に 発 信  """"




New Things 最新情報
=================================================================
= 
= 申し訳ありません。今週は、更新がありませんでした。
= 来週をお楽しみに。
= 
=================================================================




バグ・トラブル情報
??? 
??? 
??? ラスタデータを表示したときに、
??? どうもカラーパレットで指定した色と違う...
??? 
??? この原因として2つの事が考えられます。
??? 
??? (1)コントラスト自動調整が働いている。
??? (2)ピラミッドデータに異常がある。
??? 
??? 
??? (1)の場合は Group Controls ウィンドウにて、対象のラスタレイヤ
??? を Active(赤丸状態)にしてから、メニューより Layer/Controls...
??? を選択し Raster Layer Controls ウィンドウを表示します。
??? そして、Object タブの Contrast: を None にすると、コントラスト
??? 自動調整がオフになり、カラーパレットと色が一致します。
??? 
??? (2)の場合は、メインメニューより、Process/Raster/Utilities/
??? Pyramid... を選択しまして、Raster Pyramiding ウィンドウを表示
??? します。ウィンドウ内の Raster ボタンを押しまして、対象のラスタ
??? データを指定後、Compute full pyramid チェックをオンにした状態
??? で Run ボタンを押しますと、ラスタデータの表示を早くするための
??? 間引き画像データである、ピラミッドデータを再構築します。
??? 
??? 
??? どうも、イメージどおりの色合いになってくれない...
??? そんなときは上記のチェックをしてみてください。
??? 
??? 
??? 
??? 





今週の話題 on メーリングリスト
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
~ 今週はとくにありませんでした。
~ 
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 メーリングリストに参加されたい方は、
 info@opengis.co.jp までご連絡ください。
 弊社からの対応は バグ対応と、新製品告知のみですが、
 ユーザーさまの間での情報交換の場として
 多くの方に利用されています。





■■■■■■■ 現 在 の T N T m i p s の 価 格 ■■■■■■■

・TNTmips シングルライセンス 840,000円(税別・送料別)
・年間バージョンアップ(2回分) 140,000円(税別)
・年間テクニカルサポート 90,000円(税別)


普段より、元価格にあまり上乗せしていない分、
円相場に応じて価格が変動いたしますことを
ご理解よろしくお願いいたします。

わからない点などありましたら、
info@opengis.co.jp までご連絡ください。

■■■■■■■■■■■■■■■■■■■■■■■■■■■■■




############# このメールマガジンの配信について ################

このメールマガジンは、
年間テクニカル・サポートに加入されている、
ユーザーさまのみを対象としております。

基本的には、登録アカウント1つに対して
1ユーザーを配信対象としておりますが、
例えば、研究室の予算で加入された場合などは、
ご購入いただいた指導教官と、実際に TNTmips で作業を行う学生1名の、
計2ユーザーを対象として配信させていただきます。

1アカウントに対して合計3名以上の配信は行いません。
学生さん2名が同じくらい使用する場合は、
申し訳ありませんが、じゃんけん等で1人に決めてください。



また、テクニカル・サポートの期限を半年以上過ぎてしまった場合は、
メールマガジンの配信を停止させていただきます。

メールマガジンが配信されなくなった時は、
サポート期限が半年過ぎたとお考えください。



転勤、ドメイン変更などでメールアドレスが変わった場合は、
お手数ですが、弊社までご連絡ください。



テクニカル・サポート活用案内
>サービスその1:緊急なバグ情報はすぐにメール配信します。
>サービスその2:TNTmipsに関するどんな質問にも対応します。
>サービスその3:対応は最優先で行います。
>サービスその4:専用FAX質問用紙をお送りします。
>サービスその5:電話・電子メールも、もちろん最優先で対応いたします。
>サービスその6:週に1回メールマガジンを配信します。
>サービスその7:展示ブースへご来場の方には特製オリジナルグッズを進呈。
>サービスその8:サポート対応で FTP サイトを 20MB まで利用できます。
>サービスその9:現在計画中です...



「これで年間9万円は安い!!」と感じていただけるような
サービスを提供していきたいと考えております。

よろしくお願いいたします。




============ おことわり =============

※OpenGIS Mail Magazine for T.S.Uでは、
 TNTmipsに関わる新しいニュースを、
 毎週金曜日に、皆さまへご提供させていただきます。

※このメールマガジンの転載・転送はご遠慮ください。

================================