Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ModEM_Plot Depthslice #170

Open
Michelle-12v opened this issue May 21, 2022 · 0 comments
Open

ModEM_Plot Depthslice #170

Michelle-12v opened this issue May 21, 2022 · 0 comments

Comments

@Michelle-12v
Copy link

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

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant