Commit d2aa5622 authored by Martinez-Carvajal German's avatar Martinez-Carvajal German
Browse files

All these concatenating scripts were deleted.

Concatenation of raw image (when you have to
concatenate scans that had to be splitted in two
or more parts during at the laboratory) is not as
simple as one may think to do it based on computations

The best way is to do it manually, and to do it
choosing a slice that will not present noise.

I say noise, becasie whe you do several scans of
an object, at the beggining and at the end of the
scan you may find an artifact like a shadow
sweeping the image. Hopefully, this shadow
disappears at the center of the image.

In consequence, when you do a manual concatenation,
pick the most inner image you can !

Don't forget there is a script, that concatenates
all parts once you have found the indexes of the
slices where you want to concatenate
parent 01ff112f
# Connect_Image.py
def EnhanceContrast(Imp, Saturation):
Argument = "saturated=" + str(Saturation) + " normalize process_all use"
IJ.run(Imp,"Enhance Contrast...", Argument)
def Concatenate_Scans (Folder, File_Top_Image, File_Bot_Image, Saturation):
""" Connects two tomographies based on an intersection zone """
""" Top and bottom will be connected """
""" Folder is the reading and saving directory """
""" Saturation is the % of saturation in the for the Enhance Contrast routine """
""" This program must be executed in FIJI/ImageJ """
# Files Name
# Saving Folder
from os import path
assert (path.exists(File_Top_Image)) , "FILE DOES NOT EXIST: " + File_Top_Image
assert (path.exists(File_Bot_Image)) , "FILE DOES NOT EXIST: " + File_Bot_Image
assert (path.isdir(Folder)) , "Saving Folder IS NOT A FOLDER: " + folder
# Open the topmost part of the image
from ij import IJ
# New instance of ImagePlus with the topmost part of the image
print('Reading Top Image')
ImpTop = IJ.openImage(File_Top_Image)
ImpTop.show()
print(ImpTop)
# New instance of ImagePlus with the bottom part of the image
print('Reading Bot Image')
ImpBot = IJ.openImage(File_Bot_Image)
ImpBot.show()
print(ImpBot)
# Normalize by stretching
from ij import ImagePlus
print("Enhancing Contrast")
EnhanceContrast(ImpTop, Saturation = Saturation) # Saturation in %
EnhanceContrast(ImpBot, Saturation = Saturation) # Saturation in %
# FIND THE CONNECTING SLICE
# Selecting the last slice of the top part pf the image, (the last Image Processor instance)
LastSlice = ImpTop.getNSlices() # "Obiously" the number of slices is the index of the last slice
print("# of the Last Slice (Top Image) = ", LastSlice)
StackTop = ImpTop.getImageStack()
LastIp = StackTop.getProcessor(LastSlice) # Now you have the stack you can get its processor
ImpLastSlice = ImagePlus("LastSlice", LastIp)
#ImpLastSlice.show()
# Perform the difference between the last slice and all the slices of the bottom image
# Be careful, I figured out that the minus(-) operation is not defined for the object_type "Image Processor"
# You should call the Image Calculator
from ij.plugin import ImageCalculator
print("Running Image Calculator - Subtract")
Ic = ImageCalculator()
Imp_Diff = Ic.run("Subtract create stack", ImpBot, ImpLastSlice)
Imp_Diff.title = "Subtraction"
#Imp_Diff.show()
# We need a function to get the mean of an ImagePlus Object
from ij.process import ImageStatistics as IS
def getMean(imp):
""" Return MEAN for the given ImagePlus """
ip = imp.getProcessor()
global options
stats = IS.getStatistics(ip, IS.MEAN, imp.getCalibration())
return stats.mean
# The next piece of code will calculate the mean of each slice of Imp_Diff
# First we need to create an iterable made of Imp objects
print('Getting stack of the Subtraction Image')
Stack_Diff = Imp_Diff.getImageStack()
Slices = []
print('Creating Iterable with the Processors of the Subtraction Syack')
for i in xrange(1,Stack_Diff.getSize()+1): # The index must begin in 1, the last index is not taken into account, tha is the reason why we write +1
#print ("i is", i)
# get the Processor Slice at index i
Ip = Stack_Diff.getProcessor(i)
Slices.append(ImagePlus(str(i), Ip))
# Now we call the getMean function for each slice
print('Computing means of each slice of Subtraction')
MeansBottom = map(getMean, Slices)
Minimum = reduce(min, MeansBottom)
print("Minimum of the means is: ", Minimum)
# Locate the slice correspoding to the minimum difference
for i in xrange(0, len(MeansBottom)): # Means Bottom is a list, its indexing starts at 0
if Minimum == MeansBottom[i]:
#print(MeansBottom[i])
#print(i)
ConnectingSlice = i+1
print("Connecting Slice:", ConnectingSlice)
# Creating the new stack
from ij import ImageStack
StackNew = ImageStack(ImpTop.width, ImpTop.height)
print('Creating the new Stack for the concatenation')
for i in xrange(1, ImpTop.getNSlices()+1):
Ip = StackTop.getProcessor(i)
StackNew.addSlice(Ip)
StackBot = ImpBot.getImageStack()
for i in xrange(ConnectingSlice+1, ImpBot.getNSlices()+1):
Ip = StackBot.getProcessor(i)
StackNew.addSlice(Ip)
ImpNew = ImagePlus("Concatenated", StackNew)
# IJ.resetMinAndMax(ImpNew) # This will adjust the LUT of the 16-bit image without modyfing the pixel values;
# Saving
print('Saving')
from ij.io import FileSaver
fs = FileSaver(ImpNew)
filepath = Folder + "/" + "Concatenation.tif"
fs.saveAsTiff(filepath)
print(ImpNew)
ImpNew.show()
print('Finished Ok')
# MAIN
Folder = "/home/german.martinez-carvajal/Desktop/2017/Concatenation/Test"
File_Top_Image = Folder + "/35MIC_haut_sample.tif"
File_Bot_Image = Folder + "/35MIC_bas_sample.tif"
Saturation = 0.1 # in %
Concatenate_Scans (Folder, File_Top_Image, File_Bot_Image, Saturation)
# Python version of Concatenation of tomographies
import skimage.exposure
import numpy as np
import tifffile
from os import chdir
import matplotlib.pyplot as plt
def Normalization2D(Image2D, Tol):
# Stretching
print ('Stretching')
v_min, v_max = np.percentile(Image2D, (Tol, 100 - Tol)) # values for normalization defined by a tolerance in saturation
print ('Limits of Normalization are ', v_min, v_max)
Image2D = skimage.exposure.rescale_intensity(Image2D, in_range = (v_min, v_max), out_range = 'dtype')
return (Image2D)
def Normalization3D(Image3D, Tol):
# Stretching
print ('Stretching')
v_min, v_max = np.percentile(Image3D, (Tol, 100 - Tol)) # values for normalization defined by a tolerance in saturation
print ('Limits of Normalization are ', v_min, v_max)
Num_Slices = len(Image3D)
for index in range(Num_Slices):
if index % 10 == 0:
print('% progresion = ', int(index/Num_Slices*100))
Image3D[index] = skimage.exposure.rescale_intensity(Image3D[index], in_range = (v_min, v_max), out_range = 'dtype')
return (Image3D)
def sector_mask(shape,centre,radius,angle_range):
"""
CODE FROM STACK OVER FLOW
Return a boolean mask for a circular sector. The start/stop angles in
`angle_range` should be given in clockwise order.
"""
x,y = np.ogrid[:shape[0],:shape[1]]
cx,cy = centre
tmin,tmax = np.deg2rad(angle_range)
# ensure stop angle > start angle
if tmax < tmin:
tmax += 2*np.pi
# convert cartesian --> polar coordinates
r2 = (x-cx)*(x-cx) + (y-cy)*(y-cy)
theta = np.arctan2(x-cx,y-cy) - tmin
# wrap angles between 0 and 2*pi
theta %= (2*np.pi)
# circular mask
circmask = r2 <= radius*radius
# angular mask
anglemask = theta <= (tmax-tmin)
return circmask*anglemask
def plot_hist(a, file_name):
fig, ax = plt.subplots()
hist, bins_ = np.histogram(a, bins = 256)
ax.plot(bins_[1:], hist)
ax.set_xlabel('Gray value')
ax.set_ylabel('Frecuency')
fig.set_size_inches(7, 7)
fig.savefig('Hist_of_{}.png'.format(file_name[:-4]), dpi = 300, )
def Concatenate(Dir, Top, Bottom, Tol) :
"""Connects two tomographies based on an intersection zone
Top and Bottom will be connected
inputs
Top, botom: str, contain the name of the files (ex: tomographie.tif)
Dir is the reading and saving directory
Tol is the % of saturation in the for the Normalizaition routine
This program must be executed in python """
global Con_Slice, Means, Min # connecting slice, vector of mean values, min = min of vector Means
# Reading files
chdir(Dir)
print("Reading Top File")
Image1 = tifffile.imread(Top)
print("Reading Bottom File")
Image2 = tifffile.imread(Bottom)
# Histogram of original files
print("Computing Histograms of original images")
# Image1
plot_hist(Image1, Top)
# Image2
plot_hist(Image2, Bottom)
# Original dtype =
dtype_orig = Image1.dtype
# heights image1 and image2
h1 = np.shape(Image1)[0]
h2 = np.shape(Image2)[0]
# Center
# Images 1 and 2 (3D) should have the same dimensions x*y (h may be different)
x,y = np.shape(Image2)[1],np.shape(Image2)[2]
center = int(x/2), int(y/2)
radius = int(3/4*(x/2))
#Selecting circular ROI (region of interest)
mask_2D = sector_mask((x,y),center,radius,(0,360))
# Creating 3D mask - you should PROBABLY get rid of things out of your ROI, (i.e outside the cylinder in tomographies
# For image 2
mask_3D = np.zeros(Image2.shape, dtype = np.bool)
for i in range(h2):
mask_3D[i] = mask_2D
# Last slice of image 1
last_slice = Image1[-1]
# Focus on the region of interest thanks to the masks
ROI_image1 = np.copy(last_slice)
ROI_image1[~mask_2D] = 0
ROI_image2 = np.copy(Image2)
ROI_image2[~mask_3D] = 0
# Computing histograms of ROI
print("Computing Histograms of images with circular mask")
plot_hist(ROI_image1[mask_2D], '{}_1_mask.png'.format(Top[:-4]))
plot_hist(ROI_image2[mask_3D], '{}_1_mask.png'.format(Bottom[:-4]))
#Normalization of contrast
print("Enhancing contrast by normalization")
ROI_image1 = Normalization2D(ROI_image1, Tol)
#tifffile.imsave('Top_Norm', Image1)
ROI_image2 = Normalization3D(ROI_image2, Tol)
#tifffile.imsave('Bot_Norm', Image2)
# Computing histograms of normalized ROIs
plot_hist(ROI_image1[mask_2D], '{}_2_norm_mask.png'.format(Top[:-4]))
plot_hist(ROI_image2[mask_3D], '{}_2_norm_mask.png'.format(Bottom[:-4]))
print("Changing format to np.int")
# To perform a difference format must be int at least (not int8 or int16)
ROI_image1 = np.array(ROI_image1, np.int)
ROI_image2 = np.array(ROI_image2, np.int)
# Performing difference of image2 - last slice of image1
print("Computing difference")
Difference = np.zeros(ROI_image2.shape, dtype = np.int )
Difference = abs(ROI_image2-ROI_image1)
Difference = np.array(Difference, dtype = dtype_orig)
#print("Saving difference")
#tifffile.imsave('Difference_python', Difference)
# Finding connecting slice
Means = np.mean(Difference, axis = (1,2))
del Difference
Min = np.min(Means)
Con_Slice = np.where(Means == Min)
Con_Slice = np.array(Con_Slice)
assert (Con_Slice.shape) == (1,1), 'Impossible to connect - there are two or more connecting slices'
Con_Slice = int(Con_Slice)
print('Connecting slice : ' , Con_Slice)
# Creating a file that will contain the connecting slices
file_results = open ('connecting_slices_txt', mode = 'a')
file_results.write('Connecting slice for {} and {} is :\n'.format(Top[:-4],Bottom[:-4]))
file_results.write(str(Con_Slice))
file_results.write("\n")
file_results.close()
# Eliminating Intersection
# global Bottom
Image2 = Image2[Con_Slice+1:] # eliminates the intersection
Image2 = np.concatenate((Image1,Image2))
print('Saving')
File_Name3 = '{}_and_{}.tif'.format(FileTop[:-4],FileBot[:-4])
tifffile.imsave(File_Name3, Image2)
print('Image saved in: ', Dir)
#################################################################
# MAIN
#################################################################
plt.close('all')
Dir = '/home/german.martinez-carvajal/Desktop/These/Concatenate_Tomographies/Test2_exterieur_tres_bruite'
FileTop = 'A' + '.tif' # top
FileBot = 'B' + '.tif' # bottom
Tolerance = 0.1 # in %
Concatenate(Dir, FileTop, FileBot, Tol = Tolerance)
print('Finished!')
plt.close('all')
# Connect_ImageV2.py
def Concatenate_Scans (Folder, File_Top_Image, File_Bot_Image):
""" Connects two tomographies based on an intersection zone """
""" Top and bottom will be connected """
""" Folder is the reading and saving directory """
""" This program must be executed in FIJI/ImageJ """
""" Inputs
- str: Folder, File_Top_Image, File_Top_Image = complete directories of the image
IN THIS VERSION THE "TOP" AND "BOTTOM'" IMAGES ARE ALREADY NORMALIZED
(This means that the little code section for normalization in the previous version of this function was deleted)
"""
# Files Name verification
# Saving Folder verification
from os import path
assert (path.exists(File_Top_Image)) , "FILE DOES NOT EXIST: " + File_Top_Image
assert (path.exists(File_Bot_Image)) , "FILE DOES NOT EXIST: " + File_Bot_Image
assert (path.isdir(Folder)) , "Saving Folder IS NOT A FOLDER: " + folder
# Open the topmost part of the image
from ij import IJ
# New instance of ImagePlus with the topmost part of the image
print('Reading Top Image')
ImpTop = IJ.openImage(File_Top_Image)
#ImpTop.show()
print(ImpTop)
# New instance of ImagePlus with the bottom part of the image
print('Reading Bot Image')
ImpBot = IJ.openImage(File_Bot_Image)
#ImpBot.show()
print(ImpBot)
# Normalize by stretching
from ij import ImagePlus
print("Enhancing Contrast")
# FIND THE CONNECTING SLICE
# Selecting the last slice of the top part pf the image, (the last Image Processor instance)
LastSlice = ImpTop.getNSlices() # "Obiously" the number of slices is the index of the last slice
print("# of the Last Slice (Top Image) = ", LastSlice)
StackTop = ImpTop.getImageStack()
LastIp = StackTop.getProcessor(LastSlice) # Now you have the stack you can get its processor
ImpLastSlice = ImagePlus("LastSlice", LastIp)
#ImpLastSlice.show()
# Perform the difference between the last slice and all the slices of the bottom image
# Be careful, I figured out that the minus(-) operation is not defined for the object_type "Image Processor"
# You should call the Image Calculator
from ij.plugin import ImageCalculator
print("Running Image Calculator - Subtract")
Ic = ImageCalculator()
Imp_Diff = Ic.run("Subtract create stack", ImpBot, ImpLastSlice)
Imp_Diff.title = "Subtraction"
Imp_Diff.show()
# We need a function to get the mean of an ImagePlus Object
from ij.process import ImageStatistics as IS
def getMean(imp):
""" Return MEAN for the given ImagePlus """
ip = imp.getProcessor()
global options
stats = IS.getStatistics(ip, IS.MEAN, imp.getCalibration())
return stats.mean
# The next piece of code will calculate the mean of each slice of Imp_Diff
# First we need to create an iterable made of Imp objects
print('Getting stack of the Subtraction Image')
Stack_Diff = Imp_Diff.getImageStack()
Slices = []
print('Creating Iterable with the Processors of the Subtraction Syack')
for i in xrange(1,Stack_Diff.getSize()+1): # The index must begin in 1, the last index is not taken into account, tha is the reason why we write +1
#print ("i is", i)
# get the Processor Slice at index i
Ip = Stack_Diff.getProcessor(i)
Slices.append(ImagePlus(str(i), Ip))
# Now we call the getMean function for each slice
print('Computing means of each slice of Subtraction')
MeansBottom = map(getMean, Slices)
Minimum = reduce(min, MeansBottom)
print("Minimum of the means is: ", Minimum)
# Locate the slice correspoding to the minimum difference
for i in xrange(0, len(MeansBottom)): # Means Bottom is a list, its indexing starts at 0
if Minimum == MeansBottom[i]:
#print(MeansBottom[i])
#print(i)
ConnectingSlice = i+1
print("Connecting Slice:", ConnectingSlice)
# Creating the new stack
from ij import ImageStack
StackNew = ImageStack(ImpTop.width, ImpTop.height)
print('Creating the new Stack for the concatenation')
for i in xrange(1, ImpTop.getNSlices()+1):
Ip = StackTop.getProcessor(i)
StackNew.addSlice(Ip)
StackBot = ImpBot.getImageStack()
for i in xrange(ConnectingSlice+1, ImpBot.getNSlices()+1):
Ip = StackBot.getProcessor(i)
StackNew.addSlice(Ip)
ImpNew = ImagePlus("Concatenated", StackNew)
# IJ.resetMinAndMax(ImpNew) # This will adjust the LUT of the 16-bit image without modyfing the pixel values;
# Saving
print('Saving')
from ij.io import FileSaver
fs = FileSaver(ImpNew)
filepath = Folder + "/" + "Concatenation_with_Jython.tif"
fs.saveAsTiff(filepath)
print(ImpNew)
ImpNew.show()
print('Finished Ok')
# MAIN
#Folder = "/home/german.martinez-carvajal/Desktop/2018/Concatenate_Tomographies/Test"
#File_Top_Image = Folder + "/Normalization_of_sample_top.tif"
#File_Bot_Image = Folder + "/Normalization_of_sample_bottom.tif"
#Concatenate_Scans (Folder, File_Top_Image, File_Bot_Image)
Folder = "/home/german.martinez-carvajal/Desktop/2017/Tomographie/170411COP/Raw_Data/Low_resolution"
File_Top_Image = Folder + "/Normalization_of_COP01_HAUT_70MIC.tif"
File_Bot_Image = Folder + "/Normalization_of_COP01_BAS_70MIC.tif"
Concatenate_Scans (Folder, File_Top_Image, File_Bot_Image)
# -*- coding: utf-8 -*-
"""
Created on Tue Aug 6 11:09:48 2019
@author: german.martinez-carvajal
"""
from os import chdir
import tifffile
import numpy as np
def concatenate(part1, part2, saf):
my_slice = part1[-saf]
my_slice = np.array(my_slice, dtype = np.int)
difference = np.zeros(part2.shape, dtype = np.int)
for z in range(len(part2)):
A = np.array(part2[z], dtype = np.int)
difference[z] = np.abs(A - my_slice)
means = np.mean(difference, axis = (1,2))
connecting_slice = np.argmin(means)
part1_cut = part1[:-saf]
part2_cut = part2[connecting_slice:]
concatenation = np.concatenate((part1_cut, part2_cut))
difference = np.array(difference, dtype = np.uint8)
return(difference, concatenation, connecting_slice)
def main():
# slice avant la fin
saf = 130
# concatenating segmented files
path_read = "/mnt/storage/Thesis_GDMC_XRayTomography/Segmentations/Manip_Novembre"
chdir(path_read)
sample = 1
day = '00'
list_connecting_slices = []
file_part1 = 'Hash_Segm_of_MON{}_J{}_scan{}*'.format(sample, day, 1)
file_part2 = 'Hash_Segm_of_MON{}_J{}_scan{}*'.format(sample, day, 2)
print('Reading files ...', file_part1, file_part2)
part1 = tifffile.imread(file_part1)
part2 = tifffile.imread(file_part2)
print('Concatenation ...')
d,c,cs = concatenate(part1, part2, saf)
list_connecting_slices.append(cs)
print('Reading files ...')
file_part3 = 'Hash_Segm_of_MON{}_J{}_scan{}*'.format(sample, day, 3)
part3 = tifffile.imread(file_part3)
print('Concatenation ...')
d,c,cs = concatenate(c, part3, saf)
list_connecting_slices.append(cs)
print('Reading files ...')
file_part4 = 'Hash_Segm_of_MON{}_J{}_scan{}*'.format(sample, day, 4)
part4 = tifffile.imread(file_part4)
print('Concatenation ...')
d,c,cs = concatenate(c, part4, saf)
list_connecting_slices.append(cs)
print('slices_concatenation', list_connecting_slices)
print('Saving ...')
path_save = '/mnt/storage/Thesis_GDMC_XRayTomography/Segmentations/Manip_Novembre'
chdir(path_save)
tifffile.imsave('Concatenation_MON{}_J{}.tif'.format(sample, day), c)
del part1, part2, part3, part4, c
# Concatenating filtered files
chdir('/mnt/storage/Thesis_GDMC_XRayTomography/Before_Segmentation/Manip_Novembre/Aniso')
print('Reading files ...')
part1 = tifffile.imread('MON{}_J{}_scan{}_aniso.tif'.format(sample, day, 1))
part2 = tifffile.imread('MON{}_J{}_scan{}_aniso.tif'.format(sample, day, 2))
cs = list_connecting_slices[0]
print('Concatenation ...')
concatenation_filtered = np.concatenate((part1[:-saf] , part2[cs:]))
del part1, part2
print('Reading files ...')
part3 = tifffile.imread('MON{}_J{}_scan{}_aniso.tif'.format(sample, day, 3))
cs = list_connecting_slices[1]
print('Concatenation ...')
concatenation_filtered = np.concatenate((concatenation_filtered[:-saf] , part3[cs:]))
del part3
print('Reading files ...')
part4 = tifffile.imread('MON{}_J{}_scan{}_aniso.tif'.format(sample, day, 4))
cs = list_connecting_slices[2]
print('Concatenation ...')
concatenation_filtered = np.concatenate((concatenation_filtered[:-saf] , part4[cs:]))
del part4
print('Saving ...')
tifffile.imsave('Concatenation_Aniso_MON{}_J{}.tif'.format(sample, day), concatenation_filtered)
del concatenation_filtered
print('Finished !!')
main()
\ No newline at end of file
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