Homec4science
Diffusion Motion correction paper (master) Browse

R6289/divergence.py
05139327dccemaster

Last Change
Raw File

divergence.py
View Options

import matplotlib
from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
from matplotlib.figure import Figure
from matplotlib import pyplot as plt
import numpy as np
from skimage import io
import os
#paths = ['/scratch/aymanns/Data_R57C10s/2_170608_fly1_007/results/l100g0/']
#outs = ['div_2_170608_fly1_007/div_l100g0/']
#path = '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l1000g0/'
#path = '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l1000g10/'
#path = '/scratch/aymanns/Data_R57C10s/1_170426_fly1_013/results/l1000g0/'
#paths = ['/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l100g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l200g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l300g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l400g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l500g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l600g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l700g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l800g0/',
# '/scratch/aymanns/Data_R57C10s/4_170609_fly1_000/results/l900g0/']
#
#outs = ['div_4_170609_fly1_000/div_l100g0/',
# 'div_4_170609_fly1_000/div_l200g0/',
# 'div_4_170609_fly1_000/div_l300g0/',
# 'div_4_170609_fly1_000/div_l400g0/',
# 'div_4_170609_fly1_000/div_l500g0/',
# 'div_4_170609_fly1_000/div_l600g0/',
# 'div_4_170609_fly1_000/div_l700g0/',
# 'div_4_170609_fly1_000/div_l800g0/',
# 'div_4_170609_fly1_000/div_l900g0/']
#for path,out in zip(paths,outs):
for l in [100,200,300,400,500,600,700,800,900,1000]:
for g in [0,2,4,6,8,10,12,14,16,18,20]:
print('lambda={} gamma={}'.format(l,g))
path = '/scratch/aymanns/Data_R57C10s/2_170608_fly1_007/results/l{}g{}/'.format(l,g)
out = 'div_2_170608_fly1_007/results/div_l{}g{}/'.format(l,g)
images = io.imread(path+'warped1.tif')
n_frames = images.shape[0]
frobenius_output = np.zeros((n_frames,images.shape[1],images.shape[2]))
divergence_output = np.zeros((n_frames,images.shape[1],images.shape[2]))
vectors = np.zeros((n_frames,320,320,3))
#because of columns come first y corresponds to first dimension
#roi_y_min=307
#roi_y_max=364
#roi_x_min=361
#roi_x_max=422
#Entire image
roi_y_min=0
roi_y_max=images.shape[1]
roi_x_min=0
roi_x_max=images.shape[2]
x = np.linspace(1,roi_x_max-roi_x_min,roi_x_max-roi_x_min)
y = np.linspace(1,roi_y_max-roi_y_min,roi_y_max-roi_y_min)
for i in range(n_frames-1):
print('Frame {}'.format(i+1))
wx = np.genfromtxt(path+'wx_frame{}.dat'.format(i+1),dtype=np.float,delimiter=',')
wy = np.genfromtxt(path+'wy_frame{}.dat'.format(i+1),dtype=np.float,delimiter=',')
grad_wx = np.gradient(wx)
grad_wy = np.gradient(wy)
divergence = grad_wx[1]+grad_wy[0]
frobenius = np.sqrt(np.square(grad_wx[0])+np.square(grad_wx[1])+np.square(grad_wy[0])+np.square(grad_wy[1]))
#output[i,:,:images.shape[2]]=images[i+1]
#output[i,:,images.shape[2]:]=frobenius
frobenius_output[i+1,:,:]=frobenius
divergence_output[i+1,:,:]=divergence
#fig = Figure(figsize=(4, 4), dpi=80)
#canvas = FigureCanvas(fig)
#ax = fig.gca()
#ax.streamplot(x,y,wx[roi_y_min:roi_y_max,roi_x_min:roi_x_max],wy[roi_y_min:roi_y_max,roi_x_min:roi_x_max])
#ax.axis('off')
#canvas.draw() # draw the canvas, cache the renderer
#vectors[i,:,:,:] = np.reshape(np.fromstring(canvas.tostring_rgb(), dtype='uint8'),(320,320,3))
#max_frobenius = np.max(output[:,:,images.shape[2]:].flatten())
#print(max_frobenius)
#flat = images.flatten()
#idx = np.isnan(flat)
#print(np.max(flat[]))
#print(output[:,:,images.shape[2]:])
#output[:,:,images.shape[2]:]=output[:,:,images.shape[2]:]/max_frobenius*255#(2**31-1)
#print(output[:,:,images.shape[2]:])
#io.imsave('smoothness.tif',output.astype(np.float32))
os.makedirs(out,exist_ok=True)
#io.imsave(out+'roi.tif',images[:,roi_y_min:roi_y_max,roi_x_min:roi_x_max].astype(np.float32))
io.imsave(out+'frobenius.tif',frobenius_output[:,roi_y_min:roi_y_max,roi_x_min:roi_x_max].astype(np.float32))
#io.imsave(out+'vectors.tif',vectors.astype(np.uint8))
io.imsave(out+'div_out.tif',divergence_output.astype(np.float32))
#divergence_roi = divergence_output[:,roi_y_min:roi_y_max,roi_x_min:roi_x_max]
#divergence_seperated=np.zeros((n_frames,roi_y_max-roi_y_min,(roi_x_max-roi_x_min)*2))
#idx_plus = np.where(divergence_roi>0)
#divergence_seperated[idx_plus] = divergence_roi[idx_plus]
#divergence_roi[idx_plus] = 0
#divergence_seperated[:,:,roi_x_max-roi_x_min:] = divergence_roi*-1
#io.imsave(out+'divergence.tif',divergence_seperated.astype(np.float32))
#matplotlib.use('TkAgg')
#plt.figure()
#plt.hist(divergence_output.flatten())
#plt.savefig('histogram_complete.pdf')
#
#plt.figure()
#roi_y_min=307
#roi_y_max=364
#roi_x_min=361
#roi_x_max=422
#plt.hist(divergence_output[:,roi_y_min:roi_y_max,roi_x_min:roi_x_max].flatten())
#plt.savefig('histogram_roi.pdf')
c4science · Help