CA2265063A1 - Converting magnetic field depth solutions to a three-dimensional seismic data format - Google Patents

Abstract

Magnetic field depth solutions are calculated for a magnetic field survey. These depth solutions are defined spatially by x, y and z co-ordinates. A digital three dimensional zero matrix is constructed by a computer program to represent a prism of earth underlying the area of the magnetic field survey. The matrix is subdivided evenly into elements representing smaller bins or prisms. Unit values are assigned to those elements to which the depth solutions correspond spatially. The resulting matrix values therefore represent a density or frequency of occurrence per unit volume of the magnetic depth solution points. The matrix elements are then spread to surrounding elements by three-dimensional convolution methods using Gaussian wavelets in each of three directions, thereby smoothing the data. The resulting data volume, or three dimensional digital matrix is then converted to a format compatible with a computerized seismic data visualization and interpretation tool and loaded into the tool. The integrated analysis of aeromagnetic depth solution data with seismic data is thereby enabled.

Description

10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 FIELD OF THE INVENTION This invention is in the field of geological exploration and is more particularly directed to identifying sub-surface structures through the integrated analysis of magnetic field depth solutions and seismic data. BACKGROUND OF THE INVENTION In the field of geological exploration, particularly as applied in the oil and gas industry, various methods are commonly used to acquire information regarding the depth, size and orientation of geological structures and formations. Seismic methods involve the induction of acoustic energy into the earth. Reflected or refracted acoustic energy returns to the surface and is recorded. By analysis of the travel times of the energy along reflected or refracted paths, the location and character of sub-surface structures can be determined. Furthermore, three-dimensional visualisation and interpretation of this seismic data is made possible by the application of computerized seismic data interpretation tools. Another common technique of geological exploration is the analysis of measurements of the earth’s magnetic field. As is commonly known in the art, magnetization of the earth’s core and crust, as well as the interaction of the earth’s magnetic field with solar phenomena, all contribute to the magnetic field measurements of a magnetic survey. It is also well known that the non- homogeneity of the earth that is due to the geological structure of the sub- surface, is evident in the deviations of the earth’s magnetic field from the magnitudes expected assuming uniformity. Furthermore, the magnetic field of the earth’s core is commonly represented in the art by the mathematical 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 formula, the international Geomagnetic Reference Field (IGRF), that is agreed upon by the International Union of Geodesy and Geophysics (IUGG). As is well known in the art, the magnetic anomalies of the earth’s geomagnetic field can be calculated by subtracting the IGRF from the measured magnetic field and by applying corrections for diurnal variations in the field. Figure 1a illustrates a simplified cross-section of the earth in which a clipping magnetized body lies at a depth d from ground surface. Figure 1b illustrates how the presence of the aforementioned anomalous body results in a deviation of the measured magnetic field from the field expected assuming sub-surface homogeneity, the IGRF. The form of this deviation is dependent upon the shape of the body, its orientation with respect to the earth’s magnetic field, its position on the earth, the magnetic susceptibility of the body and its surroundings, and the depth of the body from ground surface. As is well known in the art, the frequencies (in the spatial domain) at which the magnetic field varies are dependent upon the sub-surface depth of the anomalous geological structures. Shallow structures produce higher frequency (in the spatial domain) variations in the magnetic field of the earth than do deeper structures. As such, the analysis of the geomagnetic field is useful in determining estimates of depth to geological structures. As is common to the art of aeromagnetic surveying, discrete measurements of the earth’s magnetic field are made across a survey area. An airplane carrying a magnetometer flies in an orderly pattern over the area and measurements of the earth’s magnetic field are taken at regular time intervals. It is common to design aeromagnetic surveys along regularly spaced, parallel flight lines across the survey area, with consideration of the 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 geology and topography specific to the area. Figure 2 shows a survey area with aeromagnetic flight lines crossing it in an orderly fashion, with arrows indicating the direction of flight during the survey. The measurements and the locations of the measurements are recorded and are later converted to grids of measurements. As is common to the art, filters may be applied to the grids of magnetic measurements in the space or frequency domain, to enhance magnetic anomalies of particular wavelengths, and maps of the filtered magnetic field are constructed to enable more insightful analysis of the earth’s magnetic field. As is -common to the art, magnetic field depth solutions may be calculated using known formulas. Calculation of magnetic depth solution points provide a means of estimating depths to anomalous magnetized bodies, by analyzing either the magnetic field profiles gathered along the flight lines or by analysis of the grid of the magnetic field over the entire survey area. Figure 1a shows a magnetic depth solution point that has been determined using a method of analyzing the geomagnetic field. The spatial coordinates of the calculated magnetic depth solution points are stored using a computer disk drive or magnetic tape. However, the prior art does not provide any fast and simple method of representing the depth solutions in a plurality of perspectives or orientations. In the prior art, interpretation of the magnetic field depth solutions has been a slow process, restricted to one or two profile directions. Furthermore, the prior art does not provide a fast and simple method to integrate magnetic field depth solutions with seismic data. 10 11 12 13 14 15 16 17 18 19 20 21 22 23 CA 02265063 l999-03- 10 It is an object of the present invention to provide a method of three- dimensional visualisation and interpretation for magnetic field depth solutions. it is another object of the present invention to provide a method for integrating magnetic field depth solution data with seismic data so they can be used simultaneously in seismic interpretation tools. Other objects and advantages of the present invention will be apparent to those of ordinary skill in the art, having reference to the following specification together with its drawings. __yMMARY QF TH§ INV§NTl0_N In accordance with the invention, magnetic depth solution points are converted to a digital format which can be integrated with a seismic visualization and interpretation tool for display. This is accomplished by practising the following steps: o storing the calculated co-ordinates of magnetic depth solution points, derived from a magnetic survey conducted over a survey area, in computer memory; o establishing in computer memory a matrix, representative of a pre- determined volume underlying the survey area and which contains a distribution of the points, said matrix being sub—divided into elements; a retrieving the co-ordinates for each point, identifying the matrix element with which it corresponds spatially and assigning a value to that element; 10 11 12 13 14 15 16 17 18 19 20 21 22 CA 02265063 l999-03- 10 o combining the assigned values for each element to create a product matrix representing the distribution of the points in the pre- determined volume; and o converting the product matrix into seismic data format and displaying it. DESCRIPTION OF THE DRAWINGS Figure 1(a) is a simplified cross-sectional view of the earth, showing an anomalous magnetized body buried at a depth d from the surface; Figure 1(b) is a comparative profile of the calculated IGRF and the magnetic field as measured across the buried magnetized body of Figure 1(a); Figure 2 is a schematic plan view illustrating a typical aeromagnetic survey showing the lines flown by the survey plane traversing the area; Figure 3 is a block diagram setting forth a preferred form of the steps of the invention; Figure 4 is a schematic perspective view illustrating a survey prism which represents the volume of earth beneath an aeromagnetic survey grid area. Smaller component prisms or bins sub-divide the greater survey prism. Flight lines are shown across the top of the survey prism; Figure 5 is a schematic perspective view of a three-dimensional zero matrix representing the prism of Figure 4; 10 11 12 13 14 15 16 17 18 19 20 21 CA 02265063 l999-03- 10 Figure 6 is a schematic perspective view of the matrix of Figure 5 comprising the inclusion of a unit value in an elemental matrix element to represent the presence of a magnetic depth solution point within the corresponding sub-volume or bin; Figure 7 illustrates a simulated example of a wavelet constructed for convolution with the digital matrix of Figure 6. The digital samples representing the wavelet in the computer convolution process are illustrated here; Figure 8 is a schematic perspective view of the matrix of Figure 6 illustrating the convolution of the matrix with a wavelet in the vertical direction; Figure 9 is a schematic perspective view of the matrix of Figure 6 illustrating the convolution of the matrix with a wavelet in the lengthwise direction; Figure 10 is a schematic perspective view of the matrix of Figure 6 illustrating the convolution of the matrix with a wavelet in the width direction; Figure 11 illustrates an example of a seismic cross-section, with an interpreted fault clearly denoted; and Figure 12 illustrates an example of a cross-section view of magnetic depth solution data, converted to SEG-Y format and displayed using a seismic interpretation tool, for comparison with the seismic data (of Figure 12) from along the same profile. 10 11 12 13 14 15 16 17 18 19 20 21 22 23 CA 02265063 l999-03- 10 QESCRIPTKEF THE PRlfil_E_RRED EMBODIMENT As will be apparent from the following description, the preferred embodiment of the present invention may be applied to any data comprised of a scattered distribution of points in three dimensional space. A particularly useful application of the preferred embodiment of the invention has been the analysis of magnetic field depth solutions. Magnetic field depth solutions are computed by analysis of geomagnetic measurement data. Geomagnetic measurement data are obtained by an aeromagnetic survey carried out in a conventional fashion. More particularly, geomagnetic measurements are recorded by a magnetometer carried in an airplane. The airplane is flown along parallel flight lines 1 over a predetermined survey area 2, as illustrated in Figure 2. A variety of algorithms such as the Werner and Euler methods are commonly used in the art, for computing magnetic depth solutions along the flight line magnetic field profiles. Although these methods are inherently different, each method computes series of co-ordinates (X, Y and Z), which are stored on a computer disk drive. These co-ordinates represent locations of anomalous magnetic sources. The method according to the preferred embodiment of the invention for analyzing magnetic depth solution data will now be described. A computer program was written for converting these magnetic depth solutions to 3D seismic data format, for visualizing and interpreting using a seismic workstation. According to the preferred embodiment of the invention, the processes of the computer program are illustrated in the flow chart of figure 3. 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 In the first step 3 of the preferred embodiment, the geometry of the survey prism is defined by the computer program. The survey prism is a volume underlying the survey area which incorporates the sedimentary and basement rock down to a depth of several kilometers. An example survey prism 4 is illustrated in figure 4. To define the geometry of the prism, the spatial co-ordinates of the corner points 5 of the survey area are input into the computer. Furthermore, the total depth 6 of the prism is entered into the computer. This rectangular prism may be constructed with any spatial orientation. However, if the flight lines 1 of the aeromagnetic survey are parallel, it is preferred that the prism be aligned with these flight lines. The survey prism is sub—divided by the computer program into a plurality of smaller rectangular prisms of equal dimensions, herein referred to as ‘bins’ 7. The bins can be of any dimensions. However, according to the preferred embodiment of the invention, the bins are defined by the computer program to have length 8 and width 9 equal to half the flight line spacing 10. According to the first step 3 of the preferred embodiment, the computer program defines the geometry of the survey prism by calculating the spatial co-ordinates of the center of each 'bin' 7 of the survey prism. The co-ordinates for the centers of all the bins are stored by the computer program in computer memory. After the bin co-ordinates are stored in the computer memory, the second step 12 of the preferred embodiment is performed. In the second step, the computer program generates a three—dimensiona| digital matrix in the computer memory. Figure 5 illustrates a three-dimensional zero—matrix 13 that was generated by the computer program to represent the survey prism 4 of 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 figure 4. The three-dimensional digital matrix is composed of a number of discrete data values. We shall refer to these data values as 'elements' 14. The computer program creates the matrix with a number of elements in length, width and depth dimensions, equal to the number of bins 7 in the aforementioned survey prism 4. According to the preferred embodiment, the computer program assigns to each matrix element a null value in this step. Each matrix element 14 will herein be identified as shown in figure 5, by trace number 15, line number 16, and depth increment number 17. Following the creation of the digital matrix 13, the third step 19 of the preferred embodiment is performed, in which the computer program retrieves the co-ordinates of the magnetic depth solutions from computer disk storage. According to the preferred embodiment, only the spatial co-ordinates of one type of magnetic depth solution are retrieved. When the desired magnetic depth solution co-ordinates have been retrieved from computer disk storage, the fourth step 20 of the preferred embodiment is performed, in which values are assigned to the matrix elements. It is contemplated that one of ordinary skill in the art having reference to this specification will be able to readily apply the present method. According to the preferred embodiment of the invention, the computer program compares each magnetic depth solution co-ordinate (retrieved in the third step 19) with the center co-ordinates of each bin (defined in the first step 3) to determine within which bin each magnetic depth solution is contained. The elements of the digital zero-matrix (constructed in the second step 12) which represent the aforementioned containing bins 7 are identified by the computer program, and a value is added to each of these matrix elements 14. 10 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 According to the preferred embodiment of the invention, a discrete value of one is added to the digital matrix element 14 to represent each magnetic depth solution. Figure 6 illustrates how a unit value 21 is added to a three dimensional matrix element 14 which represents a bin containing a magnetic depth solution point. If a plurality of magnetic depth solution points are found to be contained within a single bin, then the matrix element corresponding to that bin will equal the number of depth solutions contained therein. When a value has been assigned to the digital matrix for each of the magnetic depth solution points of the survey, the matrix values are a scalar distribution representing the density, or rather, the frequency of occurrence per unit volume of the magnetic depth solution points. It is useful for purposes of visualization, for the digital matrix elements carrying values representing magnetic depth solution points, to be surrounded by a zone of non-zero elements. This process will hereon be referred to as "spreading". The shape and volume of this zone of non-zero matrix elements is dependent upon the anticipated resolution required. According to the fourth step 20 of the preferred embodiment of the invention, "spreading" is performed by convolving the digital matrix with a wavelet in three orthogonal directions. The fifth step 22 of the preferred embodiment of the invention involves the construction of the aforementioned wavelet. Figure 7 illustrates a symmetrical wavelet 23 which has a central amplitude of one. Since the preferred embodiment of the invention is performed in digital form, the wavelet is stored in the computer memory as a set of discrete samples 24 of the wavelet's amplitude taken at equally spaced intervals 25 along the wavelet's 11 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 length. Any shaped wavelet may be used. However, it is preferable that the wavelet be all positive, and that the wavelet be of a dominant spatial frequency sufficient to provide the anticipated resolution required. Preferably, a normal distribution curve, or some similar pattern, is chosen for this purpose. Furthermore, the sampling interval 25 chosen for digital representation of the wavelet is dependent on the anticipated resolution required. After the digital wavelet 23 has been stored in computer memory, the sixth step 26 of the preferred embodiment is performed in which the computer program convolves the digital wavelet 23 with the digital matrix in three orthogonal directions. The convolution process is well known in the art of digital signal processing and in the art of seismic data analysis. It is contemplated that one of ordinary skill in the prior art will be able to readily apply the present method. For each of the orthogonal directions, the wavelet 23 is convolved with the matrix along a plurality of parallel paths, so that every element of the matrix is effectively convolved by the wavelet in three directions. The directions of convolution may be chosen arbitrarily, but the preferred directions of convolution are in the vertical direction, along the length of matrix and across the width of the matrix. Figure 8 illustrates the convolution of the three dimensional matrix with a wavelet in the vertical direction 27. According to the preferred embodiment of the invention, the wavelet is convolved down every line and trace of the matrix. After vertical convolution, the data volume is convolved with the wavelet in a direction lengthwise across the survey area. Figure 9 illustrates the convolution of the wavelet in the direction 28 lengthwise across the data volume, across every 12 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 CA 02265063 l999-03- 10 line, and at every depth increment of the matrix. Then the matrix is convolved with a wavelet in the third orthogonal direction, as is illustrated in figure 10. In this case, the wavelet is convolved in the direction 29 across every trace and at every depth increment of the matrix. The resultant data volume contains zones of non-zero matrix elements, centered on the unit value which represents the location of the magnetic depth solution point. Where magnetic depth solution points are close together, these non-zero zones constructively interfere, or rather, the values contained by their elements are combined additively. In the final step 30 of the preferred embodiment, the three dimensional digital matrix is converted by the computer program into a digital data format that is common to the art of seismic data processing and interpretation. It is contemplated that one of ordinary skill in the prior art will be able to readily apply the present method. It is preferred that the digital matrix be converted to the common SEG-Y data format. The resultant SEG-Y data set is stored by the computer program on the computer data storage device. Once converted to a seismic data format, the data is displayed using a variety of seismic data visualization, processing and interpretational computer software tools, which are commonly available for PC and UNIX system applications. The invention is exemplified by the following example. An aeromagnetic survey was flown over an area which contained a known fault. A 2D seismic line across this fault was also available. The flight lines of the aeromagnetic survey were roughly parallel and spaced at 800 meter intervals. Werner magnetic depth solutions were calculated along each 13 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 flight profile and these points were recorded by the computer according to their XYZ co-ordinates. A rectangular survey prism containing a plurality of bins was defined for the survey area by the computer program. The bins, aligned parallel to the sides of the greater rectangular prism, were 400 meters long, 400 meters wide and 20 meters thick. A three—dimensional zero matrix was created by the computer program to digitally represent all of the bins of the survey prism. The magnetic depth solution points were retrieved from computer memory and the co-ordinates of each point were analyzed sequentially. A value of 1 was added by the computer program to the elements of the three- dimensional digital matrix which corresponded to the bins which contained magnetic depth solutions. Where a plurality of depth points were confined to the same bin, the value of their corresponding digital matrix element became a summation of their contributing unit values. A normal distribution wavelet was constructed for convolution across the data volume in three orthogonal directions. The wavelet was convolved by the computer program across the three dimensional matrix, first in the vertical direction, then in the length direction, and finally in the width direction. In each of the three cases, every element of the matrix was involved in the convolution process. After convolution, the three-dimensional matrix was converted to SEG- Y format, a format common to the art of seismic data analysis. This SEG-Y data file was loaded into the computer software program SeisX, which is commonly used in the art of seismic data interpretation, for visualization and interpretation. The data volume was displayed using a variable intensity 14 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 CA 02265063 l999-03- 10 display parameter, with colors related to the values contained within the elements of the data volume. The data volume could be viewed as vertical profile sections of any orientation across the volume, as horizontal depth slices, or as three-dimensional semi-transparent images. Magnetic maps were also loaded into the computer software program, enabling a simultaneous interpretation of the magnetic field maps with the depth solution data. Furthermore, seismic data which were plotted with depth as the vertical dimension were viewed simultaneously for comparison with the magnetic data. Figure 12 is a two-dimensional seismic data image from a seismic survey acquired within the area of the example aeromagnetic survey. The vertical scale of the seismic display image is in milliseconds. The seismic image shows an intra-sedimentary fault 31 with significant throw. Figure 13 is an image, normally displayed in color, of a vertical cross—section through the example SEG-Y formatted data volume using the seismic interpretation tool, SeisX. Its vertical units are in meters. This cross—section shows a near vertical trend 32 of magnetic depth solutions, believed to be the result of a magnetized fault. This near vertical trend 32 coincides spatially with the fault predicted by the seismic data. In conclusion, the invention is characterized by the following advantages: o the quantity of paper used in the aeromagnetic data interpretation process is greatly reduced; o visualization and interpretation of the depth point solutions in three dimensions becomes possible; 15 10 11 CA 02265063 l999-03- 10 the time required for the interpretation is dramatically reduced, since interpretation of the magnetic field maps and magnetic depth solution data are simultaneous and the interpretation map is created automatically by the seismic interpretation software; simultaneous interpretation and visualization of magnetic data and seismic data is possible; and coherent trends of magnetic depth solution data across consecutive flight line profiles can be easily identified using the three- dimensional visualization and depth slice visualization techniques available using computer tools common to the art of seismic data interpretation. 16

Claims (4)

THE EMBODIMENTS OF THE INVENTION IN WHICH AN
EXCLUSIVE PROPERTY OR PRIVILEGE IS CLAIMED ARE DEFINED AS
FOLLOWS:

1. A method for converting magnetic depth solution points, derived from a magnetic survey conducted over a survey area and having spatial co-ordinates, to a digital format which can be integrated with a seismic visualization and interpretation tool for display, comprising:
storing the co-ordinates of the magnetic depth solution points in computer memory;
establishing in computer memory a matrix representative of a pre-determined volume underlying the survey area and which contains a distribution of the points, said matrix being sub-divided into elements;
retrieving the co-ordinates for each point, identifying the matrix element with which it corresponds spatially and assigning a value to that element;
combining the assigned values for each element to create a product matrix representing the distribution of the points in the pre-determined volume; and converting the product matrix into seismic data format and displaying it.

2. The method as set forth in claim 1 comprising:
smoothing the elements assigned value over adjacent matrix elements.

3. The method as set forth in claim 2 wherein:
smoothing is conducted by convolving the matrix with wavelets.

4. A method for converting magnetic depth solution points, derived from a magnetic survey conducted over a survey area and having spatial co-ordinates, to a digital format which can be integrated with a seismic visualization and interpretation tool, comprising:
storing the co-ordinates of the magnetic depth solution points in computer memory;
establishing in computer memory a matrix representative of a pre-determined prism volume underlying the survey area and which contains a distribution of the points, said matrix being sub-divided into elements representative of prism volume bins;
retrieving the co-ordinates for each point, identifying the matrix element with which it corresponds spatially and assigning a value to that element;
combining the assigned values for each element to create a product matrix representing the distribution of the points in the prism volume;
smoothing the elements assigned value over adjacent matrix elements by convolving the matrix with wavelets in three orthogonal directions; and converting the product matrix into seismic data format and displaying it.