| GSJ Open-File Report, no. 519 | December, 2009 |
Tadashi Nakatsuka (Institute of Geology and Geoinformation, GSJ, AIST)
著者は，永年にわたり空中磁気探査の研究に従事し， その中でデータ処理・解析の各種プログラムを開発してきた． それらの多くは，すでに概要を学会誌や学術講演会で発表済み [中塚・大熊，2005a, 2005b; Nakatsuka and Okuma, 2006a, 2006b] であるが，本研究資料集においてそれらのソースプログラムを公開している．
この研究資料では，次に示すように，空中磁気探査に関する
3つのグループのソフトウェアについて取り扱う．
(1) DPAM:
空中磁気探査の計測データのデータ処理プログラム，
(2) GDMP:
格子点データ操作関係プログラム，
(3) ANAM:
磁気異常データの解析処理プログラム．
このソフトウェア体系は，やや単純な機能を有する個別のプログラムを多数集めた
ライブラリの形式になっている．
従って，何らかの課題・目的に対する実際の処理プロセスは，
一連のプログラム実行処理が必要になるであろう．
このシステムで利用しているデータFormatは，我々独自のスタイルに合わせてある．
しかし，そのデザインは，ほぼ全てのデータをASCIIテキストで扱うようにしており，
ユーザがデータの必要なFormat変換を行うことは容易であろう．
この研究資料は，同名の研究資料集 no.449 [Nakatsuka, 2007] の改訂版である．主要な改訂点は，次のとおりである．
DPAM プログラム群は，以下の空中磁気探査の計測データのデータ処理プログラム
からなる．
alog2asc, xldam, xldpn, xldhg, despike,
dvcorr, ecomp, fcomp, ggrid, ggrids,
pframe, pltrk, pchkdv, pchkmag, pchkres,
pchkcomp, xslin, xslina.
プログラム alog2asc は，空中磁気探査のバイナリデータ（我々のシステム 固有）を ASCIIテキスト形式に変換する． プログラム xldam／xldpn は ASCII形式生データから 「DPAM 測線データ」を生成するが，xldam が リアルタイムGPSデータの位置データを用いるのに対し，xldpn は "PostNav" 処理したディファレンシャルGPSの位置データを用いる． また プログラム xldhg は，ヘリコプター探査の空中磁気データから類似の "HGAM 測線データ" を生成するものである． xldam／xldpn／xldhg に共通して，中塚 [2005] に基づく サブプログラム "igrf" を用いたIGRF残差の計算が行われる．
プログラム despike は，「DPAM測線データ」の磁力値データ （全磁力測定値 および IGRF残差）のスパイクノイズを除去し，プログラム dvcorr は 定点観測磁力値データを用いて日変化補正を行う． プログラム ecomp／fcomp は 「DPAM測線データ」に対して 中塚・大熊 [2005b] に基づく機体磁気補償の処理を実行するが， ecomp が 機体磁気補償の動作確認のためにテスト飛行のデータ自体に 補正処理を実行するのに対し，fcomp は そのデータを用いて 調査測線データに対する補正を実行する．
プログラム ggrid／ggrids は，「DPAM測線データ」などの ランダム点データから，Smith and Wessel [1990] の "continuous curvature splines in tension" の方法により，グリッドデータを生成する． ggrid は 磁力値（IGRF残差）のみについてグリッド化を行い， ggrids は 高度データも含めてグリッド化する．
プログラム pframe は，各種位置パラメータ指定の際の参照用に， 調査地域の座標値フレームワークを作図し，プログラム pltrk は 「DPAM測線データ」から測線位置図を作図する． プログラム pchkmag／pchkres／pchkcomp は， データ品質の点検のために磁力値データ図化表示を行うが， pchkmag は 全磁力測定値を黒線で，pchkres は IGRF残差データを 青線で図化し，pchkcomp は 機体磁気補償処理前後両方のIGRF残差データを 赤線/青線で色別して図化する．
プログラム xslin／xslina は，ランダム点データから 必要データだけ抽出し，標準測線(StdLIN)データに変換する． xslin では 種々の形式のランダム点データがサポートされるのに対し， xslina では 「DPAM測線データ」に限定して平滑化再サンプリング の機能を有する．
GDMP プログラム群は，格子点データを操作するプログラムで構成されており，
以下のプログラムがある．
sel, seldb, altx, adjlv, gadd, gsub,
gtrim, govlay, gojoin, gmerge, txproj,
gtopo, gtrf, plmap, plmapc, plmapg,
plmaps, plmapcs, shade, plmapl, plmapcl,
xplmap, xplmapc.
プログラム sel／seldb は，既存グリッドデータから 再グリッディングにより新たなグリッドデータを生成する． sel では 入力ファイル名を指定してグリッドデータを選択するのに対し， seldb では 日本空中磁気探査データベース(AMDB) [中塚ほか，2005] のデータ* をそのコード名で指定する． 標準形式磁気異常分布グリッドデータは，磁気異常値とその観測高度 の両方の情報を含むが，プログラム altx は，その標準形式グリッドデータ から観測高度のデータをグリッドデータとして取り出す． プログラム adjlv は，グリッドデータのDCレベルの調整のため， 一定値を全グリッドデータ値に加算する． プログラム gadd／gsub は 2つのグリッドデータの加算処理／ 減算処理を行って，新たなグリッドデータを作成する． プログラム gtrim は，入力グリッドデータに対して，参照グリッドデータの 未定義データ領域と同じ範囲のグリッドデータ値を未定義値に置換えて， 新たなグリッドデータとして出力する．
プログラム govlay，gojoin，gmerge は， 多数のグリッドデータを統合した新たなグリッドデータを作成する． govlay は 順次スリット入りで重ね合わせ，gojoin は 一定の遷移ゾーン幅の設定で 後から重ね合る側を優先して重ね合わせ接合し， gmerge は 両者の重複部分を遷移ゾーンに設定して重ね合わせ接合する． プログラム txproj グリッドデータを別の図法展開座標系に変換 （再グリッディング）する． プログラム gtopo は DEM40**データ から地形高度グリッドデータを作成し， プログラム gtrf は 残差グリッドデータの標準磁場モデルを更新する． （残差の再計算に相当する．）
プログラム plmap，plmapc，plmapg，plmaps， plmapcs，shade の各プログラムは，高度分布のグリッドデータ として与えられた 3D曲面の表示を行う． plmap は 線画コンター図，plmapc は カラー段彩つきコンター図， plmapg は グレイスケール段彩図，plmaps は 陰影付きコンター図， plmapcs は 陰影付きカラー段彩図，shade は 陰影図 の表現となり， いずれもA4サイズ用紙に作図する． プログラム plmapl／plmapcl は，測線位置表示つきで 線画コンター図／カラー段彩つきコンター図 をA4用紙に描く．
プログラム xplmap／xplmapc は，グリッドデータの値分布を 線画コンター図／カラー段彩つきコンター図 で表現し，各種整飾等を施した図面 （B0サイズまで）を作成する．
* AMDB のグリッドデータは，
解凍済のものを ディレクトリ /pub/AMDB/DATA/grd の下に "[コード名].gd"
のファイル名で収容しておく必要がある．
** DEM40 とは，
日本の地表高データのファイル群であり，国土地理院発行の数値地図50mメッシュ
(標高) の CD-ROM 3枚 (日本-Ⅰ，日本-Ⅱ，日本-Ⅲ) を用いて構成したものである．
元データの著作権は国土地理院にあり，ここでは DEM40データを作成する手順を
genDEM40j.html に記す．
ANAM プログラム群は，磁気異常分布データの高度リダクションおよび定量解釈
のプログラムで構成され，以下のプログラムがある．
emag, emagf, amag, amagc, cmag,
cmagf, plamag, plamagc, tmcorr, tmcfix,
lcecorr, aaptdp, galtf, galts, emeq,
ameq, ameqc, cmeq, emeqs, ameqs,
ameqsc, cmeqs, rpmeqs, edeq, adeq,
adeqc, cdeq, edeqs, adeqs, adeqsc,
cdeqs, rpdeqs, calmas.
プログラム emag，emagf，amag，amagc， cmag，cmagf は，磁化強度マッピング [Okuma et al., 1994; Nakatsuka, 1995] の処理に用いる． まず emag／emagf で 感度係数マトリクス COEF の計算を行い， 次に amag／amagc で CG法による磁化強度マッピングの逆解析処理 を行う． その結果を用いて，cmag／cmagf で指定高度面での合成磁気異常を 計算することができる． ソースモデルは，emag／cmag では 形状が矩形ブロックで 近似されるが，emagf／cmagf では 地形データの分解能に合わせて 上面の凹凸を考慮する． CG法の繰返し処理の打切りは，amag では ループ回数で制御され， amagc では 解の収束状況で判定される． プログラム plamag／plamagc は，磁化強度マッピング結果を 線画コンター図／カラー段彩つきコンター図 で表示する．その際，周辺ソース域は マスクして表示しない．
プログラム tmcorr は，観測空中磁気異常グリッドデータに対して， 地形の一様磁化を仮定してその効果を見積り，それを除去する． また，プログラム tmcfix は，指定の磁化強度の地形効果を除去する． プログラム lcecorr は，直流電車軌道ループ電流の補正処理を行い， プログラム aaptdp は，観測空中磁気異常グリッドデータに対して， 点ダイポールソースによる半自動モデリングを行うもので，いずれも 神戸-京都地域の空中磁気データの解析 [中塚ほか，2004] で用いたものである．
プログラム galtf は，(観測高度の図化表示の目的で) StdLIN測線データ から観測高度を補間してグリッドデータを作成する． また プログラム galts は，(高度リダクション処理の基準高度 に用いる目的で) StdLIN測線データから観測高度を平滑化した高度グリッドデータ を作成する．
プログラム emeq，ameq，ameqc，cmeq は， 磁気異常グリッドデータに対して，プログラム edeq，adeq， adeqc，cdeq は，磁気異常ランダム点データに対して， 等価アノマリを用いた計算法 [Nakatsuka and Okuma, 2006a, 2006b] に基づいた 高度リダクション処理に用いる． まず emeq／edeq で 感度係数マトリクス CMUP／CFUP の計算を行い， 次に ameq／adeq または ameqc／adeqc で CG法による等価アノマリ導出の逆解析計算を行い， cmeq／cdeq で 高度リダクション結果を得る． CG法の繰返し処理の打切りは， ameq／adeq では ループ回数で制御され， ameqc／adeqc では 解の収束状況で判定される．
プログラム emeqs，ameqs，ameqsc，cmeqs， rpmeqs は，磁気異常グリッドデータに対して，プログラム edeqs， adeqs，adeqsc，cdeqs，rpdeqs は， 磁気異常ランダム点データに対して，磁化方向を有する等価ソースを用いた計算法 [Nakatsuka and Okuma, 2006a, 2006b] に基づいた 高度リダクション処理に用いる． まず emeqs／edeqs で 感度係数マトリクス CMUPS／CFUPS の計算を行い，次に ameqs／adeqs または ameqsc／ adeqsc で CG法による等価ソース導出の逆解析計算を行い， cmeqs／cdeqs で 高度リダクション結果 を， rpmeqs／rpdeqs で 極磁力変換した高度リダクション結果 を得る． CG法の繰返し処理の打切りは， ameqs／adeqs では ループ回数で制御され， ameqsc／adeqsc では 解の収束状況で判定される．
上記のプログラムの開発は，GNU Compiler Collection [http://www.gnu.org/software/gcc/] のインストールされた Linux PC の上で 行っており，ほぼすべてのプログラムは，FORTRAN言語で記述されている． alog2asc といくつかのサブプログラム（FORTRANプログラムから共用的に 使用される）だけは C言語で記述されている． プログラムは，地質調査所（現 産総研地質調査総合センター）での永年の 空中磁気探査の研究の中で段階的に構築されてきたもので，実データを用いた 動作確認を経たものである． しかし，ソースコードはその後も，処理の効率化や操作性向上を志向した改良のため， さらに改訂が行われている． 従って，何がしかの気づいていないバグが残っている可能性は否定できないが， ソースコードを公開しているので，ユーザが自ら容易にバグの修正を行えるものと 考える．
これらのプログラムの中では，中塚 [2009] によるサブルーチンライブラリ 'libgm' の機能が，極めて頻繁に利用されているので，プログラムのコンパイル のためには，サブルーチンライブラリを事前にインストールしておくことが 必須である．
大部分のプログラムでは配列変数が使用され，必要な配列サイズは取り扱う 実データによって当然異なる．（十分すぎる配列サイズの定義は，正しい結果を 得る上で支障はないものの，メモリ資源の浪費になる．） そのような配列サイズの定義は，変更が容易なように，ソースコード中で しばしば parameter文 で記述されている． その定義内容は，異なるプログラム間では必ずしも一貫性が取れていない可能性 もある． ユーザは，その必要に応じて，それらパラメータを設定変更されたい．
プログラムの実行にあたっては，作業ディレクトリの設定とともに いくつかのパラメータ指定を行う必要がある． これらプログラムでは，作業ディレクトリの設定と処理パラメータ入力の支援 のための LWKDIRサブプログラム と それに付帯する 'opnpin' メカニズムの利用が 取り入れられている． その利用には，2つの利点がある． 一つは，会話形式処理のパラメータ入力において，参考情報や選択肢が表示でき プログラムの実行操作にとくに有用である点， もう一つは，全く同じプログラムが 'opnpin' メカニズムを利用した 全パラメータ事前設定の非会話形式処理にも利用できる点である． 後者のスタイルのプログラムの実行は，「Job制御ファイル」を用意して， サブルーチンライブラリ 'libgm' に付属のユーティリティプログラム job / job1 を起動することで行える． 'opnpin' メカニズムを利用した 一連のジョブの実行方法 の詳細については， 地調研究資料集 no.518 [中塚，2009] の「各種ユーティリティプログラム」 を参照されたい．
FORTRAN言語プログラムのソースファイルは，.f のファイル拡張子，
C言語プログラムのものは .c のファイル拡張子がついている．
それらプログラムを UNIX環境でコンパイルするときに利用可能なスクリプト
@mkall が，各ディレクトリの下に用意してある．
但し，サブルーチンライブラリ 'libgm' が，/home/SHARE/lib
ディレクトリ内の libgm.aファイルとしてインストールされていること，
コマンド名 fort／gcc が希望の FORTRAN言語／C言語コンパイラ
の別名（または本名）となっていること を想定している．
三つのプログラム群 DPAM・GDMP・ANAM の それぞれについて，説明書HTML が用意してあり，各個別プログラムの短い説明と， 実行時に必要なパラメータが記されている．
この資料に添付の CD-ROM には，HTMLファイル・FORTRAN言語／C言語の ソースプログラム と 非会話形式処理の場合の「Job制御ファイル」のテンプレート が，下記のツリー構造で収容してある．
<<< この資料に含まれるファイルのツリー構造 >>>
(root)
|
+-- 0519index.html 表紙ページ HTML
+-- indexj.html この資料の概要説明（この文書）
+-- dpamj.html DPAMプログラム群の説明
+-- gdmpj.html GDMPプログラム群の説明
+-- anamj.html ANAMプログラム群の説明
+-- genDEM40j.html DEM40 ファイル群を生成する方法の説明
|
+-- index.html この資料の概要説明（英語版）
+-- dpam.html DPAMプログラム群の説明（英語版）
+-- gdmp.html GDMPプログラム群の説明（英語版）
+-- anam.html ANAMプログラム群の説明（英語版）
|
+-- dpam.tgz dpamディレクトリの圧縮アーカイブ
+-- gdmp.tgz gdmpディレクトリの圧縮アーカイブ
+-- anam.tgz anamディレクトリの圧縮アーカイブ
+-- dpam_tp.tgz dpam_tpディレクトリの圧縮アーカイブ
+-- gdmp_tp.tgz gdmp_tpディレクトリの圧縮アーカイブ
+-- anam_tp.tgz anam_tpディレクトリの圧縮アーカイブ
|
+-- dpam/ (DPAMプログラム群のソースを含むディレクトリ)
| |
| +-- @mkall (全DPAMプログラムをコンパイルするスクリプト)
| +-- alog2asc.c, xldam.f, xldpn.f, xldhg.f, despike.f,
| dvcorr.f, sml5.c, ecomp.f, fcomp.f, gsurf.c,
| ggrid.f, ggrids.f, pframe.f, pltrk.f,
| pchkdv.f, pchkmag.f, pchkres.f, pchkcomp.f,
| xslin.f, xslina.f
|
+-- gdmp/ (GDMPプログラム群のソースを含むディレクトリ)
| |
| +-- @mkall (全GDMPプログラムをコンパイルするスクリプト)
| +-- sel.f, seldb.f, altx.f, adjlv.f, gadd.f, gsub.f,
| gtrim.f, govlay.f, gojoin.f, gmerge.f, txproj.f,
| gtopo.f, gtrf.f, plmap.f, plmapc.f, plmapg.f,
| plmaps.f, plmapcs.f, shade.f, plmapl.f,
| plmapcl.f, xplgobj.c, xplmap.f, xplmapc.f
|
+-- anam/ (ANAMプログラム群のソースを含むディレクトリ)
| |
| +-- @mkall (全ANAMプログラムをコンパイルするスクリプト)
| +-- emag.f, emagf.f, amag.f, amagc.f, cmag.f, cmagf.f,
| plamag.f, plamagc.f, tmcorr.f, tmcfix.f, lcecorr.f,
| aaptdp.f, galtf.f, galts.f, emeq.f, ameq.f, ameqc.f,
| cmeq.f, emeqs.f, ameqs.f, ameqsc.f, cmeqs.f,
| rpmeqs.f, edeq.f, adeq.f, adeqc.f, cdeq.f, edeqs.f,
| adeqs.f, adeqsc.f, cdeqs.f, rpdeqs.f calmas.f
|
+-- dpam_tp/ (DPAM群プログラムのテンプレートを含むディレクトリ)
| |
| +-- alog2asc.tp, xldam.tp, xldpn.tp, xldhg.tp, despike.tp,
| dvcorr.tp, ecomp.tp, fcomp.tp, ggrid.tp, ggrids.tp,
| pframe.tp, pltrk.tp, pchkdv.tp, pchkmag.tp,
| pchkres.tp, pchkcomp.tp, xslin.tp, xslina.tp
|
+-- gdmp_tp/ (GDMP群プログラムのテンプレートを含むディレクトリ)
| |
| +-- sel.tp, seldb.tp, altx.tp, adjlv.tp, gadd.tp, gsub.tp,
| gtrim.tp, govlay.tp, gojoin.tp, gmerge.tp, txproj.tp,
| gtopo.tp, gtrf.tp, plmap.tp, plmapc.tp, plmapg.tp,
| plmaps.tp, plmapcs.tp, shade.tp, plmapl.tp,
| plmapcl.tp, xplmap.tp, xplmapc.tp
|
+-- anam_tp/ (ANAM群プログラムのテンプレートを含むディレクトリ)
|
+-- emag.tp, emagf.tp, amag.tp, amagc.tp, cmag.tp, cmagf.tp,
plamag.tp, plamagc.tp, tmcorr.tp, tmcfix.tp, lcecorr.tp,
aaptdp.tp, galtf.tp, galts.tp, emeq.tp, ameq.tp,
ameqc.tp, cmeq.tp, emeqs.tp, ameqs.tp, ameqsc.tp,
cmeqs.tp, rpmeqs.tp, edeq.tp, adeq.tp, adeqc.tp,
cdeq.tp, edeqs.tp, adeqs.tp, adeqsc.tp, cdeqs.tp,
rpdeqs.tp, calmas.tp