Magnetic anomalies caused by two-dimensional structure : their computation by digital computers and their interpretation

Columbia University in the City of New York

LAMONT GEOLOGICAL OBSERVATORY

PALISADES. NEW YORK

\

MAGNETIC ANOMALIES
CAUSED BY TWO-DIMENSIONAL STRUCTURE:
THEIR COMPUTATION BY DIGITAL COMPUTERS

AND

THEIR INTERPRETATION

by

J. R. Heirtzler, G. Peter, M. Talwani and E. G. Zurflueh

Technical Report No. 6

CU'6'62 Nonr-Geology

A

October 1962

V-v.

LAMONT GEOLOGICAL OBSERVATORY
(Columbia University)
Palisades, New York

MAGNETIC ANOMALIES CAUSED BY TWO-DIMENSIONAL STRUCTURE:
THEIR COMPUTATION BY DIGITAL COMPUTERS

AND

THEIR INTERPRETATION

by

J.R. Heirtzler, G. Peter, M. Taiwan! and E.G. Zurflueh

Technical Report No. 6
CU-6-62-Nonr - Geology

October 19&2

CONTENTS

Page

CHAPTER I Introduction

by J.R. Heirtzler and E.G. Zurflueh 1-1

CHAPTER II The Mathematical Expression for the

Magnetic Anomaly over a Two-Dimensional

Body of Polygonal Cross Section

by M. Talwani and J.R. Heirtzler 2-1

A* Derivation of the Basic Formulas 2-1

B. Induced and Remanent Magnetization

and the Total Intensity Anomaly 2-10

CHAPTER III Computed Magnetic Total Intensity Curves

over Two-Dimensional Bodies

by G. Peter 3-1

A. Description of the Model 3-1

B. Discussion of the Computed Curves 3-2

CHAPTER IV A Method of Interpretation for Anomalies

of Total Magnetic Intensity caused by
Two-Dimensional Bodies

by E.G. Zurflueh 4“1

A* General Remarks 4-1

B. Mathematical Expressions 4"3

C. Description of Method 4“10

D. Results of Practical Applications 4"l5

E. Curves and Graphs l\.-2b

APPENDIX A Digital Computer Program for the

Computation of Magnetic Anomalies over

Two-Dimensional Polygonal Bodies

by M, Talwani and J.R. Heirtzler A-l

A. Framing the Problem A-l

B. Flow Diagram A-4

C. A Simple Example A-10

ACKNOWLEDGEMENTS R-l

REFERENCES R-l

i

Digitized by the Internet Archive
in 2020 with funding from
Columbia University Libraries

https://archive.org/details/magneticanomalieOOheir

Chapter I

Introduction

J.R. Heirtzler and E.G. Zurflueh

In a sense this is a "housekeeping” report in that
it brings together, in a more or less orderly fashion, for our
own use some of the analysis techniques that have been used by
various persons at the Lamont Geological Observatory, It is
felt that such material may be of interest to individuals in¬
itiating an analysis of the geomagnetic field intensity for
purposes of determining earth structure.

In the measurements of the geomagnetic field intensity
in marine areas by Lamont it has seldom been possible to take
measurements over a systematic grid and construct magnetic contour
charts. Measurements usually consist of a magnetic profile
along the track of a research vessel, these measurements being
taken with a towed magnetometer that records the total intensity
of the geomagnetic field. Over a period of time our ships and
other ships and aircraft have made sufficient measurements in
some areas to permit a general description of the magnetic a-
nomalies even though track spacing and heading prohibit contour¬
ing, Frequently such a general description can substantiate or
refute geophysical hypotheses and suggest guide-lines for future
investigations.

An attempt to deduce exact geophysical information
from a limited magnetic profile is recognized as a highly un¬
certain undertaking, especially when control of the time

1-1

1-2

variations is doubtful. Alternatively, analyzing a high densi¬
ty of profiles on a statistical basis may be of value if no
errors in basic physics are committed and the limitations of
the techniques are understood.

The structures causing magnetic anomalies can fre¬
quently be treated as cylinderical, or two-dimensional, i.e.
extending to plus and minus infinity in a direction parallel
to a coordinate axis. Two-dimensional structure seems to be
particularly common in certain marine areas, see for example
(Vacquier, et.al., 1961 ) and (Heirtzler, et.al., 1962 ). It is
instructive in such cases to determine some of the configu-
ations of magnetic materials that can produce anomalies similar
to those observed. An electronic digital computer program to
calculate the anomaly that would be caused by an assumed two-
dimensional structure has been used by this institution for
several years. It was originally written as a program for the
IBM 6£0 by M. Landisman of Lamont. It was written again for
the IBM 1620, because of the convenience of that computer. It
is the formulation of this second program that is given here
but, naturally, both programs give the same answers.

In order to give some crude feeling of how the various
geometrical factors affect the total field anomalies a brief
selection of anomalies for simple bodies is given in Chapter III.
With the large number of geometrical parameters that are important
in this calculation it is clearly impossible to give a complete
handbook of anomalies with one parameter varying at a time. The

1-3

relative importance of the parameters is discussed in Chapters
III and IV,

In Chapter IV a close examination is made of the calcu¬
lated anomalies for those bodies with vertical sides, A working
set of rules is presented to determine the configuration of such
bodies from observed profiles. In that Chapter it is also
shown how far mathematical theory can be carried without having
to make very specialized assumptions about body configuration.

In connection with the presentation of interpretational
techniques in this paper, a brief listing of some of the pertinent
literature is given here with short comments on the contents
of each publication.

Several authors give sets of model anomalies, usually
accompanied by rules for depth determination and formulas for
the particular model used. In this way Haalck (192?)* Gulatee
(1938) and Heiland (194&) compute model anomalies of the hori¬
zontal and vertical components for different bodies. Nettleton
( 1942 ) determines the vertical intensity anomaly for various
geometric configurations. Henderson and Zietz (1948) treat the
effect of point and line sources with respect to total intensity
measurements. In 1958 the same authors published a paper on
the use of magnetic doublets in the interpretation of total in¬
tensity anomalies. Vacquier, et.al. (1951) interpret total in¬
tensity maps with the help of three-dimensional (prismatic)
models and computed second vertical derivatives. Zietz and
Henderson (1954) use three-dimensional double layer model fields

i-4

for the analysis of total intensity data. Baranov (1955)
standardizes his model calculations by reducing measurements
for different inclinations to vertical field conditions.

Smellie (1956) gives total intensity depth factors for point
and line configurations of poles and dipoles. From general
mathematical considerations Smith (1959 and 1961 ) derives in¬
equalities which can be used to estimate depth and intensity
of magnetization for arbitary bodies.

Before the advent of the airborne magnetometer the
horizontal or vertical components were measured almost ex¬
clusively and the model anomalies were treated accordingly.

Since about 1946 total intensity measurements are used in
general and the newer literature usually refers to total in¬
tensity anomalies. According to potential theory it is possible
to calculate any component from measurements of any other com¬
ponent if certain conditions are fulfilled (for example see
symmetry relations on page 2-9)• A group of papers deal with
the conversion of one component into another and with the re¬
lated problems of computing derivatives and upward and downward
continuation.

Skeels (1947) was the first to explicitly point out
the general possibility of the above mentioned conversions by
surface integration. Vestine (194^) and Vestine and Davids
( 1945 ) also point out the importance of integral methods in mag¬
netic interpretation. Hughes and Pondrom (1947) describe the
calculation of the vertical intensity out of total intensity

1-5

data, Skeels and Watson (1949) treat the problem of conversions
in a more general way. Henderson and Zietz (1949a) give a
method for the computation of second vertical derivatives and
in another peper ( 1949 b) they analyze the upward continuation
of the field. Peters (1949) treats the problems of upward and
downward continuation and the calculation of derivatives. He
also derives methods for the direct determination of basement
structure. Henderson (i 960 ) relates the different operations
of conversion to one basic formula and uses expressions suita¬
ble for use on electronic computers.

Another group of papers is concerned with the inverse
problem of the calculation of magnetic anomalies caused by
bodies of any given shape. Gassmann (1951) describes a graphi¬
cal method of effecting the triple integration necessary for
the solution of this problem. Baranov (1953) works out another
graphical method with a different sequence of integrations.
Henderson and Zietz (1957) simplify the computation of a total
intensity anomaly by projecting a map of a given body in the
direction of the earth*s field and using a modified Gassmann
process for the integration. Talwani and Ewing (i 960 ) calculate
the anomaly by dividing the body into parallel and polygonal
laminae, similar to the treatment in this report. Their method
is put in a form which can be used for digital computers.

Vajk (1951) and Nettleton (1954) analyze the first
step in magnetic anomaly interpretation, namely the separation

1-6

of anomalies from an assumed regional field. A summary of the
different interpretation methods is given by Kohler (1958)*

His work also contains a rather complete list of references
which includes papers from the large Russian literature on this
subject. Lastly, there are extensive, company classified tech¬
niques used by the geophysical exploration companies.

Chapter II

The Mathematical Expression for the Magnetic Anomaly
over a Two-Dimensional Body of Polygonal Cross Section

M* Talwani and J.R. Heirtzler

A* Derivation of the Basic Formulas

It is characteristic of geophysics texts to derive the
expressions for the horizontal and vertical magnetic anomalies
by analogy with the expressions for gravitational anomalies* In
the present report the expressions for the horizontal, vertical,
and total field anomalies are derived from basic magnetostatic
theory*

Consider the volume element Ax, located in

an external magnetic field (Fig* 2-1).

Fig. 2-1

2-1

2-2

We assume that the volume element has a uniform mag¬
netization M. Let be the magnetic moment of the element*
Then

-j* * v

The magnetic potential, at the origin, due to the volume

element is

A7 ^ x v- Af ^ y -f* %

42*

where

p* = i X ^ j y -t i i :

At the origin the potential of a rod of cross section
and infinite in the y-direction is

Jo. O'*** y**+ ^j"** ^

— J Ax A 2 / *. 7

/ x 1 *-*'* /

The vertical magnetic field strength, V, is

3/*

2 *

ZX 2 ^ a - Mj(x '-Z x )

( + i L J 1

( 2 - 1 )

2-3

and the horizontal field strength, H, measured in the x direction
is

(V- -j

a z J (2-2)

The horizontal field strength measured in the y-direction is zero
since y is not contained in the expression for the potential. If
M is the induced magnetization this means, physically, that we
have not included demagnetization effects. The above represent
the magnetic anomalies caused by the rod, that is the regional
field strength is not included.

Consider the lamina shown in Fig. 2-2, infinite in the
y-direction and infinite in the positive x-direction.

M = -

d_P

= JL AX A*-

Fig. 2-2

2-4

By the use of eqs (2-1) and (2-2) one gets

\/ = 2 4%| M x zxa-H a (x*-eV

Cx*

Z A* I H a x i

t x u -(. *•>- J

(2-3)

H - 2Alj H^Zx- fc

'j ,Lx

24 %| 'V -<- M-, fc I

1 X'-H-i'- )

(2-4)

The variables x and z in (2-3) and (2-4) represent the coordinates
of the end of the lamina. If we had taken the lamina running in
the opposite direction, as indicated in Fig. 2-3

2-5

0

T

I

t

I

" s/ZZZZZZZZZZZZZZZk

Az

Y

*

Pig. 2-3

we would have obtained

v = Z&i ,

“* f’l x 2 : +• M, X

X*-h-

H =•

**"4- ^

2-6

Notice that V and H bear the opposite signs to eqs (2-3) and
(2-4)* Next consider the prism shown in Fig. 2-4, infinite in
the y direction.

Pig. 2-4

The expressions for V and H at the origin can be obtained by

2-7

the integration with respect to z in eqs (2-3) and (2-4):

(2-6)

The integration is carried out from the small to the larger 0
In eqs (2-3) and (2-4) the variables x and z represented the co
ordinates of the end of the lamina. For the prism x and z are
related by the equation for the line AD:

X = Cx , + ~ (co-t <$>) -fc j

If eqs (2-5) and (2-6) are integrated and the variables r and
used instead of x and z one obtains:

2-8

* - - - *[". f«. - •. ) p s*-** P ^~°f (** / *i ) j

M -i &,) *'—<#' - ce>j <t>tof(r L /r.)\l

SJ ( 2 - 7 )

W - 2. j ~4>[K*f(e> ir *,)sJ- < t>-n,,<t>t of (, l / r .)'l

*-M.

I c®, - 0 Jm.s t + S ^pJL.f (n/r,)}J

(2-8)

Had it been assumed that the prism was oriented as in
Pig, 2-5 the above would have been preceded by a minus sign.

0

w/e,

IX, - 2 - )

HU

Y

Fig. 2-5

2-9

In Figs. (2-4) and (2-5) the prism extended infinitely
in one of the x directions. It is clear that H and V for a body
of finite x dimensions can be treated by assuming that it consists
of two bodies (ref. Fig. 2-6). One must be careful with directions
of angular measurements and with signs in working with these
equations.

Fig. 2-6

The procedure indicated above can be extended to any
number of prismatic or polygonal bodies each with many faces and
each with its own magnetization. This derivation does not take
into account nonuniform magnetization such as would exist near
the corners of the bodies. In fact the actual geological bodies
probably do not have distinct comers and nonuniform magneti¬
zation is believed to be unimportant.

It is interesting to note the symmetrical relationships
of the factors in eqs (2-7) and (2-8). That symmetry shows that
vertical magnetization effects the horizontal anomaly to the same
extent that the transverse magnetization effects the vertical

t

2-10

anomaly and that the transverse magnetization effects the hori¬
zontal anomaly to the same extent that the vertical magnetization
effects the vertical anomaly,

B* Induced and Remanent Magnetization
and the Total Intensity Anomaly

If one wishes to utilize eqs (2-7) and (2-8) for the

calculation of the V and H anomalies due to induced magnetization

only he would write

M = k P

where k = magnetic susceptibility

F = undisturbed regional total intensity vector
and as illustrated in Pig, 2-7
Mj = k F x = k P cos I sin s
Mg = k P z = k P sin I

Pig. 2-7

2-11

I = Inclination, positive if F is below horizontal
(northern hemisphere)

s = angle of strike, measured from horizontal pro¬
jection of P and either positive or negative y-
axis. s must always be less than 180° (and thus
sin s positive) since we will specify that the co ¬
ordinate axes are oriented to make F x positive .
This orientation of the axes must be borne in mind
for the correct Interpretation of results .

If one wishes to utilize eqs (2-7) and (2-8) for the
calculation of the V and H anomalies due to a body with remanent
magnetization, and not to include induced anomalies, he would
write

M x = M rei n cos a sin b

Mj; — sin a

where

a = inclination of remanent magnetization vector below
horizontal

b = strike angle of remanent magnetization vector
measured between horizontal projection of M rem
and positive or negative y axis, less than 180°.

To apply the equations as written the axes must be
oriented to make M x positive.

If it is desired to use eqs (2-7) and (2-8) to get the
V and H anomalies due to the total magnetization (induced and

2-12

remanent) one should write

M* = “tot cos * 3in p

“y = “tot sin *

where

oc = inclination of total magnetization vector below hori¬
zontal

p = strike angle of total magnetization vector measured

between horizontal projection of M^ot and positive or
negative y axis, less than 180° To apply the equations
as written, the axes must be oriented to make positive.
Regardless of the assumptions used to determine V and
H, the anomaly in total intensity is obtained as follows:

Let T be the anomaly in total intensity, then

(P + T)2 = (F x + h)2 + (F y )2 + (f z + V)2

Since there is no anomaly in the y direction (ref, p, 2-3)*
Squaring, dropping terms that contain the square of V or H and
remembering that

Refering to Fig. 2-7

T = H cos I sin s + V sin I * (2-9)

2-13

An electronic digital computer program to calculate
V, H, and T is given in the Appendix.

Eq (2-9) illustrates that, if the anomalies are small
relative to P, then the anomaly in total intensity is the com¬
ponent of V and H in the Direction of P. See Fig. 2-8.

Pig. 2-8

The magnitude of T is composed of the projection of H
(H cos I sin s) and of V (V sin I) in the direction of P. This
geometrical picture is invalid if the anomaly is large (we
dropped squared terms in the derivation of eq 2-9) and a more
adequate vector diagram must be used.

Chapter III

Computed Magnetic Total Intensity Curves
over Two-Dimensional Bodies

G. Peter

A. Descr iption of the Model

With the aid of the computer program (see Appendix) a

few families of magnetic total intensity curves have been com¬
puted in order to show the changes in these curves, when the
ambient magnetic field directions and the geometric properties
of the model are varied. The different secs of curves have been
obtained by varying one parameter at a time. Grouping the curves
in this manner enables one to make a qualitative analysis of them
The curves are discussed in section B,

Pig. 3-1: Ground Plan and vertical
section of assumed two-dimensional model

x = line of measurement
d = depth of burial
D = vertical extent

s = strike
W = Width

The model used in these calculations is a two-dimension
al rectangular body, which has infinite length in a horizontal

3-1

3-2

direction (see fig. 3-1). It has vertical sides, except for the
group 6 curves, where the sides have ± 30 and ± 60 degree angles
from the vertical. For all calculations the model has uniform
isotropic induced magnetization of 0.001 emu, it is in the northern
hemisphere where the field strength is £0000 gammas. The pro¬
jection of the magnetic north on the profiles points toward the
right.

The total intensity values are computed for a

measurement (profile), which
the axis of the body.

The various curves

1. Variable inclination

2. Variable strike

3. Variable depth
4# Variable width

5* Variable vertical extent
6. Slanted bodies

is horizontal and at right

are grouped as follows:

I (figure 3-2)
s (figures 3-3 to 3-7)
d (figures 3-8 to 3-12)

W (figures 3-13 to 3-15)
D (figures 3-1& to 3-18)
(figures 3-19 to 3-22)

line of
angle to

These curves are only a selection for a few varied
parameters. Additional information can be found in the papers
of Haalck (1927), Gulatee (1938) and Heiland (1946).

B. Discussion of the Computed Curves

The following is a short discussion of the general
characteristics of the computed curves and of the effect of the
various magnetic and model parameters on them. Some elementary
characteristics are also mentioned. Conclusions are derived
mathematically in Chapters II and IV.

3-3

1 . Variable inclination I. (Figure 3-2).

The major features of the magnetic anomaly picture
due to induced magnetization in the northern hemisphere con¬
sist of a maximum on the south, and a minimum on the north
side of a structure. For a symmetrical model, if the strike
s = 90 °, (the strike s is defined as the angle between the
geological strike of the two-dimensional body and the magnetic
north) in the case of I = 45the size of the maximum and the
northern minimum are equal in amplitude, and the maximum and
minimum are symmetrical about the center of the model. For
smaller inclinations the northern minimum is the dominating
feature. Figure 3-2 also shows that there is a symmetry about
I = 45°; in other words for equal inclination differences from
I = 45° the amplitudes of the maximum and minimum are equal,
and their positions with the respect of the center of the body
are equidistant.

The above statements applied to a symmetrical model
are true for the negative inclinations (southern hemisphere)
except that the north and south sides of the picture are inter¬
changed.

The above mentioned symmetry around I = 45 ° exists be¬
cause of the symmetrical model used in the calculation, and is
due to the change in the relative importance of the horizontal
and vertical components of the magnetic field at I = 45°*

2. Variable strike s (Figures 3-3 to 3-7)*

Figures 3-3 to 3-7 show the effect of the strike with
various inclinations. It can be seen that with decreasing strike

3-4

angle the maximum of the total intensity curve becomes wider, it
moves toward the center of the body, and its minimum toward the
north decreases* In the s = 0° case there is no minimum (only a
maximum with all the inclinations) and the maximum is symmetrical
about the center of the body. Because of the type of the model,
at 1=0° and s = 0° there will be no anomaly observed. The type
of the model and the increasing importance of the horizontal com¬
ponent of the magnetic field toward low inclinations requires
that the strikes have a greater effect on the low inclination
curves.

d (Figures 3-8 to 3-12)
W (Figures 3-13 to 3-15)

3, Variable depth

4 »

Variable width

5* Variable vertical extent D (Figures 3-1& to 3-18)
These three parameters of the model were individually
changed, keeping all the other parameters constant, in order to
see how the shape of the total intensity curve respond to these
changes at various inclinations.

The variable depth curves show wider and smaller anoma¬
lies with increasing depth. The variable width curves cannot be
analyzed alone, because the shape of the curve depends on the W/d
ratio. It can be deduced from these curves that when the W/d ^
1, the horizontal extension (A) of the anomaly curve (the hori¬
zontal extension at maximum value, or % 0 % minimum value at
low inclinations) is quite insensitive to the width of the body.
When W/d ^ 5* the change of A almost entirely depends on the
width, and A approximately equals W.

3-5

The flattening of the maximum peak, or at low incli¬
nation of the minimum peak, is due to the effect of the lower
edge of the model; that is, the effect of the vertical extent
of the model D. On the curves near I = 45° there is a shallow¬
ing between the maximum and minimum due to this effect (see Pig*
3-17) which can clearly be seen at smaller values of D. This
flattening effect is evident when D < W. When the width is much
greater than the vertical extent, the maximums or minimums of the
anomaly curve will appear at the edges of the body, while over
the center the curve will be shallow, and nearly horizontal.

This effect can be seen also on the variable width curves at W =

15 km, when the width of the model is equal to its vertical extent.

It is evident that there is an altitude of observation
above which the effect of the lower edge cannot be observed.

6. Slanted bodies (Figures 3-19 to 3-22)

The present calculations are over four different models.
The tops are horizontal and symmetrical about the center of the
coordinate system, the sides are 30 ° and 60 ° from the vertical
toward the north, and on the other two models, toward the south.

The lower edges are horizontal also, and they are at 15 km depth.
The magnetic total intensity has been calculated for three incli¬
nations over these models, in order to show the change on these
curves due to the asymmetry of the models.

The general picture of the magnetic total intensity over
these asymmetrical bodies is the same as for symmetrical bodies:
there is a maximum to the south, and a minimum to the north of the

3-6

■

body (northern hemisphere), but now this maximum or minimum does
not depend entirely on the inclination, but depends upon the
shape of the body also. The various symmetries in the total in¬
tensity picture, which were mentioned at the group 1 calculations,
are only valid for symmetrical models. The total intensity
picture is different at all inclinations for an asymmetrical body.

Because the mass is more distributed, the amplitudes
are smaller; over the slopes of the model the anomaly becomes
wider, and over the sharp tapering edge near the observation line
(surface) the anomaly is sharp and pronounced. These effects can
be seen by comparing the two oppositely slanted models: Pig, 3-19*
Pig, 3-21, and Pig, 3-20, Fig. 3-22, All four figures show that
the shape of the upper part of the model has a dominating effect
on the anomaly curve, while the deep extension of the model only
causes the widening of the maximums or minimums over the sloping
side.

Except for the symmetry relations, the same conclusions
given in the previous sections can be drawn for these slanted
bodies at various strike angles s, and the d, W, and D changes
would also have similar effect as they had on the symmetrical

models.

3-7

3-9

3-10

3-11

3-13

3-14

3-1*

3-17

3-19

3-20

3-22

3-23

3-24

3-25

3-26

3-27

Chapter IV

A Method of Interpretation for Anomalies
of Total Magnetic Intensity caused by
Two-Dimensional Bodies

E.G. Zurflueh

A. General Remarks

There are numerous methods of interpretation for mag¬
netic data, each one having a certain range of application*

The method of Egyed (1948) which refers specifically to two-
dimensional bodies depends on the measurement of the hori¬
zontal and vertical components and can, therefore, not easily
be used for total intensity data.

Among the various existing procedures for the analysis
of total intensity anomalies, the most widely used are probably
the ones of Vacquier et.al. (1951)> Henderson and Zietz (1958)
and different mathematical methods as described by Peters (1949)*

The mathematical methods are complicated and can only
be applied to the interpretation of accurate magnetic maps.

The more graphical methods on the other hand have limitations
imposed by the geometrical configuration of bodies to be ana¬
lyzed* The method of Vacquier et.al., for instance, becomes
inaccurate when applied to bodies with a horizontal extent
equal or smaller than its depth of burial. The method of
Henderson and Zietz has been partly devised to overcome this
difficulty and gives good results for narrow bodies, while it
is, in turn, not applicable to wide bodies. Both methods re¬
late to prismatic or cylindrical bodies of more or less iso¬
metric horizontal cross section.

4-1

4-2

Since many of the important geological structures causing
magnetic anomalies are much more extended in one direction than
in others, the computations of Chapter III are expanded to
furnish a quantitative method of analysis for two-dimensional
bodies.

It is believed that the method can close a certain gap in
the range of graphical methods and it is also considered desirable
to work out a simple and uniform method which can be used for
bodies of small as well as large width compared to depth of burial.

The method has been devised for simplicity of procedure.

In order to minimize the influence of adjacent anomalies only
distances near the main extreme of the anomaly have been used for
the depth determinations. As a general rule it is only feasible
to use differences between abscissas of certain points on the
curve since all measures involving the total intensity values
are dependent on the susceptibility of the body. Prom various
distances that have been considered, the distance h between points
of half maximum (or minimum) value has been selected as the main
depth measure (see Pig. 4“2)* The computations have shown that
this distance (hereafter called half-width) depends mainly on
depth and width and is comparatively little influenced by other
factors.

The distance a has been chosen so that the ratio b = h/a
gives a clear distinction between different values of W/d. In
order to satisfy this condition two different distances a had to
be used for middle and for low and high magnetic latitudes. Por
inclinations near 0° or ± 90 ° there is a choice between both

4-3

values of a and depending on the circumstances one or the other
may give better results.

For the practical applications it is useful to keep thr
limitations of a method in mind. For this purpose the basic
assumptions that have been made are briefly summarized h^ e:

1* two-dimensional body.

2. vertical sides (this assumption implies that the
method is, in general, to be applied to anomalies
which are caused by susceptibility contrast and not
by topography).

3* vertical extent larger than depth of burial.

4. uniform magnetization.

5* no remanent magnetization (if there is a remanent

magnetization parallel to the present earth's field,
the method can be used but susceptibility estimates
will yield too high values).

B. Mathematical Expressions

The derivation of a few mathematical expressions for the
anomalies to be analyzed will be helpful in evaluating the
properties of the computed curves.

1. Formula for the anomaly curve

Assuming induced magnetization only, we obtain the following
expression for the total intensity anomaly produced by a two-
dimensional body with vertical sides and rectangular cross section
(see Fig. 4~l) by using equations (2-7)> (2-8), (2-9) and the
expressions on page 2-10:

oc/ >

(4-2)

at which

K =

and

X U F* C -L g~^s — X)

1 k F LX >^Wx ^

+ a

i

i

i

= depth of burial
= width

= vertical extent
= field point (on
the x-axis)

= origin (note that
the origin used in
this chapter is located
on the body*s axis of
symetry and not at the
field points as in
chapter II.)

section of assumed

body and coordinate system.

4-5

Using distances instead of angles and setting d = d + D for
simplification we obtain the following expressions for A and

B:

A =

B '

X -*•

sV

—— -J' — — .<4/rz4-U'

i I x j 1

A VO**, : A

o X~ v i

w- *+•

t7 A - £

/ w . ‘<V \ 1 ;

* ' V * I) J -

^ ^ ^ ^ J

— ' ->

(4-3)

r 'vc:,.

Equation (ip—ij.) can also be written as:

^ i

1 f . ^-j !

v x. yj

J

( 4 - 4 )

B

• }

+ (

x.«

4. x . ,x .. f T c i^ i- . / '?* V/M/ vyA\

* _** v ^ ^_f* - J > ' - \ A - ^ / V r ^ v p j -A- p. j

X , x. / r ^ .1 jj X'» ,,, i/' :i u-. . ,

W-;

£ y
( 4 - 4 *)

Equation (4-1) gives the resolution of the total intensity anomaly
into an even and an odd term as will be shown here:

The function A(x) represents the even term since replacing
x by -x in equation (4-3) does not change the expression. B(x)
can be written as B(x) = log f(x), where f(x) is the square root
of the fraction in equation (ip—4 a ) • If B is odd it must satisfy
the condition log f(x) = -log f(-x), and therefore, f(x) = l/f(-x).
That this condition is fulfilled can easily be seen from equation

( 4 - 4 *)•

2, Extremes positions

In order to evaluate the maxima and minima positions we
take the derivative of formula (4-1)

T» = K-jA* + K 2 B*

4-6

Taking the derivatives of A and B from equations (4-3)
and (4-4) we obtain:

d d ch d

A(v)= - z

/ W\< »x / •-c -a \ vv\<* jx ,

J +( v ^x) ^ + d + d +

(4-5)

BU) - »

Y v

y/

x*

. w

Y -v* —

_ _ 4 — [

ti W x

A +^ V T j

Y -

vY

V X

V/

* - r

( \.v

V / < t - \/4f .

d d tV < --/(4-6)

A_

In a more concise notation the derivative of the total intensity
anomaly can then be written as:

- ' X -I

. r / rr, iz' f tc^t-

i 1 <0 -

(4-7)

whe re

£>, & -.x s i' *

c 1 - 7 (^(d.D^-D.DJ),)
■=<* ?• , (* r(^r)(D 1 i) 4 i> J 'qD 1 Dj

-r 1

Expression (4-7) is of the form \ C<J *s *'

D^ J

rev )

p ^)

r ix)

where P^8)( x ) and p(5)( x ) are polynomials of the eighth and fifth
degree respectively.

Setting T* (x) = 0, p(8) ( x ) contributes the two solutions x = ± «jo
( where the anomaly curve approaches the x - axis asymptotically),
and in the following we will only consider the equation:

p(5)) x ) = K X E* + K 2 P* = 0

(4-8)

4-7

Substituting for D x through w© obtain the following ex¬
pressions for E* and F»:

= - 2WD/ •+ WD(4- <1 d (4-9)

+ -+3jw' 1 d-}w 4 d>

For the general case it is not possible to obtain an exact so¬
lution for equation (4-8)* However, since E* is an odd function
of the fifth degree and F* an even function of the fourth degree,
we can solve the problem for cases where either K-j_ or K 2 is equal
to zero*

a* Kg = 0•

According to equations (4-2) this condition is satisfied
for I = 0° or ± 90° and for s = 0°, Solving the equation E* (x) = 0
we get the following extremes positions:

This condition can also be written as tan I = ± sin s
(from (4-2). The relation can only be satisfied for + 45 < ^T?-45°«
Two extreme values are: I = 45 °» 3 = 90° and I = 0°, s = 0°.

4-8

For the latter case both K]_ and K 2 are equal to zero and no a-
nomaly will be observed.

From the equation F*(x) = 0 we obtain the extremes positions:

3• Summary

Formulas have been derived for anomaly curves over two-
dimensional, rectangular and vertical sided bodies and for the
positions of the extremes of the curves for special cases of
inclination and strike. It is not possible to obtain exact so¬
lutions for some other quantities like the half-width that are
used in the method of analysis presented here.

However, the formulas give information about the general
properties of the anomaly curves and they also allow conclusions
about the relationship between the different parameters involved.
Many of the relations mentioned here are illustrated in a quali¬
tative way by the curves shown in Chapter III of this report.

The following main conclusions can be pointed out:

The parameters k and F are merely multiplying factors and
have no influence on the shape of the anomaly curve (see equations
(4-1) and (4-2)).

The angular parameters I and s determine the relative im¬
portance of the even and odd terms of the curve. From equations
(4-2) it follows that for positive and negative inclinations of

4-9

equal magnitude the anomaly curves are symmetric with respect to
the z-axis, if the other parameters are held constant. This is
true if the origin of the coordinate system is taken over the
center of the body, and the profiles are drawn in the standard
way with the northern end to the right.

It can be shown that for s =90°, constant body parameters
and for inclinations that differ from 45° by equal and opposite
amounts, there is a central symmetry with respect to the origin.
Mathematically this can be formulated:

I = 45° ± T(x, oO = -T(-x, - oi) for s = 90 °

The two symmetry relations (see also page 2-9) reduce the
amount of computation necessary for the method of interpretation.
Since we can account for different strikes by applying a cor¬
rection to the values computed for s = 90 ° (see page 4-43)> we
can restrict ourselves to computing curves for inclinations be¬
tween 0° and + 4£°.

The first solution of equation (4-12) gives the distance
of the point with maximum value from the origin which coincides
with a zero position of the anomaly curve for the cases where
= 0, This is the distance a which is used for the interpre¬
tation at certain inclinations (see Fig. lj.-2). The formula
shows that for great width it is approximately equal to half the
body’s width and that, therefore, the main extremes are located
over the edges of the body.

With respect to the body parameters d, W and D it can be
said that, as test computations have shown, changes in D in

4-io

general affect mainly the amplitude of the anomaly, while they
have less influence on the location of the extremes and of the
zero positions. Furthermore, for physical reasons it is clear
that if D is greater than a certain value, say greater than d
or W, an increase in D will only slightly affect the anomaly
curve because the added magnetized mass is too distant from the
plane of measurement.

C. Description of Method

The quantitative analysis of total intensity anomalies de¬
scribed below is based on theoretical curves which have been com¬
puted for rectangular and vertical sided two-dimensional bodies
(see Fig. 4“l)* Sets of curves have been obtained for fixed ratios
of width to depth by varying deoth and width proportionally (Figs.
4-11 to 4-26). In addition the curves for variable strike of
Chapter III (Figs. 3-3 to 3-7) are used. The analysis can be
carried out in the following steps:

1) The inclination I of the region is determined by some
means (isogonic maps, tables or measurements)•

2) The approximate strike s of the body causing an a-
nomaly is determined from the magnetic or geological
map or from other data available.

If no information of this kind is at hand a rough
estimate of s can sometimes be obtained by comparing
the measured anomaly with the curves with variable
strike in Chapter III.

If a magnetic map is to be interpreted, a profile is

3)

4-ii

drawn through the center of an anomaly at right angles
to its strike.

If only profiles have been measured, the profile
has to be projected on the direction perpendicular to
the strike. For the practical application it is suf¬
ficient to multiply the distances h and a (see below)
by the sine of the angle between profile and strike.

4) The profile then is compared with the computed curves
and a curve that most closely resembles the actual one
is chosen, while keeping in mind the changes in shape
produced by the strike.

Since often the exact regional field is not known,
the main purpose of this comparison is to determine a
zero line for the anomaly.

5) From the measured anomaly profile the distances h and
a and the ratio b are determined from the main extreme
of the curve (see Fig. 4“2)« It is essential that the
part of the anomaly which lies between the main maximum
and main minimum be used for the determination of a.

h = x* -x . (half-width)

0.^ 0.5

a = x n - -x_ or x -x

0.5 0 m 0

according to the inclination
(see nomograms A),
b = h/a

Fig. 4-2: Total intensity anomaly and explanation

of h, a and b.

4-12

6) If the strike s differs substantially from 90° a cor¬
rection has to be applied to h and b by using tables
4-1 to 4“4 (page 4“£5)*

7) With the corrected values of h and b enter nomogram A
for the nearest inclination, (Figs, 4-28 to 4-37)*

Read a ratio of W/d for the body,

8) Nomogram A also gives an estimate of the depth of burial
d, A more accurate depth estimate can be obtained by
using nomogram B with h and the previously determined
value of W/d, The width W of the body can then be cal¬
culated from d and W/d,

9) For inclinations near the middle of the l£° intervals,
of the computed curves an interpolation is necessary.

For this purpose a W/d - ratio is determined from nomogram
A for the two adjacent inclinations. Then the average of
the two W/d - ratios is taken and a depth d determined
from nomogram B for both inclinations. The average of
these two values gives the depth,

10) The susceptibility k of the rock body can be calculated
from the maximum value of the anomaly. The model a-
nomalies have been computed for k = 10“^ emu and a total
field strength of £ 0*000 gammas. After having corrected
for the total intensity of the earth’s field in the area,
we obtain the susceptibility in electromagnetic CGS-units
by the following relationship:

4-13

k = 10*3 T max/T*max

where T max = maximum value of measured anomaly
T*max = maximum value of computed anomaly
It has to be borne in mind, however, that the vertical
extent D and the strike s have a strong influence on
the amplitude of the anomaly, and therefore, affect the
susceptibility measurements*

The interpretation procedure can be clarified by stating an ex¬
ample. The profile In Fig. 4-3 has been computed for an incli¬
nation of 26 ° and the true body parameters are indicated in the
figure*

The half-width h taken from the profile is 6*5 km, tne
distance a = Xj^-Xq is 4*2 km and the ratio b, therefore, equals
1*55* Correcting for the strike according to table 4-1 we add
l\% to the half-width and obtain h = 6.8 km. Similarly the cor¬
rected value of b is 1*44*

Entering nomogram A on page 4“£l with the adjusted values
of h and b, we estimate a W/d - ratio of 3*3* With this and
with h nomogram B on page 4-5>2 gives us the estimated depth of
the body which is approximately 2 km.

In order to complete the interpretation it is necessary
to determine the location of the body*s center. This can be
done in a more or less qualitative fashion as follows:

4—14

4-15

According to eq. (4-11) for I = 0° or ± 90° and for s
= 90 ° the anomaly is symmetric with respect to a vertical
line through the center of the body and the main extreme is
located over the center*

For the cases where the anomaly curve represents an
odd function (in the standard coordinate system, see eq* (4-
12)* The maximum and minimum are located approximately over
the edges of the body, while its center is given by the zero
position between them*

For other inclinations and strikes the center of the
body can be found between the main extreme and zero positions
by interpolation according to the computed curves*

D. Results of Practical Applications
1* Analyzed Cases

Some anomalies on Canadian aeromagnetic maps have been
interpreted with the method described above. From geological
maps it can be inferred that the chosen anomalies are caused
by rocks exposed on the surface* The depths determined,
therefore, can be directly compared with the flight altitudes*
Since no susceptibility values measured from rock samples
were immediately available, no susceptibilities have been cal¬
culated for the examples.

Fig. 4-4 4-^0 show parts of the aeromagnetic maps

with the analyzed anomalies. The chosen profiles are indicated

by the lines AB and the distances h and a are also shown in
the figures. The zero level used is given by the point 0
on the profiles. The inclination and the strike for each
case are indicated in parenthesis in the accompanying text.

a. Port Hawkesbury, N.S. (Pig. 4-4)

The ratio of W/d was found to be approximately 7 and
the depth determined is 290 m (I = 73s = 57°) • The flight
altitude is about 1000 feet or 305>m above ground. The geo¬
logical map shows an outcrop of Pre-Carboniferous granitic
rocks surrounded by Mississippian sediments at the location
of the anomaly.

b. Deskenatlata Lake South, NW Territories

For the anomaly near the Taltson River (Pig. 4"*5) a
W/d of 8 and a depth of 290m were obtained (I = 81°, s =

22 °).

The anomaly A^ in Fig. 4-6 yielded a W/d ratio of
4 and a depth of 225 m (I = 8l°, s = 13 °), the profile A 2 B 2
a W/d of 2.5 and a depth of 320m (I = 8l°, s = 10°).

The flight altitude for all three cases is 1000 feet
or 305m* The average of the three depths determined is 278 m
which differs less than 10$ from the average flight altitude.

Geologically the area is composed of Archaean and Pro¬
terozoic acid rocks.

Since 75° is the nearest inclination and because the

4-17

graphs for I = 75>° are more reliable for small depths than
the ones for I = 90 °, no interpolation was attempted and only
the graphs for I = 7£° were used for these cases.

c. Mack Township, Ontario (Pig. 4*7)

With a W/d ratio of 2 we estimate the depth to be 170m
(I = 76 °, s = 90 °)• The average flight altitude is 500 feet
or 15>3ia above ground. The rocks of the area are Huronian
quartzites with quartz diabase and quartz norite intrusions.

d. Gladstone Township, Ontario (Pig. 4“®)

The estimated depth is 130m and the W/d ratio 1 (I «s
76°, s = 55°)• The geological map indicates Huronian sedi¬
ments of the Gowganda Formation with intrusives of the same
type as found in Mack Township.

e. Bright Township, Ontario

The profile in Pig. 4~9 gives a W/d of 1.3 and a depth
of 120m (I = 7^°, s ~ 90°). The results for Fig. 4**10 are:
W/d = l/2 or smaller and d = 170m (I = 7&°, s = 90°) •

Pre-Huronian granite gneiss and related plutonic rocks
are the rock formations encountered near the anomalies.

The average value of all depth estimates for the ex¬
amined anomalies in Ontario (which are relatively close to¬
gether) is l48m* which comes very close to the average flight
altitude of 15 > 3 nu

4-18

2. Conclusions

The test cases have shown that the method gives reliable
depth estimates for anomalies that are reasonably elongated
in shape. F or carefully chosen anomalies the results often
are accurate to within 10$, although only very small depths
of burial have been considered. For greater depths it is to
be expected that the method will give even better results,
because deviations from the basic assumptions (uniform mag¬
netization etc.) will be less important for those cases.

For this method, like for others, it is of advantage
to choose isolated and simple anomalies. As a general cri¬
terion it can be said that good results can be expected, if
the shape of the anomaly profile is in its major parts com¬
parable to the theoretical curves.

For anomalies which are nearly isometric in the map
projection, the depth estimates obtained with this method
tend to become too small, because in this case the bodies
considered are not two-dimensional. However, if an anomaly
of more or less elliptic shape is about three times as long
as broad, the method can be applied successfully.

It was not possible to check the W/d ratios determined
for the actual anomalies considered above. But measurements
on computed curves have shown that they are somewhat less re¬
liable than the depth estimates.

4-19

Ctaminond Island

Morrison ■

Mclnnes Pt

Head Ba

Cove .

Ounphey
Head s'

McIntosh

Ballam

Dundee

Rear Black
/ River

Balmoral

Buchanan

FIG. 4-4= PORT HAWKESBURY, N S
MAP 238G

' *200

4—20

FIG. 4—5= DESKENATLATA LAKE SOUTH, NW TERRITORIES

MAP 107 G

4-22

FIG. 4-7' MACK TOWNSHIP, ONTARIO

4-23

FIG. 4-8= GLADSTONE TOWNSHIP, ONTARIO

2

4—24

FIG. 4-9 s

BRIGHT TOWNSHIP* ONTARIO

4-25

FIG. 4-10*

BRIGHT TOWNSHIP, ONTARIO

4-26

E. Curves and Graphs
1. Computed curves

The bodies for which the anomaly curves have been com¬
puted are of the kind described in paragraph B and the co¬
ordinate system used is essentially the same as shown in Pig
4-1• As abscissa the ratio of the horizontal distance x to
the depth of burial d for each curve has been plotted in
order to facilitate the graphical representation. The verti
cal extent D is the same for all bodies, namely 15 km. The
depth of burial d in km is indicated on each curve. In all
cases the assumed total intensity P of the undisturbed field
is 50*000 gammas and the susceptibility k is 10“3 emu.

For width to depth ratios smaller than one the curves

practically do not change their character except that their

amplitude is multiplied by a factor proportional to W/d.

That is the curves for W/d = l/2 e.g. are approximately the

same as the curves for w/d = 1 multiplied by the factor of

l/2. Therefore, only curves for W/d ^ 1 are shown here.

FIGURE 4-11

**8 1=0° S s 90°

d

4-27

FIGURE 4-12

~ » 4 1*0° S*90°

4-28

FIGURE 4-13

^ r 2 1=0° S s 90°

4-29

FIGURE 4-14

^ = I 1 = 0° S=90°

4-30

FIGURE 4-15

4-31

FIGURE 4-16

^=4 1=15° S=90°

4-32

FIGURE 4-17

4-33

FIGURE 4-18

4-34

FIGURE 4-19

^7*8 1 = 30° S=90°

a

4-35

FIGURE 4-20

4-36

FIGURE 4-21

1 = 30*

S=90*

4-37

FIGURE 4-22

= I 1 = 30° S=90°

4-38

FIGURE 4-23

4-39

FIGURE 4-24

4-40

FIGURE 4-25

4-41

FIGURE 4-26

w

d

1=45°

S=90°

4-42

2. Nomograms and Tables

The values of the tables and the nomograms have been
obtained by measuring the corresponding distances on the com¬
puted curves.

a. Nomograms

Nomogram A for each inclination shows two families of
curves. Both of them give the half-width h as a function of
the ratio b. For the solid curves the W/d ratio is held
constant, while the dashed curves are plotted for constant
depths of burial.

On nomogram B the half-width is plotted as a function
of depth for constant W/d ratios.

b. Tables

As equations (4-1) and (4-2) show, the strike influences
anomaly curves for different body parameters in a different
way. Corrections for the strike for different values of depth
and width are given in tables 4“1 to 4“4* ^ h and ^ b are

the percentages that have to be added to the measured values
of h and b in order to obtain the values which would have been
measured for a strike of 90 °*

For inclinations of less than 4£° and strikes near 0° the
positive part of the anomalies becomes larger than the negative
one, and in these cases it is more convenient to measure the
distances h and a from the maximum instead of using the minimum
(see figure 4“27)♦

4-44

4-45

h(km)

I

>b

0

2

4-46

h(km)

4-47

h(km)

A

9

e

7

6

5

A

3

2

fC

9

0

7

6

5

A

3

2

9

e

7

6

5

4

3

2

\i

4-48

m)

i=o° a + 90°

NOMOGRAM B

7 S 9 10

5 6 7 B 9 100

d (km)

FIG. 4-31

4-49

h (km)

4-50

h(km)

I = ± 15° & ± 75°

NOMOGRAM B

FIG. 4-37

4-51

h(km)

4-52

h (km)
a

| 5 - 30° a 1 60 '

NOMOGRAM B

FIG. 4-35

4-53

h(km)

4-54

NOMOGRAM

h(km)

I = ±45°

B

> d(km)

FIG. 4-37

m ^

k-H

Table 4“1

Correction for the Strike: W/d=3* d=2, D=l3

s

Ah(%)

Ab($)

Ah($)

AbW

i=i5°

6o°

+ 1.8

+ 2.8

45°

+ 4.1

+12.4

30°

+10.0

+26.9

15°

+14.4

+23.4

+11.4

+14.4

0°

- 4.2

+13.9

1=30°

6o°

+ 3.1

- 3.3

45°

+ 4.4

-10.9

- 1.1

-36.1

30°

- 1.3

-33.4

+ 4.4

-11.7

15°

- 3.4

+ 3.3

0°

-14 -7

- 7.2

i=43°

6o°

+ 3.3

+13.4

45°

0.0

+21.6

30°

- 4.0

+33.7

15°

-13.1

+33.4

0°

-l6.5

+19.7

i=6o®

6o°

- 2.7

+ 2.2

45°

- 3.8

+ 2.1

30°

- 8.3

+ 0.9

15°

-14.1

- 3.7

0°

-14.8

- 7.9

i=75°

45°

- 1.3

-10.4

15°

- 3.1

+ 2.9

0°

- 4.8

+19.7

4-56

Table 4"2

Correction for the Strike: W/d=5, d=0 # 5> D=l5

s

Ah(%)

Ab (%)

Ah(^)

AbW

1=15°

45°

+2,0

- 1.5

15°

+3.2

-12.7

+4.9

+ 2.8

0°

-1.2

+330

1=30°

45°

+4.2

-4.5

-6.5

-44*o

15°

-0.8

+10.7

0°

+5*4

+58.0

i=45°

45°

+6.4

+34-0

15°

-1.6

+57.0

0°

-3.1

+123

i=6o°

45°

+1.7

+12.6

15°

-3.2

+ 18.6

0°

-4.3

+60.5

i=75°

45°

-1.5

+ 2.0

15°

-1.9

+72.4

0°

-2.1

+172

4-57

Table 4“3

Correction for the Strike: W/d=l, d=2, D=l5

s

At W

i=*i5°

•6o°

- 2.6

45°

- 4.8

0

O

-12.6

15°

0°

-31.3

1=30°

0

0

sO

- 4-3

45“

-12.8

30°

-28.5

15°

o°

x=45°

60°

+ 7.4

45°

+16.2

O

O

+30.2

H

vn

0

+44*o

0®

+52.0

i=6o°

45°

+ 11.5

15°

+20.0

0®

+21.0

1=75*

45°

+ 1.9

15°

+ 3.0

0°

+ 5*o

-14.8

-31.0

-42.1

b(#)

-72.8

-30.0

-71.0

+ 7.2

+185

- 5.8

-19.3

-30.8

-37.9

-37.4

-14-3

-26.0

+ 8.8

+21.6

+20.1

+190

+15.0

+24*4

+55.0

+114

+306

+32.2

+99.4

+184

+ 24.5

+91.7

+154

4-58

Table 4*4

Correction

for the Strike

>: W/d=l,

d=0.5, D=i5

s

& h (%)

&b {%)

& h (%)

1=15°

45°

- $.0

-37.0

15°

-35.6

-81.7

-29.2

-75.8

0°

+ 4*6

+400

1=30°

45°

-14.0

-23.6

-34.2

-41.7

15°

+10.6

+19.6

0°

+25.0

+440

i=45°

45°

+19.2

+23.8

15°

+58*0

+188

0°

+60.7

+730

i=6o°

45°

+12.5

+23.6

15°

+22.7

+160

0°

+23.8

+396

i=75°

45°

+ 1.8

+39.4

15°

+ 5*6

+190

0°

+ 6.5

+360

Appendix

A Digital Computer Program
for the Computation of Magnetic Anomalies
over Two-Dimensional Polygonal Bodies

M. Talwani and J.R. Heirtzler

1* Framing the Problem

Rewrite eqs (2-7) and (2-8) for a prism of Fig. (2-ij.) expressing
<P in terms of distance:

x , - **. = x = - K,

J I “^X * "^11 ~ “ ^1/

5 Ls%r\ ^ IT i /

c. &

z/ */*• )

S =

IT.

6 cf > <- c>S < p > —

■^z/ * r.

*?, + *

= *>z. **/

t =

z/

* */\

/z.

’ / 2 .

A-l

A-2

This yields

V = Z (PI, o' - M t P')

(A-l)

H = t (n t p‘ *■ H t d')

with

(A-2)

*■ X

4 U

— —- JU. (r„r ,J

(A-3)

K/fc

•f- x

/v

sCo$ (r L /n )

(A-l*)

By refering back to page 2-7 , it is seen that integrations leading
to this result were carried out in a counter-clockwise direction,
that is B moving c.c.w. during the integration. For a computer
program there is an advantage in using a different labeling of
body points (see Fig. A-l):

The previous calculations can be utilized provided we replace
*9 ^7 *1 » by \r, , Q, by , and O t by €>, # These changes
give the following equations:

V ~ Z (A7 y Q - /t 1i p )

(A-5)

H - Z (PT A F -h <Q )

(A-6)

p =

it

*1,

i*t - ) +

*it *n.

toy (r u /tr { J

(A-7)

■3./

'it

O - (6>,-6>J - - *£*$(*. Jr,) (A-8)

**/t 4i/ + */t

For a given field point (point of observation or origin location)
the computer computes eqs (A-5), (A-6), and the total intensity
anomaly (eq 2-9) for each face of the body using eqs (A-7) and
(A-8) as internal subroutines# It sums the anomalies for each
consecutive face and then proceeds to do the same for other bodie
present# After the various sums are printed for a given field
point it goes to the next field point, etc#

The program has been written specifically for bodies

with induced magnetization only, but can also.be used without
modification for bodies having remanent magnetization or bodies
with both remanent and induced magnetization#

For computing convenience we can write (A-5> and (A-6)

above as:

V - fc<-oi x'/iU* s'JQ -

L (A-9)

u =r z-kFfc^i s'Jf + z';q.~\

L (A-10)

For bodies with induced magnetization only, k represents the
susceptibility and F the magnitude of the total undisturbed field,
I* =1# s’ = s as in Fig. 2-7* For bodies with remanent or mixed
magnetization kF has, of course, no significance; it is merely
retained for convenience so that one can use the same program as
for bodies with induced magnetization. In practice kF is made
to represent the magnitude of M rem for bodies with remanent mag¬
netization by giving k the dummy value of 1 and F the value for
M rem , Also I 1 = a, s = b as mentioned on p. 2-11, Similarly for
a body with mixed magnetization k is again given the dummy value
of 1 and F the value of M tot . I* = o**, s* =p as mentioned on
p, 2-12,

We note that the expression (2-9) for the evaluation of
T, the total anomaly, remains unchanged whether the magnetization
is induced, remanent, or mixed,

2, Flow Diagram

To assist in the understanding of the flow diagram and
the subsequent program the following list of symbols is defined:

A-5

Table A-l
Program Symbols

Symbol

DIP*

DIED*

D»

DD*

F*

C*

EXX, ZEE*

JTOT*

CONS*

FO*

DELX*

THETA, THETB

Meaning

Inclination or Dip, I, in degrees

Modified Inclination or Dip, I*
in degrees

Strike, s, in degrees

Modified Strike, s', in degrees

Total Magnetic Field intensity,
in gammas

Magnetic susceptibility, k, in e*m.u.

Coordinates of body comers, entered
in program in clockwise sequence in¬
cluding first twice

Total number of body points including
first twice

Constant value of z of field points,
zero for shipboard surveys, negative
for aeromagnetic surveys taking sea
level as z = 0

x value of first field point
Increment In x value of field points
Angle to corner measured from x-axis.

OMEGA

K

KTOT*

LNO*

Angular span of face from field
point,

Subscript for consecutive field
points

Total number of field points
Subscript for body number

LNOT*

Total number of bodies

H

Horizontal anomaly for single face
of body at field point, in gammas

V

Vertical anomaly for single face of
body at field point, in gammas

PSUM

The sum of P*s for the various
faces of the body

Q3 UM

The sum of Q*s for the various
faces of the body

HSUM, VSUM

Horizontal and Vertical anomalies,
in gammas, at field point for
given body

HASUM, VASUM

Horizontal and vertical anomalies,
in gammas, at field point for all
bodies

TSUM

Anomaly in total intensity, in
gammas, at field point for given body

TASUM

Anomaly in total intensity, in gammas
at field point for all bodies

•a-Data which must b© entered into computer from typewriter.
D and DIP are used in the calculation of V and H*

DD and DIPD are used to calculate TSUM and TASUM.

A-7

A-8

C FORTRAN PROGRAM FOR IBM 1620

DIMENSION FX(47) ,FZ (ij.7) ,VASUM(47) ,HASUM(47)
DIMENSION VSUM(47),HSUIvl(47) ,EXX(30),ZEE(30)
ij.00 ACCEPT,DD,DIPD,KTOT.LNOT,CONS

402 sdd=sin(. 0174533 *dd)

CDIPD=COS(.0174§33*DIFD)

SDIPD=SIN(.0174§33*DIPD)

601 ACCEPT,FO,DELX

603 DO 604 K=1,KT0T
RK=K

FX(K)=(FO-DELX)+DELX*RK

FZ(K)=CONS

VASUM(K)=0.

604 hasum(k)=o.

IO 5 ACCEPT,LNO,C,F,JTOT
ACCEPT,D,DIP

410 sd=sin(.0174533*d)

CDIP=COS (.0174533*DII>)

SDIP=SIN(.0174533*®IP)

DO 11 J=l,JTOT

411 ACCEPT,EXX(J),ZEE(J)

11 CONTINUE

DO 36 K=1,KT0T
PSUM=0.

QSUM=0.

X1=EXX(1)-FX (K)

Z1=ZEE(1)-FZ(K)

RSQI=0CL**2+ZHh*2

if (xi) no, 140,180

110 IF(Z1)120,130,130

120 THETA=ATN(Z 1 /X 1 )- 3.1415927

GO TO 200

130 THETA=ATN(Zl/Xl)+3.1415927
go to aoo

140 IF(Zl)l 50 ,l 60,170
150 THBTA=-1.5707963
GO TO 200
160 THETA=0•0
GO TO 200

170 THETApI.5707963
GO TO 200

180 THETA=*ATN (Zl/Xl)

200 J=2

201 X2=EXX(J)-FX(K)

Z2=ZEE(J)-FZ(K)

RS Q2»X2*-*2+Z2*~*2
IF(X2)210,240,280
210 IF(Z2)220,230,230
220 THETB^ATN(Z2/X2)-3.1415927
GO TO 300

230 THETB=ATN(Z2A2)+3.1415927
GO TO 300

A-9

2k0 IF(Z2)250,260,270

250 thetb=-i .5707963
GO TO 300
260 THETB= 0,0

GO TO 300

270 THETB=1 # 5707 963
GO TO 300

280 THETB=^TN(Z2/X2)

300 IP(Z1-Z2)320,31,320
31 P=0.

Q=0.

GO TO 32

320 OMEGA=THETA-THETB

IP(OMEGA)3201,3202,3202
32021P(OMEGA-3.1415927 ) 330,330,340
3201IF(OMEGA+3•1415927)34°* 330, 330
330 THETD=OMEGA
GO TO 370

340 IF(OMEGA)350,360,360
350 THETD=0MEGA+6.2831853
GO TO 370

360 thetd=omega- 6,2831853

370 XL2=X1-X2
Z21=Z2-Z1

XSQ=X12**2

ZSQfZ21**2

XZ=Z21*X12

GL=0 ,5'^LOG (RSQ2/RSQ1)

P=( (ZSQ/(XSQ+ZSQ) )*THETD)+((XZ/(XSQ+ZSQ) )**GL)
Q= (THETEHt( XZ/ (XSQ+ZSQ)))-(GL*(ZSQ/(XSQ+ZSQ)) )

32 IP(SENSESWITCH 3)33#34

33 H=2.*C-*F#( (CDIP*SD*P) + (SDIP**Q))

V=2.*CttF*( (CD I P-ttSD*«-Q) - (SDI P-«-P))

PRINT, XI, X2,Z1,Z2,XSQ,ZSQ,XZ
PRINT , RSQ1, RSQ2, THETA , THETB
PRINT,K, J, THETD,GL, P, Q,H,V

34 PSUM=PSUM+P
QSUM=QSUM+Q
X1=X2

Z1=Z2

RSQ1=RSQ2

THETA=THETB

J=J+1

JR=J-1

IP(JR-JTOT)201,202,202

202 HSUM(K)=&,*CttF*( (CDIP*SD*PSUM)+(SDIP*QSUM))
VSUM(K)=2.*C*F*( (CDIP#SD*-QSUM) -(SDIP-&PSUM))
HASUM(K)=HSUM(K)+HASUM(K)

VASUM (K) =VSUM(K) +VASUM (K)

36 CONTINUE

PRINT,LNO,C,F,D,DIP
DO 37 K=1,KT0T

A-10

TS UM= (HSUM (K) *CD I PD*SDD ) + (VSUM(K)-b-SDIPD)

37 PRINT,K,FX(K),HSUM(K) ,VSUM(K) ,TSUM
IF(LNO-LNOT )105,38,105
38 K=1

3 81 TASUlfc (HASUM (K) *CDI PD*SDD) + ( V AS UM (K) *SD I PD)
39 PRINT,K,FX(K),HASUM(K),VASUM(K),TASUM
K=K+1
KR=K-1

IF(KR-KTOT) 381,391,391
391 GO TO 400

END

A-ll

3* A Simple Example

Suppose we wish to determine the induced anomalies of
three field points for an assumed body shape as given in Pig* A-2.

Pig. A-2

Suppose the axis of the body makes an angle of 60° with magnetic
north, the inclination is 10 °, the ambient magnetic field strength
is 50,000 gammas, the susceptibility is 0.001 emu.

We have JTOT = 4, LNO = 1, LNOT = 1, CONS = 5., KTOT = 3,

P0 = 5., DELX = 3.

Information in the following form must be supplied to
the machine by typewriter

60. 10. 3 15.

5. 3.

1 .001

60. 10.

8 . 8 .

8 . 6 .

11 . 6 .

8 . 8 .

50000 . 4

A-12

With senseswitch 3 off numerical results will be printed as
follows:

(LNO)

(c)

(p)

(D)

(DIP)

1

1.0000000E-03

50000.000

60.000000

10.000000

(K)

FX(K)

HSUM(K)

VSUM(K)

TSUM(K)

1

5.0000000

11.142248

8.2548^30

10.936315

2

8.0000000

-25.503102

53.009137

-12.545855

3

11.000000

-12.179321

-43.110189

-17.873370

A-X3

4# Console Operation for Lamont IBM 1620

Switches:

Clear Memory:

Load Tape:

Load Data:

Senswitch 3 ON for detailed print out (used
for trouble shooting only)

Other senseswitches OPP
Program switches OPF
Other switches ON

Reset, Insert, 160001000000, Release, Start,
(after a few seconds) Instant Stop,

Mount object tape, set tape console switch
to Reel, Reset, Insert, 3^0000000300, Re¬
lease, Start, If machine fails to read, re¬
peat previous instructions. If functioning
properly machine will pause 2/3 way through
reading to think but will proceed to end of
tape when it will type "Load Data”,

Set typewriter tabs at approximately

16 28 45 58 70.

Reset, Insert, 4907500* Release, Start
Numbers must be typed in upper case.

Minus signs and decimal points in lower case.
Zero in either case.

A space must be typed between successive data
within a line of data,

A record mark after the last datum of a line
is optional.

A decimal point is required after floating-

Typing Error:

Stop Output:

Check Stops:

point data.

A decimal point is optional after fixed-
point data.

If incorrect data has been entered in machine
Reset, Insert, 4907500# Release, Start, re¬
enter all data.

If incorrect data has been typed but not
entered into machine Senswitch 4 ON, Release,
Start, Senswitch 4 OFF, repeat line of data.
Depress Stop. Neatness can be achieved by
depressing Stop while carriage of typewriter
is in return motion following typing of last
line of desired data. If machine is not
stopped it will automatically return to the
stage where it is again ready to accept all
data.

If a check stop occurs during data entry it
may be due to a parity error in that part of
the memory related to data storage: Reset,
Insert, 4907500, Release, Start, re-enter
all data.

If a check stop occurs at any other time or
if the above procedure fails clear memory,
re-load tape, re-enter all data.

R-l

ACKNOWLEDGEMENTS

The figures of Chapters III and IV of this report
were drawn by Jay Maxon and David Wolfe,

The aeromagnetic maps of Chapter IV-D were reproduced
with the permission of the Geological Survey of Canada,

Financial support for this work was provided by the
Office of Naval Research, Contract Nonr 266(48) and the National
Science Foundation, Grant 4l5l*

Reproduction of this document in whole, or in part,
is permitted for any purpose of the United States Government,

REFERENCES

Baranov, V, (1953)* Sur 1© calcul de 1*influence gravimeurique
des structures definies par les isobathes, Geoph,
Prosp., v. 1, no, 1, pp. 36 - 43 *

Baranov, V. (1955)s A new method of interpretation of aeromag-
netic maps; pseudo-gravimetric anomalies. Geophysics
v, 20, no. 2, pp. 379-380,

Egyed, L. (1948): The determination of an infinite inclined

dike from the results of gravity and magnetic surveys.
Geophysics v. 13, no. 3> PP* 437-442*

Gassmann, F. (1951): Graphical evaluation of the anomalies of
gravity and of the magnetic field, caused by three-

R-2

dimensional bodies, Proc. 3rd World Petrol, Congress,
Sec. 1, pp. 613-621,

Girdler, R.W,, and G. Peter (i 960 ): An example of the importance
of natural remanent magnetization in the interpretation
of magnetic anomalies, Geophys. Prosp, v. 8 , no. 3,

pp. 474-4 8 3.

Gulatee, B.L. (1938): Magnetic anomalies. Survey of India,

Prof, paper no. 29> 12 p., 11 plates.

Haalck, H. (1927): Die magnetischen Verfahren der angewandten
Geophysik. Bomtraeger, Berlin.

Heiland, C.A. (1946): Geophysical exploration. Prentice-
Hall Inc., New York.

Heirtzler, J.R., C.L. Drake and J. Hirshman ( 1962 ): Magnetic
anomalies off the east coast of North America. In
preparation.

Henderson, R.G. and I. Zietz (1946): Analysis of total mag¬
netic intensity anomalies produced by point and line
sources. Geophysics, v. 13# no. J>, pp. 426-436.

Henderson, R.G. and I. Zietz (1949&): The computation of

second vertical derivatives of geomagnetic fields.
Geophysics v. l4> no. 4> PP* 508-516.

Henderson, R.G, and I. Zietz (1949^>) : The upward continuation
of anomalies in total magnetic intensity fields.
Geophysics v. l4* no. 4* PP* 517-534*

Henderson, R.G. and I. Zietz (1957 )• Graphical calculation of

total-intensity anomalies of three-dimensional bodies
Geophysics v. 22, no. 4> PP* 887-904*

Henderson, R.G. and I. Zietz (1958)* Magnetic-doublet theory
in the analysis of total-intensity anomalies, Geol.
Surv. Bull, 1052-D, Washington,

Henderson, R.G. (i 960 ): A comprehensive system of automatic
computation in magnetic and gravity interpretation.
Geophysics v, 25, no, 3> pp. 569-585*

Hughes, D.S, and W.L. Pondrom (1947)? Computation of vertical
magnetic anomalies from total magnetic field measure
ments, Trans, Amer. Geophys, Union, v, 28, no. 2,
pp. 193-197.

Kohler, K. (1958): Grundlagen fur de Auswertung von magne-
tischen Anomalien. Preiberger Forschgshefte C4l,
Akademie-Verl., Berlin, pp. 1-128.

Nettleton, L. (194^): Gravity and magnetic calculations.
Geophysics v. 7# no. 3* PP* 293-310.

Nettleton, L. (1954) : Regionals, residuals and structures.
Geophysics v. 19# no. 1, pp. 1-22.

Peters, L.J. (1949): The direct approach to magnetic interpre
tation and its practical application. Geophysics
v. 14# no. 3, pp. 290-320.

Skeels, D.C. (1947): Ambiguity in gravity interpretation.
Geophysics v. 12, no. 1, pp. 43-58.

Skeels, D.C. and R.J. Watson (1949): Derivation of magnetic

and gravitational quantities by surface integration.
Geophysics v. l4> no. 2, pp. 133-150.

Smellie, D.W. (1956): Elementary approximations in aeromagnetic

interpretation. Geophysics v. 21, no. 4» pp. 1021-1040.
Smith, R.A. (1959 )• Some depth formulae for local magnetic
and gravity anomalies. Geop. Prosp., v. 7# no. 1 ,
pp. 55-63.

Smith, R.A. ( 1961 ): Some theorems concerning local magnetic

anomalies. Geoph. Prosp., v. 9# no. 3# PP* 399-410*
Talwani, M., and M. Ewing (i 960 ): Rapid computation of gravi¬
tational attraction of three-dimensional bodies of
arbitrary shape. Geophysics v. 25# no. 1, pp. 203-225.
Vacquier, V., N.C. Steenland, R.G. Henderson and I. Zietz (1951 )•
Interpretation of aeromagnetic maps. GSA Memoir 4?>

151 P.

Vacquier, V., A.D. Raff and R.E. Warren ( 1961 ): Horizontal

displacements in the floor of the northeastern Pacific
Ocean. GSA Bulletin v. 72, no. 8 , pp. 1251-1258.

Vajk, R. (1951 )• Regional correction of gravity data. Geo-

fisica Pura e Applicata, v. 19# fasc. 3-4* pp. 129-143*
Vestine, E.H. (1941) : On the analysis of surface magnetic

fields by integrals, part I. Terr. Mag. v. 46# no. 1,

pp. 27 -4i.

Vestine, E.H. and N. Davids (1945): Analysis and interpretation
of geomagnetic anomalies. Terr. Mag. v. 50# no. 1,

pp. I- 36 .

Zietz, I. and R.G. Henderson (1954) : Total-intensity magnetic
anomalies of three-dimensional distribution by means

of experimentally derived double layer model fields.
Science, Washington, v. 119 , no. 3088 , pp. 329-330.

List of maps used

Geological Survey of Canada (1922): Portions of the districts

of Algoma, Sudbury and Timiskaming, Ontario, map 155A.
Geological Survey of Canada (1945): Geological map of the
Dominion of Canada, map 820A.

Geological Survey of Canada (1949 )• Geological map of the
Maritime Provinces, map 910 A.

Geological Survey of Canada (1955): Port Hawkesbury, N.S.

P

Aeromagnetic series, sheet 11 map 238 G.

Geological Survey of Canada (195&): Deskenatlata Lake South,

NW Territories. Aeromagnetic series, sheet 85

15 ’

map 107G.

Province of Ontario, Department of Mines (1955): Aeromagnetic
map. Bright Township, Distr. of Algoma.

Province of Ontario, Department of Mines (1955): Aeromagnetic
map, Gladstone Township, Distr. of Algoma.

Province of Ontario, Department of Mines (1955): Aeromagnetic

map. Mack Township, Distr. of Algoma