MTgeophysics / mtpy

Python toolbox for standard Magnetotelluric (MT) data analysis
GNU General Public License v3.0
147 stars 66 forks source link

ModEM_Plot Depthslice #170

Open Michelle-12v opened 2 years ago

Michelle-12v commented 2 years ago

Dear all, I tried to use the ModEM_PlotDepthslice for the example data. Also, I tried another data. the script works well but the color scale limits for resistivity does not appear

my script is:

import os import matplotlib.pyplot as plt from matplotlib import colors import numpy as np from mtpy.modeling.modem import PlotSlices

working directory and save path

wd = r'C:/mtpywin/mtpy/examples/model_files/ModEM_2' savepath = r'C:/tmp'

model and data file names

model_fn = os.path.join(wd,'Modular_MPI_NLCG_004.rho') data_fn = os.path.join(wd,'ModEM_Data.dat') fs = 8 # fontsize on plot

create a PlotSlices object

ps = PlotSlices(model_fn=model_fn, data_fn=data_fn, save_path=wd, plot_yn='n')

create a new figure

fig = plt.figure(figsize=[6,3]) ax = plt.subplot(111)

get resistivity data along a slice defined by the stations

gd, gz, gv = ps.get_slice("STA", nsteps=200 # total number of cells across )

to get resistivity along an arbitrary profile

(in this case stations 5 to 9)

xy locations

coords = np.vstack([ps.station_east[5:10],ps.station_north[5:10]]).T

gd, gz, gv = ps.get_slice("XY", nsteps=50, coords=coords)

ci = ax.pcolormesh(gd, gz, gv, vmin=1,vmax=1e4, # min & max resistivity on colour map norm=colors.LogNorm(), # plot colours on a log scale cmap='bwr_r', # colour map rasterized=True)

set some plot parameters and add labels

ax.invert_yaxis() ax.set_aspect(1) ax.set_ylim(1e2) plt.setp(ax.get_xticklabels(),fontsize=fs) plt.setp(ax.get_yticklabels(),fontsize=fs) plt.xlabel('Distance, km',fontsize=fs) plt.ylabel('Depth, km',fontsize=fs) plt.savefig(os.path.join(savepath,'DepthSlice.png'), dpi=400) # change to your desired figure resolution