Home

User`s guide

image

Contents

1. button Once the gridding has been made the Make grid button becomes inactive and the Plot data button will show the original interpolated data To be able to re grid the original data using different step size for example one needs to apply the Reset grid button 24 Successive utilization of Plot inverse and Plot data buttons can be used to visualize the changes that are made to the data due to the frequency domain operations As an alternative the Plot diff button can be used to visualize the difference between the original interpolated data and the inverse transformed data Pressing Padding button will invalidate both the FFT and current inverse results Note Sequential Fourier operations e g first reduction to pole and then vertical gradient can be performed using the Edit Inverse Data menu item In this case the data are already defined on a regular grid as if ASC grid file had been read in and the Make grid button will be inactive The Plot data button will show the original data after previous Fourier operation and the Revert button restores the original irregularly sampled data Alternatively to make sequential Fourier operations one can save the inverse results either as DAT or ASC grid file and read the file as new input data 4 6 Frequency filtering Frequency filtering affects current FFT spectrum on which gradients components and all other processing functions are computed This means that frequency fil
2. Similarly in directional filtering the user is instructed to draw a azimuthal line along the circles to provide the upper and lower angle of the filter sector The width of the cut off range is defined separately defining a single point some distance away from the upper or lower filter sector The filter cut off range will be adjusted symmetrically at both sides Note The ring and sector values can be given in reverse order outwards vs inwards or clockwise vs counter clockwise in which case they will be automatically arranged properly 27 The interactive editing operation can be cancelled drawing outside the graph area The filter rings and sectors are shown by dotted lines on the frequency spectrum In addition the wave lengths corresponding to inner and outer ring and the limit angles of the sector are shown on the information text of the graph After the user gets accustomed with interactive filtering one should learn manual number based filtering which is more accurate than interactive filtering Manual filtering is made by providing the wave lengths in the units of data dimensions of the filter rings low pass and high pass filtering in the input dialog that appears after the filtering is initiated 4 6 4 Automatic low pass filtering The results of most frequency domain filtering operations especially the second degree gradients are strongly affected by the smoothness and continuity of the data To circumvent comp
3. amplitude and phase spectrum The discrete 2 D Fourier transform 1s computed using fast Fourier transform FFT algorithm Claerbout 1976 More information about FFT algorithm can be found from Press et al 1988 for example The Fourier transform represents a sum of sine and cosine terms of different spatial frequencies or wave numbers k and k that are based on data coverage D max x min x and D max y min y and data sampling d and d along x and y directions The highest spatial frequency is the so called Nyquist frequency e g max k 0 5 d The lowest frequency is based on the data coverage e g min k 0 5 D Considering that the reciprocal of the spatial frequency represents wave length A I k zero frequency means infinite wave length 1 e constant level of the data Because of the properties of the Fourier transform symmetry linearity shift and derivative properties several computational operations can be performed in Fourier transformed frequency kx ky domain more efficiently than in the spatial x y domain For more detailed information about Fourier transform methods in potential field analysis please see Blakely 1995 2 Installing FOURPOT The FOURPOT program can be run on a PC under Windows system and a graphics display with at least 1280x1024 pixels in size FOURPOT uses dynamic memory allocation and reserves memory for several auxiliary data matrices derived from the interpolated
4. and blanked out by white color in the map graphs See chapter 5 for more information about file formats After successful data input the program activates Plot data and Reset grid buttons and plots a contour or image map as seen in Appendix A If data were read from a DAT file it is assumed to be irregularly sampled and the Make grid button is active If data were read from an ASC grid file the Make grid button will be inactive If Reset grid button is now pressed the original ASC data will be lost and considered as irregularly spaced data read from a DAT file 20 2 D Fourier analysis version 1 3a by M Pirttijarvi 2014 File Edit Process Tools Exit Plot data Shortcuts Reset grid Upwad 2 D Fourier transform analysis Make grid L Dowiwd Original data Step x N ooo grad Directions ddd Height units oo Incl Elev deg zs Dect de 2870 2880 Dimension 52x 48 Spacing 1 000000 x 1 000000 km Field Figure 1 FOURPOT GUI after reading in the data and before gridding it 4 2 Interpolation All input data read from DAT files are assumed to be irregularly sampled Even if the original data were already evenly discretized it will be re interpolated before it s passed to FFT Therefore after reading in the data the user must first provide suitable values for the number of grid points N and M or grid stepping dx and dy in the Step x and Step y text fields and then press the Make grid button t
5. and place the grid nodes nicely and or to select a smaller area from large map The Operation can be made either manually giving the x and y coordinates of the limits or interactively approximately with the mouse In the former case the program shows an input dialog where the x and y limit coordinates of the area are given In the latter case the program enters interactive editing mode and the user draws a rectangular clip area using the mouse as described in the following After providing new coordinate limits the program redraws the data map and the gridding see chapter 4 2 of the data can be made using new data limits Note Interactive editing operations including the zooming and clipping follow the same guidelines When editing is initiated the mouse cursor becomes a cross hair cursor above the graph area and most of the GUI widgets become inactive so that the user cannot escape the editing mode accidentally A rectangle zooming clipping line frequency filtering pick profile or point direction filtering panning is then drawn with the left mouse button pressed Before the first editing operation the program shows on screen instructions In case of Zoom Clip area operation the user can revert back to the original limit values by drawing the selection rectangle outside the graph area or by giving a blank line in the manual input dialog The information text next to the graph may show the dimensions incorrectly unless
6. coordinates to get them right The correct y coordinates are computed as y y yo dy yo where yo 1s the y coordinate of the origin yllcorner and dy dy dx is the ratio of the true y step with respect to the x step The data are stored as rows from top to bottom but the header defines the x and y coordinates of the lower left corner 5 3 Regularly gridded matrix files The matrix format MAT assumes that the data that are already regularly sampled and stored without x and y coordinates The matrix format is provided for better functionality with programs such as Matlab and Maple The format of a matrix file DAT is illustrated in the example below 64 64 43 L60R 09 A3 3208709 0 43A TUE U9 The header line defines the number of data points in x and y directions NP N and MP M The following NOP NP x MP lines contain the actual data values is a single column 38 Since the matrix file does not contain coordinates it is important to know the order in which the data are stored By default the data is stored in column wise fashion from bottom to top and from left to right The origin is always located in the bottom left corner Using generic programming notation we would have do a kine do j 1 mp read or write f i j end do end do The default order can be changed so that for each column the rows will be read by giving a negative value for the number of points in x direction NP NP This works wit
7. data This may become a problem for the 32 bit version which is limited to 1 GB of continuous memory The 64 bit version can use all available system memory and is generally faster because of additional OpenMP parallelizations In case of large data arrays however the CPU becomes critical factor because the 2D FFT of matrices larger than 4096 by 4096 elements takes lots of time The program has simple graphical user interface GUI that is used to handle file input and output to initiate the computational operations and to visualize the data and the results The user interface and the data visualization are based on DISLIN graphics library http www dislin de The program requires either the 32 bit FOURPOT EXE or the 64 bit FOURPOT64 EXE executable file The 64 bit version requires the presence of the LIBIOMPS5MD DLL file for OpenMP parallelizations The distribution file FOURPOT ZIP also contains a description file _README TXT a user s manual FOURPOT_MANU PDF and some example data files DAT To install the program simply decompress it somewhere on your hard disk 2 1 Notes on version 1 3a Several small changes have been implemented in version 3a The most important changes are 1 Data files are now checked before reading and hence the header of the preformatted data files DAT does no longer contain the number of lines only the column indices 2 Interactive editing uses now rectangles and lines which are vis
8. data area and padded margins Together padding and tapering effectively remove so called Gibbs phenomenon as well as ringing and other artifacts in the inverse transformed data Padding will be discussed more in chapter 4 3 FOURPOT does not care what kind of data it operates on Strictly speaking pole reduction and pseudo gravity for example should be applied to magnetic data only For teaching and testing purposes the Data type menu item allows defining the actual data type undefined magnetic or gravity In practice this has effect only on the information text next to the graph the scaling of the radial amplitude spectrum see chapter 4 9 and the proper selection of the output model files GRABODS vs MAGBODS The Miscellaneous sub menu contains following items select GUI directed or manually a rectangular area from the data oom clip area swap between contour map and image pixel map modes Contour Image Show Hide points show hide the original data points in original data view Meters Kilometers swap dimension labels between meters and kilometers Unit conversion scale x or y distances with user defined values Data conversion scale data values with user defined value Step number distance swap gridding between steps number or step distance The minimum and maximum coordinate limits of the data matrix are based on the original data coverage The Zoom Clip area item can be used to redefine and round the coordinate limits
9. following items Median filter 3x3 smooth the data using a median of a moving 3x3 window Ridge tools pick continuous ridges fingerprint directions and remove them Profile Model tools pick a profile and create initial 2 5D prism model automatically Radial spectrum show the radial spectrum used for depth analysis Median filter 3x3 can be used to smooth gridded data and to remove outliers spikes from it Median filter can greatly improve the results of many Fourier operations which assume that the data and its derivatives are continuous The Pick ridges 3x3 and Pick ridges 5x5 items utilize an experimental algorithm that investigates if the grid point at the center of a moving 3x3 or 5x5 window represents a continuous ridge like anomaly peak Ridge is defined by three coordinate pairs the first indicates the center of the ridge the second and the third indicate where the ridge comes from and where it goes to The ridges can be straight e g W E or SE NW or bent up to 90 degrees eq W SE or SE NE Ridge points can be saved into a column formatted text file using File Save ridge data menu item The ridge data can also be saved in Atlas BNA format Hint Ridge picking is usually applied on horizontal gradient data to locate geological contacts When the operation is re applied on upward continued data one can estimate the dip directions of the contacts To do this one should save the results of the ridge picking made on horizontal gradient of or
10. for example the density contrast or magnetic susceptibility depth extent and dip angle of the geological targets giving rise to potential field anomalies The tilt gradient however contains useful information that can be used a to locate the number and the center position of the anomalous targets and b to divide the data into sub anomalies i e sections of data the inversion should use per each model prism FOURPOT uses the tilt gradient data to generate a semi automatic model for the GRABODS and MAGBODS modeling and inversion programs of M Pirttij rvi These programs which are still under development and have not been put to public distribution can be used to interpret single profile data using multiple dipping prism models 32 To create the initial model the user must first define the regional field 1 e the base anomaly for the profile This 1s made pressing Edit base button Once the interactive mode is entered the user needs to press the left mouse button over the response graph to define the knot points through which the base anomaly should go The knot points must be given in an ascending order from left to right The editing mode is ended with a single click of right mouse button If more than two base knot points are defined cubic splines are used to compute the regional trend see Fig 5 If only one or two knot points are defined the base anomaly will be a horizontal or dipping line Existing base anomaly can be removed pr
11. graph is centered so that the origin k k 0 is located at the center and the horizontal and vertical axes range between the wave numbers knx kny related to the negative and positive Nyquist frequencies fnx fry the values of which are shown in the information text next to the graph The lowest frequencies largest wave lengths widest data variations are located in the middle of the graph and the highest frequencies shortest wave lengths are located far from the origin Continuous linear features appear as linear features in the FFT spectrum as well the angle however is reversed 4 5 Inverse FFT and difference Once the FFT has been computed and the 2 D frequency spectrum is available one can perform non automatic frequency filtering and the various frequency domain processing tasks The Plot inverse button will display the current inverse results The Plot diff button will display the difference between the original gridded data and the inverse data If automatic low pass filtering is not active the Plot inverse will show the inverse of the original FFT and the difference will be very close to zero If automatic low pass filtering is active the inverse results show a smoothed version of the original map and the difference will show the high pass filtered parts of the data When processing tasks are made Process menu or quick task buttons the inverse results will be shown automatically 1 e without pressing Plot inverse
12. graph showing the inverse results or picked profile is redrawn Rotation is not made if the current view is not one of the abovementioned inverse results because they are unaffected by the rotation The Reset rotation button resets the total rotation to zero The basic gradients 1 st and 2 nd degree gradients are kept in computer memory so that derived components can be computed faster Therefore when basic gradient components is rotated all other gradient components are affected as well For example rotating d dx gradient 90 degrees makes it equal to d dy and d dy becomes equal to d dx After resetting the x and y gradients are again directed along the original x and y axes of the data In pole reduction the transformed field becomes vertical m f 0 0 1 and declination does not play any role Normally in equator reduction the transformed magnetization and field will be rotated towards north m f 0 1 0 based to the given value of declination However the direction of the transformed field reduced to equator can be rotated For example the reduced field can be computed in EW direction using an angle of 90 m f 1 0 0 When performing sun shading and generalized derivative the elevation and azimuth angles are read from the same two text fields as the inclination and declination used in pole and equator reduction and pseudo field computations Instead of editing the azimuth angle manually the Rotate button allows fast ro
13. on max horizontal gradient and principal directions shown on top of horizontal gradient 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Inverse data 1830 1830 1820 1820 gt 1810 gt 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 Dimension 64x 64 1 000 x 1 000 km X Spacing 1 000x 1 000 km Spacing Low Pass filtered 4 00 3 00 Low Pass filtered 4 00 3 00 Horizontal gradient 0 5 10 15 Max horizontal gradient Dimension 64x 64 Field Field 54
14. the Meters Kilometers menu item is utilized Although knowledge about the spatial dimensions is usually not required the pseudo gravity and pseudo magnetic field need to know the real dimensions for the amplitude of the transformed field to be more or less truthful Unit conversion can be used for example to scale all distances from miles into kilometers or from meters to kilometers Likewise Data conversion can be used for example to scale the data from mgal m into E tv s mgal km It can also be used to reverse the sign of the data e g the negative gradient of the potential In both cases the program shows an input dialog where the user provides a single character for a mathematical operation addition subtraction or x multiplication and or division and the numerical value For example string 100 will multiply all data by one hundred The operator has to be the first character on the given line An empty line or zero value cancels the operation Step number distance is a mode switch for grid sampling that changes the meaning of Step x and Step y text fields in the left control panel The original data are considered to be irregularly sampled and must me gridded and padded before FFT When Step number mode is active the number of grid nodes Nx Ny in x and y directions integer values must be given When Step distance mode is active the distance between two grid nodes dx dy in x and y direction
15. 2 000000 m Log oF Appendix D Difference between original data and inverse transformed data 2 D Fourier transform analysis Difference 1830 1820 1810 1800 1790 2840 2850 2860 2870 2880 Dimension 52x 48 Spacing 1 000000 x 1 000000 m X 0 015 0 006 0 003 0 011 Field x 1073 47 Appendix E Low pass filtered frequency spectrum 2 D Fourier transform analysis Fourier transform 0 1 0 0 Ky 0 1 0 4 0 4 0 2 0 0 0 2 0 4 Dimension 64x 64 Spacing 1562500E 01 x 1562500E 01 1 m Wavelength 2 000000 2 000000 m 2 0 0 3 1 3 3 0 Low Pass filtered 4 000000 2 998501 Appendix F Low pass filtered data and the difference with the original data high pass filtered data 2 D Fourier transform analysis Difference 2 D Fourier transform analysis Inverse data 1830 1830 gt 1810 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 X X 17 9 35 61 4 5 1 5 1 4 4 4 Field Field 48 Appendix G Low pass and direction filtered frequency spectrum 2 D Fourier transform analysis Fourier transform 0 4 0 2 0 0 0 2 0 4 Kx Dimension 4x 64 Spacing 1562500E 01 x 1562500E 01 1 m Nyquist 5000000 5000000 1 m Wave length 2 000000 2 000000 m 2 0 0 5 1 0 2 6 Low Pass filtered 4 000000 2 998501 Log oF Direction filtered 15 52443 15 33253 Appendix H Direction filtered data and difference with the original data a
16. 5 deg and Azim 10 deg shown transparently on the interpolated low pass filtered data 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X 17 9 35 61 Field 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X 17 9 35 61 Field 53 Appendix Q FOURPOT GUI Sun shading Elev 20 deg and Azim 135 deg shown transparently on the low pass filtered data using slightly altered color scale om 2 D Fourier analysis version 1 3a by M Pirttij rvi 2014 File Edit Process Tools Exit Plot data Shortcuts Reset grid __Upwad 2 D Fourier transform analysis Make gid __Downwd Inverse data Padding __PoleRed Plot FFT _EquaRed Plot inverse Zsa Plot diff __Xgrad Step x N 52 o Ygad Step Y M 48 _ Vade Directions ddad Height units oo oo e d dydy Incl Elev deg d dedy Decl Azim deg d dxdz Iv Activate rotation d dydz Rotation angle deg 45 0 HoriGr _ Roae _ otata TotalG Reset rotation TitGr A Transparent shading MaxHGr Range Theta al E SunShed Center __GenDeri_ lt gt Undo Levels ai ai Update 2870 2880 Dimension 52x 48 Spacing 1 000000 x 1 000000 km Low Pass filtered 4 000000 2 998501 Sun shading Elev amp Azimuth 20 00000 amp 135 deg Field Appendix R Ridge picking 3x3 made
17. 80 xX xX 15 6 4 13 5 as 2 6 Field Field Appendix L Horizontal gradient and total gradient analytical signal 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Inverse data 1830 1830 1820 1820 gt 1810 gt 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 X gt 4 0 5 10 15 0 6 11 16 Field Field 51 Appendix M Tilt derivative and theta map 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Inverse data 1830 1830 1820 1820 gt 1810 gt 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 X x 89 29 30 90 0 30 60 90 Field Field Appendix N Gravity potential and pseudo magnetic field k 0 01 SI Ag 0 2 g cm and To 52000 nT 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Inverse data 1830 1830 1820 1820 gt 1810 gt 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 X X 248 67 115 296 458 161 136 433 Field Field 52 Appendix O Sun shading and general derivative Elev 75 deg and Azim 10 deg 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 p4 3 0 0 8 1 5 3 7 Field 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X 1 0 0 3 0 3 1 0 Field Appendix P Sun shading and general derivative Elev 7
18. FOURPOT Potential field data processing and analysis of using 2 D Fourier transform User s guide to version 1 3a G Markku Pirttij rvi 2014 markku pirttijarvi at gmail com Table of contents jae Colne odu i iea eee ema ween eee eran rer ern ern eine enn E ct er een reer 4 Za MENS CAME FOURPO Tursin A R E eal eavieds 6 PENOS ON VERSION L daaa E N nneunea i eaaelas 6 N A a isear acta clsse oun teu aa eacear aes ceiele cee ec aoe seea etcetera 7 SE Gl SEEE Ea U Meee ee er EAE a ere e ne Mg me one ne teen tc ieee NE ener emer ete enn nN eee 7 PUASA 01s 016 ee eee Ree Cee Re ST PORTE NOE ee ee aera Peer Sere er eee T er rer ee eet 8 Fo et OCE SMH Dt 06 apenas A E nT ee oe enter ee ee ree eee Tee 11 3 3 1 Upward and downward continuation esssesseeoeeeesssssssoerressssssseeeeresssssseeeeees 1 J92 Pokandgeguator TECUCH OM x jaca EEE 13 Je ASIC SLACICIIG dsi cei E eu edenesieeues 14 3o Ae Denved CRA C IIS parasite asrrcdndie diets saute alata ste sald tee nln dio ta alate ist asta sald tedden diate eas 15 3 5 0 Potentials and pseudo fields saisis a 15 IO HUDENS T OT Mierna rr es em gat ee Pn 16 3 3 7 Sun shading and generalized derivative cc cceeeccccccccccceseeeseeececeeessaeeeseeeeeees 17 Ae MOOS MAN IM ete Recetas ahaa ata et aetna ta net cic hte ts ba came niet hes 18 rE I MMC NW act a eer E E E E A 19 4s IS i CNS OO OR ae e cat atae deaae anatase ecaeens 20 Al Readine m CNC dataserier e e E EEE 20 AEDO
19. Kilometers menu item When the model is saved the user is given an option to save the profile data as discussed in chapter 4 9 1 4 10 Radial amplitude spectrum The items in the Tools Radial spectrum sub menu become operational when the 2D Fourier spectrum is the active graph Show spectrum item shows the amplitude spectrum as a function of radial wave number k The radial ampl7itude is computed as the mean of the 2 D Fourier amplitude spectrum A FI Re F Im F 7 along rings with radius k ky k 34 centered on the origin k k 0 The amplitude spectrum has traditionally been used to estimate the depth to the top of the potential field sources This is accomplished by fitting linear lines to the decaying amplitude curve on a semi logarithmic scale An example is shown in Fig 5 Dimension 64x 64 5 Spacing 1 000 x 1 000 km Rad spectre depths Slope 2 437 km 4 Slope 0 746 km o log F Kr Figure 6 Radial amplitude spectrum with two manually fitted lines The slopes of the lines give estimates for the depth of the potential field sources The spectrum is related to the gravity data shown in the Appendices In frequency domain the Fourier transform of a potential field can be formulated as F Ce giving log F C h k Thus the depth A to the top of the source 1s related to the tangent of the line that is fitted to the linear parts of the amplitude spectrum on semi logarithmic scale The coeffici
20. are also used to compute additional gravity tensor related components namely the first and second rotational invariants see Pedersen amp Rasmussen 1990 and the second Falcon gravity gradient or curvature component 24 gyy 2 14 3 3 4 Derived gradients The derived gradient operations include horizontal gradient H total gradient A sometimes also called analytical signal tilt gradient 7 maximum horizontal gradient amplitude MH also called total horizontal derivative horizontal gradient of tilt gradient S and theta map TM Blakely 1995 and Cooper amp Cowan 2006 du 2 du H IE E 8 du du du 2 A z E 3 dU 1 ae T tan TEES 10 dU MH Re tant 11 d T d T Palas age 74 TM tan H 13 The derived gradients are used to enhance small amplitude features in the data to delineate and locate geological targets and contacts and derive information about structural geology Examples of derived gradient operations are shown in Appendix L and M 3 3 5 Potentials and pseudo fields Potential fields are defined as the negative gradient of corresponding scalar potential F V dd dx d dy d dz In Fourier domain the gradients of the field were computed by simple multiplication Eqs 5 7 Thus the potential can be obtained from the vertical component of the field by dividing its Fourier transform by the radial wave number k In case of the magnetic potent
21. axis on the left Processed data e g vertical gradient in Fig 3 are shown against the auxiliary y axis on the right The profile output thus allows detailed comparison of the effects of various processing operations If the user changes the view form profile graph to 2D map 1 e using Plot data Plot inverse or Plot FFT buttons one must use the Show prof button to see the profile response again 31 JUSIPeIS ZP AP 20 Unprocessed data 30 0 5 10 15 20 Distance km Figure 4 Example of the profile response The solid and dashed curve correspond to the unprocessed frequency filtered data left y axis and the first vertical derivative right y axis interpolated along the profile Once the profile has been generated the data along it can be saved into column formatted text file using File Save profile data menu item The file will contain lots of information columns 1 4 contain xyz coordinates and running profile coordinate distance columns 5 12 contain the base level discussed later unprocessed data the first xyz gradients horizontal gradient total gradient tilt gradient and the last column contains whatever processed data is currently active When profile data is saved the program asks if the data are saved with or without the GRABODS MAGBODS header information discussed later 4 9 2 Prism model generation Once profile data has been saved to disk one could start model based inversion to estimate
22. d the smooth bell shaped filter function of c low pass d high pass and e band pass filters 4 6 1 Low pass and high pass filtering In low pass filtering all data outside the outer ring high frequencies are nulled and all data inside the inner ring are preserved see Fig 2a and 2c On the contrary in high pass filtering all data inside the inner ring low frequencies are nulled and all data outside the outer ring are preserved see Fig 2d The radii of the rings are defined as wave length 1 which is the inverse of radial frequency k ky thy where k and ky are the axes of the 2 D frequency graph Note The x and y axes of the 2 D FFT plot are not scaled If sampling dx and dy 1s different in the x and y directions the Nyquist frequencies 1 e the minimum and maximum values of the axes are different In the present version of FOURPOT the frequency filtering uses the same wave length both in x and y directions The lengths of the k and k axes however are the same As a consequence circular filter rings will appear as ellipses in the 2 D FFT plot 26 4 6 2 Directional filtering Directional filtering is equal to the fk filtering commonly used in seismic and GPR data processing Elongated linear features in the data appear as linear features also in the FFT spectrum The slopes are reciprocal so that horizontal features parallel to x axis in the data are shown close to the vertical k axis in FFT spectrum Dir
23. direction of the inducing magnetic field The behavior of Earth s magnetic field can be estimated using the IGRF model international geomagnetic reference field Pole reduction transforms magnetic field data as if it were measured on a magnetic pole where inclination is vertical This helps estimating the true dip and strike directions of the targets Since pole reduction is an unstable operation at low latitudes Equator reduction can be performed instead These operations involve a phase transformation and they can be formulated as Blakely 1995 OlmO Om Of F U F U x 4 Where and are functions of the direction vectors of the original magnetization m m my m and external magnetic field f f f f and O m and O s are functions of the direction vectors of the transformed magnetization and field m and f Unless remanent magnetization is taken into consideration m f and m f Note The inclination and the declination angles in degrees of the original magnetic field and magnetization are read from the ncl and Decl text fields in the left control panel Inclination J 90 is taken from horizontal plane and it is positive downwards northern hemisphere and negative upwards southern hemisphere Declination D 180 is taken from the direction of positive y axis geographical north and it is positive in clock wise orientation See chapter 4 7 for information about rotating the equator reduc
24. ectional filter removes linear features from the data by nulling frequency content inside given sector see Fig 2b The directional filter is defined by two angles that correspond either a to the upper and lower limit of the sector interactive mode or b to the center and the width of the sector manual mode The angles are defined in degrees in counter clockwise direction taken from the positive k axis Because of the symmetry of the frequency spectrum the filter will be automatically mirrored around the origin Thus the user needs to define a single sector either above the k axis or left to the k axis In addition the width of the cut off range must be defined by providing a third angle 4 6 3 Interactive and manual filtering Frequency filtering can be done either a interactively with the computer mouse or b manually by providing the numerical values of the filter cut off ranges The selection is made when the corresponding menu item is chosen in the Edit menu When interactive mode is chosen most of the other program widgets become inactive and the mouse cursor changes from normal arrow into a cross hair above the graph area note only above the graph area The 2 D frequency spectrum is shown together with auxiliary polar coordinate grid The user is given instructions to draw a radial line along away from the origin The distance from the origin determines the radii of the inner and outer filter rings in terms of wave lengths
25. ed data 13 3 3 3 Basic gradients Vertical dU dz and horizontal gradients dU dx and dU dy are computed automatically after FFT and kept in compute memory so that other derived gradient components can be computed faster Also second degree gradients are computed automatically and kept in memory so that they can be rotated horizontally faster Here x and y directions represents the horizontal axis East and vertical axis North of the mapped data The vertical gradients dU dz for magnetic data and d U dz for gravity data are particularly useful in sharpening the anomalies below the anomalous sources The less common gradient components d U dxdy d U dxdz d U dydz correspond to the off diagonal components of the gravity tensor 8xy 2x and gyz provided that they are computed from the gravity potential Examples of gradient operations are shown in Appendices J and K In frequency domain the Fourier transforms of the basic gradients F d U dx are computed by multiplying the Fourier transform of the field F U with wave numbers kx ky and k The gradient operations can be formulated as Blakely 1995 F ik F U F ikyF U and F 2 k F U 5 r 4 k F U r Z k F U and r k 2F U 6 and r k kyF U F z k k F U and r ee ky k F U 7 where n 1 or 2 is the degree of the gradient and i is the imaginary number i 1 The above mentioned basic gradients
26. effective Because the scale widgets give slightly different results for different graphs one can reset the color mapping parameters to their defaults using the Reset quick button Please be reminded that FOURPOT was not designed to be used as a plotting program The user is advised to save the mapped results into text files and use more advanced programs like Golden Software s Surfer for publication quality plotting 4 9 Profile Model tools The six push buttons below the Levels scale widget at the bottom of the left control panel are used a to see the results of various filtering operations as a line graph along a straight user defined profile and b to create an initial model based on 3 D dipping prisms for the interpretation of the potential field data along the picked profile These tasks are also available through the Tools Profile Model tools sub menu 4 9 1 Profile plotting Pick prof button is used to define the x and y coordinates of the profile start and end This can be done either interactively with the mouse or more precisely giving the coordinates in an input dialog After applying the Pick prof button the gridded and frequency filtered data is interpolated bilinear method along the profile and shown in the graph area An example of the profile response is shown in Fig 4 Unless some processing operation has already been made the profile graph shows only the unprocessed frequency filtered inverse data against the y
27. ent C depends on the selected data type Edit Data type for generic undefined data C 1 for gravity data C 1 k kyk and for magnetic data C 1 k k Depending on the selected data type and the data the tangents give different depth estimates that may or may not be realistic In such cases one should prefer to the un normalized spectrum C 1 Please be noted that the method used for depth determination has not been verified Refer to Bhattacharyya 1967 and Ruotoistenmaki 1987 for more information about the use of radial amplitude spectrum in depth determination Interactive line fitting is initiated using the Add new line item in the Tools Radial spectrum sub menu The normal mouse pointer arrow will change into a cross hair cursor over the graph area and most of the GUI widgets become inactive The user should then press the left 35 mouse button at the starting point of the line over the radial amplitude spectrum While keeping the left mouse button pressed and moving the mouse a line is being updated and the user should release the left mouse button at the location that defines the location and angle of the line for the best The graph will be redrawn with the new line on it and the slope of the line will be displayed on the information text on the right side of the graph Multiple lines can be added to the same graph The Delete last line item removes the most recent line Note Because of the cut off range the radial s
28. essing the left mouse button once outside to the left or to the right the graph area New knots can be added to existing base anomaly pressing the left mouse button once either above or below the profile graph and then once at the new knot location over the graph Likewise knots can be removed from the base anomaly pressing the left mouse button once over the graphs close to the selected point and then once above or below the response graph Z 38 4 T 18 id ge D 2 S 2 22 Z E 42 62 82 af 2 F 38 7 5 oF nn O 4 a 3 0 5 10 15 20 Distance km Figure 5 Example of the prism model generated using tilt gradient data The solid and dashed curves are the unprocessed frequency filtered and the tilt gradient interpolated along the profile The dotted line is the base level fitted using splines through the three base level knots The dotted horizontal line at the top of prisms is the width of the sub areas to be used in the GRABODS inversion algorithm 33 Once the base anomaly has been defined the Make mod button for the automatic prism model becomes active The number and location of the prism models is based on the positive peaks of the tilt gradient data The tilt gradient varies between 90 and 90 degrees Usually positive potential field anomalies are reflected as positive values of tilt gradient and its sign becomes negative between two anomalies The only input parameter needed for model genera
29. fuilters based on the local phase Computers amp Geosciences 1585 1591 Cooper G R J amp Cowan D R 2011 A generalized derivative operator for potential field data Geoph Prospecting 59 188 194 Pedersen L B and Rasmussen T M 1990 The gradient tensor of potential field anomalies some implications on data collection and data processing of maps Geophysics 55 1558 1566 Press W H Flannery B P Teukolsky S A amp Vetterling W T 1988 Numerical recipes in C The art of scientific computing Cambridge Univ Press 735 p Ruotoistenm ki T 1987 Estimation of depth to potential field sources using the Fourier amplitude spectrum Bulletin 340 Geological Survey of Finland PhD thesis 44 8 Terms of use and disclaimer The FOURPOT software is free for personal and academic use It should not be used for commercial purposes without permission and must not be sold or distributed onwards for profit If you decide to publish results computed with FOURPOT please make a reference to this user manual and FOURPOT homepage at https wiki oulu fi x OoU7AQ If you find FOURPOT useful please send me e mail The FOURPOT software is provided as is The author MP and the University of Oulu disclaim all warranties expressed or implied with regard to this software In no event shall the author or the University of Oulu be liable for any indirect or consequential damages or any damages whatsoever resulting from los
30. h data stored using following do loop do j 1 mp d i Lanp read or write f 1 J end do end do Furthermore if MP lt 0 the data are stored in multiple columns If NP gt 0 and MP lt 0 each line is considered to be one y columns of the mapped data as in the default case If NP lt 0 and MP lt Q the lines are considered to be data rows as in the ASC file Using generic programming notation we would have in this case do 3 ah read LEI Ty 2 l np end do Note Because the matrix format does not contain coordinates the sampling frequencies and corresponding wave lengths shown by FOURPOT are determined by an input dialog window that appears on the screen after opening the MAT file The input dialog defines the x and y coordinates of the origin and the x and y step between the matrix elements If the x and or y steps are negative the axis gets reversed Moreover the data are mirrored in x or y direction if the last two parameters in the dialog are set to one instead of zero Because normally the origin is considered to be in the lower left corner the negative steps and mirroring allow redefining the order in which the data are stored 39 5 4 Graph parameters Several graph parameters and text strings can be changed by manually editing the FOURPOT DIS file with a normal text editor This allows translating the graphs to suit the needs of the user bit better The format of the FOURPOT DIS file is important If it becomes co
31. he screen 19 4 Using the program At start up FOURPOT reads graph parameters and some additional settings from the FOURPOT DIS file see chapter 5 4 for more information If this file does not exist default parameters will be assigned and the file is created automatically The program then builds up the GUI shown in Fig 1 None of the program controls widgets however are active before data has been read in The data processing and analysis is performed in few successive steps that are discussed next The preliminary step is to check that the input data file is in correct format so that it can be read in see chapter 5 on file formats After reading in the data one should apply a Make grid b Padding and c Plot FFT buttons to perform interpolation and padding and to compute the 2 D FFT and its inverse Only after the FFT has been computed one can perform frequency filtering Edit menu and other frequency domain data processing Process menu and take a look at the inverse Plot inverse and difference Plot difference 4 1 Reading in the data The first thing to do after starting up the program is to read in the data using the Read data DAT Read grid ASC or Read matrix MAT menu items Note that Geosoft XYZ formatted files can be read as normal DAT files provided that the file header is suitable Also note any missing data values inside ASC files are replaced for data processing by the median of the special points
32. ial the direction of the inducing field must be taken into account Furthermore the gravity and magnetic potentials and hence the fields of an anomalous target with the same dimensions but different petrophysical properties density and magnetic susceptibility are interrelated via Poisson s relation Thus it is possible to compute so called 15 pseudo gravity field from the measured anomalous magnetic data and anomalous pseudo magnetic field from gravity data The pseudo gravity and magnetic fields gps and AT and the gravity and magnetic potentials and V are computed from following equations Blakely 1995 F 9ps x F AT and F ATps F On 0 x Flg 04 Fld 7 x F g and F V Pa x F AT 15 where AT is the anomalous magnetic field reduced to the pole g is the vertical component of the gravity field and A is a function that depends on the density contrast and the intensity of magnetization which in turn is a function of magnetic susceptibility and intensity of the inducing magnetic field The values of inclination and declination defined in GUI window are used for the direction of magnetization and magnetic field The program asks for the values of density contrast g cm susceptibility SI and intensity of the inducing magnetic field nT Examples of gravity potential and pseudo magnetic field are shown in Appendix N 3 3 6 Hilbert transform Hilbert transform is a phase shift operati
33. ible during editing and thus help user in the drawing task 3 An error that created unnecessary missing data points was fixed in data interpolation 3 Menu items 3 1 File menu The File menu contains following items Read data dat read in irregularly spaced data in default data format s Read grid asc read in regularly spaced data in ASC ascii grid format read in regularly spaced data in matrix format k save the data of the current view in general xyz column format Save current data dat Save current data grid asc Save current data matrix mat Save ridge data Save profile data Save model file Read BNA map Read DIS file Save graph as PS Save graph as EPS Save graph as PDF Save graph as WMF Save graph as PNG Save graph as GIF save the data of the current view in a ASC grid format save the data of the current view in a matrix format save the discrete data points derived from ridge analysis save the data interpolated along a user defined profile save the GRABODS MAGBODS model file read an overlay map file in Atlas BNA file format read new graph parameters from a DIS file save the current graph in Adobe s Postscript format save the current graph in Adobe s Encapsulated PS format save the current graph in Adobe s PDF format save the current graph in Windows metafile format save the current graph in Portable Network Graphics format save the current graph in GIF file format These
34. iginal data in BNA format perform upward continuation apply Inverse data menu item perform the ridge picking on the upward continued data and read in the previous BNA file for comparison Pick directions is somewhat similar operation that computes the orientation or principal axis directions of the data based on the x and y gradients of the data Directions are saved as angles in degrees taken from the horizontal x axis east towards positive y axis north Both 18 ridge picking and computation of principal directions query for a threshold value in percent that will ignore smaller features when increased The items in the Profile Model tools and Radial spectrum sub menus are discussed in chapters 4 9 and 4 10 3 5 Exit menu Exit menu has two items The first one restarts the whole program GUI window so that it fits better either normal 4 3 or widescreen displays When changing to widescreen mode the program asks for a scaling value for the aspect ratio Value 1 0 is used for 4 3 displays and values between 0 65 0 85 suit widescreen displays OK to Exit item saves FOURPOT DIS file automatically and closes the program without any additional warnings Therefore the user should take care to save the results before exiting Errors that are encountered before the GUI starts up can be found from the FOURPOT ERR file When operating in GUI mode run time errors arising from improper parameter values for example are shown on t
35. iple columns of which only one will be read Note Since version 1 3a the number of lines does not need to be defined in the header line M If the first line starts with a comment character slash fence exclamation or lower case l or upper case L the column indices ICO1 ICO2 ICO3 can be given interactively using an input dialog that appears Because FOURPOT ignores lines starting with these special characters FOURPOT can read Geosoft XYZ files X YZ An example below illustrates the generic column file format EEE MM Line 1 0 00 0 00 0 8665106E 01 0 33489 T0200 Q00 O 855865LE 01 0 44134 20 400 O00 NII 0020 SRT OE 0265474 37 5 2 ASC grid files The ArcInfo ASCII grid file ASC or ESRI grid format is supported because it is quite common and Golden Software Surfer supports it directly The example below illustrates the ASC format Please see the Wikipedia article http en wikipedia org wiki Esri_grid for the details of the file format ncols 141 nrows 241 xllcorner 3097500 yllcorner 6597500 cellsize 5000 nodata value 1 7014100000000001E 038 94 30492Z0061652583 S149 64149697181 1 49 685123416514031 6 3991202604341 52 4 6 Note The ASC format requires square grid cells dx dy However FOURPOT saves the data in ASC file even if the x and y steps are not equal In this case cellsize in the header corresponds to grid step in x direction and the user needs to manipulate the y
36. irregularly sampled In the latter case the data will be interpolated on a regular grid first The processing operations are performed in frequency domain accomplished by a 2 D Fourier transform The frequency domain operations include Frequency filtering low pass high pass amp angular direction filtering Upward and downward continuation Reduction to magnetic pole and equator with rotation st vertical and horizontal gradients d dz d dx amp d dy 2 nd degree gradients d dx d dy d dz d dxdy d dxdz amp d dydz Rotational invariants of gravity tensor and Falcon curvature response Horizontal gradient and Total gradient Tilt gradient Theta map and Maximum horizontal gradient ee ee a ae ee er Pseudo magnetic field and pseudo gravimetric field p lt Gravity and magnetic potential 11 Hilbert transform 1D or 2D 12 Sun shading and Generalized derivative with rotation Considering a continuous function f x with continuous first derivatives the Fourier transform F k and the inverse Fourier transform can be written as Blakely 1995 F k eO x dx and f x 2 e F k_ dk Aa 1 Considering discrete 2 D data fx f x y the governing equations become N M E 1 N IM 1 Dae Le ee oie fri and fri ae D gt 2 k d 1 NM n 0m 0 where N and M are the number of data values in x and y directions The transform is complex which means that it consists of both
37. lO enea S S ZZ A a Padd SAIN CMAP UL IN Sele dahed cn dual ws ual A 22 AA Fast Fourier transform seses a a a A T 24 2 vee FEET ando eren CC ecan A atates datiesceneesetjaaters 24 ARR EE E a EE A O EEE A ES staan E T ST EE 25 4 6 1 Low pass and high pass filtering xcs tisae cab hexssstecichiickie veer celehouseteena toed 26 A OZ LIPS COOMA TUSTIN sate ssuteit esi annkineudtins weaantensaudadh waAraneiien dite mianatebanandidhateAeneniser ZI 46 3 Interac ve an iManUal Tenne eienn a a 2T 4 6 4 Automatic low pass filtering cceccccccccccssssssseeceeeeeeaaeeeseeecceeeeeeaeeeseseeeeeeeeeaas 28 4 7 Rotati n ANG Transparent SHACING sac acted masivassvaueetamduiwassdademabivesawedimdulenstiewetts 29 4 8 Color mapping range center and levels ccccssiicecassieivasiosedtsrwevsaci eecataveiereeeders 30 AO Prone Model TOONS rastaeenes cozeraitetidawan lt dtareieentavens a a a a 31 AIT POG DION aie a ten aseae ataaeestieeeeel 31 A 2 Prisi model Jeneral o sno i 32 AO Radiailamphitude Se C HUM siie a E E S JRE TOMA esn a E ee eee eee ere ere J 1 Column formatted data TES sccccansveincdaeuceussaunadsceted a E E IASC PATE ee ne ene Re te eR nen re ea es ony ta oe ey RC re Soe Reoulanly eridded Mati k TICS sce scr E a eewies ied enetuteien Japa Paramete S cca oieks ety n eine hares eaten pe donate ear aaeale ee Aces Addition inora Oie a glestausniante seasaecensanioel Mestaneniet oF sR CPO ICE oren T T S lt Terms Of USE ANd GISCLAIMER secssse
38. ls Range and Center scale widgets can be used to change the color mapping in contour and image maps Levels widget is used to define the number of contour levels default 1s 21 in contour maps or the refinement of the data pixels in image maps In the latter case the refinement varies from 1 original data pixels are shown to 10 depending on the value of Levels 5 50 If the number of data is large the refinement is reduced to draw the graphs faster By default the color scale e g rainbow or reverse rainbow is evenly distributed between the minimum and the maximum data values of the current map Range scale widget changes the width of the color scale between 5 and 125 of the original minimum and maximum values Decreasing the value is useful if the maximum data value of the map is larger than the mean data level due to outlier for example Increasing its value is useful if the maps get saturated at their minimum or maximum limit values Center scale widget changes the middle point of the color scale between 85 and 85 of the current data range Increasing the center value will emphasize small data values and increases the amount of colors at the beginning part of the color scale Decreasing the center value will emphasize large data values and increases the amount of colors at the end of the color scale 30 Note After changing any of the abovementioned scales the user must press the Update quick button to make the changes
39. menu items bring up the standard file selection dialog of the operating system that is used to provide the filename for I O operations All data files DAT ASC MAT are in text format see chapter 5 DISLIN cannot be used for direct output to a printer but the graphs can be saved into image files PS EPS PDF WMF PNG GIF for printing The graphs are saved as they appear on the screen in A4 size and landscape mode see Appendix A The preferred output format is PDF or PS because their resolution is much better than that of GIF and WMF files The PDF format is particularly useful because Adobe s Acrobat Reader can be used to take bitmap snapshots at desired resolution 3 2 Edit menu The Edit menu contains following items Inverse gt Data set current inverse results as new initial data grid Subtract inverse subtract the inverse results from the original data Add inverse add the inverse results to the original data change the color scale rainbow gray scale etc Color scale d Padding mode change the method used to pad and taper data Data type define type of potential field data gravity vs magnetic Miscellaneous a sub menu for some additional settings Lepegess 23 perform GUI guided or manual low pass filtering High pass filter perform GUI guided or manual high pass filtering Direction filter perform GUI guided or manual direction filtering Inverse gt data menu item makes the current inverse re
40. nd generalized derivative Sun shading and generalized derivative are Fourier operations that can be used to enhance small and linear features in the data This is accomplished by altering the direction of illumination or gradient direction defined by two angles elevation from the vertical axis zenith and azimuth from the positive x axis east They are defined by equations Cooper amp Cowan 2011 1 p sin8 q cos cos p sin Pr opta aad CDE J acot coset psie 19 J1 p2 q2 1 p q y where po cos tan g go sin tang p dU dx and q dU dy and and are angles that define the azimuth and elevation of the illumination or gradient direction and A x y is the total gradient defined in Eq 9 Examples of gravity sun shading and generalized derivative are shown in Appendix O Note The elevation and the azimuth angles in degrees of the sun shading and generalized derivative are read from the same two text field as the inclination and declination Elevation 0 180 is taken from the zenith out of the plane of the mapped data and azimuth 180 is taken from positive x axis east and it is positive in clock wise orientation See chapter 4 7 for more information about rotating the azimuth direction 17 Undo processing menu item can be used to revert back to the unprocessed status Possible low pass high pass or directional filtering operations will not be unmade however 3 4 Tools menu The Tools menu contains
41. ng the median of few selected data points instead of zeros b Linear extrapolation which will extend with the outmost data point c Mean value based which will use the nearest points inside some increasing search radius and d Gradient based which uses the mean value and the derivative of the field to extend the data to the padding zone The Shift yes no option is used to define if padding is added only to the top and to the right of original data area or if padding is made all around the data area which looks like the original data area were shifted Note The best results are usually obtained using gradient based padding with shift The other padding options are preserved primarily for testing and teaching purposes so that the user can see the effect and importance of padding and tapering 23 4 4 Fast Fourier transform After sampling interpolation and padding the Plot FFT button must be pressed to compute the FFT and to see the 2 D frequency spectrum see Appendix C When Plot FFT button is pressed if valid inverse results do not exist the automatic low pass filtering if activated is also made and the inverse FFT is computed and inverse results and the difference between original interpolated data and inverse transformed data are computed The true frequency spectrum F is complex F Re F i Im The FFT graph however shows the 10 base logarithm logio of the amplitude spectrum A F Re F Im F 7 a The
42. ntinuation is equivalent to low pass filtering since it removes high frequency contents of the data On the contrary downward continuation shifts the data below the plane of measurements e g inside the earth It is used to enhance the high frequency content of the data and to estimate the depth to the top of the targets The upward and downward continuation can be formulated as Blakely 1995 F U F U x e 47 and F U4 F U x et47 3 where F U F U and F U4 are the Fourier transforms of the potential field U upward continued field U and downward continued field U4 Az gt O is the elevation difference and k ky k 1 is the radial wave number The transformed fields are then computed by taking the inverse Fourier transform of F U and F U zg Downward continuation is not a stable operation If the plane of continuation is located below the actual potential field sources the results become erratic Laplace s equation does not hold anymore Examples of upward and downward continuation are shown in Appendix I 12 Note The height difference or elevation Az used in upward or downward continuation is read from the Height text field in the control panel left to the graph area The value of height difference is always positive 3 3 2 Pole and equator reduction Unlike the gravity field the static magnetic field of a symmetric body e g vertical prism exhibits non symmetric anomaly shape because of the inclined
43. o Teslas of the magnetic field 4 the density contrast in grams per cubic centimeter and 5 susceptibility dimensionless in SI units The extra parameters are reserved for future use e The following block of text lines define various text labels used in the graphs The maximum length of each text label is 70 characters e The last block of text lines define the labels of the frequency domain operations shown in the graph max 70 characters The format of the DIS file is likely to change in the future The character caret is used to define superscripts exponents the _ character underscore is used to generate subscripts and the character dollar is used to move the text back to the baseline for example Log_10 F 42 6 Additional information I started writing FOURPOT when I worked at the Geological Survey of Finland for the 3 D crustal model project funded by the Academy of Finland The original idea was to utilize the depth analysis methods of Ruotoistenm ki 1987 After I become a lecturer of applied geophysics at the University of Oulu in 2005 I started to modify FOURPOT for educational purposes Since then I ve used it in teaching gravity and magnetic data processing Every now and then I ve added new features and fixed bugs to make it more and more useful The Fourier transform is based on the FFT algorithm FORK by Jon Claerbout 1976 The bi linear interpolation algorithm is adapted from Numerical Reci
44. o perform the required interpolation The grid spacing dx and dy does not need to be equal in x and y directions After interpolation Make grid button becomes inactive and Plot data button is used to visualize the original automatically interpolated data To define new grid sampling one must first revert back to original automatically computed sampling using the Reset grid button edit new values for the Step x and Step y text fields and press the newly active Make grid button Normally Step x and Step y text fields show the number of nodes N and M Only when looking at the original data Plot data button and when Step distance option Edit Miscellaneous 1s active they can show the grid node step distances dx and dy Note The interpolation discussed above is performed using the DISLIN subroutine GETMAT see DISLIN user manual This subroutine works best if the data are already regularly sampled and the new grid coincides with the original grid If the grid sampling is too dense compared to data sampling the interpolation does not give good results and a blank white map may be shown In this case larger step values should be provided Alternatively highly irregular data should be interpolated on a regular grid using more advanced interpolation algorithm e g minimum curvature of third party software e g Golden Software Surfer 4 3 Padding and tapering Padding and tapering are essential parts of successful Fourier transfo
45. on which is related to the concept of analytic signal It is usually performed on magnetic data because the analytic signal does not depend on the direction of magnetization Hilbert transform is closely related to pole reduction although it is much easier to compute For a 3 D case analytic signal is defined as a vector field Blakely 1995 dU dU dU Ay yz nito tiaz 17 where the real and imaginary parts are a Hilbert transform pair The norm A IAI is equal to the total gradient in Eq 9 According to Blakely 1995 the Fourier transform of 1 D analytic signal can be computed as F A F U 1 sgn k 18 16 where sgn k 1 if kK lt O and sgn k 1 if k gt 0 In practice the 2 D Hilbert transform is accomplished by multiplying the spectrum with two and nulling all items in the lower left quarter of the centered amplitude spectrum kx lt 0 and ky lt 0 When initiated the program asks if the transform is made for 2 D map data or 1 D trace data In the latter case the program also asks if the transform is made along x axis half spectrum k lt 0 is nulled or y axis half spectrum ky lt 0 is nulled Please note that in FOURPOT the results are returned as the real part of the inverse FFT In case of Hilbert transform however the program will ask the user if the results are given as the absolute of the complex value which will yield the so called envelope function of the original function 3 3 7 Sun shading a
46. pectrum of low pass and high pass filtered data may look quite different from the spectrum of unfiltered data The horizontal axis starts from zero and ends either at the minimum Nyquist frequency or at the outer ring of the low pass filtering In other words fully low pass filtered part of the spectrum is not shown at all In this case a vertical dotted line indicates the position of inner cut off rings of low pass filter or high pass filter if applied 36 5 File formats Data can be read and saved in three file formats 1 column formatted data file DAT and 2 ArcInfo ASCII grid file format ASC and 3 gridded data matrix MAT When saving data into file the contents of the current graph will be saved Thus either a the original or reassigned data Inverse data or b the interpolated data together with the inverse transformed and processed data and their difference or c the padded and tapered data or d the centered Fourier transform will be saved In the b case the program asks if the data will be interpolated bilinear interpolation at original data locations or not 5 1 Column formatted data files The format of a DAT file is illustrated below O00 0200 H 866SD06RFOL 0233469 LOe00 O00 OseS5e65LE701L OAAS ZU 200 D300 SOLEO OS ETOL GGA 74 The first line defines the index numbers of the columns that refer to the x and y coordinates ICO1 ICO2 and the data ICO3 The file can contain mult
47. pes Press et al 1988 Most of the Fourier operations have been documented in Blakely 1995 and Cooper and Cowan 2011 The principal directions are computed using RIDGEORIENT algorithm for fingerprint detection by Raymond Thai amp Peter Kovesi http www csse uwa edu au pk research matlabfns FingerPrints ridgeorient m Thanks to Richard Stuart for his comments on smooth frequency filtering FOURPOT is written in Fortran 90 and compiled with Intel Visual Fortran 15 The graphical user interface is based on the DISLIN graphics library version 10 4 http www dislin de by Helmut Michels In principle FOURPOT could be compiled and run on other operating systems Linux Mac OSX without major modifications provided that DISLIN is also installed At the moment however the source code is not made available and I do not provide active support for the software Nonetheless if you find the results erroneous or if you have suggestions for improvements please inform me 43 7 References Bhattacharyya B K 1967 Some general properties of potential fields in space and frequency domain a review Geoexploration 5 127 143 Blakely R J 1995 Potential theory in gravity and magnetic applications Cambridge Univ Press 441 p Claerbout J 1976 Fundamentals of geophysical data processing With applications to petroleum prospecting McGraw Hill Book Co Cooper G R J amp Cowan D R 2006 Enhancing potential field data using
48. rm processing The Padding button automatically adds extra columns and rows around the interpolated data matrix so that the grid dimensions N and M will become even powers of two e g 64 128 256 512 etc that is required by the FFT algorithm The padding can be performed with or without tapering Tapering means that artificial values are added to the padded area so that the continuity of the data is preserved at the edges See Fig 2 and Appendix B for examples of padding a normal padding b shifted padding b padding tapering padded values tapered values padded values original data original data area area original data area padded values tapered values 0 0 N N 00 N 00 N Figure 2 Schematic view of a padding without shift b padding with shift and c shifted padding with tapering The new dimensions N and M are powers of two Sometimes the original grid dimension may already be so close to some power of two that the tapering will not ensure enough data continuity To further increase the padded dimension from 64 to 128 for example one needs to change the current values of Step X or Step Y to some number between 65 and 128 and then press the Padding button again The dimensions of the grid can be decreased only using re interpolation Make grid button Padding and tapering depend on the selection made with the Padding mode item in Edit menu The padding modes are a Zeros which actually uses for paddi
49. rrupted FOURPOT is likely to crash when started In this case one should delete the file and a new one with default values will be generated automatically next time the program is started The file format is illustrated below Fourpot ver 1 30 parameter file 36 oe 20 O 0 1 0 2 1 1 il il 0 0 SoU BOO Oe90 O90 Ose O00 Gago U0 74 7 Sef DILS67 U U05000 U0 LOGO OO USU 2 D Fourier transform analysis Original data Interpolated data Padded data Fourier transform Inverse data Difference Radial ave spectre Ky Log 10 F Kr Log F leg Fy C Magnetic data Gravity data Undefined data Dimension Spacing Nyquist Wave length Rad spectre depths Tensor LOtatLon Rotation Declin rotation Elev amp Azimuth Distance Depth Low Pass filtered High Pass filtered Direction filtered Upward continued Downward continued Reduced to pole 40 dF dz gradient dF dx gradient dF dy gradient a 25F dz 2S gradient AS 25F dx 2S gradient G Z2SP ay 25 Gradient d 2SF dxdy gradient d 2SF dydz gradient d 2SF dydz gradient Horizontal gradient Total gradient Tilt gradient Pseudo gravity Pseudo magnetic Gravity potential Magnetic potential Theta map ROts invariant l Rots invariant 2 Reduced to equator d ASP Lay 29F dn 2n EJ dR 2p 72 Maxe hOra zontal gradient Hilbert transform Sun shading General derivative Unprocessed data e The first line is used as a header to identif
50. s floating point values must be given 10 The remaining items in the Edit menu are related to frequency filtering operations that can be performed either interactively with the mouse or manually by providing the numeric values via GUI Low pass and high pass filtering are the two most typical operations Considering 2 D FFT low pass filtering removes nullifies values outside certain radius from the origin of the amplitude spectrum see Appendices C E and F On the contrary high pass filtering removes low frequency data around and near to the origin of the spectrum Combined use of both low and high pass filtering is equivalent to band pass filtering The direction filtering can be used to remove linear features in the original data see Appendices G and H The frequency filtering will be discussed more in the chapter 4 6 3 3 Process menu The Process menu contains the items for inverse processing tasks Most of these operations can also be accessed using a quick tool button in the control panel next to the graph area Upward continuation performs upward continuation of potential field data Downward continuation performs downward continuation Pole reduction performs reduction to magnetic pole Equator reduction performs reduction to magnetic equator vertical gradient of the data dF dz gradient dF dx gradient horizontal x gradient of the data dF dy gradient horizontal y gradient of the data aaa second vertical gradien
51. s image maps 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Difference 1830 1830 1820 1810 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 2880 X Dimension 52x 48 X Dimension 52x 48 Spacing 1 000000 x 1 000000 m Spacing 1 000000 x 1 000000 m Low Pass filtered 4 000000 2 998501 Low Pass filtered 4 000000 2 998501 Direction filtered 15 52443 15 33253 Direction filtered 15 52443 15 33253 Unprocessed data Difference Field 49 Appendix I Upward and downward continued data h 2 km 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X 27 11 29 47 Field 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X aiti 227 56 139 Field Appendix J First vertical and horizontal x gradient 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 4 15 5 4 14 Field 2 D Fourier transform analysis Inverse data 1830 1820 gt 1810 1800 1790 2840 2850 2860 2870 2880 X 13 5 2 10 Field 50 Appendix K Second vertical and horizontal xy gradient 2 D Fourier transform analysis 2 D Fourier transform analysis Inverse data Inverse data 1830 1830 1820 1820 gt 1810 gt 1810 1800 1800 1790 1790 2840 2850 2860 2870 2880 2840 2850 2860 2870 28
52. s of use data or profits arising out of or in connection with the use or performance of this software 45 Appendix A FOURPOT GUI after data input and default interpolation _ w 2 D Fourier analysis version 1 3a by M Pirttij rvi 2014 SE File Edit Process Tools Exit Plot data Shortcuts Reset grid __Upward 2 D Fourier transform analysis T _Dewriwa Interpolated data Padding pa Plot inverse Zm p 1830 a r sovm p Directions G dds Height units ooo e d dud Incl Elev deg 23 i ded 1820 Decl Azim deg 24s ded Rotation angle deg 45 0 E ansparent ng a ea Colors Macr 1810 Range tat ET Lj Tej rer 1 00 Cae Center xs m maj o Levels lod 1800 2840 2850 2860 2870 2880 Dimension 52x 48 Spacing 1 000000 x 1 000000 km l Field Appendix B Data map after padding with shift and tapering The dotted rectangle depicts the area of the original data 2 D Fourier transform analysis Interpolated data 1840 1830 1820 gt 1810 1800 1790 1780 2840 2850 2860 2870 2880 2890 Generic data Dimension 64x 64 Spacing 1 000x 1 000 km 46 Appendix C 2 D frequency amplitude spectrum 2 D Fourier transform analysis Fourier transform 0 4 0 3 0 2 0 1 0 1 0 4 0 2 0 0 0 2 0 4 Dimension 64x 64 1562500E 01 1 m Spacing 1562500E 01 x Wavelength 2 000000
53. scrcseenntenddoee tan teen siadensaebeddans E E T SSS Appendices A R Keywords potential fields Fourier transformation frequency filtering data processing List of Appendices Appendix A FOURPOT GUI after data input and default interpolation Appendix B Data map after padding with shift and tapering Appendix C 2 D frequency amplitude spectrum Appendix D Difference between original data and inverse transformed data Appendix E Low pass filtered frequency spectrum Appendix F Low pass filtered data and difference with the original data Appendix G Direction filtered frequency spectrum Appendix H Direction filtered data and difference with the original data Appendix I Upward and downward continued data Appendix J First vertical and horizontal x gradient Appendix K Second vertical and horizontal xy gradient Appendix L Horizontal gradient and total gradient Appendix M Tilt derivative and theta map Appendix N Gravity potential and pseudo magnetic field Appendix O Sun shading and general derivative Appendix P Transparency amp sun shading and general derivative Appendix Q Sun shading shown transparently on the low pass filtered data Appendix R Ridge picking and principal directions 1 Introduction FOURPOT is a computer program made for processing and analysis of two dimensional 2 D potential field data arising from geophysical gravity and magnetic field measurements in particular The data can be regularly or
54. sults the new interpolated input data for successive frequency domain operations For example tilt gradient can be computed after upward continuation or reduction to the pole without saving the UC or RTP results in a new data file Subtract inverse item allows computing the residual field as the difference between the original interpolated data and the upward continued or low pass filtered data Similarly Add inverse item can be used to add the current inverse results to the original data Although this item is not often used it could be used to enhance high frequency content of the data for example The items in Color scale sub menu are quite self explaining they change the color scale of the image and contour maps Comparable colors item imposes the same limit values on the color scales of the original and the processed data and thus makes their comparison easier Padding mode sub menu is used to demonstrate different padding and tapering modes Padding is an operation where artificial or null values are added around the original data area see Appendix B The purpose of padding is to extend the data N and M to an even power of two e g 64 128 256 512 etc as required by the FFT algorithm Tapering is an operation where the padded values are made such that they prevent rapid amplitude changes at the borders of the data area This means that the data and their derivative are made more or less continuous across the original
55. t of the data d F dz gradient second horizontal x gradient of the data d F dx gradient a second horizontal y gradient of the data d F dy gradient Fe fdidy ieakent xy gradient of the data dF dxdz gradient XZ ST adient of the data d F dydz gradient yz gradient of the data d F dx d F dy 2 tensor component measured by Falcon gravity gradient system L n first rotational invariant gravity tensor Rot invariant 2 second rotational invariant gravity tensor 11 Hori gradient horizontal gradient Total gradient total gradient analytic signal Tilt gradient tilt gradient Max hori grad Tose maximum horizontal gradient amplitude EE pseudo gravity field from magnetic data Pseudo gravity pseudo magnetic field from gravity data imine mee gravity potential from gravity field Bouguer data Gravity potential magnetic potential from total field TMI data Magnetic potential Se Hilbert transform phase shift Sun shading sun shading to be used with transparency General derivative generalized derivative operator from G F Cooper i reverts back to unprocessed but frequency filtered data Undo processing 3 3 1 Upward and downward continuation Upward continuation is an operation that shifts the data by a constant height level above the surface of the earth or the plane of measurements It is used to estimate the large scale or regional low frequency or long wave length trends of the data To some degree upward co
56. tation of the direction of sun shading and generalized derivative Resetting makes the azimuth zero along x axis Most of the other gradient operators as well as sun shading and generalized derivative are more informative when displayed using greyscale colors Activating the Transparent shading check box allows displaying the inverse results as greyscale colors transparently below the original unprocessed but gridded data that uses the current user selected color scale When 29 activated the program asks the user for the strength of transparency Transparency ranges between O fully transparent and 1 fully opaque Negative value of transparency will apply inverse greyscale Examples of transparent color shading are shown in Appendix P for sun shading and generalized derivative Together with transparent shading the rotation can be used to find the direction or directions that enhance wanted features in the data Remember however that some gradient components and derived gradients are invariant to horizontal rotation Also note that some inverse results pole and equator reduction potentials and pseudo fields and Hilbert transform are not shown above the original data transparently Instead the inverse results are shown transparently on top of themselves Although this may sound unnecessary it often improves the visibility of data maxima or minima depending on the current color scale 4 8 Color mapping range center and leve
57. tered data are always passed to further processing operations without the use of nverse data menu item As a matter of fact low pass filtering is often an essential preliminary step for successful operation of the gradient filters Frequency filtering means that some parts of the spectrum are nulled Rapid changes in the spectrum however can cause oscillations in the inverse transformed data This is known as the Gibbs phenomenon Therefore instead of using rectangular or box car shaped filter functions bell shaped smoothing is applied over the cut off range of the filter In low pass and high pass filtering this means in practice that instead of a single filter ring the user provides the radii of the inner and the outer ring of the cut off range Between the two rings half sine or half cosine function is used to define the smoothing In directional filtering the user provides the direction angles of the upper and lower limit of the notch filters and the width angle of the cut off range The filter rings and sectors and the cut off ranges are illustrated in Fig 3 25 a 2 D low pass filter b 2 D directional filter mailed cut off range as cut off values ranges preserved values upper and lower sector range 228 d high pass filter e band pass filter o w Ly Pi x See od k G ess PELI Pid MOTELI Figure 3 Schematic view of a cut off rings and b sectors of 2 D filters an
58. tion is a threshold value which can be used to ignore small anomalies The threshold value depends on the tilt gradient data but usually the default value of 0 0 is a safe choice Reducing the threshold below zero allows smaller anomalies to be added to model generation Vice versa increasing the value rejects small anomalies and builds model only below the strongest anomalies When the model has been generated a vertical cross section of the model will appear below the response graph see Fig 5 The model generation uses the minima around the positive peaks of the tilt gradient data to define the sub areas for the inversion dotted lines through the top of the prisms in Fig 4 The width of the sub areas is then used to define the initial value for the thickness of the prism models and the anti symmetry of the sub areas with respect to the center of the bodies is used to define the initial dip angle The strike extent width of the prism models in a direction perpendicular to the profile and depth extent are based on length of the profile The density or magnetic susceptibility is based on the amplitude of the potential field data at the center location The model can be saved into MAGBODS or GRABODS input file format using File Save model file menu item For this purpose the data type needs to be set accordingly using Edit Data type sub menu The user should also take care to define the spatial dimension correctly using Miscellaneous Meters
59. utational artifacts in the inverse transformed data it is often necessary to perform low pass filtering before any other processing operations Automatic LP filter item enables disables automatic low pass filtering To activate automatic low pass filtering the user needs to provide the radii of the filter rings manually In this case the inner and outer radii are defined as a frequencies normalized by the maximum Nyquist frequency 1 e a value between 0 and 1 in both x and y directions 1 e not as wave lengths Thus low pass filter with default values 0 5 and 0 67 removes frequencies greater than two thirds of Nyquist frequency 0 67 0 5 dx or wave lengths less than 3 dx also in y direction if dx lt dy Default values are used if both values are zero Note Automatic low pass filtering is disabled if the input values are omitted totally 1 e a blank line is given as input 28 4 7 Rotation and transparent shading Horizontal rotation around z axis can be applied to basic gradient components other than first and second vertical derivative as well as reduction to equator sun shading and generalized derivative operations to enhance linear features of varying direction in the potential field data The Activate rotation check box must first be activated and a non zero value provided for the Rotation angle text field The Rotate button is then used to make one rotation step which 1s added to the current value of total rotation and the
60. y the file and version number e The 1 st relevant line defines three character height values 1 main title 2 axes labels and 3 the plot information text The other two values are reserved for future use e The 2 nd line defines some integer valued options 1 the screen mode 0 1 between the normal 4 3 aspect ratio and wide screen mode 2 the color scale 0 6 3 padding mode 0 3 4 tapering shift modes 0 1 5 contour vs image map 0 1 6 automatic low pass filtering 0 1 and 7 step number vs step distance 0 1 The extra parameters are reserved for future use e The 3 rd line defines the x horizontal and y vertical position of the origin of the main graph in pixels from the bottom left corner of the page integers and the length of the x and y axis relative to the size of the origin subtracted width and height of the plot area The total size of the plot area is always 2970x2100 pixels landscape A4 The fifth parameter defines the aspect ratio of the graph area for widescreen mode Aspect ratio 1 refers to normal 4 3 screen and values less to zero indicate widescreen displays The remaining parameters are reserved for future use e The 4 th line defines the initial values for the magnetic and gravity field components used in pole reduction and pseudo field computations These parameters include 1 the inclination and 2 declination in degrees from horizontal plane and true north and 3 41 intensity in nan

Download Pdf Manuals

image

Related Search

Related Contents

Bedienungsanleitung  TAFCO WINDOWS NU2-150S-W Installation Guide  INSTALLATION MANUAL MANUAL DE INSTALACIÓN MANUEL D  the attached document  Samsung WA9011B0 คู่มือการใช้งาน  ASUS G10AC PL8185 User's Manual  環境報告書 2001  

Copyright © All rights reserved.
Failed to retrieve file