<jupyter_start><jupyter_text>This script converts the continuous mseed files of the days & the stations required to be tested into numpy arrays The output of this script will be numpy array of each station separately for one day log.  Each subsequent 30-second window starts 10 seconds after the start time of the preceding time-window. This makes each time sample represented three times.<jupyter_code>import obspy, os, shutil
import matplotlib.pyplot as plt
import numpy as np
import datetime
full_date='jan.2018.24'
print('The full date of the day is: ',full_date)
print(' ')
print(' ')
print('Defining some Functions...')
print(' ')
print(' ')
def unify_start_time(st):
    st.sort(['starttime'])       #sort function works in place
    start_time=st[-1].stats['starttime']
    st.trim(start_time)
    st.sort(['station'])
def unify_end_time(st):
    st.sort(['endtime'])
    end_time=st[0].stats['endtime']
    st.trim(endtime=end_time)
    st.sort(['station'])
def check_channel_arrangement(st):
    for num,tr in enumerate(st):
        if num%3==0:
            if not tr.stats.channel=='HH1':
                print(tr.stats.station,tr.stats.channel,tr.stats.starttime)
                print('number in stream: ', num)
                print('Check channel arrangement!!!')
        elif num%3==1:
            if not tr.stats.channel=='HH2':
                print(tr.stats.station,tr.stats.channel,tr.stats.starttime)
                print('number in stream: ', num)
                print('Check channel arrangement!!!')
        elif num%3==2:
            if not tr.stats.channel=='HHZ':
                print(tr.stats.station,tr.stats.channel,tr.stats.starttime)
                print('number in stream: ', num)
                print('Check channel arrangement!!!')
def check_borehole(st):
    if len(st)%12==0:

    
        for x in range(0,len(st),12):
            z=[y.stats.station[:-1] for y in st[x:x+12]]
            t=[str(y.stats.starttime) for y in st[x:x+12]]
            
            if len(set(z))>1:
                print
                print(st[x:x+12])
                print('Check borehole arrangement, sorting probelm')
            if len(set(t))>1:
                print(x)
                print('Check borehole arrangement, starttime problem')
    else:
        print('They are not divisible by 12, check input files')
print('The current path directory is : ',os.getcwd())

folder='Dec_Jan_2018_records'
print(' ')
print(' ')
print('The folder where the mseed data exists is : ', folder)
print(' ')
print(' ')
month=datetime.datetime.strptime(full_date,'%b.%Y.%d').month
#print(month)
year=datetime.datetime.strptime(full_date,'%b.%Y.%d').year
#print(year)
day=datetime.datetime.strptime(full_date,'%b.%Y.%d').day
#print(day)
month_str=full_date.split('.')[0]
#print(month_str)

st=obspy.Stream()
for root,dirs,files in os.walk('../Dec_Jan_2018_records',topdown=True):
    #print(root)
    try:
        st+=obspy.read(os.path.join(root,'*{}*'.format(full_date)))
    except:
        pass

print(' ')
print(' ')
print(' The date of the day is: ')
print(' ')
print(' ')
print('year: ',year,' -- month: ',month_str,' -- day: ',day)
print(' ')
print(' ')
print('Stations loaded..  ')
print(' ')
print(' ')
print(st.__str__(extended=True))
#obspy.UTCDateTime(year,month,day,23,59,58)+datetime.timedelta(days=10)<jupyter_output><empty_output><jupyter_text>day=1st=obspy.read('../May_2019_mseed_records/may_{}/*'.format(day))stst.sort(['station','channel'])print(st.__str__(extended=True))<jupyter_code>x=list()
for tr in st:
    if tr.stats.starttime>(obspy.UTCDateTime(year,month,day,23,59,58)-datetime.timedelta(days=1)) or tr.stats.endtime<obspy.UTCDateTime(year,month,day,0,1,0)+datetime.timedelta(days=1):
        #print(tr.stats.station)
        x.append(tr.stats.station[:-1])
x=list(set(x))
print('Stations that do not have continuous records are: ')
print(x)
for bore in x:
    for tr in st:
        if bore==tr.stats.station[:-1]:
            st.remove(tr)

    
st.sort(['station','channel'])
print('Trimming remaining records to day limits ')
st.trim(obspy.UTCDateTime(year,month,day,23,59,58.5)-datetime.timedelta(days=1),
        obspy.UTCDateTime(year,month,day,0,0,1.5)+datetime.timedelta(days=1))


st.sort(['station','channel','starttime'])
print('Data after preparation...  ')
print('  ')
print(st.__str__(extended=True))


check_borehole(st)
check_channel_arrangement(st)
borehole={}
for x in range(0,len(st),12):
    borehole[st[x].stats.station[:-1]]=st[x:x+12]
(borehole)
print('Stations that will be processed: ')
for key in borehole.keys():
    print(key)
# THIS IS IF YOU WANT EACH STATION IN A SEPARATE MATRIX
print(' ')
print('Starting processing.. ')

min_freq=5.0
max_freq=25.0
corner=4
phase=False

window_length=30
taper_length=window_length*0.05

twl=window_length+taper_length*2
stp=10

npts_of_a_trace=3001
print('tapper_length ', taper_length)
print('output window length', window_length )
print('total sliding window length ', twl)




arr = np.empty((0,4,npts_of_a_trace,3))
test_dict={}

#time_matrix=[]
#test_data_dict={}
#test_data_dict['starttime']=list(=)

for key,st in borehole.items():
    print(key)

    arr = np.empty((0,4,npts_of_a_trace,3),dtype='float32')
    for st2 in st.slide(window_length=twl, step=stp,include_partial_windows=False):

        st_slide=st2.copy()

        st_slide.detrend('demean')
        st_slide.detrend('linear')
        st_slide.taper(type='hann',max_percentage=0.05,max_length=taper_length)
        st_slide.filter('bandpass', freqmin= min_freq ,freqmax=max_freq, corners=corner, zerophase=phase)
        st_slide.decimate(2)
        st_slide.trim(starttime=st_slide[0].stats.starttime+taper_length,
                      endtime=st_slide[0].stats.endtime-taper_length)
        #print(set([str(x.stats.starttime) for x in st_slide]))
        #time_matrix.append(st_slide[0].stats.starttime)
        tr_data=[tr.data for tr in st_slide]
        #print(len(tr_data))
        for x in range(0,len(tr_data)):
            if tr_data[x].shape[0]!=1:
                tr_data[x]=tr_data[x].reshape((1,1,tr_data[x].shape[0],1))
            else:

                print('error in tr_data[x].shape[0]', tr_data[x],tr_data[x].shape)
        thre_chan=[]
        for x in range(0,len(tr_data),3):
            thre_chan.append(np.concatenate((tr_data[x],tr_data[x+1],tr_data[x+2]),axis=3))
        #print('length of three channels: ',len(thre_chan))                
        for x in range(0,len(thre_chan)):

            abs_max=np.max(np.absolute(thre_chan[x]))
            thre_chan[x]= thre_chan[x]/abs_max
            if np.max(np.absolute(thre_chan[x]))!=1:
                print('normalization error')
        borehole_matrix=[]

        for x in range(0,len(thre_chan),4):
            borehole_matrix.append(np.concatenate((thre_chan[x],thre_chan[x+1],thre_chan[x+2],thre_chan[x+3]),axis=1))

        #print('length of the borehole: ',len(borehole_matrix))  
        if borehole_matrix[0].shape[2]==npts_of_a_trace:
            arr=np.concatenate((arr,borehole_matrix[0]),axis=0)
            if arr.shape[0]%1000==0:
                print('array shape update',arr.shape)



            #test_data_dict['starttime'].append([trace.stats.starttime for trace in st_slide])

        else:
            print('This shape is not compatible: ',borehole_matrix[0].shape[2])
    test_dict['{}'.format(key)]=arr
    if not os.path.isdir('./cont_data_numpy_records/{}_{}'.format(month_str,day)):
        os.mkdir('./cont_data_numpy_records/{}_{}'.format(month_str,day))
    np.save('./cont_data_numpy_records/{}_{}/{}_{}_{}_30_sec'.format(month_str,day,month_str,day,key),test_dict[key])
    del test_dict[key]
test_dict<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.16024314536989137, 'fraction_numerical': 0.020822555613036732, 'mean_word_length': 3.5595518867924527, 'pattern_counts': {'":': 0, '<': 10, '<?xml version=': 0, '>': 12, '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': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '306720', 'license': 'MIT', 'n_tokens_mistral': 2797, 'n_tokens_neox': 2698, 'n_words': 474, 'path': 'scripts/2_continuous_data_conversion_into_numpy/cont_data_preparation.ipynb', 'repo_name': 'AhmedShaheen01/GroningenNet', 'size': '7732'}