.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/generated/scripts_synth_data/single_survey_line.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_generated_scripts_synth_data_single_survey_line.py: Invert a single survey line for a thin plate ============================================ .. raw:: html .. GENERATED FROM PYTHON SOURCE LINES 13-18 |Open In Colab| .. |Open In Colab| image:: https://img.shields.io/badge/open%20in-Colab-b5e2fa?logo=googlecolab&style=flat-square&color=ffd670 :target: https://colab.research.google.com/github/inlab-geo/cofi-examples/blob/main/examples/airborne_em/airborne_em_single_line_transmitters.ipynb .. GENERATED FROM PYTHON SOURCE LINES 21-40 .. raw:: html .. raw:: html .. If you are running this notebook locally, make sure you’ve followed `steps here `__ to set up the environment. (This `environment.yml `__ file specifies a list of packages required to run the notebooks) .. GENERATED FROM PYTHON SOURCE LINES 40-52 .. code-block:: Python # -------------------------------------------------------- # # # # Uncomment below to set up environment on "colab" # # # # -------------------------------------------------------- # # !pip install -U cofi # !pip install git+https://github.com/JuergHauser/PyP223.git # !git clone https://github.com/inlab-geo/cofi-examples.git # %cd cofi-examples/examples/vtem_ma .. GENERATED FROM PYTHON SOURCE LINES 54-59 .. code-block:: Python # If this notebook is run locally PyP223 needs to be installed separately by uncommenting the following line, # that is by removing the # and the white space between it and the exclamation mark. # !pip install git+https://github.com/JuergHauser/PyP223.git .. GENERATED FROM PYTHON SOURCE LINES 61-84 .. code-block:: Python import numpy import matplotlib.pyplot as plt import matplotlib.ticker as ticker from matplotlib.lines import Line2D import cofi from vtem_max_forward_lib import ( problem_setup, system_spec, survey_setup, true_model, ForwardWrapper, plot_transient, plot_predicted_profile, plot_observed_profile, plot_plate_faces, plot_plate_faces_single ) numpy.random.seed(42) .. GENERATED FROM PYTHON SOURCE LINES 89-100 Background ---------- This example inverts a single survey line of VTEM max data using the inline and vertical component for a thin plate target. The dip azimuth of the plate is in the direction of the survey line and we only seek to recover the plate dip, the plate width, that is the extent of the plate from the plate reference point in the direction of dip and the location of the plate along the survey line. The forward problem and model setup is described `here <./thin_plate_inversion.ipynb>`__. .. GENERATED FROM PYTHON SOURCE LINES 103-114 Problem definition ------------------ Here we create a survey that consists of 12 fiducials forming a survey line in the direction of the dip azimuth of the thin plate. This allows us to infer more than just a single plate parameter. We are esstially inferring a 2D model and can only hope to directly recover parameters of the thin plate that map dirrectly into the the 2D mdoel. That is we can constrain its depth, position along the flight direction, dip and width, and they become our declared model paramters that are exposed to CoFI. .. GENERATED FROM PYTHON SOURCE LINES 114-135 .. code-block:: Python tx_min = 115 tx_max = 281 tx_interval = 15 n_transmitters = (tx_max - tx_min - 1) // tx_interval + 1 tx = numpy.arange(tx_min, tx_max, tx_interval) survey_setup = { "tx": tx, # transmitter easting/x-position "ty": numpy.array([100]*n_transmitters), # transmitter northing/y-position "tz": numpy.array([50]*n_transmitters), # transmitter height/z-position "tazi": numpy.deg2rad(numpy.array([90]*n_transmitters)), # transmitter azimuth "tincl": numpy.deg2rad(numpy.array([6]*n_transmitters)), # transmitter inclination "rx": tx, # receiver easting/x-position "ry": numpy.array([100]*n_transmitters), # receiver northing/y-position "rz": numpy.array([50]*n_transmitters), # receiver height/z-position "trdx": numpy.array([0]*n_transmitters), # transmitter receiver separation inline "trdy": numpy.array([0]*n_transmitters), # transmitter receiver separation crossline "trdz": numpy.array([0]*n_transmitters), # transmitter receiver separation vertical } .. GENERATED FROM PYTHON SOURCE LINES 137-141 .. code-block:: Python forward = ForwardWrapper(true_model, problem_setup, system_spec, survey_setup, ["pdip", "peast", "ptop", "pwdth2"]) .. rst-class:: sphx-glr-script-out .. code-block:: none ['pdip', 'peast', 'ptop', 'pwdth2'] .. GENERATED FROM PYTHON SOURCE LINES 143-147 .. code-block:: Python # check the order of parameters in a model vector forward.params_to_invert .. rst-class:: sphx-glr-script-out .. code-block:: none ['pdip', 'peast', 'ptop', 'pwdth2'] .. GENERATED FROM PYTHON SOURCE LINES 149-152 .. code-block:: Python true_param_value = numpy.array([60, 175, 30, 90]) .. GENERATED FROM PYTHON SOURCE LINES 157-160 True model ~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 160-177 .. code-block:: Python _, axes = plt.subplots(2, 2) axes[1,1].axis("off") plot_plate_faces( "plate_true", forward, true_param_value, axes[0,0], axes[0,1], axes[1,0], color="purple", label="True model" ) plt.tight_layout() point = Line2D([0], [0], label='Fiducial', marker='o', markersize=5, markeredgecolor='orange', markerfacecolor='orange', linestyle='') handles, labels = axes[1,0].get_legend_handles_labels() handles.extend([point]) axes[1,0].legend(handles=handles,bbox_to_anchor=(1.04, 0), loc="lower left") .. image-sg:: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_001.png :alt: single survey line :srcset: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 182-185 Generate synthetic data ~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 185-197 .. code-block:: Python # The data absolute_noise= 0.05 # create data and ad a realisation of the noise data_pred_true = forward(true_param_value) data_obs = data_pred_true + numpy.random.randn(len(data_pred_true))*absolute_noise # define data covariance matrix sigma=absolute_noise Cdinv=numpy.identity(len(data_obs))*(1.0/(sigma*sigma)) .. GENERATED FROM PYTHON SOURCE LINES 202-205 Parameter estimation -------------------- .. GENERATED FROM PYTHON SOURCE LINES 208-210 **Initialise a model for inversion** .. GENERATED FROM PYTHON SOURCE LINES 210-213 .. code-block:: Python init_param_value = numpy.array([50, 155, 35, 70]) .. GENERATED FROM PYTHON SOURCE LINES 218-220 **Define helper functions for CoFI** .. GENERATED FROM PYTHON SOURCE LINES 220-247 .. code-block:: Python def my_objective(model): dpred = forward(model) residual = dpred - data_obs return residual.T @ Cdinv @ residual def my_gradient(model): dpred = forward(model) jacobian = forward.jacobian(model, relative_step=0.1) residual = dpred - data_obs return jacobian.T @ Cdinv @ residual def my_hessian(model): jacobian = forward.jacobian(model) return jacobian.T @ Cdinv @ jacobian class PerIterationCallbackFunction: def __init__(self): self.x = None self.i = 0 def __call__(self, xk): print(f"Iteration #{self.i+1}") print(f" objective value: {my_problem.objective(xk)}") self.x = xk self.i += 1 .. GENERATED FROM PYTHON SOURCE LINES 252-254 **Define CoFI problem** .. GENERATED FROM PYTHON SOURCE LINES 254-261 .. code-block:: Python my_problem = cofi.BaseProblem() my_problem.set_objective(my_objective) my_problem.set_gradient(my_gradient) my_problem.set_hessian(my_hessian) my_problem.set_initial_model(init_param_value) .. GENERATED FROM PYTHON SOURCE LINES 266-268 **Define CoFI options** .. GENERATED FROM PYTHON SOURCE LINES 268-273 .. code-block:: Python my_options = cofi.InversionOptions() my_options.set_tool("scipy.optimize.minimize") my_options.set_params(method="Newton-CG", options={"maxiter": 10},callback=PerIterationCallbackFunction()) .. GENERATED FROM PYTHON SOURCE LINES 278-280 **Run CoFI inversion** .. GENERATED FROM PYTHON SOURCE LINES 280-285 .. code-block:: Python my_inversion = cofi.Inversion(my_problem, my_options) my_result = my_inversion.run() print(my_result.model) .. rst-class:: sphx-glr-script-out .. code-block:: none Iteration #1 objective value: 12965.94204951452 Iteration #2 objective value: 12507.956940690212 Iteration #3 objective value: 3547.6083610260457 Iteration #4 objective value: 2914.0976629697516 Iteration #5 objective value: 2760.5961826863527 Iteration #6 objective value: 2699.7730944662535 Iteration #7 objective value: 2695.053833161359 [ 57.51701555 169.1248781 16.63134734 63.14639711] .. GENERATED FROM PYTHON SOURCE LINES 290-293 Plotting -------- .. GENERATED FROM PYTHON SOURCE LINES 296-299 Data - Fiducials ~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 299-310 .. code-block:: Python _, (ax1, ax2) = plt.subplots(1, 2) plot_transient(true_param_value, forward, "Data from true model", ax1, ax2, color="purple") plot_transient(init_param_value, forward, "Data from init model", ax1, ax2, color="green", linestyle=":") plot_transient(my_result.model, forward, "Data from MAP model", ax1, ax2, color="red", linestyle="-.") ax1.legend(loc="upper center") ax2.legend(loc="upper center") ax1.set_title("Vertical") ax2.set_title("Inline") plt.tight_layout() .. image-sg:: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_002.png :alt: Vertical, Inline :srcset: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 315-318 Data - Profiles ~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 318-338 .. code-block:: Python idx_to_plot = numpy.arange(8, 20) # subset of numpy.arange(0, 44) _, (ax1 ,ax2) = plt.subplots(2, 1, figsize=(12,6)) x_labels = ax1.get_xticks() ax1.xaxis.set_major_formatter(ticker.FormatStrFormatter('%5f')) plot_predicted_profile(true_param_value, forward, "Data from true model", gate_idx=idx_to_plot, ax=ax1, cmp='vertical', color="purple") plot_predicted_profile(init_param_value, forward, "Data from starting model", gate_idx=idx_to_plot, ax=ax1,cmp='vertical', color="green", linestyle=":") plot_predicted_profile(my_result.model, forward, "Data from MAP model", gate_idx=idx_to_plot, ax=ax1, cmp='vertical',color="red", linestyle="-.") ax1.legend(bbox_to_anchor=(1.04, 0), loc="lower left") plot_predicted_profile(true_param_value, forward, "Data from true model", gate_idx=idx_to_plot, ax=ax2, cmp='inline', color="purple") plot_predicted_profile(init_param_value, forward, "Data from starting model", gate_idx=idx_to_plot, ax=ax2,cmp='inline', color="green", linestyle=":") plot_predicted_profile(my_result.model, forward, "Data from MAP model", gate_idx=idx_to_plot, ax=ax2, cmp='inline',color="red", linestyle="-.") ax2.legend(bbox_to_anchor=(1.04, 0), loc="lower left") plt.tight_layout() .. image-sg:: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_003.png :alt: single survey line :srcset: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 343-346 Model ~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 346-370 .. code-block:: Python _, axes = plt.subplots(2, 2) axes[1,1].axis("off") plot_plate_faces( "plate_true", forward, true_param_value, axes[0,0], axes[0,1], axes[1,0], color="purple", label="True model" ) plot_plate_faces( "plate_init", forward, init_param_value, axes[0,0], axes[0,1], axes[1,0], color="green", label="Starting model" ) plot_plate_faces( "plate_inverted", forward, my_result.model, axes[0,0], axes[0,1], axes[1,0], color="red", label="MAP solution", linestyle="dotted" ) plt.tight_layout() point = Line2D([0], [0], label='Fiducial', marker='o', markersize=5, markeredgecolor='orange', markerfacecolor='orange', linestyle='') handles, labels = axes[1,0].get_legend_handles_labels() handles.extend([point]) axes[1,0].legend(handles=handles,bbox_to_anchor=(1.04, 0), loc="lower left") .. image-sg:: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_004.png :alt: single survey line :srcset: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 372-379 .. code-block:: Python _, ax = plt.subplots(1, 1) plot_plate_faces_single("plate_true", "xz", forward, true_param_value, ax, color="purple", label="True model") plot_plate_faces_single("plate_init", "xz", forward, init_param_value, ax, color="green", label="Starting model") plot_plate_faces_single("plate_inverted", "xz", forward, my_result.model, ax, color="red", label="MAP solution", linestyle="dotted") ax.legend(); .. image-sg:: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_005.png :alt: single survey line :srcset: /examples/generated/scripts_synth_data/images/sphx_glr_single_survey_line_005.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 384-397 -------------- Watermark ========= .. raw:: html .. raw:: html .. GENERATED FROM PYTHON SOURCE LINES 397-403 .. code-block:: Python watermark_list = ["cofi", "numpy", "scipy", "matplotlib"] for pkg in watermark_list: pkg_var = __import__(pkg) print(pkg, getattr(pkg_var, "__version__")) .. rst-class:: sphx-glr-script-out .. code-block:: none cofi 0.2.11 numpy 2.3.5 scipy 1.17.0 matplotlib 3.10.8 .. GENERATED FROM PYTHON SOURCE LINES 404-404 sphinx_gallery_thumbnail_number = -1 .. rst-class:: sphx-glr-timing **Total running time of the script:** (1 minutes 48.656 seconds) .. _sphx_glr_download_examples_generated_scripts_synth_data_single_survey_line.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: single_survey_line.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: single_survey_line.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: single_survey_line.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_