Laboratoire RDP » Biophysics Team » FLORIVAR

#!/usr/bin/env python

# before lounching the executable activate the timagetk conda environement

# conda activate florivar

# general imports :

from sys import path, argv

import os

import time

import numpy as np

from skimage.filters import threshold_otsu

from titk_tools.io import imread, imsave, SpatialImage

#from timagetk.plugins import linear_filtering

from scipy.ndimage import gaussian_filter

# multiprocessing

from multiprocessing import Pool

nproc = 20 # number of processors

print('=======================================')

def gauss_titk(img, radius):

	'''

	gaussian filter with radius in real units

	- if radius=a number : isotropic

	- if radius= a list of numbers : anisotropic

	'''

	img_filtered = linear_filtering(img, std_dev=radius, method='gaussian_smoothing')#, real=True)

	return img_filtered

def gauss_nd(img, radius):

	'''

	gaussian filter with radius in real units

	- if radius=a number : isotropic

	- if radius= a list of numbers : anisotropic

	'''

	img_filtered = gaussian_filter(img, radius)

	return img_filtered

def iterated_gauss_2d(img):

	vs=img.voxelsize

	rad=2*vs[2]/vs[0]

	#rad=10

	radius = [rad, rad, 0]

	img_filtered=gaussian_filter(img, radius)

	rad=rad/2

	while rad>1:

		print(rad)

		radius=[rad,rad,0]

		img_filtered = gaussian_filter(img_filtered, radius)

		rad=rad/2

	return SpatialImage(img_filtered, dtype='uint8', voxelsize=vs, origin=[0, 0, 0])

indir = argv[1]

T = int(argv[2])

newdir=indir+'_normalised-OtsuT'+str(T)

os.system("mkdir -p "+newdir)

print('Normalising images in the folder:',indir)

print('Common Otsu threshold for normalisation:',T)

print('Results will be written in the folder:',newdir)

pathdir=indir

listfiles=os.listdir(pathdir)

def treat_file(i):

	filename = listfiles[i]

	extension = filename.split('.')[-1]

	if extension=='tif':

		print('--------------------------')

		print(filename)

		im = imread(pathdir+'/'+filename)

		# filter for computing Otsu

		#imf = iterated_gauss_2d(im)

		imf=im

		if len(np.unique(imf))>0:#20:

			iotsu = threshold_otsu(imf)

			# normalise the original image

			imin = np.min(im)+0.0

			imax= np.max(im)+0.0

			a = T/(iotsu-imin)

			b = - imin*a

			transformed_im = SpatialImage(a*im+b, dtype='uint16', voxelsize=im.voxelsize, origin=[0, 0, 0])

			transformed_im[transformed_im>255]=255 		

			transformed_im = SpatialImage(transformed_im, dtype='uint8', voxelsize=im.voxelsize, origin=[0, 0, 0])

			tiotsu = threshold_otsu(transformed_im)

			imf = iterated_gauss_2d(transformed_im)

			tiotsuf = threshold_otsu(imf)

			print('   otsu=',iotsu,'  otsu normalised=',tiotsu,'  otsu normalised filtered=',tiotsuf)

			newname = 'norm'+str(T)+'-'+filename

			newnamef = 'filtered-'+newname

			#imsave(newdir+'/'+newname,transformed_im)

			imsave(newdir+'/'+newnamef,imf)

	return 0

inputs = range(len(listfiles))

for i in range(len(inputs)):

	treat_file(i)

'''

pool = Pool(processes=nproc)

pool.map(treat_file, inputs)

pool.close()

pool.join()

'''