%
%
%
%define time
Fs = 11025;
dt = 1/Fs;
t = 0:dt:10;

%define AP 1
AP1amp = 0.8;
AP1dur1 = .001;
AP1dur2 = .002;
AP1prob = .001;
%make single AP
%by concantenating two 1/2 sin cycles
oneT = 0:dt:(AP1dur1+AP1dur2);
oneAP1 = AP1amp*sin(2*pi/(2*AP1dur1)*oneT).*(oneT<=AP1dur1) + ...
    AP1amp*0.5*sin(2*pi/(2*AP1dur2)*(oneT+AP1dur1)).*(oneT>AP1dur1) ;

%define AP 2
AP2amp = 0.4;
AP2dur1 = .0015;
AP2dur2 = .003;
AP2prob = .001;
%make single AP2
%by concantenating two 1/2 sin cycles
oneT = 0:dt:(AP2dur1+AP2dur2);
oneAP2 = AP2amp*sin(2*pi/(2*AP2dur1)*oneT).*(oneT<=AP2dur1) + ...
    AP2amp*0.5*sin(2*pi/(2*AP2dur2)*(oneT+AP2dur1)).*(oneT>AP2dur1) ;

%make a bunch of AP1s
locationAP1 = rand(size(t))<AP1prob; %poisson
%locationAP1 = mod([1:length(t)],457)==0; %uniform
vAP = conv(locationAP1,oneAP1);
vAP = vAP(1:length(t)); %truncate (because of convolution)

%make and add a bunch of AP2s
locationAP2 = rand(size(t))<AP2prob; %poisson
vAP2 = conv(locationAP2,oneAP2);
vAP = vAP + vAP2(1:length(t)); %truncate and add

%make gaussian noise plus 60 Hz 120 vlot power line interference
thermalNoiseAmp = .05;
thermalNoiseCutoff = 3000 ; 
[b,a] = butter(2,thermalNoiseCutoff/(Fs/2));
lineNoiseAmp = .05;
vnoise = thermalNoiseAmp*filter(b,a,randn(size(t))) + lineNoiseAmp*sin(60*2*pi*t);

%now add niose to APs and 
vAPnoise =  vAP + vnoise;

% nplot=2
% figure(10);clf;
% 
% subplot(nplot,1,1)
% plot(t,vAP)
% title('poisson APs')
% set(gca,'ylim',[-1 1.5])
% 
% subplot(nplot,1,2)
% plot(t,vAPnoise)
% title('poisson APs with noise')
% set(gca,'ylim',[-1 1.5])

sound([vAPnoise',vAP'], Fs)


%