scientific_comp_projects/CODE/[python]thesis_old_scripts/alone_contour.py

185 lines
6.6 KiB
Python

# -*- coding: utf-8 -*-
"""
Created on Thu Jun 28 15:30:23 2018
@author: Armando
"""
#import os
import ReadIM
import numpy as np
import matplotlib.pyplot as plt
#import matplotlib.animation as animation
#import glob
from matplotlib.colors import LogNorm
#from matplotlib.animation import FFMpegWriter
#from skimage.filters.rank import mean_bilateral
#from skimage.morphology import disk
#import scipy.misc
from skimage import exposure
#Courbe fit values - as well as the maximum and minimum values for the inverse fit function
global a, b, c, d
a = 0.316617
b = 1.26288
c = -7.70744
d = 0.6722
y_min = 0.36
y_max = 0.98
#Numerical to physical equivalences
nappe_laser = 0.250 #Thickness of laser en mm
dimension_pixel = 0.026244 #Equivalence mm/px
#Conversion to SI units
dimension_pixel_m = dimension_pixel * 10**(-3) #Equivalence mm/px
nappe_laser_m = nappe_laser * 10**(-3) #Thickness of laser en m
surface_pixel = dimension_pixel_m**2 #Surface equivalence of a pixel (m^2)
#The working directory, where all the images are stored
working_directory = 'D:\\python_processing\\imagesV2\\eauXG\\'
#help to accelerate script by not showing images
#Plotting options
global interpolation, vmax_value, vmin_value, colormap
colormap = plt.get_cmap('viridis')
vmin_value = 0.001 #Minimum value for the pixel value of a processed image
vmax_value = 0.6 #Maximum value for the pixel value of a processed image
interpolation = 'bicubic' #Output image interpolation (see matplotlib interpolations for more information)
origin = []
#Fit-inverse function
def inverse_tanhfit(y):
return (1/2*np.log((a+y-d)/(a-y+d))-c)/(b)
cut_image_y0 = 0 #The borders that will be cropped from the image [y0,yf,x0,xf]
cut_image_yf = 0 #xf and yf, are the values to be rested from the total images(hence if 0 the max values are used)
cut_image_x0 = 0
cut_image_xf = 0
radius = 0
r=[]
#The experiment and image to process (only one of each will be used since it is a test)
working_dir = 'C:\\Users\\Armando\\Desktop\\conf\\'
image_to_process = 'B00210.im7'
vbuff, vatts = ReadIM.extra.get_Buffer_andAttributeList(working_dir+ image_to_process)
v_array, vbuff = ReadIM.extra.buffer_as_array(vbuff)
original_array = v_array[1][cut_image_y0:len(v_array[0])-cut_image_yf,cut_image_x0:len(v_array[0][0])-cut_image_xf]
original_array = np.array(original_array,dtype = float)
modif_array = np.array(original_array,dtype = float)
original_array = exposure.rescale_intensity(original_array)
fig, ax = plt.subplots()
#
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
plt.rcParams['figure.figsize'] = (10.0, 6.0)
plt.rcParams['font.size'] = 24
plt.rcParams['lines.linewidth'] = 2
plt.rcParams['lines.markersize'] = 24
plt.rcParams['xtick.direction'] = 'out'
plt.rcParams['ytick.direction'] = 'out'
plt.rcParams['legend.fontsize'] = 24
plt.rcParams['axes.titlesize'] = 24
plt.rcParams['axes.labelsize'] = 24
#
vmin = 0.05
vmax = 0.8
v =[0-17.8423, 33.59-17.8423, 44.88-6.72752, 0-6.72752]
im = ax.imshow(original_array, cmap=colormap, origin='upper', extent = v, norm = LogNorm(vmin , vmax))
plt.ylabel('Z (mm)')
plt.axis('scaled')
#ax.set_yticks([0, 44.88,10])
plt.xlabel('X (mm)')
cbar = fig.colorbar(im, ticks=[vmin, vmax])
cbar.ax.set_yticklabels(['Low', 'High']) # horizontal colorbar
plt.tight_layout()
plt.show()
#%%
vbuff, vatts = ReadIM.extra.get_Buffer_andAttributeList(working_dir+ image_to_process)
v_array, vbuff = ReadIM.extra.buffer_as_array(vbuff)
original_array = v_array[0][cut_image_y0:len(v_array[0])-cut_image_yf,cut_image_x0:len(v_array[0][0])-cut_image_xf]
original_array = np.array(original_array,dtype = float)
modif_array = np.array(original_array,dtype = float)
#original_array = exposure.rescale_intensity(original_array)
for y in range(len(modif_array)-1):
for x in range(len(modif_array[0])-1):
if modif_array[y][x] >= y_max:
modif_array[y][x] = 9
elif modif_array[y][x] <= y_min:
modif_array[y][x] = 9
else:
modif_array[y][x] = inverse_tanhfit(modif_array[y][x])
#The CO2 concentration is then obtained from the pH value
#modif_array = 207608319.9386*np.power(modif_array,-9.7399) #Concentration CO2 en g/m^3 (equal to mg/l)(from measures made by Tom LACASSAGNE)
modif_array = 115467070469*np.exp(-4.63*modif_array) #Concentration CO2 en g/m^3 (equal to mg/l)(theoretical values)
modif_array[modif_array<8e-13] = 0
modif_array[modif_array>1] = 0
#The values under 0.1 and over 40 (g/m^3) are considered noise and hence turned to 0
cut_image1 = modif_array[0:len(modif_array),0:len(modif_array[0])//2]
cut_image2 = np.fliplr(modif_array[0:len(modif_array),len(modif_array[0])//2:len(modif_array[0])-1])
#Then we take the mean value for both images.
modif_array2 = np.mean(np.array([cut_image1,cut_image2]), axis=0,dtype=np.float64)
#This section, creates the r-vector (which have the same size as cutted images)
if r == []:
for k in range(len(modif_array2[0])):
radius += dimension_pixel_m
r.append(radius)
#Multiplication of the image [C] by the {r} radius vector
modif_array2 = (modif_array2) @ r
#Final calculation to obtain the mass of CO2 in the image
mass_co2 = (surface_pixel * np.pi * 2 * modif_array2.sum())/nappe_laser_m #The mass is obtained in (g)
#%%This is the second section that will be used to plot the data obtained
modif_array2 = modif_array.copy()
#modif_array = exposure.rescale_intensity(modif_array)
fig, ax = plt.subplots()
#
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
plt.rcParams['figure.figsize'] = (10.0, 6.0)
plt.rcParams['font.size'] = 24
plt.rcParams['lines.linewidth'] = 2
plt.rcParams['lines.markersize'] = 24
plt.rcParams['xtick.direction'] = 'out'
plt.rcParams['ytick.direction'] = 'out'
plt.rcParams['legend.fontsize'] = 24
plt.rcParams['axes.titlesize'] = 24
plt.rcParams['axes.labelsize'] = 24
#
v =[0-5.95617, 11.3637-5.95617, 27.976, 0]
im = ax.imshow(modif_array, cmap=colormap, origin='upper',
norm = LogNorm(vmin = 0.00001, vmax = 1), extent = v)
plt.ylabel('Z (mm)')
plt.axis('scaled')
#ax.set_yticks([0, 44.88,10])
plt.xlabel('X (mm)')
cbar = fig.colorbar(im, ticks=[0.00001, 1] )
cbar.ax.set_yticklabels(['Low', 'High']) # horizontal colorbar
plt.tight_layout()
plt.show()