Commit c251b4f0 authored by Jirka's avatar Jirka

class, & S-T-A experiments

parent f4fe5f85
This diff is collapsed.
This diff is collapsed.
......@@ -5,6 +5,7 @@ import glob
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import multiprocessing as mproc
from scipy import ndimage
from skimage import io, draw, transform, filters
# from PIL import Image
......@@ -20,15 +21,15 @@ else:
# PATH_DATA_SYNTH = '/datagrid/Medical/microscopy/drosophila_segmOvary/'
DEFAULT_PATH_DATA = '/datagrid/temporary/Medical/'
DEFAULT_DIR_APD = 'atomicPatternDictionary_v2'
DEFAULT_DIR_APD = 'atomicPatternDictionary_vx'
DEFAULT_PATH_APD = os.path.join(DEFAULT_PATH_DATA, DEFAULT_DIR_APD)
# DEFAULT_IM_SIZE = (512, 512)
# DEFAULT_IM_SIZE = (256, 256)
DEFAULT_IM_SIZE = (128, 128)
# DEFAULT_IM_SIZE = (64, 64)
DEFAULT_NB_PTNS = 25
# DEFAULT_IM_SIZE = (128, 128)
DEFAULT_IM_SIZE = (64, 64)
DEFAULT_NB_PTNS = 3
# DEFAULT_NB_SPLS = 1500
DEFAULT_NB_SPLS = 1200
DEFAULT_NB_SPLS = 200
DEFAULT_DIR_DICT = 'dictionary'
DEFAULT_IM_POSIX = '.png'
DEFAULT_IM_PATTERN = 'pattern_{:03d}'
......@@ -153,21 +154,20 @@ def dictionary_generate_atlas(p_out=DEFAULT_PATH_APD, nb=DEFAULT_NB_PTNS,
plt.show()
atlas_final = np.zeros(im_size, dtype=np.uint8)
# export to dictionary
idx, l_imgs = 0, []
logger.debug('... post-processing over generated patterns: {}'.format(
l_imgs = []
logger.info('... post-processing over generated patterns: {}'.format(
np.unique(atlas).tolist()))
for id in np.unique(atlas):
for i, id in enumerate(np.unique(atlas)[1:]):
im = np.zeros(im_size, dtype=np.uint8)
# create pattern
im[atlas== (id+1)] = 1
im[atlas == (id)] = 1
# remove all smaller unconnected elements
im = extract_image_largest_element(im)
if np.sum(im) == 0: continue
l_imgs.append(im)
export_image(out_dir, im, idx, DEFAULT_IM_PATTERN)
export_image(out_dir, im, i, DEFAULT_IM_PATTERN)
# add them to the final arlas
atlas_final[im==1] = id
idx += 1
atlas_final[im==1] = i + 1
export_image(out_dir, atlas_final, 'atlas')
return l_imgs
......@@ -262,8 +262,10 @@ def export_image(outDir, im, im_name, n_template=DEFAULT_IM_SEGM):
if not type(im_name)==str:
im_name = n_template.format(im_name)
p_img = os.path.join(outDir, im_name + DEFAULT_IM_POSIX)
logger.debug(' .. saving image "{}...{}"'.format(p_img[:25], p_img[-25:]))
io.imsave(p_img, im/float(np.max(im))*1.)
im_norm = im / (np.max(im) * 1.)
logger.debug(' .. saving image {} to "{}...{}"'.format(im_norm.shape,
p_img[:25], p_img[-25:]))
io.imsave(p_img, im_norm)
# Image.fromarray(im).save(pImg)
return p_img
......@@ -315,7 +317,7 @@ def image_transform_binary2prob(im, coef=0.1):
:param coef: float, influence hoe strict the boundary between F-B is
:return: np.array<w, h> float image
"""
logger.debug('... transform binary image to probability')
logger.info('... transform binary image to probability')
im_dist = ndimage.distance_transform_edt(im)
im_dist -= ndimage.distance_transform_edt(1-im)
im_prob = 1. / (1. + np.exp(-coef * im_dist))
......@@ -333,7 +335,7 @@ def add_image_prob_noise(im, ration=0.1):
:param ration: float (0, 1) means 0 = no noise
:return: np.array<w, h> float image
"""
logger.debug('... add smooth noise to a probability image')
logger.info('... add smooth noise to a probability image')
rnd = 2* (np.random.random(im.shape) -0.5)
rnd[abs(rnd) > ration] = 0
im_noise = np.abs(im - rnd)
......@@ -363,7 +365,8 @@ def dataset_prob_construct(imgs, out_dir, coef=0.5):
def dataset_load_images(name=DEFAULT_DATASET, path_base=DEFAULT_PATH_APD,
im_ptn='', im_posix=DEFAULT_IM_POSIX, nb_spls=None):
im_ptn='*', im_posix=DEFAULT_IM_POSIX, nb_spls=None,
nb_jobs=1):
""" load complete dataset or just a subset
:param name: str, name od particular dataset
......@@ -375,17 +378,34 @@ def dataset_load_images(name=DEFAULT_DATASET, path_base=DEFAULT_PATH_APD,
"""
p_dir = os.path.join(path_base, name)
logger.info('loading folder ({}) <- "{}"'.format(os.path.exists(p_dir), p_dir))
p_imgs = glob.glob(os.path.join(p_dir, im_ptn + '*' + im_posix))
p_search = os.path.join(p_dir, im_ptn + im_posix)
logger.debug('image search "{}"'.format(p_search))
p_imgs = glob.glob(p_search)
logger.debug('number spls {} in dataset "{}"'.format(len(p_imgs), name))
p_imgs = sorted(p_imgs)[:nb_spls]
imgs = []
for i, pIm in enumerate(p_imgs):
im = io.imread(pIm)
imgs.append(im / float(np.max(im)))
if nb_jobs > 1:
logger.debug('running in {} threads...'.format(nb_jobs))
mproc_pool = mproc.Pool(nb_jobs)
imgs = mproc_pool.map(load_image, p_imgs)
mproc_pool.close()
mproc_pool.join()
else:
logger.debug('running in single threads...')
imgs = []
for i, p_im in enumerate(p_imgs):
imgs.append(load_image(p_im))
im_names = [os.path.splitext(os.path.basename(p))[0] for p in p_imgs]
return imgs, im_names
def load_image(p_img):
im = io.imread(p_img)
im = im / float(np.max(im))
return im
def dataset_load_weights(path_base=DEFAULT_PATH_APD, n_file=DEFAULT_NAME_WEIGHTS):
""" loading all true wieghts for given dataset
......@@ -405,7 +425,7 @@ def dataset_load_weights(path_base=DEFAULT_PATH_APD, n_file=DEFAULT_NAME_WEIGHTS
def dataset_create_atlas(name='dictionary', path_base=DEFAULT_PATH_APD,
im_ptn='pattern_'):
im_ptn='pattern_*'):
""" load all independent patterns and compose them into single m-label atlas
:param name: str, name of dataset
......@@ -421,9 +441,10 @@ def dataset_create_atlas(name='dictionary', path_base=DEFAULT_PATH_APD,
return np.array(atlas, dtype=np.uint8)
def dataset_export_images(p_out, imgs, names=None):
def dataset_export_images(p_out, imgs, names=None, nb_jobs=1):
""" export complete dataset
:param parallel: bool
:param p_out: str
:param imgs: [np.array<w, h>]
:param names: [str] or None (use indexes)
......@@ -433,11 +454,29 @@ def dataset_export_images(p_out, imgs, names=None):
logger.debug('export {} images into "{}"'.format(len(imgs), p_out))
if names is None:
names = range(len(imgs))
for i, im in enumerate(imgs):
export_image(p_out, im, names[i])
if nb_jobs > 1:
logger.debug('running in {} threads...'.format(nb_jobs))
sp = [(p_out, im, names[i]) for i, im in enumerate(imgs)]
mproc_pool = mproc.Pool(nb_jobs)
mproc_pool.map(export_image_sp, sp)
mproc_pool.close()
mproc_pool.join()
else:
logger.debug('running in single threads...')
for i, im in enumerate(imgs):
export_image(p_out, im, names[i])
p_npz = os.path.join(p_out, 'input_images.npz')
np.savez(open(p_npz, 'w'), imgs)
return
def export_image_sp(sp):
p_out, im, name = sp
export_image(p_out, im, name)
return None
def dataset_convert_nifti(path_in, path_out, posix=DEFAULT_IM_POSIX):
""" having a datset of png images conver them into nifti images
......
......@@ -8,6 +8,8 @@ run experiments with Atomic Learning Pattern Encoding
# matplotlib.use('Agg')
import os
import gc
import time
import numpy as np
import generate_dataset as gen_data
import multiprocessing as mproc
......@@ -22,16 +24,19 @@ if jirka:
else:
DEFAULT_PATH_DATA = '/datagrid/Medical/microscopy/drosophila/real_segmentations'
REAL_DATASET_NAME = '108_genes_expression'
# BINARY_DATASET_POSIX = '_binary'
DEFAULT_NB_THREADS = int(mproc.cpu_count() * 0.8)
# REAL_DATASET_NAME = '108_genes_expression'
REAL_DATASET_NAME = '1000_imgs'
IMAGE_PATTERN = '*_exp'
DEFAULT_PARAMS = {
'computer': os.uname(),
'in_path': DEFAULT_PATH_DATA,
'dataset': REAL_DATASET_NAME,
'out_path': DEFAULT_PATH_DATA,
}
IMAGE_BINARY_THRESHOLD = 0.5
BINARY_POSIX = '_binary'
NB_THREADS = int(mproc.cpu_count() * 0.7)
IMAGE_BINARY_THRESHOLD = 0.95
def extend_images(imgs):
......@@ -42,11 +47,13 @@ def extend_images(imgs):
:return: [np.array<w, h>]
"""
im_sizes = [im.shape for im in imgs]
logger.debug('different image sizes: {}'.format(set(im_sizes)))
if len(set(im_sizes)) == 1:
return imgs
w_max = max([w for w, h in im_sizes])
h_max = max([h for w, h in im_sizes])
m_size = (w_max, h_max)
logger.debug('max image size: {}'.format(m_size))
for i, im in enumerate(imgs):
if im.shape == m_size:
continue
......@@ -57,27 +64,35 @@ def extend_images(imgs):
return imgs
def find_borders(im_sum):
def find_borders(im_mean, v_lim=0.):
""" find from sided the rows and col where the image set is zero
:param im_sum: np.array<w, h>
:return:
"""
border = {}
for i in range(len(im_sum.shape)):
idx = np.nonzero(im_sum.sum(axis=i))[0]
border[i] = min(idx) - 1, max(idx) + 1
logger.debug('crop value threshold: {}'.format(v_lim))
for i in range(len(im_mean.shape)):
vals = im_mean.mean(axis=i)
vals[vals < v_lim] = 0
idx = np.nonzero(vals)[0]
if len(idx) > 0:
border[i] = min(idx), max(idx)
else:
logger.warning('nothing to chose in axis: {}'.format(i))
border[i] = 0, -1
return border
def crop_images(imgs):
def crop_images(imgs, v_lim=0.):
""" try to cut out image rows and colums that are useless
:param imgs: [np.array<w, h>]
:return: [np.array<w, h>]
"""
im_sum = np.array(imgs).sum(axis=0)
border = find_borders(im_sum)
logger.info('perform image crop...')
im_mean = np.array(imgs).mean(axis=0)
border = find_borders(im_mean, v_lim)
logger.debug('image crop limits: {}'.format(border))
imgs_crop = [None] * len(imgs)
logger.debug('crop all {} images'.format(len(imgs)))
......@@ -101,21 +116,22 @@ def threshold_image_otsu(img):
def threshold_image_adapt(img):
img_th = filters.threshold_adaptive(img, 120)
selem = morphology.disk(15)
img_th = filters.threshold_adaptive(img, 450)
selem = morphology.disk(12)
img_th = morphology.opening(img_th, selem)
return img_th
def threshold_images(imgs, fn_th):
def threshold_images(imgs, fn_th, nb_jobs=NB_THREADS):
""" threshold images by specific level
:param imgs: [np.array<w, h>]
:param th: float
:return: [np.array<w, h>]
"""
logger.debug('binary images...')
mproc_pool = mproc.Pool(DEFAULT_NB_THREADS)
logger.info('binary images by "{}"...'.format(fn_th.__name__))
logger.debug('running in {} threads...'.format(nb_jobs))
mproc_pool = mproc.Pool(nb_jobs)
imgs_th = mproc_pool.map(fn_th, imgs)
mproc_pool.close()
mproc_pool.join()
......@@ -127,26 +143,30 @@ def binarise_main(params=DEFAULT_PARAMS):
:return:
"""
imgs, names = gen_data.dataset_load_images(params['dataset'], params['in_path'])
imgs, names = gen_data.dataset_load_images(params['dataset'], params['in_path'],
im_ptn=IMAGE_PATTERN)
logger.debug('loaded {} images of size {}'.format(len(imgs), imgs[0].shape))
imgs = extend_images(imgs)
imgs = crop_images(imgs)
imgs = crop_images(imgs, 1e-3)
gc.collect(), time.sleep(1)
p_export = os.path.join(params['in_path'], params['dataset'] + '_binary')
p_export = os.path.join(params['in_path'], params['dataset'] + BINARY_POSIX)
if not os.path.exists(p_export):
os.mkdir(p_export)
imgs_th = threshold_images(imgs, threshold_image_fix)
p_out = os.path.join(p_export, 'binary-fix')
gen_data.dataset_export_images(p_out, imgs_th, names)
p_out = os.path.join(p_export, 'binary-fix_{}'.format(IMAGE_BINARY_THRESHOLD))
gen_data.dataset_export_images(p_out, imgs_th, names, nb_jobs=NB_THREADS)
gc.collect(), time.sleep(1)
imgs_th = threshold_images(imgs, threshold_image_otsu)
p_out = os.path.join(p_export, 'binary-otsu')
gen_data.dataset_export_images(p_out, imgs_th, names)
# imgs_th = threshold_images(imgs, threshold_image_otsu)
# p_out = os.path.join(p_export, 'binary-otsu')
# gen_data.dataset_export_images(p_out, imgs_th, names, nb_jobs=NB_THREADS)
# gc.collect(), time.sleep(1)
imgs_th = threshold_images(imgs, threshold_image_adapt)
p_out = os.path.join(p_export, 'binary-adapt')
gen_data.dataset_export_images(p_out, imgs_th, names)
gen_data.dataset_export_images(p_out, imgs_th, names, nb_jobs=NB_THREADS)
return
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment