diff --git a/Python_script/PostPlot.py b/Python_script/PostPlot.py index 115dbffbe601fb9d5fda88b5564cb565201ec72b..4a64f8547a5a85e5a55d52d50a425893e0ecc410 100644 --- a/Python_script/PostPlot.py +++ b/Python_script/PostPlot.py @@ -196,7 +196,7 @@ if __name__ == '__main__': # plot regression coefficients plot_obj.edit_annotation_in_plot(annotate_string = reg_coeff) - plot_obj.edit_annotation_in_plot(annotate_string = reg_coeff) + plot_obj.edit_annotation_in_plot(annotate_string = reg_coeff) if plot_trace_curvefit == True: diff --git a/Python_script/curvefit_and_correlation.py b/Python_script/curvefit_and_correlation.py new file mode 100644 index 0000000000000000000000000000000000000000..9952f9221418c6a2c6928b8dade1b965165f0ccc --- /dev/null +++ b/Python_script/curvefit_and_correlation.py @@ -0,0 +1,254 @@ +# -*- coding: utf-8 -*- +""" +Created on Thu Jul 27 08:53:36 2023 + +@author: pawelzik +""" +import numpy as np +from scipy.optimize import curve_fit + + + +class reg_and_corr: + + def __init__(self, post_plot_obj, legend_loc = 'upper left', \ + legend_bbox_to_anchor = (1.09, 1)): + + self.post_plot_obj = post_plot_obj + self.K_phases = None + self.phase_t0 = None + self.legend_loc = legend_loc + self.legend_bbox_to_anchor = legend_bbox_to_anchor + + + def calc_correlation(self, func_1 , func_2 , state = 'False'): + # calulate correlation between two parameter from data frame + corr_coeff= np.corrcoef(func_1, func_2)[0][1] + + annotate_string = "$R_{\phi(S_{21}),Temp_{DUT}}$ = %.3f" % corr_coeff + + return annotate_string + + def calc_regression_coeff_phase_S21(self, data_frame, time_unit = 'min', state = False): + + phases = data_frame.S21_PHASE.tolist() + self.phase_t0 = phases[0] + + self.K_phases = np.median(data_frame.S21_PHASE) + + # time stamp of index = 0 is the start point of time axis, Michael + # substract all other timestamps with time stamp of index zero to get time scale in + # seconds, Michael + + time_sec = data_frame.TIMESTAMP - data_frame.TIMESTAMP[0] + + # set scaling of time axis depending on the time unit given in parameter list, Michael + # default unit is minutes, Michael + time_vals = self.scaling_time_axes(time_sec, time_unit) + + func, popt = self.choose_fit(time_vals, data_frame.S21_PHASE) + + annotate_string_fit = self.plot_fitted_func_param(func, *popt, time_unit = time_unit) + + return annotate_string_fit + + def plot_phase_curve_fit(self, data_frame, postplot_obj, time_unit = 'min', state = False): + + if state == True: + + + phases = data_frame.S21_PHASE.tolist() + self.phase_t0 = phases[0] + + self.K_phases = np.median(data_frame.S21_PHASE) + + # time stamp of index = 0 is the start point of time axis, Michael + # substract all other timestamps with time stamp of index zero to get time scale in + # seconds, Michael + + time_sec = data_frame.TIMESTAMP - data_frame.TIMESTAMP[0] + + # set scaling of time axis depending on the time unit given in parameter list, Michael + # default unit is minutes, Michael + time_vals = self.scaling_time_axes(time_sec, time_unit) + + # call chose fit function + func, popt = self.choose_fit(time_vals, data_frame.S21_PHASE) + + # genrate trace for curvefit + fit_phaeses = func(time_vals, *popt) + + # determine max, min of measured phases + min_phases, max_phases = self.get_extended_min_max(data_frame.S21_PHASE) + + # determine max, min of curvefit + min_fit_phases, max_fit_phases = self.get_extended_min_max(fit_phaeses) + + # determine min of measure phases and min of curvefit + minimum = min(min_phases, min_fit_phases) + + # determine max of measured phases and max of curvefit + maximum = max(max_phases, max_fit_phases) + + # set axis for phaseplot to min, max values + postplot_obj.ax1[0].set_ylim(minimum, maximum) + + # set color of trace for curvefit result to green (trace is visible) + postplot_obj.path_collection_fit.set_color('green') + + # set label of trace curvefit results to name of fit function name + postplot_obj.path_collection_fit.set_label('fitted with\n' + func.__name__) + + + + postplot_obj.path_collection_fit.set_offsets(np.c_[time_vals, func(time_vals, *popt)]) + + else: + # set color of trace for curvefit result to None (trace is invisible) + postplot_obj.path_collection_fit.set_color('None') + # set label of trace curvefit results to '' (no label) + postplot_obj.path_collection_fit.set_label('') + + # get legend handles and labels y-axes from subplot magnitude, phase + handles_phase, labels_phase = postplot_obj.ax1[0].get_legend_handles_labels() + handles_mag, labels_mag = postplot_obj.magnitude_axis.get_legend_handles_labels() + handles_equi0, labels_eqi0 = postplot_obj.equi_axis0.get_legend_handles_labels() + + handles = handles_phase + handles_mag + handles_equi0 + labels = labels_phase + labels_mag + labels_eqi0 + + # update legend subplot phase, magnitude + postplot_obj.ax1[0].legend(handles, labels, loc = postplot_obj.legend_loc, \ + bbox_to_anchor = postplot_obj.legend_bbox_to_anchor) + + # refresh plot window + postplot_obj.fig.canvas.draw() + postplot_obj.fig.canvas.flush_events() + + + + +# step response of PT1 with initial behavior , Michael + def PT1(self, x, K, T1): + return K * (1 - np.exp(-x/T1)) + self.phase_t0 * np.exp(-x/T1) + + # step response of PT2 (D > 1) with initial behavior + def aperiodicPT2(self, x, K, T1, T2,): + return K + self.phase_t0* np.exp(-x/T1)/(T2-T1) - K *T1* np.exp(-x/T1)/(T2-T1) \ + - self.phase_t0 * np.exp(-x/T2)/(T2-T1) + K* T2* np.exp(-x/T2)/(T2-T1) + + # step response of PT2 (D = 1) with initial behavior, Michael + def aper_borderPT2(self, x, K, T): + return K + self.phase_t0*np.exp(-x/T) - K * np.exp(-x/T) + self.phase_t0*np.exp(-x/T)*x/T + \ + - K * np.exp(-x/T) * x/T + + # step response of PT2 (0 < D < 1) with initial behavior, Michael + def periodicPT2(self, x, D, T): + return self.K_phases + (self.phase_t0 -self.K_phases)* \ + np.exp(-D*x/T)*np.cos(np.sqrt(1-np.square(D))*x/T) + \ + (self.phase_t0 - self.K_phases) * D * \ + np.exp(-D*x/T)*np.sin(np.sqrt(1-np.square(D))*x/T)/np.sqrt(1-np.square(D)) + + # step response series of PT2 ( D = 1) and PT1 with initial behavior, Michael + def aper_borderPT2_PT1(self, x, a, b, c, d, T1, T2): + return a + b*np.exp(-x/T1)+c*x*np.exp(-x/T1) + d*np.exp(-x/T2) + + # method for cuvefit, Michael + # two different algorithms for curvefit are use depending on used step response for fit, Michael + def curvefit(self, func, time, phases): + + if func.__name__ == 'periodicPT2': + popt, pcov = curve_fit(func, time, phases, method = 'trf', \ + bounds = ([0.001, 0.001], [0.99, np.inf])) + else: + popt, pcov = curve_fit(func, time, phases, method = 'lm') + + return popt, pcov + + # method adds the parameter of the fit with the lowest error to annotation string, Michael + + def plot_fitted_func_param(self, func, *popt, time_unit): + + if func.__name__ == 'PT1': + K = popt[0] + T1 = popt[1] + annotate_string = "K: %.2f\n$T_1$: %.3f %s" %(K, T1, time_unit) + elif func.__name__ == 'aperiodicPT2': + K = popt[0] + T1 = popt[1] + T2 = popt[2] + annotate_string = "K: %.2f\n$T_1$: %.3f %s\n$T_2$: %.3f %s" % (K, T1, time_unit, T2, \ + time_unit) + elif func.__name__ == 'aper_borderPT2': + K = popt[0] + T = popt[1] + annotate_string = "K: %.2f\nT: %.3f %s" %(K, T, time_unit) + elif func.__name__ == 'periodicPT2': + K = self.K_phases + D = popt[0] + T = popt[1] + T_2perc = 4*T/D + annotate_string = "K: %.2f\nD: %.3f\nT: %.3f %s\n$T_E$(2%%): %.2f %s" \ + % (K, D, T, time_unit, T_2perc, time_unit) + elif func.__name__ == 'aper_borderPT2_PT1': + T1 = popt[4] + T2 = popt[5] + annotate_string = "$T_1$: %.3f %s\n$T_2$: %.3f %s" % (T1, time_unit, T2, time_unit) + else: + # default value + annotate_string ='' + return annotate_string + + # calculates the fit for all step responses and calulates the lowest error, Michael + # parameter of the plot with the lowest error is returned by this function, Michael + # if curve fit fails for one stept response it takes the next one and plots in command + # window the functions + # which failes, Michael + def choose_fit(self, time, phases): + popt = [] + perr = [] + used_functions = [] + functions = [self.PT1, self.aper_borderPT2, self.aperiodicPT2, self.periodicPT2, \ + self.aper_borderPT2_PT1] + for func in functions: + try: + # call curvefit method and put the fitted parameter for step response function + # to a list, Michael + popt.append(self.curvefit(func, time, phases)[0]) + # calculation of fitting error + perr.append(np.sqrt(np.square(np.subtract(phases, func(time, *popt[-1]))).mean())) + used_functions.append(func) + except RuntimeError: + print('curvefit of %s fails' %func.__name__) + # get position of parameters in list which has the lowest calculated error, Michael + pos = perr.index(min(perr)) + # return function name and parameter with the lowest calculated error during curvefit, Michael + return used_functions[pos], popt[pos] + + + + # scaling time axes depending on the time unit + def scaling_time_axes(self, time_sec, time_unit): + if time_unit == 'min': + time_vals = time_sec/60 + + elif time_unit == 'hours': + time_vals = time_sec/3600 + + elif time_unit == 'sec': + time_vals = time_sec + + else: + time_vals = time_sec/60 + + return time_vals + + # add 5 % of the distance between min and max to the range + @staticmethod + def get_extended_min_max(array): + distance = array.max() - array.min() + if distance == 0.: + distance = 1 + return array.min()-0.05*distance, array.max()+0.05*distance + + \ No newline at end of file diff --git a/Python_script/the.pdf b/Python_script/the.pdf index ac2e8eb8f746bcf84571f39dc69f2410d5f8c126..efd92b1502d1ef54657b7af20c41b5fa49f809df 100644 Binary files a/Python_script/the.pdf and b/Python_script/the.pdf differ