""" script reads preconditioing data (single or multiple files), picks all peaks, finds % relaxation after 10s and after 1000 cycles, saves figure and saves outputs to csv file
# Input needed : at least one pc data file, area of cross section"""
import sys
import matplotlib.pyplot as plt
import peakutils
import numpy as np
import csv
import statistics
import BasicUtilities
def USAGE(argv):
print
print ('USAGE: ') + argv[0] + ' <at least one data file required>'
print
def main(argv):
if len(argv)<2:
USAGE(argv)
sys.exit(1)
num = len(sys.argv[1:]) # count number of arguments after the python script name, allows any number of preconditioing datafiles
percent_relax_overall = []
percent_relax_early = []
output_filename = argv[1][0:20]
input_tags = ""
for j in range (1,num):
input_tags = input_tags + argv[j][22] # grabbing file tags
time = [[] for x in range (num)]
load = [[] for x in range (num)]
disp = [[] for x in range (num)]
for i in range(1,num):
# extract data from each file
f_name = argv[i]
infile = open(f_name,'r')
for line in infile:
line = line.rstrip('\n')
line = line.split()
line = map(float,line)
time[i-1].append(line[0])
disp[i-1].append(line[2])
load[i-1].append(line[1])
infile.close()
load[i-1] = [0.0098 * r for r in load[i-1]]
# plot
plt.figure(1)
plt.plot(time[i-1],load[i-1],label=('test'+str(i)))
plt.legend(loc=1)
plt.legend(loc=1)
plt.xlabel('time, s')
plt.ylabel('load,N')
""" grab peaks, filter out false peaks, plot true peaks, peak at 10s and first and last peaks. """
# find last peak, may not be real peak
k = -1
while load[i-1][k] < load[i-1][k-1] :
k = k -1
peak_end = load[i-1][k]
# find first peak, may not be real peak
m = 0
while load[i-1][m] < load[i-1][m+1] :
m = m+1
peak_one = load[i-1][m]
threshold = peak_end / peak_one # possibly lower than optimum threshold, may get several false positives
print "threshold =", threshold
plt.figure(2)
plt.plot(time[i-1][k],load[i-1][k],'bx')
plt.plot(time[i-1][m],load[i-1][m],'bx')
# find all peaks
time_p = []
load_p = []
load[i-1] = np.array(load[i-1])
time[i-1] = np.array(time[i-1])
# find indices of all peaks, default min dist = 0.45, thres = threshold, FMC-L,cc = thres = 0.2
indexes = peakutils.peak.indexes(load[i-1], thres=threshold, min_dist=0.45)
time[i-1] = np.array(time[i-1]) # since list can't be indexed using another list, convert to array
load[i-1] = np.array(load[i-1])
time_p = time[i-1][indexes]
load_p = load[i-1][indexes]
plt.plot(time_p,load_p, 'r*')
# remove false peaks
print 'Removing false peaks'
time_p = list(time_p)
load_p = list(load_p)
time_p_actual = []
load_p_actual = []
index_peak = load_p.index(max(load_p)) # index of peak
del(load_p[0:index_peak])
del(time_p[0:index_peak])
del_index = []
for j in range (0, len(time_p)-1):
if time_p[j+1]-time_p[j] >= 0.4: # time between peaks> 0.5s (.4 works for most, FMC-M,unc = 0.23, FMC-M,cc = 0.23, FMC-L,cc = 0.16)
time_p_actual.append(time_p[j])
load_p_actual.append(load_p[j])
elif time_p[j]-time_p[j-1] >= 0.4: # (.4 works for most (FMC-M, cc = 0.35, FMC-L,cc = 0.16 )
time_p_actual.append(time_p[j])
load_p_actual.append(load_p[j])
del_index.append(j+1)
time_p_actual = [d for e, d in enumerate(time_p_actual) if e not in del_index]
load_p_actual = [d for e, d in enumerate(load_p_actual) if e not in del_index]
# time_p_actual.append(time_p[-1]) # grab last peak
# load_p_actual.append(load_p[-1])
plt.figure(3)
plt.plot(time_p_actual,load_p_actual,'.', label=('test'+str(i)))
percent_relax_overall.append((((load_p_actual[0] - load_p_actual[-1]))/load_p_actual[0]) * 100 )
print '% relaxation =', percent_relax_overall
# % change in the first 10s
def find_nearest(array, value):
array = np.asarray(array)
idx = (np.abs(array - value)).argmin()
return idx
idx = find_nearest(time_p_actual, 10)
time_init = time_p_actual[idx]
load_init = load_p_actual[idx]
plt.plot(time_init, load_init, "y*", ms= 15)
plt.legend(loc=1)
plt.legend(loc=1)
plt.xlabel('time, s')
plt.ylabel('load,N')
plt.figure(4)
stress_p_actual = []
stress_p_actual = [ (c / float(argv[-1])) for c in load_p_actual]
stress_init = stress_p_actual[idx]
plt.plot(time_p_actual, stress_p_actual,'.' , label=('test'+str(i)))
plt.plot(time_init, stress_init, "y*", ms= 15)
plt.legend(loc=1)
plt.legend(loc=1)
plt.xlabel('time, s')
plt.ylabel('stress,MPa')
percent_relax_early.append((((load_p_actual[0] - load_p_actual[idx]))/load_p_actual[0]) * 100)
print '% relaxation first 10s =', percent_relax_early
# save fig
pc_Plot = output_filename + '_' + input_tags + '_pc_repeat'+ '.svg'
plt.savefig(pc_Plot)
plt.show()
# calculate mean, SD and COV
mean_overall_relaxation = statistics.mean(percent_relax_overall)
SD_overall_relaxation = statistics.stdev(percent_relax_overall)
cov_overall_relaxation = (SD_overall_relaxation / mean_overall_relaxation) * 100
mean_10s_relaxation = statistics.mean(percent_relax_early)
SD_10s_relaxation = statistics.stdev(percent_relax_early)
cov_10s_relaxation = (SD_10s_relaxation / mean_10s_relaxation) * 100
# write to output file created by stress relaxation data analysis file , save relaxation %
output_file = output_filename + '_' + input_tags + '_output_pc_repeat' # repeatability output file
out_datafile_id = open(output_file + '.csv', 'w')
out_writer = csv.writer(out_datafile_id)
out_writer.writerow(['overall_%_relaxation',percent_relax_overall])
out_writer.writerow(['%_relaxation_first_10s',percent_relax_early])
out_writer.writerow(['mean_overall_relaxation_for_test_set',mean_overall_relaxation])
out_writer.writerow(['SD_overall_relaxation_for_test_set',SD_overall_relaxation])
out_writer.writerow(['COV_overall_relaxation_for_test_set',cov_overall_relaxation])
out_writer.writerow(['mean_10s_relaxation_for_test_set',mean_10s_relaxation])
out_writer.writerow(['SD_10s_relaxation_for_test_set',SD_10s_relaxation])
out_writer.writerow(['COV_10s_relaxation_for_test_set',cov_10s_relaxation])
# plot % change vs time
# plt.figure(4)
# a = time_p_actual[-1]
# plt.plot([10,10,10], percent_relax_early, 'ro')
# plt.plot([a,a,a],percent_relax_overall,'bo')
if __name__=="__main__":
main(sys.argv)