new scripts

#!/usr/bin/env python3

'''

/**

 * @file imagej_tiff_saver.py

 * @brief save tiffs for imagej (1.52d+) - with stacks and hyperstacks

 * @par <b>License</b>:

 *  This program is free software: you can redistribute it and/or modify

 *  it under the terms of the GNU General Public License as published by

 *  the Free Software Foundation, either version 3 of the License, or

 *  (at your option) any later version.

 *

 *  This program is distributed in the hope that it will be useful,

 *  but WITHOUT ANY WARRANTY; without even the implied warranty of

 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the

 *  GNU General Public License for more details.

 *

 *  You should have received a copy of the GNU General Public License

 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.

*/

'''

__copyright__ = "Copyright 2018, Elphel, Inc."

__license__   = "GPL-3.0+"

__email__     = "oleg@elphel.com"

from PIL import Image, TiffImagePlugin

import numpy as np

import math

def __get_IJ_IFD(t,z,c):

  ifd = TiffImagePlugin.ImageFileDirectory_v2()

  #is_hyperstack = 'true' if len(shape)>1 else 'false'

  #if (len(shape)>0):

  return ifd

def save(path,images,force_stack=False,force_hyperstack=False):

  # Got images, analyze shape:

  #   - possible formats (c == depth):

  #     -- (t,z,h,w,c)

  #     -- (t,h,w,c), t or z does not matter

  #     -- (h,w,c)

  #     -- (h,w)

  # save single layer image

  if   len(images.shape)==2:

    # (h,w) -> (h,w,c=1)

    images = images[:,:,np.newaxis]

  elif len(images.shape)>2:

    # do nothing

    pass

  else:

    # (w,) -> (h=1,w,c=1)

    images = images[np.newaxis,:,np.newaxis]

  #imlist = np.ravel(Image.fromarray(images))

  t,z,h,w,c  = images.shape

  c_axis = len(images.shape)-1

  channels = np.squeeze(np.split(images,c,axis=c_axis))

  split_channels = np.concatenate(channels,axis=-3)

  images_flat = np.reshape(split_channels,(-1,h,w))

  imlist = []

  for i in range(images_flat.shape[0]):

    imlist.append(Image.fromarray(images_flat[i]))

  imlist[0].save(path,save_all=True,append_images=imlist[1:],tiffinfo=__get_IJ_IFD(t,z,c))

# Testing

if __name__ == "__main__":

  def hamming_window(x,N):

    y = 0.54 - 0.46*math.cos(2*math.pi*x/(N-1))

    return y

  hw = hamming_window

  NT = 5

  NC = 2

  NZ = 3

  NX = 512

  NY = 512

  images = np.empty((NT,NZ,NY,NX,NC))

  import time

  print(str(time.time())+": Generating test images")

  for t in range(NT):

    for z in range(NZ):

      for c in range(NC):

        images[t,z,:,:,c] = np.array([[(255-t*25)*hw(i,512)*hw(j,512) for i in range(NX)] for j in range(NY)],np.float32)

  print(str(time.time())+": Test images generated")

  print("Images shape: "+str(images.shape))

  v = save("result_2.tiff",images)

Open all tiffs in a folder, combine a single tiff from randomly selected

tiles from originals

'''

import tensorflow as tf

#import tensorflow.contrib.slim as slim

from PIL import Image

import time

sys.exit()

#http://stackoverflow.com/questions/287871/print-in-terminal-with-colors-using-python

class bcolors:

    HEADER = '\033[95m'

except IndexError:

  src = "."

print("Importing TensorCrawl")

print_time()

import tensorflow as tf

import tensorflow.contrib.slim as slim

print("TensorCrawl imported")

print_time()

IS_TEST = False

# BEGIN IF IS_TEST

if not IS_TEST:

  tlist = glob.glob(src+"/*.tiff")

print("Found "+str(len(tlist))+" preprocessed tiff files:")

  print("\n".join(tlist))

  print("Found "+str(len(tlist))+" preprocessed tiff files:")

  print_time()

  ''' WARNING, assuming:

        - timestamps and part of names match

  packed_tiles = np.dstack((packed_tiles,values[:,:,0]))

  print("Packed (81x4 -> 1x(25*4+1)) tiled input shape: "+str(packed_tiles.shape))

  print("Values shape "+str(values.shape))

  print_time()

else:

  print("Init test data")

  ptab_name = "tile_packing_table.xml"

  pt = pile.PackingTable(ptab_name,LAYERS_OF_INTEREST).lut

  # 9x9 2 layers, no neighbors

  l = np.zeros((9,9))

  for y,x in itertools.product(range(l.shape[0]),range(l.shape[1])):

    l[y,x] = 9*y + x

  l_value = np.array([2.54,3.54,0.5])

  #print(l)

  l1 = l

  l2 = l*2

  l3 = l*3

  l4 = l*4

  ls = np.dstack((l1,l2,l3,l4))

  #print(ls.shape)

  l_packed_pre = pile.pack_tile(ls,pt)

  #print(l_packed_pre.shape)

  #print(l_packed_pre)

  l_packed = np.hstack((l_packed_pre,l_value[0]))

  #print(l_packed.shape)

  #print(l_packed)

  # use l_packed

  packed_tiles = np.empty([1,1,l_packed.shape[0]])

  values =       np.empty([1,1,2])

  print(packed_tiles.shape)

  print(values.shape)

  packed_tiles[0,0] = l_packed

  values[0,0] = l_value[1:3]

  print(packed_tiles[0,0])

  print(values[0,0])

# END IF IS_TEST

#print("CHECKPOINTE")

#for i in range(tiles.shape[0]):

#  for j in range(tiles.shape[1]):

#    nn_input = pile.get_tile_with_neighbors(tiles,i,j,RADIUS)

#print_time()

result_dir = './result/'

checkpoint_dir = './result/'

save_freq = 500

def lrelu(x):

def network(input):

  fc1 = slim.fully_connected(input,42,activation_fn=lrelu,scope='g_fc1')

  fc2 = slim.fully_connected(fc1,  21,activation_fn=lrelu,scope='g_fc2')

  fc3 = slim.fully_connected(fc2,   1,activation_fn=lrelu,scope='g_fc3')

  fc1 = slim.fully_connected(input,101,activation_fn=lrelu,scope='g_fc1')

  fc2 = slim.fully_connected(fc1,  101,activation_fn=lrelu,scope='g_fc2')

  fc3 = slim.fully_connected(fc2,  101,activation_fn=lrelu,scope='g_fc3')

  fc4 = slim.fully_connected(fc3,  101,activation_fn=lrelu,scope='g_fc4')

  fc5 = slim.fully_connected(fc4,    2,activation_fn=lrelu,scope='g_fc5')

  return fc3

  return fc5

sess = tf.session()

sess = tf.Session()

in_tile = tf.placeholder(tf.float32,[None,None,101])

gt      = tf.placeholder(tf.float32,[None,None,2])

in_tile = tf.placeholder(tf.float32,[None,101])

gt      = tf.placeholder(tf.float32,[None,2])

#cf_cutoff = tf.constant(tf.float32,[None,1])

out = network(in_tile)

G_loss = tf.reduce_mean(tf.abs(out[:,0]-gt[:,0]))

# min cutoff

cf_cutoff = 0.173303

cf_w = tf.pow(tf.maximum(gt[:,1]-cf_cutoff,0.0),1)

cf_wsum = tf.reduce_sum(cf_w)

cf_w_norm = cf_w/cf_wsum

#out_cf = out[:,1]

G_loss = tf.reduce_mean(tf.abs(out[:,0]-cf_w_norm*gt[:,0]))

t_vars=tf.trainable_variables()

lr=tf.placeholder(tf.float32)

G_opt=tf.train.AdamOptimizer(learning_rate=lr).minimize(G_loss,var_list=[var for var in t_vars if var.name.startswith('g_')])

saver=tf.train.Saver()

for folder in allfolders:

  lastepoch = np.maximum(lastepoch, int(folder[-4:]))

#g_loss = np.zeros((,1))

g_loss = np.zeros((packed_tiles.shape[0]*packed_tiles.shape[1],1))

learning_rate = 1e-4

for epoch in range(lastepoch,4001):

print(bcolors.HEADER+"Last Epoch = "+str(lastepoch)+bcolors.ENDC)

for epoch in range(lastepoch,1):

#for epoch in range(lastepoch,4001):

  if os.path.isdir("result/%04d"%epoch):

    continue

  cnt=0

  if epoch > 2000:

    learning_rate = 1e-5

  for ind in np.random.permutation(tiles.shape[0]*tiles.shape[1]):

  for ind in np.random.permutation(packed_tiles.shape[0]*packed_tiles.shape[1]):

    input_patch = tiles[i,j]

    gt_patch = values[i,j,1:2]

    #print("Iteration "+str(cnt))

    st=time.time()

    cnt+=1

    i = int(ind/packed_tiles.shape[1])

    j = ind%packed_tiles.shape[1]

    #input_patch = tiles[i,j]

    input_patch = np.empty((1,packed_tiles.shape[2]))

    input_patch[0] = packed_tiles[i,j]

    gt_patch = np.empty((1,2))

    if not IS_TEST:

      gt_patch[0] = values[i,j,1:3]

    else:

      gt_patch[0] = values[i,j]

    #print(input_patch)

    #print(gt_patch)

    gt_patch[gt_patch==-256] = np.nan

    skip_iteration = False

    # if nan skip run!

    if np.isnan(np.sum(gt_patch[0])):

      skip_iteration = True

    if np.isnan(np.sum(input_patch[0])):

      skip_iteration = True

    if skip_iteration:

      #print(bcolors.WARNING+"Found NaN, skipping iteration for tile "+str(i)+","+str(j)+bcolors.ENDC)

      pass

    else:

      _,G_current,output = sess.run([G_opt,G_loss,out],feed_dict={in_tile:input_patch,gt:gt_patch,lr:learning_rate})

      g_loss[ind]=G_current

      print("%d %d Loss=%.3f CurrentLoss=%.3f Time=%.3f"%(epoch,cnt,np.mean(g_loss[np.where(g_loss)]),G_current,time.time()-st))

    if epoch%save_freq==0:

      if not os.path.isdir(result_dir + '%04d'%epoch):

        os.makedirs(result_dir + '%04d'%epoch)

  saver.save(sess, checkpoint_dir + 'model.ckpt')

print_time()

print(bcolors.OKGREEN+"time: "+str(time.time())+bcolors.ENDC)

#!/usr/bin/env python3

__copyright__ = "Copyright 2018, Elphel, Inc."

__license__   = "GPL-3.0+"

__email__     = "oleg@elphel.com"

'''

Open all tiffs in a folder, combine a single tiff from randomly selected

tiles from originals

'''

from PIL import Image

import os

import sys

import glob

import imagej_tiff as ijt

import pack_tile as pile

import numpy as np

import itertools

import time

#http://stackoverflow.com/questions/287871/print-in-terminal-with-colors-using-python

class bcolors:

    HEADER = '\033[95m'

    OKBLUE = '\033[94m'

    OKGREEN = '\033[92m'

    WARNING = '\033[38;5;214m'

    FAIL = '\033[91m'

    ENDC = '\033[0m'

    BOLD = '\033[1m'

    BOLDWHITE = '\033[1;37m'

    UNDERLINE = '\033[4m'

def print_time():

  print(bcolors.BOLDWHITE+"time: "+str(time.time())+bcolors.ENDC)

# USAGE: python3 test_3.py some-path

VALUES_LAYER_NAME = 'other'

LAYERS_OF_INTEREST = ['diagm-pair', 'diago-pair', 'hor-pairs', 'vert-pairs']

RADIUS = 1

try:

  src = sys.argv[1]

except IndexError:

  src = "."

print("Importing TensorCrawl")

print_time()

import tensorflow as tf

import tensorflow.contrib.slim as slim

print("TensorCrawl imported")

print_time()

result_dir = './result/'

checkpoint_dir = './result/'

save_freq = 500

def lrelu(x):

    return tf.maximum(x*0.2,x)

def network(input):

  fc1 = slim.fully_connected(input,101,activation_fn=lrelu,scope='g_fc1')

  fc2 = slim.fully_connected(fc1,  101,activation_fn=lrelu,scope='g_fc2')

  fc3 = slim.fully_connected(fc2,  101,activation_fn=lrelu,scope='g_fc3')

  fc4 = slim.fully_connected(fc3,  101,activation_fn=lrelu,scope='g_fc4')

  fc5 = slim.fully_connected(fc4,    2,activation_fn=lrelu,scope='g_fc5')

  return fc5

sess = tf.Session()

in_tile = tf.placeholder(tf.float32,[None,101])

gt      = tf.placeholder(tf.float32,[None,2])

out = network(in_tile)

#G_loss = tf.reduce_mean(tf.abs(out[:,0]-gt[:,0]))

#t_vars=tf.trainable_variables()

#lr=tf.placeholder(tf.float32)

#G_opt=tf.train.AdamOptimizer(learning_rate=lr).minimize(G_loss,var_list=[var for var in t_vars if #var.name.startswith('g_')])

saver=tf.train.Saver()

sess.run(tf.global_variables_initializer())

ckpt=tf.train.get_checkpoint_state(checkpoint_dir)

if ckpt:

  print('loaded '+ckpt.model_checkpoint_path)

  saver.restore(sess,ckpt.model_checkpoint_path)

# do not need output for now

#if not os.path.isdir(result_dir + 'final/'):

#    os.makedirs(result_dir + 'final/')

tlist = glob.glob(src+"/*.tiff")

print("\n".join(tlist))

print("Found "+str(len(tlist))+" preprocessed tiff files:")

print_time()

''' WARNING, assuming:

      - timestamps and part of names match

      - layer order and names are identical

'''

# Now PROCESS

for item in tlist:

  print(bcolors.OKGREEN+"Processing "+item+bcolors.ENDC)

  # open the first one to get dimensions and other info

  tiff = ijt.imagej_tiff(item)

  # shape as tiles? make a copy or make writeable

  # (242, 324, 9, 9, 5)

  # get labels

  labels = tiff.labels.copy()

  labels.remove(VALUES_LAYER_NAME)

  print("Image data layers:  "+str(labels))

  print("Layers of interest: "+str(LAYERS_OF_INTEREST))

  print("Values layer: "+str([VALUES_LAYER_NAME]))

  # create copies

  tiles  = np.copy(tiff.getstack(labels,shape_as_tiles=True))

  # need values? Well, just to compare

  values = np.copy(tiff.getvalues(label=VALUES_LAYER_NAME))

  #gt = values[:,:,1:3]

  print("Mixed tiled input data shape: "+str(tiles.shape))

  #print_time()

  # now pack from 9x9 to 1x25

  # tiles and values

  # Parse packing table

  # packing table name

  ptab_name = "tile_packing_table.xml"

  ptab = pile.PackingTable(ptab_name,LAYERS_OF_INTEREST).lut

  # might not need it because going to loop through anyway

  packed_tiles = np.array([[pile.pack_tile(tiles[i,j],ptab) for j in range(tiles.shape[1])] for i in range(tiles.shape[0])])

  packed_tiles = np.dstack((packed_tiles,values[:,:,0]))

  # flatten

  packed_tiles_flat = packed_tiles.reshape(-1, packed_tiles.shape[-1])

  values_flat       = values.reshape(-1, values.shape[-1])

  print("Packed (81x4 -> 1x(25*4+1)) tiled input shape: "+str(packed_tiles_flat.shape))

  print("Values shape "+str(values_flat.shape))

  print_time()

  # now run prediction

  output = sess.run(out,feed_dict={in_tile:packed_tiles_flat})

  print("Output shape: "+str(output.shape))

  # so, let's print

  for i in range(output.shape[0]):

    p  = output[i,0]

    pc = output[i,1]

    fv = values_flat[i,0]

    gt = values_flat[i,1]

    cf = values_flat[i,2]

    vstring = "["+"{0:.2f}".format(fv)+", "+"{0:.2f}".format(gt)+", "+"{0:.2f}".format(cf)+"]"

    pstring = "["+"{0:.2f}".format(p)+", "+"{0:.2f}".format(pc)+"]"

    if not np.isnan(p):

      outstring = "i: "+str(i)+"    Values:  "+vstring+"    Prediction:  "+pstring

      if cf<0.5:

        print(outstring)

      else:

        print(bcolors.WARNING+outstring+bcolors.ENDC)

    #else:

    #  print("i: "+str(i)+" NaNs")

# check histogram:

#import matplotlib.pyplot as plt

#x = np.histogram(values_flat[:,2])

#plt.hist(x, bins=100)

#plt.ylabel('Confidence')

print_time()

print(bcolors.OKGREEN+"time: "+str(time.time())+bcolors.ENDC)