fitt
  • AI Chat
  • Code
  • Report
  • Beta
    Spinner
    import numpy as np
    import pandas as pd
    from scipy.optimize import curve_fit
    import matplotlib.pyplot as plt
    
    
    df = pd.DataFrame({'energy': [81, 121.78, 244.31, 276.398, 302.85, 344.31, 356.01, 778.9, 963.39, 1112.07, 1408.01],
                       'efficiency_meas_ring1': [0.0031, 0.0044, 0.0056, 0.0052, 0.0050, 0.0046, 0.0045, 0.0030, 0.0026, 0.0024, 0.0021], 
                       'efficiency_meas_ring2': [0.0052, 0.0105, 0.0102, 0.0095, 0.0091, 0.0085, 0.0083, 0.0056, 0.0050, 0.0046, 0.0039], 
                      'yerr_1': [0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1], 
                      'yerr_2': [0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1]})
    
    def eff_fit(x, a, b, c, d, e, f, g):
        return (a + b*np.log(x) + c*np.log(x)**2 + d*np.log(x)**3 + e*np.log(x)**4 + f*np.log(x)**5 + g*np.log(x)**6)/x
        
    xdata = df['energy']
    ydata = df['efficiency_meas_ring1']*(1/0.012)
    ydata_2 = df['efficiency_meas_ring2']*(1/0.012)
    
    sigma_ring1 = df['yerr_1']*0.25
    sigma_ring2 = df['yerr_2']*0.25
    
    # Plot the actual data
    plt.plot(xdata, ydata, ".", label="Data_ring1");
    plt.plot(xdata, ydata_2, "*", label="Data_ring2");
    
    # The actual curve fitting happens here
    optimizedParameters1, pcov = curve_fit(eff_fit, xdata, ydata);
    optimizedParameters2, pcov = curve_fit(eff_fit, xdata, ydata_2);
    
    # Use the optimized parameters to plot the best fit
    x = np.linspace(80,1500,10000)
    plt.plot(x, eff_fit(x, *optimizedParameters1), label="fit_ring1");
    plt.plot(x, eff_fit(x, *optimizedParameters2), label="fit_ring2");
    
    plt.errorbar(xdata, ydata, sigma_ring1, fmt = 'o',color = 'orange', 
                ecolor = 'lightgreen', elinewidth = 2, capsize=5)
    plt.errorbar(xdata, ydata_2, sigma_ring2, fmt = 'o',color = 'orange', 
                ecolor = 'lightgreen', elinewidth = 2, capsize=5)
    
    
    # Show the graph
    plt.xlabel('E[keV]')
    plt.ylabel('efficiency [a.u]')
    plt.legend();
    plt.show();