diff --git a/XAS.py b/XAS.py index 8ca0b7cc683e59ca3c0a9bdc6127e60c70fb2ca4..2ff3c914555f848828a7075431e5ee6c8857af03 100644 --- a/XAS.py +++ b/XAS.py @@ -14,11 +14,15 @@ import matplotlib.gridspec as gridspec import matplotlib.pyplot as plt import re -def absorption(T, Io): - """ Compute the absorption A = -ln(T/Io) +def absorption(T, Io, fluorescence=False): + """ Compute the absorption A = -ln(T/Io) (or A = T/Io + for fluorescence) + Inputs: T: 1-D transmission value array of length N Io: 1-D Io monitor value array of length N + fluorescence: boolean, if False, compute A as + negative log, if True, compute A as ratio Output: a structured array with: @@ -61,14 +65,19 @@ def absorption(T, Io): p = np.corrcoef(T, Io)[0,1] - muA = -np.log(muT/muIo) - # from error propagation for correlated data - sigmaA = (np.sqrt((sigmaT/muT)**2 + (sigmaIo/muIo)**2 - - 2*p*sigmaIo*sigmaT/(muIo*muT))) + if fluorescence: + muA = muT/muIo + sigmaA = np.abs(muA)*(np.sqrt((sigmaT/muT)**2 + + (sigmaIo/muIo)**2 - 2*p*sigmaIo*sigmaT/(muIo*muT))) + else: + muA = -np.log(muT/muIo) + # from error propagation for correlated data + sigmaA = (np.sqrt((sigmaT/muT)**2 + (sigmaIo/muIo)**2 + - 2*p*sigmaIo*sigmaT/(muIo*muT))) - # direct calculation - #mask = (Io != 0) - #sigmaA = np.nanstd(-np.log(T[mask]/Io[mask])) + # direct calculation + #mask = (Io != 0) + #sigmaA = np.nanstd(-np.log(T[mask]/Io[mask])) return np.array([(muA, sigmaA, weights, muT, sigmaT, muIo, sigmaIo, p, counts)], dtype=fdtype) @@ -108,7 +117,8 @@ def binning(x, data, func, bins=100, bin_length=None): return bins, res -def xas(nrun, bins=None, Iokey='SCS_SA3', Itkey='MCP3apd', nrjkey='nrj', Iooffset=0, plot=False): +def xas(nrun, bins=None, Iokey='SCS_SA3', Itkey='MCP3apd', nrjkey='nrj', Iooffset=0, plot=False, + fluorescence = False): """ Compute the XAS spectra from a xarray nrun. Inputs: @@ -120,6 +130,8 @@ def xas(nrun, bins=None, Iokey='SCS_SA3', Itkey='MCP3apd', nrjkey='nrj', Iooffse NRJkey: string for the nrj fields, typically 'nrj' Iooffset: offset to apply on Io plot: boolean, displays a XAS spectrum if True + fluorescnce: boolean, if True, absorption is the ratio, + if False, absorption is negative log Outputs: a dictionnary containing: @@ -157,9 +169,9 @@ def xas(nrun, bins=None, Iokey='SCS_SA3', Itkey='MCP3apd', nrjkey='nrj', Iooffse if len(data) == 0: - return absorption([], []) + return absorption([], [], fluorescence) else: - return absorption(It_sign*data['It'], Io_sign*data['Io']) + return absorption(It_sign*data['It'], Io_sign*data['Io'], fluorescence) if bins is None: num_bins = 80 @@ -183,7 +195,10 @@ def xas(nrun, bins=None, Iokey='SCS_SA3', Itkey='MCP3apd', nrjkey='nrj', Iooffse gs = gridspec.GridSpec(2,1,height_ratios=[4,1]) ax1 = plt.subplot(gs[0]) ax1.plot(bins_c, muA, color='C1', label=r'$\sigma$') - ax1.set_ylabel('XAS') + if fluorescence: + ax1.set_ylabel('XAS (fluorescence)') + else: + ax1.set_ylabel('XAS (-log)') ax1.set_xlabel('Energy (eV)') ax1.legend() ax1_twin = ax1.twinx()