| GSJ Open-File Report, no. 519 | December, 2009 |
Tadashi Nakatsuka (Institute of Geology and Geoinformation, GSJ, AIST)
I have been engaged in the research on aeromagnetic survey, and developed various kinds of programs for data processing, analysis and interpretation. The outlines of them have already been reported on the scientific journals and at the scientific meetings [Nakatsuka and Okuma, 2005a, 2005b, 2006a, 2006b; etc.]. Their source codes are open to public through the GSJ Open-file Reports.
In this report, I deal with three groups of programs for the aeromagnetic
survey, as follows.
(1) DPAM:
programs for aeromagnetic survey data processing,
(2) GDMP:
programs for grid data manipulation,
(3) ANAM:
programs for analysis of magnetic anomaly data.
This software system is in a style of program libraries, consisting of many
individual programs of rather simple functions.
Then the actual process to some target/objective would require a series of
program executions.
The data formats used in this system are customized for our own original
style.
However, as for its design, all data are handled in forms of ASCII text,
so that the user can easily convert data to/from another style of formats.
This report is a revised version for the previous GSJ Open-File Report no.449 [Nakatsuka, 2007]. The points of revision are:
DPAM group consists of programs for aeromagnetic survey data
processing as follows.
alog2asc, xldam, xldpn, xldhg, despike,
dvcorr, ecomp, fcomp, ggrid, ggrids,
pframe, pltrk, pchkdv, pchkmag, pchkres,
pchkcomp, xslin and xslina.
The program alog2asc converts aeromagnetic binary raw data (specialized for our survey system) into ASCII file format, the programs xldam and xldpn generate line data (DPAM line data) file from aeromagnetic ASCII raw data, making use of Real-time GPS and PostNav processed DGPS position data, respectively, and the program xldhg generates similar line data (HGAM line data) file from helicopter-borne gradiometer aeromagnetic data acquired. In addition, xldam, xldpn and xldhg also calculate the IGRF residuals, using the subprogram "igrf" by Nakatsuka [2005].
The program despike eliminates spike noises in magnetic field values (total magnetic field and IGRF residual) of DPAM line data, the program dvcorr performs the diurnal correction of magnetic field data making use of Ground-Station magnetic field data, and the programs ecomp and fcomp execute the compensation of aircraft's magnetic effect for compensation flight DPAM line data and for survey flight DPAM line data, respectively, using the method described by Nakatsuka and Okuma [2005b].
The programs ggrid and ggrids generate grid data from random points data (i.e., DPAM line data or equivalent) by using the method of "continuous curvature splines in tension" as developed by Smith and Wessel [1990]. Only the magnetic field (IGRF residual) data are converted into grid by ggrid, while ggrids applies the method also on the altitude data.
The program pframe draws an illustration of coordinates framework of survey area, for easier setting of grid location parameters, and the program pltrk plots trackline paths from DPAM line data. The programs pchkmag, pchkres and pchkcomp are used to plot magnetic field data in order to inspect the data quality. Total magnetic field data is plotted in black by pchkmag, and IGRF residual data is plotted in blue by pchkres, while pchkcomp plots IGRF residual magnetic data both before and after the compensation of aircraft's magnetic noise in red and blue colors.
The programs xslin and xslina extract simple standard line (StdLIN) data from random point data. Various formats of random point (line) data are supported by xslin as input data, while xslina deals with only DPAM line data but with the function of averaged re-sampling.
GDMP group consists of programs for grid data manipulation:
sel, seldb, altx, adjlv, gadd, gsub,
gtrim, govlay, gojoin, gmerge, txproj,
gtopo, gtrf, plmap, plmapc, plmapg,
plmaps, plmapcs, shade, plmapl, plmapcl,
xplmap and xplmapc.
The programs sel and seldb generate new grid data from an existing grid data with re-gridding. For sel program, the existing grid data should be assigned as the input filename parameter, while AMDB grid data* (Nakatsuka et al., [2005]) is selected by the code-name input for seldb. The standard format magnetic anomaly grid data includes both information of magnetic anomalies and the altitudes of observation surface. The program altx extracts the altitude information into a standard grid data format. The program adjlv adds a constant value to all grid data in order to adjust DC level, the program gadd adds 2 sets of grid data, and the program gsub subtracts the 2nd set grid data values from the 1st set grids, as generating a new grid data. The program gtrim replaces values of input grid data with NULL (undefined) data value for the range where the values of reference grid data are NULL data value, as generating a new grid data.
The program govlay overlays multiple grid data sequentially with placing slit zone, the program gojoin joins multiple grid data sequentially with placing transient zone surrounding the overlaid data, and the program gmerge merges multiple grid data sequentially with placing transient zone along overlapping zone, to create a new grid data. The program txproj translates grid data into another map projection with re-gridding, the program gtopo generates topographic altitudes grid data from DEM40** (Digital Elevation Model 40m) data, and the program gtrf translates IGRF residual grid data into another reference model.
The programs plmap, plmapc, plmapg, plmaps, plmapcs and shade produce illustrations of 3D surface given as altitude distribution of grid data, by line contour map (plmap), color-graded contour map (plmapc), gray-scale grading map (plmapg), shaded-relief contour map (plmaps), shaded-relief color-grading map (plmapcs) and shaded-relief map (shade), all on an A4 sheet. The program plmapl generates line contour map with drawing trackline paths, and the program plmapcl generates color-graded contour map with trackline paths, both on an A4 sheet.
The programs xplmap and xplmapc create line contour map and color-graded contour map of grid data, respectively, with placing various supplementary items on a selectable size of sheet up to B0.
* AMDB grid files should be
recovered from TGZ file and stored in the directory /pub/AMDB/DATA/grd
with filenames as "code-name.gd".
** DEM40 is a group of files
of terrain elevation data in Japan, extracted from Digital Maps (50m mesh
elevation, Japan-I, Japan-II and Japan-III) by the Geographical Survey
Institute, Japan (GSI).
The source data are subject to GSI's copyright. The procedure how to generate
DEM40 files from GSI's Digital Maps is presented in
genDEM40j.html (only in Japanese).
ANAM group consists of programs for the altitude reduction and quantitative
interpretation of magnetic anomaly data:
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 and calmas.
The programs emag, emagf, amag, amagc, cmag and cmagf are used for the magnetization intensity mapping [Okuma et al., 1994; Nakatsuka, 1995]. First, emag / emagf calculates the contribution coefficients matrix COEF, then amag / amagc executes inversion process of magnetization intensity mapping by the method of conjugate gradients (CG), and cmag / cmagf may be used to calculate the synthetic magnetic anomalies on a desired altitude surface by using the result of amag / amagc. The source model is approximated into rectangular blocks in emag / cmag, while the surface undulation is considered in emagf / cmagf. The iteration of CG method is controlled by the number of loop count (amag) or by convergency tolerance (amagc). The programs plamag and plamagc plot the result of the magnetization intensity mapping in the form of line contour map and color-graded contour map, respectively, with masking the surrounding source region.
The program tmcorr executes the correction of topographic effect on the magnetic anomaly grid data, assuming a uniform terrain magnetization, while the program tmcfix executes the correction of topographic effect with the assigned value of terrain magnetization. The program lcecorr estimates the railway loop-current effect in the observed magnetic anomaly grid data, and correct for it. The program aaptdp is for the semi-automatic modeling by point-dipole sources for observed magnetic anomaly grid data. These processes were used in the data interpretation of aeromagnetic data in the Kobe-Kyoto area [Nakatsuka et al., 2004].
The program galtf interpolates observation altitudes of StdLIN data into grid data (for the purpose of illustrating observation altitudes), and the program galts generates grid data of smoothed observation altitude from StdLIN data (for the purpose of using as the reference in altitude reduction).
The programs emeq, ameq, ameqc and cmeq are for the altitude reduction from grid data, and the programs edeq, adeq, adeqc and cdeq are for the altitude reduction from random point data, by equivalent anomaly method [Nakatsuka and Okuma, 2006a, 2006b]. First, emeq or edeq calculates the contribution coefficients matrix CMUP or CFUP, and ameq / ameqc or adeq / adeqc executes inversion process of equivalent anomaly derivation by CG method. The iteration of CG method is controlled by the number of loop count in ameq / adeq, or by convergency tolerance in ameqc / adeqc. Then the reduction calculation is performed by cmeq or cdeq program.
The programs emeqs, ameqs, ameqsc, cmeqs and rpmeqs are for the altitude reduction from grid data, and the programs edeqs, adeqs, adeqsc, cdeqs and rpdeqs are for the altitude reduction from random point data, by equivalent source magnetization method [Nakatsuka and Okuma, 2006a, 2006b]. First, emeqs or edeqs calculates the contribution coefficients matrix CMUPS or CFUPS, and ameqs / ameqsc or adeqs / adeqsc executes inversion process of deriving equivalent source magnetization by CG method. The iteration of CG method is controlled by the number of loop count in ameqs / adeqs, or by convergency tolerance in ameqsc / adeqsc. Then the reduction calculation is performed by cmeqs or cdeqs program, and the reduction-to-pole anomaly on that surface is calculated by rpmeqs or rpdeqs program.
The program calmas calculates theoretical magnetic anomaly distribution on the specified observation surface, caused by a simple source model.
I have developed these programs on a Linux machine, on which the GNU Compiler Collection [http://www.gnu.org/software/gcc/] is installed, and all programs are coded in FORTRAN language, except alog2asc program and a few subprograms (for the common use from FORTRAN programs) are coded in C language. The programs were developed step by step during the years of research activities on aeromagnetic survey at the Geological Survey of Japan, and they were verified to function well by some practical data. The program codes, however, have been revised further for the improvement on reasonable data handling or easier operation. So there is a little possibility that some unaware bugs are left. But, I convince the user will easily correct bugs, as the source codes are now presented.
In these programs, the subroutine library 'libgm' by Nakatsuka [2009] is utilized quite frequently. So the installation of the subroutine library 'libgm' in advance is required to compile these programs.
Most programs use some array variables, and the necessary size (dimension) of the array is dependent on the actual data to be dealt with. (Although too large size definition does not harm the correct result, the memory resource may be wasted.) The definitions of such dimensions are often given by using 'parameter' statement in source codes, to enable easier adjustment. Actually the definitions of such parameters may not be consistent among programs even in the same group. The user is requested to customize them according to his / her requirement.
When executing these programs, a few to several or more parameters have to be specified, along with selecting working directory path. In these programs, the LWKDIR (Assist to set-up working dir. and process parameters, etc.) function with 'opnpin' mechanism of the subroutine library 'libgm' mentioned above is implemented. The use of this function provides us two merits. The one is that the parameter input is prompted with informative message or valid choice, which is very useful to learn how to operate the program. The other is that the same program is used to the non-interactive process for predetermined parameters with the help of 'opnpin' mechanism. The latter style of execution of program is realized by creating 'Job Control File' and submitting job / job1 command of the Utility program appended to the subroutine library. Refer to the section of Job Control File, "Misc. Utility Program" of GSJ Open-file Report no.518 [Nakatsuka, 2009], with respect to the detail of executing a series of job steps with the help of 'opnpin' mechanism.
Source programs coded in FORTRAN language have the filenames with .f
(specifier) extensions, and those in C language with .c (specifier)
extensions.
To compile these programs under UNIX environment, @mkall script in each
directory can be used, provided that the subroutine library 'libgm' by
Nakatsuka [2005] is installed as the archive file libgm.a under the
directory of /home/SHARE/lib, and the command names
'fort' and 'gcc' is aliased to (or the name of) the desired
FORTRAN and C compilers.
There are also HTML documentations for each of three program groups, DPAM, GDMP and ANAM. The explanation of the function of individual program is given briefly, and the parameters required to execute the program are described.
In the CD-ROM appended to this report, HTML files, FORTRAN and C source files, and templates of Job Control files for non-interactive execution are included, in the following tree structure.
<<< Tree structure of files included in this report >>>
(root)
|
+-- 0519index.html Cover page HTML
+-- index.html Introduction HTML in English (This document)
+-- dpam.html English Documantation for DPAM programs
+-- gdmp.html English Documantation for GDMP programs
+-- anam.html English Documantation for ANAM programs
|
+-- indexj.html Introduction HTML in Japanese
+-- dpamj.html Japanese Documantation for DPAM programs
+-- gdmpj.html Japanese Documantation for GDMP programs
+-- anamj.html Japanese Documantation for ANAM programs
+-- genDEM40j.html (in Japanese) How to generate DEM40 file
|
+-- dpam.tgz Gzipped tar archive of the directory 'dpam'
+-- gdmp.tgz Gzipped tar archive of the directory 'gdmp'
+-- anam.tgz Gzipped tar archive of the directory 'anam'
+-- dpam_tp.tgz Gzipped tar archive of the directory 'dpam_tp'
+-- gdmp_tp.tgz Gzipped tar archive of the directory 'gdmp_tp'
+-- anam_tp.tgz Gzipped tar archive of the directory 'anam_tp'
|
+-- dpam/ (Directory containing DPAM source programs)
| |
| +-- @mkall (Script to compile all DPAM programs below)
| +-- 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/ (Directory containing GDMP source programs)
| |
| +-- @mkall (Script to compile all GDMP programs below)
| +-- 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/ (Directory containing ANAM source programs)
| |
| +-- @mkall (Script to compile all ANAM programs below)
| +-- 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/ (Directory containing templates for DPAM programs)
| |
| +-- 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/ (Directory containing templates for GDMP programs)
| |
| +-- 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/ (Directory containing templates for ANAM programs)
|
+-- 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