流域推定プログラム
【課題】流域に関する知識のないユーザによる正確な流域推定を可能にする。
【解決手段】地図上の地点を指示する入力を受け付ける手順と、前記指示された地点を含む前記地図上の格子を、処理対象である第１格子として選択する手順と、前記第１格子を含む領域内の格子の標高値をメモリに格納する手順と、前記領域内の窪地を検出する手順と、前記標高値を変更する手順と、前記領域内の格子の標高値に基づいて、前記領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する手順と、前記第２格子を、処理対象である第１格子として新たに選択する手順と、前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する手順と、を実行するプログラム。
【解決手段】地図上の地点を指示する入力を受け付ける手順と、前記指示された地点を含む前記地図上の格子を、処理対象である第１格子として選択する手順と、前記第１格子を含む領域内の格子の標高値をメモリに格納する手順と、前記領域内の窪地を検出する手順と、前記標高値を変更する手順と、前記領域内の格子の標高値に基づいて、前記領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する手順と、前記第２格子を、処理対象である第１格子として新たに選択する手順と、前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する手順と、を実行するプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本願明細書で開示される技術は、コンピュータを使って流域を推定する方法に関する。本願明細書で開示される技術は特に、１以上の地点の位置を入力データとし、それぞれの地点の流域を数値地形モデル（ｄｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅ，ＤＥＭ）に基づいて推定するコンピュータプログラムに関する。
【背景技術】
【０００２】
流域は、河川、ダム、下水などの水を管理する土木構造物にとって非常に重要な情報である。多くの場合、土木構造物に入ってくる水の起源は降雨である。つまり、降雨が地表面に落ち、地表面を流れ下って構造物に入る。したがって、構造物に至る全ての流れをさかのぼれば、降雨が落ちた地表面に達する。この地表面を構造物の流域という。
【０００３】
近年では、地表面の標高は電子データとして得られる。このような標高の電子データは、一般的にＤｉｇｉｔａｌ Ｅｌｅｖａｔｉｏｎ Ｍｏｄｅｌ （ＤＥＭ）と呼ばれる。ＤＥＭでは、ある領域の地表面が格子（ｇｒｉｄ又はｃｅｌｌ）と呼ばれる小領域に分割され、それぞれの格子に標高値が与えられる。多くの場合、ＤＥＭの格子は正方形である。
【０００４】
このＤＥＭを使い、土木構造物の流域を推定するコンピュータプログラムが作られている。通常、地表面に落ちた降雨は最も急な下り勾配（ｇｒａｄｉｅｎｔ）の方向に流れ下る。したがって、ＤＥＭの格子ごとの勾配を求めれば、流域を推定できる。勾配から流域を推定する方法は、非特許文献１に詳しい。
【０００５】
このような流域を推定するコンピュータプログラムに共通の課題は、ＤＥＭの誤差である。勾配はベクトル解析における勾配演算によって算出される。勾配演算は元データの誤差に敏感である。したがって、ＤＥＭデータから流域を推定する方法は誤差に敏感である。
【０００６】
図１８は、ＤＥＭの誤差の影響を説明するための図である。
【０００７】
グラフ１８００は、ある断面に沿った標高である。格子１８０１から格子１８０９までの格子に標高が与えられ、勾配が算出されている。勾配を矢印（たとえば１８１１）で示す。グラフ１８００では、格子１８０６の標高が、隣接する格子１８０５及び１８０７の標高よりも低い。下り勾配の方向（すなわち図の矢印の方向）と水が流れる方向とが一致すると仮定すると、格子１８０１から格子１８０７までに降った雨は格子１８０６に流れ出て、格子１８０６からは水が流れ出ないことになる。また、格子１８０９の流域は格子１８０７及び１８０８となる。この結果と、地表面を流れる水は最終的に標高の低い場所に至るという私たちの常識とは矛盾する。この矛盾の原因は、格子１８０６近傍の標高データの誤差にあると推測される。
【０００８】
誤差の修正方法として、いくつかの方法が提案されている。最も一般的な修正方法は、窪地埋め（ｓｉｎｋ ｆｉｌｌｉｎｇ）である。窪地とは１以上の格子からなる、周囲に比べ相対的に標高が低い領域である。窪地埋めの方法は非特許文献１に詳しく説明されている。この方法では、ある領域内の全ての格子のうち、ある格子の標高が、隣接するすべての格子の標高と同じ又は低い場合、ごく小さな値をその格子の標高に加える。この処理を繰り返すことによって修正が完了する。このように誤差が修正されたＤＥＭを元に推定した流域と、私たちの常識とは一致する。
【０００９】
図１８のグラフ１８２０は、グラフ１８００に示したＤＥＭデータをこの方法によって修正した結果である。格子１８２１から格子１８２９は、それぞれ、修正前のグラフ１８００の格子１８０１から格子１８０９に対応する。窪地（すなわち格子１８２５及び１８２６）が埋められた結果、格子１８２９の流域が格子１８２１から格子１８２９までとなった。
【００１０】
さらに別の誤差修正方法として、特許文献１に開示された「Ｍｅｔｈｏｄ ｆｏｒ ｆｌｏｗ ａｎａｌｙｓｉｓ ｕｔｉｌｉｚｉｎｇ ｃｅｌｌ ｅｌｅｖａｔｉｏｎ ｔｏｐｏｇｒａｐｈｙ」が挙げられる。この方法では、あらかじめ領域の内部が複数のｍａｊｏｒ ｃｅｌｌに分割される。さらに、ｍａｊｏｒ ｃｅｌｌの内部が、少なくとも横方向に３分割、縦方向に３分割することによって、少なくとも９個のｍｉｎｏｒ ｃｅｌｌに分割される。あるｍｉｎｏｒ ｃｅｌｌにおける勾配を算出するときは、そのｍｉｎｏｒ ｃｅｌｌを含むｍａｊｏｒ ｃｅｌｌと、さらにそのｍａｊｏｒ ｃｅｌｌに隣接するｍａｊｏｒ ｃｅｌｌについてＤＥＭデータのエラーが修正される。その後、これらのｍａｊｏｒ ｃｅｌｌに含まれるｍｉｎｏｒ ｃｅｌｌについて勾配が再度算出される。
【００１１】
非特許文献１では処理の対象とする領域全体に誤差修正処理をする必要があったが、特許文献１に開示された方法では処理の対象とする領域を複数に分割し、分割された領域ごとに勾配を求めることが可能になる。
【００１２】
いずれの誤差修正方法であっても、ＤＥＭデータの誤差を修正した後に流域を推定する。つまり、まずＤＥＭデータを修正する。次にＤＥＭデータから勾配データを算出する。この勾配データを使って上流方向にたどると、流域が求まる。
【特許文献１】米国特許出願公開第２００５／００４９８０２号明細書
【非特許文献１】Ｃｏｓｔａ−Ｃａｂｒａｌ，Ｍ．Ｃ．，ａｎｄ Ｓ．Ｊ．Ｂｕｒｇｅｓ（１９９４），Ｄｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅｌ ｎｅｔｗｏｒｋｓ（ＤＥＭＯＮ）： Ａ ｍｏｄｅｌ ｏｆ ｆｌｏｗ ｏｖｅｒ ｈｉｌｌｓｌｏｐｅｓ ｆｏｒ ｃｏｍｐｕｔａｔｉｏｎ ｏｆ ｃｏｎｔｒｉｂｕｔｉｎｇ ａｎｄ ｄｉｓｐｅｒｓａｌ ａｒｅａｓ，Ｗａｔｅｒ Ｒｅｓｏｕｒ．Ｒｅｓ．，３０（６），１６８１−１６９２．
【発明の開示】
【発明が解決しようとする課題】
【００１３】
従来の流域推定プログラムは、ＤＥＭデータの誤差を修正した後に、流域を推定する。したがって、流域を推定する前に、修正すべきＤＥＭデータの範囲を知っておく必要がある。ところが流域を推定したい場合、修正すべきＤＥＭデータの範囲は流域である。つまり従来の流域推定プログラムのユーザは、そのプログラムの実行結果である流域を、プログラムを実行する前に知っておく必要がある。したがって、従来の流域推定プログラムを使えるユーザは、流域をあらかじめ知っているユーザに限られていた。流域を知らないユーザが従来の流域推定プログラムを使用すると、ＤＥＭデータを適切に修正できないために、流域を正確に推定することが困難であった。
【００１４】
本発明はこの問題に鑑みてなされたものである。本発明の目的は、流域に関する事前知識のないユーザが使用しても、流域を正確に推定することができる流域推定コンピュータプログラムを提供することである。
【課題を解決するための手段】
【００１５】
本願で開示する代表的な発明は、計算機を制御するプログラムであって、前記計算機は、前記プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備え、前記プログラムは、地図上の地点を指示する入力を受け付ける第１手順と、前記地図を互いに直交して設定された第１軸及び第２軸に沿ってそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、前記第１領域内の窪地を検出する第４手順と、前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を前記プロセッサに実行させることを特徴とする。
【発明の効果】
【００１６】
本発明の一実施形態によれば、流域に関する事前知識を要求しない、流域推定コンピュータプログラムが提供される。コンピュータが使えるユーザであれば、このプログラムを使用することによって、専門家と同等の品質の流域推定をすることができる。
【発明を実施するための最良の形態】
【００１７】
最初に、本発明の実施形態の概要を説明する。本実施形態では、ＤＥＭデータの誤差を修正しながら流域を推定する。また、この処理をコンピュータで効率的に行うために、部分領域（ｓｕｂｄｏｍａｉｎ）を使用する。部分領域とは、複数の格子からなる領域である。本発明の処理の概要は以下の通りである。まず、現在流域とわかっている領域を含む部分領域との距離が所定の値より近い部分領域のＤＥＭデータを読み込む。次に、読み込んだＤＥＭデータを修正する。ＤＥＭデータの修正とは窪地を埋める処理である。この窪地を埋める処理を、窪地がなくなるまで繰り返す。ただし、窪地が所定の体積又は面積よりも大きい場合、その窪地は埋めない。また、すでにＤＥＭデータを読み込んである部分領域よりも外側に窪地が拡がる可能性がある場合、その窪地が拡がる可能性がある部分領域のＤＥＭデータが読み込まれるまでその窪地を埋める処理を保留する。続いて、現在流域とわかっている領域の１格子外側の格子が、流域か否かを判定する。その格子が流域であると判定された場合、その格子を流域に含める。この一連の処理を繰り返し行う。
【００１８】
以下、本発明の実施形態を添付図面に基づいて説明する。
【００１９】
図１は、本発明の実施形態におけるコンピュータシステムの構成を示す概略ブロック図である。
【００２０】
本実施の形態のコンピュータシステムは、ＰＣ（Ｐｅｒｓｏｎａｌ Ｃｏｍｐｕｔｅｒ）１００によって構成される。ＰＣ１００は、入出力部１１０及びＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ）１９０を備える。入出力部１１０は、キーボード１１１及びびマウス１１２などからなる入力部と、ディスプレイ１１３などからなる出力部とによって構成される。ユーザは、入力部から命令を入力し、出力部から結果を確認することができる。
【００２１】
さらにＰＣ１００は、メモリ１２０及びＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）１３０を備える。メモリ１２０には流域推定プログラム１４１が展開される。この流域推定プログラム１４１はＣＰＵへの命令によって構成される。ＣＰＵ１３０は、流域推定プログラム１４１の命令に従って、計算、ＨＤＤ１９０へのアクセス、及び入出力部１１０との情報の交換などを行う。以下の説明において流域推定プログラム１４１が実行する処理は、実際には、流域推定プログラム１４１の命令に従うＣＰＵ１３０によって実行される。
【００２２】
ＨＤＤ１９０には地図データ１９１及び標高データ１９２が格納される。流域推定プログラム１４１はこれらのデータを読み込むことができる。さらにＨＤＤ１９０には流向データ１９３を格納することができる。流域推定プログラム１４１はこのデータを書き込むことができる。なお、ＨＤＤ１９０は、ＰＣ１００に内蔵されてもよいし、ＰＣ１００の外部に接続されてもよいし、ＰＣ１００からネットワーク（図示省略）経由でアクセスできる別の装置に接続されてもよい。
【００２３】
地図データ１９１は、ユーザに場所を理解させるために必要な情報を含んでいれば、どのような種類のものであってもよい。地図データ１９１に適している情報として、道路、川、山頂などの位置、等高線、地名、衛星から撮影した地表面の画像、及び緯度経度を示す線などが挙げられる。本実施形態では、川の位置、山頂の位置、及び等高線を含む地図データ１９１が使用される。
【００２４】
標高データ１９２は、一般的にＤｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅｌ（ＤＥＭ）と呼ばれるデータである。すなわち、標高データ１９２は、地表面に区画された格子とその格子における地表面の標高とを記載した情報である。例えば、標高データ１９２は、地表面に等間隔で並べた正方形の格子（ｓｑｕａｒｅ ｇｒｉｄ）、あるいは不規則三角形ネットワーク（ｔｒｉａｎｇｕｌａｒ ｉｒｒｅｇｕｌａｒ ｎｅｔｗｏｒｋ）で表現されている格子における標高の値を含む。
【００２５】
流向データ１９３は、流域推定プログラム１４１によって得られた結果に関する情報である。この流向データ１９３については、後述する図６で詳細に説明する。
【００２６】
図２は、本発明の実施形態のＰＣ１００によって実行される処理を示すシーケンス図である。
【００２７】
以下、図２にしたがってシステムの処理シーケンスを説明する。
【００２８】
ＰＣ１００で流域推定プログラム１４１が起動すると、ＣＰＵ１３０は、流域推定プログラム１４１の命令に基づいて、ＨＤＤ１９０から地図データ１９１を読み出す（ステップ２０１）。以下、図２においてＣＰＵ１３０が実行する処理は、流域推定プログラム１４１の命令に基づいて実行される。
【００２９】
データを読み出したＣＰＵ１３０は、ウィンドウ３００からなるＧＵＩ（Ｇｒａｐｈｉｃａｌ Ｕｓｅｒ Ｉｎｔｅｒｆａｃｅ）をディスプレイ１１３に出力する（ステップ２０２）。
【００３０】
図３は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの入力待ちの状態を示す説明図である。
【００３１】
図３のＧＵＩは、ディスプレイ１１３の所定の領域に表示されるウィンドウ３００からなり、キーボード１１１及びマウス１１２による入力を受け付ける。
【００３２】
ウィンドウ３００は、少なくとも地図３１０を表示する。ユーザは、所定の操作をすることによって、地図３１０として表示される地域を所望の地域に変更することができる。
【００３３】
地図３１０の縮尺は、流出現象の空間スケールと同程度であることが望ましい。したがって典型的な流出現象であれば、１／２５，０００かそれより粗い縮尺が適当である。
【００３４】
本実施形態の地図３１０は、等高線３１１、川３１２及び山頂３１３の位置を示す情報を含む。
【００３５】
さらにウィンドウ３００は少なくとも、ユーザに１以上の地点を指定させる機能、及び、ユーザにその指定を確定させる機能を有する。後述するように、１以上の地点を指定するためにマウスポインタ３０１が、およびその指定を確定するためにＯＫボタン３０２が使用される。
【００３６】
このＧＵＩを用いてユーザは、興味がある地点を１点以上指定する（ステップ２０３）。興味がある地点を以後、下流端（ｔｅｒｍｉｎａｌ）と呼ぶ。例えば、ユーザは、ある地点に流入する水の流域を知りたい場合、その地点を下流端として指定する。本実施形態では、ユーザがマウス１１２を動かすことによってマウスポインタ３０１を地図３１０上で動かし、マウス１１２をクリックすることによって１地点を指定する。なお、指定した地点にアイコン３２１を表示することが望ましい。これによって、ユーザの操作結果をユーザにわかりやすく伝えることができる。
【００３７】
ユーザは、興味ある全ての地点を入力した後、ＯＫボタン３０２をクリックする。ＯＫボタン３０２がクリックされると、ＣＰＵ１３０はステップ２０３を終了する。図３は、ユーザが川付近の１地点をアイコン３２１によって下流端として指定した状態を示している。したがって、ＣＰＵ１３０はこれからこの１地点の上流域を解析する。
【００３８】
次に流域推定プログラム１４１は、流向算出処理（ステップ２０４）を繰り返し実行する。ステップ２０４の詳細は後述する図１１で説明する。
【００３９】
ステップ２０４は、標高データ１９２読み出し処理（ステップ２０５）及び処理済み領域表示処理（２０６）を含む。ステップ２０５では、ＣＰＵ１３０がＨＤＤ１９０から標高データ１９２を読み込む。ステップ２０６では、ＣＰＵ１３０がウィンドウ３００に処理済み領域を表示する。処理済み領域とは、流向算出処理が終了した領域である。
【００４０】
図４は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの解析中の状態を示す第１の説明図である。
【００４１】
図４の例では、処理済み領域を示す矩形４１３が地図３１０上に表示される。なお、ＣＰＵ１３０が流域推定プログラム１４１を処理中であることをマウスポインタのアイコン４０１によって表すことが望ましい。
【００４２】
図５は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの解析中の状態を示す第２の説明図である。
【００４３】
図５は、図４の状態のＧＵＩが表示された時点からしばらく時間が経った時点で表示される。つまり、図４の処理済み領域（矩形４１３で囲まれた領域）よりも図５の処理済み領域（矩形５１３で囲まれた領域）の方が広い。このように、ステップ２０４が繰り返し実行されている間に処理済み領域が拡大し、そのことがユーザに提示される。
【００４４】
ステップ２０４が終了した時点で、処理済み領域に含まれる各格子の流向が確定する。
【００４５】
次に、ＣＰＵ１３０は、確定した流向を示す情報をＨＤＤ１９０に流向データ１９３として書き込む（ステップ２０８）。
【００４６】
さらに、ＣＰＵ１３０は、図６に示すように、確定した流向をウィンドウ３００に表示する（ステップ２０９）。
【００４７】
図６は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの流域推定済みの状態を示す説明図である。
【００４８】
流向を示す線（流向ベクトル；ｆｌｏｗ ｄｉｒｅｃｔｉｏｎ ｖｅｃｔｏｒ）６１４が地図３１０上に表示される。それぞれの流向ベクトル６１４は、下流端（すなわちアイコン３２１によって指示された地点）に達することが確認できる。また、図５と図６を比較すると、図５の処理済み領域５１３の内に、流向ベクトル６１４が生成されることが確認できる。この流向ベクトル６１４が表示された領域が、下流端に流れ込む水の流域である。
【００４９】
図６には、ステップ２０９において流向を示す情報が表示される例を示した。しかし、ある格子が流域に含まれると判定された場合、その格子が流域に含まれることを示す情報のみが地図３１０上に表示されてもよい。少なくとも、ステップ２０９において、格子を識別する情報（例えば、格子の位置を示す情報）、及び、その格子が流域に含まれることを示す情報が出力される。
【００５０】
以上、図２にしたがってＰＣ１００によって実行される処理を示すシーケンスについて説明した。流域推定プログラム１４１は、空間をコンピュータで表現するために変数及び格子を利用する。また、流域推定プログラム１４１は、格子をコンピュータで効率的に扱うために部分領域を利用する。そこで、流域推定プログラム１４１が利用する格子、部分領域及び変数を説明する。
【００５１】
図７は、本発明の実施形態の流域推定プログラム１４１が利用する格子及び部分領域を示す説明図である。
【００５２】
本実施形態では、空間の東西方向（７０１）をｘ軸、南北方向（７０２）をｙ軸とする直交座標系（Ｃａｒｔｅｓｉａｎ ｃｏｏｒｄｉｎａｔｅ ｓｙｓｔｅｍ）を利用する。さらにこの空間をｘ軸及びｙ軸に平行な格子線（図中の細線、例えば格子線７０３）で分割する。格子線のｘ軸方向の間隔（７０４）をｄｘ、ｙ軸方向の間隔（７０５）をｄｙとする。格子線で囲まれた正方形（又は長方形）の領域が格子である。
【００５３】
格子は、その格子の位置を示す情報（例えば、座標（ｉ，ｊ））によって識別されてもよい。地点（ｘ，ｙ）を含む格子（ｉ，ｊ）は以下の式で求められる。
【００５４】
i=floor(x/dx), j=floor(y/dy) ・・・（１）
ただし、ｆｌｏｏｒ（）は、括弧内の値の小数部分を切り捨てることによって、括弧内の値の整数部分のみを取得する関数である。
【００５５】
さらに本実施形態では、空間をｘ軸及びｙ軸に平行な部分領域線（図中の太線、例えば部分領域線７０６）で分割する。部分領域線のｘ軸方向の間隔（７０７）をｓｉ＊ｄｘ、ｙ軸方向の間隔（７０８）をｓｊ＊ｄｙとする。ただし、ｓｉ及びｓｊは１以上の整数である。それぞれの格子は１つの部分領域に属することになる。
【００５６】
格子（ｉ，ｊ）が属する部分領域（ｂｉ，ｂｊ）は以下の式で求められる。
【００５７】
bi=floor(i/si), bj=floor(j/sj) ・・・（２）
部分領域（ｂｉ，ｂｊ）に属する格子は以下の式を満たす格子（ｉ，ｊ）である。
【００５８】
bi*si<=i<(bi+1)*si, bj*sj<=j<(bj+1)*sj ・・・（３）
それぞれの格子には、図８に示す変数が割り当てられる。
【００５９】
図８は、本発明の実施形態の格子（ｉ，ｊ）に割り当てられる変数の説明図である。
【００６０】
各格子には、少なくとも図８に示した情報、すなわち、標高（８１０）、最大勾配方向配列（８２０）、流向候補配列（８３０）及び流域フラグ（８５０）が割り当てられる。言い換えると、格子ごとに、図８に示す情報が生成され、メモリ１２０に保持される。
【００６１】
説明のために、それぞれの変数に名前をつける。標高（８１０）は変数Ｇｌｅｖ、最大勾配方向配列（８２０）は配列ＤｏｗｎＤｉｒ０、流向候補配列（８３０）は配列ＤｏｗｎＤｉｒ１、流域フラグ（８５０）は変数ＩｓＢａｓｉｎにそれぞれ格納される。
【００６２】
なお、上記のような変数が格納される代わりに、図８で示した各変数と同等の情報がコンピュータに格納されてもよい。たとえば、変数に別の名前をつけてもよい。
【００６３】
標高（８１０）は、格子の標高値を格納するための変数である。
【００６４】
格子ごとの標高値に基づいて、格子ごとの勾配を算出することができる。
【００６５】
図９は、本発明の実施形態において算出される勾配の説明図である。
【００６６】
格子（ｉ，ｊ）における勾配は、格子（ｉ，ｊ）の標高値とその格子に隣接する格子における標高値との差を格子間の距離で割った値である。１格子に隣接する格子は８ある。したがって、格子（ｉ，ｊ）における勾配の方向は８ある。すなわち、東、北東、北、北西、西、南西、南、南東である。勾配の方向をそれぞれ方向９０１から方向９０８の矢印で示す。説明のために、格子９００における８方向の勾配を、変数Ｇｒａｄ１からＧｒａｄ８と記載する。
【００６７】
変数Ｇｒａｄ１からＧｒａｄ８の値は、以下式によって求められる。
【００６８】
Grad1=(Glev[i,j]-Glev[i-1,j+1])/hypot(dx,dy) ・・・（４−１）
Grad2=(Glev[i,j]-Glev[i,j+1])/dy ・・・（４−２）
Grad3=(Glev[i,j]-Glev[i+1,j+1])/hypot(dx,dy) ・・・（４−３）
Grad4=(Glev[i,j]-Glev[i-1,j])/dx ・・・（４−４）
Grad5=(Glev[i,j]-Glev[i+1,j])/dx ・・・（４−５）
Grad6=(Glev[i,j]-Glev[i-1,j-1])/hypot(dx,dy) ・・・（４−６）
Grad7=(Glev[i,j]-Glev[i,j-1])/dy ・・・（４−７）
Grad8=(Glev[i,j]-Glev[i+1,j-1])/hypot(dx,dy) ・・・（４−８）
ただし、Ｇｌｅｖ［ｉ，ｊ］は格子（ｉ，ｊ）における勾配、ｈｙｐｏｔ（ｘ，ｙ）は原点（０，０）と点（ｘ，ｙ）との距離を返す関数である。
【００６９】
上記の式によって算出された変数Ｇｒａｄ１からＧｒａｄ８の値は、格子（ｉ，ｊ）から、その格子の周囲に隣接する８個の格子へ向かう方向の下り勾配を示す。例えば、格子（ｉ，ｊ）から格子（ｉ−１，ｊ＋１）に向かう方向の勾配が下りである（すなわち、格子（ｉ，ｊ）の標高より格子（ｉ−１，ｊ＋１）の標高が低い）場合、変数Ｇｒａｄ１の値は正となる。そして、その下り勾配が急である（すなわち、格子（ｉ，ｊ）の標高と格子（ｉ−１，ｊ＋１）の標高との差が大きい）ほど、変数Ｇｒａｄ１の値は大きくなる。
【００７０】
なお、本実施形態では、格子と、その格子に隣接する８格子との間の勾配が計算される。しかし、格子と、その格子の東西方向及び南北方向に隣接する４格子との間の勾配が計算されてもよい。その場合は、上記８方向の勾配変数のうちＧｒａｄ２、Ｇｒａｄ４、Ｇｒａｄ５及びＧｒａｄ７のみが使用される。ただし、勾配が計算される隣接格子の数が増えるほど、精度のよい結果が期待できる。
【００７１】
図８に示した最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）は、これら８方向の変数（Ｇｒａｄ１からＧｒａｄ８）のうち、値が最大である変数に対応する方向を保持する配列である。すなわち、配列ＤｏｗｎＤｉｒ０に格納できる要素は、東、北東、北、北西、西、南西、南及び南東の８方向のうち少なくとも一つを識別するための情報であり、要素の数（配列サイズ）は０から８のいずれかである。この配列の初期値は、要素数０を示す状態であるｅｍｐｔｙである。一つの格子について計算されたＧｒａｄ１からＧｒａｄ８のうち、値が最大であるものが複数存在する場合、複数の方向を識別する値が最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に保持される。
【００７２】
流向候補配列ＤｏｗｎＤｉｒ１（８３０）は、最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に格納された要素のうち、さらに別の基準を満たす要素のみを保持する配列である。この配列の初期値は、要素数０を示す状態であるｅｍｐｔｙである。流向候補配列ＤｏｗｎＤｉｒ１（８３０）に保持される要素については、後に図１５を参照して詳述する。
【００７３】
流域フラグ（８５０）は、その格子が流域に含まれるか否かを判別した結果を格納する変数である。この変数の値は、ｔｒｕｅ又はｆａｌｓｅのいずれかとなる。ｔｒｕｅは格子が流域に含まれることを、ｆａｌｓｅは格子が流域に含まれないことを示す。この変数の初期値が、ｆａｌｓｅである。
【００７４】
それぞれの格子に対応する、図８に示す変数が、メモリ１２０に展開される。格子の数は、極めて多くなる場合がある。その場合、メモリ１２０に展開される変数が非常に多くなる。変数の展開を効率的に行うためには、変数の実装をコンピュータの構造に適した設計にする必要がある。本発明では、格子をコンピュータで効率的に扱うために、部分領域にゴーストセル（ｇｈｏｓｔ ｃｅｌｌ）を追加する。
【００７５】
図１０は、本発明の実施形態のメモリ１２０に展開された状態の部分領域を説明するための図である。
【００７６】
部分領域には、ゴーストセルが付加される。図においてゴーストセルは斜線でハッチングされた格子である。一方、ドットで網掛けされた格子は、部分領域に含まれる格子（すなわち内部セル）である。
【００７７】
ある部分領域に付加されたゴーストセルは、その部分領域に隣接する部分領域に属する内部セルのコピーである。図１０には、（ｂｉ，ｂｊ）＝（０，０）を満たす部分領域１０１０及び（ｂｉ，ｂｊ）＝（１，０）を満たす部分領域１０２０が書かれている。部分領域１０２０に付加されたゴーストセル１０２１は、部分領域１０１０に含まれる内部セル１０１１のコピーである。同様に、部分領域１０１０に付加されたゴーストセル１０１２は、部分領域１０２０に含まれる内部セル１０２２のコピーである。このように、全てのゴーストセルは対応する内部セルのコピーである。
【００７８】
部分領域（ｂｉ，ｂｊ）に含まれる内部セルは、前述したとおり（３）式を満たす格子である。
【００７９】
部分領域（ｂｉ，ｂｊ）に付加されたゴーストセルは、（３）式を満たさず、かつ以下の（５）式を満たす格子である。
【００８０】
bi*si-mi<=I<(bi+1)*si+mi, bj*sj-mj<=j<(bj+1)*sj+mj ・・・（５）
ｍｉ及びｍｊは１以上の整数である。これらの数は、部分領域の内部セルにおける変数を求めるために参照される必要がある格子が、ゴーストセルとしてその部分領域に付加されるように定められる。本実施形態においては、ｍｉ＝１，ｍｊ＝１とした。これは、ある格子の勾配を計算するために、その格子に隣接する８格子の標高値を参照する必要があるためである。ゴーストセルの標高値を参照しなければ、部分領域の端に位置する内部セルの勾配を算出することができない。しかし、ｍｉ＝１，ｍｊ＝１とすれば、一つの部分領域の内部セル及びその部分領域に付加されたゴーストセルの標高値に基づいて、その部分領域の全内部セルの勾配を計算することができる。
【００８１】
このように、部分領域にゴーストセルが付加される場合、部分領域の内部セルを処理するために必要な情報が全てその部分領域及びゴーストセルに含まれる。このことは、複数の部分領域を並列処理することを可能にする。したがって、ゴーストセルが付加された部分領域を利用することで、プログラムの効率が向上することが期待できる。
【００８２】
以上、図７から図１０を使って、流域推定プログラム１４１の処理で使用される格子、部分領域、及び変数を説明した。続いて、流域推定プログラム１４１における流向算出処理（ステップ２０４）の詳細を説明する。
【００８３】
図１１は、本発明の実施形態の流域推定プログラム１４１の流向算出処理（ステップ２０４）を示すフローチャートである。
【００８４】
まず、ステップ１１０１において、ＣＰＵ１３０は、下流端を含む格子を現在の検査対象として指定する。下流端は地図３１０に入力される（図３参照）。ＣＰＵ１３０は、下流端の位置情報を地図３１０における座標系から直交座標系に変換する。続いてＣＰＵ１３０は、（１）式によって直交座標（ｘ，ｙ）を格子（ｉ，ｊ）に変換する。ＣＰＵ１３０は、この格子（ｉ，ｊ）に関する情報を現在の検査対象としてメモリ１２０に記憶する。現在の検査対象を記憶するための配列を、ここではＣｕｒＣｅｌｌｓと記載する。
【００８５】
ステップ１１０２において、ＣＰＵ１３０は、現在の検査対象の格子を含む部分領域を読み込む。このとき、図１０に示すように、ゴーストセルも読み込まれ、部分領域に付加される。このステップの詳細を後述する図１２及び図１３で説明する。
【００８６】
ステップ１１０３において、ＣＰＵ１３０は、処理済み領域をウィンドウ３００に出力する（図４、図５参照）。この処理は、すでに図２で説明したステップ２０６に相当する。
【００８７】
ステップ１１０４において、ＣＰＵ１３０は、部分領域に付加されたゴーストセルの同期を行う。このステップ１１０４の詳細を後述する図１６で説明する。
【００８８】
ステップ１１０５において、ＣＰＵ１３０は、現在の検査対象の上流の格子を次の検査対象として記憶する。具体的には、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌｓに格納された格子に隣接し、かつその格子の上流に位置する格子を次の検査対象としてメモリ１２０に記憶する。次の検査対象を記憶するための配列を、ここではＮｅｘｔＣｅｌｌｓと記載する。このステップ１１０５の詳細を後述する図１５で説明する。
【００８９】
ステップ１１０６において、ＣＰＵ１３０は、次の検査対象を現在の検査対象として新たに指定する。具体的には、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌｓを空にした後、配列ＮｅｘｔＣｅｌｌｓの全ての要素を配列ＣｕｒＣｅｌｌｓに格納する。さらに、ＣＰＵ１３０は、配列ＮｅｘｔＣｅｌｌｓを空にする。
【００９０】
ステップ１１０７において、ＣＰＵ１３０は、現在の検査対象の要素数が０より大か否かを判定する。具体的には、配列ＣｕｒＣｅｌｌｓの要素数と０とが比較される。
【００９１】
ステップ１１０７の結果が真である場合、現在の検査対象の格子の上流の格子が存在する可能性がある。この場合、処理はステップ１１０２に戻る。一方、ステップ１１０７の結果が偽である場合、前回のステップ１１０５において、配列ＮｅｘｔＣｅｌｌｓに追加された格子が存在しない。すなわち、現在までに検出された流域の上流の格子が存在しない。この場合、上流の探索が終了し、処理は図２に示したステップ２０８に進む。
【００９２】
このように、ステップ１１０２からステップ１１０７の処理によって、現在の検査対象の格子に隣接し、かつ上流にある格子が検査対象として指定される。また、ステップ１１０５において現在の検査対象の格子の上流の格子が存在する限りステップ１１０２からステップ１１０７が繰り返される。したがって、図１１に示したステップによって、一つの流域に含まれる下流端から上流までに位置する全ての格子が、必ず一度検査対象として指定される。
【００９３】
図１２は、本発明の実施形態の部分領域読み込み処理（ステップ１１０２）の詳細を示す第１のフローチャートである。
【００９４】
ステップ１２０１において、ＣＰＵ１３０は、現在の検査対象として指定された全ての格子についてステップ１２０２からステップ１２０８までの処理を繰り返す。具体的には、ＣＰＵ１３０は、現在の検査対象として指定された格子のうち一つを選択し、選択された格子についてステップ１２０２から１２０８までの処理を実行する。この選択された格子は、図１２の説明において対象格子と記載される。現在の検査対象として指定された全ての格子についてステップ１２０２から１２０８までの処理が終了するまで、上記の選択及び処理が繰り返される。
【００９５】
ステップ１２０２において、ＣＰＵ１３０は、対象格子を含む部分領域の識別子（部分領域ＩＤ）を算出する。部分領域ＩＤとして、式（２）によって算出される部分領域の座標（ｂｉ，ｂｊ）が使用される。格子（ｉ，ｊ）を含む部分領域座標（ｂｉ，ｂｊ）の算出方法は式（２）に示した。算出された部分領域ＩＤは、図１２の説明において対象部分領域ＩＤと記載される。対象部分領域ＩＤによって識別される部分領域は、図１２の説明において対象部分領域と記載される。
【００９６】
ステップ１２０３において、ＣＰＵ１３０は、対象部分領域ＩＤが配列Ｄｏｍａｉｎに含まれているか否かを判定する。この配列Ｄｏｍａｉｎは、ステップ１２０４以下のステップが実行された部分領域ＩＤを記憶するための配列であり、要素数は０以上である。ステップ１２０３の結果が真であれば、対象部分領域について、既にステップ１２０４以降の処理が実行されている。この場合、処理はステップ１２０９に進む。ステップ１２０３の結果が偽であれば、処理はステップ１２０４以下に進む。
【００９７】
ステップ１２０４において、ＣＰＵ１３０は、対象部分領域ＩＤを配列Ｄｏｍａｉｎに追加する。
【００９８】
ステップ１２０５において、ＣＰＵ１３０は、対象部分領域から所定の距離内にある部分領域をメモリ１２０に展開する。したがって、少なくとも対象部分領域と、対象部分領域に隣接する部分領域がメモリ１２０に展開される。具体的には、それらの部分領域に属する内部セル及びそれらに付加されたゴーストセルのデータ（図８参照）がメモリ１２０に格納される。対象部分領域に隣接する部分領域とは、東西方向３列、南北方向３列の９個の部分領域である。それらの９個の部分領域の中心が対象部分領域である。
【００９９】
なお、後述するように、窪地を埋める処理によって、格子の標高値が変更される場合がある。しかし、対象部分領域から上記の所定の距離より遠くにある窪地を埋めても、対象部分領域の標高は変更されない。これは、窪地を埋める処理において、窪地の面積が所定の値を超えた場合には窪地の標高値が変更されないためである。この処理については後述する（図１４参照）。
【０１００】
上記の所定の距離は、その土地にあった値に設定されることが望ましい。たとえば降雨量が多く土地の浸食の激しい日本では、この所定の距離を１ｋｍ程度にするのが適当と考えられる。逆に、降雨量が少なく土地の浸食が緩やかな米国中西部では、この所定の距離を１０ｋｍ程度にするのが適当と考えられる。部分領域の詳細は上記図１０にて説明した。
【０１０１】
ステップ１２０６において、ＣＰＵ１３０は、標高データ１９２を読み出す。この処理はすでに図２でステップ２０５として説明した。ＣＰＵ１３０は、読み出した標高データのうち、ステップ１２０５において展開された部分領域に属する格子の標高データを、その格子の標高Ｇｌｅｖ（８１０）に代入する。
【０１０２】
ステップ１２０７において、ＣＰＵ１３０は、最大勾配方向を算出する。勾配についてはすでに図９で説明した。各方向の勾配は、式（４−１）から式（４−８）によって算出される。ＣＰＵ１３０は、算出した８方向の勾配の最大値に対応する勾配方向を最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に格納する。
【０１０３】
ステップ１２０８において、ＣＰＵ１３０は、最大勾配の値が０以下である格子を初期値として、窪地データを生成する。さらに、ＣＰＵ１３０は、生成した窪地データを配列Ｓｉｎｋｓに格納する。ここで、配列Ｓｉｎｋｓは窪地データを記憶するための配列である。
【０１０４】
図１７は、本発明の実施形態の窪地データの詳細を示す説明図である。
【０１０５】
窪地データは、現在の標高（１７１０）、新しい標高（１７２０）、格子配列（１７３０）、隣接格子配列（１７４０）及び窪地埋め終了フラグ（１７５０）からなる。説明のため、これらをそれぞれＧｌｅｖＣｕｒ、ＧｌｅｖＮｅｗ、Ｃｅｌｌｓ、Ｎｅｉｂｓ及びＩｓＦｉｎｉｓｈと記載する。
【０１０６】
格子配列Ｃｅｌｌｓは、窪地に含まれる格子を識別する情報を格納する配列である。
【０１０７】
隣接格子配列Ｎｅｉｂｓは、窪地に含まれる格子に隣接する格子を識別する情報を格納する配列である。より詳細には、隣接格子配列Ｎｅｉｂｓには、窪地に含まれる各格子に隣接する全ての格子を識別する情報のうち、格子配列Ｃｅｌｌｓに含まれないものが格納される。
【０１０８】
格子配列Ｃｅｌｌｓ及び隣接格子配列Ｎｅｉｂｓの初期値は、要素のないｅｍｐｔｙ状態である。この窪地データの初期化は、１つの格子ｔｈｅＣｅｌｌを与えることによって実行される。ここで与えられるｔｈｅＣｅｌｌは、窪地の底に相当する格子である。
【０１０９】
現在の標高ＧｌｅｖＣｕｒは、窪地内の格子の標高値の候補を格納する変数である。現在の標高ＧｌｅｖＣｕｒの初期値は、その格子ｔｈｅＣｅｌｌの標高である。
【０１１０】
新しい標高ＧｌｅｖＮｅｗは、窪地内の格子の標高値の候補となる可能性がある値を一時的に格納する変数である。
【０１１１】
格子配列Ｃｅｌｌｓには、初期値として、その格子ｔｈｅＣｅｌｌが格納される。
【０１１２】
窪地埋め終了フラグＩｓＦｉｎｉｓｈの初期値は、窪地を埋める処理（図１４参照）が終了していないことを示す偽（ｆａｌｓｅ）である。この値は、窪地を埋める処理が終了すると、真（ｔｒｕｅ）に更新される。
【０１１３】
図１７に示すデータの一組が、配列Ｓｉｎｋｓの一要素に相当する。
【０１１４】
再び図１２に戻り説明する。ステップ１２０９において、ステップ１２０１のループが完了するまで処理はステップ１２０１に戻る。ステップ１２０１のループが完了した場合、処理は図１３に示したステップ１３０１に進む。
【０１１５】
図１２の例では、ステップ１２０５において、対象格子を含む部分領域及びその周囲の一つ以上の部分領域のデータがメモリ１２０に格納される。しかし、図１２は一例であり、メモリ１２０に格納されるデータ量がより多く（又はより少なく）なるように制御されてもよい。
【０１１６】
少なくとも、対象格子に流入する水の流向を算出するためには、対象格子の周囲に隣接する８個の格子の最大勾配方向を算出する必要がある。対象格子が格子（ｉ，ｊ）である場合、対象格子の周囲に隣接する８個の格子とは、ｘ＝ｉ−１、ｉ又はｉ＋１、ｙ＝ｊ−１、ｊ又はｊ＋１を満たす９個の格子（ｘ、ｙ）のうち、対象格子を除く８個の格子である。
【０１１７】
対象格子の周囲に隣接する８個の格子の最大勾配方向を算出するためには、対象格子の標高値、及び、対象格子の周囲に隣接する８個の格子の標高値に加えて、対象格子の周囲に隣接する８個の格子のさらに周囲に隣接する１６個の格子の標高値を取得する必要がある。対象格子が格子（ｉ，ｊ）である場合、上記の１６個の格子は、ｘ＝ｉ−２、ｉ−１、ｉ、ｉ＋１又はｉ＋２、ｙ＝ｊ−２、ｊ−１、ｊ、ｊ＋１又はｊ＋２を満たす２５個の格子（ｘ、ｙ）のうち、対象格子及び対象格子に隣接する８個の格子を除く１６個の格子である。
【０１１８】
このため、ステップ１２０５において、少なくとも、対象格子、対象格子の周囲に隣接する８個の格子、及び、それらの８個の格子の周囲に隣接する１６個の格子からなる２５個の格子のデータがメモリ１２０に格納される必要がある。
【０１１９】
図１３は、本発明の実施形態の部分領域読み込み処理（ステップ１１０２）の詳細を示す第２のフローチャートである。
【０１２０】
ステップ１３０１において、ＣＰＵ１３０は、配列Ｓｉｎｋｓ内の窪地データについてステップ１３０２からステップ１３０５までの処理を繰り返す。
【０１２１】
ステップ１３０２において、ＣＰＵ１３０は、窪地を埋める処理を実行する。この処理の詳細を後述する図１４にて説明する。窪地を埋める処理によって、窪地を埋める処理が終了したことを示す窪地埋め終了フラグＩｓＦｉｎｉｓｈ（１７５０）が変化する。
【０１２２】
ステップ１３０４において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈが真であるか否かを判定する。窪地埋め終了フラグＩｓＦｉｎｉｓｈが真の場合、窪地を埋める処理が終了している。この場合、処理はステップ１３０５に進む。一方、窪地埋め終了フラグＩｓＦｉｎｉｓｈが偽の場合、窪地を埋める処理が終了していない。この場合、処理はステップ１３０６に進む。
【０１２３】
ステップ１３０５において、ＣＰＵ１３０は、窪地を埋める処理が終了した窪地データを配列Ｓｉｎｋから削除する。
【０１２４】
ステップ１３０６において、ＣＰＵ１３０は、ステップ１３０１のループが完了するまでステップ１３０１に戻る。
【０１２５】
ステップ１３０７において、ＣＰＵ１３０は、ステップ１２０７と同じ処理（最大勾配方向の算出）を再び実行する。このステップ１３０７でステップ１２０７と同じ処理を再び行う理由は、ステップ１３０２の窪地埋め処理によって格子の標高値が変更されたためである。すなわち、ステップ１３０７において、変更された標高値を用いて最大勾配方向が算出される。
【０１２６】
図１４は、本発明の実施形態の窪地を埋める処理（ステップ１３０２）の詳細を示すフローチャートである。
【０１２７】
まず、ステップ１４０１において、ＣＰＵ１３０は、配列Ｃｅｌｌｓ内の格子に隣接する格子を配列Ｎｅｉｂｓに格納する。
【０１２８】
ステップ１４０２において、ＣＰＵ１３０は、窪地を埋めるために必要な標高データがメモリに展開されているか否かを判定する。具体的には、ＣＰＵ１３０は、配列Ｎｅｉｂｓ内の全ての格子がメモリ１２０に展開されているか否かを判定する。これは、配列Ｎｅｉｂｓ内の格子（ｉ，ｊ）を含む部分領域（ｂｉ，ｂｊ）を（２）式によって求め、その値が配列Ｄｏｍａｉｎに含まれているか否かを判定することによって実行される。
【０１２９】
ステップ１４０２の判定の結果が真であれば、この処理をさらに進めるためステップ１４０３に進む。判定の結果が偽であれば、窪地埋め処理を実行するために必要な標高データの少なくとも一部がまだメモリ１２０に展開されていない。この場合、この窪地埋め処理を保留するため、処理は図１３のステップ１３０４に進む。
【０１３０】
ステップ１４０３において、ＣＰＵ１３０は、配列Ｎｅｉｂｓ内の格子で最低の標高値を変数ＧｌｅｖＮｅｗに格納する。
【０１３１】
ステップ１４０４において、ＣＰＵ１３０は、窪地を埋める処理が終了したか否かを判定する。これは、条件式ＧｌｅｖＮｅｗ＞ＧｌｅｖＣｕｒによって判定される。この条件式が真であれば、窪地を埋める処理は終了していないためステップ１４０８に進む。偽であれば窪地を埋める処理は終了したためステップ１４０５に進む。
【０１３２】
ステップ１４０５において、ＣＰＵ１３０は、配列Ｃｅｌｌｓ内の格子の標高をＧｌｅｖＣｕｒに書き換える。
【０１３３】
ステップ１４０６において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈをｔｒｕｅにし、この窪地埋め処理を終了するため図１３のステップ１３０４に進む。
【０１３４】
ステップ１４０８において、ＣＰＵ１３０は、現在の標高ＧｌｅｖＣｕｒに新しい標高ＧｌｅｖＮｅｗの値を代入する。
【０１３５】
ステップ１４０９において、ＣＰＵ１３０は、変数ＧｌｅｖＣｕｒの値と等しい標高の格子を配列Ｃｅｌｌｓに追加する。
【０１３６】
ステップ１４１０において、ＣＰＵ１３０は、窪地の大きさ（例えば面積）が所定の基準値よりも小さいか否かを判定する。窪地の面積が所定の基準値よりも大きい場合、その窪地は、標高の誤差に起因する窪地ではなく、実際にその場所に存在する窪地であるとみなされる。この基準値は、その土地にあった値に設定するとよい。たとえば降雨量が多く土地の浸食の激しい日本ではこの基準値を１平方キロメートル程度にするのが適当と考えられる。逆に、降雨量が少なく土地の浸食が緩やかな米国中西部では、この基準値を１０平方キロメートル程度にするのが適当と考えられる。
【０１３７】
このステップ１４１０は、実際には、Ｃｅｌｌｓの要素数が、上記の所定の基準値を格子の面積で割った値より小さいか否かを判定することによって実行されてもよい。あるいは、ステップ１４１０の判定基準として窪地の面積ではなく、窪地の体積が用いられてもよい。ただし、判定基準として面積を使った方が、ＣＰＵ１３０の負荷は軽くなる。
【０１３８】
窪地の体積は、窪地に含まれる各格子の深さを格子の面積に乗算した値を合計することによって算出される。各格子の深さは、各格子の現在の標高ＧｌｅｖＣｕｒからその格子の標高Ｇｌｅｖを減算することによって算出される。
【０１３９】
このステップ１４１０の判別結果が真であれば、窪地は地形データの誤差に起因すると見なされるため、処理はステップ１４０１に戻る。ステップ１４１０の判別結果が偽であれば、窪地は実際の地形に起因すると見なされるため、処理はステップ１４１１に進む。
【０１４０】
ステップ１４１１において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈをｔｒｕｅにし、この窪地埋め処理を終了するため図１３のステップ１３０４に進む。
【０１４１】
上記のように、ステップ１４０４の判定結果が偽（すなわちＮｏ）であった場合、窪地内の格子の標高が変更され、確定する。ステップ１４１０の判定結果が偽であった場合、窪地内の格子の標高が変更されないことが確定する。
【０１４２】
一方、ステップ１４０２の判定結果が偽であった場合、窪地埋め処理が保留される。すなわち、この場合、窪地内の格子の標高がまだ確定しない。その後、図１２のステップ１２０５において、新たな部分領域がメモリ１２０に展開されると、ステップ１４０２の判定結果が真（すなわちＹｅｓ）となり、窪地内の格子の標高を確定させるための処理が続行される。
【０１４３】
しかし、ステップ１３０７が実行される時点で、最大勾配方向が算出される格子及びその格子に隣接する格子の標高が確定している必要がある。言い換えると、ステップ１３０１から１３０６までのループが終了した時点で、少なくとも、現在の検査対象の格子を含む部分領域内の全ての格子の標高が確定している必要がある。この確定のために、ステップ１２０５において、ステップ１４１０の基準値が示す面積より十分に広い範囲の部分領域がメモリ１２０に展開される必要がある。
【０１４４】
具体的には、ステップ１２０５において、対象部分領域内の格子から半径ｒ（ｍ）以内の全ての格子がメモリ１２０に展開される必要がある。ここで、ステップ１４１０の基準値がＡ（平方ｍ）、一格子が１辺ｄｘ（ｍ）の正方形である(すなわちｄｘとｄｙが等しい)場合、半径ｒは、次の式（６）を満たす必要がある。
【０１４５】
ｒ≧Ａ／ｄｘ・・・（６）
上記の条件が満たされる場合、対象部分領域内の格子を含み、かつ、面積が基準値Ａ以内である窪地は、それがどのような形状であっても、ステップ１２０５において展開された部分領域の中に含まれる。したがって、この場合、対象部分領域の少なくとも一部の格子が窪地に含まれていても、その対象部分領域内の全ての格子の標高は、ステップ１３０７が実行される時点で全て確定している。
【０１４６】
図１５は、本発明の実施形態の上流探索処理（ステップ１１０５）の詳細を示すフローチャートである。
【０１４７】
ステップ１５０１において、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌ内の格子について以降の処理を繰り返す。図１５において、ステップ１５０１から１５０７までのループの処理対象の格子を、Ｃｅｌｌ０と記載する。
【０１４８】
ステップ１５０２において、ＣＰＵ１３０は、Ｃｅｌｌ０に隣接する８個の格子について以降の処理を繰り返す。図１５において、ステップ１５０２から１５０６までのループの処理対象の格子を、Ｃｅｌｌ１と記載する。したがって、Ｃｅｌｌ０とＣｅｌｌ１とは隣接している。
【０１４９】
ステップ１５０３において、ＣＰＵ１３０は、Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０があり、かつ、Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）がＦａｌｓｅであるか否かを判定する。Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０がある場合、そのＣｅｌｌ１から流出した水がＣｅｌｌ０に流入すると推定される。Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）がＦａｌｓｅである場合、そのＣｅｌｌ１がどの流域に属するかまだ確定していない。
【０１５０】
ステップ１５０３の判定の結果が真の場合、Ｃｅｌｌ１とＣｅｌｌ０とが一つの流域に含まれる。この場合、処理はステップ１５０４に進む。
【０１５１】
一方、判定の結果が偽の場合（すなわち、Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０がないか、又は、Ｃｅｌｌ０の流域フラグ（ＩｓＢａｓｉｎ）がＴｒｕｅである場合）、処理はステップ１５０６に進む。
【０１５２】
ステップ１５０４において、ＣＰＵ１３０は、Ｃｅｌｌ１の流向候補配列（ＤｏｗｎＤｉｒ１）にＣｅｌｌ０の方向を追加する。たとえば、Ｃｅｌｌ０がＣｅｌｌ１の北西方向にあれば、その方向が流向の候補として配列に追加される。
【０１５３】
ステップ１５０５において、ＣＰＵ１３０は、Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）を真にする。
【０１５４】
ステップ１５０６において、ＣＰＵ１３０は、ステップ１５０２のループが完了するまでステップ１５０２に戻る。
【０１５５】
ステップ１５０２からステップ１５０６までのループが完了した時点で、Ｃｅｌｌ１の流向候補配列（ＤｏｗｎＤｉｒ１）に含まれる流向候補が確定する。例えば、ステップ１５０４において、Ｃｅｌｌ０がＣｅｌｌ１の北西方向にある場合、その方向（北西方向）がＣｅｌｌ１の流向の候補として追加される。そして、ループが完了すると、この流向が確定する。流向の確定は、以下のようにすればよい。流向候補配列ＤｏｗｎＤｉｒ１に１つの方向が格納されている場合、その方向を流向として確定する。流向候補配列ＤｏｗｎＤｉｒ１に２以上の方向が格納されている場合、あらかじめ決めた方法にしたがって1つの方向を選択する。選択する方法は、乱数により決めてもいいし、あらかじめ決めた優先度（たとえば、東、北東、北、北西、西、南西、南、南東の順に優先して選択する）に従って決めてもよい。Ｃｅｌｌ１の水がＣｅｌｌ０に流出する場合、Ｃｅｌｌ０とＣｅｌｌ１とが一つの流域に含まれると判定される。このＣｅｌｌ１は、Ｃｅｌｌ０の次の検査対象として選択される。すなわち、このＣｅｌｌ１は、配列ＮｅｘｔＣｅｌｌｓに追加され、図１１のステップ１１０６において、現在の検査対象として新たに選択される。
【０１５６】
ステップ１５０７において、ＣＰＵ１３０は、ステップ１５０１のループが完了するまでステップ１５０１に戻る。
【０１５７】
図１６は、本発明の実施形態の部分領域に付加されたゴーストセルの同期処理（ステップ１１０４）の詳細を示すフローチャートである。
【０１５８】
ステップ１６０１において、ＣＰＵ１３０は、配列Ｄｏｍａｉｎ内の部分領域について以降の処理を繰り返す。
【０１５９】
ステップ１６０２において、ＣＰＵ１３０は、部分領域に付加されたすべてのゴーストセルについて以降の処理を繰り返す。
【０１６０】
ステップ１６０３において、ＣＰＵ１３０は、ゴーストセルに対応する内部セルの変数値をコピーする。例えば、図１０のゴーストセル１０１２についてステップ１６０３が実行される場合、ゴーストセル１０１２に対応する内部セル１０２２の変数値（すなわち、図８に示す各変数の値）がゴーストセル１０１２にコピーされる。その結果、内部セル１０２２の標高値が窪地埋め処理によって変更された場合、その変更がゴーストセル１０１２に反映される。
【０１６１】
ステップ１６０４において、ＣＰＵ１３０は、ステップ１６０２のループが完了するまでステップ１６０２に戻る。
【０１６２】
ステップ１６０５において、ＣＰＵ１３０は、ステップ１６０１のループが完了するまでステップ１６０１に戻る。
【０１６３】
以上述べたとおり、本発明の実施形態によれば、流域に関する事前知識を要求としない、流域推定コンピュータプログラムを提供することができる。このプログラムを使用することによって、コンピュータが使えるユーザであれば、専門家と同等の品質の流域推定が可能になる。
【図面の簡単な説明】
【０１６４】
【図１】本発明の実施形態におけるコンピュータシステムの構成を示す概略ブロック図である。
【図２】本発明の実施形態のＰＣによって実行される処理を示すシーケンス図である。
【図３】本発明の実施形態の流域推定プログラムが提供するＧＵＩの入力待ちの状態を示す説明図である。
【図４】本発明の実施形態の流域推定プログラムが提供するＧＵＩの解析中の状態を示す第１の説明図である。
【図５】本発明の実施形態の流域推定プログラムが提供するＧＵＩの解析中の状態を示す第２の説明図である。
【図６】本発明の実施形態の流域推定プログラムが提供するＧＵＩの流域推定済みの状態を示す説明図である。
【図７】本発明の実施形態の流域推定プログラムが利用する格子及び部分領域を示す説明図である。
【図８】本発明の実施形態の格子に割り当てられる変数の説明図である。
【図９】本発明の実施形態において算出される勾配の説明図である。
【図１０】本発明の実施形態のメモリに展開された状態の部分領域を説明するための図である。
【図１１】本発明の実施形態の流域推定プログラムの流向算出処理を示すフローチャートである。
【図１２】本発明の実施形態の部分領域読み込み処理の詳細を示す第１のフローチャートである。
【図１３】本発明の実施形態の部分領域読み込み処理の詳細を示す第２のフローチャートである。
【図１４】本発明の実施形態の窪地を埋める処理の詳細を示すフローチャートである。
【図１５】本発明の実施形態の上流探索処理の詳細を示すフローチャートである。
【図１６】本発明の実施形態の部分領域に付加されたゴーストセルの同期処理の詳細を示すフローチャートである。
【図１７】本発明の実施形態の窪地データの詳細を示す説明図である。
【図１８】ＤＥＭの誤差の影響を説明するための図である。
【符号の説明】
【０１６５】
１００ ＰＣ
３１０ 地図
３２１ 下流端を示すアイコン
４１３、５１３ 処理済み領域
６１４ 流向ベクトル
７０３ 格子線
７０６ 部分領域線
１０１０、１０２０ 部分領域
【技術分野】
【０００１】
本願明細書で開示される技術は、コンピュータを使って流域を推定する方法に関する。本願明細書で開示される技術は特に、１以上の地点の位置を入力データとし、それぞれの地点の流域を数値地形モデル（ｄｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅ，ＤＥＭ）に基づいて推定するコンピュータプログラムに関する。
【背景技術】
【０００２】
流域は、河川、ダム、下水などの水を管理する土木構造物にとって非常に重要な情報である。多くの場合、土木構造物に入ってくる水の起源は降雨である。つまり、降雨が地表面に落ち、地表面を流れ下って構造物に入る。したがって、構造物に至る全ての流れをさかのぼれば、降雨が落ちた地表面に達する。この地表面を構造物の流域という。
【０００３】
近年では、地表面の標高は電子データとして得られる。このような標高の電子データは、一般的にＤｉｇｉｔａｌ Ｅｌｅｖａｔｉｏｎ Ｍｏｄｅｌ （ＤＥＭ）と呼ばれる。ＤＥＭでは、ある領域の地表面が格子（ｇｒｉｄ又はｃｅｌｌ）と呼ばれる小領域に分割され、それぞれの格子に標高値が与えられる。多くの場合、ＤＥＭの格子は正方形である。
【０００４】
このＤＥＭを使い、土木構造物の流域を推定するコンピュータプログラムが作られている。通常、地表面に落ちた降雨は最も急な下り勾配（ｇｒａｄｉｅｎｔ）の方向に流れ下る。したがって、ＤＥＭの格子ごとの勾配を求めれば、流域を推定できる。勾配から流域を推定する方法は、非特許文献１に詳しい。
【０００５】
このような流域を推定するコンピュータプログラムに共通の課題は、ＤＥＭの誤差である。勾配はベクトル解析における勾配演算によって算出される。勾配演算は元データの誤差に敏感である。したがって、ＤＥＭデータから流域を推定する方法は誤差に敏感である。
【０００６】
図１８は、ＤＥＭの誤差の影響を説明するための図である。
【０００７】
グラフ１８００は、ある断面に沿った標高である。格子１８０１から格子１８０９までの格子に標高が与えられ、勾配が算出されている。勾配を矢印（たとえば１８１１）で示す。グラフ１８００では、格子１８０６の標高が、隣接する格子１８０５及び１８０７の標高よりも低い。下り勾配の方向（すなわち図の矢印の方向）と水が流れる方向とが一致すると仮定すると、格子１８０１から格子１８０７までに降った雨は格子１８０６に流れ出て、格子１８０６からは水が流れ出ないことになる。また、格子１８０９の流域は格子１８０７及び１８０８となる。この結果と、地表面を流れる水は最終的に標高の低い場所に至るという私たちの常識とは矛盾する。この矛盾の原因は、格子１８０６近傍の標高データの誤差にあると推測される。
【０００８】
誤差の修正方法として、いくつかの方法が提案されている。最も一般的な修正方法は、窪地埋め（ｓｉｎｋ ｆｉｌｌｉｎｇ）である。窪地とは１以上の格子からなる、周囲に比べ相対的に標高が低い領域である。窪地埋めの方法は非特許文献１に詳しく説明されている。この方法では、ある領域内の全ての格子のうち、ある格子の標高が、隣接するすべての格子の標高と同じ又は低い場合、ごく小さな値をその格子の標高に加える。この処理を繰り返すことによって修正が完了する。このように誤差が修正されたＤＥＭを元に推定した流域と、私たちの常識とは一致する。
【０００９】
図１８のグラフ１８２０は、グラフ１８００に示したＤＥＭデータをこの方法によって修正した結果である。格子１８２１から格子１８２９は、それぞれ、修正前のグラフ１８００の格子１８０１から格子１８０９に対応する。窪地（すなわち格子１８２５及び１８２６）が埋められた結果、格子１８２９の流域が格子１８２１から格子１８２９までとなった。
【００１０】
さらに別の誤差修正方法として、特許文献１に開示された「Ｍｅｔｈｏｄ ｆｏｒ ｆｌｏｗ ａｎａｌｙｓｉｓ ｕｔｉｌｉｚｉｎｇ ｃｅｌｌ ｅｌｅｖａｔｉｏｎ ｔｏｐｏｇｒａｐｈｙ」が挙げられる。この方法では、あらかじめ領域の内部が複数のｍａｊｏｒ ｃｅｌｌに分割される。さらに、ｍａｊｏｒ ｃｅｌｌの内部が、少なくとも横方向に３分割、縦方向に３分割することによって、少なくとも９個のｍｉｎｏｒ ｃｅｌｌに分割される。あるｍｉｎｏｒ ｃｅｌｌにおける勾配を算出するときは、そのｍｉｎｏｒ ｃｅｌｌを含むｍａｊｏｒ ｃｅｌｌと、さらにそのｍａｊｏｒ ｃｅｌｌに隣接するｍａｊｏｒ ｃｅｌｌについてＤＥＭデータのエラーが修正される。その後、これらのｍａｊｏｒ ｃｅｌｌに含まれるｍｉｎｏｒ ｃｅｌｌについて勾配が再度算出される。
【００１１】
非特許文献１では処理の対象とする領域全体に誤差修正処理をする必要があったが、特許文献１に開示された方法では処理の対象とする領域を複数に分割し、分割された領域ごとに勾配を求めることが可能になる。
【００１２】
いずれの誤差修正方法であっても、ＤＥＭデータの誤差を修正した後に流域を推定する。つまり、まずＤＥＭデータを修正する。次にＤＥＭデータから勾配データを算出する。この勾配データを使って上流方向にたどると、流域が求まる。
【特許文献１】米国特許出願公開第２００５／００４９８０２号明細書
【非特許文献１】Ｃｏｓｔａ−Ｃａｂｒａｌ，Ｍ．Ｃ．，ａｎｄ Ｓ．Ｊ．Ｂｕｒｇｅｓ（１９９４），Ｄｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅｌ ｎｅｔｗｏｒｋｓ（ＤＥＭＯＮ）： Ａ ｍｏｄｅｌ ｏｆ ｆｌｏｗ ｏｖｅｒ ｈｉｌｌｓｌｏｐｅｓ ｆｏｒ ｃｏｍｐｕｔａｔｉｏｎ ｏｆ ｃｏｎｔｒｉｂｕｔｉｎｇ ａｎｄ ｄｉｓｐｅｒｓａｌ ａｒｅａｓ，Ｗａｔｅｒ Ｒｅｓｏｕｒ．Ｒｅｓ．，３０（６），１６８１−１６９２．
【発明の開示】
【発明が解決しようとする課題】
【００１３】
従来の流域推定プログラムは、ＤＥＭデータの誤差を修正した後に、流域を推定する。したがって、流域を推定する前に、修正すべきＤＥＭデータの範囲を知っておく必要がある。ところが流域を推定したい場合、修正すべきＤＥＭデータの範囲は流域である。つまり従来の流域推定プログラムのユーザは、そのプログラムの実行結果である流域を、プログラムを実行する前に知っておく必要がある。したがって、従来の流域推定プログラムを使えるユーザは、流域をあらかじめ知っているユーザに限られていた。流域を知らないユーザが従来の流域推定プログラムを使用すると、ＤＥＭデータを適切に修正できないために、流域を正確に推定することが困難であった。
【００１４】
本発明はこの問題に鑑みてなされたものである。本発明の目的は、流域に関する事前知識のないユーザが使用しても、流域を正確に推定することができる流域推定コンピュータプログラムを提供することである。
【課題を解決するための手段】
【００１５】
本願で開示する代表的な発明は、計算機を制御するプログラムであって、前記計算機は、前記プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備え、前記プログラムは、地図上の地点を指示する入力を受け付ける第１手順と、前記地図を互いに直交して設定された第１軸及び第２軸に沿ってそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、前記第１領域内の窪地を検出する第４手順と、前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を前記プロセッサに実行させることを特徴とする。
【発明の効果】
【００１６】
本発明の一実施形態によれば、流域に関する事前知識を要求しない、流域推定コンピュータプログラムが提供される。コンピュータが使えるユーザであれば、このプログラムを使用することによって、専門家と同等の品質の流域推定をすることができる。
【発明を実施するための最良の形態】
【００１７】
最初に、本発明の実施形態の概要を説明する。本実施形態では、ＤＥＭデータの誤差を修正しながら流域を推定する。また、この処理をコンピュータで効率的に行うために、部分領域（ｓｕｂｄｏｍａｉｎ）を使用する。部分領域とは、複数の格子からなる領域である。本発明の処理の概要は以下の通りである。まず、現在流域とわかっている領域を含む部分領域との距離が所定の値より近い部分領域のＤＥＭデータを読み込む。次に、読み込んだＤＥＭデータを修正する。ＤＥＭデータの修正とは窪地を埋める処理である。この窪地を埋める処理を、窪地がなくなるまで繰り返す。ただし、窪地が所定の体積又は面積よりも大きい場合、その窪地は埋めない。また、すでにＤＥＭデータを読み込んである部分領域よりも外側に窪地が拡がる可能性がある場合、その窪地が拡がる可能性がある部分領域のＤＥＭデータが読み込まれるまでその窪地を埋める処理を保留する。続いて、現在流域とわかっている領域の１格子外側の格子が、流域か否かを判定する。その格子が流域であると判定された場合、その格子を流域に含める。この一連の処理を繰り返し行う。
【００１８】
以下、本発明の実施形態を添付図面に基づいて説明する。
【００１９】
図１は、本発明の実施形態におけるコンピュータシステムの構成を示す概略ブロック図である。
【００２０】
本実施の形態のコンピュータシステムは、ＰＣ（Ｐｅｒｓｏｎａｌ Ｃｏｍｐｕｔｅｒ）１００によって構成される。ＰＣ１００は、入出力部１１０及びＨＤＤ（Ｈａｒｄ Ｄｉｓｋ Ｄｒｉｖｅ）１９０を備える。入出力部１１０は、キーボード１１１及びびマウス１１２などからなる入力部と、ディスプレイ１１３などからなる出力部とによって構成される。ユーザは、入力部から命令を入力し、出力部から結果を確認することができる。
【００２１】
さらにＰＣ１００は、メモリ１２０及びＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）１３０を備える。メモリ１２０には流域推定プログラム１４１が展開される。この流域推定プログラム１４１はＣＰＵへの命令によって構成される。ＣＰＵ１３０は、流域推定プログラム１４１の命令に従って、計算、ＨＤＤ１９０へのアクセス、及び入出力部１１０との情報の交換などを行う。以下の説明において流域推定プログラム１４１が実行する処理は、実際には、流域推定プログラム１４１の命令に従うＣＰＵ１３０によって実行される。
【００２２】
ＨＤＤ１９０には地図データ１９１及び標高データ１９２が格納される。流域推定プログラム１４１はこれらのデータを読み込むことができる。さらにＨＤＤ１９０には流向データ１９３を格納することができる。流域推定プログラム１４１はこのデータを書き込むことができる。なお、ＨＤＤ１９０は、ＰＣ１００に内蔵されてもよいし、ＰＣ１００の外部に接続されてもよいし、ＰＣ１００からネットワーク（図示省略）経由でアクセスできる別の装置に接続されてもよい。
【００２３】
地図データ１９１は、ユーザに場所を理解させるために必要な情報を含んでいれば、どのような種類のものであってもよい。地図データ１９１に適している情報として、道路、川、山頂などの位置、等高線、地名、衛星から撮影した地表面の画像、及び緯度経度を示す線などが挙げられる。本実施形態では、川の位置、山頂の位置、及び等高線を含む地図データ１９１が使用される。
【００２４】
標高データ１９２は、一般的にＤｉｇｉｔａｌ ｅｌｅｖａｔｉｏｎ ｍｏｄｅｌ（ＤＥＭ）と呼ばれるデータである。すなわち、標高データ１９２は、地表面に区画された格子とその格子における地表面の標高とを記載した情報である。例えば、標高データ１９２は、地表面に等間隔で並べた正方形の格子（ｓｑｕａｒｅ ｇｒｉｄ）、あるいは不規則三角形ネットワーク（ｔｒｉａｎｇｕｌａｒ ｉｒｒｅｇｕｌａｒ ｎｅｔｗｏｒｋ）で表現されている格子における標高の値を含む。
【００２５】
流向データ１９３は、流域推定プログラム１４１によって得られた結果に関する情報である。この流向データ１９３については、後述する図６で詳細に説明する。
【００２６】
図２は、本発明の実施形態のＰＣ１００によって実行される処理を示すシーケンス図である。
【００２７】
以下、図２にしたがってシステムの処理シーケンスを説明する。
【００２８】
ＰＣ１００で流域推定プログラム１４１が起動すると、ＣＰＵ１３０は、流域推定プログラム１４１の命令に基づいて、ＨＤＤ１９０から地図データ１９１を読み出す（ステップ２０１）。以下、図２においてＣＰＵ１３０が実行する処理は、流域推定プログラム１４１の命令に基づいて実行される。
【００２９】
データを読み出したＣＰＵ１３０は、ウィンドウ３００からなるＧＵＩ（Ｇｒａｐｈｉｃａｌ Ｕｓｅｒ Ｉｎｔｅｒｆａｃｅ）をディスプレイ１１３に出力する（ステップ２０２）。
【００３０】
図３は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの入力待ちの状態を示す説明図である。
【００３１】
図３のＧＵＩは、ディスプレイ１１３の所定の領域に表示されるウィンドウ３００からなり、キーボード１１１及びマウス１１２による入力を受け付ける。
【００３２】
ウィンドウ３００は、少なくとも地図３１０を表示する。ユーザは、所定の操作をすることによって、地図３１０として表示される地域を所望の地域に変更することができる。
【００３３】
地図３１０の縮尺は、流出現象の空間スケールと同程度であることが望ましい。したがって典型的な流出現象であれば、１／２５，０００かそれより粗い縮尺が適当である。
【００３４】
本実施形態の地図３１０は、等高線３１１、川３１２及び山頂３１３の位置を示す情報を含む。
【００３５】
さらにウィンドウ３００は少なくとも、ユーザに１以上の地点を指定させる機能、及び、ユーザにその指定を確定させる機能を有する。後述するように、１以上の地点を指定するためにマウスポインタ３０１が、およびその指定を確定するためにＯＫボタン３０２が使用される。
【００３６】
このＧＵＩを用いてユーザは、興味がある地点を１点以上指定する（ステップ２０３）。興味がある地点を以後、下流端（ｔｅｒｍｉｎａｌ）と呼ぶ。例えば、ユーザは、ある地点に流入する水の流域を知りたい場合、その地点を下流端として指定する。本実施形態では、ユーザがマウス１１２を動かすことによってマウスポインタ３０１を地図３１０上で動かし、マウス１１２をクリックすることによって１地点を指定する。なお、指定した地点にアイコン３２１を表示することが望ましい。これによって、ユーザの操作結果をユーザにわかりやすく伝えることができる。
【００３７】
ユーザは、興味ある全ての地点を入力した後、ＯＫボタン３０２をクリックする。ＯＫボタン３０２がクリックされると、ＣＰＵ１３０はステップ２０３を終了する。図３は、ユーザが川付近の１地点をアイコン３２１によって下流端として指定した状態を示している。したがって、ＣＰＵ１３０はこれからこの１地点の上流域を解析する。
【００３８】
次に流域推定プログラム１４１は、流向算出処理（ステップ２０４）を繰り返し実行する。ステップ２０４の詳細は後述する図１１で説明する。
【００３９】
ステップ２０４は、標高データ１９２読み出し処理（ステップ２０５）及び処理済み領域表示処理（２０６）を含む。ステップ２０５では、ＣＰＵ１３０がＨＤＤ１９０から標高データ１９２を読み込む。ステップ２０６では、ＣＰＵ１３０がウィンドウ３００に処理済み領域を表示する。処理済み領域とは、流向算出処理が終了した領域である。
【００４０】
図４は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの解析中の状態を示す第１の説明図である。
【００４１】
図４の例では、処理済み領域を示す矩形４１３が地図３１０上に表示される。なお、ＣＰＵ１３０が流域推定プログラム１４１を処理中であることをマウスポインタのアイコン４０１によって表すことが望ましい。
【００４２】
図５は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの解析中の状態を示す第２の説明図である。
【００４３】
図５は、図４の状態のＧＵＩが表示された時点からしばらく時間が経った時点で表示される。つまり、図４の処理済み領域（矩形４１３で囲まれた領域）よりも図５の処理済み領域（矩形５１３で囲まれた領域）の方が広い。このように、ステップ２０４が繰り返し実行されている間に処理済み領域が拡大し、そのことがユーザに提示される。
【００４４】
ステップ２０４が終了した時点で、処理済み領域に含まれる各格子の流向が確定する。
【００４５】
次に、ＣＰＵ１３０は、確定した流向を示す情報をＨＤＤ１９０に流向データ１９３として書き込む（ステップ２０８）。
【００４６】
さらに、ＣＰＵ１３０は、図６に示すように、確定した流向をウィンドウ３００に表示する（ステップ２０９）。
【００４７】
図６は、本発明の実施形態の流域推定プログラム１４１が提供するＧＵＩの流域推定済みの状態を示す説明図である。
【００４８】
流向を示す線（流向ベクトル；ｆｌｏｗ ｄｉｒｅｃｔｉｏｎ ｖｅｃｔｏｒ）６１４が地図３１０上に表示される。それぞれの流向ベクトル６１４は、下流端（すなわちアイコン３２１によって指示された地点）に達することが確認できる。また、図５と図６を比較すると、図５の処理済み領域５１３の内に、流向ベクトル６１４が生成されることが確認できる。この流向ベクトル６１４が表示された領域が、下流端に流れ込む水の流域である。
【００４９】
図６には、ステップ２０９において流向を示す情報が表示される例を示した。しかし、ある格子が流域に含まれると判定された場合、その格子が流域に含まれることを示す情報のみが地図３１０上に表示されてもよい。少なくとも、ステップ２０９において、格子を識別する情報（例えば、格子の位置を示す情報）、及び、その格子が流域に含まれることを示す情報が出力される。
【００５０】
以上、図２にしたがってＰＣ１００によって実行される処理を示すシーケンスについて説明した。流域推定プログラム１４１は、空間をコンピュータで表現するために変数及び格子を利用する。また、流域推定プログラム１４１は、格子をコンピュータで効率的に扱うために部分領域を利用する。そこで、流域推定プログラム１４１が利用する格子、部分領域及び変数を説明する。
【００５１】
図７は、本発明の実施形態の流域推定プログラム１４１が利用する格子及び部分領域を示す説明図である。
【００５２】
本実施形態では、空間の東西方向（７０１）をｘ軸、南北方向（７０２）をｙ軸とする直交座標系（Ｃａｒｔｅｓｉａｎ ｃｏｏｒｄｉｎａｔｅ ｓｙｓｔｅｍ）を利用する。さらにこの空間をｘ軸及びｙ軸に平行な格子線（図中の細線、例えば格子線７０３）で分割する。格子線のｘ軸方向の間隔（７０４）をｄｘ、ｙ軸方向の間隔（７０５）をｄｙとする。格子線で囲まれた正方形（又は長方形）の領域が格子である。
【００５３】
格子は、その格子の位置を示す情報（例えば、座標（ｉ，ｊ））によって識別されてもよい。地点（ｘ，ｙ）を含む格子（ｉ，ｊ）は以下の式で求められる。
【００５４】
i=floor(x/dx), j=floor(y/dy) ・・・（１）
ただし、ｆｌｏｏｒ（）は、括弧内の値の小数部分を切り捨てることによって、括弧内の値の整数部分のみを取得する関数である。
【００５５】
さらに本実施形態では、空間をｘ軸及びｙ軸に平行な部分領域線（図中の太線、例えば部分領域線７０６）で分割する。部分領域線のｘ軸方向の間隔（７０７）をｓｉ＊ｄｘ、ｙ軸方向の間隔（７０８）をｓｊ＊ｄｙとする。ただし、ｓｉ及びｓｊは１以上の整数である。それぞれの格子は１つの部分領域に属することになる。
【００５６】
格子（ｉ，ｊ）が属する部分領域（ｂｉ，ｂｊ）は以下の式で求められる。
【００５７】
bi=floor(i/si), bj=floor(j/sj) ・・・（２）
部分領域（ｂｉ，ｂｊ）に属する格子は以下の式を満たす格子（ｉ，ｊ）である。
【００５８】
bi*si<=i<(bi+1)*si, bj*sj<=j<(bj+1)*sj ・・・（３）
それぞれの格子には、図８に示す変数が割り当てられる。
【００５９】
図８は、本発明の実施形態の格子（ｉ，ｊ）に割り当てられる変数の説明図である。
【００６０】
各格子には、少なくとも図８に示した情報、すなわち、標高（８１０）、最大勾配方向配列（８２０）、流向候補配列（８３０）及び流域フラグ（８５０）が割り当てられる。言い換えると、格子ごとに、図８に示す情報が生成され、メモリ１２０に保持される。
【００６１】
説明のために、それぞれの変数に名前をつける。標高（８１０）は変数Ｇｌｅｖ、最大勾配方向配列（８２０）は配列ＤｏｗｎＤｉｒ０、流向候補配列（８３０）は配列ＤｏｗｎＤｉｒ１、流域フラグ（８５０）は変数ＩｓＢａｓｉｎにそれぞれ格納される。
【００６２】
なお、上記のような変数が格納される代わりに、図８で示した各変数と同等の情報がコンピュータに格納されてもよい。たとえば、変数に別の名前をつけてもよい。
【００６３】
標高（８１０）は、格子の標高値を格納するための変数である。
【００６４】
格子ごとの標高値に基づいて、格子ごとの勾配を算出することができる。
【００６５】
図９は、本発明の実施形態において算出される勾配の説明図である。
【００６６】
格子（ｉ，ｊ）における勾配は、格子（ｉ，ｊ）の標高値とその格子に隣接する格子における標高値との差を格子間の距離で割った値である。１格子に隣接する格子は８ある。したがって、格子（ｉ，ｊ）における勾配の方向は８ある。すなわち、東、北東、北、北西、西、南西、南、南東である。勾配の方向をそれぞれ方向９０１から方向９０８の矢印で示す。説明のために、格子９００における８方向の勾配を、変数Ｇｒａｄ１からＧｒａｄ８と記載する。
【００６７】
変数Ｇｒａｄ１からＧｒａｄ８の値は、以下式によって求められる。
【００６８】
Grad1=(Glev[i,j]-Glev[i-1,j+1])/hypot(dx,dy) ・・・（４−１）
Grad2=(Glev[i,j]-Glev[i,j+1])/dy ・・・（４−２）
Grad3=(Glev[i,j]-Glev[i+1,j+1])/hypot(dx,dy) ・・・（４−３）
Grad4=(Glev[i,j]-Glev[i-1,j])/dx ・・・（４−４）
Grad5=(Glev[i,j]-Glev[i+1,j])/dx ・・・（４−５）
Grad6=(Glev[i,j]-Glev[i-1,j-1])/hypot(dx,dy) ・・・（４−６）
Grad7=(Glev[i,j]-Glev[i,j-1])/dy ・・・（４−７）
Grad8=(Glev[i,j]-Glev[i+1,j-1])/hypot(dx,dy) ・・・（４−８）
ただし、Ｇｌｅｖ［ｉ，ｊ］は格子（ｉ，ｊ）における勾配、ｈｙｐｏｔ（ｘ，ｙ）は原点（０，０）と点（ｘ，ｙ）との距離を返す関数である。
【００６９】
上記の式によって算出された変数Ｇｒａｄ１からＧｒａｄ８の値は、格子（ｉ，ｊ）から、その格子の周囲に隣接する８個の格子へ向かう方向の下り勾配を示す。例えば、格子（ｉ，ｊ）から格子（ｉ−１，ｊ＋１）に向かう方向の勾配が下りである（すなわち、格子（ｉ，ｊ）の標高より格子（ｉ−１，ｊ＋１）の標高が低い）場合、変数Ｇｒａｄ１の値は正となる。そして、その下り勾配が急である（すなわち、格子（ｉ，ｊ）の標高と格子（ｉ−１，ｊ＋１）の標高との差が大きい）ほど、変数Ｇｒａｄ１の値は大きくなる。
【００７０】
なお、本実施形態では、格子と、その格子に隣接する８格子との間の勾配が計算される。しかし、格子と、その格子の東西方向及び南北方向に隣接する４格子との間の勾配が計算されてもよい。その場合は、上記８方向の勾配変数のうちＧｒａｄ２、Ｇｒａｄ４、Ｇｒａｄ５及びＧｒａｄ７のみが使用される。ただし、勾配が計算される隣接格子の数が増えるほど、精度のよい結果が期待できる。
【００７１】
図８に示した最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）は、これら８方向の変数（Ｇｒａｄ１からＧｒａｄ８）のうち、値が最大である変数に対応する方向を保持する配列である。すなわち、配列ＤｏｗｎＤｉｒ０に格納できる要素は、東、北東、北、北西、西、南西、南及び南東の８方向のうち少なくとも一つを識別するための情報であり、要素の数（配列サイズ）は０から８のいずれかである。この配列の初期値は、要素数０を示す状態であるｅｍｐｔｙである。一つの格子について計算されたＧｒａｄ１からＧｒａｄ８のうち、値が最大であるものが複数存在する場合、複数の方向を識別する値が最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に保持される。
【００７２】
流向候補配列ＤｏｗｎＤｉｒ１（８３０）は、最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に格納された要素のうち、さらに別の基準を満たす要素のみを保持する配列である。この配列の初期値は、要素数０を示す状態であるｅｍｐｔｙである。流向候補配列ＤｏｗｎＤｉｒ１（８３０）に保持される要素については、後に図１５を参照して詳述する。
【００７３】
流域フラグ（８５０）は、その格子が流域に含まれるか否かを判別した結果を格納する変数である。この変数の値は、ｔｒｕｅ又はｆａｌｓｅのいずれかとなる。ｔｒｕｅは格子が流域に含まれることを、ｆａｌｓｅは格子が流域に含まれないことを示す。この変数の初期値が、ｆａｌｓｅである。
【００７４】
それぞれの格子に対応する、図８に示す変数が、メモリ１２０に展開される。格子の数は、極めて多くなる場合がある。その場合、メモリ１２０に展開される変数が非常に多くなる。変数の展開を効率的に行うためには、変数の実装をコンピュータの構造に適した設計にする必要がある。本発明では、格子をコンピュータで効率的に扱うために、部分領域にゴーストセル（ｇｈｏｓｔ ｃｅｌｌ）を追加する。
【００７５】
図１０は、本発明の実施形態のメモリ１２０に展開された状態の部分領域を説明するための図である。
【００７６】
部分領域には、ゴーストセルが付加される。図においてゴーストセルは斜線でハッチングされた格子である。一方、ドットで網掛けされた格子は、部分領域に含まれる格子（すなわち内部セル）である。
【００７７】
ある部分領域に付加されたゴーストセルは、その部分領域に隣接する部分領域に属する内部セルのコピーである。図１０には、（ｂｉ，ｂｊ）＝（０，０）を満たす部分領域１０１０及び（ｂｉ，ｂｊ）＝（１，０）を満たす部分領域１０２０が書かれている。部分領域１０２０に付加されたゴーストセル１０２１は、部分領域１０１０に含まれる内部セル１０１１のコピーである。同様に、部分領域１０１０に付加されたゴーストセル１０１２は、部分領域１０２０に含まれる内部セル１０２２のコピーである。このように、全てのゴーストセルは対応する内部セルのコピーである。
【００７８】
部分領域（ｂｉ，ｂｊ）に含まれる内部セルは、前述したとおり（３）式を満たす格子である。
【００７９】
部分領域（ｂｉ，ｂｊ）に付加されたゴーストセルは、（３）式を満たさず、かつ以下の（５）式を満たす格子である。
【００８０】
bi*si-mi<=I<(bi+1)*si+mi, bj*sj-mj<=j<(bj+1)*sj+mj ・・・（５）
ｍｉ及びｍｊは１以上の整数である。これらの数は、部分領域の内部セルにおける変数を求めるために参照される必要がある格子が、ゴーストセルとしてその部分領域に付加されるように定められる。本実施形態においては、ｍｉ＝１，ｍｊ＝１とした。これは、ある格子の勾配を計算するために、その格子に隣接する８格子の標高値を参照する必要があるためである。ゴーストセルの標高値を参照しなければ、部分領域の端に位置する内部セルの勾配を算出することができない。しかし、ｍｉ＝１，ｍｊ＝１とすれば、一つの部分領域の内部セル及びその部分領域に付加されたゴーストセルの標高値に基づいて、その部分領域の全内部セルの勾配を計算することができる。
【００８１】
このように、部分領域にゴーストセルが付加される場合、部分領域の内部セルを処理するために必要な情報が全てその部分領域及びゴーストセルに含まれる。このことは、複数の部分領域を並列処理することを可能にする。したがって、ゴーストセルが付加された部分領域を利用することで、プログラムの効率が向上することが期待できる。
【００８２】
以上、図７から図１０を使って、流域推定プログラム１４１の処理で使用される格子、部分領域、及び変数を説明した。続いて、流域推定プログラム１４１における流向算出処理（ステップ２０４）の詳細を説明する。
【００８３】
図１１は、本発明の実施形態の流域推定プログラム１４１の流向算出処理（ステップ２０４）を示すフローチャートである。
【００８４】
まず、ステップ１１０１において、ＣＰＵ１３０は、下流端を含む格子を現在の検査対象として指定する。下流端は地図３１０に入力される（図３参照）。ＣＰＵ１３０は、下流端の位置情報を地図３１０における座標系から直交座標系に変換する。続いてＣＰＵ１３０は、（１）式によって直交座標（ｘ，ｙ）を格子（ｉ，ｊ）に変換する。ＣＰＵ１３０は、この格子（ｉ，ｊ）に関する情報を現在の検査対象としてメモリ１２０に記憶する。現在の検査対象を記憶するための配列を、ここではＣｕｒＣｅｌｌｓと記載する。
【００８５】
ステップ１１０２において、ＣＰＵ１３０は、現在の検査対象の格子を含む部分領域を読み込む。このとき、図１０に示すように、ゴーストセルも読み込まれ、部分領域に付加される。このステップの詳細を後述する図１２及び図１３で説明する。
【００８６】
ステップ１１０３において、ＣＰＵ１３０は、処理済み領域をウィンドウ３００に出力する（図４、図５参照）。この処理は、すでに図２で説明したステップ２０６に相当する。
【００８７】
ステップ１１０４において、ＣＰＵ１３０は、部分領域に付加されたゴーストセルの同期を行う。このステップ１１０４の詳細を後述する図１６で説明する。
【００８８】
ステップ１１０５において、ＣＰＵ１３０は、現在の検査対象の上流の格子を次の検査対象として記憶する。具体的には、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌｓに格納された格子に隣接し、かつその格子の上流に位置する格子を次の検査対象としてメモリ１２０に記憶する。次の検査対象を記憶するための配列を、ここではＮｅｘｔＣｅｌｌｓと記載する。このステップ１１０５の詳細を後述する図１５で説明する。
【００８９】
ステップ１１０６において、ＣＰＵ１３０は、次の検査対象を現在の検査対象として新たに指定する。具体的には、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌｓを空にした後、配列ＮｅｘｔＣｅｌｌｓの全ての要素を配列ＣｕｒＣｅｌｌｓに格納する。さらに、ＣＰＵ１３０は、配列ＮｅｘｔＣｅｌｌｓを空にする。
【００９０】
ステップ１１０７において、ＣＰＵ１３０は、現在の検査対象の要素数が０より大か否かを判定する。具体的には、配列ＣｕｒＣｅｌｌｓの要素数と０とが比較される。
【００９１】
ステップ１１０７の結果が真である場合、現在の検査対象の格子の上流の格子が存在する可能性がある。この場合、処理はステップ１１０２に戻る。一方、ステップ１１０７の結果が偽である場合、前回のステップ１１０５において、配列ＮｅｘｔＣｅｌｌｓに追加された格子が存在しない。すなわち、現在までに検出された流域の上流の格子が存在しない。この場合、上流の探索が終了し、処理は図２に示したステップ２０８に進む。
【００９２】
このように、ステップ１１０２からステップ１１０７の処理によって、現在の検査対象の格子に隣接し、かつ上流にある格子が検査対象として指定される。また、ステップ１１０５において現在の検査対象の格子の上流の格子が存在する限りステップ１１０２からステップ１１０７が繰り返される。したがって、図１１に示したステップによって、一つの流域に含まれる下流端から上流までに位置する全ての格子が、必ず一度検査対象として指定される。
【００９３】
図１２は、本発明の実施形態の部分領域読み込み処理（ステップ１１０２）の詳細を示す第１のフローチャートである。
【００９４】
ステップ１２０１において、ＣＰＵ１３０は、現在の検査対象として指定された全ての格子についてステップ１２０２からステップ１２０８までの処理を繰り返す。具体的には、ＣＰＵ１３０は、現在の検査対象として指定された格子のうち一つを選択し、選択された格子についてステップ１２０２から１２０８までの処理を実行する。この選択された格子は、図１２の説明において対象格子と記載される。現在の検査対象として指定された全ての格子についてステップ１２０２から１２０８までの処理が終了するまで、上記の選択及び処理が繰り返される。
【００９５】
ステップ１２０２において、ＣＰＵ１３０は、対象格子を含む部分領域の識別子（部分領域ＩＤ）を算出する。部分領域ＩＤとして、式（２）によって算出される部分領域の座標（ｂｉ，ｂｊ）が使用される。格子（ｉ，ｊ）を含む部分領域座標（ｂｉ，ｂｊ）の算出方法は式（２）に示した。算出された部分領域ＩＤは、図１２の説明において対象部分領域ＩＤと記載される。対象部分領域ＩＤによって識別される部分領域は、図１２の説明において対象部分領域と記載される。
【００９６】
ステップ１２０３において、ＣＰＵ１３０は、対象部分領域ＩＤが配列Ｄｏｍａｉｎに含まれているか否かを判定する。この配列Ｄｏｍａｉｎは、ステップ１２０４以下のステップが実行された部分領域ＩＤを記憶するための配列であり、要素数は０以上である。ステップ１２０３の結果が真であれば、対象部分領域について、既にステップ１２０４以降の処理が実行されている。この場合、処理はステップ１２０９に進む。ステップ１２０３の結果が偽であれば、処理はステップ１２０４以下に進む。
【００９７】
ステップ１２０４において、ＣＰＵ１３０は、対象部分領域ＩＤを配列Ｄｏｍａｉｎに追加する。
【００９８】
ステップ１２０５において、ＣＰＵ１３０は、対象部分領域から所定の距離内にある部分領域をメモリ１２０に展開する。したがって、少なくとも対象部分領域と、対象部分領域に隣接する部分領域がメモリ１２０に展開される。具体的には、それらの部分領域に属する内部セル及びそれらに付加されたゴーストセルのデータ（図８参照）がメモリ１２０に格納される。対象部分領域に隣接する部分領域とは、東西方向３列、南北方向３列の９個の部分領域である。それらの９個の部分領域の中心が対象部分領域である。
【００９９】
なお、後述するように、窪地を埋める処理によって、格子の標高値が変更される場合がある。しかし、対象部分領域から上記の所定の距離より遠くにある窪地を埋めても、対象部分領域の標高は変更されない。これは、窪地を埋める処理において、窪地の面積が所定の値を超えた場合には窪地の標高値が変更されないためである。この処理については後述する（図１４参照）。
【０１００】
上記の所定の距離は、その土地にあった値に設定されることが望ましい。たとえば降雨量が多く土地の浸食の激しい日本では、この所定の距離を１ｋｍ程度にするのが適当と考えられる。逆に、降雨量が少なく土地の浸食が緩やかな米国中西部では、この所定の距離を１０ｋｍ程度にするのが適当と考えられる。部分領域の詳細は上記図１０にて説明した。
【０１０１】
ステップ１２０６において、ＣＰＵ１３０は、標高データ１９２を読み出す。この処理はすでに図２でステップ２０５として説明した。ＣＰＵ１３０は、読み出した標高データのうち、ステップ１２０５において展開された部分領域に属する格子の標高データを、その格子の標高Ｇｌｅｖ（８１０）に代入する。
【０１０２】
ステップ１２０７において、ＣＰＵ１３０は、最大勾配方向を算出する。勾配についてはすでに図９で説明した。各方向の勾配は、式（４−１）から式（４−８）によって算出される。ＣＰＵ１３０は、算出した８方向の勾配の最大値に対応する勾配方向を最大勾配方向配列ＤｏｗｎＤｉｒ０（８２０）に格納する。
【０１０３】
ステップ１２０８において、ＣＰＵ１３０は、最大勾配の値が０以下である格子を初期値として、窪地データを生成する。さらに、ＣＰＵ１３０は、生成した窪地データを配列Ｓｉｎｋｓに格納する。ここで、配列Ｓｉｎｋｓは窪地データを記憶するための配列である。
【０１０４】
図１７は、本発明の実施形態の窪地データの詳細を示す説明図である。
【０１０５】
窪地データは、現在の標高（１７１０）、新しい標高（１７２０）、格子配列（１７３０）、隣接格子配列（１７４０）及び窪地埋め終了フラグ（１７５０）からなる。説明のため、これらをそれぞれＧｌｅｖＣｕｒ、ＧｌｅｖＮｅｗ、Ｃｅｌｌｓ、Ｎｅｉｂｓ及びＩｓＦｉｎｉｓｈと記載する。
【０１０６】
格子配列Ｃｅｌｌｓは、窪地に含まれる格子を識別する情報を格納する配列である。
【０１０７】
隣接格子配列Ｎｅｉｂｓは、窪地に含まれる格子に隣接する格子を識別する情報を格納する配列である。より詳細には、隣接格子配列Ｎｅｉｂｓには、窪地に含まれる各格子に隣接する全ての格子を識別する情報のうち、格子配列Ｃｅｌｌｓに含まれないものが格納される。
【０１０８】
格子配列Ｃｅｌｌｓ及び隣接格子配列Ｎｅｉｂｓの初期値は、要素のないｅｍｐｔｙ状態である。この窪地データの初期化は、１つの格子ｔｈｅＣｅｌｌを与えることによって実行される。ここで与えられるｔｈｅＣｅｌｌは、窪地の底に相当する格子である。
【０１０９】
現在の標高ＧｌｅｖＣｕｒは、窪地内の格子の標高値の候補を格納する変数である。現在の標高ＧｌｅｖＣｕｒの初期値は、その格子ｔｈｅＣｅｌｌの標高である。
【０１１０】
新しい標高ＧｌｅｖＮｅｗは、窪地内の格子の標高値の候補となる可能性がある値を一時的に格納する変数である。
【０１１１】
格子配列Ｃｅｌｌｓには、初期値として、その格子ｔｈｅＣｅｌｌが格納される。
【０１１２】
窪地埋め終了フラグＩｓＦｉｎｉｓｈの初期値は、窪地を埋める処理（図１４参照）が終了していないことを示す偽（ｆａｌｓｅ）である。この値は、窪地を埋める処理が終了すると、真（ｔｒｕｅ）に更新される。
【０１１３】
図１７に示すデータの一組が、配列Ｓｉｎｋｓの一要素に相当する。
【０１１４】
再び図１２に戻り説明する。ステップ１２０９において、ステップ１２０１のループが完了するまで処理はステップ１２０１に戻る。ステップ１２０１のループが完了した場合、処理は図１３に示したステップ１３０１に進む。
【０１１５】
図１２の例では、ステップ１２０５において、対象格子を含む部分領域及びその周囲の一つ以上の部分領域のデータがメモリ１２０に格納される。しかし、図１２は一例であり、メモリ１２０に格納されるデータ量がより多く（又はより少なく）なるように制御されてもよい。
【０１１６】
少なくとも、対象格子に流入する水の流向を算出するためには、対象格子の周囲に隣接する８個の格子の最大勾配方向を算出する必要がある。対象格子が格子（ｉ，ｊ）である場合、対象格子の周囲に隣接する８個の格子とは、ｘ＝ｉ−１、ｉ又はｉ＋１、ｙ＝ｊ−１、ｊ又はｊ＋１を満たす９個の格子（ｘ、ｙ）のうち、対象格子を除く８個の格子である。
【０１１７】
対象格子の周囲に隣接する８個の格子の最大勾配方向を算出するためには、対象格子の標高値、及び、対象格子の周囲に隣接する８個の格子の標高値に加えて、対象格子の周囲に隣接する８個の格子のさらに周囲に隣接する１６個の格子の標高値を取得する必要がある。対象格子が格子（ｉ，ｊ）である場合、上記の１６個の格子は、ｘ＝ｉ−２、ｉ−１、ｉ、ｉ＋１又はｉ＋２、ｙ＝ｊ−２、ｊ−１、ｊ、ｊ＋１又はｊ＋２を満たす２５個の格子（ｘ、ｙ）のうち、対象格子及び対象格子に隣接する８個の格子を除く１６個の格子である。
【０１１８】
このため、ステップ１２０５において、少なくとも、対象格子、対象格子の周囲に隣接する８個の格子、及び、それらの８個の格子の周囲に隣接する１６個の格子からなる２５個の格子のデータがメモリ１２０に格納される必要がある。
【０１１９】
図１３は、本発明の実施形態の部分領域読み込み処理（ステップ１１０２）の詳細を示す第２のフローチャートである。
【０１２０】
ステップ１３０１において、ＣＰＵ１３０は、配列Ｓｉｎｋｓ内の窪地データについてステップ１３０２からステップ１３０５までの処理を繰り返す。
【０１２１】
ステップ１３０２において、ＣＰＵ１３０は、窪地を埋める処理を実行する。この処理の詳細を後述する図１４にて説明する。窪地を埋める処理によって、窪地を埋める処理が終了したことを示す窪地埋め終了フラグＩｓＦｉｎｉｓｈ（１７５０）が変化する。
【０１２２】
ステップ１３０４において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈが真であるか否かを判定する。窪地埋め終了フラグＩｓＦｉｎｉｓｈが真の場合、窪地を埋める処理が終了している。この場合、処理はステップ１３０５に進む。一方、窪地埋め終了フラグＩｓＦｉｎｉｓｈが偽の場合、窪地を埋める処理が終了していない。この場合、処理はステップ１３０６に進む。
【０１２３】
ステップ１３０５において、ＣＰＵ１３０は、窪地を埋める処理が終了した窪地データを配列Ｓｉｎｋから削除する。
【０１２４】
ステップ１３０６において、ＣＰＵ１３０は、ステップ１３０１のループが完了するまでステップ１３０１に戻る。
【０１２５】
ステップ１３０７において、ＣＰＵ１３０は、ステップ１２０７と同じ処理（最大勾配方向の算出）を再び実行する。このステップ１３０７でステップ１２０７と同じ処理を再び行う理由は、ステップ１３０２の窪地埋め処理によって格子の標高値が変更されたためである。すなわち、ステップ１３０７において、変更された標高値を用いて最大勾配方向が算出される。
【０１２６】
図１４は、本発明の実施形態の窪地を埋める処理（ステップ１３０２）の詳細を示すフローチャートである。
【０１２７】
まず、ステップ１４０１において、ＣＰＵ１３０は、配列Ｃｅｌｌｓ内の格子に隣接する格子を配列Ｎｅｉｂｓに格納する。
【０１２８】
ステップ１４０２において、ＣＰＵ１３０は、窪地を埋めるために必要な標高データがメモリに展開されているか否かを判定する。具体的には、ＣＰＵ１３０は、配列Ｎｅｉｂｓ内の全ての格子がメモリ１２０に展開されているか否かを判定する。これは、配列Ｎｅｉｂｓ内の格子（ｉ，ｊ）を含む部分領域（ｂｉ，ｂｊ）を（２）式によって求め、その値が配列Ｄｏｍａｉｎに含まれているか否かを判定することによって実行される。
【０１２９】
ステップ１４０２の判定の結果が真であれば、この処理をさらに進めるためステップ１４０３に進む。判定の結果が偽であれば、窪地埋め処理を実行するために必要な標高データの少なくとも一部がまだメモリ１２０に展開されていない。この場合、この窪地埋め処理を保留するため、処理は図１３のステップ１３０４に進む。
【０１３０】
ステップ１４０３において、ＣＰＵ１３０は、配列Ｎｅｉｂｓ内の格子で最低の標高値を変数ＧｌｅｖＮｅｗに格納する。
【０１３１】
ステップ１４０４において、ＣＰＵ１３０は、窪地を埋める処理が終了したか否かを判定する。これは、条件式ＧｌｅｖＮｅｗ＞ＧｌｅｖＣｕｒによって判定される。この条件式が真であれば、窪地を埋める処理は終了していないためステップ１４０８に進む。偽であれば窪地を埋める処理は終了したためステップ１４０５に進む。
【０１３２】
ステップ１４０５において、ＣＰＵ１３０は、配列Ｃｅｌｌｓ内の格子の標高をＧｌｅｖＣｕｒに書き換える。
【０１３３】
ステップ１４０６において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈをｔｒｕｅにし、この窪地埋め処理を終了するため図１３のステップ１３０４に進む。
【０１３４】
ステップ１４０８において、ＣＰＵ１３０は、現在の標高ＧｌｅｖＣｕｒに新しい標高ＧｌｅｖＮｅｗの値を代入する。
【０１３５】
ステップ１４０９において、ＣＰＵ１３０は、変数ＧｌｅｖＣｕｒの値と等しい標高の格子を配列Ｃｅｌｌｓに追加する。
【０１３６】
ステップ１４１０において、ＣＰＵ１３０は、窪地の大きさ（例えば面積）が所定の基準値よりも小さいか否かを判定する。窪地の面積が所定の基準値よりも大きい場合、その窪地は、標高の誤差に起因する窪地ではなく、実際にその場所に存在する窪地であるとみなされる。この基準値は、その土地にあった値に設定するとよい。たとえば降雨量が多く土地の浸食の激しい日本ではこの基準値を１平方キロメートル程度にするのが適当と考えられる。逆に、降雨量が少なく土地の浸食が緩やかな米国中西部では、この基準値を１０平方キロメートル程度にするのが適当と考えられる。
【０１３７】
このステップ１４１０は、実際には、Ｃｅｌｌｓの要素数が、上記の所定の基準値を格子の面積で割った値より小さいか否かを判定することによって実行されてもよい。あるいは、ステップ１４１０の判定基準として窪地の面積ではなく、窪地の体積が用いられてもよい。ただし、判定基準として面積を使った方が、ＣＰＵ１３０の負荷は軽くなる。
【０１３８】
窪地の体積は、窪地に含まれる各格子の深さを格子の面積に乗算した値を合計することによって算出される。各格子の深さは、各格子の現在の標高ＧｌｅｖＣｕｒからその格子の標高Ｇｌｅｖを減算することによって算出される。
【０１３９】
このステップ１４１０の判別結果が真であれば、窪地は地形データの誤差に起因すると見なされるため、処理はステップ１４０１に戻る。ステップ１４１０の判別結果が偽であれば、窪地は実際の地形に起因すると見なされるため、処理はステップ１４１１に進む。
【０１４０】
ステップ１４１１において、ＣＰＵ１３０は、窪地埋め終了フラグＩｓＦｉｎｉｓｈをｔｒｕｅにし、この窪地埋め処理を終了するため図１３のステップ１３０４に進む。
【０１４１】
上記のように、ステップ１４０４の判定結果が偽（すなわちＮｏ）であった場合、窪地内の格子の標高が変更され、確定する。ステップ１４１０の判定結果が偽であった場合、窪地内の格子の標高が変更されないことが確定する。
【０１４２】
一方、ステップ１４０２の判定結果が偽であった場合、窪地埋め処理が保留される。すなわち、この場合、窪地内の格子の標高がまだ確定しない。その後、図１２のステップ１２０５において、新たな部分領域がメモリ１２０に展開されると、ステップ１４０２の判定結果が真（すなわちＹｅｓ）となり、窪地内の格子の標高を確定させるための処理が続行される。
【０１４３】
しかし、ステップ１３０７が実行される時点で、最大勾配方向が算出される格子及びその格子に隣接する格子の標高が確定している必要がある。言い換えると、ステップ１３０１から１３０６までのループが終了した時点で、少なくとも、現在の検査対象の格子を含む部分領域内の全ての格子の標高が確定している必要がある。この確定のために、ステップ１２０５において、ステップ１４１０の基準値が示す面積より十分に広い範囲の部分領域がメモリ１２０に展開される必要がある。
【０１４４】
具体的には、ステップ１２０５において、対象部分領域内の格子から半径ｒ（ｍ）以内の全ての格子がメモリ１２０に展開される必要がある。ここで、ステップ１４１０の基準値がＡ（平方ｍ）、一格子が１辺ｄｘ（ｍ）の正方形である(すなわちｄｘとｄｙが等しい)場合、半径ｒは、次の式（６）を満たす必要がある。
【０１４５】
ｒ≧Ａ／ｄｘ・・・（６）
上記の条件が満たされる場合、対象部分領域内の格子を含み、かつ、面積が基準値Ａ以内である窪地は、それがどのような形状であっても、ステップ１２０５において展開された部分領域の中に含まれる。したがって、この場合、対象部分領域の少なくとも一部の格子が窪地に含まれていても、その対象部分領域内の全ての格子の標高は、ステップ１３０７が実行される時点で全て確定している。
【０１４６】
図１５は、本発明の実施形態の上流探索処理（ステップ１１０５）の詳細を示すフローチャートである。
【０１４７】
ステップ１５０１において、ＣＰＵ１３０は、配列ＣｕｒＣｅｌｌ内の格子について以降の処理を繰り返す。図１５において、ステップ１５０１から１５０７までのループの処理対象の格子を、Ｃｅｌｌ０と記載する。
【０１４８】
ステップ１５０２において、ＣＰＵ１３０は、Ｃｅｌｌ０に隣接する８個の格子について以降の処理を繰り返す。図１５において、ステップ１５０２から１５０６までのループの処理対象の格子を、Ｃｅｌｌ１と記載する。したがって、Ｃｅｌｌ０とＣｅｌｌ１とは隣接している。
【０１４９】
ステップ１５０３において、ＣＰＵ１３０は、Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０があり、かつ、Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）がＦａｌｓｅであるか否かを判定する。Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０がある場合、そのＣｅｌｌ１から流出した水がＣｅｌｌ０に流入すると推定される。Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）がＦａｌｓｅである場合、そのＣｅｌｌ１がどの流域に属するかまだ確定していない。
【０１５０】
ステップ１５０３の判定の結果が真の場合、Ｃｅｌｌ１とＣｅｌｌ０とが一つの流域に含まれる。この場合、処理はステップ１５０４に進む。
【０１５１】
一方、判定の結果が偽の場合（すなわち、Ｃｅｌｌ１の最大勾配方向にＣｅｌｌ０がないか、又は、Ｃｅｌｌ０の流域フラグ（ＩｓＢａｓｉｎ）がＴｒｕｅである場合）、処理はステップ１５０６に進む。
【０１５２】
ステップ１５０４において、ＣＰＵ１３０は、Ｃｅｌｌ１の流向候補配列（ＤｏｗｎＤｉｒ１）にＣｅｌｌ０の方向を追加する。たとえば、Ｃｅｌｌ０がＣｅｌｌ１の北西方向にあれば、その方向が流向の候補として配列に追加される。
【０１５３】
ステップ１５０５において、ＣＰＵ１３０は、Ｃｅｌｌ１の流域フラグ（ＩｓＢａｓｉｎ）を真にする。
【０１５４】
ステップ１５０６において、ＣＰＵ１３０は、ステップ１５０２のループが完了するまでステップ１５０２に戻る。
【０１５５】
ステップ１５０２からステップ１５０６までのループが完了した時点で、Ｃｅｌｌ１の流向候補配列（ＤｏｗｎＤｉｒ１）に含まれる流向候補が確定する。例えば、ステップ１５０４において、Ｃｅｌｌ０がＣｅｌｌ１の北西方向にある場合、その方向（北西方向）がＣｅｌｌ１の流向の候補として追加される。そして、ループが完了すると、この流向が確定する。流向の確定は、以下のようにすればよい。流向候補配列ＤｏｗｎＤｉｒ１に１つの方向が格納されている場合、その方向を流向として確定する。流向候補配列ＤｏｗｎＤｉｒ１に２以上の方向が格納されている場合、あらかじめ決めた方法にしたがって1つの方向を選択する。選択する方法は、乱数により決めてもいいし、あらかじめ決めた優先度（たとえば、東、北東、北、北西、西、南西、南、南東の順に優先して選択する）に従って決めてもよい。Ｃｅｌｌ１の水がＣｅｌｌ０に流出する場合、Ｃｅｌｌ０とＣｅｌｌ１とが一つの流域に含まれると判定される。このＣｅｌｌ１は、Ｃｅｌｌ０の次の検査対象として選択される。すなわち、このＣｅｌｌ１は、配列ＮｅｘｔＣｅｌｌｓに追加され、図１１のステップ１１０６において、現在の検査対象として新たに選択される。
【０１５６】
ステップ１５０７において、ＣＰＵ１３０は、ステップ１５０１のループが完了するまでステップ１５０１に戻る。
【０１５７】
図１６は、本発明の実施形態の部分領域に付加されたゴーストセルの同期処理（ステップ１１０４）の詳細を示すフローチャートである。
【０１５８】
ステップ１６０１において、ＣＰＵ１３０は、配列Ｄｏｍａｉｎ内の部分領域について以降の処理を繰り返す。
【０１５９】
ステップ１６０２において、ＣＰＵ１３０は、部分領域に付加されたすべてのゴーストセルについて以降の処理を繰り返す。
【０１６０】
ステップ１６０３において、ＣＰＵ１３０は、ゴーストセルに対応する内部セルの変数値をコピーする。例えば、図１０のゴーストセル１０１２についてステップ１６０３が実行される場合、ゴーストセル１０１２に対応する内部セル１０２２の変数値（すなわち、図８に示す各変数の値）がゴーストセル１０１２にコピーされる。その結果、内部セル１０２２の標高値が窪地埋め処理によって変更された場合、その変更がゴーストセル１０１２に反映される。
【０１６１】
ステップ１６０４において、ＣＰＵ１３０は、ステップ１６０２のループが完了するまでステップ１６０２に戻る。
【０１６２】
ステップ１６０５において、ＣＰＵ１３０は、ステップ１６０１のループが完了するまでステップ１６０１に戻る。
【０１６３】
以上述べたとおり、本発明の実施形態によれば、流域に関する事前知識を要求としない、流域推定コンピュータプログラムを提供することができる。このプログラムを使用することによって、コンピュータが使えるユーザであれば、専門家と同等の品質の流域推定が可能になる。
【図面の簡単な説明】
【０１６４】
【図１】本発明の実施形態におけるコンピュータシステムの構成を示す概略ブロック図である。
【図２】本発明の実施形態のＰＣによって実行される処理を示すシーケンス図である。
【図３】本発明の実施形態の流域推定プログラムが提供するＧＵＩの入力待ちの状態を示す説明図である。
【図４】本発明の実施形態の流域推定プログラムが提供するＧＵＩの解析中の状態を示す第１の説明図である。
【図５】本発明の実施形態の流域推定プログラムが提供するＧＵＩの解析中の状態を示す第２の説明図である。
【図６】本発明の実施形態の流域推定プログラムが提供するＧＵＩの流域推定済みの状態を示す説明図である。
【図７】本発明の実施形態の流域推定プログラムが利用する格子及び部分領域を示す説明図である。
【図８】本発明の実施形態の格子に割り当てられる変数の説明図である。
【図９】本発明の実施形態において算出される勾配の説明図である。
【図１０】本発明の実施形態のメモリに展開された状態の部分領域を説明するための図である。
【図１１】本発明の実施形態の流域推定プログラムの流向算出処理を示すフローチャートである。
【図１２】本発明の実施形態の部分領域読み込み処理の詳細を示す第１のフローチャートである。
【図１３】本発明の実施形態の部分領域読み込み処理の詳細を示す第２のフローチャートである。
【図１４】本発明の実施形態の窪地を埋める処理の詳細を示すフローチャートである。
【図１５】本発明の実施形態の上流探索処理の詳細を示すフローチャートである。
【図１６】本発明の実施形態の部分領域に付加されたゴーストセルの同期処理の詳細を示すフローチャートである。
【図１７】本発明の実施形態の窪地データの詳細を示す説明図である。
【図１８】ＤＥＭの誤差の影響を説明するための図である。
【符号の説明】
【０１６５】
１００ ＰＣ
３１０ 地図
３２１ 下流端を示すアイコン
４１３、５１３ 処理済み領域
６１４ 流向ベクトル
７０３ 格子線
７０６ 部分領域線
１０１０、１０２０ 部分領域
【特許請求の範囲】
【請求項１】
計算機を制御するプログラムであって、
前記計算機は、前記プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備え、
前記プログラムは、
地図上の地点を指示する入力を受け付ける第１手順と、
前記地図を互いに直交して設定された第１軸及び第２軸に平行にそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、
前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、
前記第１領域内の窪地を検出する第４手順と、
前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、
前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、
前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を前記プロセッサに実行させることを特徴とするプログラム。
【請求項２】
前記第１領域は、少なくとも、前記第１格子の周囲に隣接する８個の前記第２格子、及び、前記８個の第２格子の周囲に隣接する１６個の格子を含むことを特徴とする請求項１に記載のプログラム。
【請求項３】
前記第４手順は、前記第１領域内の格子の標高値が、その格子に隣接する８個の格子のうち最も標高が低い格子の標高値以下である場合、その格子を含む窪地を検出する手順を含み、
前記第５手順は、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子が、前記第１領域に含まれ、かつ、前記検出された窪地の大きさが所定の閾値を超えない場合に、前記検出された窪地に含まれる格子の標高値を変更すると判定する手順と、
前記検出された窪地に含まれる格子の標高値を変更すると判定された場合、前記検出された窪地に含まれる全ての格子の標高値を、前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値に変更する手順と、を含むことを特徴とする請求項１に記載のプログラム。
【請求項４】
前記第５手順は、さらに、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値を、前記窪地に含まれる全ての格子の標高値の候補として前記メモリに保持する手順と、
前記検出された窪地の大きさが前記所定の閾値を超えた場合、前記標高値の候補を前記メモリから削除する手順と、を含み、
前記格子の標高値の候補が保持されている場合、前記第４手順において、前記保持されている標高値の候補が前記格子の標高値として使用されることを特徴とする請求項３に記載のプログラム。
【請求項５】
前記検出された窪地の大きさは、前記検出された窪地の面積又は前記検出された窪地の体積であることを特徴とする請求項３に記載のプログラム。
【請求項６】
前記第６手順は、
前記第２格子の標高値から、前記第２格子に隣接する８個の各格子の標高値を減算した値を、前記第２格子と前記第２格子に隣接する各格子との間の距離で除算することによって、前記第２格子から前記第２格子に隣接する各格子への下り勾配の値を算出する手順と、
前記第２格子から前記第１格子への下り勾配の値が、前記第２格子に隣接する８個の各格子への下り勾配の値の最大値である場合、前記第２格子から流出する水が前記第１格子へ流れると判定する手順と、を含むことを特徴とする請求項１に記載のプログラム。
【請求項７】
前記第９手順は、前記第１格子として新たに選択された格子から流出する水の流向を示す情報を出力する手順を含むことを特徴とする請求項１に記載のプログラム。
【請求項８】
前記プログラムは、前記第８手順において処理対象である第１格子が新たに選択された場合、さらに前記第３手順から前記第８手順を前記プロセッサに実行させることを特徴とする請求項１に記載のプログラム。
【請求項９】
プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備える計算機であって、
前記プロセッサは、
地図上の地点を指示する入力を受け付ける第１手順と、
前記地図を互いに直交して設定された第１軸及び第２軸に平行にそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、
前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、
前記第１領域内の窪地を検出する第４手順と、
前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、
前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、
前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を実行することを特徴とする計算機。
【請求項１０】
前記第１領域は、少なくとも、前記第１格子の周囲に隣接する８個の前記第２格子、及び、前記８個の第２格子の周囲に隣接する１６個の格子を含むことを特徴とする請求項９に記載の計算機。
【請求項１１】
前記第４手順は、前記第１領域内の格子の標高値が、その格子に隣接する８個の格子のうち最も標高が低い格子の標高値以下である場合、その格子を含む窪地を検出する手順を含み、
前記第５手順は、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子が、前記第１領域に含まれ、かつ、前記検出された窪地の大きさが所定の閾値を超えない場合に、前記検出された窪地に含まれる格子の標高値を変更すると判定する手順と、
前記検出された窪地に含まれる格子の標高値を変更すると判定された場合、前記検出された窪地に含まれる全ての格子の標高値を、前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値に変更する手順と、を含むことを特徴とする請求項９に記載の計算機。
【請求項１２】
前記第５手順は、さらに、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全て格子の標高値のうち最も低い値を、前記窪地に含まれる全ての格子の標高値の候補として前記メモリに保持する手順と、
前記検出された窪地の大きさが前記所定の閾値を超えた場合、前記標高値の候補を前記メモリから削除する手順と、を含み、
前記格子の標高値の候補が保持されている場合、前記第４手順において、前記保持されている標高値の候補が前記格子の標高値として使用されることを特徴とする請求項１１に記載の計算機。
【請求項１３】
前記検出された窪地の大きさは、前記検出された窪地の面積又は前記検出された窪地の体積であることを特徴とする請求項１１に記載の計算機。
【請求項１４】
前記第６手順は、
前記第２格子の標高値から、前記第２格子に隣接する８個の各格子の標高値を減算した値を、前記第２格子と前記第２格子に隣接する各格子との間の距離で除算することによって、前記第２格子から前記第２格子に隣接する各格子への下り勾配の値を算出する手順と、
前記第２格子から前記第１格子への下り勾配の値が、前記第２格子に隣接する８個の各格子への下り勾配の値の最大値である場合、前記第２格子から流出する水が前記第１格子へ流れると判定する手順と、を含むことを特徴とする請求項９に記載の計算機。
【請求項１５】
前記第９手順は、前記第１格子として新たに選択された格子から流出する水の流向を示す情報を出力する手順を含むことを特徴とする請求項９に記載の計算機。
【請求項１６】
前記プロセッサは、前記第８手順において処理対象である第１格子が新たに選択された場合、さらに前記第３手順から前記第８手順を実行することを特徴とする請求項９に記載の計算機。
【請求項１】
計算機を制御するプログラムであって、
前記計算機は、前記プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備え、
前記プログラムは、
地図上の地点を指示する入力を受け付ける第１手順と、
前記地図を互いに直交して設定された第１軸及び第２軸に平行にそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、
前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、
前記第１領域内の窪地を検出する第４手順と、
前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、
前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、
前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を前記プロセッサに実行させることを特徴とするプログラム。
【請求項２】
前記第１領域は、少なくとも、前記第１格子の周囲に隣接する８個の前記第２格子、及び、前記８個の第２格子の周囲に隣接する１６個の格子を含むことを特徴とする請求項１に記載のプログラム。
【請求項３】
前記第４手順は、前記第１領域内の格子の標高値が、その格子に隣接する８個の格子のうち最も標高が低い格子の標高値以下である場合、その格子を含む窪地を検出する手順を含み、
前記第５手順は、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子が、前記第１領域に含まれ、かつ、前記検出された窪地の大きさが所定の閾値を超えない場合に、前記検出された窪地に含まれる格子の標高値を変更すると判定する手順と、
前記検出された窪地に含まれる格子の標高値を変更すると判定された場合、前記検出された窪地に含まれる全ての格子の標高値を、前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値に変更する手順と、を含むことを特徴とする請求項１に記載のプログラム。
【請求項４】
前記第５手順は、さらに、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値を、前記窪地に含まれる全ての格子の標高値の候補として前記メモリに保持する手順と、
前記検出された窪地の大きさが前記所定の閾値を超えた場合、前記標高値の候補を前記メモリから削除する手順と、を含み、
前記格子の標高値の候補が保持されている場合、前記第４手順において、前記保持されている標高値の候補が前記格子の標高値として使用されることを特徴とする請求項３に記載のプログラム。
【請求項５】
前記検出された窪地の大きさは、前記検出された窪地の面積又は前記検出された窪地の体積であることを特徴とする請求項３に記載のプログラム。
【請求項６】
前記第６手順は、
前記第２格子の標高値から、前記第２格子に隣接する８個の各格子の標高値を減算した値を、前記第２格子と前記第２格子に隣接する各格子との間の距離で除算することによって、前記第２格子から前記第２格子に隣接する各格子への下り勾配の値を算出する手順と、
前記第２格子から前記第１格子への下り勾配の値が、前記第２格子に隣接する８個の各格子への下り勾配の値の最大値である場合、前記第２格子から流出する水が前記第１格子へ流れると判定する手順と、を含むことを特徴とする請求項１に記載のプログラム。
【請求項７】
前記第９手順は、前記第１格子として新たに選択された格子から流出する水の流向を示す情報を出力する手順を含むことを特徴とする請求項１に記載のプログラム。
【請求項８】
前記プログラムは、前記第８手順において処理対象である第１格子が新たに選択された場合、さらに前記第３手順から前記第８手順を前記プロセッサに実行させることを特徴とする請求項１に記載のプログラム。
【請求項９】
プログラムが格納されるメモリと、前記メモリに格納された前記プログラムを実行するプロセッサと、を備える計算機であって、
前記プロセッサは、
地図上の地点を指示する入力を受け付ける第１手順と、
前記地図を互いに直交して設定された第１軸及び第２軸に平行にそれぞれ等間隔に区切る線によって囲まれた格子のうち、前記入力によって指定された地点を含む格子を、処理対象である第１格子として選択する第２手順と、
前記第１格子を含む第１領域内の格子に対して定められた標高値を前記メモリに格納する第３手順と、
前記第１領域内の窪地を検出する第４手順と、
前記検出された窪地に含まれる格子の標高値を、前記窪地に含まれる格子に隣接する一つの格子の標高値に変更する第５手順と、
前記第１領域内の格子の標高値に基づいて、前記第１領域内の格子のうち、前記第１格子に隣接する第２格子から流出する水の流向を判定する第６手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子と前記第１格子とが一つの流域に含まれると判定する第７手順と、
前記第２格子から流出する水が前記第１格子に流れると判定された場合、前記第２格子を、処理対象である第１格子として新たに選択する第８手順と、
前記第１格子として新たに選択された格子が前記一つの流域に含まれることを示す情報を出力する第９手順と、を実行することを特徴とする計算機。
【請求項１０】
前記第１領域は、少なくとも、前記第１格子の周囲に隣接する８個の前記第２格子、及び、前記８個の第２格子の周囲に隣接する１６個の格子を含むことを特徴とする請求項９に記載の計算機。
【請求項１１】
前記第４手順は、前記第１領域内の格子の標高値が、その格子に隣接する８個の格子のうち最も標高が低い格子の標高値以下である場合、その格子を含む窪地を検出する手順を含み、
前記第５手順は、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子が、前記第１領域に含まれ、かつ、前記検出された窪地の大きさが所定の閾値を超えない場合に、前記検出された窪地に含まれる格子の標高値を変更すると判定する手順と、
前記検出された窪地に含まれる格子の標高値を変更すると判定された場合、前記検出された窪地に含まれる全ての格子の標高値を、前記検出された窪地に含まれる全ての格子の周囲に隣接する全ての格子の標高値のうち最も低い値に変更する手順と、を含むことを特徴とする請求項９に記載の計算機。
【請求項１２】
前記第５手順は、さらに、
前記検出された窪地に含まれる全ての格子の周囲に隣接する全て格子の標高値のうち最も低い値を、前記窪地に含まれる全ての格子の標高値の候補として前記メモリに保持する手順と、
前記検出された窪地の大きさが前記所定の閾値を超えた場合、前記標高値の候補を前記メモリから削除する手順と、を含み、
前記格子の標高値の候補が保持されている場合、前記第４手順において、前記保持されている標高値の候補が前記格子の標高値として使用されることを特徴とする請求項１１に記載の計算機。
【請求項１３】
前記検出された窪地の大きさは、前記検出された窪地の面積又は前記検出された窪地の体積であることを特徴とする請求項１１に記載の計算機。
【請求項１４】
前記第６手順は、
前記第２格子の標高値から、前記第２格子に隣接する８個の各格子の標高値を減算した値を、前記第２格子と前記第２格子に隣接する各格子との間の距離で除算することによって、前記第２格子から前記第２格子に隣接する各格子への下り勾配の値を算出する手順と、
前記第２格子から前記第１格子への下り勾配の値が、前記第２格子に隣接する８個の各格子への下り勾配の値の最大値である場合、前記第２格子から流出する水が前記第１格子へ流れると判定する手順と、を含むことを特徴とする請求項９に記載の計算機。
【請求項１５】
前記第９手順は、前記第１格子として新たに選択された格子から流出する水の流向を示す情報を出力する手順を含むことを特徴とする請求項９に記載の計算機。
【請求項１６】
前記プロセッサは、前記第８手順において処理対象である第１格子が新たに選択された場合、さらに前記第３手順から前記第８手順を実行することを特徴とする請求項９に記載の計算機。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図１３】
【図１４】
【図１５】
【図１６】
【図１７】
【図１８】
【公開番号】特開２００９−１６３６２１（Ｐ２００９−１６３６２１Ａ）
【公開日】平成２１年７月２３日（２００９．７．２３）
【国際特許分類】
【出願番号】特願２００８−２２８４（Ｐ２００８−２２８４）
【出願日】平成２０年１月９日（２００８．１．９）
【出願人】（０００２３３０４４）株式会社日立エンジニアリング・アンド・サービス (276)
【公開日】平成２１年７月２３日（２００９．７．２３）
【国際特許分類】
【出願日】平成２０年１月９日（２００８．１．９）
【出願人】（０００２３３０４４）株式会社日立エンジニアリング・アンド・サービス (276)
[ Back to top ]