From 9c237c64156bf8adc3b66e8413ce0d50d124e606 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 09:52:45 +0100 Subject: [PATCH 01/10] use balance for all diverging colormaps, update cpt limits and colorbar settings --- .../block_averaged_sources.py | 18 ++--- .../equivalent_sources/cartesian.py | 17 +++-- .../equivalent_sources/gradient_boosted.py | 16 ++-- .../equivalent_sources/spherical.py | 18 +++-- doc/gallery_src/forward/point_gravity.py | 6 +- doc/gallery_src/forward/prism_layer.py | 2 +- .../forward/prisms_topo_gravity.py | 10 ++- doc/gallery_src/forward/tesseroid_layer.py | 10 ++- doc/gallery_src/gravity_disturbance.py | 9 ++- .../gravity_disturbance_topofree.py | 9 ++- doc/gallery_src/isostatic_moho_airy.py | 7 +- .../transformations/reduction_to_pole.py | 10 +-- doc/gallery_src/transformations/tga.py | 12 +-- doc/gallery_src/transformations/tilt.py | 16 ++-- .../transformations/upward_continuation.py | 9 ++- .../transformations/upward_derivative.py | 8 +- .../equivalent_sources/block-averaged-eqs.rst | 4 +- .../eq-sources-spherical.rst | 2 +- .../eqs-parameters-estimation.rst | 2 +- .../gradient-boosted-eqs.rst | 2 +- doc/user_guide/equivalent_sources/index.rst | 4 +- doc/user_guide/forward_modelling/point.rst | 6 +- .../forward_modelling/tesseroid.rst | 1 + doc/user_guide/gravity_disturbance.rst | 2 +- doc/user_guide/igrf.rst | 76 ++++++++++--------- doc/user_guide/topographic_correction.rst | 6 +- 26 files changed, 161 insertions(+), 121 deletions(-) diff --git a/doc/gallery_src/equivalent_sources/block_averaged_sources.py b/doc/gallery_src/equivalent_sources/block_averaged_sources.py index 0a50e41aa..150630a39 100644 --- a/doc/gallery_src/equivalent_sources/block_averaged_sources.py +++ b/doc/gallery_src/equivalent_sources/block_averaged_sources.py @@ -41,6 +41,7 @@ """ import ensaio +import numpy as np import pandas as pd import pygmt import pyproj @@ -106,14 +107,14 @@ title = "Observed magnetic anomaly data" +# Get the 99.9th percentile of the absolute value of the point data to use as color +# scale limits +cpt_lim = np.quantile(np.abs(data.total_field_anomaly_nt), 0.999) + # Make colormap of data -# Get the 95 percentile of the maximum absolute value between the original and -# gridded data so we can use the same color scale for both plots and have 0 -# centered at the white color. -maxabs = vd.maxabs(data.total_field_anomaly_nt, grid.magnetic_anomaly.values) * 0.95 pygmt.makecpt( - cmap="vik", - series=(-maxabs, maxabs), + cmap="balance+h0", + series=[-cpt_lim, cpt_lim], background=True, ) @@ -129,8 +130,7 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a400f100", "x+lnT"]) - +fig.colorbar(cmap=True, frame=["x+lnT"], position="+e") fig.shift_origin(xshift=fig_width + 1) title = "Gridded and upward-continued" @@ -142,6 +142,6 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a400f100", "x+lnT"]) +fig.colorbar(cmap=True, frame=["x+lnT"], position="+e") fig.show() diff --git a/doc/gallery_src/equivalent_sources/cartesian.py b/doc/gallery_src/equivalent_sources/cartesian.py index ad67084df..c789ffcc1 100644 --- a/doc/gallery_src/equivalent_sources/cartesian.py +++ b/doc/gallery_src/equivalent_sources/cartesian.py @@ -31,6 +31,7 @@ """ import ensaio +import numpy as np import pandas as pd import pygmt import pyproj @@ -95,14 +96,14 @@ title = "Observed magnetic anomaly data" +# Get the 99.9th percentile of the absolute value of the point data to use as color +# scale limits +cpt_lim = np.quantile(np.abs(data.total_field_anomaly_nt), 0.999) + # Make colormap of data -# Get the 95 percentile of the maximum absolute value between the original and -# gridded data so we can use the same color scale for both plots and have 0 -# centered at the white color. -maxabs = vd.maxabs(data.total_field_anomaly_nt, grid.magnetic_anomaly.values) * 0.95 pygmt.makecpt( - cmap="vik", - series=(-maxabs, maxabs), + cmap="balance+h0", + series=[-cpt_lim, cpt_lim], background=True, ) @@ -118,7 +119,7 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a400f100", "x+lnT"]) +fig.colorbar(cmap=True, frame=["x+lnT"], position="+e") fig.shift_origin(xshift=fig_width + 1) @@ -131,6 +132,6 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a400f100", "x+lnT"]) +fig.colorbar(cmap=True, frame=["x+lnT"], position="+e") fig.show() diff --git a/doc/gallery_src/equivalent_sources/gradient_boosted.py b/doc/gallery_src/equivalent_sources/gradient_boosted.py index 72b5881f6..1a1260d46 100644 --- a/doc/gallery_src/equivalent_sources/gradient_boosted.py +++ b/doc/gallery_src/equivalent_sources/gradient_boosted.py @@ -31,6 +31,7 @@ import boule as bl import ensaio +import numpy as np import pandas as pd import pygmt import pyproj @@ -115,13 +116,14 @@ title = "Observed gravity disturbance data" +# Get the 99.9th percentile of the absolute value of the point data to use as color +# scale limits +cpt_lim = np.quantile(np.abs(data.gravity_disturbance), 0.999) + # Make colormap of data pygmt.makecpt( - cmap="vik", - series=( - -data.gravity_disturbance.quantile(0.99), - data.gravity_disturbance.quantile(0.99), - ), + cmap="balance+h0", + series=[-cpt_lim, cpt_lim], background=True, ) @@ -137,7 +139,7 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a50f25", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.shift_origin(xshift=fig_width + 1) @@ -150,6 +152,6 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a50f25", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.show() diff --git a/doc/gallery_src/equivalent_sources/spherical.py b/doc/gallery_src/equivalent_sources/spherical.py index d41b43906..255329625 100644 --- a/doc/gallery_src/equivalent_sources/spherical.py +++ b/doc/gallery_src/equivalent_sources/spherical.py @@ -80,17 +80,19 @@ # Plot observed and gridded gravity disturbance fig = pygmt.Figure() +# Get the 99.9th percentile of the absolute value of the point data to use as color +# scale limits +cpt_lim = np.quantile(np.abs(gravity_disturbance), 0.999) + # Make colormap of data -# Get the 90% of the maximum absolute value between the original and gridded -# data so we can use the same color scale for both plots and have 0 centered -# at the white color. -maxabs = vd.maxabs(gravity_disturbance, grid.gravity_disturbance.values) * 0.90 pygmt.makecpt( - cmap="vik", - series=(-maxabs, maxabs), + cmap="balance+h0", + series=[-cpt_lim, cpt_lim], background=True, ) +title = "Observed gravity disturbance data" + fig.plot( projection="M10c", region=region, @@ -102,7 +104,7 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a100f50", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.shift_origin(xshift="w+3c") @@ -112,6 +114,6 @@ cmap=True, ) -fig.colorbar(cmap=True, frame=["a100f50", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.show() diff --git a/doc/gallery_src/forward/point_gravity.py b/doc/gallery_src/forward/point_gravity.py index 38c4a2214..b822188b5 100644 --- a/doc/gallery_src/forward/point_gravity.py +++ b/doc/gallery_src/forward/point_gravity.py @@ -56,7 +56,7 @@ maxabs = vd.maxabs(gravity) * 0.80 -pygmt.makecpt(cmap="vik", series=(-maxabs, maxabs, 0.3)) +pygmt.makecpt(cmap="vik", series=[-maxabs, maxabs, 0.3], background=True) with pygmt.config(FONT_TITLE="16p"): fig.grdimage( @@ -70,6 +70,8 @@ fig.plot(x=easting, y=northing, style="c0.2c", fill="grey") -fig.colorbar(cmap=True, position="JMR", frame=["a.6f.2", "x+lmGal"]) +fig.colorbar(cmap=True, position="JMR+e", frame=["x+lmGal"]) + +fig.legend() fig.show() diff --git a/doc/gallery_src/forward/prism_layer.py b/doc/gallery_src/forward/prism_layer.py index bec45c826..4cf887b73 100644 --- a/doc/gallery_src/forward/prism_layer.py +++ b/doc/gallery_src/forward/prism_layer.py @@ -62,6 +62,6 @@ cmap="viridis", ) -fig.colorbar(cmap=True, position="JMR", frame=["a2f1", "x+lmGal"]) +fig.colorbar(cmap=True, position="JMR", frame=["x+lmGal"]) fig.show() diff --git a/doc/gallery_src/forward/prisms_topo_gravity.py b/doc/gallery_src/forward/prisms_topo_gravity.py index d71d7f651..de4046264 100644 --- a/doc/gallery_src/forward/prisms_topo_gravity.py +++ b/doc/gallery_src/forward/prisms_topo_gravity.py @@ -87,15 +87,21 @@ title = "Gravitational acceleration of the topography" +# Get the max absolute value to use as color scale limits +cpt_lims = vd.maxabs(grid.gravity) + +# Make colormap of data +pygmt.makecpt(cmap="balance+h0", series=[-cpt_lims, cpt_lims]) + with pygmt.config(FONT_TITLE="14p"): fig.grdimage( region=xy_region, projection=fig_proj, grid=grid.gravity, frame=["ag", f"+t{title}"], - cmap="vik", + cmap=True, ) -fig.colorbar(cmap=True, frame=["a100f50", "x+lmGal"]) +fig.colorbar(cmap=True, position="JMR", frame=["a100f50", "x+lmGal"]) fig.show() diff --git a/doc/gallery_src/forward/tesseroid_layer.py b/doc/gallery_src/forward/tesseroid_layer.py index 4d26564df..5e42984c6 100644 --- a/doc/gallery_src/forward/tesseroid_layer.py +++ b/doc/gallery_src/forward/tesseroid_layer.py @@ -55,13 +55,19 @@ # Plot gravity field fig = pygmt.Figure() -maxabs = vd.maxabs(gravity.g_z) -pygmt.makecpt(cmap="polar", series=(-maxabs, maxabs)) + +# Get the max absolute value to use as color scale limits +cpt_lims = vd.maxabs(gravity.g_z) + +# Make colormap of data +pygmt.makecpt(cmap="balance+h0", series=[-cpt_lims, cpt_lims]) fig.grdimage( gravity.g_z, projection="M15c", nan_transparent=True, + cmap=True, ) + fig.basemap(frame=True) fig.colorbar(frame='af+l"Gravity [mGal]"', position="JCR") fig.coast(shorelines="0.5p,black", borders=["1/0.5p,black"]) diff --git a/doc/gallery_src/gravity_disturbance.py b/doc/gallery_src/gravity_disturbance.py index 398771d00..df48156ae 100644 --- a/doc/gallery_src/gravity_disturbance.py +++ b/doc/gallery_src/gravity_disturbance.py @@ -18,6 +18,7 @@ import boule as bl import ensaio +import numpy as np import pygmt import xarray as xr @@ -36,7 +37,11 @@ # Make a plot of data using PyGMT fig = pygmt.Figure() -pygmt.grd2cpt(grid=disturbance, cmap="polar", continuous=True) +# Get the 99.9th percentile of the absolute value to use as color scale limits +cpt_lims = np.quantile(np.abs(disturbance), 0.999) + +# Make colormap of data +pygmt.makecpt(cmap="balance+h0", series=[-cpt_lims, cpt_lims], background=True) title = "Gravity disturbance of the Earth" @@ -50,6 +55,6 @@ fig.coast(shorelines="0.5p,black", resolution="crude") -fig.colorbar(cmap=True, frame=["a100f50", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.show() diff --git a/doc/gallery_src/gravity_disturbance_topofree.py b/doc/gallery_src/gravity_disturbance_topofree.py index 07601d7c3..94e353090 100644 --- a/doc/gallery_src/gravity_disturbance_topofree.py +++ b/doc/gallery_src/gravity_disturbance_topofree.py @@ -28,6 +28,7 @@ import boule as bl import ensaio +import numpy as np import pygmt import xarray as xr @@ -62,7 +63,11 @@ # Make a plot of data using PyGMT fig = pygmt.Figure() -pygmt.grd2cpt(grid=disturbance_topofree, cmap="vik+h0", continuous=True) +# Get the 99th percentile of the absolute value to use as color scale limits +cpt_lims = np.quantile(np.abs(disturbance_topofree), 0.99) + +# Make colormap of data +pygmt.makecpt(cmap="balance+h0", series=[-cpt_lims, cpt_lims], background=True) title = "Topography-free (Bouguer) gravity disturbance of the Earth" @@ -77,6 +82,6 @@ fig.coast(shorelines="0.5p,black", resolution="crude") -fig.colorbar(cmap=True, frame=["a200f50", "x+lmGal"]) +fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") fig.show() diff --git a/doc/gallery_src/isostatic_moho_airy.py b/doc/gallery_src/isostatic_moho_airy.py index cc1e0180a..92e6590b1 100644 --- a/doc/gallery_src/isostatic_moho_airy.py +++ b/doc/gallery_src/isostatic_moho_airy.py @@ -56,10 +56,11 @@ print("\nMoho depth grid:") print(moho) -# Draw the maps +# Make a plot of data using PyGMT fig = pygmt.Figure() -pygmt.grd2cpt(grid=moho, cmap="viridis", reverse=True, continuous=True) +# Make colormap of data +pygmt.makecpt(cmap="viridis", series=[moho.to_numpy().min(), moho.to_numpy().max()]) title = "Airy isostatic Moho depth of Africa" @@ -73,6 +74,6 @@ fig.coast(shorelines="0.5p,black", resolution="crude") -fig.colorbar(cmap=True, frame=["a10000f2500", "x+lmeters"]) +fig.colorbar(cmap=True, frame=["x+lmeters"]) fig.show() diff --git a/doc/gallery_src/transformations/reduction_to_pole.py b/doc/gallery_src/transformations/reduction_to_pole.py index 8a785a36b..d0e5944d2 100644 --- a/doc/gallery_src/transformations/reduction_to_pole.py +++ b/doc/gallery_src/transformations/reduction_to_pole.py @@ -52,9 +52,9 @@ # Plot original magnetic anomaly and the reduced to the pole fig = pygmt.Figure() with fig.subplot(nrows=1, ncols=2, figsize=("28c", "15c"), sharey="l"): - # Make colormap for both plots (saturate it a little bit) - scale = 0.5 * vd.maxabs(magnetic_grid, rtp_grid) - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + # Make colormap for both plots + cpt_lim = 0.5 * vd.maxabs(magnetic_grid, rtp_grid) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) with fig.set_panel(panel=0): # Plot magnetic anomaly grid fig.grdimage( @@ -65,7 +65,7 @@ # Add colorbar fig.colorbar( frame='af+l"Magnetic anomaly [nT]"', - position="JBC+h+o0/1c+e", + position="JBC+h+o0/1c+ef", ) with fig.set_panel(panel=1): # Plot upward reduced to the pole grid @@ -73,6 +73,6 @@ # Add colorbar fig.colorbar( frame='af+l"Reduced to the pole [nT]"', - position="JBC+h+o0/1c+e", + position="JBC+h+o0/1c+ef", ) fig.show() diff --git a/doc/gallery_src/transformations/tga.py b/doc/gallery_src/transformations/tga.py index 21c710e5a..d5fe8396e 100644 --- a/doc/gallery_src/transformations/tga.py +++ b/doc/gallery_src/transformations/tga.py @@ -45,8 +45,8 @@ with fig.subplot(nrows=1, ncols=2, figsize=("28c", "15c"), sharey="l"): with fig.set_panel(panel=0): # Make colormap of data - scale = 2500 - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + cpt_lim = vd.maxabs(magnetic_grid) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot magnetic anomaly grid fig.grdimage( grid=magnetic_grid, @@ -56,17 +56,17 @@ # Add colorbar fig.colorbar( frame='af+l"Magnetic anomaly [nT]"', - position="JBC+h+o0/1c+e", + position="JBC+h+o0/1c", ) with fig.set_panel(panel=1): # Make colormap for total gradient amplitude (saturate it a little bit) - scale = 0.6 * vd.maxabs(tga) - pygmt.makecpt(cmap="polar+h", series=[0, scale], background=True) + cpt_lim = 0.6 * vd.maxabs(tga) + pygmt.makecpt(cmap="balance+h0", series=[0, cpt_lim], background=True) # Plot total gradient amplitude fig.grdimage(grid=tga, projection="X?", cmap=True) # Add colorbar fig.colorbar( frame='af+l"Total Gradient Amplitude [nT/m]"', - position="JBC+h+o0/1c+e", + position="JBC+h+o0/1c+ef", ) fig.show() diff --git a/doc/gallery_src/transformations/tilt.py b/doc/gallery_src/transformations/tilt.py index 9f9c1b579..e62234dd7 100644 --- a/doc/gallery_src/transformations/tilt.py +++ b/doc/gallery_src/transformations/tilt.py @@ -79,10 +79,10 @@ sharey="l", margins=["1c", "1c"], ): - scale = 0.5 * vd.maxabs(magnetic_grid, rtp_grid) + cpt_lim = 0.5 * vd.maxabs(magnetic_grid, rtp_grid) with fig.set_panel(panel=0): # Make colormap of data - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot magnetic anomaly grid fig.grdimage( grid=magnetic_grid, @@ -92,7 +92,7 @@ ) with fig.set_panel(panel=1): # Make colormap of data - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot reduced to the pole magnetic anomaly grid fig.grdimage( grid=rtp_grid, @@ -104,13 +104,13 @@ label = "nT" fig.colorbar( frame=f"af+l{label}", - position="JMR+o1/-0.25c+e", + position="JMR+o1/-0.25c+ef", ) - scale = 0.6 * vd.maxabs(tilt_grid, tilt_rtp_grid) + cpt_lim = vd.maxabs(tilt_grid, tilt_rtp_grid) with fig.set_panel(panel=2): # Make colormap for tilt (saturate it a little bit) - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot tilt fig.grdimage( grid=tilt_grid, @@ -120,7 +120,7 @@ ) with fig.set_panel(panel=3): # Make colormap for tilt rtp (saturate it a little bit) - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot tilt fig.grdimage( grid=tilt_rtp_grid, @@ -132,6 +132,6 @@ label = "rad" fig.colorbar( frame=f"af+l{label}", - position="JMR+o1/-0.25c+e", + position="JMR+o1/-0.25c", ) fig.show() diff --git a/doc/gallery_src/transformations/upward_continuation.py b/doc/gallery_src/transformations/upward_continuation.py index db6ce9829..edd76893d 100644 --- a/doc/gallery_src/transformations/upward_continuation.py +++ b/doc/gallery_src/transformations/upward_continuation.py @@ -11,6 +11,7 @@ import ensaio import pygmt +import verde as vd import xarray as xr import xrft @@ -44,9 +45,9 @@ # Plot original magnetic anomaly and the upward continued grid fig = pygmt.Figure() with fig.subplot(nrows=1, ncols=2, figsize=("28c", "15c"), sharey="l"): - # Make colormap for both plots data - scale = 2500 - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + # Make colormap based on original data + cpt_lim = vd.maxabs(magnetic_grid) * 0.6 + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) with fig.set_panel(panel=0): # Plot magnetic anomaly grid fig.grdimage( @@ -65,6 +66,6 @@ # Add colorbar fig.colorbar( frame='af+l"Upward continued to 1000m [nT]"', - position="JBC+h+o0/1c+e", + position="JBC+h+o0/1c", ) fig.show() diff --git a/doc/gallery_src/transformations/upward_derivative.py b/doc/gallery_src/transformations/upward_derivative.py index 2d48b89db..30c74618a 100644 --- a/doc/gallery_src/transformations/upward_derivative.py +++ b/doc/gallery_src/transformations/upward_derivative.py @@ -46,8 +46,8 @@ with fig.subplot(nrows=1, ncols=2, figsize=("28c", "15c"), sharey="l"): with fig.set_panel(panel=0): # Make colormap of data - scale = 2500 - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + cpt_lim = vd.maxabs(magnetic_grid) * 0.6 + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot magnetic anomaly grid fig.grdimage( grid=magnetic_grid, @@ -61,8 +61,8 @@ ) with fig.set_panel(panel=1): # Make colormap for upward derivative (saturate it a little bit) - scale = 0.6 * vd.maxabs(deriv_upward) - pygmt.makecpt(cmap="polar+h", series=[-scale, scale], background=True) + cpt_lim = vd.maxabs(deriv_upward) * 0.6 + pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot upward derivative fig.grdimage(grid=deriv_upward, projection="X?", cmap=True) # Add colorbar diff --git a/doc/user_guide/equivalent_sources/block-averaged-eqs.rst b/doc/user_guide/equivalent_sources/block-averaged-eqs.rst index c27e08dbb..77c5f28f4 100644 --- a/doc/user_guide/equivalent_sources/block-averaged-eqs.rst +++ b/doc/user_guide/equivalent_sources/block-averaged-eqs.rst @@ -94,7 +94,7 @@ And project the geographic coordinates to plain Cartesian ones: title = "Observed total-field magnetic anomaly" pygmt.makecpt( - cmap="polar+h0", + cmap="balance+h0", series=(-maxabs, maxabs), background=True, ) @@ -182,7 +182,7 @@ we are efectivelly upward continuing the data. title = "Observed magnetic anomaly data" pygmt.makecpt( - cmap="polar+h0", + cmap="balance+h0", series=(-maxabs, maxabs), background=True) diff --git a/doc/user_guide/equivalent_sources/eq-sources-spherical.rst b/doc/user_guide/equivalent_sources/eq-sources-spherical.rst index 9d98993a3..32b91a654 100644 --- a/doc/user_guide/equivalent_sources/eq-sources-spherical.rst +++ b/doc/user_guide/equivalent_sources/eq-sources-spherical.rst @@ -167,7 +167,7 @@ Lets plot it: fig = pygmt.Figure() # Make colormap of data - pygmt.makecpt(cmap="polar+h0",series=(-maxabs, maxabs,)) + pygmt.makecpt(cmap="balance+h0",series=(-maxabs, maxabs,)) title = "Block-median reduced gravity disturbance" fig.plot( diff --git a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst index 22b5c2c57..615da657f 100644 --- a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst +++ b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst @@ -216,7 +216,7 @@ Lets plot it: fig = pygmt.Figure() # Make colormap of data - pygmt.makecpt(cmap="polar+h0",series=(-maxabs, maxabs,)) + pygmt.makecpt(cmap="balance+h0",series=(-maxabs, maxabs,)) title = "Gravity disturbance with first guess" diff --git a/doc/user_guide/equivalent_sources/gradient-boosted-eqs.rst b/doc/user_guide/equivalent_sources/gradient-boosted-eqs.rst index 5cb53f00b..2aa6a2784 100644 --- a/doc/user_guide/equivalent_sources/gradient-boosted-eqs.rst +++ b/doc/user_guide/equivalent_sources/gradient-boosted-eqs.rst @@ -173,7 +173,7 @@ And plot it: fig = pygmt.Figure() # Make colormap of data - pygmt.makecpt(cmap="polar+h0",series=(-maxabs, maxabs,)) + pygmt.makecpt(cmap="balance+h0",series=(-maxabs, maxabs,)) title = "Observed gravity disturbance data" with pygmt.config(FONT_TITLE="14p"): diff --git a/doc/user_guide/equivalent_sources/index.rst b/doc/user_guide/equivalent_sources/index.rst index 266a9e3e7..bbdec13d6 100644 --- a/doc/user_guide/equivalent_sources/index.rst +++ b/doc/user_guide/equivalent_sources/index.rst @@ -143,7 +143,7 @@ And plot it: fig_proj = f"x1:{fig_ratio}" fig = pygmt.Figure() - pygmt.makecpt(cmap="polar+h0", series=[-maxabs, maxabs]) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs]) title="Predicted gravity disturbance" with pygmt.config(FONT_TITLE="14p"): fig.plot( @@ -204,7 +204,7 @@ And plot it maxabs = vd.maxabs(grid.gravity_disturbance) fig = pygmt.Figure() - pygmt.makecpt(cmap="polar+h0", series=[-maxabs, maxabs]) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) fig.grdimage( frame=['af', 'WSen'], grid=grid.gravity_disturbance, diff --git a/doc/user_guide/forward_modelling/point.rst b/doc/user_guide/forward_modelling/point.rst index 004058c0d..8f6999474 100644 --- a/doc/user_guide/forward_modelling/point.rst +++ b/doc/user_guide/forward_modelling/point.rst @@ -113,7 +113,7 @@ Lets plot this gravitational field: fig = pygmt.Figure() maxabs = vd.maxabs(g_z) - pygmt.makecpt(cmap="polar", series=(-maxabs, maxabs), background=True) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) fig.grdimage( region=(-250, 1250, -250, 1250), @@ -239,8 +239,8 @@ Lets plot these results using :mod:`pygmt`: fig = pygmt.Figure() title = "Gravitational acceleration (downward)" - maxabs = vd.maxabs(g_z)*.95 - pygmt.makecpt(cmap="polar", series=(-maxabs, maxabs), background=True) + maxabs = vd.maxabs(g_z) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) fig.grdimage( region=(-72, -68, -46, -42), diff --git a/doc/user_guide/forward_modelling/tesseroid.rst b/doc/user_guide/forward_modelling/tesseroid.rst index 21c4a5c47..4db314fcc 100644 --- a/doc/user_guide/forward_modelling/tesseroid.rst +++ b/doc/user_guide/forward_modelling/tesseroid.rst @@ -111,6 +111,7 @@ And finally plot the computed gravitational field .. jupyter-execute:: import pygmt + grid = vd.make_xarray_grid( coordinates, gravity, data_names="gravity", extra_coords_names="extra") diff --git a/doc/user_guide/gravity_disturbance.rst b/doc/user_guide/gravity_disturbance.rst index a27da1c5e..6aef78354 100644 --- a/doc/user_guide/gravity_disturbance.rst +++ b/doc/user_guide/gravity_disturbance.rst @@ -143,7 +143,7 @@ And plot it: maxabs = vd.maxabs(gravity_disturbance) fig = pygmt.Figure() - pygmt.makecpt(series=[-maxabs, maxabs], cmap="polar+h") + pygmt.makecpt(series=[-maxabs, maxabs], cmap="balance+h0") fig.grdimage( gravity_disturbance, projection="W20c", diff --git a/doc/user_guide/igrf.rst b/doc/user_guide/igrf.rst index 9eff60d82..eef99ee2f 100644 --- a/doc/user_guide/igrf.rst +++ b/doc/user_guide/igrf.rst @@ -15,6 +15,7 @@ generation IGRF field with :class:`harmonica.IGRF14`. Here's how it works. import pygmt import boule as bl import harmonica as hm + import verde as vd All of the functionality is wrapped in the :class:`~harmonica.IGRF14` class. When creating an instance of it, we need to provide the date on which we want @@ -147,24 +148,31 @@ We can plot the three components using :mod:`pygmt` on a nice map: .. jupyter-execute:: fig = pygmt.Figure() + for c in ["b_east", "b_north", "b_up"]: - fig.grdimage( - grid[c], cmap="polar+h", projection="W15c", region="g", - ) - fig.coast(shorelines=True) - fig.colorbar( - position="JMR+ml+o0.5c", - frame=[ + max_abs = vd.maxabs(grid[c]) + pygmt.makecpt( + cmap="balance+h0", + series=[-max_abs, max_abs], + background=True, + ) + fig.grdimage( + grid[c], cmap=True, projection="W15c", region="g", + ) + fig.coast(shorelines=True) + fig.colorbar( + position="JMR+ml+o0.5c", + frame=[ "a10000", f"x+l{grid[c].attrs['long_name']}", f"y+l{grid[c].attrs['units']}", - ] - ) - if c == "b_east": - fig.basemap(frame=["a", f"+t{grid.attrs['title']}"]) - else: - fig.basemap(frame="a") - fig.shift_origin(yshift="-h-0.5c") + ] + ) + if c == "b_east": + fig.basemap(frame=["a", f"+t{grid.attrs['title']}"]) + else: + fig.basemap(frame="a") + fig.shift_origin(yshift="-h-0.5c") fig.show() The grid spacing was calculated automatically to match the maximum degree of @@ -201,33 +209,33 @@ with PyGMT the same way we did the vector components: fig = pygmt.Figure() cmaps = { - "intensity": "viridis", - "inclination": "polar+h", - "declination": "polar+h", + "intensity": "viridis", + "inclination": "balance+h0", + "declination": "balance+h0", } cb_annot = { - "intensity": "a10000", - "inclination": "a20+u\\260", - "declination": "a40+u\\260", + "intensity": "a10000", + "inclination": "a20+u\\260", + "declination": "a40+u\\260", } for c in ["intensity", "inclination", "declination"]: - fig.grdimage( - grid[c], cmap=cmaps[c], projection="W0/15c", - ) - fig.coast(shorelines=True) - fig.colorbar( - position="JMR+ml+o0.5c", - frame=[ + fig.grdimage( + grid[c], cmap=cmaps[c], projection="W0/15c", + ) + fig.coast(shorelines=True) + fig.colorbar( + position="JMR+ml+o0.5c", + frame=[ cb_annot[c], f"x+l{grid[c].attrs['long_name']}", f"y+l{grid[c].attrs['units']}", - ] - ) - if c == "intensity": - fig.basemap(frame=["a", f"+t{grid.attrs['title']}"]) - else: - fig.basemap(frame="a") - fig.shift_origin(yshift="-h-0.5c") + ] + ) + if c == "intensity": + fig.basemap(frame=["a", f"+t{grid.attrs['title']}"]) + else: + fig.basemap(frame="a") + fig.shift_origin(yshift="-h-0.5c") fig.show() We can clearly see the `South Atlantic Magnetic Anomaly diff --git a/doc/user_guide/topographic_correction.rst b/doc/user_guide/topographic_correction.rst index ec7b60c07..0e1fb0080 100644 --- a/doc/user_guide/topographic_correction.rst +++ b/doc/user_guide/topographic_correction.rst @@ -72,7 +72,7 @@ And plot it: maxabs = vd.maxabs(data.gravity_disturbance_mgal) fig = pygmt.Figure() - pygmt.makecpt(cmap="polar+h0", series=[-maxabs, maxabs]) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs]) fig.plot( x=data.longitude, y=data.latitude, @@ -125,7 +125,7 @@ We can now compute the Bouguer disturbance and plot it: maxabs = vd.maxabs(bouguer_disturbance) fig = pygmt.Figure() - pygmt.makecpt(cmap="polar+h0", series=[-maxabs, maxabs]) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs]) fig.plot( x=data.longitude, y=data.latitude, @@ -245,7 +245,7 @@ And plot it: maxabs = vd.maxabs(topo_free_disturbance) fig = pygmt.Figure() - pygmt.makecpt(cmap="polar+h0", series=[-maxabs, maxabs]) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs]) fig.plot( x=data.longitude, y=data.latitude, From 15c844a9d13fbc4267f75bb460863c8af4c6f9fc Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 09:54:37 +0100 Subject: [PATCH 02/10] use pygmt instead of matplotlib for documentation --- doc/user_guide/forward_modelling/dipole.rst | 70 ++++- doc/user_guide/forward_modelling/prism.rst | 289 +++++++++++++------ doc/user_guide/transformations.rst | 300 +++++++++++++------- 3 files changed, 471 insertions(+), 188 deletions(-) diff --git a/doc/user_guide/forward_modelling/dipole.rst b/doc/user_guide/forward_modelling/dipole.rst index e7a5e9a6d..3d42046ff 100644 --- a/doc/user_guide/forward_modelling/dipole.rst +++ b/doc/user_guide/forward_modelling/dipole.rst @@ -92,20 +92,70 @@ every observation point: b_e, b_n, b_u = hm.dipole_magnetic(coordinates, dipoles, magnetic_moments, field="b") +We can reshape the results into variables of an dataset: + +.. jupyter-execute:: + + grid = vd.make_xarray_grid( + coordinates, + data=(b_e, b_n, b_u), + data_names=["b_e", "b_n", "b_u"], + extra_coords_names="extra" + ) + .. jupyter-execute:: + :hide-code: + + import pygmt + + # Needed so that displaying works on jupyter-sphinx and sphinx-gallery at + # the same time. Using PYGMT_USE_EXTERNAL_DISPLAY="false" in the Makefile + # for sphinx-gallery to work means that fig.show won't display anything here + # either. + pygmt.set_display(method="notebook") - import matplotlib.pyplot as plt - fig, axes = plt.subplots(nrows=1, ncols=3, sharey=True, figsize=(12, 8)) +.. jupyter-execute:: + + import pygmt + + fig = pygmt.Figure() + + maxabs = vd.maxabs(grid.b_e) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_e, + projection="X10c", + cmap=True, + frame=["WSne+tEasting component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.b_n) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_n, + projection="X10c", + cmap=True, + frame=["wSne+tNorthing component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.b_u) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_u, + projection="X10c", + cmap=True, + frame=["wSnE+tUpward component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') - fields = {"b_e": b_e, "b_n": b_n, "b_u": b_u} - for field, ax in zip(fields, axes): - tmp = ax.pcolormesh(coordinates[0], coordinates[1], fields[field]) - ax.set_aspect("equal") - ax.set_title(field) - ax.ticklabel_format(style="sci", scilimits=(0, 0), axis="both") - plt.colorbar(tmp, ax=ax, orientation="horizontal", label="nT", pad=0.008) - plt.show() + fig.show() ---- diff --git a/doc/user_guide/forward_modelling/prism.rst b/doc/user_guide/forward_modelling/prism.rst index 1794defb0..08a75112f 100644 --- a/doc/user_guide/forward_modelling/prism.rst +++ b/doc/user_guide/forward_modelling/prism.rst @@ -90,54 +90,167 @@ point: for field in fields: results[field] = hm.prism_gravity(coordinates, prism, density, field=field) +We can reshape the results into variables of an dataset: + +.. jupyter-execute:: + + grid = vd.make_xarray_grid( + coordinates, + tuple(results.values()), + data_names=results.keys(), + extra_coords_names="extra", + ) + print(grid) + Plot the results: +.. jupyter-execute:: + :hide-code: + + import pygmt + + # Needed so that displaying works on jupyter-sphinx and sphinx-gallery at + # the same time. Using PYGMT_USE_EXTERNAL_DISPLAY="false" in the Makefile + # for sphinx-gallery to work means that fig.show won't display anything here + # either. + pygmt.set_display(method="notebook") + .. jupyter-execute:: - import matplotlib.pyplot as plt + import pygmt - plt.pcolormesh(coordinates[0], coordinates[1], results["potential"]) - plt.gca().set_aspect("equal") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(label="J/kg") - plt.show() + fig = pygmt.Figure() + + fig.grdimage( + projection="X10c", + grid=grid.potential, + frame=["a", "x+leasting (m)", "y+lnorthing (m)"], + cmap="viridis", + ) + + fig.colorbar(cmap=True, position="JMR", frame=["x+lJ kg@+-1@+"]) + fig.show() .. jupyter-execute:: - fig, axes = plt.subplots(nrows=1, ncols=3, sharey=True, figsize=(12, 8)) + fig = pygmt.Figure() + + maxabs = vd.maxabs(grid.g_e) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_e, + projection="X10c", + cmap=True, + frame=["WSne+tEasting component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') - for field, ax in zip(("g_e", "g_n", "g_z"), axes): - tmp = ax.pcolormesh(coordinates[0], coordinates[1], results[field]) - ax.set_aspect("equal") - ax.set_title(field) - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, ax=ax, label="mGal", orientation="horizontal", pad=0.08) - plt.show() + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.g_n) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_n, + projection="X10c", + cmap=True, + frame=["wSne+tNorthing component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.g_z) + pygmt.makecpt(cmap="balance+h0", series=[0, maxabs], background=True) + fig.grdimage( + grid=grid.g_z, + projection="X10c", + cmap=True, + frame=["wSnE+tDownward component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.show() .. jupyter-execute:: - fig, axes = plt.subplots(nrows=1, ncols=3, sharey=True, figsize=(12, 8)) + fig = pygmt.Figure() - for field, ax in zip(("g_ee", "g_nn", "g_zz"), axes): - tmp = ax.pcolormesh(coordinates[0], coordinates[1], results[field]) - ax.set_aspect("equal") - ax.set_title(field) - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, ax=ax, label="Eotvos", orientation="horizontal", pad=0.08) - plt.show() + maxabs = vd.maxabs(grid.g_ee) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_ee, + projection="X10c", + cmap=True, + frame=["WSne+tEasting-easting tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.g_nn) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_nn, + projection="X10c", + cmap=True, + frame=["wSne+tNorthing-northing tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.g_zz) + pygmt.makecpt(cmap="balance+h0", series=[0, maxabs], background=True) + fig.grdimage( + grid=grid.g_zz, + projection="X10c", + cmap=True, + frame=["wSnE+tDownward-downward tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.show() .. jupyter-execute:: - fig, axes = plt.subplots(nrows=1, ncols=3, sharey=True, figsize=(12, 8)) + fig = pygmt.Figure() + + maxabs = vd.maxabs(grid.g_en) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_en, + projection="X10c", + cmap=True, + frame=["WSne+tEasting-northing tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.shift_origin(xshift="11c") - for field, ax in zip(("g_en", "g_ez", "g_nz"), axes): - tmp = ax.pcolormesh(coordinates[0], coordinates[1], results[field]) - ax.set_aspect("equal") - ax.set_title(field) - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, ax=ax, label="Eotvos", orientation="horizontal", pad=0.08) - plt.show() + maxabs = vd.maxabs(grid.g_ez) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_ez, + projection="X10c", + cmap=True, + frame=["wSne+tEasting-downward tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.g_nz) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.g_nz, + projection="X10c", + cmap=True, + frame=["wSnE+tNorthing-downward tensor", "xa","ya"], + ) + fig.colorbar(frame='+lEotvos') + + fig.show() Passing multiple prisms @@ -186,20 +299,7 @@ generated by the whole set of prisms on every computation point: Lets plot this gravitational field: .. jupyter-execute:: - :hide-code: - - import pygmt - # Needed so that displaying works on jupyter-sphinx and sphinx-gallery at - # the same time. Using PYGMT_USE_EXTERNAL_DISPLAY="false" in the Makefile - # for sphinx-gallery to work means that fig.show won't display anything here - # either. - pygmt.set_display(method="notebook") - - -.. jupyter-execute:: - - import pygmt grid = vd.make_xarray_grid( coordinates, g_z, data_names="g_z", extra_coords_names="extra") fig = pygmt.Figure() @@ -259,40 +359,57 @@ points by choosing ``field="b"``: b_e, b_n, b_u = hm.prism_magnetic(coordinates, prisms, magnetization, field="b") + +We can reshape the results into variables of an dataset: + .. jupyter-execute:: - fig, axes = plt.subplots(nrows=1, ncols=3, sharey=True, figsize=(12, 6)) - - for ax, mag_component, title in zip(axes, (b_e, b_n, b_u), ("Be", "Bn", "Bu")): - maxabs = vd.maxabs(mag_component) - tmp = ax.pcolormesh( - coordinates[0], - coordinates[1], - mag_component, - vmin=-maxabs, - vmax=maxabs, - cmap="RdBu_r", - ) - ax.contour( - coordinates[0], - coordinates[1], - mag_component, - colors="k", - linewidths=0.5, - ) - ax.set_title(title) - ax.set_aspect("equal") - plt.colorbar( - tmp, - ax=ax, - orientation="horizontal", - label="nT", - pad=0.08, - aspect=42, - shrink=0.8, - ) - - plt.show() + grid = vd.make_xarray_grid( + coordinates, + data=(b_e, b_n, b_u), + data_names=["b_e", "b_n", "b_u"], + extra_coords_names="extra" + ) + +.. jupyter-execute:: + + fig = pygmt.Figure() + + maxabs = vd.maxabs(grid.b_e) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_e, + projection="X10c", + cmap=True, + frame=["WSne+tEasting component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.b_n) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_n, + projection="X10c", + cmap=True, + frame=["wSne+tNorthing component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.shift_origin(xshift="11c") + + maxabs = vd.maxabs(grid.b_u) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_u, + projection="X10c", + cmap=True, + frame=["wSnE+tUpward component", "xa","ya"], + ) + fig.colorbar(frame='+lnT') + + fig.show() Alternatively, we can compute just a single component by choosing ``field`` to @@ -318,16 +435,22 @@ generated by these two prisms: .. jupyter-execute:: - maxabs = vd.maxabs(b_u) + fig = pygmt.Figure() + + grid = vd.make_xarray_grid( + coordinates, b_u, data_names="b_u", extra_coords_names="extra") - tmp = plt.pcolormesh( - coordinates[0], coordinates[1], b_u, vmin=-maxabs, vmax=maxabs, cmap="RdBu_r" + maxabs = vd.maxabs(grid.b_u) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + fig.grdimage( + grid=grid.b_u, + projection="X10c", + cmap=True, + frame=["WSne+tUpward component", "xa","ya"], ) - plt.contour(coordinates[0], coordinates[1], b_u, colors="k", linewidths=0.5) - plt.title("Bu") - plt.gca().set_aspect("equal") - plt.colorbar(tmp, label="nT", pad=0.03, aspect=42, shrink=0.8) - plt.show() + fig.colorbar(frame='+lnT') + + fig.show() .. _prism_layer: diff --git a/doc/user_guide/transformations.rst b/doc/user_guide/transformations.rst index ce9ac9734..769b717d5 100644 --- a/doc/user_guide/transformations.rst +++ b/doc/user_guide/transformations.rst @@ -25,15 +25,39 @@ We can load the data file using :mod:`xarray`: And plot it: .. jupyter-execute:: + :hide-code: - import matplotlib.pyplot as plt + import pygmt - tmp = magnetic_grid.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Magnetic anomaly grid") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT") - plt.show() + # Needed so that displaying works on jupyter-sphinx and sphinx-gallery at + # the same time. Using PYGMT_USE_EXTERNAL_DISPLAY="false" in the Makefile + # for sphinx-gallery to work means that fig.show won't display anything here + # either. + pygmt.set_display(method="notebook") + + +.. jupyter-execute:: + + import pygmt + import verde as vd + + fig = pygmt.Figure() + + maxabs = vd.maxabs(magnetic_grid) * .6 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + magnetic_grid, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tMagnetic anomaly"] + ) + fig.colorbar( + position="JCB+e", + frame=["af", 'x+lnT'], + ) + + fig.show() .. seealso:: @@ -68,12 +92,23 @@ needed by the :func:`xrft.pad` function): .. jupyter-execute:: - tmp = magnetic_grid_padded.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Padded magnetic anomaly grid") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT") - plt.show() + fig = pygmt.Figure() + + maxabs = vd.maxabs(magnetic_grid_padded) * .6 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + magnetic_grid_padded, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tPadded magnetic anomaly"] + ) + fig.colorbar( + position="JCB+e", + frame=["af", 'x+lnT'], + ) + + fig.show() Now that we have the padded grid, we can apply any grid transformation. @@ -103,12 +138,23 @@ And plot it: .. jupyter-execute:: - tmp = deriv_upward.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Upward derivative of the magnetic anomaly") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT/m") - plt.show() + fig = pygmt.Figure() + + maxabs = vd.maxabs(deriv_upward) * .5 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + deriv_upward, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tUpward derivative of the magnetic anomaly"] + ) + fig.colorbar( + position="JCB+e", + frame=["af", 'x+lnT/m'], + ) + + fig.show() Horizontal derivatives @@ -132,24 +178,35 @@ And plot them: .. jupyter-execute:: - fig, (ax1, ax2) = plt.subplots( - nrows=1, ncols=2, sharey=True, figsize=(12, 8) + fig = pygmt.Figure() + + maxabs = vd.maxabs(deriv_easting, deriv_northing) * .4 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + deriv_easting, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tEasting derivative of the magnetic anomaly"] ) - cbar_kwargs=dict( - label="nT/m", orientation="horizontal", shrink=0.8, pad=0.08, aspect=42 + fig.shift_origin(xshift="16c") + + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + deriv_northing, + projection="X15c", + cmap=True, + frame=["af", "wESn+tNorthing derivative of the magnetic anomaly"] ) - kwargs = dict(center=0, cmap="seismic", cbar_kwargs=cbar_kwargs) - tmp = deriv_easting.plot(ax=ax1, **kwargs) - tmp = deriv_northing.plot(ax=ax2, **kwargs) + fig.colorbar( + position="JBC+h+e+o-8c/1c+w15c/.8c", + frame=["af", 'x+lnT/m'], + ) - ax1.set_title("Easting derivative of the magnetic anomaly") - ax2.set_title("Northing derivative of the magnetic anomaly") - for ax in (ax1, ax2): - ax.set_aspect("equal") - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - plt.show() + fig.show() By default, these two functions compute the horizontal derivatives using central finite differences methods. We can choose to use either the finite @@ -172,25 +229,35 @@ frequency domain: .. jupyter-execute:: - fig, (ax1, ax2) = plt.subplots( - nrows=1, ncols=2, sharey=True, figsize=(12, 8) - ) + fig = pygmt.Figure() + + maxabs = vd.maxabs(deriv_easting, deriv_northing) * .4 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) - cbar_kwargs=dict( - label="nT/m", orientation="horizontal", shrink=0.8, pad=0.08, aspect=42 + fig.grdimage( + deriv_easting, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tEasting derivative of the magnetic anomaly"] ) - kwargs = dict(center=0, cmap="seismic", cbar_kwargs=cbar_kwargs) - tmp = deriv_easting.plot(ax=ax1, **kwargs) - tmp = deriv_northing.plot(ax=ax2, **kwargs) + fig.shift_origin(xshift="16c") + + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) - ax1.set_title("Easting derivative of the magnetic anomaly") - ax2.set_title("Northing derivative of the magnetic anomaly") - for ax in (ax1, ax2): - ax.set_aspect("equal") - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - plt.show() + fig.grdimage( + deriv_northing, + projection="X15c", + cmap=True, + frame=["af", "wESn+tNorthing derivative of the magnetic anomaly"] + ) + + fig.colorbar( + position="JBC+h+e+o-8c/1c+w15c/.8c", + frame=["af", 'x+lnT/m'], + ) + fig.show() .. important:: @@ -228,12 +295,23 @@ And plot it: .. jupyter-execute:: - tmp = upward_continued.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Upward continued magnetic anomaly to 1000m") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT") - plt.show() + fig = pygmt.Figure() + + maxabs = vd.maxabs(upward_continued) + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + upward_continued, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tUpward continued magnetic anomaly at 1000m"] + ) + fig.colorbar( + position="JCB", + frame=["af", 'x+lnT'], + ) + + fig.show() Reduction to the pole @@ -280,12 +358,23 @@ And plot it: .. jupyter-execute:: - tmp = rtp_grid.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Magnetic anomaly reduced to the pole") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT") - plt.show() + fig = pygmt.Figure() + + maxabs = vd.maxabs(rtp_grid) * .8 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + rtp_grid, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tReduced to the pole magnetic anomaly"] + ) + fig.colorbar( + position="JCB+e", + frame=["af", 'x+lnT'], + ) + + fig.show() If on the other hand we have any knowledge about the orientation of the magnetization vector of the sources, we can specify the @@ -309,12 +398,23 @@ magnetization vector of the sources, we can specify the .. jupyter-execute:: - tmp = rtp_grid.plot(cmap="seismic", center=0, add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Reduced to the pole with remanence") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT") - plt.show() + fig = pygmt.Figure() + + maxabs = vd.maxabs(rtp_grid) * .8 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + rtp_grid, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tReduced to the pole with remanence"] + ) + fig.colorbar( + position="JCB+e", + frame=["af", 'x+lnT'], + ) + + fig.show() Gaussian filters @@ -367,34 +467,35 @@ Let's plot the results side by side: .. jupyter-execute:: - import verde as vd + fig = pygmt.Figure() - fig, (ax1, ax2) = plt.subplots( - nrows=1, ncols=2, sharey=True, figsize=(12, 8) + maxabs = vd.maxabs(magnetic_low_freqs, magnetic_high_freqs) * .6 + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + magnetic_low_freqs, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tLow-pass filtered magnetic anomaly"] ) - maxabs = vd.maxabs(magnetic_low_freqs, magnetic_high_freqs) - kwargs = dict(cmap="seismic", vmin=-maxabs, vmax=maxabs, add_colorbar=False) - - tmp = magnetic_low_freqs.plot(ax=ax1, **kwargs) - tmp = magnetic_high_freqs.plot(ax=ax2, **kwargs) - - ax1.set_title("Magnetic anomaly after low-pass filter") - ax2.set_title("Magnetic anomaly after high-pass filter") - for ax in (ax1, ax2): - ax.set_aspect("equal") - ax.ticklabel_format(style="sci", scilimits=(0, 0)) - - plt.colorbar( - tmp, - ax=[ax1, ax2], - label="nT", - orientation="horizontal", - aspect=42, - shrink=0.8, - pad=0.08, + fig.shift_origin(xshift="16c") + + pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) + + fig.grdimage( + magnetic_high_freqs, + projection="X15c", + cmap=True, + frame=["af", "wESn+tHigh-pass filtered magnetic anomaly"] ) - plt.show() + + fig.colorbar( + position="JBC+h+e+o-8c/1c+w15c/.8c", + frame=["af", 'x+lnT'], + ) + + fig.show() Total gradient amplitude @@ -433,14 +534,23 @@ And plot it: .. jupyter-execute:: - import verde as vd + fig = pygmt.Figure() + + maxabs = vd.maxabs(tga_grid) + pygmt.makecpt(cmap="viridis", series=[0, maxabs], background=True) + + fig.grdimage( + tga_grid, + projection="X15c", + cmap=True, + frame=["af", "WeSn+tTotal gradient amplitude of the magnetic anomaly"] + ) + fig.colorbar( + position="JCB", + frame=["af", 'x+lnT/m'], + ) - tmp = tga_grid.plot(cmap="viridis", add_colorbar=False) - plt.gca().set_aspect("equal") - plt.title("Total gradient amplitude of the magnetic anomaly") - plt.gca().ticklabel_format(style="sci", scilimits=(0, 0)) - plt.colorbar(tmp, label="nT/m") - plt.show() + fig.show() ---- From bf00da8eb1bebbbe90f12f61d4fc577c444497d6 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 09:55:31 +0100 Subject: [PATCH 03/10] wrap returns in print in .rst file --- .../equivalent_sources/eqs-parameters-estimation.rst | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst index 615da657f..e7f9fcb3d 100644 --- a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst +++ b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst @@ -39,7 +39,7 @@ Lets start by loading some sample gravity data: fname = ensaio.fetch_bushveld_gravity(version=1) data = pd.read_csv(fname) - data + print(data) And project their coordinates using a Mercator projection: @@ -80,7 +80,7 @@ value of the score obtained after each cross validation. data.gravity_disturbance_mgal, ) ) - score_first_guess + print(score_first_guess) The resulting score corresponds to the R^2. It represents how well the equivalent sources can reproduce the variation of our data. As closer it gets @@ -145,7 +145,7 @@ And now we can actually ran one cross validation for each pair of parameters: ) ) scores.append(score) - scores + print(scores) Once every score has been computed, we can obtain the best score and the corresponding parameters that generate it: From 38afc87a8f0b068c8f45cad7b8bb348cf6094913 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 09:56:00 +0100 Subject: [PATCH 04/10] add hidden pygmt import to igrf.rst --- doc/user_guide/igrf.rst | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/doc/user_guide/igrf.rst b/doc/user_guide/igrf.rst index eef99ee2f..ce0623ee0 100644 --- a/doc/user_guide/igrf.rst +++ b/doc/user_guide/igrf.rst @@ -8,6 +8,17 @@ magnetic field from 1900 until 5 years after the latest release (based on predictions of the secular variation). Harmonica allows calculating the 14th generation IGRF field with :class:`harmonica.IGRF14`. Here's how it works. +.. jupyter-execute:: + :hide-code: + + import pygmt + + # Needed so that displaying works on jupyter-sphinx and sphinx-gallery at + # the same time. Using PYGMT_USE_EXTERNAL_DISPLAY="false" in the Makefile + # for sphinx-gallery to work means that fig.show won't display anything here + # either. + pygmt.set_display(method="notebook") + .. jupyter-execute:: import datetime From e6d05c8020693617f1d5c43b231249625ff36c03 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 10:00:07 +0100 Subject: [PATCH 05/10] remove quotes from and update and standardize figure titles and labels --- doc/gallery_src/equivalent_sources/spherical.py | 8 +++++--- doc/gallery_src/forward/prism_layer.py | 4 ++-- doc/gallery_src/forward/tesseroid.py | 2 +- doc/gallery_src/forward/tesseroid_layer.py | 6 +++++- .../forward/tesseroid_variable_density.py | 2 +- doc/gallery_src/transformations/reduction_to_pole.py | 12 +++++++++--- doc/gallery_src/transformations/tga.py | 12 +++++++++--- doc/gallery_src/transformations/tilt.py | 6 ++---- .../transformations/upward_continuation.py | 12 +++++++++--- doc/gallery_src/transformations/upward_derivative.py | 12 +++++++++--- doc/user_guide/forward_modelling/point.rst | 3 +-- doc/user_guide/gravity_disturbance.rst | 6 +++--- doc/user_guide/topographic_correction.rst | 2 +- 13 files changed, 57 insertions(+), 30 deletions(-) diff --git a/doc/gallery_src/equivalent_sources/spherical.py b/doc/gallery_src/equivalent_sources/spherical.py index 255329625..e940f271b 100644 --- a/doc/gallery_src/equivalent_sources/spherical.py +++ b/doc/gallery_src/equivalent_sources/spherical.py @@ -94,9 +94,9 @@ title = "Observed gravity disturbance data" fig.plot( - projection="M10c", + projection="M12c", region=region, - frame=["WSne", "xa5", "ya4"], + frame=[f"WSne+t{title}", "xa5", "ya4"], x=longitude, y=latitude, fill=gravity_disturbance, @@ -108,8 +108,10 @@ fig.shift_origin(xshift="w+3c") +title = "Gridded and upward-continued" + fig.grdimage( - frame=["ESnw", "xa5", "ya4"], + frame=[f"ESnw+t{title}", "xa5", "ya4"], grid=grid.gravity_disturbance, cmap=True, ) diff --git a/doc/gallery_src/forward/prism_layer.py b/doc/gallery_src/forward/prism_layer.py index 4cf887b73..22a39ba9c 100644 --- a/doc/gallery_src/forward/prism_layer.py +++ b/doc/gallery_src/forward/prism_layer.py @@ -51,14 +51,14 @@ # Plot gravity field fig = pygmt.Figure() -title = "Gravitational acceleration of a layer of prisms" +title = "Gravitational acceleration of topography with prisms" with pygmt.config(FONT_TITLE="14p"): fig.grdimage( region=region, projection="X10c/10c", grid=grid.gravity, - frame=["a", f"+t{title}", 'x+l"easting (m)"', 'y+l"northing (m)"'], + frame=["a", f"+t{title}", "x+leasting (m)", "y+lnorthing (m)"], cmap="viridis", ) diff --git a/doc/gallery_src/forward/tesseroid.py b/doc/gallery_src/forward/tesseroid.py index 9e805598b..65f240740 100644 --- a/doc/gallery_src/forward/tesseroid.py +++ b/doc/gallery_src/forward/tesseroid.py @@ -50,7 +50,7 @@ # Plot the gravitational field fig = pygmt.Figure() -title = "Downward component of gravitational acceleration" +title = "Gravitational acceleration of a tesseroid" with pygmt.config(FONT_TITLE="16p"): fig.grdimage( diff --git a/doc/gallery_src/forward/tesseroid_layer.py b/doc/gallery_src/forward/tesseroid_layer.py index 5e42984c6..a760cbcb4 100644 --- a/doc/gallery_src/forward/tesseroid_layer.py +++ b/doc/gallery_src/forward/tesseroid_layer.py @@ -61,14 +61,18 @@ # Make colormap of data pygmt.makecpt(cmap="balance+h0", series=[-cpt_lims, cpt_lims]) + +title = "Gravitational acceleration of topography with tesseroids" + fig.grdimage( gravity.g_z, + frame=f"+t{title}", projection="M15c", nan_transparent=True, cmap=True, ) fig.basemap(frame=True) -fig.colorbar(frame='af+l"Gravity [mGal]"', position="JCR") +fig.colorbar(frame="af+lmGal", position="JCR") fig.coast(shorelines="0.5p,black", borders=["1/0.5p,black"]) fig.show() diff --git a/doc/gallery_src/forward/tesseroid_variable_density.py b/doc/gallery_src/forward/tesseroid_variable_density.py index 16be3382d..7225c1e64 100644 --- a/doc/gallery_src/forward/tesseroid_variable_density.py +++ b/doc/gallery_src/forward/tesseroid_variable_density.py @@ -78,7 +78,7 @@ def density(radius): # Plot the gravitational field fig = pygmt.Figure() -title = "Downward component of gravitational acceleration" +title = "Gravitational acceleration of variable density tesseroids" with pygmt.config(FONT_TITLE="16p"): fig.grdimage( diff --git a/doc/gallery_src/transformations/reduction_to_pole.py b/doc/gallery_src/transformations/reduction_to_pole.py index d0e5944d2..8666b984f 100644 --- a/doc/gallery_src/transformations/reduction_to_pole.py +++ b/doc/gallery_src/transformations/reduction_to_pole.py @@ -61,18 +61,24 @@ grid=magnetic_grid, projection="X?", cmap=True, + frame="+tMagnetic anomaly" ) # Add colorbar fig.colorbar( - frame='af+l"Magnetic anomaly [nT]"', + frame="af+lnT", position="JBC+h+o0/1c+ef", ) with fig.set_panel(panel=1): # Plot upward reduced to the pole grid - fig.grdimage(grid=rtp_grid, projection="X?", cmap=True) + fig.grdimage( + grid=rtp_grid, + projection="X?", + cmap=True, + frame="+tReduced to the pole", + ) # Add colorbar fig.colorbar( - frame='af+l"Reduced to the pole [nT]"', + frame="af+lnT", position="JBC+h+o0/1c+ef", ) fig.show() diff --git a/doc/gallery_src/transformations/tga.py b/doc/gallery_src/transformations/tga.py index d5fe8396e..049eae1df 100644 --- a/doc/gallery_src/transformations/tga.py +++ b/doc/gallery_src/transformations/tga.py @@ -52,10 +52,11 @@ grid=magnetic_grid, projection="X?", cmap=True, + frame="+tMagnetic Anomaly", ) # Add colorbar fig.colorbar( - frame='af+l"Magnetic anomaly [nT]"', + frame="af+lnT", position="JBC+h+o0/1c", ) with fig.set_panel(panel=1): @@ -63,10 +64,15 @@ cpt_lim = 0.6 * vd.maxabs(tga) pygmt.makecpt(cmap="balance+h0", series=[0, cpt_lim], background=True) # Plot total gradient amplitude - fig.grdimage(grid=tga, projection="X?", cmap=True) + fig.grdimage( + grid=tga, + projection="X?", + cmap=True, + frame="+tTotal Gradient Amplitude", + ) # Add colorbar fig.colorbar( - frame='af+l"Total Gradient Amplitude [nT/m]"', + frame="af+lnT/m", position="JBC+h+o0/1c+ef", ) fig.show() diff --git a/doc/gallery_src/transformations/tilt.py b/doc/gallery_src/transformations/tilt.py index e62234dd7..96dad14df 100644 --- a/doc/gallery_src/transformations/tilt.py +++ b/doc/gallery_src/transformations/tilt.py @@ -101,9 +101,8 @@ frame=["a", "+tReduced to the pole (RTP)"], ) # Add colorbar - label = "nT" fig.colorbar( - frame=f"af+l{label}", + frame="af+lnT", position="JMR+o1/-0.25c+ef", ) @@ -129,9 +128,8 @@ frame=["a", "+tTilt of RTP grid"], ) # Add colorbar - label = "rad" fig.colorbar( - frame=f"af+l{label}", + frame="af+lradians", position="JMR+o1/-0.25c", ) fig.show() diff --git a/doc/gallery_src/transformations/upward_continuation.py b/doc/gallery_src/transformations/upward_continuation.py index edd76893d..84ca04f2c 100644 --- a/doc/gallery_src/transformations/upward_continuation.py +++ b/doc/gallery_src/transformations/upward_continuation.py @@ -54,18 +54,24 @@ grid=magnetic_grid, projection="X?", cmap=True, + frame="+tMagnetic Anomaly at 500m", ) # Add colorbar fig.colorbar( - frame='af+l"Magnetic anomaly at 500m [nT]"', + frame="af+lnT", position="JBC+h+o0/1c+e", ) with fig.set_panel(panel=1): # Plot upward continued grid - fig.grdimage(grid=upward_continued, projection="X?", cmap=True) + fig.grdimage( + grid=upward_continued, + projection="X?", + frame="+tUpward continued to 1000m", + cmap=True, + ) # Add colorbar fig.colorbar( - frame='af+l"Upward continued to 1000m [nT]"', + frame="af+lnT", position="JBC+h+o0/1c", ) fig.show() diff --git a/doc/gallery_src/transformations/upward_derivative.py b/doc/gallery_src/transformations/upward_derivative.py index 30c74618a..a77856dfb 100644 --- a/doc/gallery_src/transformations/upward_derivative.py +++ b/doc/gallery_src/transformations/upward_derivative.py @@ -53,10 +53,11 @@ grid=magnetic_grid, projection="X?", cmap=True, + frame="+tMagnetic Anomaly", ) # Add colorbar fig.colorbar( - frame='af+l"Magnetic anomaly [nT]"', + frame="af+lnT", position="JBC+h+o0/1c+e", ) with fig.set_panel(panel=1): @@ -64,10 +65,15 @@ cpt_lim = vd.maxabs(deriv_upward) * 0.6 pygmt.makecpt(cmap="balance+h0", series=[-cpt_lim, cpt_lim], background=True) # Plot upward derivative - fig.grdimage(grid=deriv_upward, projection="X?", cmap=True) + fig.grdimage( + grid=deriv_upward, + projection="X?", + cmap=True, + frame="+tUpward Derivative", + ) # Add colorbar fig.colorbar( - frame='af+l"Upward derivative [nT/m]"', + frame="af+lnT/m", position="JBC+h+o0/1c+e", ) fig.show() diff --git a/doc/user_guide/forward_modelling/point.rst b/doc/user_guide/forward_modelling/point.rst index 8f6999474..cf03b5178 100644 --- a/doc/user_guide/forward_modelling/point.rst +++ b/doc/user_guide/forward_modelling/point.rst @@ -238,7 +238,6 @@ Lets plot these results using :mod:`pygmt`: coordinates_spherical, g_z, data_names="g_z", extra_coords_names="extra") fig = pygmt.Figure() - title = "Gravitational acceleration (downward)" maxabs = vd.maxabs(g_z) pygmt.makecpt(cmap="balance+h0", series=[-maxabs, maxabs], background=True) @@ -246,7 +245,7 @@ Lets plot these results using :mod:`pygmt`: region=(-72, -68, -46, -42), projection="M10c", grid=grid.g_z, - frame=[f"WSne+t{title}", "x", "y"], + frame=["a", "x", "y"], cmap=True,) fig.colorbar(cmap=True, position="JMR", frame=["a0.000000005", "x+lmGal"]) diff --git a/doc/user_guide/gravity_disturbance.rst b/doc/user_guide/gravity_disturbance.rst index 6aef78354..a8130a56d 100644 --- a/doc/user_guide/gravity_disturbance.rst +++ b/doc/user_guide/gravity_disturbance.rst @@ -89,7 +89,7 @@ Lets plot it: fig.basemap(frame=["af", "WEsn"]) fig.colorbar( position="JCB+w10c", - frame=["af", 'y+l"mGal"', 'x+l"observed gravity"'], + frame=["af", 'y+lmGal', 'x+lObserved gravity'], ) fig.coast(shorelines=True, resolution="c", area_thresh=1e4) fig.show() @@ -122,7 +122,7 @@ And plot it: fig.basemap(frame=["af", "WEsn"]) fig.colorbar( position="JCB+w10c", - frame=["af", 'y+l"mGal"', 'x+l"normal gravity"'], + frame=["af", 'y+lmGal', 'x+lNormal gravity'], ) fig.coast(shorelines=True, resolution="c", area_thresh=1e4) fig.show() @@ -153,7 +153,7 @@ And plot it: fig.basemap(frame=["af", "WEsn"]) fig.colorbar( position="JCB+w10c", - frame=["af", 'y+l"mGal"', 'x+l"gravity disturbance"'], + frame=["af", 'y+lmGal', 'x+lGravity disturbance'], ) fig.coast(shorelines=True, resolution="c", area_thresh=1e4) fig.show() diff --git a/doc/user_guide/topographic_correction.rst b/doc/user_guide/topographic_correction.rst index 0e1fb0080..7e980e9db 100644 --- a/doc/user_guide/topographic_correction.rst +++ b/doc/user_guide/topographic_correction.rst @@ -82,7 +82,7 @@ And plot it: projection="M15c", frame=['ag', 'WSen'], ) - fig.colorbar(cmap=True, frame=["a50f25", "x+lgravity disturbance", "y+lmGal"]) + fig.colorbar(cmap=True, frame=["a50f25", "x+lGravity disturbance", "y+lmGal"]) fig.show() From 725169b57c3693b52164f5ace669c849d126a465 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 10:00:47 +0100 Subject: [PATCH 06/10] standardize forward calculation example colorbars --- doc/gallery_src/forward/tesseroid.py | 2 +- doc/gallery_src/forward/tesseroid_variable_density.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/gallery_src/forward/tesseroid.py b/doc/gallery_src/forward/tesseroid.py index 65f240740..d18a50fc7 100644 --- a/doc/gallery_src/forward/tesseroid.py +++ b/doc/gallery_src/forward/tesseroid.py @@ -61,7 +61,7 @@ cmap="viridis", ) -fig.colorbar(cmap=True, frame=["a200f50", "x+lmGal"]) +fig.colorbar(cmap=True, position="JMR", frame=["a200f50", "x+lmGal"]) fig.coast(shorelines="1p,black") diff --git a/doc/gallery_src/forward/tesseroid_variable_density.py b/doc/gallery_src/forward/tesseroid_variable_density.py index 7225c1e64..3142e3f0f 100644 --- a/doc/gallery_src/forward/tesseroid_variable_density.py +++ b/doc/gallery_src/forward/tesseroid_variable_density.py @@ -89,7 +89,7 @@ def density(radius): cmap="viridis", ) -fig.colorbar(cmap=True, frame=["a200f50", "x+lmGal"]) +fig.colorbar(cmap=True, position="JMR", frame=["a200f50", "x+lmGal"]) fig.coast(shorelines="1p,black") From 1694600215117b946d341566586c7d4b4349706d Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 10:01:03 +0100 Subject: [PATCH 07/10] add labels to point masses in plot --- doc/gallery_src/forward/point_gravity.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/doc/gallery_src/forward/point_gravity.py b/doc/gallery_src/forward/point_gravity.py index b822188b5..9c5be7555 100644 --- a/doc/gallery_src/forward/point_gravity.py +++ b/doc/gallery_src/forward/point_gravity.py @@ -68,7 +68,14 @@ cmap=True, ) -fig.plot(x=easting, y=northing, style="c0.2c", fill="grey") +fig.plot( + x=easting, + y=northing, + style="c0.15c", + fill="white", + pen="1p,black", + label="Point masses", +) fig.colorbar(cmap=True, position="JMR+e", frame=["x+lmGal"]) From 0bfeb31912efa183e23de99270f01607b0ba3f0f Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 10:01:18 +0100 Subject: [PATCH 08/10] make nans transparent in plot --- doc/gallery_src/equivalent_sources/spherical.py | 1 + 1 file changed, 1 insertion(+) diff --git a/doc/gallery_src/equivalent_sources/spherical.py b/doc/gallery_src/equivalent_sources/spherical.py index e940f271b..81551a8fe 100644 --- a/doc/gallery_src/equivalent_sources/spherical.py +++ b/doc/gallery_src/equivalent_sources/spherical.py @@ -114,6 +114,7 @@ frame=[f"ESnw+t{title}", "xa5", "ya4"], grid=grid.gravity_disturbance, cmap=True, + nan_transparent=True, ) fig.colorbar(cmap=True, frame=["x+lmGal"], position="+e") From 1cf65809e43bdce0352eef4c5651f9fd0060a614 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Wed, 19 Nov 2025 14:17:10 +0100 Subject: [PATCH 09/10] format fix --- doc/gallery_src/transformations/reduction_to_pole.py | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/doc/gallery_src/transformations/reduction_to_pole.py b/doc/gallery_src/transformations/reduction_to_pole.py index 8666b984f..0741eb28c 100644 --- a/doc/gallery_src/transformations/reduction_to_pole.py +++ b/doc/gallery_src/transformations/reduction_to_pole.py @@ -58,10 +58,7 @@ with fig.set_panel(panel=0): # Plot magnetic anomaly grid fig.grdimage( - grid=magnetic_grid, - projection="X?", - cmap=True, - frame="+tMagnetic anomaly" + grid=magnetic_grid, projection="X?", cmap=True, frame="+tMagnetic anomaly" ) # Add colorbar fig.colorbar( From b839ef1fa87f88bd05f02cd319e4e45a138e5505 Mon Sep 17 00:00:00 2001 From: mdtanker Date: Fri, 28 Nov 2025 12:32:05 +0100 Subject: [PATCH 10/10] remove print statement to keep html rendering of dataframe --- doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst index e7f9fcb3d..2053c2dc5 100644 --- a/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst +++ b/doc/user_guide/equivalent_sources/eqs-parameters-estimation.rst @@ -39,7 +39,7 @@ Lets start by loading some sample gravity data: fname = ensaio.fetch_bushveld_gravity(version=1) data = pd.read_csv(fname) - print(data) + data And project their coordinates using a Mercator projection: