In [1]:

    

data_dir = 'C:/Data/Antonio/Philip/081114Patch clamp/nanorods 630/fov1/114640_take1 100Hz'


    

In [2]:

    

from __future__ import division
%matplotlib inline
import matplotlib.pyplot as plt
import seaborn as sns
from IPython.display import display
import numpy as np
import scipy.ndimage as ndi
sns.set_context(rc={'lines.markeredgewidth': 1})  # workaround for a bug in matplotlib 1.4.2


    

In [3]:

    

from IPython.core.display import HTML
HTML(open("./styles/custom2.css", "r").read())


    

    
Out[3]:

In [4]:

    

from figscroller import FrameScroller, HorizontalScroller
from patchclamp import PatchDataset
import timetraces as tt


    

In [5]:

    

data = PatchDataset(data_dir)


    

In [6]:

    

data.params['Square frame rate (Hz)']


    

    
Out[6]:





400.0

In [7]:

    

plt.plot(data.time[:100], data.voltage[:100], label='Voltage')
plt.ylabel('Voltage (V)')
plt.twinx()
plt.plot(data.time[:100], data.current[:100], label='Current',
         color=sns.color_palette()[1])
plt.ylabel('Current (pA)')
plt.grid(False)
plt.xlabel('Time (s)')


    

    
Out[7]:





<matplotlib.text.Text at 0x1938f630>






    

In [8]:

    

sns.set_style('dark')


    

In [9]:

    

fig, ax = plt.subplots(figsize=(10, 6))
im = ax.imshow(data.video.mean(0), cmap='cubehelix')
plt.colorbar(im);
ax.set_title('Mean frame')


    

    
Out[9]:





<matplotlib.text.Text at 0x18e62908>






    

In [10]:

    

mvideo = data.video.mean(0)
mvideo.shape, mvideo.dtype


    

    
Out[10]:





((192L, 272L), dtype('float64'))

In [11]:

    

mvideo_highpass = mvideo - ndi.gaussian_filter(mvideo, sigma=8)
th = mvideo_highpass > 0.4*mvideo_highpass.max()


    

In [12]:

    

fig, ax = plt.subplots(1, 2, figsize=(20, 6))
im = ax[0].imshow(mvideo_highpass, cmap='cubehelix')
#plt.colorbar(im);
im = ax[1].imshow(th, cmap='cubehelix')
ax[0].set_title('High-pass Mean frame')


    

    
Out[12]:





<matplotlib.text.Text at 0x1b6fa588>






    

In [13]:

    

video_highpass = data.video - ndi.gaussian_filter(data.video.astype('float'), sigma=(6, 6, 10))
# video_highpass.shape, video_highpass.dtype

# fig, ax = plt.subplots(figsize=(10, 6))
# im = ax.imshow(video_highpass.mean(0), cmap='cubehelix')
# plt.colorbar(im);
# ax.set_title('Mean frame after HPF')


    

In [14]:

    

video_highpass.shape


    

    
Out[14]:





(1996L, 192L, 272L)

In [15]:

    

data.video.mean(), video_highpass.mean()


    

    
Out[15]:





(115.14464272448082, -9.3453705835467859e-14)

In [16]:

    

video_highpass_tm = video_highpass - video_highpass.mean(axis=0)
video_highpass_tm.shape, video_highpass_tm.dtype


    

    
Out[16]:





((1996L, 192L, 272L), dtype('float64'))

In [17]:

    

video_highpass_on = np.zeros(data.video.shape[1:])
for row in range(data.video.shape[1]):
    for col in range(data.video.shape[2]):
        on_periods = video_highpass_tm[:, row, col] > 0
        video_highpass_on[row, col] = video_highpass_tm[on_periods, row, col].mean()
video_highpass_on.shape, video_highpass_on.dtype


    

    
Out[17]:





((192L, 272L), dtype('float64'))

In [18]:

    

hp_frame = ndi.gaussian_filter(video_highpass_on, sigma=0.3)


    

In [19]:

    

fig, ax = plt.subplots(figsize=(18, 11))
im = ax.imshow(video_highpass.mean(0), cmap='cubehelix', interpolation='none')
plt.colorbar(im);
ax.set_title('Mean frame after HPF');


    

In [20]:

    

fig, ax = plt.subplots(figsize=(18, 11))
im = ax.imshow(video_highpass_on, cmap='cubehelix', interpolation='none')
plt.colorbar(im);
ax.set_title('Mean frame after HPF');


    

In [21]:

    

# %matplotlib qt
# sns.set_style('dark')


    

In [22]:

    

# fig, ax = plt.subplots(figsize=(14, 9))
# im = ax.imshow(data.video.mean(0), cmap='cubehelix', interpolation='none')
# plt.colorbar(im)

# scroller = FrameScroller(fig, im, data.video)
# points = plt.ginput(10, timeout=0)
# points


    

In [23]:

    

# fig, ax = plt.subplots(figsize=(14, 9))
# im = ax.imshow(frame, cmap='cubehelix', interpolation='none')
# plt.colorbar(im)
# points = plt.ginput(10, timeout=0)
# points


    

In [24]:

    

# fig, ax = plt.subplots(figsize=(14, 9))
# im = ax.imshow(frame, cmap='cubehelix')
# plt.colorbar(im)
# #points = plt.ginput(10, timeout=0)
# #points


    

In [25]:

    

frame = data.video.mean(0)


    

In [26]:

    

%matplotlib inline
sns.set_style('dark')


    

In [27]:

    

points_patched = \
    [(129.84400198396881, 60.998118646576899),
     (154.63469015102521, 86.020150644805398),
     (144.95690265150074, 73.968566211435302),
     (144.04390383079087, 71.959968805873615),
     (88.898775059915508, 109.94071974740365),
     (117.93213755848896, 126.92249781260698),
     (113.91494274736559, 41.100608665880458),
]    
nrows = 2
npoints = len(points_patched)
even_npoints = npoints if npoints % 2 == 0 else npoints + 1
ncols = even_npoints // 2

fig, ax = plt.subplots(nrows, ncols, figsize=(3.5*ncols, 3.5*nrows))
axr = ax.ravel()
for ip, point in enumerate(points_patched):
    rpoint = np.round(point)
    roi = tt.get_roi_square(rpoint, pad=4)
    axr[ip].imshow(hp_frame[roi], cmap='cubehelix', interpolation='none',
                   extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    axr[ip].plot(point[0], point[1], 'r+')
    axr[ip].annotate('%d' % ip, point, color='red', fontsize=14)


    

In [28]:

    

%matplotlib inline
sns.set_style('dark')


    

In [29]:

    

points_unpatched = \
    [(226.12893521850839, 166.78365250387941),
     (240.92794772704292, 168.97609880144006),
     (226.67704679289855, 159.84090589493726),
     (114.13147018478396, 26.101681743736179),
     (116.8, 11.0),
     (107.00601971771177, 27.92872032503675),
     (234, 181.6),
    ]
    
nrows = 2
npoints = len(points_unpatched)
even_npoints = npoints if npoints % 2 == 0 else npoints + 1
ncols = even_npoints // 2

fig, ax = plt.subplots(nrows, ncols, figsize=(3.5*ncols, 3.5*nrows))
axr = ax.ravel()
for ip, point in enumerate(points_unpatched):
    rpoint = np.round(point)
    roi = tt.get_roi_square(rpoint, pad=4)
    axr[ip].imshow(hp_frame[roi], cmap='cubehelix', interpolation='none',
                   extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    axr[ip].plot(point[0], point[1], 'r+')
    axr[ip].annotate('%d' % ip, point, color='red', fontsize=14)


    

In [30]:

    

fig, ax = plt.subplots(figsize=(18, 11))
for ip, point in enumerate(points_patched):
    plt.plot(point[0], point[1], 'r+')
    plt.annotate('%d' % ip, point, color='white', fontsize=12)
im = ax.imshow(data.video.mean(0), cmap='cubehelix', interpolation='none')
plt.colorbar(im);
ax.set_title('QD on patched cell')


    

    
Out[30]:





<matplotlib.text.Text at 0x5e4784e0>






    

In [31]:

    

fig, ax = plt.subplots(figsize=(18, 11))
for ip, point in enumerate(points_unpatched):
    plt.plot(point[0], point[1], 'r+')
    plt.annotate('%d' % ip, point, color='white', fontsize=12)
im = ax.imshow(data.video.mean(0), cmap='cubehelix')
plt.colorbar(im);
ax.set_title('Unpatched QD');


    

In [32]:

    

sns.set_style('darkgrid')


    

In [33]:

    

shape2d = (16, 16)
radii = [1, 1.6, 2, 2.2, 2.5]
point = (6, 8)

fig, ax = plt.subplots(1, len(radii), figsize=(3*len(radii), 3))
for i, clip_radius in enumerate(radii):
    mask = tt.get_roi_circle(point, clip_radius, shape2d)
    a = np.zeros(shape2d)
    a[mask] = 1
    im = ax[i].imshow(a, interpolation='none', cmap='gray', alpha=0.3)
    ax[i].axvline(point[0])
    ax[i].axhline(point[1])
    ax[i].set_title('Radius = %.2f, # Pixels = %d' % (clip_radius, len(mask[0])))
    
point = (point[0]+0.5, point[1]+0.5)
fig, ax = plt.subplots(1, len(radii), figsize=(3*len(radii), 3))
for i, clip_radius in enumerate(radii):
    mask = tt.get_roi_circle(point, clip_radius, shape2d)
    a = np.zeros(shape2d)
    a[mask] = 1
    im = ax[i].imshow(a, interpolation='none', cmap='gray', alpha=0.3)
    ax[i].axvline(point[0])
    ax[i].axhline(point[1])
    ax[i].set_title('Radius = %.2f, # Pixels = %d' % (clip_radius, len(mask[0])))

point = (point[0]-0.4, point[1]+0.4)
fig, ax = plt.subplots(1, len(radii), figsize=(3*len(radii), 3))
for i, clip_radius in enumerate(radii):
    mask = tt.get_roi_circle(point, clip_radius, shape2d)
    a = np.zeros(shape2d)
    a[mask] = 1
    im = ax[i].imshow(a, interpolation='none', cmap='gray', alpha=0.3)
    ax[i].axvline(point[0])
    ax[i].axhline(point[1])
    ax[i].set_title('Radius = %.2f, # Pixels = %d' % (clip_radius, len(mask[0])))


    

In [34]:

    

sns.set_style('dark')


    

In [35]:

    

shape2d = (11, 11)
clip_radius = 1.8

nrows = 2
npoints = len(points_patched)
even_npoints = npoints if npoints % 2 == 0 else npoints + 1
ncols = even_npoints // 2

fig, ax = plt.subplots(nrows, ncols, figsize=(3.5*ncols, 3.5*nrows))
axr = ax.ravel()
for ip, point in enumerate(points_patched):
    rpoint = np.round(point)
    roi = tt.get_roi_square(rpoint, pad=(shape2d[0]-1)//2)
    axr[ip].imshow(hp_frame[roi], cmap='cubehelix', interpolation='none',
                   extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    axr[ip].plot(point[0], point[1], 'r+')
    axr[ip].annotate('%d' % ip, point, color='red', fontsize=14)
    
    mask = tt.get_roi_circle(point, clip_radius, frame.shape)
    a = np.ones(frame.shape)
    a[mask] = 0
    im = axr[ip].imshow(a[roi], interpolation='none', cmap='gray', vmin=0.5, alpha=0.7,
                        extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    im.cmap.set_under(alpha=0)
    axr[ip].set_title('Radius = %.2f, # Pixels = %d' % (clip_radius, len(mask[0])))


    

In [36]:

    

shape2d = (11, 11)
clip_radius = 1.8

nrows = 2
npoints = len(points_unpatched)
even_npoints = npoints if npoints % 2 == 0 else npoints + 1
ncols = even_npoints // 2

fig, ax = plt.subplots(nrows, ncols, figsize=(3.5*ncols, 3.5*nrows))
axr = ax.ravel()
for ip, point in enumerate(points_unpatched):
    rpoint = np.round(point)
    roi = tt.get_roi_square(rpoint, pad=(shape2d[0]-1)//2)
    axr[ip].imshow(hp_frame[roi], cmap='cubehelix', interpolation='none',
                   extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    axr[ip].plot(point[0], point[1], 'r+')
    axr[ip].annotate('%d' % ip, point, color='red', fontsize=14)
    
    mask = tt.get_roi_circle(point, clip_radius, frame.shape)
    a = np.ones(frame.shape)
    a[mask] = 0
    im = axr[ip].imshow(a[roi], interpolation='none', cmap='gray', vmin=0.5, alpha=0.7,
                        extent=(roi[1].start-0.5, roi[1].stop-0.5, roi[0].stop-0.5, roi[0].start-0.5))
    im.cmap.set_under(alpha=0)
    axr[ip].set_title('Radius = %.2f, # Pixels = %d' % (clip_radius, len(mask[0])))


    

In [37]:

    

%matplotlib inline


    

In [38]:

    

sns.set_style('darkgrid')


    

In [39]:

    

clip_radius = 1.8
vert_offset = 25


    

In [40]:

    

fig, ax = plt.subplots(1, 1, figsize=(20, 1.5*len(points_patched)), sharex=True)
fig.subplots_adjust(left=0.04, right=0.96, hspace=0)
for i, point in enumerate(points_patched):
    timetrace = tt.get_timetrace(data.video, point, clip_radius=clip_radius, detrend_sigma=250)
    ftimetrace = ndi.filters.gaussian_filter1d(timetrace, 10)
    vpos = vert_offset*i
    
    l, = ax.plot(data.time, timetrace + vpos, 
                 label='QD %d (%.1f, %.1f)' % (i, point[0], point[1]))
    ax.plot(data.time, ftimetrace + vpos, color='k')
#ax.set_xlim(0, 2)
ax.set_ylim(-20)
plt.legend(ncol=5)


    

    
Out[40]:





<matplotlib.legend.Legend at 0x5eca9f98>






    

In [41]:

    

colors = sns.color_palette('deep', max(len(points_patched), len(points_unpatched)))


    

In [42]:

    

clip_radius = 1.8
detrend_sigma = 300
lowpass_sigma = 10
thresholds = [-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1]

points_sel = points_patched
fig, ax = plt.subplots(len(points_sel), 1, figsize=(20, 2*len(points_sel)), sharex=True)
fig.tight_layout()
fig.suptitle('Patched')
fig.subplots_adjust(left=0.04, right=0.96, hspace=0.05)
for i, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    ftimetrace = ndi.filters.gaussian_filter1d(timetrace, 10)
    on_periods = tt.get_on_periods_slices(timetrace, thresholds[i], lowpass_sigma=lowpass_sigma, align=4)
    
    ax[i].plot(data.time, timetrace, color=colors[i], alpha=0.3,
               label='QD %d (%.1f, %.1f)' % (i, point[0], point[1]))
    for on_period in on_periods:
        ax[i].plot(data.time[on_period], timetrace[on_period], color=colors[i], alpha=1)
    #ax[i].plot(data.time[time_on], trace_on, '.', color=colors[i], alpha=1)
    ax[i].plot(data.time, ftimetrace, color='k')
    ax[i].axhline(thresholds[i], color='k', ls='--', alpha=0.7)
    ax[i].locator_params(axis='y', tight=True, nbins=2)
    ax[i].legend()


    

In [43]:

    

points_sel = points_unpatched

fig, ax = plt.subplots(len(points_sel), 1, figsize=(20, 2*len(points_sel)), sharex=True)
fig.tight_layout()
fig.suptitle('Unpatched')
fig.subplots_adjust(left=0.04, right=0.96, hspace=0.05)
for i, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    ftimetrace = ndi.filters.gaussian_filter1d(timetrace, 10)
    on_periods = tt.get_on_periods_slices(timetrace, thresholds[i], lowpass_sigma=lowpass_sigma, align=4)
    
    ax[i].plot(data.time, timetrace, color=colors[i], alpha=0.3,
               label='QD %d (%.1f, %.1f)' % (i, point[0], point[1]))
    for on_period in on_periods:
        ax[i].plot(data.time[on_period], timetrace[on_period], color=colors[i], alpha=1)
    #ax[i].plot(data.time[time_on], trace_on, '.', color=colors[i], alpha=1)
    ax[i].plot(data.time, ftimetrace, color='k')
    ax[i].axhline(thresholds[i], color='k', ls='--', alpha=0.7)
    ax[i].locator_params(axis='y', tight=True, nbins=2)
    ax[i].legend()


    

In [44]:

    

%matplotlib inline


    

In [45]:

    

sns.set_style('darkgrid')


    

In [46]:

    

point = points_patched[1]

fig, ax = plt.subplots(2, 1, figsize=(22, 8), sharex=True)
timetrace1 = tt.get_timetrace(data.video, point, clip_radius=clip_radius, detrend_sigma=None)
trend = ndi.filters.gaussian_filter1d(timetrace1, detrend_sigma)
ftimetrace = ndi.filters.gaussian_filter1d(timetrace1 -  trend, lowpass_sigma)
ax[0].plot(data.time, timetrace1, label='QD (%.1f, %.1f)' % point);
ax[0].plot(data.time, trend, color='k');
ax[1].plot(data.time, timetrace1 - trend, label='QD (%.1f, %.1f)' % point);
ax[1].plot(data.time, ftimetrace, color='k');
ax[0].set_title('Difference between background corrected and raw timetrace');


    

In [47]:

    

f0 = 100
sim_signal = 1 + np.cos(2*np.pi*f0*data.time - np.pi/4)


    

In [48]:

    

points_sel = points_patched
clip_radius = 1.8
detrend_sigma = 300

mdouble_diff, mdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
sdouble_diff, sdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace += 0.5*sim_signal
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    double_diff = tt.edge_diff_sum(timetrace_on, offset=0)
    double_diffo = tt.edge_diff_sum(timetrace_on, offset=1)
    
    mdouble_diff[ip], mdouble_diffo[ip] = double_diff.mean(), double_diffo.mean()
    sdouble_diff[ip], sdouble_diffo[ip] = double_diff.std(), double_diffo.std()

fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(mdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[0].plot(mdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[0].set_title('Mean of (sum of falling and rising edges)')
ax[1].plot(sdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[1].plot(sdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[1].set_title('Std.Dev. of (sum of falling and rising edges)')
ax[2].plot(mdouble_diff/sdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[2].plot(mdouble_diffo/sdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[2].set_title('SNR of (sum of falling and rising edges)')
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Patched');
plt.close(fig)
fig_patched = fig


    

In [49]:

    

points_sel = points_unpatched

mdouble_diff, mdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
sdouble_diff, sdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    double_diff = tt.edge_diff_sum(timetrace_on, offset=0)
    double_diffo = tt.edge_diff_sum(timetrace_on, offset=1)
    
    mdouble_diff[ip], mdouble_diffo[ip] = double_diff.mean(), double_diffo.mean()
    sdouble_diff[ip], sdouble_diffo[ip] = double_diff.std(), double_diffo.std()

fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(mdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[0].plot(mdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[0].set_title('Mean of (sum of falling and rising edges)')
ax[1].plot(sdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[1].plot(sdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[1].set_title('Std.Dev. of (sum of falling and rising edges)')
ax[2].plot(mdouble_diff/sdouble_diff, '-o', mew=0, ms=6, label='2-edges sum')
ax[2].plot(mdouble_diffo/sdouble_diffo, '-o', mew=0, ms=6, label='2-edges sum out-of-phase')
ax[2].set_title('SNR of (sum of falling and rising edges)')
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Unpatched');
plt.close(fig)
fig_unpatched = fig


    

In [50]:

    

display(fig_patched)
display(fig_unpatched)


    

In [51]:

    

%matplotlib inline


    

In [52]:

    

sns.set_style('darkgrid')


    

In [53]:

    

points_sel = points_patched
ncycles = 9
running_avg = True
clip_radius = 1.8
detrend_sigma = 300

n_mdouble_diff, n_mdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
n_sdouble_diff, n_sdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace += 0.5*sim_signal
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    double_diff = tt.edge_diff_sum(timetrace_on)
    double_diffo = tt.edge_diff_sum(timetrace_on, offset=1)
    n_double_diff = tt.avg_nsamples(double_diff, num_samples=ncycles, running_avg=running_avg)
    n_double_diffo = tt.avg_nsamples(double_diffo, num_samples=ncycles, running_avg=running_avg)
    
    n_mdouble_diff[ip], n_mdouble_diffo[ip] = n_double_diff.mean(), n_double_diffo.mean()
    n_sdouble_diff[ip], n_sdouble_diffo[ip] = n_double_diff.std(),  n_double_diffo.std()

name = '%d-cycles avg of 2-edges sum' % ncycles 
fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(n_mdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[0].plot(n_mdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[0].set_title('Mean of (%s)' % name)
ax[1].plot(n_sdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[1].plot(n_sdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[1].set_title('Std.Dev. of (%s)' % name)
ax[2].plot(n_mdouble_diff/n_sdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[2].plot(n_mdouble_diffo/n_sdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[2].set_title('SNR of (%s)' % name)
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Patched, %d-cycles average' % ncycles);
plt.close(fig)
fig_patched = fig


    

In [54]:

    

points_sel = points_unpatched

n_mdouble_diff, n_mdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
n_sdouble_diff, n_sdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    double_diff = tt.edge_diff_sum(timetrace_on)
    double_diffo = tt.edge_diff_sum(timetrace_on, offset=1)
    n_double_diff = tt.avg_nsamples(double_diff, num_samples=ncycles, running_avg=running_avg)
    n_double_diffo = tt.avg_nsamples(double_diffo, num_samples=ncycles, running_avg=running_avg)
    
    n_mdouble_diff[ip], n_mdouble_diffo[ip] = n_double_diff.mean(), n_double_diffo.mean()
    n_sdouble_diff[ip], n_sdouble_diffo[ip] = n_double_diff.std(),  n_double_diffo.std()

name = '%d-cycles avg of 2-edges sum' % ncycles 
fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(n_mdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[0].plot(n_mdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[0].set_title('Mean of (%s)' % name)
ax[1].plot(n_sdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[1].plot(n_sdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[1].set_title('Std.Dev. of (%s)' % name)
ax[2].plot(n_mdouble_diff/n_sdouble_diff, '-o', mew=0, ms=6, label='in-phase')
ax[2].plot(n_mdouble_diffo/n_sdouble_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[2].set_title('SNR of (%s)' % name)
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Unpatched, %d-cycles average' % ncycles);
plt.close(fig)
fig_unpatched = fig


    

In [55]:

    

display(fig_patched)
display(fig_unpatched)


    

In [56]:

    

points_sel = points_patched
ncycles = 18
running_avg = True
clip_radius = 1.8
detrend_sigma = 300

n_mean_diff, n_mean_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
n_std_diff, n_std_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    all_diff = tt.edge_diff_all(timetrace_on)
    all_diffo = tt.edge_diff_all(timetrace_on, offset=1)
    n_all_diff = tt.avg_nsamples(all_diff, num_samples=ncycles, running_avg=running_avg)
    n_all_diffo = tt.avg_nsamples(all_diffo, num_samples=ncycles, running_avg=running_avg)

    n_mean_diff[ip], n_mean_diffo[ip] = n_all_diff.mean(), n_all_diffo.mean()
    n_std_diff[ip], n_std_diffo[ip] = n_all_diff.std(), n_all_diffo.std()

name = '%d-cycles avg of alternated diff' % ncycles 
fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(n_mean_diff, '-o', mew=0, ms=6, label='in-phase')
ax[0].plot(n_mean_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[0].set_title('Mean of (%s)' % name)
ax[1].plot(n_std_diff, '-o', mew=0, ms=6, label='in-phase')
ax[1].plot(n_std_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[1].set_title('Std.Dev. of (%s)' % name)
ax[2].plot(n_mean_diff/n_std_diff, '-o', mew=0, ms=6, label='in-phase')
ax[2].plot(n_mean_diffo/n_std_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[2].set_title('SNR of (%s)' % name)
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Patched');
plt.close(fig)
fig_patched = fig;


    

In [57]:

    

points_sel = points_unpatched

n_mean_diff, n_mean_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
n_std_diff, n_std_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace_on = timetrace
    time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    all_diff = tt.edge_diff_all(timetrace_on)
    all_diffo = tt.edge_diff_all(timetrace_on, offset=1)
    n_all_diff = tt.avg_nsamples(all_diff, num_samples=ncycles, running_avg=running_avg)
    n_all_diffo = tt.avg_nsamples(all_diffo, num_samples=ncycles, running_avg=running_avg)

    n_mean_diff[ip], n_mean_diffo[ip] = n_all_diff.mean(), n_all_diffo.mean()
    n_std_diff[ip], n_std_diffo[ip] = n_all_diff.std(), n_all_diffo.std()


fig, ax = plt.subplots(1, 3, figsize=(18, 4))
ax[0].plot(n_mean_diff, '-o', mew=0, ms=6, label='in-phase')
ax[0].plot(n_mean_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[0].set_title('Mean of (%s)' % name)
ax[1].plot(n_std_diff, '-o', mew=0, ms=6, label='in-phase')
ax[1].plot(n_std_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[1].set_title('Std.Dev. of (%s)' % name)
ax[2].plot(n_mean_diff/n_std_diff, '-o', mew=0, ms=6, label='in-phase')
ax[2].plot(n_mean_diffo/n_std_diffo, '-o', mew=0, ms=6, label='out-of-phase')
ax[2].set_title('SNR of (%s)' % name)
ax[1].set_xlabel('Quantum Dot')
ax[2].legend(ncol=2);
fig.suptitle('Unpatched');
plt.close(fig)
fig_unpatched = fig;


    

In [58]:

    

display(fig_patched)
display(fig_unpatched)


    

In [59]:

    

def t_run(sample_diff, num_samples, running_avg):
    if running_avg:
        return np.arange(0, sample_diff.size)*2 + 0.5 + num_samples
    else:
        return np.arange(0, sample_diff.size)*2*num_samples + num_samples


    

In [60]:

    

points_sel = points_patched
ncycles = 16
running_avg = True
clip_radius = 1.8
detrend_sigma = 300

n_alt_diff, n_alt_diffo = [], []
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    #timetrace += 0.5*sim_signal
    timetrace_on = timetrace
    #time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    n_alt_diff.append(
        tt.avg_nsamples(tt.edge_diff_all(timetrace, offset=0), 
                        num_samples=ncycles, running_avg=running_avg))
    n_alt_diffo.append(
        tt.avg_nsamples(tt.edge_diff_all(timetrace, offset=1), 
                        num_samples=ncycles, running_avg=running_avg))
        
alpha = 0.8
name = '%d-cycles avg of 2-edges sum' % ncycles
fig, ax = plt.subplots(len(points_sel), 1, 
                       figsize=(18, 4*len(points_sel)),
                       sharex=False)
fig.tight_layout()
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    signal_clue = n_alt_diff[ip]**2
    signal_clueo = n_alt_diffo[ip]**2
    ax[ip].plot(t_run(signal_clue, ncycles, running_avg), 
                signal_clue, label='QD %d: signal^2 (%d cycles) in-phase' % (ip, ncycles), alpha=alpha)
    ax[ip].plot(t_run(signal_clueo, ncycles, running_avg), 
                signal_clueo, label='QD %d: signal^2 (%d cycles) out-of-phase' % (ip, ncycles), alpha=alpha)

    ax2 = ax[ip].twinx()
    ax2.grid(False)
    ax2.plot(timetrace, 'k', alpha=0.3)
    ax2.set_ylim(-40)

    ax[0].set_title('Patched: %d-cycles average' % ncycles)
    ax[ip].legend(loc='upper left')
    #ax[ip].set_ylim(-0, 20)
plt.close(fig)
fig_patched = fig;


    

In [61]:

    

points_sel = points_unpatched

n_alt_diff, n_alt_diffo = [], []
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    timetrace_on = timetrace
    #time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

    n_alt_diff.append(
        tt.avg_nsamples(tt.edge_diff_all(timetrace, offset=0), 
                        num_samples=ncycles, running_avg=running_avg))
    n_alt_diffo.append(
        tt.avg_nsamples(tt.edge_diff_all(timetrace, offset=1), 
                        num_samples=ncycles, running_avg=running_avg))
        

fig, ax = plt.subplots(len(points_sel), 1, 
                       figsize=(18, 4*len(points_sel)),
                       sharex=False)
fig.tight_layout()
for ip, point in enumerate(points_sel):
    timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
    signal_clue = n_alt_diff[ip]**2 
    signal_clueo = n_alt_diffo[ip]**2
    ax[ip].plot(t_run(signal_clue, ncycles, running_avg), 
                signal_clue, label='QD %d: signal^2 (%d cycles) in-phase' % (ip, ncycles), alpha=alpha)
    ax[ip].plot(t_run(signal_clueo, ncycles, running_avg), 
                signal_clueo, label='QD %d: signal^2 (%d cycles) out-of-phase' % (ip, ncycles), alpha=alpha)
    
    ax2 = ax[ip].twinx()
    ax2.grid(False)
    ax2.plot(timetrace, 'k', alpha=0.3)
    ax2.set_ylim(-40)

    ax[0].set_title('Unpatched: %d-cycles average' % ncycles)
    ax[ip].legend(loc='upper left')
    #ax[ip].set_ylim(-0, 20)
plt.close(fig)
fig_unpatched = fig;


    

In [62]:

    

fig_patched


    

    
Out[62]:

In [63]:

    

fig_unpatched


    

    
Out[63]:

In [64]:

    

points_sel = points_patched
multi_ncycles = [1, 4, 9, 14]
running_avg = True
clip_radius = 1.8
detrend_sigma = 300

mdouble_diff, mdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
sdouble_diff, sdouble_diffo = np.zeros(len(points_sel)), np.zeros(len(points_sel))
n_double_diff, n_double_diffo = {n: [] for n in multi_ncycles}, {n: [] for n in multi_ncycles}
for ncycles in multi_ncycles:
    for ip, point in enumerate(points_sel):
        timetrace = tt.get_timetrace(data.video, point=point, clip_radius=clip_radius, detrend_sigma=detrend_sigma)
        #timetrace += 0.5*sim_signal
        #timetrace_on = timetrace
        time_on, timetrace_on = tt.get_on_periods_timetrace(timetrace, thresholds[ip], lowpass_sigma=10, align=4)

        double_diff = tt.edge_diff_sum(timetrace_on)
        double_diffo = tt.edge_diff_sum(timetrace_on, offset=1)
        n_double_diff[ncycles].append(tt.avg_nsamples(double_diff, num_samples=ncycles, running_avg=running_avg))
        n_double_diffo[ncycles].append(tt.avg_nsamples(double_diffo, num_samples=ncycles, running_avg=running_avg))

alpha = 0.5
nbins = 25
name = '%d-cycles avg of 2-edges sum' % ncycles
fig, ax = plt.subplots(len(points_sel), len(multi_ncycles), 
                       figsize=(4*len(multi_ncycles), 3*len(points_sel)),
                       sharex=False)
fig.tight_layout()
for ip, point in enumerate(points_sel):
    for cidx, ncycles in enumerate(multi_ncycles):
        sample, sampleo = n_double_diff[ncycles][ip], n_double_diffo[ncycles][ip]
        
        bin_abs_max = np.max([np.abs(sample).max(), np.abs(sampleo).max()])
        bins = np.linspace(-bin_abs_max*1.1, bin_abs_max*1.1, nbins+1)
        centers = bins[:-1] + 0.5*(bins[1] - bins[0])

        ax[ip, cidx].hist(n_double_diff[ncycles][ip], bins=bins, label='QD %d in-phase' % ip, alpha=alpha)
        ax[ip, cidx].hist(n_double_diffo[ncycles][ip], bins=bins, label='QD %d out-of-phase' % ip, alpha=alpha)
        ax[ip, cidx].axvline(n_double_diff[ncycles][ip].mean(), color=colors[0])
        ax[ip, cidx].axvline(n_double_diffo[ncycles][ip].mean(), color=colors[1])
        
        ax[0, cidx].set_title('Patched: %d-cycles average' % ncycles)
    ax[ip, 0].legend(loc='upper left')
plt.close(fig)
fig_patched = fig;