% time/freq reassignemnt
clear all

Fs = 22050 ;
tmax = 1.0 ;
t = 0:1/Fs:tmax ;

%first FM signal
FM1_center_freq = 2000 ;
FM1_depth = 0.1 ;
FM1_mod_freq = 25 ;
FM1_amp = 0.5 ;
AM1_freq = 5 ;
FM1_mod = FM1_depth * cumsum(cos(2*pi*FM1_mod_freq*t)) ;
y = FM1_amp * cos(2*pi*FM1_center_freq*t + FM1_mod) .* ...
cos(2*pi*AM1_freq*t).^2 .* (cos(2*pi*AM1_freq*t)>0) ;

%second FM signal
FM2_center_freq = 1500 ;
FM2_depth = 0.05 ;
FM2_mod_freq = 20 ;
FM2_amp = 0.25 ;
AM2_freq = 8 ;
FM2_mod = FM2_depth * cumsum(cos(2*pi*FM2_mod_freq*t)) ;
y = y + FM2_amp * cos(2*pi*FM2_center_freq*t + FM2_mod) .* ...
cos(2*pi*AM2_freq*t).^2 .* (cos(2*pi*AM2_freq*t)>0) ;

%additive, low passed, gaussian noise
noise = randn(size(t));
[b,a] = butter(4,0.2);
noise = filter(b,a,noise) ;
%y = y + 0.001*noise ;

%impulse-like noise burst
imp_phase = pi/2 ;
imp_amp = cos(2*pi*AM1_freq*t+imp_phase).^100 ...
.* (cos(2*pi*AM1_freq*t+imp_phase)>0) ;
%y = y + 0.5 * imp_amp .* noise ;

% linear chirp
y = y + .3*chirp(t, 1000, tmax, 3000) ;

%pname = 'C:\Documents and Settings\Bruce\My Documents\Octave\RA\' ;
%wavwrite(y, Fs, 16, [pname,'test2.wav']);

nfft = 256 ;
overlap = fix(0.8 * nfft) ;
ampthreshold = -40 ;
plot_freq_range = 4000 ;
% syntax is METHOD_MARKERSIZE
% 'color_5' means use color with markersize 5
% method options are 'fast', 'color'
render_method = 'fast_1_dot' ; 

%
figure(1); clf; 
subplot(2,1,1)
plot(t,y);
set(gca,'ylim', [-1,1]) ;
%set(gca,'xlim', [0 .05]) ;

subplot(2,1,2)
specgram(y, nfft, Fs, [], overlap);
set(gca,'ylim', [0, plot_freq_range]) ;
colormap gray
drawnow

figure(2);clf;
raspecgram(y, nfft, Fs, [], overlap, ampthreshold, render_method);
set(gca,'ylim', [0, plot_freq_range]) ;