python应用 曲线拟合04

  1 import xlrd
  2 import xlwt
  3 import numpy as np
  4 import matplotlib.pyplot as plt
  5 
  6 
  7 def pol2(x, a, b, c):
  8     x = np.array(x)
  9     return a + b*x + c*pow(x, 2)
 10 
 11 
 12 def pol3(x, a, b, c, d):
 13     x = np.array(x)
 14     return a + b*x + c*pow(x, 2) + d*pow(x, 3)
 15 
 16 
 17 def pol4(x, a, b, c, d, e):
 18     x = np.array(x)
 19     return a + b*x + c*pow(x, 2) + d*pow(x, 3) + e*pow(x, 4)
 20 
 21 
 22 def pol5(x, a, b, c, d, e, f):
 23     x = np.array(x)
 24     return a + b*x + c*pow(x, 2) + d*pow(x, 3) + e*pow(x, 4) + f*pow(x, 5)
 25 
 26 
 27 def pol6(x, a, b, c, d, e, f, g):
 28     x = np.array(x)
 29     return a + b*x + c*pow(x, 2) + d*pow(x, 3) + e*pow(x, 4) + f*pow(x, 5) + 30         g*pow(x, 6)
 31 
 32 
 33 def main():
 34     data_file = xlrd.open_workbook(‘./data_78As.xlsx‘)
 35     sheet_proton = data_file.sheet_by_name(‘proton‘)
 36     sheet_neutron = data_file.sheet_by_name(‘neutron‘)
 37 
 38     write_excl = xlwt.Workbook(encoding=‘utf-8‘)
 39     write_sheet_proton = write_excl.add_sheet(‘proton‘)
 40     write_sheet_neutron = write_excl.add_sheet(‘neutron‘)
 41 
 42     x_proton = sheet_proton.col_values(0)
 43     x_proton_add = x_proton[:]
 44     x_proton_add.insert(0, 0.23)
 45     x_proton_add.insert(0, 0.22)
 46     x_proton_add.insert(0, 0.21)
 47     x_proton_add.insert(0, 0.20)
 48     x_proton_add.append(0.36)
 49     x_proton_add.append(0.37)
 50     x_proton_add.append(0.38)
 51     x_proton_add.append(0.39)
 52     x_proton_add.append(0.40)
 53 
 54     j = 0
 55     for temp in x_proton_add:
 56         write_sheet_proton.write(j, 0, temp)
 57         j = j + 1
 58 
 59     for i in range(18):
 60         y = sheet_proton.col_values(i+1)
 61         plt.scatter(x_proton, y)
 62         print(y)
 63         #  popt, pcov = curve_fit(pol5, x_proton, y, [1, 1, 1, 1, 1, 1])
 64         #  print(popt)
 65 
 66         z = np.polyfit(x_proton, y, 3)
 67 
 68         y_fit = pol3(x_proton_add, z[3], z[2], z[1], z[0])
 69         plt.plot(x_proton_add, y_fit, color=‘green‘, linewidth=1.0)
 70         plt.show()
 71         j = 0
 72         for temp in y_fit:
 73             print(x_proton_add[j])
 74             if x_proton_add[j] <= 0.23:
 75                 write_sheet_proton.write(j, i+1, temp)
 76                 j = j + 1
 77             elif x_proton_add[j] > 0.23 and x_proton_add[j] < 0.35:
 78                 write_sheet_proton.write(j, i+1, y[j-4])
 79                 j = j + 1
 80             else:
 81                 write_sheet_proton.write(j, i+1, temp)
 82                 j = j + 1
 83 
 84     plt.clf()
 85     x_neutron = sheet_neutron.col_values(0)
 86     x_neutron_add = x_neutron[:]
 87     x_neutron_add.insert(0, 0.23)
 88     x_neutron_add.insert(0, 0.22)
 89     x_neutron_add.insert(0, 0.21)
 90     x_neutron_add.insert(0, 0.20)
 91     x_neutron_add.append(0.36)
 92     x_neutron_add.append(0.37)
 93     x_neutron_add.append(0.38)
 94     x_neutron_add.append(0.39)
 95     x_neutron_add.append(0.40)
 96 
 97     j = 0
 98     for temp in x_neutron_add:
 99         write_sheet_neutron.write(j, 0, temp)
100         j = j + 1
101 
102     for i in range(18):
103         y = sheet_neutron.col_values(i+1)
104         plt.scatter(x_neutron, y)
105         print(y)
106         #  popt, pcov = curve_fit(pol5, x_proton, y, [1, 1, 1, 1, 1, 1])
107         #  print(popt)
108 
109         z = np.polyfit(x_neutron, y, 3)
110 
111         y_fit = pol3(x_neutron_add, z[3], z[2], z[1], z[0])
112         plt.plot(x_neutron_add, y_fit, color=‘green‘, linewidth=1.0)
113         plt.show()
114         j = 0
115         for temp in y_fit:
116             if x_neutron_add[j] <= 0.23:
117                 write_sheet_neutron.write(j, i+1, temp)
118                 j = j + 1
119             elif x_neutron_add[j] > 0.23 and x_neutron_add[j] < 0.35:
120                 write_sheet_neutron.write(j, i+1, y[j-4])
121                 j = j + 1
122             else:
123                 write_sheet_neutron.write(j, i+1, temp)
124                 j = j + 1
125 
126     write_excl.save("result_78As.xls")
127 
128 
129 if __name__ == ‘__main__‘:
130     main()