In [1]:
import math
import matplotlib.pyplot as plt
import numpy as np

In [39]:
### file EDI

fDir = 'C:\\Users\\NSS06\\Downloads\\Edi_File\\model ori\\AZ_054.edi'

In [40]:
def readData(component, fList, nFreq):
    readEdi = fList.readline().split()                                             # read first for detect number of column
#     print(readEdi)
    nRow = math.ceil(int(nFreq)/(len(readEdi)))                                    # calculate nRow                                      
    
    for i in range(0, len((readEdi))):
            component.append(float(readEdi[i]))
            
    for j in range(1, nRow):                                                       # read block data
        readEdi = fList.readline().split()
#         print(readEdi)
        for i in range(0, len((readEdi))):
            component.append(float(readEdi[i]))

    return component

In [41]:
fList = open(fDir, 'r')
readEdi = ' '

print(fDir)
while not readEdi == '>END':                                                                    # read until eof
    readEdi = fList.readline()
#     print(readEdi[0:5])

    if readEdi[0:5] == 'REFLO':                                                                # read station
        refLoc = str(readEdi[8:12])
    if readEdi[0:5] == 'REFLA':                                                                # read latitude
        refLat = readEdi[7:]
    if readEdi[0:5] == 'REFLO':                                                                # read longitude
        refLong = readEdi[8:]
    if readEdi[0:5] == 'REFEL':                                                                # read elevation
        refEle = readEdi[8:]
    if readEdi[0:5] == 'NFREQ':                                                                 # read number of freq
        nFreq = readEdi[6:]
        
    if readEdi[0:5] == '>FREQ': #______________________________________________________________ FREQ
        freq = []
        freq = readData(freq, fList, nFreq)
        
    if readEdi[0:5] == '>ZROT': #______________________________________________________________ ZROT
        zrot = []
        zrot = readData(zrot, fList, nFreq)
    
    if readEdi[0:5] == '>ZXXR': #______________________________________________________________ ZXX
        zxxr = []
        zxxr = readData(zxxr, fList, nFreq)
    
    if readEdi[0:5] == '>ZXXI': 
        zxxi = []
        zxxi = readData(zxxi, fList, nFreq)
    
    if readEdi[0:5] == '>ZXX.': 
        zxxv = []
        zxxv = readData(zxxv, fList, nFreq)
    
    if readEdi[0:5] == '>ZXYR': #______________________________________________________ _______ ZXY
        zxyr = []
        zxyr = readData(zxyr, fList, nFreq)
    
    if readEdi[0:5] == '>ZXYI': 
        zxyi = []
        zxyi = readData(zxyi, fList, nFreq)
    
    if readEdi[0:5] == '>ZXY.': 
        zxyv = []
        zxyv = readData(zxyv, fList, nFreq)
    
    if readEdi[0:5] == '>ZYXR': #______________________________________________________________ ZYX
        zyxr = []
        zyxr = readData(zyxr, fList, nFreq)
    
    if readEdi[0:5] == '>ZYXI': 
        zyxi = []
        zyxi = readData(zyxi, fList, nFreq)
    
    if readEdi[0:5] == '>ZYX.': 
        zyxv = []
        zyxv = readData(zyxv, fList, nFreq)
    
    if readEdi[0:5] == '>ZYYR': #______________________________________________________________ ZYY
        zyyr = []
        zyyr = readData(zyyr, fList, nFreq)
    
    if readEdi[0:5] == '>ZYYI': 
        zyyi = []
        zyyi = readData(zyyi, fList, nFreq)
    
    if readEdi[0:5] == '>ZYY.': 
        zyyv = []
        zyyv = readData(zyyv, fList, nFreq)
    
    if readEdi[0:5] == '>TROT': #______________________________________________________________ TROT
        trot = []
        trot = readData(trot, fList, nFreq)
        
    if readEdi[0:5] == '>TXR.': #______________________________________________________________ TX
        txr = []
        txr = readData(txr, fList, nFreq)
    
    if readEdi[0:5] == '>TXI.': 
        txi = []
        txi = readData(txi, fList, nFreq)
    
    if readEdi[0:5] == '>TXVA': 
        txv = []
        txv = readData(txv, fList, nFreq)
    
    if readEdi[0:5] == '>TYR.': #______________________________________________________________ TY
        tyr = []
        tyr = readData(tyr, fList, nFreq)
    
    if readEdi[0:5] == '>TYI.': 
        tyi = []
        tyi = readData(tyi, fList, nFreq)
    
    if readEdi[0:5] == '>TYVA': 
        tyv = []
        tyv = readData(tyv, fList, nFreq)


C:\Users\NSS06\Downloads\Edi_File\model ori\AZ_054.edi

In [69]:
plt.figure(1, figsize=(20, 30))

plt.subplot(4,1,1)
plt.semilogx(freq,zxxr,'ro', label='zxyr')
plt.semilogx(freq,zxxi,'bo', label='zxyi')
plt.semilogx(freq,zxxv,'go', label='zxyv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('ZXXR')
plt.legend()

# plt.subplot(4,3,2)
# plt.semilogx(freq,zxxi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZXXI')

# plt.subplot(4,3,3)
# plt.semilogx(freq,zxxv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZXXV')

plt.subplot(4,1,2)
plt.loglog(freq,zxyr,'ro', label='zxyr')
plt.loglog(freq,zxyi,'bo', label='zxyi')
plt.loglog(freq,zxyv,'go', label='zxyv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('ZXYR')
plt.legend()

# plt.subplot(4,3,5)
# plt.semilogx(freq,zxyi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZXYI')

# plt.subplot(4,3,6)
# plt.semilogx(freq,zxyv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZXYV')

plt.subplot(4,1,3)
plt.semilogx(freq,zyxr,'ro', label='zyxr')
plt.semilogx(freq,zyxi,'bo', label='zyxi')
plt.semilogx(freq,zyxv,'go', label='zyxv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('ZYXR')
plt.legend()

# plt.subplot(4,3,8)
# plt.semilogx(freq,zyxi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZYXI')

# plt.subplot(4,3,9)
# plt.semilogx(freq,zyxv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZYXV')

plt.subplot(4,1,4)
plt.semilogx(freq,zyyr,'ro', label='zyyr')
plt.semilogx(freq,zyyi,'bo', label='zyyi')
plt.semilogx(freq,zyyv,'go', label='zyyv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('ZYYR')
plt.legend()

# plt.subplot(4,3,11)
# plt.semilogx(freq,zyyi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZYYI')

# plt.subplot(4,3,12)
# plt.semilogx(freq,zyyv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('ZYYV')

plt.show()



In [ ]:


In [68]:
plt.figure(2, figsize=(20, 15))

plt.subplot(2,1,1)
plt.semilogx(freq,txr,'ro', label='txr')
plt.semilogx(freq,txi,'bo', label='txi')
plt.semilogx(freq,txv,'go', label='txv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('TXR')
plt.legend()

# plt.subplot(2,3,2)
# plt.semilogx(freq,txi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('TXI')

# plt.subplot(2,3,3)
# plt.semilogx(freq,txv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('TXV')

plt.subplot(2,1,2)
plt.semilogx(freq,tyr,'ro', label='tyr')
plt.semilogx(freq,tyi,'bo', label='tyi')
plt.semilogx(freq,tyv,'go', label='tyv')
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
# plt.ylabel('Probability')
plt.title('TYR')
plt.legend()

# plt.subplot(2,3,5)
# plt.semilogx(freq,tyi,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('TYI')

# plt.subplot(2,3,6)
# plt.semilogx(freq,tyv,'bo')
# plt.gca().invert_xaxis()
# plt.xlabel('Freq (Hz)')
# # plt.ylabel('Probability')
# plt.title('TYV')

plt.show()



In [56]:
RHOxy = []
RHOyx = []
PHASExy = []
PHASEyx = []

for i in range(0, int(nFreq)):  # komponen 0-60
    RHOxy.append((0.2*float(1/freq[i]))*(np.abs(complex((zxyr[i]), (zxyi[i]))))**2)
    RHOyx.append((0.2*float(1/freq[i]))*(np.abs(complex((zyxr[i]), (zyxi[i]))))**2)
    PHASExy.append(np.degrees(np.arctan(((float(zxyi[i]) / float(zxyr[i]))))))
    PHASEyx.append(np.degrees(np.arctan(((float(zyxi[i]) / float(zyxr[i]))))))

In [72]:
plt.figure(3, figsize=(20, 20))

plt.subplot(2,1,1)
plt.loglog(freq,RHOxy,'ro',label='xy')
plt.loglog(freq,RHOyx,'bo',label='yx')
plt.xlim([0.001, 100])
plt.ylim([1, 1000])
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
plt.ylabel('Ohm.meter')
plt.title('RHO')
plt.legend()

plt.subplot(2,1,2)
plt.semilogx(freq,PHASExy,'ro', label='xy')
plt.semilogx(freq,PHASEyx,'bo', label='yx')
plt.xlim([0.001, 100])
plt.ylim([0, 90])
# plt.ylim([0, 180])
# plt.ylim([0, 270])
# plt.ylim([-90, 90])
# plt.ylim([-180, 180])
# plt.ylim([-270, 270])
plt.gca().invert_xaxis()
plt.xlabel('Freq (Hz)')
plt.ylabel('Degree')
plt.title('PHASE')
plt.legend()

plt.show()



In [ ]:


In [ ]:


In [ ]: