Stex6299

CreaDati - Risoluzione - AutoGUI

May 1st, 2021
1,019
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. # -*- coding: utf-8 -*-
  2. """
  3. Created on Sun Apr 25 10:44:22 2021
  4.  
  5. Autore: Stefano Carsi
  6. Data ultima modifica: 202104271819
  7.  
  8.    
  9. - Cosa ho fatto:
  10.    *
  11.    
  12. - Cosa devo fare:
  13.    *
  14.    
  15. - Per relazione:
  16.    *
  17.    
  18.    
  19.  
  20. - Info per l'uso:
  21.    > Serve un pizzico di manualità
  22.    
  23.    > Runnare la cella, aggiustare le cose e provare a far convergere
  24.      il plot
  25.      
  26.    > Solo dopo che si è soddisfatti, runnare la cella sottostante per memo-
  27.      rizzarsi i parametri ottimizzati dal fit
  28.      
  29.    > Infine esportare il file con tutti i dati
  30.    
  31.    > Le figure dei fit, se serve, meglio esportarle a mano
  32.    
  33.    > Nota bene: il vettore degli starting point non può essere vuoto!
  34.  
  35.  
  36.  
  37.    
  38. """
  39.  
  40.  
  41. #INIZIO STANDARD
  42. #%matplotlib auto
  43. import matplotlib.pyplot as plt
  44. import numpy as np
  45. from scipy.optimize import curve_fit
  46.  
  47. from matplotlib import use as mplUse
  48. mplUse("Qt5Agg")
  49. #%matplotlib inline
  50.  
  51.  
  52. import os
  53. path = r"F:\Corsi\Laboratorio 4 - Fisica nucleare e subnucleare\3 - Spettroscopia gamma\PITONE"
  54. os.chdir(path)
  55.  
  56. Export = False
  57.  
  58.  
  59.  
  60.  
  61. #%% Carico i dati
  62.  
  63. listaMatrici = []
  64. listaNum = []
  65.  
  66. for e in os.scandir(r"..\datadir_voltagescan_2021"):
  67.     if e.name.split(".")[1] != "dat": continue
  68.  
  69.     tmp = np.loadtxt(e.path)
  70.     listaMatrici.append(tmp.copy())
  71.     listaNum.append(int(e.name.split(".")[0][-2:]))
  72.    
  73. #listaNum = np.array(listaNum, dtype = np.int8)    
  74. zipIt = zip(listaNum, range(len(listaNum)))
  75. myDict = dict(zipIt)
  76.    
  77.  
  78.  
  79.  
  80. #%% Inizio - Modelli di fitting
  81.  
  82. from dataphile.statistics.regression.modeling import Parameter, Model, CompositeModel, AutoGUI
  83.  
  84.  
  85. def gauss(x,a,mu,sigma):
  86.     return a/(sigma*np.sqrt(2*np.pi))* np.exp(-(x-mu)**2 / (2*sigma**2))
  87.  
  88. def line(x, m, q):
  89.     return m*x+q
  90.  
  91. def exp(x,a,k):
  92.     return a*np.exp(-k*x)
  93.  
  94.  
  95.  
  96. #%% 52 - funziona, ma non usato
  97. numRun = 52
  98. dati = listaMatrici[myDict[numRun]][:,1]
  99. dati = dati[dati<4095]
  100.  
  101. h, bins = np.histogram(dati, bins=100)
  102. binc = bins[:-1] + (bins[1] - bins[0])
  103.    
  104.    
  105. fig = plt.figure()
  106. fig.set_size_inches(12,5)
  107. fig.suptitle("Titolo", fontsize=16)
  108.  
  109. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  110. ax.set_yscale("log")
  111. ax.grid(True)
  112. ax.set_xlabel("ADC")
  113. ax.set_ylabel("Entries")
  114.  
  115.  
  116.  
  117. # Plotto tutti i dati
  118. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  119. # Plotto i dati da fittare
  120.  
  121. cond = (binc>1273) & (binc<1560)
  122. xVect = binc[cond]
  123. yVect = h[cond]
  124.  
  125. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  126.                       ds="steps-mid", c="tab:green")
  127.  
  128.  
  129.  
  130.  
  131. model = CompositeModel(
  132.     Model(gauss,
  133.           Parameter(value=443026, bounds=(3e5, 3e6),      label='A'),
  134.           Parameter(value=1434.95,   bounds=(0, 4095),   label='mu'),
  135.           Parameter(value=42.4055,   bounds=(0, 300), label='sigma'),
  136.           label='gaussiana'),
  137.     Model(line,
  138.           Parameter(value=-2.54899, bounds=(-10, 10),  label='m'),
  139.           Parameter(value=4205.7, bounds=(-100, 10000), label='q'),
  140.           label='retta'),
  141.     label='gaussian_peaks')
  142.  
  143.  
  144. xsample = np.linspace(0, 4095, 1500)
  145. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', label='model')
  146. ax.legend();
  147.  
  148.  
  149. fig
  150.  
  151.  
  152. ax.set_position([0.15, 0.30, 0.84, 0.56])
  153. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.12], figure=fig,
  154.               slider_options={'color': 'steelblue'}, data=(xVect, yVect));
  155.  
  156.  
  157. A_Riassunto = model.summary()
  158.  
  159.  
  160.  
  161. #%% 52 Ufficiale (old)
  162.  
  163. numRun = 52
  164. dati = listaMatrici[myDict[numRun]][:,1]
  165. dati = dati[dati<4095]
  166.  
  167. h, bins = np.histogram(dati, bins=150)
  168. binc = bins[:-1] + (bins[1] - bins[0])
  169.    
  170.    
  171. fig = plt.figure()
  172. fig.set_size_inches(12,5)
  173. fig.suptitle("226Ra - 600V - Gain: 1000 ", fontsize=16)
  174.  
  175. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  176. #ax.set_yscale("log")
  177. ax.grid(True)
  178. ax.set_xlabel("ADC", fontsize=14)
  179. ax.set_ylabel("Entries", fontsize=14)
  180.  
  181.  
  182.  
  183. # Plotto tutti i dati
  184. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  185. # Plotto i dati da fittare
  186.  
  187. cond = (binc<1560) & (binc>243)
  188. xVect = binc[cond]
  189. yVect = h[cond]
  190.  
  191. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  192.                       ds="steps-mid", c="tab:green")
  193.  
  194.  
  195.  
  196.  
  197. model = CompositeModel(
  198.     Model(gauss,
  199.           Parameter(value=463743, bounds=(1e5, 3e6),      label='A'),
  200.           Parameter(value=304.195,   bounds=(0, 4095),   label='mu'),
  201.           Parameter(value=16.2228,   bounds=(0, 100), label='sigma'),
  202.           label='gaussiana 1'),
  203.     Model(gauss,
  204.           Parameter(value=128113, bounds=(1e5, 3e6),      label='A'),
  205.           Parameter(value=536.718,   bounds=(0, 4095),   label='mu'),
  206.           Parameter(value=21.1004,   bounds=(0, 100), label='sigma'),
  207.           label='gaussiana 2'),
  208.     Model(gauss,
  209.           Parameter(value=134248, bounds=(1e5, 3e6),      label='A'),
  210.           Parameter(value=655.77,   bounds=(0, 4095),   label='mu'),
  211.           Parameter(value=21.771,   bounds=(0, 100), label='sigma'),
  212.           label='gaussiana 3'),
  213.     Model(gauss,
  214.           Parameter(value=296409, bounds=(1e5, 3e6),      label='A'),
  215.           Parameter(value=768.839,   bounds=(0, 4095),   label='mu'),
  216.           Parameter(value=25.0682,   bounds=(0, 100), label='sigma'),
  217.           label='gaussiana 4'),
  218.     Model(gauss,
  219.           Parameter(value=474508, bounds=(1e5, 3e6),      label='A'),
  220.           Parameter(value=888.935,   bounds=(0, 4095),   label='mu'),
  221.           Parameter(value=28.5505,   bounds=(0, 100), label='sigma'),
  222.           label='gaussiana 5'),
  223.     Model(gauss,
  224.           Parameter(value=276598, bounds=(1e5, 3e6),      label='A'),
  225.           Parameter(value=1427.27,   bounds=(0, 4095),   label='mu'),
  226.           Parameter(value=40.4266,   bounds=(0, 100), label='sigma'),
  227.           label='gaussiana 6'),
  228.  
  229. #    Model(line,
  230. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  231. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  232. #          label='retta'),
  233.     Model(exp,
  234.           Parameter(value=9142.08, bounds=(0, 1e5),  label='a'),
  235.           Parameter(value=0.00208337, bounds=(0, 1e-2), label='k'),
  236.           label='exp'),
  237.     label='gaussian_peaks')
  238.  
  239.  
  240. xsample = np.linspace(0, 4095, 1500)
  241. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  242.                        label='Curva fittata')
  243. ax.legend();
  244.  
  245.  
  246. fig
  247.  
  248.  
  249. ax.set_position([0.15, 0.30, 0.84, 0.56])
  250. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.12], figure=fig,
  251.               slider_options={'color': 'steelblue'}, data=(xVect, yVect));
  252.  
  253. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  254. # è facilmente castabile in un array di numpy
  255. Riassunto_52 = model.summary()
  256. Riassunto52 = np.array(Riassunto_52)
  257.  
  258. if Export: pass
  259.     #plt.savefig(r'../RELAZIONE/FIGURE/risoluzione-52.eps', format='eps')
  260. plt.show()
  261.  
  262.  
  263.  
  264. #%% 52 nuovo ufficiale OK
  265. # Nota, ho alzato il bin150->200, ora converge
  266.  
  267. # Carico i dati
  268. numRun = 52
  269. dati = listaMatrici[myDict[numRun]][:,1]
  270. dati = dati[dati<4095]
  271.  
  272.  
  273. # Costruisco l'istogramma e trovo il centro
  274. h, bins = np.histogram(dati, bins=200)
  275. binc = bins[:-1] + (bins[1] - bins[0])
  276.    
  277.    
  278. # Creo figura e assi su cui disegnare
  279. fig = plt.figure()
  280. fig.set_size_inches(12,5)
  281. fig.suptitle("226Ra - 600V - Gain: 1000 ", fontsize=16)
  282.  
  283. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  284. ax.set_yscale("log")
  285. ax.grid(True)
  286. ax.set_xlabel("ADC", fontsize=14)
  287. ax.set_ylabel("Entries", fontsize=14)
  288.  
  289.  
  290.  
  291. # Plotto tutti i dati
  292. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  293.  
  294.  
  295. # Seleziono e Plotto i dati da fittare
  296. cond = (binc<2633) & (binc>243)
  297. xVect = binc[cond]
  298. yVect = h[cond]
  299.  
  300. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  301.                       ds="steps-mid", c="tab:green")
  302.  
  303.  
  304.    
  305. # Starting points: inseriti copiando quelli di convergenza
  306. startPars = np.array([
  307.        3.14020187e+05, 2.98589322e+02, 1.43257226e+01,
  308.        9.41411082e+04, 5.34214983e+02, 2.02581202e+01,
  309.        1.04459955e+05, 6.51821651e+02, 2.23514168e+01,
  310.        
  311.        2.30961972e+05, 7.65474019e+02, 2.49634128e+01,
  312.        3.65195534e+05, 8.86502718e+02, 2.82249031e+01,
  313.        2.14699579e+05, 1.42224403e+03, 4.12398651e+01,
  314.        
  315.        1.13048417e+05, 1.76150564e+03, 1.49847505e+02,
  316.        6.24648382e+04, 2.11018949e+03, 1.12551885e+02,
  317.        9.82593192e+04, 2.47588684e+03, 1.06711679e+02,
  318.        
  319.        7.81582354e+03, 2.32745796e-03])
  320.  
  321.  
  322. # Movello di fitting
  323. model = CompositeModel(
  324.     Model(gauss,
  325.           Parameter(value=startPars[0], bounds=(1e5, 3e6),      label='A'),
  326.           Parameter(value=startPars[1],   bounds=(0, 4095),   label='mu'),
  327.           Parameter(value=startPars[2],   bounds=(0, 100), label='sigma'),
  328.           label='gaussiana 1'),
  329.     Model(gauss,
  330.           Parameter(value=startPars[3], bounds=(1e5, 3e6),      label='A'),
  331.           Parameter(value=startPars[4],   bounds=(0, 4095),   label='mu'),
  332.           Parameter(value=startPars[5],   bounds=(0, 100), label='sigma'),
  333.           label='gaussiana 2'),
  334.     Model(gauss,
  335.           Parameter(value=startPars[6], bounds=(1e5, 3e6),      label='A'),
  336.           Parameter(value=startPars[7],   bounds=(0, 4095),   label='mu'),
  337.           Parameter(value=startPars[8],   bounds=(0, 100), label='sigma'),
  338.           label='gaussiana 3'),
  339.    
  340.     Model(gauss,
  341.           Parameter(value=startPars[9], bounds=(1e5, 3e6),      label='A'),
  342.           Parameter(value=startPars[10],   bounds=(0, 4095),   label='mu'),
  343.           Parameter(value=startPars[11],   bounds=(0, 100), label='sigma'),
  344.           label='gaussiana 4'),
  345.     Model(gauss,
  346.           Parameter(value=startPars[12], bounds=(1e5, 3e6),      label='A'),
  347.           Parameter(value=startPars[13],   bounds=(0, 4095),   label='mu'),
  348.           Parameter(value=startPars[14],   bounds=(0, 100), label='sigma'),
  349.           label='gaussiana 5'),
  350.     Model(gauss,
  351.           Parameter(value=startPars[15], bounds=(1e5, 3e6),      label='A'),
  352.           Parameter(value=startPars[16],   bounds=(0, 4095),   label='mu'),
  353.           Parameter(value=startPars[17],   bounds=(0, 100), label='sigma'),
  354.           label='gaussiana 6'),
  355.    
  356.     Model(gauss,
  357.           Parameter(value=startPars[18], bounds=(1e4, 3e6),      label='A'),
  358.           Parameter(value=startPars[19],   bounds=(0, 4095),   label='mu'),
  359.           Parameter(value=startPars[20],   bounds=(0, 100), label='sigma'),
  360.           label='gaussiana 7'),
  361.     Model(gauss,
  362.           Parameter(value=startPars[21], bounds=(1e4, 3e6),      label='A'),
  363.           Parameter(value=startPars[22],   bounds=(0, 4095),   label='mu'),
  364.           Parameter(value=startPars[23],   bounds=(0, 100), label='sigma'),
  365.           label='gaussiana 8'),
  366.     Model(gauss,
  367.           Parameter(value=startPars[24], bounds=(1e4, 3e6),      label='A'),
  368.           Parameter(value=startPars[25],   bounds=(0, 4095),   label='mu'),
  369.           Parameter(value=startPars[26],   bounds=(0, 100), label='sigma'),
  370.           label='gaussiana 9'),
  371.  
  372. #    Model(line,
  373. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  374. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  375. #          label='retta'),
  376.     Model(exp,
  377.           Parameter(value=startPars[27], bounds=(0, 1e5),  label='a'),
  378.           Parameter(value=startPars[28], bounds=(0, 1e-2), label='k'),
  379.           label='exp'),
  380.     label='gaussian_peaks')
  381.  
  382.  
  383. # Disegno la curva da muovere
  384. xsample = np.linspace(0, 4095, 1500)
  385. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  386.                        label='Curva fittata')
  387. ax.legend();
  388.  
  389.  
  390. # Mostro la figura
  391. fig
  392.  
  393.  
  394. # Ridimensiono gli assi per fare spazio ai cursori
  395. ax.set_position([0.10, 0.30, 0.84, 0.62])
  396.  
  397. # Chiamo ciò che fa la magia
  398. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.16], figure=fig,
  399.               slider_options={'color': 'steelblue'}, data=(xVect, yVect),
  400.               sigma=np.sqrt(yVect));
  401.  
  402. #%%
  403. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  404. # è facilmente castabile in un array di numpy
  405. Riassunto_52 = model.summary()
  406. Riassunto52 = np.array(Riassunto_52)
  407.  
  408. if Export: pass
  409.     #plt.savefig(r'../RELAZIONE/FIGURE/risoluzione-52.eps', format='eps')
  410. #plt.show()
  411.  
  412.  
  413.  
  414.  
  415.  
  416.  
  417. #%% 53 - Fatto OK
  418. # Nota, ho alzato il bin150->200, ora converge
  419.  
  420. # Carico i dati
  421. numRun = 53
  422. dati = listaMatrici[myDict[numRun]][:,1]
  423. dati = dati[dati<4095]
  424.  
  425.  
  426. # Costruisco l'istogramma e trovo il centro
  427. h, bins = np.histogram(dati, bins=200)
  428. binc = bins[:-1] + (bins[1] - bins[0])
  429.    
  430.    
  431. # Creo figura e assi su cui disegnare
  432. fig = plt.figure()
  433. fig.set_size_inches(12,5)
  434. fig.suptitle("226Ra - 650V - Gain: 500 ", fontsize=16)
  435.  
  436. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  437. ax.set_yscale("log")
  438. ax.grid(True)
  439. ax.set_xlabel("ADC", fontsize=14)
  440. ax.set_ylabel("Entries", fontsize=14)
  441.  
  442.  
  443.  
  444. # Plotto tutti i dati
  445. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  446.  
  447.  
  448. # Seleziono e Plotto i dati da fittare
  449. cond = (binc<2360) & (binc>220)
  450. xVect = binc[cond]
  451. yVect = h[cond]
  452.  
  453. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  454.                       ds="steps-mid", c="tab:green")
  455.  
  456.  
  457.  
  458. # Starting points: inseriti copiando quelli di convergenza
  459. startPars = np.array([
  460.        3.42507160e+05, 2.80817395e+02, 1.24616954e+01,
  461.        1.03919529e+05, 4.90186785e+02, 2.03356879e+01,
  462.        1.18326472e+05, 5.99110882e+02, 2.07618554e+01,
  463.        
  464.        2.37704890e+05, 7.02525659e+02, 2.16424925e+01,
  465.        3.85261776e+05, 8.10998553e+02, 2.55724833e+01,
  466.        2.21297335e+05, 1.29723065e+03, 3.59272711e+01,
  467.        
  468.        1.12844278e+05, 1.60399837e+03, 1.34508650e+02,
  469.        7.36554529e+04, 1.92969659e+03, 1.30558647e+02,
  470.        8.67953207e+04, 2.24743304e+03, 7.91527781e+01,
  471.        
  472.        8.00390614e+03, 2.48022403e-03])
  473.  
  474.  
  475. # Movello di fitting
  476. model = CompositeModel(
  477.     Model(gauss,
  478.           Parameter(value=startPars[0], bounds=(1e5, 3e6),      label='A'),
  479.           Parameter(value=startPars[1],   bounds=(0, 4095),   label='mu'),
  480.           Parameter(value=startPars[2],   bounds=(0, 100), label='sigma'),
  481.           label='gaussiana 1'),
  482.     Model(gauss,
  483.           Parameter(value=startPars[3], bounds=(1e5, 3e6),      label='A'),
  484.           Parameter(value=startPars[4],   bounds=(0, 4095),   label='mu'),
  485.           Parameter(value=startPars[5],   bounds=(0, 100), label='sigma'),
  486.           label='gaussiana 2'),
  487.     Model(gauss,
  488.           Parameter(value=startPars[6], bounds=(1e5, 3e6),      label='A'),
  489.           Parameter(value=startPars[7],   bounds=(0, 4095),   label='mu'),
  490.           Parameter(value=startPars[8],   bounds=(0, 100), label='sigma'),
  491.           label='gaussiana 3'),
  492.    
  493.     Model(gauss,
  494.           Parameter(value=startPars[9], bounds=(1e5, 3e6),      label='A'),
  495.           Parameter(value=startPars[10],   bounds=(0, 4095),   label='mu'),
  496.           Parameter(value=startPars[11],   bounds=(0, 100), label='sigma'),
  497.           label='gaussiana 4'),
  498.     Model(gauss,
  499.           Parameter(value=startPars[12], bounds=(1e5, 3e6),      label='A'),
  500.           Parameter(value=startPars[13],   bounds=(0, 4095),   label='mu'),
  501.           Parameter(value=startPars[14],   bounds=(0, 100), label='sigma'),
  502.           label='gaussiana 5'),
  503.     Model(gauss,
  504.           Parameter(value=startPars[15], bounds=(1e5, 3e6),      label='A'),
  505.           Parameter(value=startPars[16],   bounds=(0, 4095),   label='mu'),
  506.           Parameter(value=startPars[17],   bounds=(0, 100), label='sigma'),
  507.           label='gaussiana 6'),
  508.    
  509.     Model(gauss,
  510.           Parameter(value=startPars[18], bounds=(1e4, 3e6),      label='A'),
  511.           Parameter(value=startPars[19],   bounds=(0, 4095),   label='mu'),
  512.           Parameter(value=startPars[20],   bounds=(0, 100), label='sigma'),
  513.           label='gaussiana 7'),
  514.     Model(gauss,
  515.           Parameter(value=startPars[21], bounds=(1e4, 3e6),      label='A'),
  516.           Parameter(value=startPars[22],   bounds=(0, 4095),   label='mu'),
  517.           Parameter(value=startPars[23],   bounds=(0, 100), label='sigma'),
  518.           label='gaussiana 8'),
  519.     Model(gauss,
  520.           Parameter(value=startPars[24], bounds=(1e4, 3e6),      label='A'),
  521.           Parameter(value=startPars[25],   bounds=(0, 4095),   label='mu'),
  522.           Parameter(value=startPars[26],   bounds=(0, 100), label='sigma'),
  523.           label='gaussiana 9'),
  524.  
  525. #    Model(line,
  526. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  527. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  528. #          label='retta'),
  529.     Model(exp,
  530.           Parameter(value=startPars[27], bounds=(0, 1e5),  label='a'),
  531.           Parameter(value=startPars[28], bounds=(0, 1e-2), label='k'),
  532.           label='exp'),
  533.     label='gaussian_peaks')
  534.  
  535.  
  536. # Disegno la curva da muovere
  537. xsample = np.linspace(0, 4095, 1500)
  538. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  539.                        label='Curva fittata')
  540. ax.legend();
  541.  
  542.  
  543. # Mostro la figura
  544. fig
  545.  
  546.  
  547. # Ridimensiono gli assi per fare spazio ai cursori
  548. ax.set_position([0.10, 0.30, 0.84, 0.62])
  549.  
  550. # Chiamo ciò che fa la magia
  551. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.16], figure=fig,
  552.               slider_options={'color': 'steelblue'}, data=(xVect, yVect),
  553.               sigma=np.sqrt(yVect));
  554.  
  555. #%%
  556. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  557. # è facilmente castabile in un array di numpy
  558. Riassunto_53 = model.summary()
  559. Riassunto53 = np.array(Riassunto_53)
  560.  
  561.  
  562.  
  563.  
  564.  
  565.  
  566. #%% 54 - Fatto
  567. # I parametri soo fatti sui 200 bin
  568.  
  569. # Carico i dati
  570. numRun = 54
  571. dati = listaMatrici[myDict[numRun]][:,1]
  572. dati = dati[dati<4095]
  573.  
  574.  
  575. # Costruisco l'istogramma e trovo il centro
  576. h, bins = np.histogram(dati, bins=200)
  577. binc = bins[:-1] + (bins[1] - bins[0])
  578.    
  579.    
  580. # Creo figura e assi su cui disegnare
  581. fig = plt.figure()
  582. fig.set_size_inches(12,5)
  583. fig.suptitle("226Ra - 700V - Gain: 200 ", fontsize=16)
  584.  
  585. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  586. ax.set_yscale("log")
  587. ax.grid(True)
  588. ax.set_xlabel("ADC", fontsize=14)
  589. ax.set_ylabel("Entries", fontsize=14)
  590.  
  591.  
  592.  
  593. # Plotto tutti i dati
  594. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  595.  
  596.  
  597. # Seleziono e Plotto i dati da fittare
  598. cond = (binc<1705) & (binc>166)
  599. xVect = binc[cond]
  600. yVect = h[cond]
  601.  
  602. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  603.                       ds="steps-mid", c="tab:green")
  604.  
  605.    
  606.  
  607. # Starting points: inseriti copiando quelli di convergenza
  608. startPars = np.array([
  609.        4.62943107e+05, 2.26198954e+02, 1.42027020e+01,
  610.        4.30064153e+05, 3.51438064e+02, 6.29233950e+01,
  611.        1.10434659e+05, 4.57722754e+02, 1.32710721e+01,
  612.        
  613.        2.80277831e+05, 5.27133855e+02, 1.74496915e+01,
  614.        4.17788000e+05, 6.04375499e+02, 1.88606343e+01,
  615.        2.22374631e+05, 9.48896242e+02, 2.56032992e+01,
  616.        
  617.        7.15420610e+04, 1.16368668e+03, 7.36926220e+01,
  618.        6.46185002e+04, 1.38351571e+03, 8.84068589e+01,
  619.        8.02610692e+04, 1.61800918e+03, 5.50063439e+01,
  620.        
  621.        6.00904050e+03, 2.55138106e-03])
  622.  
  623.  
  624. # Movello di fitting
  625. model = CompositeModel(
  626.     Model(gauss,
  627.           Parameter(value=startPars[0], bounds=(1e5, 3e6),      label='A'),
  628.           Parameter(value=startPars[1],   bounds=(0, 4095),   label='mu'),
  629.           Parameter(value=startPars[2],   bounds=(0, 100), label='sigma'),
  630.           label='gaussiana 1'),
  631.     Model(gauss,
  632.           Parameter(value=startPars[3], bounds=(1e5, 3e6),      label='A'),
  633.           Parameter(value=startPars[4],   bounds=(0, 4095),   label='mu'),
  634.           Parameter(value=startPars[5],   bounds=(0, 100), label='sigma'),
  635.           label='gaussiana 2'),
  636.     Model(gauss,
  637.           Parameter(value=startPars[6], bounds=(1e5, 3e6),      label='A'),
  638.           Parameter(value=startPars[7],   bounds=(0, 4095),   label='mu'),
  639.           Parameter(value=startPars[8],   bounds=(0, 100), label='sigma'),
  640.           label='gaussiana 3'),
  641.    
  642.     Model(gauss,
  643.           Parameter(value=startPars[9], bounds=(1e5, 3e6),      label='A'),
  644.           Parameter(value=startPars[10],   bounds=(0, 4095),   label='mu'),
  645.           Parameter(value=startPars[11],   bounds=(0, 100), label='sigma'),
  646.           label='gaussiana 4'),
  647.     Model(gauss,
  648.           Parameter(value=startPars[12], bounds=(1e5, 3e6),      label='A'),
  649.           Parameter(value=startPars[13],   bounds=(0, 4095),   label='mu'),
  650.           Parameter(value=startPars[14],   bounds=(0, 100), label='sigma'),
  651.           label='gaussiana 5'),
  652.     Model(gauss,
  653.           Parameter(value=startPars[15], bounds=(1e5, 3e6),      label='A'),
  654.           Parameter(value=startPars[16],   bounds=(0, 4095),   label='mu'),
  655.           Parameter(value=startPars[17],   bounds=(0, 100), label='sigma'),
  656.           label='gaussiana 6'),
  657.    
  658.     Model(gauss,
  659.           Parameter(value=startPars[18], bounds=(1e4, 3e6),      label='A'),
  660.           Parameter(value=startPars[19],   bounds=(0, 4095),   label='mu'),
  661.           Parameter(value=startPars[20],   bounds=(0, 100), label='sigma'),
  662.           label='gaussiana 7'),
  663.     Model(gauss,
  664.           Parameter(value=startPars[21], bounds=(1e4, 3e6),      label='A'),
  665.           Parameter(value=startPars[22],   bounds=(0, 4095),   label='mu'),
  666.           Parameter(value=startPars[23],   bounds=(0, 100), label='sigma'),
  667.           label='gaussiana 8'),
  668.     Model(gauss,
  669.           Parameter(value=startPars[24], bounds=(1e4, 3e6),      label='A'),
  670.           Parameter(value=startPars[25],   bounds=(0, 4095),   label='mu'),
  671.           Parameter(value=startPars[26],   bounds=(0, 100), label='sigma'),
  672.           label='gaussiana 9'),
  673.  
  674. #    Model(line,
  675. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  676. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  677. #          label='retta'),
  678.     Model(exp,
  679.           Parameter(value=startPars[27], bounds=(0, 1e5),  label='a'),
  680.           Parameter(value=startPars[28], bounds=(0, 1e-2), label='k'),
  681.           label='exp'),
  682.     label='gaussian_peaks')
  683.  
  684.  
  685. # Disegno la curva da muovere
  686. xsample = np.linspace(0, 4095, 1500)
  687. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  688.                        label='Curva fittata')
  689. ax.legend();
  690.  
  691.  
  692. # Mostro la figura
  693. fig
  694.  
  695.  
  696. # Ridimensiono gli assi per fare spazio ai cursori
  697. ax.set_position([0.10, 0.30, 0.84, 0.62])
  698.  
  699. # Chiamo ciò che fa la magia
  700. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.16], figure=fig,
  701.               slider_options={'color': 'steelblue'}, data=(xVect, yVect),
  702.               sigma=np.sqrt(yVect));
  703.  
  704. #%%
  705. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  706. # è facilmente castabile in un array di numpy
  707. Riassunto_54 = model.summary()
  708. Riassunto54 = np.array(Riassunto_54)
  709.  
  710.  
  711.  
  712.  
  713.  
  714.  
  715. #%% 55 - Fatto OK
  716. # Nota, ho alzato il bin150->200, ora converge
  717.  
  718. # Carico i dati
  719. numRun = 55
  720. dati = listaMatrici[myDict[numRun]][:,1]
  721. dati = dati[dati<4095]
  722.  
  723.  
  724. # Costruisco l'istogramma e trovo il centro
  725. h, bins = np.histogram(dati, bins=200)
  726. binc = bins[:-1] + (bins[1] - bins[0])
  727.    
  728.    
  729. # Creo figura e assi su cui disegnare
  730. fig = plt.figure()
  731. fig.set_size_inches(12,5)
  732. fig.suptitle("226Ra - 750V - Gain: 200 ", fontsize=16)
  733.  
  734. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  735. ax.set_yscale("log")
  736. ax.grid(True)
  737. ax.set_xlabel("ADC", fontsize=14)
  738. ax.set_ylabel("Entries", fontsize=14)
  739.  
  740.  
  741.  
  742. # Plotto tutti i dati
  743. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  744.  
  745.  
  746. # Seleziono e Plotto i dati da fittare
  747. cond = (binc<2803) & (binc>243)
  748. xVect = binc[cond]
  749. yVect = h[cond]
  750.  
  751. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  752.                       ds="steps-mid", c="tab:green")
  753.  
  754.  
  755.  
  756. # Starting points: inseriti copiando quelli di convergenza
  757. startPars = np.array([
  758.        3.27007341e+05, 3.13976549e+02, 1.51769090e+01,
  759.        1.10288407e+05, 5.66753954e+02, 2.90578040e+01,
  760.        1.16014769e+05, 6.98831507e+02, 2.36358556e+01,
  761.        
  762.        2.39377275e+05, 8.23293555e+02, 2.62860390e+01,
  763.        3.79262612e+05, 9.53982266e+02, 2.88812872e+01,
  764.        2.18900196e+05, 1.53509437e+03, 4.27045957e+01,
  765.        
  766.        1.04806678e+05, 1.88950564e+03, 1.57225121e+02,
  767.        7.48170187e+04, 2.27747309e+03, 1.51051278e+02,
  768.        8.90622408e+04, 2.65905490e+03, 9.61389581e+01,
  769.        
  770.        6.20122878e+03, 2.03132197e-03])
  771.  
  772.  
  773. # Movello di fitting
  774. model = CompositeModel(
  775.     Model(gauss,
  776.           Parameter(value=startPars[0], bounds=(1e5, 3e6),      label='A'),
  777.           Parameter(value=startPars[1],   bounds=(0, 4095),   label='mu'),
  778.           Parameter(value=startPars[2],   bounds=(0, 100), label='sigma'),
  779.           label='gaussiana 1'),
  780.     Model(gauss,
  781.           Parameter(value=startPars[3], bounds=(1e5, 3e6),      label='A'),
  782.           Parameter(value=startPars[4],   bounds=(0, 4095),   label='mu'),
  783.           Parameter(value=startPars[5],   bounds=(0, 100), label='sigma'),
  784.           label='gaussiana 2'),
  785.     Model(gauss,
  786.           Parameter(value=startPars[6], bounds=(1e5, 3e6),      label='A'),
  787.           Parameter(value=startPars[7],   bounds=(0, 4095),   label='mu'),
  788.           Parameter(value=startPars[8],   bounds=(0, 100), label='sigma'),
  789.           label='gaussiana 3'),
  790.    
  791.     Model(gauss,
  792.           Parameter(value=startPars[9], bounds=(1e5, 3e6),      label='A'),
  793.           Parameter(value=startPars[10],   bounds=(0, 4095),   label='mu'),
  794.           Parameter(value=startPars[11],   bounds=(0, 100), label='sigma'),
  795.           label='gaussiana 4'),
  796.     Model(gauss,
  797.           Parameter(value=startPars[12], bounds=(1e5, 3e6),      label='A'),
  798.           Parameter(value=startPars[13],   bounds=(0, 4095),   label='mu'),
  799.           Parameter(value=startPars[14],   bounds=(0, 100), label='sigma'),
  800.           label='gaussiana 5'),
  801.     Model(gauss,
  802.           Parameter(value=startPars[15], bounds=(1e5, 3e6),      label='A'),
  803.           Parameter(value=startPars[16],   bounds=(0, 4095),   label='mu'),
  804.           Parameter(value=startPars[17],   bounds=(0, 100), label='sigma'),
  805.           label='gaussiana 6'),
  806.    
  807.     Model(gauss,
  808.           Parameter(value=startPars[18], bounds=(1e4, 3e6),      label='A'),
  809.           Parameter(value=startPars[19],   bounds=(0, 4095),   label='mu'),
  810.           Parameter(value=startPars[20],   bounds=(0, 100), label='sigma'),
  811.           label='gaussiana 7'),
  812.     Model(gauss,
  813.           Parameter(value=startPars[21], bounds=(1e4, 3e6),      label='A'),
  814.           Parameter(value=startPars[22],   bounds=(0, 4095),   label='mu'),
  815.           Parameter(value=startPars[23],   bounds=(0, 100), label='sigma'),
  816.           label='gaussiana 8'),
  817.     Model(gauss,
  818.           Parameter(value=startPars[24], bounds=(1e4, 3e6),      label='A'),
  819.           Parameter(value=startPars[25],   bounds=(0, 4095),   label='mu'),
  820.           Parameter(value=startPars[26],   bounds=(0, 100), label='sigma'),
  821.           label='gaussiana 9'),
  822.  
  823. #    Model(line,
  824. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  825. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  826. #          label='retta'),
  827.     Model(exp,
  828.           Parameter(value=startPars[27], bounds=(0, 1e5),  label='a'),
  829.           Parameter(value=startPars[28], bounds=(0, 1e-2), label='k'),
  830.           label='exp'),
  831.     label='gaussian_peaks')
  832.  
  833.  
  834. # Disegno la curva da muovere
  835. xsample = np.linspace(0, 4095, 1500)
  836. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  837.                        label='Curva fittata')
  838. ax.legend();
  839.  
  840.  
  841. # Mostro la figura
  842. fig
  843.  
  844.  
  845. # Ridimensiono gli assi per fare spazio ai cursori
  846. ax.set_position([0.10, 0.30, 0.84, 0.62])
  847.  
  848. # Chiamo ciò che fa la magia
  849. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.16], figure=fig,
  850.               slider_options={'color': 'steelblue'}, data=(xVect, yVect),
  851.               sigma=np.sqrt(yVect));
  852.  
  853. #%%
  854. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  855. # è facilmente castabile in un array di numpy
  856. Riassunto_55 = model.summary()
  857. Riassunto55 = np.array(Riassunto_55)
  858.  
  859.  
  860.  
  861.  
  862.  
  863.  
  864. #%% 56 - Fatto OK
  865. # Nota, ho alzato il bin150->200, ora converge
  866.  
  867. # Carico i dati
  868. numRun = 56
  869. dati = listaMatrici[myDict[numRun]][:,1]
  870. dati = dati[dati<4095]
  871.  
  872.  
  873. # Costruisco l'istogramma e trovo il centro
  874. h, bins = np.histogram(dati, bins=200)
  875. binc = bins[:-1] + (bins[1] - bins[0])
  876.    
  877.    
  878. # Creo figura e assi su cui disegnare
  879. fig = plt.figure()
  880. fig.set_size_inches(12,5)
  881. fig.suptitle("226Ra - 800V - Gain: 100 ", fontsize=16)
  882.  
  883. ax = fig.add_axes([0.15, 0.14, 0.84, 0.70])
  884. ax.set_yscale("log")
  885. ax.grid(True)
  886. ax.set_xlabel("ADC", fontsize=14)
  887. ax.set_ylabel("Entries", fontsize=14)
  888.  
  889.  
  890.  
  891. # Plotto tutti i dati
  892. ax.plot(binc, h, color='black', lw=1, label='Dati completi', ds="steps-mid")
  893.  
  894.  
  895. # Seleziono e Plotto i dati da fittare
  896. cond = (binc<2243) & (binc>213)
  897. xVect = binc[cond]
  898. yVect = h[cond]
  899.  
  900. data_graph, = ax.plot(xVect, yVect, lw=1, label='Dati per fit',\
  901.                       ds="steps-mid", c="tab:green")
  902.  
  903.  
  904.  
  905. # Starting points: inseriti copiando quelli di convergenza
  906. startPars = np.array([
  907.        3.41057648e+05, 2.74009366e+02, 1.21576501e+01,
  908.        9.65121696e+04, 4.77133343e+02, 1.87943137e+01,
  909.        1.23293444e+05, 5.84143244e+02, 1.98466608e+01,
  910.        
  911.        2.36748059e+05, 6.82714027e+02, 2.20560395e+01,
  912.        3.80932808e+05, 7.88170951e+02, 2.36810707e+01,
  913.        2.29361831e+05, 1.25312868e+03, 3.60282529e+01,
  914.        
  915.        9.47249003e+04, 1.53170346e+03, 1.09290881e+02,
  916.        8.10669238e+04, 1.83305278e+03, 1.16567200e+02,
  917.        8.98130756e+04, 2.13912681e+03, 7.72823389e+01,
  918.        
  919.        8.22077398e+03, 2.52499607e-03])
  920.  
  921.  
  922. # Movello di fitting
  923. model = CompositeModel(
  924.     Model(gauss,
  925.           Parameter(value=startPars[0], bounds=(1e5, 3e6),      label='A'),
  926.           Parameter(value=startPars[1],   bounds=(0, 4095),   label='mu'),
  927.           Parameter(value=startPars[2],   bounds=(0, 100), label='sigma'),
  928.           label='gaussiana 1'),
  929.     Model(gauss,
  930.           Parameter(value=startPars[3], bounds=(1e5, 3e6),      label='A'),
  931.           Parameter(value=startPars[4],   bounds=(0, 4095),   label='mu'),
  932.           Parameter(value=startPars[5],   bounds=(0, 100), label='sigma'),
  933.           label='gaussiana 2'),
  934.     Model(gauss,
  935.           Parameter(value=startPars[6], bounds=(1e5, 3e6),      label='A'),
  936.           Parameter(value=startPars[7],   bounds=(0, 4095),   label='mu'),
  937.           Parameter(value=startPars[8],   bounds=(0, 100), label='sigma'),
  938.           label='gaussiana 3'),
  939.    
  940.     Model(gauss,
  941.           Parameter(value=startPars[9], bounds=(1e5, 3e6),      label='A'),
  942.           Parameter(value=startPars[10],   bounds=(0, 4095),   label='mu'),
  943.           Parameter(value=startPars[11],   bounds=(0, 100), label='sigma'),
  944.           label='gaussiana 4'),
  945.     Model(gauss,
  946.           Parameter(value=startPars[12], bounds=(1e5, 3e6),      label='A'),
  947.           Parameter(value=startPars[13],   bounds=(0, 4095),   label='mu'),
  948.           Parameter(value=startPars[14],   bounds=(0, 100), label='sigma'),
  949.           label='gaussiana 5'),
  950.     Model(gauss,
  951.           Parameter(value=startPars[15], bounds=(1e5, 3e6),      label='A'),
  952.           Parameter(value=startPars[16],   bounds=(0, 4095),   label='mu'),
  953.           Parameter(value=startPars[17],   bounds=(0, 100), label='sigma'),
  954.           label='gaussiana 6'),
  955.    
  956.     Model(gauss,
  957.           Parameter(value=startPars[18], bounds=(1e4, 3e6),      label='A'),
  958.           Parameter(value=startPars[19],   bounds=(0, 4095),   label='mu'),
  959.           Parameter(value=startPars[20],   bounds=(0, 100), label='sigma'),
  960.           label='gaussiana 7'),
  961.     Model(gauss,
  962.           Parameter(value=startPars[21], bounds=(1e4, 3e6),      label='A'),
  963.           Parameter(value=startPars[22],   bounds=(0, 4095),   label='mu'),
  964.           Parameter(value=startPars[23],   bounds=(0, 100), label='sigma'),
  965.           label='gaussiana 8'),
  966.     Model(gauss,
  967.           Parameter(value=startPars[24], bounds=(1e4, 3e6),      label='A'),
  968.           Parameter(value=startPars[25],   bounds=(0, 4095),   label='mu'),
  969.           Parameter(value=startPars[26],   bounds=(0, 100), label='sigma'),
  970.           label='gaussiana 9'),
  971.  
  972. #    Model(line,
  973. #          Parameter(value=-4.55851, bounds=(-10, 10),  label='m'),
  974. #          Parameter(value=5881.62, bounds=(-100, 10000), label='q'),
  975. #          label='retta'),
  976.     Model(exp,
  977.           Parameter(value=startPars[27], bounds=(0, 1e5),  label='a'),
  978.           Parameter(value=startPars[28], bounds=(0, 1e-2), label='k'),
  979.           label='exp'),
  980.     label='gaussian_peaks')
  981.  
  982.  
  983. # Disegno la curva da muovere
  984. xsample = np.linspace(0, 4095, 1500)
  985. model_graph, = ax.plot(xsample, model(xsample), color='steelblue', \
  986.                        label='Curva fittata')
  987. ax.legend();
  988.  
  989.  
  990. # Mostro la figura
  991. fig
  992.  
  993.  
  994. # Ridimensiono gli assi per fare spazio ai cursori
  995. ax.set_position([0.10, 0.30, 0.84, 0.62])
  996.  
  997. # Chiamo ciò che fa la magia
  998. gui = AutoGUI(model, [model_graph], bbox=[0.20, 0.07, 0.75, 0.16], figure=fig,
  999.               slider_options={'color': 'steelblue'}, data=(xVect, yVect),
  1000.               sigma=np.sqrt(yVect));
  1001.  
  1002. #%%
  1003. # Eseguire a mano dopo aver premuto il pulsantone ed ispezionare a mano
  1004. # è facilmente castabile in un array di numpy
  1005. Riassunto_56 = model.summary()
  1006. Riassunto56 = np.array(Riassunto_56)
  1007.  
  1008.  
  1009.  
  1010. #%% Esporto i dati
  1011. np.savez("Dati", rias52 = Riassunto52, rias53 = Riassunto53, \
  1012.                  rias54 = Riassunto54, rias55 = Riassunto55, \
  1013.                  rias56 = Riassunto56, )
  1014.  
  1015.  
RAW Paste Data

Adblocker detected! Please consider disabling it...

We've detected AdBlock Plus or some other adblocking software preventing Pastebin.com from fully loading.

We don't have any obnoxious sound, or popup ads, we actively block these annoying types of ads!

Please add Pastebin.com to your ad blocker whitelist or disable your adblocking software.

×