In [2]:
import os, sys, glob, re
import datetime as dt
import numpy as np
from calendar import monthrange
from matplotlib.dates import date2num, num2date
# import h5py
from sklearn import decomposition
import nimfa

sys.path.insert(0,'/home/wu-jung/code_git/mi-instrument')
from mi.instrument.kut.ek60.ooicore.zplsc_b import *
# from concat_raw import get_num_days_pings, get_data_from_h5
# from echogram_decomp import find_nearest_time_idx,reshape_into_3freq,reshape_into_1freq,\
#     sep_into_freq,plot_decomp_v,plot_decomp_transform
from db_diff import *

import matplotlib.pyplot as plt
# from modest_image import imshow
# from mpl_toolkits.axes_grid1 import make_axes_locatable


/home/wu-jung/anaconda3/envs/py27/lib/python2.7/site-packages/matplotlib/__init__.py:1401: UserWarning:  This call to matplotlib.use() has no effect
because the backend has already been chosen;
matplotlib.use() must be called *before* pylab, matplotlib.pyplot,
or matplotlib.backends is imported for the first time.

  warnings.warn(_use_error_msg)

In [1]:
%matplotlib inline

In [3]:
date_start = dt.date(2015,8,20)
date_end = dt.date(2015,9,20)
date_wanted = [(date_start+dt.timedelta(days=xx)).strftime('%Y%m%d') for xx in range((date_end-date_start).days + 1)]

In [16]:
data_path = '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/'
save_path = '/media/wu-jung/wjlee_apl_2/ooi_zplsc_new/'
save_fname = '%s-%s.h5' % (date_start.strftime('%Y%m%d'),date_end.strftime('%Y%m%d'))

In [5]:
# Set param
ping_bin_range = 40
depth_bin_range = 10
tvg_correction_factor = 2

hour_all = range(24)
min_all = range(20)
sec_all = range(0,60,5)
ping_per_day = len(hour_all)*len(min_all)*len(sec_all)

In [6]:
# Set up array
particle_data, data_times, power_data_dict, freq, bin_size, config_header, config_transducer = \
    parse_echogram_file(glob.glob(os.path.join(data_path,'OOI-D%s*.raw' %date_wanted[0]))[0])


2017-10-14 19:34:10,809 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150820-T000000.raw'

In [8]:
Sv_raw_all = np.ma.empty((len(power_data_dict),power_data_dict[1].shape[0],\
                          len(date_wanted)*ping_per_day))
Sv_corr_all = np.ma.empty((len(power_data_dict),power_data_dict[1].shape[0],\
                           len(date_wanted)*ping_per_day))
Sv_noise_all = np.ma.empty((len(power_data_dict),power_data_dict[1].shape[0],len(date_wanted)))
time_all = np.ma.empty((len(power_data_dict),len(date_wanted)*ping_per_day))

In [ ]:
for iD,dd in zip(range(len(date_wanted)),date_wanted):

    # load files and get calibration params
    fname = glob.glob(os.path.join(data_path,'OOI-D%s*.raw' %dd))[0]
    particle_data, data_times, power_data, freq, bin_size, config_header, config_transducer = \
        parse_echogram_file(fname)
    cal_params = get_cal_params(power_data,particle_data,config_header,config_transducer)

    # swap sequence of 120 kHz and 38 kHz cal_params and data
    cal_params = [cal_params[fn] for fn in [1,0,2]]
    power_data = get_power_data_mtx(power_data,freq)

    # clean data
    num_freq = len(power_data)
    Sv_raw = np.ma.empty((power_data.shape))
    Sv_corr = np.ma.empty((power_data.shape))
    Sv_noise = np.ma.empty((power_data.shape[0:2]))
    for fn in range(num_freq):
        noise_est,ping_bin_num = get_noise(power_data[fn,:,:],bin_size,ping_bin_range,depth_bin_range)
        Sv_raw[fn,:,:],Sv_corr[fn,:,:],Sv_noise[fn,:] = \
                remove_noise(power_data[fn,:,:],cal_params[fn],noise_est.min(),ping_bin_range,tvg_correction_factor)

    # set up indexing to get wanted pings
    dd = dt.datetime.strptime(dd,'%Y%m%d')
    time_wanted = [dt.datetime(dd.year,dd.month,dd.day,hh,mm,ss) for hh in hour_all for mm in min_all for ss in sec_all]
    idx_wanted = [find_nearest_time_idx(data_times,tt,2) for tt in time_wanted]
    notnanidx = np.argwhere(~np.isnan(idx_wanted)).flatten()
    notnanidx_in_all = np.array(idx_wanted)[notnanidx].astype(int)

    # save wanted pings
    idx_save = iD*ping_per_day + notnanidx
    idx_save_mask = iD*ping_per_day + np.argwhere(np.isnan(idx_wanted))
    Sv_raw_all[:,:,idx_save] = Sv_raw[:,:,notnanidx_in_all]
    Sv_raw_all[:,:,idx_save_mask] = np.ma.masked
    Sv_corr_all[:,:,idx_save] = Sv_corr[:,:,notnanidx_in_all]
    Sv_corr_all[:,:,idx_save_mask] = np.ma.masked
    Sv_noise_all[:,:,iD] = Sv_noise


2017-10-14 07:59:59,681 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150820-T000000.raw'
2017-10-14 08:00:23,878 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150821-T000000.raw'
2017-10-14 08:00:47,481 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150822-T000000.raw'
2017-10-14 08:01:11,608 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150823-T000000.raw'
2017-10-14 08:01:35,506 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150824-T000000.raw'
2017-10-14 08:01:59,508 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150825-T000000.raw'
2017-10-14 08:02:23,467 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150826-T000000.raw'
2017-10-14 08:02:47,091 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150827-T000000.raw'
2017-10-14 08:03:10,991 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150828-T000000.raw'
2017-10-14 08:03:34,634 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150829-T000000.raw'
2017-10-14 08:03:58,013 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150830-T000000.raw'
2017-10-14 08:04:21,578 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150831-T000000.raw'
2017-10-14 08:04:50,611 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150901-T000000.raw'
2017-10-14 08:05:15,260 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150902-T000000.raw'
2017-10-14 08:05:39,400 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150903-T000000.raw'
2017-10-14 08:06:03,302 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150904-T000000.raw'
2017-10-14 08:06:27,248 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150905-T000000.raw'
2017-10-14 08:06:51,176 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150906-T000000.raw'
2017-10-14 08:07:15,087 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150907-T000000.raw'
2017-10-14 08:07:39,305 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150908-T000000.raw'
2017-10-14 08:08:03,214 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150909-T000000.raw'
2017-10-14 08:08:27,584 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150910-T000000.raw'
2017-10-14 08:08:51,996 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150911-T000000.raw'
2017-10-14 08:09:16,534 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150912-T000000.raw'
2017-10-14 08:09:40,562 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150913-T000000.raw'
2017-10-14 08:10:04,609 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150914-T000000.raw'
2017-10-14 08:10:28,550 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150915-T000000.raw'
2017-10-14 08:10:52,570 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150916-T000000.raw'
2017-10-14 08:11:16,783 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150917-T000000.raw'
2017-10-14 08:11:40,633 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150918-T000000.raw'
2017-10-14 08:12:04,319 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150919-T000000.raw'
2017-10-14 08:12:28,035 INFO     mi.instrument.kut.ek60.ooicore.zplsc_b Begin processing echogram data: '/media/wu-jung/wjlee_apl_2/ooi_zplsc_600m/OOI-D20150920-T000000.raw'

In [103]:
MVBS = get_MVBS(Sv_corr_all,bin_size,ping_bin_range,depth_bin_range=2.5)
MVBS_5 = get_MVBS(Sv_corr_all,bin_size,ping_bin_range,depth_bin_range=5)
fig,ax = plt.subplots(3,1,figsize=(15,9))
ax[0].imshow(Sv_corr_all[0,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)
ax[1].imshow(MVBS[0,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)
ax[2].imshow(MVBS_5[0,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)


Out[103]:
<matplotlib.image.AxesImage at 0x7f53d56b7450>

In [107]:
fig,ax = plt.subplots(3,1,figsize=(15,9))
ax[0].imshow(Sv_corr_all[1,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)
ax[1].imshow(MVBS[1,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)
ax[2].imshow(MVBS_5[1,:,:],aspect='auto',cmap=e_cmap,norm=e_norm)


Out[107]:
<matplotlib.image.AxesImage at 0x7f53d4c51c50>

In [ ]:


In [ ]:


In [ ]:


In [ ]:


In [ ]:

Freq-differencing


In [104]:
Sv_1 = MVBS[2,:,:]
Sv_2 = MVBS[0,:,:]

yes_1 = ~np.isnan(Sv_1)
yes_2 = ~np.isnan(Sv_2)
Sv_diff_12 = Sv_1 - Sv_2
Sv_diff_12[yes_1 & ~yes_2] = np.inf
Sv_diff_12[~yes_1 & yes_2] = -np.inf

idx_fish = (np.isneginf(Sv_diff_12) | (Sv_diff_12<-2)) & (Sv_diff_12>-16)
idx_fish_2 = (np.isneginf(Sv_diff_12) | (Sv_diff_12<-2))
idx_zoop = np.isposinf(Sv_diff_12) | (Sv_diff_12>-2)
idx_other = Sv_diff_12<-16

In [105]:
fig,(ax0,ax1,ax2) = plt.subplots(ncols=3,figsize=(22,5))
im0 = ax0.imshow(np.ma.masked_where(~idx_fish,MVBS[0,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im0,ax=ax0)
ax0.set_title('Fish region, 38 kHz')
im1 = ax1.imshow(np.ma.masked_where(~idx_fish,MVBS[1,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im1,ax=ax1)
ax1.set_title('Fish region, 120 kHz')
im2 = ax2.imshow(np.ma.masked_where(~idx_fish,MVBS[2,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im2,ax=ax2)
ax2.set_title('Fish region, 200 kHz')


Out[105]:
<matplotlib.text.Text at 0x7f53d55185d0>

In [106]:
fig,(ax0,ax1,ax2) = plt.subplots(ncols=3,figsize=(22,5))
im0 = ax0.imshow(np.ma.masked_where(~idx_zoop,MVBS[0,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im0,ax=ax0)
ax0.set_title('Zoop region, 38 kHz')
im1 = ax1.imshow(np.ma.masked_where(~idx_zoop,MVBS[1,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im1,ax=ax1)
ax1.set_title('Zoop region, 120 kHz')
im2 = ax2.imshow(np.ma.masked_where(~idx_zoop,MVBS[2,:,:]),aspect='auto',cmap=e_cmap,norm=e_norm)
fig.colorbar(im2,ax=ax2)
ax2.set_title('Zoop region, 200 kHz')


Out[106]:
<matplotlib.text.Text at 0x7f53d5046050>

In [ ]:


In [ ]:


In [ ]: