<jupyter_start><jupyter_text>Get Real Event Data<jupyter_code>from gnam.events.gevents import gevents as gevents
from gnam.events.gstations import gstations as gstations
import matplotlib.pyplot as plt
import numpy as np
import pickle

# Unpickle
f = open('./model_bbox.pickle', 'rb')
gf_bbox = pickle.load(f)
f.close()
print('bbox after pickle :\n',gf_bbox)

f = open('./enorth.pickle', 'rb')
enorth = pickle.load(f)
f.close()

f = open('./gevents_test2.pickle', 'rb')
dill_snorth = pickle.load(f)
f.close()


#create streams
e_st3 = dill_snorth.getStreamZ(3,3)
e_st4 = dill_snorth.getStreamZ(3,4)

e_st = e_st3 + e_st4

#e_st.filter('bandpass',freqmin=1.0,freqmax=18.0,corners=4)
#e_st.plot()<jupyter_output><empty_output><jupyter_text>Better plotting of traces<jupyter_code>from gnam.events.gevents import gevents as gevents
from gnam.events.gstations import gstations as gstations
import matplotlib.pyplot as plt
import numpy as np
import pickle
import copy

dc_st = copy.deepcopy(e_st) 
dc_st2 = copy.deepcopy(e_st) 
stlen = len(e_st)
for i in range(stlen):
    '''
    fig, ax = plt.subplots(3,figsize=(15,4.8))
    #e_st[i].plot()
    #e_st[i].filter('bandpass',freqmin=1.0,freqmax=18.0,corners=4)
    #e_st[i].plot()
    '''
    print('delta: ', e_st[i].stats.delta)
    print('sampr: ', e_st[i].stats.sampling_rate)
    print('startt:', e_st[i].stats.starttime)
    print('endt:  ', e_st[i].stats.endtime)
    print('len(data):', len(e_st[i].data))
    raw = e_st[i].data
    t = np.arange(len(raw))
    dc_st[i].filter('bandpass',freqmin=1.0,freqmax=18.0,corners=4)
    dc_st2[i].filter('bandpass',freqmin=1.0,freqmax=18.0,corners=4,zerophase=True)
    filt = dc_st[i].data
    filt2 = dc_st2[i].data
    '''
    ax[0].plot(t, raw)
    ax[1].plot(t, filt)
    ax[2].plot(t, filt2)
    plt.show()
    '''

'''
fig, ax = plt.subplots()

signal = s[0].data
t = np.arange(signal.size, dtype = 'float')/200.0

ax.plot(t, signal)
'''<jupyter_output><empty_output><jupyter_text>Read SPECFEM Syn Data and Store into Stream<jupyter_code>import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import fnmatch
import os
import numpy as np
import pandas as pd

#fdir = 'semd_data'
#fdir = 'dyncm_data'
fdir_nQ = 'syn_nQ_data'
fdir_wQ = 'syn_wQ_data'
file_list_nQ = []
file_list_wQ = []
for file in os.listdir(fdir_nQ):
    if fnmatch.fnmatch(file, '*semd'):
        file_list_nQ.append(file)
        
for file in os.listdir(fdir_wQ):
    if fnmatch.fnmatch(file, '*semd'):
        file_list_wQ.append(file)
        

evnt3 = enorth.getIncCatalog()[3]
etime = evnt3.origins[0].time
print('etime:',etime)

syn_st_nQ = Stream()
syn_st_wQ = Stream()
for l in [[file_list_nQ,fdir_nQ,syn_st_nQ],[file_list_wQ,fdir_wQ,syn_st_wQ]]:
    file_list = l[0]
    fdir      = l[1]
    syn_st    = l[2]
    for f in file_list:
        fqn = fdir + '/' + f
        df = pd.io.parsers.read_csv(fqn,sep="\s+",header=None, usecols=[0,1])
        data = df[[1]].to_numpy().astype(np.float32).flatten()
        fhdr = f.split('.')
        # Fill header attributes
        stats = {'network': fhdr[0], 'station': fhdr[1], 'location': '',
                 'channel': fhdr[2], 'npts': len(data), 'delta': 0.001}
        # set current time
        #stats['starttime'] = UTCDateTime()
        syntime = df[[0]].to_numpy().astype(np.float64).flatten()
        #print('syntime:\n',syntime)
        print('starttime:        ',etime)
        stats['starttime'] = etime + syntime[0]
        print('starttime shifted:',stats['starttime'])
        st = Stream([Trace(data=data, header=stats)])
        print('Stream Before Resample:\n', st)
        print('endtime shifted:  ',st[0].stats.endtime)
        print('sampling:',st[0].stats.sampling_rate)
        st[0].resample(200)
        print('resampling:',st[0].stats.sampling_rate)
        print('Stream After Resample:\n', st)
        print('deltal:  ',st[0].stats.delta)
        print('resamp starttime:  ',st[0].stats.starttime)
        print('resamp endtime:    ',st[0].stats.endtime)
        print('resamp len(data)   ',len(st[0].data))
        endtime = st[0].stats.endtime
        st[0] = st[0].slice(etime,endtime)
        print('sliced starttime:  ',st[0].stats.starttime)
        print('sliced endtime:    ',st[0].stats.endtime)
        print('sliced len(data)   ',len(st[0].data))
        syn_st += st
    
orig_gpu_st_nQ = copy.deepcopy(syn_st_nQ)
orig_gpu_st_wQ = copy.deepcopy(syn_st_wQ)
#syn_st.filter('bandpass',freqmin=1.0,freqmax=18.0,corners=4)
#syn_st.plot()<jupyter_output><empty_output><jupyter_text>Register Syn and Real Data<jupyter_code>aresame = len(syn_st_nQ) == len(e_st)
print(aresame)
aresame = len(syn_st_wQ) == len(e_st)
print(aresame)

# sort stations
sort_syn_st_nQ = Stream()
sort_syn_st_wQ = Stream()
for i in range(len(e_st)):
    for j in range(len(syn_st_nQ)):
        if syn_st_nQ[j].stats.station == e_st[i].stats.station:
            sort_syn_st_nQ += syn_st_nQ[j]
            break
            
for i in range(len(e_st)):
    for j in range(len(syn_st_wQ)):
        if syn_st_wQ[j].stats.station == e_st[i].stats.station:
            sort_syn_st_wQ += syn_st_wQ[j]
            break

syn_st_nQ = sort_syn_st_nQ
syn_st_wQ = sort_syn_st_wQ
for i in range(len(e_st)):
    for syn_st in [syn_st_nQ,syn_st_wQ]:
        print('syn.station,real.station = (%s,%s)' %(syn_st[i].stats.station,e_st[i].stats.station))
        print('syn.starttime  = %s' %(syn_st[i].stats.starttime))
        print('real.starttime = %s' %(e_st[i].stats.starttime))
        print('syn.endtime    = %s' %(syn_st[i].stats.endtime))
        print('real.endtime   = %s' %(e_st[i].stats.endtime))
        print('syn.delta      = %s' %(syn_st[i].stats.delta))
        print('real.delta     = %s' %(e_st[i].stats.delta))
        print('syn.len        = %s' %(len(syn_st[i].data)))
        print('real.len       = %s' %(len(e_st[i].data)))
        start_t = e_st[i].stats.starttime
        end_t   = syn_st[i].stats.endtime
        e_st[i] = e_st[i].slice(etime,endtime)
        print('sliced syn.starttime  = %s' %(syn_st[i].stats.starttime))
        print('sliced real.starttime = %s' %(e_st[i].stats.starttime))
        print('sliced syn.endtime    = %s' %(syn_st[i].stats.endtime))
        print('sliced real.endtime   = %s' %(e_st[i].stats.endtime))
        print('sliced syn.delta      = %s' %(syn_st[i].stats.delta))
        print('sliced real.delta     = %s' %(e_st[i].stats.delta))
        print('sliced syn.len        = %s' %(len(syn_st[i].data)))
        print('sliced real.len       = %s' %(len(e_st[i].data)))<jupyter_output><empty_output><jupyter_text>Plot and Compare Syn vs. Real Data<jupyter_code>import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import matplotlib.pyplot as plt
import numpy as np
import copy 

dc_real = copy.deepcopy(e_st)
dc_syn_nQ  = copy.deepcopy(syn_st_nQ)
dc_syn_wQ  = copy.deepcopy(syn_st_wQ)
plt_scale = int(1*len(dc_real))
#plt_scale = int(3*len(dc_real))
#plt_scale = int(2*len(dc_real))
#fig, ax = plt.subplots(2,figsize=(14,4.8))
#fig, ax = plt.subplots(plt_scale,figsize=(14,2*plt_scale))
fig, ax = plt.subplots(plt_scale,figsize=(14,3*plt_scale))
fig.tight_layout()
for i in range(len(dc_real)):
    dc_real[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    #dc_syn[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    dc_syn_nQ[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    dc_syn_wQ[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    real_data = dc_real[i].data
    #syn_data  = dc_syn[i].data
    syn_data_nQ  = dc_syn_nQ[i].data
    syn_data_wQ  = dc_syn_wQ[i].data
    #star_t = dc_real[i].stats.starttime
    dt     = dc_real[i].stats.delta
    nt     = len(real_data)
    t      = dt*np.arange(nt)
    #ax[i*3].plot(t, real_data,c='black')
    #ax[i*3+1].plot(t, syn_data_nQ,c='blue')
    #ax[i*3+2].plot(t, syn_data_wQ,c='orange')
    ax[i].plot(t, real_data,c='black')
    ax[i].plot(t, syn_data_nQ,c='blue')
    ax[i].plot(t, syn_data_wQ,c='orange')
    #ax[i].plot(t, syn_data_wQ,c='orange',linestyle='dashed')
    rstation = dc_real[i].stats.station
    #sstation = dc_syn[i].stats.station
    sstation_nQ = dc_syn_nQ[i].stats.station
    sstation_wQ = dc_syn_wQ[i].stats.station
    overlay_title = 'Real(black), Syn_nQ(blue), Syn_wQ(orange):' + str(rstation)
    '''
    r_title_str = 'Real Station:' + str(rstation)
    s_title_str_nQ = 'Syn_nQ Station:' + str(sstation_nQ)
    s_title_str_wQ = 'Syn_wQ Station:' + str(sstation_wQ)
    ax[i*3].set_title(r_title_str)
    ax[i*3+1].set_title(s_title_str_nQ)
    ax[i*3+2].set_title(s_title_str_wQ)
    '''
    ax[i].set_title(overlay_title)
plt.show()
'''

#dc_syn_nQ.write('nam_nQ.sac', format='SAC') 
#dc_syn_wQ.write('nam_wQ.sac', format='SAC') 
dc_syn_nQ.write('./output/nam_nQ.semd', format='SACXY') 
dc_syn_wQ.write('./output/nam_wQ.semd', format='SACXY') 
'''
fig.savefig("overlay_syn_noQ_wQ2.pdf", bbox_inches='tight')
#fig.savefig("compare_real_syn_noQ_wQ2.pdf", bbox_inches='tight')<jupyter_output><empty_output><jupyter_text>Forcing Function Plot<jupyter_code>import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import matplotlib.pyplot as plt
import numpy as np
import copy 
import pandas as pd

fqn = './plot_source_time_function.txt'
df = pd.io.parsers.read_csv(fqn,sep="\s+",header=None, usecols=[0,1])
print('df:\n',df)
fdata = df[[1]].to_numpy().astype(np.float32).flatten()
tdata = df[[0]].to_numpy().astype(np.float32).flatten() + 0.01
print('fdata:\n',fdata[:15])
print('tdata:\n',tdata[:15])


fig, ax = plt.subplots()
ax.plot(tdata[:30],fdata[:30])
plt.show()<jupyter_output><empty_output><jupyter_text>Read CPU Syn and Create Stream<jupyter_code>import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import fnmatch
import os
import numpy as np
import pandas as pd

fdir = 'cpu_data'
file_list = []
for file in os.listdir(fdir):
    if fnmatch.fnmatch(file, '*semd'):
        file_list.append(file)
        

evnt3 = enorth.getIncCatalog()[3]
etime = evnt3.origins[0].time
print('etime:',etime)

cpu_st = Stream()
for f in file_list:
    fqn = fdir + '/' + f
    df = pd.io.parsers.read_csv(fqn,sep="\s+",header=None, usecols=[0,1])
    data = df[[1]].to_numpy().astype(np.float32).flatten()
    fhdr = f.split('.')
    # Fill header attributes
    stats = {'network': fhdr[0], 'station': fhdr[1], 'location': '',
             'channel': fhdr[2], 'npts': len(data), 'delta': 0.001}
    # set current time
    #stats['starttime'] = UTCDateTime()
    syntime = df[[0]].to_numpy().astype(np.float64).flatten()
    #print('syntime:\n',syntime)
    print('starttime:        ',etime)
    stats['starttime'] = etime + syntime[0]
    print('starttime shifted:',stats['starttime'])
    st = Stream([Trace(data=data, header=stats)])
    print('endtime shifted:  ',st[0].stats.endtime)
    print('sampling:',st[0].stats.sampling_rate)
    st[0].resample(200)
    print('sampling:',st[0].stats.sampling_rate)
    print('deltal:  ',st[0].stats.delta)
    print('resamp starttime:  ',st[0].stats.starttime)
    print('resamp endtime:    ',st[0].stats.endtime)
    print('resamp len(data)   ',len(st[0].data))
    endtime = st[0].stats.endtime
    st[0] = st[0].slice(etime,endtime)
    print('sliced starttime:  ',st[0].stats.starttime)
    print('sliced endtime:    ',st[0].stats.endtime)
    print('sliced len(data)   ',len(st[0].data))
    cpu_st += st

orig_cpu_st = copy.deepcopy(cpu_st)<jupyter_output><empty_output><jupyter_text>Register CPU and Real Data<jupyter_code>aresame = len(cpu_st) == len(e_st)
print(aresame)

# sort stations
sort_cpu_st = Stream()
for i in range(len(e_st)):
    for j in range(len(cpu_st)):
        if cpu_st[j].stats.station == e_st[i].stats.station:
            sort_cpu_st += cpu_st[j]
            break

cpu_st = sort_cpu_st
for i in range(len(e_st)):
    print('cpu.station,real.station = (%s,%s)' %(cpu_st[i].stats.station,e_st[i].stats.station))
    print('cpu.starttime  = %s' %(cpu_st[i].stats.starttime))
    print('real.starttime = %s' %(e_st[i].stats.starttime))
    print('cpu.endtime    = %s' %(cpu_st[i].stats.endtime))
    print('real.endtime   = %s' %(e_st[i].stats.endtime))
    print('cpu.delta      = %s' %(cpu_st[i].stats.delta))
    print('real.delta     = %s' %(e_st[i].stats.delta))
    print('cpu.len        = %s' %(len(cpu_st[i].data)))
    print('real.len       = %s' %(len(e_st[i].data)))
    start_t = e_st[i].stats.starttime
    end_t   = cpu_st[i].stats.endtime
    e_st[i] = e_st[i].slice(etime,endtime)
    print('sliced cpu.starttime  = %s' %(cpu_st[i].stats.starttime))
    print('sliced real.starttime = %s' %(e_st[i].stats.starttime))
    print('sliced cpu.endtime    = %s' %(cpu_st[i].stats.endtime))
    print('sliced real.endtime   = %s' %(e_st[i].stats.endtime))
    print('sliced cpu.delta      = %s' %(cpu_st[i].stats.delta))
    print('sliced real.delta     = %s' %(e_st[i].stats.delta))
    print('sliced cpu.len        = %s' %(len(cpu_st[i].data)))
    print('sliced real.len       = %s' %(len(e_st[i].data)))<jupyter_output><empty_output><jupyter_text>Plot and Compare CPU-Syn vs Real Data<jupyter_code>import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import matplotlib.pyplot as plt
import numpy as np
import copy 

dc_real = copy.deepcopy(e_st)
dc_cpu  = copy.deepcopy(cpu_st)
plt_scale = int(2*len(dc_real))
#fig, ax = plt.subplots(2,figsize=(14,4.8))
fig, ax = plt.subplots(plt_scale,figsize=(14,2*plt_scale))
fig.tight_layout()
for i in range(len(dc_real)):
    dc_real[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    dc_cpu[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    real_data = dc_real[i].data
    cpu_data  = dc_cpu[i].data
    #star_t = dc_real[i].stats.starttime
    dt     = dc_real[i].stats.delta
    nt     = len(real_data)
    t      = dt*np.arange(nt)
    ax[i*2].plot(t, real_data,c='black')
    ax[i*2+1].plot(t, cpu_data,c='blue')
    rstation = dc_real[i].stats.station
    sstation = dc_cpu[i].stats.station
    r_title_str = 'Real Station:' + str(rstation)
    s_title_str = 'cpu  Station:' + str(sstation)
    ax[i*2].set_title(r_title_str)
    ax[i*2+1].set_title(s_title_str)
plt.show()<jupyter_output><empty_output><jupyter_text>CPU vs GPU data<jupyter_code>aresame = len(orig_cpu_st) == len(orig_gpu_st)
print(aresame)

# sort stations
gpu_cpu_st = Stream()
for i in range(len(orig_gpu_st)):
    for j in range(len(orig_cpu_st)):
        if orig_cpu_st[j].stats.station == orig_gpu_st[i].stats.station:
            gpu_cpu_st += orig_cpu_st[j]
            break

orig_cpu_st = gpu_cpu_st
for i in range(len(orig_gpu_st)):
    print('cpu.station,gpu.station = (%s,%s)' %(orig_cpu_st[i].stats.station,orig_gpu_st[i].stats.station))
    print('cpu.starttime  = %s' %(orig_cpu_st[i].stats.starttime))
    print('gpu.starttime = %s' %(orig_gpu_st[i].stats.starttime))
    print('cpu.endtime    = %s' %(orig_cpu_st[i].stats.endtime))
    print('gpu.endtime   = %s' %(orig_gpu_st[i].stats.endtime))
    print('cpu.delta      = %s' %(orig_cpu_st[i].stats.delta))
    print('gpu.delta     = %s' %(orig_gpu_st[i].stats.delta))
    print('cpu.len        = %s' %(len(orig_cpu_st[i].data)))
    print('gpu.len       = %s' %(len(orig_gpu_st[i].data)))
    start_t = orig_gpu_st[i].stats.starttime
    end_t   = orig_cpu_st[i].stats.endtime
    orig_gpu_st[i] = orig_gpu_st[i].slice(etime,endtime)
    print('sliced cpu.starttime  = %s' %(orig_cpu_st[i].stats.starttime))
    print('sliced gpu.starttime = %s' %(orig_gpu_st[i].stats.starttime))
    print('sliced cpu.endtime    = %s' %(orig_cpu_st[i].stats.endtime))
    print('sliced gpu.endtime   = %s' %(orig_gpu_st[i].stats.endtime))
    print('sliced cpu.delta      = %s' %(orig_cpu_st[i].stats.delta))
    print('sliced gpu.delta     = %s' %(orig_gpu_st[i].stats.delta))
    print('sliced cpu.len        = %s' %(len(orig_cpu_st[i].data)))
    print('sliced gpu.len       = %s' %(len(orig_gpu_st[i].data)))
import obspy
from obspy import Stream
from obspy import Trace
from obspy import UTCDateTime
import matplotlib.pyplot as plt
import numpy as np
import copy 

dc_real = copy.deepcopy(orig_gpu_st)
dc_cpu  = copy.deepcopy(orig_cpu_st)
plt_scale = int(2*len(dc_real))
#fig, ax = plt.subplots(2,figsize=(14,4.8))
fig, ax = plt.subplots(plt_scale,figsize=(14,2*plt_scale))
fig.tight_layout()
for i in range(len(dc_real)):
    dc_real[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    dc_cpu[i].filter('bandpass',freqmin=2.0,freqmax=6.0,corners=4,zerophase=True)
    real_data = dc_real[i].data
    cpu_data  = dc_cpu[i].data
    #star_t = dc_real[i].stats.starttime
    dt     = dc_real[i].stats.delta
    nt     = len(real_data)
    t      = dt*np.arange(nt)
    ax[i*2].plot(t, real_data,c='black')
    ax[i*2+1].plot(t, cpu_data,c='blue')
    rstation = dc_real[i].stats.station
    sstation = dc_cpu[i].stats.station
    r_title_str = 'gpu Station:' + str(rstation)
    s_title_str = 'cpu  Station:' + str(sstation)
    ax[i*2].set_title(r_title_str)
    ax[i*2+1].set_title(s_title_str)
plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.17347587402484832, 'fraction_numerical': 0.014388904940768565, 'mean_word_length': 4.007523148148148, 'pattern_counts': {'":': 0, '<': 41, '<?xml version=': 0, '>': 41, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '415834', 'license': 'Apache-2.0', 'n_tokens_mistral': 7312, 'n_tokens_neox': 7201, 'n_words': 1083, 'path': 'notebooks/Compare Synthetic and Real Data.ipynb', 'repo_name': 'code-cullison/gnam', 'size': '17305'}