1GRAVFFT(1) GMT GRAVFFT(1)
2
6 gravfft - Compute gravitational attraction of 3-D surfaces in the
7 wavenumber (or frequency) domain
8
10 gravfft ingrid [ ingrid2 ] -Goutfile [ -Cn/wavelength/mean_depth/tbw
11 ] [ -Ddensity|rhogrid ] [ -En_terms ] [ -F[f[+]|g|v|n|e] ] [
12 -Iw|b|c|t |k ] [ -Nparams ] [ -Q ] [ -Tte/rl/rm/rw[/ri][+m] ] [
13 -V[level] ] [ -Wwd] [ -Zzm[zl] ] [ -fg ]
14 Note: No space is allowed between the option flag and the associated
16 arguments.
17
19 gravfft can be used into three main modes. Mode 1: Simply compute the
20 geopotential due to the surface given in the topo.grd file. Requires a
21 density contrast (-D) and possibly a different observation level (-W).
22 It will take the 2-D forward FFT of the grid and use the full Parker's
23 method up to the chosen terms. Mode 2: Compute the geopotential
24 response due to flexure of the topography file. It will take the 2-D
25 forward FFT of the grid and use the full Parker's method applied to
26 the chosen isostatic model. The available models are the "loading from
27 top", or elastic plate model, and the "loading from below" which
28 accounts for the plate's response to a sub-surface load (appropriate
29 for hot spot modeling - if you believe them). In both cases, the model
30 parameters are set with -T and -Z options. Mode 3: compute the admit‐
31 tance or coherence between two grids. The output is the average in the
32 radial direction. Optionally, the model admittance may also be calcu‐
33 lated. The horizontal dimensions of the grdfiles are assumed to be in
34 meters. Geographical grids may be used by specifying the -fg option
35 that scales degrees to meters. If you have grids with dimensions in km,
36 you could change this to meters using grdedit or scale the output with
37 grdmath. Given the number of choices this program offers, is difficult
38 to state what are options and what are required arguments. It depends
39 on what you are doing; see the examples for further guidance.
40
42 ingrid 2-D binary grid file to be operated on. (See GRID FILE FORMATS
43 below). For cross-spectral operations, also give the second
44 grid file ingrd2.
45 -Goutfile
47 Specify the name of the output grid file or the 1-D spectrum ta‐
48 ble (see -E). (See GRID FILE FORMATS below).
49
51 -Cn/wavelength/mean_depth/tbw
52 Compute only the theoretical admittance curves of the selected
53 model and exit. n and wavelength are used to compute (n * wave‐
54 length) the total profile length in meters. mean_depth is the
55 mean water depth. Append dataflags (one or two) of tbw in any
56 order. t = use "from top" model, b = use "from below" model.
57 Optionally specify w to write wavelength instead of frequency.
58 -Ddensity|rhogrid
60 Sets density contrast across surface. Used, for example, to com‐
61 pute the gravity attraction of the water layer that can later be
62 combined with the free-air anomaly to get the Bouguer anomaly.
63 In this case do not use -T. It also implicitly sets -N+h.
64 Alternatively, specify a co-registered grid with density con‐
65 trasts if a variable density contrast is required.
66 -En_terms
68 Number of terms used in Parker expansion (limit is 10, otherwise
69 terms depending on n will blow out the program) [Default = 3]
70 -F[f[+]|g|v|n|e]
72 Specify desired geopotential field: compute geoid rather than
73 gravity
74 f = Free-air anomalies (mGal) [Default]. Append + to add in
75 the slab implied when removing the mean value from the topog‐
76 raphy. This requires zero topography to mean no mass anom‐
77 aly.
78 g = Geoid anomalies (m).
80 v = Vertical Gravity Gradient (VGG; 1 Eotvos = 0.1 mGal/km).
82 e = East deflections of the vertical (micro-radian).
84 n = North deflections of the vertical (micro-radian).
86 -Iw|b|c|t |k
88 Use ingrd2 and ingrd1 (a grid with topography/bathymetry) to
89 estimate admittance|coherence and write it to stdout (-G ignored
90 if set). This grid should contain gravity or geoid for the same
91 region of ingrd1. Default computes admittance. Output contains 3
92 or 4 columns. Frequency (wavelength), admittance (coherence) one
93 sigma error bar and, optionally, a theoretical admittance.
94 Append dataflags (one to three) from w|b|c|t. w writes wave‐
95 length instead of wavenumber, k selects km for wavelength unit
96 [m], c computes coherence instead of admittance, b writes a
97 fourth column with "loading from below" theoretical admittance,
98 and t writes a fourth column with "elastic plate" theoretical
99 admittance.
100 -N[a|f|m|r|s|nx/ny][+a|[+d|h|l][+e|n|m][+twidth][+v][+w[suffix]][+z[p]]
102 Choose or inquire about suitable grid dimensions for FFT and set
103 optional parameters. Control the FFT dimension:
104 -Na lets the FFT select dimensions yielding the most accurate
105 result.
106 -Nf will force the FFT to use the actual dimensions of the
108 data.
109 -Nm lets the FFT select dimensions using the least work mem‐
111 ory.
112 -Nr lets the FFT select dimensions yielding the most rapid
114 calculation.
115 -Ns will present a list of optional dimensions, then exit.
117 -Nnx/ny will do FFT on array size nx/ny (must be >= grid file
119 size). Default chooses dimensions >= data which optimize
120 speed and accuracy of FFT. If FFT dimensions > grid file
121 dimensions, data are extended and tapered to zero.
122 Control detrending of data: Append modifiers for removing a lin‐
124 ear trend:
125 +d: Detrend data, i.e. remove best-fitting linear trend
126 [Default].
127 +a: Only remove mean value.
129 +h: Only remove mid value, i.e. 0.5 * (max + min).
131 +l: Leave data alone.
133 Control extension and tapering of data: Use modifiers to control
135 how the extension and tapering are to be performed:
136 +e extends the grid by imposing edge-point symmetry
137 [Default],
138 +m extends the grid by imposing edge mirror symmetry
140 +n turns off data extension.
142 Tapering is performed from the data edge to the FFT grid edge
144 [100%]. Change this percentage via +twidth. When +n is in
145 effect, the tapering is applied instead to the data margins
146 as no extension is available [0%].
147 Control messages being reported: +v will report suitable
149 dimensions during processing.
150 Control writing of temporary results: For detailed investigation
152 you can write the intermediate grid being passed to the forward
153 FFT; this is likely to have been detrended, extended by
154 point-symmetry along all edges, and tapered. Append +w[suffix]
155 from which output file name(s) will be created (i.e.,
156 ingrid_prefix.ext) [tapered], where ext is your file extension.
157 Finally, you may save the complex grid produced by the forward
158 FFT by appending +z. By default we write the real and imaginary
159 components to ingrid_real.ext and ingrid_imag.ext. Append p to
160 save instead the polar form of magnitude and phase to files
161 ingrid_mag.ext and ingrid_phase.ext.
162 -Q Writes out a grid with the flexural topography (with z positive
164 up) whose average depth was set by -Zzm and model parameters by
165 -T (and output by -G). That is the "gravimetric Moho". -Q
166 implicitly sets -N+h
167 -S Computes predicted gravity or geoid grid due to a subplate load
169 produced by the current bathymetry and the theoretical model.
170 The necessary parameters are set within -T and -Z options. The
171 number of powers in Parker expansion is restricted to 1. See an
172 example further down.
173 -Tte/rl/rm/rw[/ri][+m]
175 Compute the isostatic compensation from the topography load
176 (input grid file) on an elastic plate of thickness te. Also
177 append densities for load, mantle, water and infill in SI units.
178 If ri is not provided it defaults to rl. Give average mantle
179 depth via -Z. If the elastic thickness is > 1e10 it will be
180 interpreted as the flexural rigidity (by default it is computed
181 from te and Young modulus). Optionally, append +m to write a
182 grid with the Moho's geopotential effect (see -F) from model
183 selected by -T. If te = 0 then the Airy response is returned.
184 -T+m implicitly sets -N+h
185 -Wwd Set water depth (or observation height) relative to topography
187 [0]. Append k to indicate km.
188 -Zzm[zl]
190 Moho [and swell] average compensation depths (in meters positive
191 dows -- the depth). For the "load from top" model you only have
192 to provide zm, but for the "loading from below" don't forget zl.
193 -V[level] (more ...)
195 Select verbosity level [c].
196 -fg Geographic grids (dimensions of longitude, latitude) will be
198 converted to meters via a "Flat Earth" approximation using the
199 current ellipsoid parameters.
200 -^ or just -
202 Print a short message about the syntax of the command, then
203 exits (NOTE: on Windows just use -).
204 -+ or just +
206 Print an extensive usage (help) message, including the explana‐
207 tion of any module-specific option (but not the GMT common
208 options), then exits.
209 -? or no arguments
211 Print a complete usage (help) message, including the explanation
212 of all options, then exits.
213
215 By default GMT writes out grid as single precision floats in a
216 COARDS-complaint netCDF file format. However, GMT is able to produce
217 grid files in many other commonly used grid file formats and also
218 facilitates so called "packing" of grids, writing out floating point
219 data as 1- or 2-byte integers. (more ...)
220
222 If the grid does not have meter as the horizontal unit, append +uunit
223 to the input file name to convert from the specified unit to meter. If
224 your grid is geographic, convert distances to meters by supplying -fg
225 instead.
226
228 netCDF COARDS grids will automatically be recognized as geographic. For
229 other grids geographical grids were you want to convert degrees into
230 meters, select -fg. If the data are close to either pole, you should
231 consider projecting the grid file onto a rectangular coordinate system
232 using grdproject.
233
235 The FFT solution to elastic plate flexure requires the infill density
236 to equal the load density. This is typically only true directly
237 beneath the load; beyond the load the infill tends to be lower-density
238 sediments or even water (or air). Wessel [2001] proposed an approxima‐
239 tion that allows for the specification of an infill density different
240 from the load density while still allowing for an FFT solution. Basi‐
241 cally, the plate flexure is solved for using the infill density as the
242 effective load density but the amplitudes are adjusted by a factor A =
243 sqrt ((rm - ri)/(rm - rl)), which is the theoretical difference in
244 amplitude due to a point load using the two different load densities.
245 The approximation is very good but breaks down for large loads on weak
246 plates, a fairy uncommon situation.
247
249 To compute the effect of the water layer above the bat.grd bathymetry
250 using 2700 and 1035 for the densities of crust and water and writing
251 the result on water_g.grd (computing up to the fourth power of bathyme‐
252 try in Parker expansion):
253 gmt gravfft bat.grd -D1665 -Gwater_g.grd -E4
255 Now subtract it from your free-air anomaly faa.grd and you will get the
257 Bouguer anomaly. You may wonder why we are subtracting and not adding.
258 After all the Bouguer anomaly pretends to correct the mass deficiency
259 presented by the water layer, so we should add because water is less
260 dense than the rocks below. The answer relies on the way gravity
261 effects are computed by the Parker's method and practical aspects of
262 using the FFT.
263 gmt grdmath faa.grd water_g.grd SUB = bouguer.grd
265 Want an MBA anomaly? Well compute the crust mantle contribution and add
267 it to the sea-bottom anomaly. Assuming a 6 km thick crust of density
268 2700 and a mantle with 3300 density we could repeat the command used to
269 compute the water layer anomaly, using 600 (3300 - 2700) as the density
270 contrast. But we now have a problem because we need to know the mean
271 Moho depth. That is when the scale/offset that can be appended to the
272 grid's name comes in hand. Notice that we didn't need to do that before
273 because mean water depth was computed directly from data (notice also
274 the negative sign of the offset due to the fact that z is positive up):
275 gmt gravfft bat.grd=nf/1/-6000 -D600 -Gmoho_g.grd
277 Now, subtract it from the Bouguer to obtain the MBA anomaly. That is:
279 gmt grdmath bouguer.grd moho_g.grd SUB = mba.grd
281 To compute the Moho gravity effect of an elastic plate bat.grd with Te
283 = 7 km, density of 2700, over a mantle of density 3300, at an average
284 depth of 9 km
285 gmt gravfft bat.grd -Gelastic.grd -T7000/2700/3300/1035+m -Z9000
287 If you add now the sea-bottom and Moho's effects, you will get the full
289 gravity response of your isostatic model. We will use here only the
290 first term in Parker expansion.
291 gmt gravfft bat.grd -D1665 -Gwater_g.grd -E1
293 gmt gravfft bat.grd -Gelastic.grd -T7000/2700/3300/1035+m -Z9000 -E1
294 gmt grdmath water_g.grd elastic.grd ADD = model.grd
295 The same result can be obtained directly by the next command. However,
297 PAY ATTENTION to the following. I don't yet know if it's because of a
298 bug or due to some limitation, but the fact is that the following and
299 the previous commands only give the same result if -E1 is used. For
300 higher powers of bathymetry in Parker expansion, only the above example
301 seams to give the correct result.
302 gmt gravfft bat.grd -Gmodel.grd -T7000/2700/3300/1035 -Z9000 -E1
304 And what would be the geoid anomaly produced by a load at 50 km depth,
306 below a region whose bathymetry is given by bat.grd, a Moho at 9 km
307 depth and the same densities as before?
308 gmt gravfft topo.grd -Gswell_geoid.grd -T7000/2700/3300/1035 -Fg -Z9000/50000 -S -E1
310 To compute the admittance between the topo.grd bathymetry and faa.grd
312 free-air anomaly grid using the elastic plate model of a crust of 6 km
313 mean thickness with 10 km effective elastic thickness in a region of 3
314 km mean water depth:
315 gmt gravfft topo.grd faa.grd -It -T10000/2700/3300/1035 -Z9000
317 To compute the admittance between the topo.grd bathymetry and geoid.grd
319 geoid grid with the "loading from below" (LFB) model with the same as
320 above and sub-surface load at 40 km, but assuming now the grids are in
321 geographic and we want wavelengths instead of frequency:
322 gmt gravfft topo.grd geoid.grd -Ibw -T10000/2700/3300/1035 -Z9000/40000 -fg
324 To compute the gravity theoretical admittance of a LFB along a 2000 km
326 long profile using the same parameters as above
327 gmt gravfft -C400/5000/3000/b -T10000/2700/3300/1035 -Z9000/40000
329
331 Luis, J.F. and M.C. Neves. 2006, The isostatic compensation of the
332 Azores Plateau: a 3D admittance and coherence analysis. J. Geothermal
333 Volc. Res. Volume 156, Issues 1-2, Pages 10-22,
334 http://dx.doi.org/10.1016/j.jvolgeores.2006.03.010
335 Parker, R. L., 1972, The rapid calculation of potential anomalies, Geo‐
337 phys. J., 31, 447-455.
338 Wessel. P., 2001, Global distribution of seamounts inferred from grid‐
340 ded Geosat/ERS-1 altimetry, J. Geophys. Res., 106(B9), 19,431-19,441,
341 http://dx.doi.org/10.1029/2000JB000083
342
344 gmt, grdfft, grdmath, grdproject
345
347 2019, P. Wessel, W. H. F. Smith, R. Scharroo, J. Luis, and F. Wobbe
3485.4.5 Feb 24, 2019 GRAVFFT(1)