Hi Thibaut, here is the code. the observation data format is like below(I can’t upload the data from here).
10m multinodes wenner
10.00
1
573
1
0
15 10 31.4
30 20 31.0
45 30 30.1
60 40 27.6
75 50 24.8
90 60 22.5
105 70 21.9
120 80 19.4
135 90 16.8
25 10 31.2
40 20 31.6
35 10 29.7
50 20 30.8
65 30 29.9
80 40 29.1
95 50 26.0
155 90 18.6
45 10 29.0
60 20 29.8
120 60 22.1
55 10 27.8
70 20 29.0
85 30 29.9
100 40 28.9
.
.
.
If you need this data from me, let me know, I’ll send it to you through other means.
read in observation data
contents = np.genfromtxt(data_filename, skip_header = 1, delimiter=’ \n’, dtype=np.str)
n_sources = int(contents[2].split()[0]) # it is at the 3rd line fist element
n_topo = int(contents[579].split()[0]) # the 579th line 1st element
x_locations = np.zeros(n_sources) # zero vector
a_spacings = np.zeros(n_sources)
apparent_resistivity_values = np.zeros(n_sources)
a_locations = np.zeros(n_sources)
b_locations = np.zeros(n_sources) # zero vector
m_locations = np.zeros(n_sources)
n_locations = np.zeros(n_sources)
observed_data =np.zeros(n_sources)
Horizontal_xs = np.zeros(n_topo)
Vertical_ys = np.zeros(n_topo)
content_index = 4
loop over sources at certain range in this case, from 5th row to 577th row
for i in range(0,573): # range(4,576)
start by reading in the source info
content_index = content_index + 1 # read the next line
x_location, a_spacing, apparent_resistivity_value = contents[content_index].split() # this is a string(read the 5th line values and assign them to parameters)
convert the strings to a int for locations, ‘a’ spacing and float for apparent resisitivty values
x_locations[i] = int(x_location) # mid point location for wenner geometry
a_spacings[i] = int(a_spacing)
apparent_resistivity_values[i] = float(apparent_resistivity_value)
a_electrodes = x_locations - (a_spacings/2 + a_spacings)
b_electrodes = x_locations + (a_spacings/2 + a_spacings)
m_electrodes = x_locations - (a_spacings/2)
n_electrodes = x_locations + (a_spacings/2)
d_obs = apparent_resistivity_values
convert to the UBC format of wenner
a_electrodes = np.vstack([a_electrodes, np.zeros_like(a_electrodes), np.zeros_like(a_electrodes)]).T
b_electrodes = np.vstack([b_electrodes, np.zeros_like(b_electrodes), np.zeros_like(b_electrodes)]).T
m_electrodes = np.vstack([m_electrodes, np.zeros_like(m_electrodes), np.zeros_like(m_electrodes)]).T
n_electrodes = np.vstack([n_electrodes, np.zeros_like(n_electrodes), np.zeros_like(n_electrodes)]).T
Define survey
unique_tx, k = np.unique(np.c_[a_electrodes, b_electrodes], axis=0, return_index=True)
k= np.r_[k,len(a_electrodes)+1] # k = 572, len(a_electrodes) = 573
source_list = []
for ii in range(0,n_sources):
# MN electrode locations for receivers. Each is an (N, 3) numpy array
M_locations = m_electrodes[k[ii]:k[ii+1],:]
N_locations = n_electrodes[k[ii]:k[ii+1],:]
receiver_list = [dc.receivers.Dipole(M_locations, N_locations)]
# AB electrode locations for source. Each is a (1, 3) numpy array
A_location = a_electrodes[k[ii],:]
B_location = b_electrodes[k[ii],:]
source_list.append(dc.sources.Dipole(receiver_list, A_location, B_location))
survey = dc.survey.Survey_ky(source_list)
Define the a data object. Uncertainties are added later
dc_data = Data(survey=survey,dobs= d_obs)
plot psuedosection
mpl.rcParams.update({“font.size”: 12})
fig = plt.figure(figsize=(12, 5))
ax1 = fig.add_axes([0.05, 0.05, 0.8, 0.9])
dc.utils.plot_pseudoSection(
dc_data,
ax=ax1,
survey_type=“wenner”,
data_type=“appResistivity”,
space_type=“half-space”,
scale=“log”,
pcolorOpts={“cmap”: “viridis”},
)
ax1.set_title(“Apparent Resistivity [S/m]”)
plt.show()
Thanks!
