In [1]:

    

%matplotlib inline


    

    




/usr/local/python-2.7.6/lib/python2.7/site-packages/matplotlib-1.5.1-py2.7-linux-x86_64.egg/matplotlib/font_manager.py:273: UserWarning: Matplotlib is building the font cache using fc-list. This may take a moment.
  warnings.warn('Matplotlib is building the font cache using fc-list. This may take a moment.')

In [2]:

    

from numpy import array
import matplotlib.pyplot as plt
import seaborn as sns


    

    




/usr/local/python-2.7.6/lib/python2.7/site-packages/matplotlib-1.5.1-py2.7-linux-x86_64.egg/matplotlib/__init__.py:872: UserWarning: axes.color_cycle is deprecated and replaced with axes.prop_cycle; please use the latter.
  warnings.warn(self.msg_depr % (key, alt_key))

In [3]:

    

from showit import image, tile
import matplotlib.animation as animation


    

In [4]:

    

import thunder


    

In [5]:

    

import json


    

In [6]:

    

from os.path import join, exists
from os import mkdir, makedirs


    

In [7]:

    

directory = '/tier2/freeman/Nick/mVR/sessions'


    

In [8]:

    

key = '000133'
prefix = 'trial'
path = directory + '/' + key
print path
print exists(path)


    

    




/tier2/freeman/Nick/mVR/sessions/000133
True

In [9]:

    

savepath = path + '/movies'
if not exists(savepath):
    makedirs(savepath)


    

In [52]:

    

Writer = animation.writers['ffmpeg']
writer = Writer(fps=15, metadata=dict(artist='Me'), bitrate=40000)


    

In [11]:

    

data = thunder.images.frombinary(path + '/images', engine=sc)


    

In [12]:

    

summarypath = path + '/summary'
with open(summarypath + '/meta.json') as infile:
    meta = json.load(infile)


    

In [13]:

    

meta


    

    
Out[13]:





{u'averaging': u'_Inf_',
 u'depths': 280,
 u'nchannels': 1,
 u'nplanes': 1,
 u'nrois': 1,
 u'power': 22,
 u'rois': [{u'center': [-0.006609804447150028, 4.996003611375001e-16],
   u'depths': 0.0,
   u'npixels': [1792.0, 1766.0],
   u'size': [4.262350232069999, 4.200000000525001]}],
 u'shape': [1183, 1766, 1792],
 u'volumeRate': 1.817867823}

In [14]:

    

pixelsPerUmOrig = meta['rois'][0]['npixels'][0]/meta['rois'][0]['size'][0]/1000
aspect = meta['rois'][0]['size'][0]/meta['rois'][0]['size'][1]
volumeRate = meta['volumeRate']
scaleBar = 100*round(meta['rois'][0]['size'][0]*1.8)


    

In [53]:

    

fSpace = 1
dsTime = 1
dsSpace = 1
frameRate = volumeRate/dsTime
pixelsPerUm = pixelsPerUmOrig/dsSpace


    

In [54]:

    

from skimage.restoration import denoise_bilateral
from skimage import img_as_uint


    

In [55]:

    

def bilateral_filter_2D(img):
    return map(lambda x: img_as_uint(denoise_bilateral(x, win_size=fSpace, sigma_range=0.1, sigma_spatial=50)), img)


    

In [56]:

    

if fSpace > 1:
    #filtered = data.median_filter((1, fSpace, fSpace))
    filtered = data.clip(min=0).map(bilateral_filter_2D)
else:
    filtered = data


    

In [57]:

    

if dsSpace > 1:
    filtered = filtered.subsample((1, dsSpace, dsSpace))


    

In [58]:

    

from numpy import convolve, ones


    

In [59]:

    

def smooth(s, dsTime):
    r = convolve(s, ones(dsTime)/dsTime, 'same').astype('int16')
    return r[::dsTime]


    

In [60]:

    

if dsTime > 1:
    smoothed = filtered.map_as_series(lambda x: smooth(x, dsTime))
else:
    smoothed = filtered


    

In [61]:

    

result = smoothed.toarray()


    

In [319]:

    

from numpy import percentile


    

In [320]:

    

def dff(F):
    #F0 = F.mean()
    F0 = percentile(F,20)
    offset = 20
    return (F - F0)/(F0+offset)


    

In [321]:

    

normed = smoothed.map_as_series(dff)
#normedF = normed.clip(min=0).map(lambda x: denoise_bilateral(x, win_size=fSpace, sigma_range=0.1, sigma_spatial=20))


    

In [322]:

    

result = normed.toarray()


    

In [62]:

    

plane = 0
if len(result.shape) == 3:
    plane = 0
    movie = result
else:
    movie = result[:,plane,:,:]
means = movie.mean()
label = 'movie'
clim = means #3.0/4 #means


    

In [63]:

    

frame = movie.shape[0]/2
img = movie[frame,:,:]

fig = plt.figure(figsize=[12, 12.0*img.shape[0]/img.shape[1]])
fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=None, hspace=None)
ax = plt.gca()
im = image(img, clim=(0, clim*8), ax=ax)

time = ax.text(.97*img.shape[1], .04*img.shape[0], '%.1f s' % (frame/frameRate), color='white', fontsize=22, ha='right', fontdict={'family': 'monospace'});
ax.plot([.04*img.shape[1], .04*img.shape[1]+scaleBar*pixelsPerUm], [.94*img.shape[0], .94*img.shape[0]], 'w');
sclae = ax.text(.04*img.shape[1]+scaleBar*pixelsPerUm/2, .97*img.shape[0], '%d um' % scaleBar, color='white', fontsize=22, ha='center', fontdict={'family': 'monospace'});
plt.xlim([0, img.shape[1]]);
plt.ylim([img.shape[0], 0]);


    

In [64]:

    

nframes = movie.shape[0]-dsTime

def update(frame):
    im.set_array(movie[dsTime/2+frame,:,:])
    time.set_text('%.1f s' % ((dsTime/2+frame)/frameRate))

ani = animation.FuncAnimation(fig, update, nframes, blit=False)
ani.save(join(savepath, '%s2_ds%d_fs%d_dt%d-%d.mp4' % (label, dsSpace, fSpace, dsTime, plane)), writer=writer)


    

In [32]:

    

movie = result
means = movie.mean(axis=(0,2,3))
label = 'tiled'
clim = means #[3.0/4, 3.0/4, 3.0/4, 3.0/4] #means


    

In [36]:

    

fig = plt.figure(figsize=[12, 12.0*movie.shape[2]/movie.shape[3]])
fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=None, hspace=None)
frame = movie.shape[0]/2

axes, g = tile(movie[frame].transpose(1,2,0), clim=[(0, x*6) for x in clim], axis=2, fig=fig)

img = movie[frame,0]
time = g[1].text(.92*img.shape[1], .07*img.shape[0], '%.1f s' % (frame/frameRate), color='white', fontsize=22, ha='right', fontdict={'family': 'monospace'});
g[2].plot([.05*img.shape[1], .05*img.shape[1]+scaleBar*pixelsPerUm], [.92*img.shape[0], .92*img.shape[0]], 'w');
g[2].text(.05*img.shape[1]+scaleBar*pixelsPerUm/2, .975*img.shape[0], '%d um' % scaleBar, color='white', fontsize=22, ha='center', fontdict={'family': 'monospace'});
for i in range(4):
    g[i].set_xlim([0, img.shape[1]]);
    g[i].set_ylim([img.shape[0], 0]);


    

In [38]:

    

nframes = movie.shape[0]-dsTime

def update(frame):
    [axes[i].set_array(movie[frame+dsTime/2, i, :, :]) for i in range(movie.shape[1])]
    time.set_text('%.1f s' % ((frame+dsTime/2)/frameRate))

ani = animation.FuncAnimation(fig, update, nframes, blit=False)
ani.save(join(savepath, '%s_ds%d_fs%d_dt%d-tiled.mp4' % (label, dsSpace, fSpace, dsTime)), writer=writer)


    

In [ ]:

    




    

In [ ]:

    




    

In [ ]:

    




    

In [ ]: