SURFACE

surface − 調整可能な張力による連続曲面グリッド化アルゴリズム

surface [ xyzfile ] −Goutputfile.grd −Ixinc[unit][=|+][/yinc[unit][=|+]] −Rwest/east/south/north[r] [ −Aaspect_ratio ] [ −Cconvergence_limit ] [ −H[i][nrec] ] [ −Lllower ] [ −Luupper ] [ −Nmax_iterations ] [ −Q ] [ −Ssearch_radius[m] ] [ −Ttension_factor[i|b] ] [ −V ] [ −Zover-relaxation_factor ] [ −:[i|o] ] [ −bi[s|S|d|D][ncol] ] [ −fcolinfo ]

surface はランダムに分布する (x,y,z) の 3 つの値を標 準 入 力[ ま た は xyzfile ]から読み、次の式を解くことにより格子点値 z(x,y) のバイナリグリ ッドファイルを作成する。

(1 - T) * L (L (z)) + T * L (z) = 0

ここで T は 0 から 1 までの値をとる張力ファクタであり、 L はラプラス 作 用素を示す。 T = 0 は SuperMISP と ISM パッケージと同等の「最小曲率」解 を与える。最小曲率は望まない振動や偽の極大・極小の原因となることがあ る の で( Smith and Wessel, 1990 を参照)、これらの影響を防止したければ T > 0 を用いると良い。実験的には、 T ~ 0.25 が通常はポテンシャル場のデー タ に は 良 いようであり、急峻な地形データに対しては T をもっと大きく(T ~ 0.35)した方が良いようである。 T = 1 は調和曲面(支配データ点以外の場所で 極大・極小ができない)を与える。空間的なエイリアシングの防止や余分なデー タの消去のために、データを blockmean 、 blockmedian 又は blockmode によ り 前処理することを推奨する。解の上限や下限を強制したい場合は、これらを 固定値、値を持ったグリッドファイル、又は単に入力データの値の最大/最小値 といった形で入れることができる。

xyzfile

(x,y,z) のデータ値を持つ 3 列のアスキーファイル [又はバイナリ、 −b 参照]。ファイルを指定しないときは、 surface は標準入力から読 み込む。

−G

出力されるグリッドファイルの名前。出力はバイナリの 2 次元の .grd ファイル。

−I

x_inc [及び随意で y_inc ] は格子間隔である。分を示すには m を、秒を示すには c を付けること。もし代わりに e, k, i, n のいず れかの単位が付け加えられれば、それぞれ刻み幅はメートル、キロ メ ー トル、マイル、海里のいずれかと見なされ、その領域の中央の緯度 における経度の等価な度単位に変換される（変換方法は ELLIPSOID に 依 存 す る ）。もし /y_inc が与えられて 0 に設定されていれば、 /x_inc と同じ値に再設定され、そうでない場合は緯度の度に変換され る 。 もし = が付加されていたら、対応する最大の x (東) 又は y ( 北) が与えられた刻み幅にちょうど一致するように微調整される [ デ フォルトでは刻み幅は与えられた領域に合うように微調整される]。最 後に、刻み幅を与える代わりに、与えられた刻み幅に + を付加するこ と により、必要な格子点の数を指定することができる。刻み幅は格子 点の数と領域から再計算される。再計算される刻み幅は格子保存グ リ ッ ドか画素保存グリッドのどちらを選ぶかに依存する。詳細は付録 B を参照すること。

−R

xmin, xmax, ymin, ymax は対象とする領域を指定する。地理的領 域 に対しては、これらの境界は西、東、南、北に対応し、小数値を含 む度又は [+-]dd:mm[:ss.xxx][W|E|S|N] のフォーマットで指定できる 。東西南北の代わりに地図座標の左下と右上を使うときは r を付ける こと。 2 つの略号 −Rg −Rd は全球領域(それぞれ経度 0/360 又 は -180/+180 と緯度 -90/+90 )を表す。暦時間座標系に対しては相対時 間(選択された TIME_EPOCH に対して相対的で、選択された TIME_UNIT を単位とする。 t を −JX|x に付けること)又は [date]T[clock] 形式 の絶対時刻(T を −JX|x に付けること)のどちらかを与えることができ る。少なくとも date と clock のどちらかが存在していなければなら ない。 T は常に必要である。 date の文字列は [-]yyyy[-mm[-dd]] ( グレゴリオ暦)又は yyyy[-Www[-d]] (ISO 週暦)の形式でなければなら ない。 clock の文字列は hh:mm:ss[.xxx] の形式でなければならない 。 区切り文字を利用する場合にはその型と位置が指定どおりでなけれ ばならない(ただし入出力と描画フォーマットは変更可能である)。

−A

アスペクト比。もし望むならば、グリッドの非等方性を方程式 に 加 え る こ と が で きる。 aspect_ratio を入れると、 dy = dx / aspect_ratio によりグリッドの大きさと関係付けられる [デフォルト = 1 では等方的グリッドを仮定]。

−C

収束限界。反復計算では、全てのグリッドの値の変化の絶対値の 最大値が convergence_limit 未満になったとき収束したと仮定する ( 単 位は z と同じ) [デフォルトでは入力データの典型的な勾配の 0.1 パーセントにスケールされる]。

−H

入力ファイルにヘッダ行がある時に使う。 ヘ ッ ダ 行 の 数 は .gmtdefaults4 ファイルを編集することにより変更できる。このオプ ションが使われた場合、 GMT のデフォルトではヘッダ行は 1 行で あ る。入力データだけにヘッダ行が必要な場合は −Hi を使うこと [デフ ォルトでは入力データにヘッダ行があればそれも出力される]。バイナ リデータには使わないこと。

−L

出 力 される解に制限を加える。 llower で下限を設定する。 lower には下限値を持つグリッドファイルの名前、固定値、 d を付け て入力の最小値を設定、又は u を付けて非制限[デフォルト]とする。 uupper で上限を設定し、上限値を持つグリッドファイルの名前、固定 値、 d を付けて入力の最小値を設定、又は u を付けて非制限[デフォ ルト]とする。

−N

反復計算の回数。反復は convergence_limit に達したとき、又は 反 復回数が max_iterations に達したときに停止する [デフォルトは 250 ]。

−Q

多くの数から合成された最大公約数となるグリッドの大きさを 提 案 する。これにより、解にいくつかの中間段階を使うことによって、 実行時間が短くなり、より良い結果を得ることができる。 −Q によ っ て 提案されたサイズは、 −R や −I を変更することによって実現でき る。 surface の出力に対して grdsample 又は grdcut を用いるこ と により必要な −R や −I に戻すことができる。

−S

探索半径。 x,y データと同じ単位で search_radius を設定する こと。分を示すには m を付けること。これは反復の最初の計算の前に グ リッドを初期化するために使われる。グリッド格子が素数で領域分 割が不可能な場合でもない限り、これに時間を費やす価値はない [ デ フォルト = 0.0 で探索をしない]。

−T

張力ファクタ。これらは 0 と 1 の間でなければならない。張力 は内部の解(上の方程式で不要な振幅を防止する)と境界条件(解を端に 近づくにつれて平らにする)で用いられる。両方にゼロを使うと自由端 を持つ最小曲率解になる。すなわち自然 3 次スプラインである。内部 の張力を設定するには −Ttension_factori を、境界の張力を設定する には −Ttension_factorb を使うこと。もし i 又は b を付けなければ 、両方とも同じ値に設定される。 [デフォルトは両方とも = 0 で最小 曲率解を与える]。

−V

冗長モードを選択する。標準エラー出力に経過報告を送る [デ フ ォ ルトでは「黙って」実行する]。 −Vl は反復毎に収束状況を報告す る。 −V は各領域のグリッドが収束した後にだけ報告する。

−Z

過大緩和因子。このパラメータは収束を速めるために使い、 1 と 2 の間の値をとる。 1 のときは方程式をその通りに反復計算し、常に 安定収束を保証する。それより大きな値では収束する間の増分を過 大 評 価するので、より早く解に到達するが不安定になることがある。こ のファクタに大きな値を用いるときは、 −Vl オプションで毎回の反復 計算を監視すると良い [デフォルト = 1.4 で早く収束しほとんどの場 合に安定である]。

−:

入出力の(経度、緯度)と(緯度、経度)を入れ替える [デフォル ト は(経度、緯度)]。 i を付けると入力のみが、 o を付けると出力のみ が入れ替わる [デフォルトでは両方入れ替わる]。

−bi

バイナリ入力を選択する。単精度に対しては s を付けること [ デ フォルトは d (倍精度)]。大文字の S (又は D) にするとバイトス ワップを行う。随意で、バイナリファイルの列数がこのプログラム で 必 要な列数を越えるときは、バイナリファイルの列数 ncol を付ける こと [デフォルトは 3 ]。

−f

入出力列の特殊フォーマットを使う(時間又は地 理 デ ー タ) 。 i(nput) 又は o(utput) を指定すること [デフォルトでは入力と出力 の両方]。コンマ区切りで 1 列以上(又は列の範囲)を与えること。 各 列 又 は 列 の範囲の項目に対して T (絶対暦時刻)、 t (選択された TIME_EPOCH に対する相対時間)、 x (経度)、 y (緯度)、 f ( 小 数) の いずれかを付けること。 −f[i|o]g は −f[i|o]0x,1y (地理座標)の 省略形である。

入力データの精度に関わらず、グリッドファイルを GMT のプログラムは 4 バ イ ト浮動小数点配列として内部に値を持っている。これはメモリの節約のため であり、更には全ての実数は 4 バイト浮動小数点値を用いて格納することが可 能だからである。より高精度のデータ(すなわち倍精度値)は一度 GMT がグリッ ド処理をして新しいグリッドを出力すると精度を失うだろう。データ処理の 時 に 精度の損失を制限するにはデータを処理する前に規格化することを常に考慮 すべきである。

hawaii_5x5.xyg のアスキーデータの重力をブロック平均して経度 5 分、緯 度 5 分にグリッド化するため、 tension_factor = 0.25 、 convergence_limit = 0.1 ミリガルを用い、途中計算を監視し、結果を hawaii_grd.grd というフ ァ イルに出力する。

surface hawaii_5x5.xyg −R198/208/18/25 −I5m −Ghawaii_grd.grd −T0.25 −C0.1 −VL

surface は各格子点について複数のデータ が あ る と 、 先 に blockmean, blockmedian, 又は blockmode を実行するように提案する。 blockm* を実行し てもまだこのメッセージが出力されるときは、格子間隔が blockm* が用いる出 力 フォーマットにおいて更に多くの桁を必要とするほど小さいことを意味する ことが多い。 blockm* を実行する前に .gmtdefaults4 ファイル の D_FORMAT パ ラメータを変更して桁数を増やすなり、単精度または倍精度の記憶領域を用 いたバイナリ入/出力を選択するなりした方が良い。

blockmean(l), blockmedian(l), blockmode(l), GMT(l), nearneighbor(l), triangulate(l)

Smith, W. H. F, and P. Wessel, 1990, Gridding with continuous curvature splines in tension, Geophysics, 55, 293−305.