Commit 823d0384 authored by Michael Rudolf's avatar Michael Rudolf
Browse files

Major changes and performance improvements.

- Changed from Pandas dataframes to numpy arrays or native python types
- Replaced savgol-filter with scipy.signal.savgol_filter (250% faster)
- Improved peak-detection (2000% faster)
- Optimized the mutual linear regression with list comprehensions
- Added a bootstrap linear regression function
- Replaced hardcoded axis scaling by ratio (x1.125 instead of +2000)
- Removed all r'strings' for a more consistent layout, using \\ for TeX
- The timeseries plot is now flexible for experiments with less runs
- The saving functions were adjusted to match the new data types
- Some variables have been renamed
- Cleaned and reformatted code to match pep8 style
parent c869694e
This diff is collapsed.
......@@ -2,14 +2,14 @@
"""
Created on Mon Jul 23 14:32:17 2018
@author: Michael Warsitzka
@author: Michael Warsitzka, Michael Rudolf
"""
# *************************************************************************
# %%===========================IMPORT==========================================
import numpy as np
import pandas as pd
# import pandas as pd
import os
import RST_Func
import RST_Func as rfnc
# %%==========================NAMES============================================
projectname = 'example'
......@@ -36,80 +36,48 @@ Vars = {'A': 0.022619, # area of shear zone (= surface area of lid) (m^2)
'stressnoise': 0.025,
'smoothwindow': 51} # define window for smooth shear stress curve
# %%============================VARIABLES======================================
N = len(os.listdir(path_in))
len_df = np.zeros(N)
TS = pd.DataFrame()
Sigma = np.zeros(N)
Sigma_avg = np.zeros(N)
Peakstress = np.zeros((N, 2))
Dynstress = np.zeros((N, 2))
Reactstress = np.zeros((N, 2))
Dilation = np.zeros(N)
Weak = np.zeros(N)
Weak_p = np.zeros(N)
# %%=====================READ, CONVERT & EVALUATE DATA=========================
for i in range(0, N):
current_set = RST_Func.convert(path_in, os.listdir(path_in)[i], Vars)
Sigma_avg[i] = round(np.mean(current_set['normalstress'])/Vars['prec'])*Vars['prec'] # only for the file name
len_df[i] = len(current_set['time'])
TS = pd.concat([TS, pd.Series((current_set.shearforce).values)],
ignore_index=True, axis=1)
if rst == 1:
RST_var = RST_Func.eval_shearstress(current_set, Vars)
Sigma[i] = RST_var[0]
Peakstress[i, :] = RST_var[1]
Dynstress[i, :] = RST_var[2]
Reactstress[i, :] = RST_var[3]
Dilation[i] = RST_var[4]
Weak[i] = RST_var[5]
Weak_p[i] = RST_var[6]
else:
RST_var = []
# position of longest measurement
pos_dfmax = int(np.where(len_df == len_df.max())[0][0])
time_max = RST_Func.convert(path_in, os.listdir(path_in)[pos_dfmax], Vars)['time']
# Generate a list containing all files
file_list = [f for f in os.listdir(path_in)]
# data frame of sorted time series for pdf plot:
# sort according to max values
TS_sort = TS.iloc[:, TS.max().sort_values(ascending=True).index]
# sort Sigma avg according to columns
Sigma_avg_sort = Sigma_avg[TS_sort.columns]
# Data is stored as a list of dictionaries containing the data
exp_data = [rfnc.convert(path_in, f, Vars) for f in file_list]
TS_sort = pd.concat([time_max, TS_sort], axis=1)
# combine ts data and new header
TS_sort = pd.DataFrame(TS_sort.values,
columns=['Time [s]'] + [str(i) for i in Sigma_avg_sort])
# Evaluate all data sets and store results in a list
eval_data = [rfnc.eval_shearstress(cur_dat, Vars) for cur_dat in exp_data]
# data frame of time series for export to txt file:
TS = pd.concat([time_max, TS], axis=1) # add time column to data frame
# combine ts data and new header
TS = pd.DataFrame(TS.values,
columns=['Time [s]'] + [str(i) for i in Sigma_avg])
# Create lists of parameters for analysis
normal_stress = [ev[0] for ev in eval_data]
peak_stress = [ev[1][1] for ev in eval_data]
dyn_stress = [ev[2][1] for ev in eval_data]
stat_stress = [ev[3][1] for ev in eval_data]
dilation = [ev[4] for ev in eval_data]
weak = [ev[5] for ev in eval_data]
weak_p = [ev[6] for ev in eval_data]
# %%==================ANALYSIS OF FRICTION COEFFICIENTS AND COHESIONS==========
if rst == 1:
peakstress = [Peakstress[i][1] for i in range(0, N)]
dynstress = [Dynstress[i][1] for i in range(0, N)]
reactstress = [Reactstress[i][1] for i in range(0, N)]
# peak_stress = [peak_stress[i][1] for i in range(0, N)]
# dyn_stress = [dyn_stress[i][1] for i in range(0, N)]
# stat_stress = [stat_stress[i][1] for i in range(0, N)]
# %%===========Analysis 1: Mutual linear regression========================
peakfric_mut, peakdata_mut = RST_Func.rst_analmut(Sigma_avg, peakstress)
dynfric_mut, dyndata_mut = RST_Func.rst_analmut(Sigma_avg, dynstress)
reactfric_mut, reactdata_mut = RST_Func.rst_analmut(Sigma_avg, reactstress)
peakfric_mut, peakdata_mut = rfnc.rst_analmut(normal_stress, peak_stress)
dynfric_mut, dyndata_mut = rfnc.rst_analmut(normal_stress, dyn_stress)
reactfric_mut, reactdata_mut = rfnc.rst_analmut(normal_stress, stat_stress)
# %%===========ANALYSIS 2: Standard regression of all data pairs ==========
peakfric_std = RST_Func.rst_analstd(Sigma_avg, peakstress)
dynfric_std = RST_Func.rst_analstd(Sigma_avg, dynstress)
reactfric_std = RST_Func.rst_analstd(Sigma_avg, reactstress)
peakfric_std = rfnc.rst_analstd(normal_stress, peak_stress)
dynfric_std = rfnc.rst_analstd(normal_stress, dyn_stress)
reactfric_std = rfnc.rst_analstd(normal_stress, stat_stress)
# %%==========MERGE DATA===================================================
Strength = (Sigma_avg, peakstress, dynstress, reactstress, Weak, Weak_p)
Fric_mut = (peakfric_mut, dynfric_mut, reactfric_mut)
Data_mut = (peakdata_mut, dyndata_mut, reactdata_mut)
Fric_std = (peakfric_std, dynfric_std, reactfric_std)
strength = (normal_stress, peak_stress, dyn_stress,
stat_stress, weak, weak_p)
fric_mut = (peakfric_mut, dynfric_mut, reactfric_mut)
data_mut = (peakdata_mut, dyndata_mut, reactdata_mut)
fric_std = (peakfric_std, dynfric_std, reactfric_std)
print('==================================================================')
print('RST analysis successful')
print('==================================================================')
......@@ -120,27 +88,27 @@ else:
# %%====================PLOT DATA==============================================
if rst == 1:
RST_Func.plotstd(path_out, projectname, Strength, Fric_std)
RST_Func.plothist(path_out, projectname, Strength, Data_mut)
rfnc.plotstd(path_out, projectname, strength, fric_std)
rfnc.plothist(path_out, projectname, strength, data_mut)
print('==================================================================')
print('>>> Fiction data plotted')
print('==================================================================')
if plot_ts == 1:
RST_Func.plotts(path_out, projectname, TS_sort, Sigma_avg_sort, Vars)
rfnc.plotts(path_out, projectname, exp_data, normal_stress)
print('==================================================================')
print('>>> Time series data plotted')
print('==================================================================')
# %%====================Save DATA==============================================
if rst == 1:
RST_Func.saveStrength(path_out, projectname, Strength)
RST_Func.saveFric(path_out, projectname, Fric_mut, Fric_std)
# RST_Funct.RST_plotstd(path_out, name_rst2, RST2)
rfnc.saveStrength(path_out, projectname, strength)
rfnc.saveFric(path_out, projectname, fric_mut, fric_std)
# rfnct.RST_plotstd(path_out, name_rst2, RST2)
print('==================================================================')
print('>>> Friction data saved')
print('==================================================================')
if save_ts == 1:
RST_Func.saveTS(path_out, projectname, TS)
rfnc.saveTS(path_out, projectname, exp_data)
print('==================================================================')
print('>>> Time series data saved')
print('==================================================================')
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