Chcę dodać losowy szum do około 100 bin sygnału, który symuluję w Pythonie - żeby był bardziej realistyczny.
Na podstawowym poziomie, moją pierwszą myślą było przejście do segmentu segment po segmencie i po prostu wygenerowanie liczby losowej z określonego zakresu i dodanie lub odjęcie tej liczby od sygnału.
Miałem nadzieję (ponieważ to jest Python), że może być bardziej inteligentny sposób na zrobienie tego za pomocą numpy lub czegoś podobnego. (Przypuszczam, że najlepiej byłoby, gdyby liczba pobrana z rozkładu Gaussa i dodana do każdego pojemnika byłaby lepsza).
Z góry dziękuję za wszelkie odpowiedzi.
Jestem dopiero na etapie planowania kodu, więc nie mam nic do pokazania. Właśnie pomyślałem, że może istnieć bardziej wyrafinowany sposób generowania hałasu.
Pod względem wydajności, gdybym miał 10 pojemników o następujących wartościach:
Pojemnik 1: 1 Pojemnik 2: 4 Pojemnik 3: 9 Pojemnik 4:16 Pojemnik 5:25 Pojemnik 6:25 Pojemnik 7:16 Pojemnik 8: 9 Pojemnik 9: 4 Pojemnik 10: 1
Zastanawiałem się tylko, czy istnieje wstępnie zdefiniowana funkcja, która mogłaby dodać szum, aby dać mi coś takiego:
Pojemnik 1: 1,13 Pojemnik 2: 4,21 Pojemnik 3: 8,79 Pojemnik 4: 16,08 Pojemnik 5: 24,97 Pojemnik 6: 25,14 Pojemnik 7: 16,22 Pojemnik 8: 8,90 Pojemnik 9: 4,02 Pojemnik 10: 0,91
Jeśli nie, po prostu przejdę bin-by-bin i dodam do każdego liczbę wybraną z rozkładu Gaussa.
Dziękuję Ci.
Właściwie to sygnał z radioteleskopu, który symuluję. Chcę móc ostatecznie wybrać stosunek sygnału do szumu mojej symulacji.
Odpowiedzi:
Możesz wygenerować tablicę szumów i dodać ją do swojego sygnału
import numpy as np noise = np.random.normal(0,1,100) # 0 is the mean of the normal distribution you are choosing from # 1 is the standard deviation of the normal distribution # 100 is the number of elements you get in array noiseźródło
... a dla tych, którzy - jak ja - są na bardzo wczesnym etapie uczenia się odrętwienia,
import numpy as np pure = np.linspace(-1, 1, 100) noise = np.random.normal(0, 1, 100) signal = pure + noiseźródło
Dla tych, którzy próbują nawiązać połączenie między SNR z normalną zmienną losową generowaną przez numpy:
[1]
, gdzie należy pamiętać, że P to moc średnia .
Lub w dB:
[2]
W tym przypadku mamy już sygnał i chcemy wygenerować szum, aby uzyskać pożądany SNR.
Chociaż szum może mieć różne smaki w zależności od tego, co modelujesz, dobrym początkiem (szczególnie w przypadku tego przykładu radioteleskopu) jest Addytywny Biały Szum Gaussowski (AWGN) . Jak stwierdzono w poprzednich odpowiedziach, aby modelować AWGN, należy dodać zerową średnią gaussowską zmienną losową do oryginalnego sygnału. Wariancja tej zmiennej losowej wpłynie na średnią moc szumów.
Dla losowej zmiennej Gaussa X, średnia moc
, znana również jako drugi moment , wynosi 
[3]
Więc dla białego szumu,
a średnia moc jest równa wariancji 
.
Podczas modelowania tego w Pythonie możesz:
1. Obliczyć wariancję na podstawie pożądanego SNR i zestawu istniejących pomiarów, co zadziała, jeśli spodziewasz się, że pomiary będą miały dość spójne wartości amplitudy.
2. Alternatywnie, możesz ustawić moc szumu na znanym poziomie, aby dopasować coś takiego jak szum odbiornika. Szum odbiornika można zmierzyć, kierując teleskop w wolną przestrzeń i obliczając średnią moc.
Tak czy inaczej, ważne jest, aby upewnić się, że dodajesz szum do sygnału i bierzesz średnie w przestrzeni liniowej, a nie w jednostkach dB.
Oto kod generujący sygnał i wykreślający napięcie, moc w watach i moc w dB:
# Signal Generation # matplotlib inline import numpy as np import matplotlib.pyplot as plt t = np.linspace(1, 100, 1000) x_volts = 10*np.sin(t/(2*np.pi)) plt.subplot(3,1,1) plt.plot(t, x_volts) plt.title('Signal') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() x_watts = x_volts ** 2 plt.subplot(3,1,2) plt.plot(t, x_watts) plt.title('Signal Power') plt.ylabel('Power (W)') plt.xlabel('Time (s)') plt.show() x_db = 10 * np.log10(x_watts) plt.subplot(3,1,3) plt.plot(t, x_db) plt.title('Signal Power in dB') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()Oto przykład dodawania AWGN na podstawie pożądanego SNR:
# Adding noise using target SNR # Set a target SNR target_snr_db = 20 # Calculate signal power and convert to dB sig_avg_watts = np.mean(x_watts) sig_avg_db = 10 * np.log10(sig_avg_watts) # Calculate noise according to [2] then convert to watts noise_avg_db = sig_avg_db - target_snr_db noise_avg_watts = 10 ** (noise_avg_db / 10) # Generate an sample of white noise mean_noise = 0 noise_volts = np.random.normal(mean_noise, np.sqrt(noise_avg_watts), len(x_watts)) # Noise up the original signal y_volts = x_volts + noise_volts # Plot signal with noise plt.subplot(2,1,1) plt.plot(t, y_volts) plt.title('Signal with noise') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() # Plot in dB y_watts = y_volts ** 2 y_db = 10 * np.log10(y_watts) plt.subplot(2,1,2) plt.plot(t, 10* np.log10(y_volts**2)) plt.title('Signal with noise (dB)') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()A oto przykład dodawania AWGN na podstawie znanej mocy szumów:
# Adding noise using a target noise power # Set a target channel noise power to something very noisy target_noise_db = 10 # Convert to linear Watt units target_noise_watts = 10 ** (target_noise_db / 10) # Generate noise samples mean_noise = 0 noise_volts = np.random.normal(mean_noise, np.sqrt(target_noise_watts), len(x_watts)) # Noise up the original signal (again) and plot y_volts = x_volts + noise_volts # Plot signal with noise plt.subplot(2,1,1) plt.plot(t, y_volts) plt.title('Signal with noise') plt.ylabel('Voltage (V)') plt.xlabel('Time (s)') plt.show() # Plot in dB y_watts = y_volts ** 2 y_db = 10 * np.log10(y_watts) plt.subplot(2,1,2) plt.plot(t, 10* np.log10(y_volts**2)) plt.title('Signal with noise') plt.ylabel('Power (dB)') plt.xlabel('Time (s)') plt.show()źródło
Dla tych, którzy chcą dodać szum do wielowymiarowego zestawu danych załadowanego w ramce danych pandy, a nawet numpy ndarray, oto przykład:
import pandas as pd # create a sample dataset with dimension (2,2) # in your case you need to replace this with # clean_signal = pd.read_csv("your_data.csv") clean_signal = pd.DataFrame([[1,2],[3,4]], columns=list('AB'), dtype=float) print(clean_signal) """ print output: A B 0 1.0 2.0 1 3.0 4.0 """ import numpy as np mu, sigma = 0, 0.1 # creating a noise with the same dimension as the dataset (2,2) noise = np.random.normal(mu, sigma, [2,2]) print(noise) """ print output: array([[-0.11114313, 0.25927152], [ 0.06701506, -0.09364186]]) """ signal = clean_signal + noise print(signal) """ print output: A B 0 0.888857 2.259272 1 3.067015 3.906358 """źródło
W prawdziwym życiu chcesz zasymulować sygnał z białym szumem. Powinieneś dodać do swojego sygnału losowe punkty, które mają normalny rozkład Gaussa. Jeśli mówimy o urządzeniu, które ma czułość podaną w jednostkach / SQRT (Hz), to musisz obliczyć z niego odchylenie standardowe punktów. Tutaj podaję funkcję "white_noise", która robi to za Ciebie, reszta kodu to demonstracja i sprawdzam, czy robi to, co powinna.
%matplotlib inline import numpy as np import matplotlib.pyplot as plt from scipy import signal """ parameters: rhp - spectral noise density unit/SQRT(Hz) sr - sample rate n - no of points mu - mean value, optional returns: n points of noise signal with spectral noise density of rho """ def white_noise(rho, sr, n, mu=0): sigma = rho * np.sqrt(sr/2) noise = np.random.normal(mu, sigma, n) return noise rho = 1 sr = 1000 n = 1000 period = n/sr time = np.linspace(0, period, n) signal_pure = 100*np.sin(2*np.pi*13*time) noise = white_noise(rho, sr, n) signal_with_noise = signal_pure + noise f, psd = signal.periodogram(signal_with_noise, sr) print("Mean spectral noise density = ",np.sqrt(np.mean(psd[50:])), "arb.u/SQRT(Hz)") plt.plot(time, signal_with_noise) plt.plot(time, signal_pure) plt.xlabel("time (s)") plt.ylabel("signal (arb.u.)") plt.show() plt.semilogy(f[1:], np.sqrt(psd[1:])) plt.xlabel("frequency (Hz)") plt.ylabel("psd (arb.u./SQRT(Hz))") #plt.axvline(13, ls="dashed", color="g") plt.axhline(rho, ls="dashed", color="r") plt.show()źródło
Świetne odpowiedzi powyżej. Niedawno miałem potrzebę generowania symulowanych danych i właśnie z tego skorzystałem. Dzielenie się w razie potrzeby również innym,
import logging __name__ = "DataSimulator" logging.basicConfig(level=logging.INFO) logger = logging.getLogger(__name__) import numpy as np import pandas as pd def generate_simulated_data(add_anomalies:bool=True, random_state:int=42): rnd_state = np.random.RandomState(random_state) time = np.linspace(0, 200, num=2000) pure = 20*np.sin(time/(2*np.pi)) # concatenate on the second axis; this will allow us to mix different data # distribution data = np.c_[pure] mu = np.mean(data) sd = np.std(data) logger.info(f"Data shape : {data.shape}. mu: {mu} with sd: {sd}") data_df = pd.DataFrame(data, columns=['Value']) data_df['Index'] = data_df.index.values # Adding gaussian jitter jitter = 0.3*rnd_state.normal(mu, sd, size=data_df.shape[0]) data_df['with_jitter'] = data_df['Value'] + jitter index_further_away = None if add_anomalies: # As per the 68-95-99.7 rule(also known as the empirical rule) mu+-2*sd # covers 95.4% of the dataset. # Since, anomalies are considered to be rare and typically within the # 5-10% of the data; this filtering # technique might work #for us(https://en.wikipedia.org/wiki/68%E2%80%9395%E2%80%9399.7_rule) indexes_furhter_away = np.where(np.abs(data_df['with_jitter']) > (mu + 2*sd))[0] logger.info(f"Number of points further away : {len(indexes_furhter_away)}. Indexes: {indexes_furhter_away}") # Generate a point uniformly and embed it into the dataset random = rnd_state.uniform(0, 5, 1) data_df.loc[indexes_furhter_away, 'with_jitter'] += random*data_df.loc[indexes_furhter_away, 'with_jitter'] return data_df, indexes_furhter_awayźródło
AWGN Podobny do funkcji Matlab
def awgn(sinal): regsnr=54 sigpower=sum([math.pow(abs(sinal[i]),2) for i in range(len(sinal))]) sigpower=sigpower/len(sinal) noisepower=sigpower/(math.pow(10,regsnr/10)) noise=math.sqrt(noisepower)*(np.random.uniform(-1,1,size=len(sinal))) return noiseźródło