SR = 22050; % sampling rate
F0 = 50; Period = (1/F0)*SR;
NN = 3*SR;  % length of 3 seconds
% Make impulse train
stim = zeros(3*SR,1);
indices = round(Period:Period:(NN-Period)); Nindices = length(indices);
% Add random offsets
Rindices = round(indices + 0.35*Period*random("uniform", -1, 1, [1 Nindices]));
stim(Rindices) = 1;
% 100 msec onset and offset ramp
LL = round(0.1*SR); 
upramp = ((1:LL)/LL)'; downramp = ((LL:-1:1)/LL)';
stim(1:LL) = stim(1:LL).*upramp;
stim((NN-LL+1):NN) = stim((NN-LL+1):NN).*downramp;
% Low-pass filter, roll-off at 2000 Hz
  fc = 2000;
[b a] = butter(15, 2*fc/SR);
Xstim = filtfilt(b,a,stim);
wavwrite(Xstim, SR, "fry35_50.wav");

