Adding and generation of shot noise to a gaussian chirp signal.

Hi,

How can i adding shot noise to the gaussian chirp signal?

signal codes :

clc; clear; close all;

%--------------------------------------------------------------

fs = 1e3; dt = 1/fs; t = 0:dt:22.5e-3; % signal evaluation time

f1 = 29.25e6; f2 = 30.75e6;

c = (f2-f1)/(t(end)-t(1));

n = (-511:512)'; % Timeline

sigmas = exp(log(2):.3:log(200)); % Scales

w = exp(-.5*(n.^2)*sigmas.^(-2)); % Window Function

% Gaussian Chirp

chirps = w.*cos(n*sigmas.^(-1)*2*pi*5 + (c/2)*(n/sigmas').^2);

x = abs(chirps(:, 15)); % base signal

figure; plot(x); title('bidirectional gaussian chirp signal')

xlabel('Time(s)'); ylabel('Amplitude'); axis tight

%--------------------------------------------------------------

%adding shot noise

y = []; num_itter = 30; % number of itterations

for i = 1:num_itter

% append

y = [y; x];

end

Y = zeros(size(t)); Y(1:length(y)) = y;

figure; plot(t, abs(Y)); title('the final signal') axis tight

%//////////////////////////////////////////////////////////////////////////////////////////////////////

i write another code for generating this chirp gaussian signal but i can add shot noise to it:

clear all; clc

%--------------------------------------------------------------

% chirp signal generation

fs=1e3; t=0:1/fs:1;

f1=0; f2=20;

a=(f2-f1)/(t(end)-t(1)); f=f1+a*t;

x=cos(2*pi*f.*t);

subplot(2,1,1); plot(t,x); title(['down chirp signal']); xlabel('Time(s)'); ylabel('Amplitude');

% gaussian pulse generation

sigma=0.1; t1=-0.5:1/fs:0.5;

variance=sigma^2; w=1/(sqrt(2*pi*variance))*(exp(-t1.^2/(2*variance)));

subplot(2,1,2) plot(t1,w,'b'); title(['Gaussian Pulse \sigma=', num2str(sigma),'s']);

xlabel('Time(s)'); ylabel('Amplitude');

% gaussian chirp pulse generation

figure; xw = w(:) .* x(:);

plot(t,abs(xw)) xlim([0 1])\ xlabel('Time(s)'); ylabel('Amplitude');

title(['gaussian down chirp signal']);

%adding shot noise

y = []; num_itter = 32; % number of itterations

for i = 1:num_itter

% append

y = [y; xw];

end

T=0:1/fs:32.1; Y = zeros(size(T)); Y(1:length(y)) = y;

figure; plot(T,abs( Y));

title('the transmitting signal') xlabel('Time(s)'); ylabel('Amplitude');

axis tight

please help me!!!

thanks,Majid.

Answer 1

[Mathieu NOE]

hello

minor bug fixed in the last lines (see comments ) - it's the time vector definition (length) was incorrect

%--------------------------------------------------------------
% chirp signal generation
fs=1e3;  t=0:1/fs:1; 
f1=0; f2=20;
a=(f2-f1)/(t(end)-t(1));   f=f1+a*t;
x=cos(2*pi*f.*t);
subplot(2,1,1); plot(t,x); title(['down chirp signal']); xlabel('Time(s)'); ylabel('Amplitude');
% gaussian pulse generation
sigma=0.1; t1=-0.5:1/fs:0.5;
variance=sigma^2; w=1/(sqrt(2*pi*variance))*(exp(-t1.^2/(2*variance)));
subplot(2,1,2) ;plot(t1,w,'b'); title(['Gaussian Pulse \sigma=', num2str(sigma),'s']);
xlabel('Time(s)'); ylabel('Amplitude');
% gaussian chirp pulse generation
figure; xw = w(:) .* x(:);
plot(t,abs(xw));
xlim([0 1]); xlabel('Time(s)'); ylabel('Amplitude');
title(['gaussian down chirp signal']);
%adding shot noise
y = [];  num_itter = 32; % number of itterations
for i = 1:num_itter
    % append
    y = [y; xw]; 
end
T=(0:length(y)-1) *1/fs;                    % time axis fixed
% Y = zeros(size(T));  Y(1:length(y)) = y;  % no need for this line 
figure; plot(T,abs(y));
title('the transmitting signal') ;
xlabel('Time(s)'); ylabel('Amplitude');
axis tight

[majid]

thanks, Mathieu NOE

Do you have any idea about adding shot noise to this signal?

[Mathieu NOE]

hello

sorry but I don't understand what do you mean by "shot noise "

the code you created here has (to me) no special meaning except replicating num_itter times the signal (y)

%adding shot noise
y = [];  num_itter = 30; % number of itterations
for i = 1:num_itter
    % append
    y = [y; x]; 
end
Y = zeros(size(t));  Y(1:length(y)) = y; 
figure; plot(t, abs(Y)); title('the final signal') axis tight

so what is "shot noise" ?

[majid]

%adding shot noise >>>>>(I mean i want to add shot noise code here.)

Shot noise or Poisson noise is a type of noise which can be modeled by a Poisson process. In electronics shot noise originates from the discrete nature of electric charge. Shot noise also occurs in photon counting in optical devices, where shot noise is associated with the particle nature of light.

yes,this code is about replicating signal.

y = []; num_itter = 30; % number of itterations

for i = 1:num_itter

% append

y = [y; x];

end

Y = zeros(size(t)); Y(1:length(y)) = y;

figure; plot(t, abs(Y)); title('the final signal') axis tight

[Mathieu NOE]

hello

I found this Fex submission

Shot Noise Simulation - File Exchange - MATLAB Central (mathworks.com)

I slightly modified the code to not use exprnd which requires the Statistics and Machine Learning Toolbox.

% Simulate Shot Noise
% Here we generate a Poisson Impulse Process (a kind of white noise)
clear all; close all; clc

%%
mu = 1; % average wait time between events
N = 50; % number of events along the simulation
t_tick = 1e-3; % time tick for plotting/simulating
% time_between = exprnd(mu, N, 1); % random waiting times between events (exponential pdf)
time_between = -mu*log(rand(N, 1)); % random waiting times between events (exponential pdf)

time = cumsum( time_between ); % actual time of events is the cumulative sum of waiting times
time = round( time/t_tick )*t_tick; % round
time_data1 = [time'; ones(1,N)]'; % a matrix containing event times + '1' per each event
time_data2 = [ (0:t_tick:max(time) )'    zeros( 1+max(time)/t_tick, 1 ) ]; % a matrix containing non-event time ticks + '0' per each time tick. This is used for plotting
time_data3=[time_data1; time_data2]; % concatenate both matrices
time_data4=sortrows(time_data3,1); % sort according to time (column 1)
plot(time_data4(:,1), time_data4(:,2) )
ylim( [-0.1 1.1]); grid
ylabel('Events'); title(['Average=' num2str(1/mu) ' events/second']);

%% Now plot fft of the generated Poisson Impulse Process (white noise)
Y=abs(fft(time_data4(:,2)));
f=0:1/max(time):0.5/t_tick;
figure
plot( f,Y(1:length(f)) )
ylabel('psd'); xlabel('frequency [Hz]'); title('PSD');xlim([-10 500]); ylim([-1 54])

maybe you can add this portion of code to yours

[majid]

thank you!

I saw these codes before but I cant add it to my signal.