.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "content/examples/07-nsem/plot_fwd_nsem_MTTipper3D.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_content_examples_07-nsem_plot_fwd_nsem_MTTipper3D.py: MT: 3D: Forward =============== Forward model 3D MT data. Test script to use SimPEG.NSEM platform to forward model impedance and tipper synthetic data. .. GENERATED FROM PYTHON SOURCE LINES 10-131 .. rst-class:: sphx-glr-horizontal * .. image-sg:: /content/examples/07-nsem/images/sphx_glr_plot_fwd_nsem_MTTipper3D_001.png :alt: Slice 15 :srcset: /content/examples/07-nsem/images/sphx_glr_plot_fwd_nsem_MTTipper3D_001.png :class: sphx-glr-multi-img * .. image-sg:: /content/examples/07-nsem/images/sphx_glr_plot_fwd_nsem_MTTipper3D_002.png :alt: plot fwd nsem MTTipper3D :srcset: /content/examples/07-nsem/images/sphx_glr_plot_fwd_nsem_MTTipper3D_002.png :class: sphx-glr-multi-img .. code-block:: Python import discretize from SimPEG.electromagnetics import natural_source as NSEM from SimPEG import utils import numpy as np import matplotlib.pyplot as plt try: from pymatsolver import Pardiso as Solver except ImportError: from SimPEG import Solver def run(plotIt=True): """ MT: 3D: Forward =============== Forward model 3D MT data. """ # Make a mesh M = discretize.TensorMesh( [ [(100, 9, -1.5), (100.0, 13), (100, 9, 1.5)], [(100, 9, -1.5), (100.0, 13), (100, 9, 1.5)], [(50, 10, -1.6), (50.0, 10), (50, 6, 2)], ], x0=["C", "C", -14926.8217], ) # Setup the model conds = [1, 1e-2] sig = utils.model_builder.create_block_in_wholespace( M.gridCC, [-100, -100, -350], [100, 100, -150], conds ) sig[M.gridCC[:, 2] > 0] = 1e-8 sig[M.gridCC[:, 2] < -1000] = 1e-1 sigBG = np.zeros(M.nC) + conds[1] sigBG[M.gridCC[:, 2] > 0] = 1e-8 if plotIt: collect_obj = M.plot_slice(np.log10(sig), grid=True, normal="X")[0] plt.colorbar(collect_obj) # Setup the the survey object # Receiver locations rx_x, rx_y = np.meshgrid(np.arange(-600, 601, 100), np.arange(-600, 601, 100)) rx_loc = np.hstack( (utils.mkvc(rx_x, 2), utils.mkvc(rx_y, 2), np.zeros((np.prod(rx_x.shape), 1))) ) # Make a receiver list receiver_list = [] for rx_orientation in ["xx", "xy", "yx", "yy"]: receiver_list.append(NSEM.Rx.PointNaturalSource(rx_loc, rx_orientation, "real")) receiver_list.append(NSEM.Rx.PointNaturalSource(rx_loc, rx_orientation, "imag")) for rx_orientation in ["zx", "zy"]: receiver_list.append(NSEM.Rx.Point3DTipper(rx_loc, rx_orientation, "real")) receiver_list.append(NSEM.Rx.Point3DTipper(rx_loc, rx_orientation, "imag")) # Source list source_list = [ NSEM.Src.PlanewaveXYPrimary(receiver_list, freq) for freq in np.logspace(4, -2, 13) ] # Survey MT survey = NSEM.Survey(source_list) # Setup the problem object problem = NSEM.Simulation3DPrimarySecondary( M, survey=survey, solver=Solver, sigma=sig, sigmaPrimary=sigBG, forward_only=True, ) # Calculate the data # data = problem.make_synthetic_data(relative_error=0.1, add_noise=True) data = NSEM.Data(survey=survey, dobs=problem.dpred()) # Add standard deviation to the data - 10% relative error and 0 floor data.relative_error = 0.1 data.noise_floor = 0.0 # Add plots if plotIt: # Plot the data # On and off diagonal (on left and right axis, respectively) fig, axes = plt.subplots(2, 1, figsize=(7, 5)) plt.subplots_adjust(right=0.8) [(ax.invert_xaxis(), ax.set_xscale("log")) for ax in axes] ax_r, ax_p = axes ax_r.set_yscale("log") ax_r.set_ylabel("Apparent resistivity [xy-yx]") ax_r_on = ax_r.twinx() ax_r_on.set_yscale("log") ax_r_on.set_ylabel("Apparent resistivity [xx-yy]") ax_p.set_ylabel("Apparent phase") ax_p.set_xlabel("Frequency [Hz]") # Start plotting ax_r = data.plot_app_res( np.array([-200, 0]), components=["xy", "yx"], ax=ax_r, errorbars=True ) ax_r_on = data.plot_app_res( np.array([-200, 0]), components=["xx", "yy"], ax=ax_r_on, errorbars=True ) ax_p = data.plot_app_phs( np.array([-200, 0]), components=["xx", "xy", "yx", "yy"], ax=ax_p, errorbars=True, ) ax_p.legend(bbox_to_anchor=(1.05, 1), loc=2) if __name__ == "__main__": do_plots = True run(do_plots) if do_plots: plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** (3 minutes 11.869 seconds) **Estimated memory usage:** 1017 MB .. _sphx_glr_download_content_examples_07-nsem_plot_fwd_nsem_MTTipper3D.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_fwd_nsem_MTTipper3D.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_fwd_nsem_MTTipper3D.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_