In [5]:
0.03*1264
Out[5]:
In [6]:
%pylab inline
In [7]:
import pandas as pd
In [8]:
data = pd.read_excel("../data/Kimberlite-2015-07-17.xls")
In [9]:
data_active = data.loc[np.logical_and((data['Facies'] == 'XVK')|(data['Facies'] == 'PK')|(data['Facies'] == 'HK')|(data['Facies'] == 'VK'), data.notnull()['Rinf']==True)][["0LabID (PCG)", "Facies", "Peregrine ID", "(Latitude)", "(Longitude)", "Depth (m)","Mag Susc [SI]","Resistivity [Ohm.m]","Sat Geometric Dens [g/cc]","Chargeability [ms]","Rinf","Rh","Qh","Ph","Rm","Qm","Pm","Rl","Ql","Pl","Re","Qe","Pe-f","Pe-i"]]
In [10]:
data_active
Out[10]:
In [11]:
PK = data.loc[np.logical_and((data['Facies'] == 'PK'), data.notnull()['Rinf']==True)][["0LabID (PCG)", "Facies", "Peregrine ID", "(Latitude)", "(Longitude)", "Depth (m)","Mag Susc [SI]","Resistivity [Ohm.m]","Sat Geometric Dens [g/cc]","Chargeability [ms]","Rinf","Rm","Qm","Pm","Rl","Ql","Pl"]]
HK = data.loc[np.logical_and((data['Facies'] == 'HK'), data.notnull()['Rinf']==True)][["0LabID (PCG)", "Facies", "Peregrine ID", "(Latitude)", "(Longitude)", "Depth (m)","Mag Susc [SI]","Resistivity [Ohm.m]","Sat Geometric Dens [g/cc]","Chargeability [ms]","Rinf","Rm","Qm","Pm","Rl","Ql","Pl"]]
In [12]:
HK
Out[12]:
In [13]:
PK
Out[13]:
In [27]:
def CPEfun(Rx, Qx, px, freq):
out = np.zeros_like(freq, dtype=complex128)
out = 1./(1./Rx + Qx*(np.pi*2*freq*1j)**px)
return out
def CPEfunElec(Rx, Qx, pex, pix, freq):
out = np.zeros_like(freq, dtype=complex128)
out = 1./(1./Rx + (1j)**pix*Qx*(np.pi*2*freq)**pex)
return out
def CPEfunSeries(Rx, Qx, px, freq):
out = np.zeros_like(freq, dtype=complex128)
out = Rx + 1./(Qx*(np.pi*2*freq*1j)**px)
return out
f0peak = lambda R, Q, P: (R*Q)**(-1./P)/np.pi/2.
taupeak = lambda R, Q, P: (R*Q)**(1./P)
In [28]:
def TKCColeColeElec(frequency, PID, data):
Rh, Qh, Ph = data[data["Peregrine ID"]==PID]['Rh'].values[0], data[data["Peregrine ID"]==PID]['Qh'].values[0], data[data["Peregrine ID"]==PID]['Ph'].values[0]
Rm, Qm, Pm = data[data["Peregrine ID"]==PID]['Rm'].values[0], data[data["Peregrine ID"]==PID]['Qm'].values[0], data[data["Peregrine ID"]==PID]['Pm'].values[0]
Rl, Ql, Pl = data[data["Peregrine ID"]==PID]['Rl'].values[0], data[data["Peregrine ID"]==PID]['Ql'].values[0], data[data["Peregrine ID"]==PID]['Pl'].values[0]
Re, Qe, Pe_f, Pe_i = data[data["Peregrine ID"]==PID]['Re'].values[0], data[data["Peregrine ID"]==PID]['Qe'].values[0], data[data["Peregrine ID"]==PID]['Pe-f'].values[0],data[data["Peregrine ID"]==PID]['Pe-i'].values[0]
rhoinf = data[data["Peregrine ID"]==PID]['Rinf'].values[0]
fpeakm = f0peak(Rm, Qm, Pm)
fpeakl = f0peak(Rl, Ql, Pl)
rhoh = CPEfun(Rh, Qh, Ph, frequency)
rhom = CPEfun(Rm, Qm, Pm, frequency)
rhol = CPEfun(Rl, Ql, Pl, frequency)
rhoe = CPEfunElec(Re, Qe, Pe_f, Pe_i, frequency)
rho = rhoh+rhom+rhol+rhoinf+rhoe
m = (rho.real[0]-rho.real[-1])/rho.real[0]
fig, ax = plt.subplots(3, 2, figsize = (15, 9))
ax[2,0].semilogx(frequency, rho.real, 'k-', lw=2)
ax[0,0].semilogx(frequency, rhom.real+rhoinf, 'b-', lw=2)
ax[1,0].semilogx(frequency, rhol.real+rhoinf, 'g-', lw=2)
ax2 = ax[0,0].twinx()
ax3 = ax[1,0].twinx()
ax4 = ax[2,0].twinx()
ax4.semilogx(frequency, (rho.imag), 'k--', lw=2)
ax2.semilogx(frequency, (rhom.imag), 'b--', lw=2)
ax3.semilogx(frequency, (rhol.imag), 'g--', lw=2)
ax2.semilogx(np.r_[fpeakm, fpeakm], np.r_[rhom.imag.min(), rhom.imag.max()], 'b', lw=2)
ax3.semilogx(np.r_[fpeakl, fpeakl], np.r_[rhol.imag.min(), rhol.imag.max()], 'g', lw=2)
ax2.invert_yaxis()
ax3.invert_yaxis()
ax4.invert_yaxis()
ax[0,0].grid(True)
ax[1,0].grid(True)
ax[2,0].grid(True)
ax[2,1].plot(rho.real, rho.imag, 'k-')
ax[0,1].plot((rhom+rhoinf).real, (rhom+rhoinf).imag, 'b-')
ax[1,1].plot((rhol+rhoinf).real, (rhol+rhoinf).imag, 'g-')
# ax[2,0].plot(temp[:,0], temp[:,4], 'k.')
# ax[2,1].plot(temp[:,4], temp[:,5], 'k.')
# ax4.plot(temp[:,0], temp[:,5], 'k.')
ax[0,1].invert_yaxis()
ax[0,1].grid(True)
ax[1,1].invert_yaxis()
ax[1,1].grid(True)
ax[2,1].invert_yaxis()
ax[2,1].grid(True)
# print data[data["0LabID (PCG)"]==PID]['Facies'].values[0], PID
print data[data["Peregrine ID"]==PID]['Facies'].values[0], PID
print "Rinf = ", rhoinf
print "Chargeability = ", m
print "Taum = ", 1./fpeakm
print "Taul = ", 1./fpeakl
return
In [29]:
PK
Out[29]:
In [30]:
# frequency = np.logspace(-2, 6, 211)
# TKCColeColeElec(frequency, "K1P-0825", data_active)
In [31]:
frequency = np.logspace(-2, 6, 211)
TKCColeColeElec(frequency, "K1P-0807", data_active)
In [32]:
frequency = np.logspace(-2, 6, 211)
TKCColeColeElec(frequency, "K2P-0017", data_active)
In [39]:
def TKCColeCole(frequency, PID, data):
Rm, Qm, Pm = data[data["Peregrine ID"]==PID]['Rm'].values[0], data[data["Peregrine ID"]==PID]['Qm'].values[0], data[data["Peregrine ID"]==PID]['Pm'].values[0]
Rl, Ql, Pl = data[data["Peregrine ID"]==PID]['Rl'].values[0], data[data["Peregrine ID"]==PID]['Ql'].values[0], data[data["Peregrine ID"]==PID]['Pl'].values[0]
rhoinf = data[data["Peregrine ID"]==PID]['Rinf'].values[0]
fpeakm = f0peak(Rm, Qm, Pm)
fpeakl = f0peak(Rl, Ql, Pl)
rhom = CPEfun(Rm, Qm, Pm, frequency)
rhol = CPEfun(Rl, Ql, Pl, frequency)
rho = rhom+rhol+rhoinf
m = (rho.real[0]-rho.real[-1])/rho.real[0]
fig, ax = plt.subplots(3, 2, figsize = (15, 9))
ax[2,0].semilogx(frequency, rho.real, 'k-', lw=2)
ax[0,0].semilogx(frequency, rhom.real+rhoinf, 'b-', lw=2)
ax[1,0].semilogx(frequency, rhol.real+rhoinf, 'g-', lw=2)
ax2 = ax[0,0].twinx()
ax3 = ax[1,0].twinx()
ax4 = ax[2,0].twinx()
ax4.semilogx(frequency, (rho.imag), 'k--', lw=2)
ax2.semilogx(frequency, (rhom.imag), 'b--', lw=2)
ax3.semilogx(frequency, (rhol.imag), 'g--', lw=2)
ax2.semilogx(np.r_[fpeakm, fpeakm], np.r_[rhom.imag.min(), rhom.imag.max()], 'b', lw=2)
ax3.semilogx(np.r_[fpeakl, fpeakl], np.r_[rhol.imag.min(), rhol.imag.max()], 'g', lw=2)
ax2.invert_yaxis()
ax3.invert_yaxis()
ax4.invert_yaxis()
ax[0,0].grid(True)
ax[1,0].grid(True)
ax[2,0].grid(True)
ax[2,1].plot(rho.real, rho.imag, 'k-')
ax[0,1].plot((rhom+rhoinf).real, (rhom+rhoinf).imag, 'b-')
ax[1,1].plot((rhol+rhoinf).real, (rhol+rhoinf).imag, 'g-')
ax[0,1].invert_yaxis()
ax[0,1].grid(True)
ax[1,1].invert_yaxis()
ax[1,1].grid(True)
ax[2,1].invert_yaxis()
ax[2,1].grid(True)
print data[data["Peregrine ID"]==PID]['Facies'].values[0], PID
print "Rinf = ", rhoinf
print "Chargeability = ", m
print "Taum = ", 1./fpeakm
print "Taul = ", 1./fpeakl
return
In [40]:
RlHK = data_active[data_active['Facies']=='HK']["Rl"].values[:]
QlHK = data_active[data_active['Facies']=='HK']["Ql"].values[:]
PlHK = data_active[data_active['Facies']=='HK']["Pl"].values[:]
fpeakHK = f0peak(RlHK, QlHK, PlHK)
In [41]:
RlPK = data_active[data_active['Facies']=='PK']["Rl"].values[:]
QlPK = data_active[data_active['Facies']=='PK']["Ql"].values[:]
PlPK = data_active[data_active['Facies']=='PK']["Pl"].values[:]
fpeakPK = f0peak(RlPK, QlPK, PlPK)
In [42]:
RlXVK = data_active[data_active['Facies']=='XVK']["Rl"].values[:]
QlXVK = data_active[data_active['Facies']=='XVK']["Ql"].values[:]
PlXVK = data_active[data_active['Facies']=='XVK']["Pl"].values[:]
fpeakXVK = f0peak(RlXVK, QlXVK, PlXVK)
In [43]:
RlVK = data_active[data_active['Facies']=='VK']["Rl"].values[:]
QlVK = data_active[data_active['Facies']=='VK']["Ql"].values[:]
PlVK = data_active[data_active['Facies']=='VK']["Pl"].values[:]
fpeakVK = f0peak(RlVK, QlVK, PlVK)
In [ ]:
In [44]:
TKCColeCole(frequency, "K1P-0535", data_active)
In [45]:
TKCColeCole(frequency, "K1P-0540", data_active)
In [47]:
TKCColeCole(frequency, "K1P-0591", data_active)
In [48]:
TKCColeCole(frequency, "K1P-0593", data_active)
In [50]:
TKCColeCole(frequency, "K1P-0595", data_active)
In [51]:
data_active[data_active['Facies']=='PK']
Out[51]:
In [52]:
TKCColeCole(frequency, "K1P-0807", data_active)
In [54]:
TKCColeCole(frequency, "K1P-0825", data_active)
In [55]:
data_active[data_active['Facies']=='XVK']
Out[55]:
In [56]:
TKCColeCole(frequency, "K2P-0017", data_active)
In [57]:
TKCColeCole(frequency, "K2P-0020", data_active)
In [58]:
TKCColeCole(frequency, "K2P-0024", data_active)
In [59]:
TKCColeCole(frequency, "K2P-0031", data_active)
In [60]:
TKCColeCole(frequency, "K2P-0077", data_active)
In [61]:
TKCColeCole(frequency, "K2P-0147", data_active)
In [62]:
TKCColeCole(frequency, "K2P-0157", data_active)
In [293]:
# Ro Rh Qh Ch Ph Fh pRh pQh pCh Rm Qm Cm Pm Fm pRm pQm pCm Rl Ql Cl Pl Fl pRl pQl pCl Re Qe Pe-f Pe-i
svec = ['h', 'm', 'l']
PK = {}
for s in svec:
datatypeIP = ["Depth (m)", "Resistivity [Ohm.m]", "Rinf", "Ro", "R"+s, "Q"+s, "P"+s, "F"+s, "pR"+s, "pQ"+s, "pC"+s]
# PK[s]=data.loc[data['Peregrine ID'] == 'K1P-0807'][datatypeIP]
PK[s]=data.loc[data['Peregrine ID'] == 'K1P-0825'][datatypeIP]
In [294]:
PK['l']
Out[294]:
In [295]:
PK['m']
Out[295]:
In [73]:
PK['h']
Out[73]:
In [74]:
PK['l'][datatypeIP[4:]].values[0,:]
Out[74]:
In [75]:
frequency = np.logspace(-6, 8, 211)
In [77]:
In [ ]:
In [81]:
In [83]:
print Pm, Pl
print 1./fpeakm, 1./fpeakl
print rho.real[0], rho.real[-1], m
In [67]:
rhoh = CPEfun(Rh, Qh, Ph, frequency)
In [68]:
matplotlib.rcParams.update({'font.size': 14})
In [69]:
In [156]:
test = data.loc[np.logical_and((data['Facies'] == 'XVK')|(data['Facies'] == 'PK')|(data['Facies'] == 'HK')|(data['Facies'] == 'VK'), data.notnull()['Rinf']==True)][["Facies", "(Latitude)", "(Longitude)", "Depth (m)","Mag Susc [SI]","Resistivity [Ohm.m]","Sat Geometric Dens [g/cc]"]]
In [157]:
test[test["Facies"]=="HK"]
Out[157]:
In [147]:
plt.semilogy(test[test["Facies"]=="XVK"]["Mag Susc [SI]"], test[test["Facies"]=="XVK"]["Resistivity [Ohm.m]"], 'ko')
plt.semilogy(test[test["Facies"]=="PK"]["Mag Susc [SI]"], test[test["Facies"]=="PK"]["Resistivity [Ohm.m]"], 'bo')
plt.semilogy(test[test["Facies"]=="HK"]["Mag Susc [SI]"], test[test["Facies"]=="HK"]["Resistivity [Ohm.m]"], 'ro')
plt.semilogy(test[test["Facies"]=="VK"]["Mag Susc [SI]"], test[test["Facies"]=="VK"]["Resistivity [Ohm.m]"], 'go')
plt.grid(True)
plt.legend(("XVK", "PK", "HK", "VK"), loc=4)
plt.xlabel("Susceptibility (SI)")
plt.ylabel("Resistivity (ohm-m)")
Out[147]:
In [150]:
plt.plot(test[test["Facies"]=="XVK"]["Mag Susc [SI]"], test[test["Facies"]=="XVK"]["Sat Geometric Dens [g/cc]"], 'ko')
plt.plot(test[test["Facies"]=="PK"]["Mag Susc [SI]"], test[test["Facies"]=="PK"]["Sat Geometric Dens [g/cc]"], 'bo')
plt.plot(test[test["Facies"]=="HK"]["Mag Susc [SI]"], test[test["Facies"]=="HK"]["Sat Geometric Dens [g/cc]"], 'ro')
plt.plot(test[test["Facies"]=="VK"]["Mag Susc [SI]"], test[test["Facies"]=="VK"]["Sat Geometric Dens [g/cc]"], 'go')
plt.grid(True)
plt.legend(("XVK", "PK", "HK", "VK"), loc=4)
plt.xlabel("Susceptibility (SI)")
plt.ylabel("Density (g/cc)")
Out[150]:
In [151]:
test
Out[151]:
In [115]:
plt.plot(test["Depth (m)"].values, test["Resistivity [Ohm.m]"].values, 'ko')
Out[115]:
In [96]:
print(data.loc[np.logical_and(data['Facies'] == 'XVK', data.notnull()['Rinf']==True)][["Depth (m)", "Peregrine ID", "Resistivity [Ohm.m]", "Rinf", "Chargeability [ms]"]])
In [97]:
print(data.loc[np.logical_and(data['Facies'] == 'VK', data.notnull()['Rinf']==True)][["Depth (m)", "Peregrine ID", "Resistivity [Ohm.m]", "Rinf", "Chargeability [ms]"]])
In [ ]:
In [ ]: