To receive notifications about scheduled maintenance, please subscribe to the mailing-list gitlab-operations@sympa.ethz.ch. You can subscribe to the mailing-list at https://sympa.ethz.ch

Commit 985fec4c authored by felikskiszkurno's avatar felikskiszkurno
Browse files

Date: 15.01.2021

Initial commit.
Begin of work on tests with scikit classificators.
parents
# SlopeStability
** Author: Feliks Kiszkurno **
## Description:
Scripts and tools developed for the purpose of my Master Thesis
## Structure
### Main script
This is where most of the work happens.
### Modules
Functions and classes, that shouldn't be put in the main files. Assigned to different modules based on their purpose.
#### slopestabilitytools
All tools, that are not directly related to modeling, inversion, or classification and evaluation of the data.
\ No newline at end of file
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Fri Jan 8 10:29:00 2021
@author: felikskrno
"""
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import slopestabilitytools
import slopestabilityML
import slostabcreatedata
import pygimli as pg
import pygimli.meshtools as mt
import pygimli.physics.ert as ert
is_success = slopestabilitytools.folder_structure.create_folder_structure()
number_of_tests = 10
rho_spread_factor = 1.5
rho_max = 150
layers_min = 1
layers_max = 3
min_depth = 1
max_depth = 10
tests_horizontal = slopestabilitytools.model_params(number_of_tests,
rho_spread_factor, rho_max,
layers_min, layers_max,
min_depth, max_depth)
for test_name in tests_horizontal['names'].values():
slostabcreatedata.create_data()
world_boundary_v = [-200, 0] # [right, left border] relatively to the middle
world_boundary_h = [200, -100] # [top, bottom border]
test_results = {}
test_results_grid = {}
test_input = {}
test_name = 'hor_1'
#tests_horizontal['layers_pos'][test_name] = [-5]
# INPUT MODEL - SUBSURFACE START #
world = mt.createWorld(start=world_boundary_v, end=world_boundary_h, layers=tests_horizontal['layers_pos'][test_name])#,
#marker=np.linspace(1, tests_horizontal['layer_n']['hor_1'],
# tests_horizontal['layer_n']['hor_1']))
geometry = world # +block
measurement_scheme = ert.createERTData(elecs=np.linspace(start=-45, stop=45, num=91), schemeName='dd')
for electrode in measurement_scheme.sensors():
geometry.createNode(electrode)
geometry.createNode(electrode - [0, 0.1]) # What does it do?
mesh = mt.createMesh(geometry, quality=34) # , area=2)#
resistivity_map = tests_horizontal['rho_values'][test_name] # [0]
#resistivity_map[0] = [1, 50.0]
#resistivity_map[1] = [2, 150.0]
input_model = pg.solver.parseMapToCellArray(resistivity_map, mesh) # rename to input_mesh
# INPUT MODEL - SUBSURFACE MODEL END ###
# SIMULATE ERT MEASUREMENT - START ###
mesh_pd = [] # add new mesh
data = ert.simulate(mesh, scheme=measurement_scheme, res=resistivity_map, noiseLevel=1, noiseAbs=1e-6, seed=1337)
data.remove(data['rhoa'] < 0)
# SIMULATE ERT MEASUREMENT - END ###
ert_manager = ert.ERTManager(sr=False, useBert=True, verbose=True, debug=False)
# RUN INVERSION #
k0 = pg.physics.ert.createGeometricFactors(data)
model_inverted = ert_manager.invert(data=data, lam=20, paraDX=0.25, paraMaxCellSize=5, paraDepth=10, quality=34,
zPower=0.4)
result = ert_manager.inv.model
result_array = result.array()
input_model2 = pg.interpolate(srcMesh=mesh, inVec=input_model, destPos=ert_manager.paraDomain.cellCenters())
input_model2_array = input_model2.array()
experiment_results = pd.DataFrame(data={'X': ert_manager.paraDomain.cellCenters().array()[:, 0],
'Y': ert_manager.paraDomain.cellCenters().array()[:, 1],
'Z': ert_manager.paraDomain.cellCenters().array()[:, 2],
'INM': input_model2_array,
'RES': result_array})
test_results[test_name] = experiment_results
fig, ax = plt.subplots(3)
fig.suptitle(test_name)
pg.show(ert_manager.paraDomain, input_model2, ax=ax[0])
ax[0].set_title("Model on the output mesh")
pg.show(ert_manager.paraDomain, result, ax=ax[1])
ax[1].set_title("Inverted")
pg.show(ert_manager.paraDomain, result - input_model2, ax=ax[2])
ax[2].set_title("Diff")
fig.savefig('results/figs/hor_{}_results.eps'.format(test_name))
fig.savefig('results/figs/hor_{}_results.png'.format(test_name))
fig_input, ax_input = plt.subplots(1)
pg.show(mesh, input_model, ax=ax_input)
ax_input.set_title('1 Geometry of the model')
fig_input.savefig('results/figs/hor_{}_input.eps'.format(test_name))
fig_input.savefig('results/figs/hor_{}_input.png'.format(test_name))
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Fri Jan 8 10:29:00 2021
@author: felikskrno
"""
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import slopestabilitytools
import slopestabilityML
import slostabcreatedata
import pygimli as pg
import pygimli.meshtools as mt
import pygimli.physics.ert as ert
is_success = slopestabilitytools.folder_structure.create_folder_structure()
number_of_tests = 10
rho_spread_factor = 1.5
rho_max = 150
layers_min = 1
layers_max = 3
min_depth = 1
max_depth = 10
tests_horizontal = slopestabilitytools.model_params(number_of_tests,
rho_spread_factor, rho_max,
layers_min, layers_max,
min_depth, max_depth)
world_boundary_v = [-200, 0] # [right, left border] relatively to the middle
world_boundary_h = [200, -100] # [top, bottom border]
test_results = {}
test_results_grid = {}
test_input = {}
test_name = 'hor_1'
#tests_horizontal['layers_pos'][test_name] = [-5]
# INPUT MODEL - SUBSURFACE START #
world = mt.createWorld(start=world_boundary_v, end=world_boundary_h, layers=tests_horizontal['layers_pos'][test_name])#,
#marker=np.linspace(1, tests_horizontal['layer_n']['hor_1'],
# tests_horizontal['layer_n']['hor_1']))
geometry = world # +block
measurement_scheme = ert.createERTData(elecs=np.linspace(start=-45, stop=45, num=91), schemeName='dd')
for electrode in measurement_scheme.sensors():
geometry.createNode(electrode)
geometry.createNode(electrode - [0, 0.1]) # What does it do?
mesh = mt.createMesh(geometry, quality=34) # , area=2)#
resistivity_map = tests_horizontal['rho_values'][test_name] # [0]
#resistivity_map[0] = [1, 50.0]
#resistivity_map[1] = [2, 150.0]
input_model = pg.solver.parseMapToCellArray(resistivity_map, mesh) # rename to input_mesh
# INPUT MODEL - SUBSURFACE MODEL END ###
# SIMULATE ERT MEASUREMENT - START ###
mesh_pd = [] # add new mesh
data = ert.simulate(mesh, scheme=measurement_scheme, res=resistivity_map, noiseLevel=1, noiseAbs=1e-6, seed=1337)
data.remove(data['rhoa'] < 0)
# SIMULATE ERT MEASUREMENT - END ###
ert_manager = ert.ERTManager(sr=False, useBert=True, verbose=True, debug=False)
# RUN INVERSION #
k0 = pg.physics.ert.createGeometricFactors(data)
model_inverted = ert_manager.invert(data=data, lam=20, paraDX=0.25, paraMaxCellSize=5, paraDepth=10, quality=34,
zPower=0.4)
result = ert_manager.inv.model
result_array = result.array()
input_model2 = pg.interpolate(srcMesh=mesh, inVec=input_model, destPos=ert_manager.paraDomain.cellCenters())
input_model2_array = input_model2.array()
experiment_results = pd.DataFrame(data={'X': ert_manager.paraDomain.cellCenters().array()[:, 0],
'Y': ert_manager.paraDomain.cellCenters().array()[:, 1],
'Z': ert_manager.paraDomain.cellCenters().array()[:, 2],
'INM': input_model2_array,
'RES': result_array})
test_results[test_name] = experiment_results
fig, ax = plt.subplots(3)
fig.suptitle(test_name)
pg.show(ert_manager.paraDomain, input_model2, ax=ax[0])
ax[0].set_title("Model on the output mesh")
pg.show(ert_manager.paraDomain, result, ax=ax[1])
ax[1].set_title("Inverted")
pg.show(ert_manager.paraDomain, result - input_model2, ax=ax[2])
ax[2].set_title("Diff")
fig.savefig('results/figs/hor_{}_results.eps'.format(test_name))
fig.savefig('results/figs/hor_{}_results.png'.format(test_name))
fig_input, ax_input = plt.subplots(1)
pg.show(mesh, input_model, ax=ax_input)
ax_input.set_title('1 Geometry of the model')
fig_input.savefig('results/figs/hor_{}_input.eps'.format(test_name))
fig_input.savefig('results/figs/hor_{}_input.png'.format(test_name))
This diff is collapsed.
numpy~=1.18.5
matplotlib~=3.3.2
pandas~=1.1.5
pygimli~=1.1.0
\ No newline at end of file
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
""" This module is created to contain all tools, that are not directly related to inversion or processing
List:
model_params: generates parameters used to develop models
directory_structure: create directory structure to contain figures and other files
"""
from .model_params import model_params
from .plot_and_save import plot_and_save
from .folder_structure.create_folder_structure import create_folder_structure
from .folder_structure.create_folder_for_test import create_folder_for_test
from .folder_structure.check_create_folder import check_create_folder
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
from .create_folder_for_test import create_folder_for_test
from .create_folder_structure import create_folder_structure
from .check_create_folder import check_create_folder
\ No newline at end of file
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
import os
def check_create_folder(folder_path):
worked = True
if os.path.isdir(folder_path):
print("Folder for figures exists!")
else:
print(folder_path)
os.mkdir(folder_path)
print("Created folder for figures!")
return worked
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
import os
from .check_create_folder import check_create_folder
def create_folder_for_test(test_name):
folder_path = os.getcwd() + '/results/results/%'.format(test_name)
is_success = check_create_folder(folder_path)
return is_success
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
import os
from .check_create_folder import check_create_folder
def create_folder_structure():
is_success = True
# Folder for figures
folder_path = os.getcwd()+'/results/'
is_success = check_create_folder(folder_path)
folder_path = os.getcwd() + '/results/figures/'
is_success = check_create_folder(folder_path)
# Folder for results
folder_path = os.getcwd()+'/results/'
is_success = check_create_folder(folder_path)
folder_path = os.getcwd() + '/results/results/'
is_success = check_create_folder(folder_path)
return is_success
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Fri Jan 8 10:29:00 2021
@author: felikskrno
"""
import numpy as np
from numpy import random
random.seed(999)
# number_of_tests = 1
# rho_spread_factor = 1.5
# max_rho = 150
# layers_min = 1
# layers_max = 3
# world_boundary_v = [-100, 0] # [right, left border] relatively to the middle
# world_boundary_h = [200, -200] # [top, bottom border]
def model_params(n_of_tests, rho_spread, rho_max, layers_n_min, layers_n_max, depth_min, depth_max):
test_names = {}
test_n_layers = {}
test_rho = {}
test_layers_pos = {}
for test_id in range(n_of_tests):
test_names[test_id] = 'hor_{}'.format(str(test_id + 1))
test_n_layers[test_names[test_id]] = np.random.randint(layers_n_min, layers_n_max)
rho_temp = []
rho_used = []
layer_pos_temp = []
layer_pos_used = []
layer = 0
while layer < test_n_layers[test_names[test_id]]:
new_layer_pos = random.randint(depth_min, depth_max)
if len(layer_pos_temp) == 0:
layer_pos_temp.append(new_layer_pos)
layer_pos_used.append(new_layer_pos)
else:
if new_layer_pos not in layer_pos_used:
new_layer_pos = random.randint(depth_min, depth_max)
layer_pos_temp.append(new_layer_pos)
layer_pos_used.append(new_layer_pos)
else:
while new_layer_pos in layer_pos_used:
new_layer_pos = random.randint(depth_min, depth_max)
layer_pos_temp.append(new_layer_pos)
layer_pos_used.append(new_layer_pos)
layer += 1
layer = 0
while layer < test_n_layers[test_names[test_id]] + 1:
new_rho = int(random.rand(1)[0] * rho_max)
if len(rho_temp) == 0:
#rho_temp.append([layer + 1, new_rho])
rho_used.append(new_rho)
else:
if new_rho not in rho_used:
new_rho = int(random.rand(1)[0] * rho_max)
#rho_temp.append([layer + 1, new_rho])
rho_used.append(new_rho)
else:
while new_rho in rho_used:
new_rho = int(random.rand(1)[0] * rho_max)
#rho_temp.append([layer + 1, new_rho])
rho_used.append(new_rho)
layer += 1
rho_temp = np.sort(rho_used)
rho_final = []
layer_id = 1
for rho in rho_temp:
print(rho)
rho_final.append([layer_id, rho])
layer_id += 1
test_layers_pos[test_names[test_id]] = -1*(np.sort(np.array(layer_pos_temp)))#np.flip(np.sort(np.array(layer_pos_temp)))
test_rho[test_names[test_id]] = rho_final#np.sort(rho_temp)
rho_temp = []
rho_used = []
layer_pos = []
layer_pos_temp = []
layer_pos_used = []
tests_horizontal = {'names': test_names,
'layer_n': test_n_layers,
'rho_values': test_rho,
'layers_pos': test_layers_pos}
tests_horizontal2 = {}
for test_name in test_names:
return tests_horizontal
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on 15.01.2021
@author: Feliks Kiszkurno
"""
def plot_and_save():
return
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on
@author:
"""
def create_data(tests_horizontal):
world_boundary_v = [-200, 0] # [right, left border] relatively to the middle
world_boundary_h = [200, -100] # [top, bottom border]
test_results = {}
test_results_grid = {}
test_input = {}
test_name = 'hor_1'
for test_name in tests_horizontal['']
# tests_horizontal['layers_pos'][test_name] = [-5]
# INPUT MODEL - SUBSURFACE START #
world = mt.createWorld(start=world_boundary_v, end=world_boundary_h,
layers=tests_horizontal['layers_pos'][test_name]) # ,
# marker=np.linspace(1, tests_horizontal['layer_n']['hor_1'],
# tests_horizontal['layer_n']['hor_1']))
geometry = world # +block
measurement_scheme = ert.createERTData(elecs=np.linspace(start=-45, stop=45, num=91), schemeName='dd')
for electrode in measurement_scheme.sensors():
geometry.createNode(electrode)
geometry.createNode(electrode - [0, 0.1]) # What does it do?
mesh = mt.createMesh(geometry, quality=34) # , area=2)#
resistivity_map = tests_horizontal['rho_values'][test_name] # [0]
# resistivity_map[0] = [1, 50.0]
# resistivity_map[1] = [2, 150.0]
input_model = pg.solver.parseMapToCellArray(resistivity_map, mesh) # rename to input_mesh
# INPUT MODEL - SUBSURFACE MODEL END ###
# SIMULATE ERT MEASUREMENT - START ###
mesh_pd = [] # add new mesh
data = ert.simulate(mesh, scheme=measurement_scheme, res=resistivity_map, noiseLevel=1, noiseAbs=1e-6, seed=1337)
data.remove(data['rhoa'] < 0)
# SIMULATE ERT MEASUREMENT - END ###
ert_manager = ert.ERTManager(sr=False, useBert=True, verbose=True, debug=False)