diff options
author | Eugeniy E. Mikhailov <evgmik@gmail.com> | 2024-07-13 21:04:38 -0400 |
---|---|---|
committer | Eugeniy E. Mikhailov <evgmik@gmail.com> | 2024-07-13 21:04:38 -0400 |
commit | fa086b54e06684dc8e666df6245629e896b215fd (patch) | |
tree | c197839cba4af4690a65cd00c141a7836ce0810e | |
parent | 705cc2eb112ab5e4c8329ab45fa2c89593048f2b (diff) | |
download | qolab-fa086b54e06684dc8e666df6245629e896b215fd.tar.gz qolab-fa086b54e06684dc8e666df6245629e896b215fd.zip |
black formatter
-rw-r--r-- | qolab/math/spectral_utils.py | 26 |
1 files changed, 14 insertions, 12 deletions
diff --git a/qolab/math/spectral_utils.py b/qolab/math/spectral_utils.py index ad5ce6e..797267e 100644 --- a/qolab/math/spectral_utils.py +++ b/qolab/math/spectral_utils.py @@ -1,6 +1,7 @@ import numpy as np from scipy.fft import rfft, rfftfreq + def spectrum(t, y): """ Calculate spectrum of real value signal with stripped away zero frequency component. @@ -8,23 +9,25 @@ def spectrum(t, y): rate and time interval. """ N = t.size - dt=np.mean(np.diff(t)) - freq=rfftfreq(N, dt) + dt = np.mean(np.diff(t)) + freq = rfftfreq(N, dt) # y= y - np.mean(y) yf = rfft(y) - yf *= 2/N; # produce normalized amplitudes - return freq[1:], yf[1:]; # strip of boring freq=0 + yf *= 2 / N # produce normalized amplitudes + return freq[1:], yf[1:] # strip of boring freq=0 + -def noise_density_spectrum(t,y): +def noise_density_spectrum(t, y): """ Calculate noise amplitude spectral density (ASD), the end results has unitis of y/sqrt(Hz) i.e. it does sqrt(PSD) where PSD is powerd spectrum density. Preserves the density independent of sampling rate and time interval. """ freq, yf = spectrum(t, y) - yf = yf*np.sqrt(t[-1]-t[0]) # scales with 1/sqrt(RBW) + yf = yf * np.sqrt(t[-1] - t[0]) # scales with 1/sqrt(RBW) return freq, yf + def noise_spectrum_smooth(fr, Ampl, Nbins=100): """ Smooth amplitude spectrum, especially at high frequency end. @@ -35,16 +38,15 @@ def noise_spectrum_smooth(fr, Ampl, Nbins=100): """ frEdges = np.logspace(np.log10(fr[0]), np.log10(fr[-1]), Nbins) - frCenter = np.zeros(frEdges.size-1) - power = np.zeros(frEdges.size-1) + frCenter = np.zeros(frEdges.size - 1) + power = np.zeros(frEdges.size - 1) for i, (frStart, frEnd) in enumerate(zip(frEdges[0:-1], frEdges[1:])): # print (f"{i=} {frStart=} {frEnd}") ind = (frStart <= fr) & (fr <= frEnd) - frCenter[i] = np.mean( fr[ind] ) - power [i] = np.mean( np.power( np.abs(Ampl[ind]),2) ) + frCenter[i] = np.mean(fr[ind]) + power[i] = np.mean(np.power(np.abs(Ampl[ind]), 2)) ind = np.logical_not(np.isnan(frCenter)) frCenter = frCenter[ind] power = power[ind] # print(f'{frCenter=} {power=}') - return frCenter, np.sqrt( power ) - + return frCenter, np.sqrt(power) |