clear all;


%% T1. Convolution

% Params:
freq=5;                  % frequency of the sine wave (Hz)
sampRate=512;            % sampling rate

% Centered time vector (kernel length)
time = -2:1/sampRate:2;
% make a wavelet
c = 5;
s = c/(2*pi*freq);
cmw = exp(1i*2*pi*freq*time) .* exp(-time.^2./(2*s^2));

% Load the sample data:
load("FCz_sample.mat")

% Convolution and FFT parameters:
kernel = length(time);
htime = (length(time)-1)/2;  % half of the wavelet timing
ndata  = length(signal);
nconv  = (ndata + kernel) - 1;

% FFTs of kernel and signal
cx = fft(cmw,nconv);
dx = fft(signal,nconv);
% Normalize kernel (optional)
cx = cx./max(cx);

% Convolution
rc = cx .* dx;
as = ifft(rc,nconv);
% % Cut off wings
as = as(htime+1:end-htime);
% % Take a magnitude of the complex signal
as = abs(as).^2;

%% T2

% Load the sample data:
load("FCz_sample.mat")

% Wavelet parameters
sampRate = 512;       
mtime = -2:1/sampRate:2;
htime = (length(mtime)-1)/2;
num_frex    = 30;
min_freq    = 1;
max_freq    = 30;
frex        = logspace(log10(min_freq),log10(max_freq),num_frex); 
rcycles     = [3 10];
ncycles     = logspace(log10(rcycles(1)),log10(rcycles(end)),num_frex);

% Convolution and FFT parameters

% Create a set of wavelets
% preallocate arrays
cmwFFT  = zeros(length(frex),nconv);
for i = 1:length(frex)
%...
cmwFFT(i,:) = cmw./max(cmw);
end

% Convolution across frequencies
tf_result = zeros(length(frex),ndata);
for j = 1:length(frex)
%...
tf_result(j,:) = tf_power;
end


% Plot the result
timevec = linspace(1,ndata,ndata);
figure; clf;set(gcf, 'color','w', 'WindowState', 'maximize');
contourf(timevec,frex,tf_result,96,'linecolor','none','fill','on'),hold on 
set(gca,'ylim',ylim,'xlim',xlim,'clim',[0 20],'xtick',[1,1000,2000,3000,4000,ndata],'FontSize',35,...
    'LineWidth',3) 
    box off;colormap (jet)