Commit 54baef77 authored by Marina Kreme's avatar Marina Kreme
Browse files

load mask

parent 46411446
Pipeline #6048 passed with stage
in 4 minutes and 50 seconds
......@@ -4,13 +4,9 @@ clc; clear; close all;
% The mask used for figure 1 is the one that has been computed with
% the time-frequency parameters generated with a Gauss window of 256 time samples.
%% Load time-frequency mask from .mat file
mask_gauss = load('mask_with_gauss256_win.mat');
mask_hann = load('mask_with_hann512_win.mat');
%% Load time-frequency mask from .mat file
%save masks in a structure
mask_dict.('Gauss256') = mask_gauss;
mask_dict.('Hann512') = mask_hann;
mask_gauss = load('mask.mat');
%% directory for figures
......@@ -64,69 +60,63 @@ end
%% Compute the evd of the multiplier by using the DGT and IDGT
% generated from the Gauss window. The calculation of the evd is
% done for each of the two Gauss and Hann masks.
% generated from the Gauss window
%%
nb_eigvalues=4000;
masks = fieldnames(mask_dict);
evdn_gauss ={};
sig_len = param_dict.Gauss256.signal_params.sig_len;
for k=1: length(masks)
disp(k);
gab_mul = gen_gabmul_operator(param_dict.('Gauss256').dgt, ...,
param_dict.('Gauss256').idgt, mask_dict.(masks{k}).mask);
tic;
q_mat = randomized_range_finder(gab_mul, sig_len, nb_eigvalues);
t1= toc;
fprintf(" runtimes for Q: %f\n", t1);
tic;
evdn = EVD_nystrom(gab_mul, q_mat);
t2 = toc;
evdn_gauss.(masks{k}) = evdn;
fprintf(" runtimes for nystrom - gauss: %f\n", t2);
end
sig_len = param_dict.Gauss256.signal_params.sig_len; % signal length
gab_mul = gen_gabmul_operator(param_dict.('Gauss256').dgt, ...,
param_dict.('Gauss256').idgt, mask_gauss.mask);
tic;
q_mat = randomized_range_finder(gab_mul, sig_len, nb_eigvalues);
t1= toc;
fprintf(" runtimes for Q: %f\n", t1);
tic;
evdn_gauss = EVD_nystrom(gab_mul, q_mat);
t2 = toc;
fprintf(" runtimes for nystrom - gauss: %f\n", t2);
%% Compute the evd of the multiplier by using the DGT and IDGT
% generated from the Hanns window. The calculation of the evd is
% done for each of the two Gauss and Hann masks.
% generated from the Hann window. The calculation of the evd is
% done with the Gauss mask
%%
evdn_hann = {};
%%
for k=1: length(masks)
gab_mul = gen_gabmul_operator(param_dict.('Hann512').dgt, ...,
param_dict.('Hann512').idgt, mask_dict.(masks{k}).mask);
tic;
q_mat = randomized_range_finder(gab_mul, sig_len, nb_eigvalues);
t1= toc;
fprintf(" runtimes for Q: %f\n", t1);
tic;
evdn = EVD_nystrom(gab_mul, q_mat);
t2 = toc;
evdn_hann.((masks{k})) = evdn;
fprintf(" runtimes for nystrom-hann: %f\n", t2);
end
% gabor multiplier
gab_mul = gen_gabmul_operator(param_dict.('Hann512').dgt, ...,
param_dict.('Hann512').idgt, mask_gauss.mask);
%rand-evd
tic;
q_mat = randomized_range_finder(gab_mul, sig_len, nb_eigvalues); %rand-evd
t1= toc;
fprintf(" runtimes for Q: %f\n", t1);
tic;
evdn_hann = EVD_nystrom(gab_mul, q_mat);
t2 = toc;
fprintf(" runtimes for nystrom-hann: %f\n", t2);
%% eigenvalues
%%
figure;
D1 = diag(evdn_hann.Gauss256.D);
D4= diag(evdn_gauss.Gauss256.D);
semilogy(D4,'b', 'Linewidth',4);
hold on;
figure;
D1 = diag(evdn_hann.D);
D2= diag(evdn_gauss.D);
semilogy(D2,'b', 'Linewidth',4);
hold on;
semilogy(D1, 'r','Linewidth',4);
semilogy(D4(1),'bo', 'Linewidth', 15,'MarkerSize',10)
semilogy(D2(1),'bo', 'Linewidth', 15,'MarkerSize',10)
semilogy(D1(1),'r^','Linewidth', 15,'MarkerSize',10)
semilogy(3539,D4(3539),'bo','Linewidth', 15,'MarkerSize',10)
semilogy(3408,D1(3408),'r^','Linewidth', 15,'MarkerSize',10)
semilogy(3199,D2(3199),'bo','Linewidth', 15,'MarkerSize',10)
semilogy(3072,D1(3072),'r^','Linewidth', 15,'MarkerSize',10)
grid on;
......@@ -136,8 +126,8 @@ ylabel({' Eigenvalues $\sigma[k]$'},'Interpreter','latex');
set(gca, 'FontSize', 20, 'fontName','Times');
l = legend('Gauss','Hann',...,
'$\sigma[1]$ = 1, $\sigma[3539] = 2.699 \times 10^{-4}$',...,
'$\sigma[1]$ =1, $\sigma[3408] = 1.258 \times 10^{-4}$',...,
'$\sigma[1]$ = 1, $\sigma[3072] = 3.914 \times 10^{-4}$',...,
'$\sigma[1]$ =1, $\sigma[3199] = 3.865 \times 10^{-4}$',...,
'Location','southwest');
set(l, 'interpreter', 'latex')
......@@ -146,34 +136,34 @@ saveas(gcf,fullfile(fig_dir, 'eigenvalues_gauss_hann.png'));
%% eigenvectors
eigs_gauss = evdn_gauss.Gauss256.U;
eigs_gauss = evdn_gauss.U;
figure;
figure;
set(gcf,'position',[1, 1 1000 800]);
subplot(221);
plot_spectrogram(eigs_gauss(:,147), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,1), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(222);
plot_spectrogram(eigs_gauss(:,86), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,100), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(223);
plot_spectrogram(eigs_gauss(:,3039), param_dict.('Gauss256').dgt_params,...,
subplot(223);
plot_spectrogram(eigs_gauss(:,3199), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(224);
plot_spectrogram(eigs_gauss(:,3046), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,3072), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
ylabel('Frequency (kHz)')
yticks([0,1000,2000,3000,4000]);
......@@ -194,7 +184,7 @@ saveas(gcf,fullfile(fig_dir, 'mask_gauss.pdf'));
%%
save('evdn_gauss_hann.mat','evdn_hann','evdn_gauss','param_dict','mask_dict');
save('evdn_gauss_hann.mat','evdn_hann','evdn_gauss','param_dict','mask_gauss');
......
clc; clear; close all;
% The script permet d'tuider l'estimation du rang par Halko vs eigs
%% Repertoires pour les figures
% The script compares the rank of the Gabor multiplier estimated
% by rand-evd (see [1]) and eigs (see [2]).
%[1] *Finding structure with randomness: Probabilistic algorithms for`
% constructing approximate matrix decompositions*, by N. Halko, P.-G. Martinsson,
% and J. A. Tropp,
%[2]*A Krylov-Schur Algorithm for Large Eigenproblems*, by Stewart, G.W.
%% figures directory
pwd;
pathname ='figures_JSTSP';
if ~exist('figures_JSTSP','dir')
mkdir('figures_JSTSP');
pathname ='fig_rank_randEvd_vs_eigs';
if ~exist(pathname,'dir')
mkdir(pathname);
end
addpath('figures_JSTSP')
addpath(pathname)
%% Parametres du signal
%% Signal_params
resampling_fs=8000; % frequence d'echantillonage
sig_len = 8192; % longueur du signal
fs=8000;
sig_len = 16384;
signal_params = generate_signal_parameters(fs, sig_len);
tolerance = 1e-6; % parametres pour Halko
proba = 1-1e-4;
r = compute_r(sig_len, sig_len, proba);
%% Generation des parametres de la DGT
%% Rand-evd params
tolerance = 1e-6; %
proba = 1-1e-4;
r = compute_r(sig_len, sig_len, proba);
%% DGT params and operators
t_lim = [0.4, 0.6];
f_lims = 0.2:0.05:0.6;
win_type = 'gauss';
approx_win_len = 64; % changer sinon trop long.
win_len = 256;
params = get_params(win_len, win_type);
dgt_params = generate_dgt_parameters(win_type, approx_win_len);
signal_params = generate_signal_parameters(resampling_fs, sig_len);
[dgt, idgt] = get_stft_operators(dgt_params, signal_params);
hop = params.hop;
nbins =params.nbins;
dgt_params = generate_dgt_parameters(win_type, win_len, hop, nbins, sig_len);
[dgt, idgt] = get_stft_operators(dgt_params, signal_params);
%%
mask_area_list = zeros(length(f_lims),1);
......@@ -43,12 +55,12 @@ t_eigs = zeros(length(f_lims),1);
s_vec_list = cell(length(f_lims),1);
seuil = 10^(-14); % pour la EVD via eigs
thres = 10^(-14); % using for evd with eigs
%%
for k =1:length(f_lims)
fprintf("Je suis a l'iteration numero %.f patiente\n",k);
fprintf("first iteration %.f please wait\n",k);
% Mask Generation
%%
f_lim =[0.1, f_lims(k)];
mask = generate_rectangular_mask(dgt_params.nbins,dgt_params.hop, signal_params.sig_len, t_lim, f_lim);
......@@ -61,16 +73,17 @@ for k =1:length(f_lims)
end
figure(k);
plot_mask(mask, dgt_params.nbins, dgt_params.hop, signal_params.fs);
%% Gabor multiplier
gab_mul = gen_gabmul_operator(dgt, idgt, mask);
%% EVD via Halko
%Gabor multiplier
gab_mul = gen_gabmul_operator(dgt, idgt, mask);
%rand-evd
tic;
q_mat = adaptative_randomized_range_finder(gab_mul, sig_len, tolerance, r);
t_arrf(k) = toc;
ranks_arrf(k)= size(q_mat,2);
%% EVD via eigs
eigs
tic;
[u_mat,s] = eigs(gab_mul, signal_params.sig_len, signal_params.sig_len);
t_eigs(k) = toc;
......@@ -80,9 +93,7 @@ for k =1:length(f_lims)
ranks_eigs(k) = length(s_vec(s_vec > seuil));
end
%%
save('rank_estimation.mat','ranks_arrf','ranks_eigs', 'mask_area_list',...,
't_arrf','t_eigs', 's_vec_list');
%% plot figures
......@@ -97,3 +108,7 @@ legend('Rand-EVD', 'eigs', 'Location','northwest');
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname,'rank_estimation_gauss.png'));
saveas(gcf,fullfile(pathname,'rank_estimation_gauss.fig'));
%%
save('rank_estimation.mat','ranks_arrf','ranks_eigs', 'mask_area_list',...,
't_arrf','t_eigs', 's_vec_list');
clc; clear; close all;
%%
%%
load('evdn_gauss_hann.mat');
......@@ -12,18 +12,20 @@ end
addpath(fig_dir)
%%
%%
figure;
D1 = diag(evdn_hann.Gauss256.D);
D4= diag(evdn_gauss.Gauss256.D);
semilogy(D4,'b', 'Linewidth',4);
hold on;
figure;
D1 = diag(evdn_hann.D);
D2= diag(evdn_gauss.D);
semilogy(D2,'b', 'Linewidth',4);
hold on;
semilogy(D1, 'r','Linewidth',4);
semilogy(D4(1),'bo', 'Linewidth', 15,'MarkerSize',10)
semilogy(D2(1),'bo', 'Linewidth', 15,'MarkerSize',10)
semilogy(D1(1),'r^','Linewidth', 15,'MarkerSize',10)
semilogy(3539,D4(3539),'bo','Linewidth', 15,'MarkerSize',10)
semilogy(3408,D1(3408),'r^','Linewidth', 15,'MarkerSize',10)
semilogy(3199,D2(3199),'bo','Linewidth', 15,'MarkerSize',10)
semilogy(3072,D1(3072),'r^','Linewidth', 15,'MarkerSize',10)
grid on;
......@@ -33,40 +35,45 @@ ylabel({' Eigenvalues $\sigma[k]$'},'Interpreter','latex');
set(gca, 'FontSize', 20, 'fontName','Times');
l = legend('Gauss','Hann',...,
'$\sigma[1]$ = 1, $\sigma[3539] = 2.699 \times 10^{-4}$',...,
'$\sigma[1]$ =1, $\sigma[3408] = 1.258 \times 10^{-4}$',...,
'$\sigma[1]$ = 1, $\sigma[3072] = 3.914 \times 10^{-4}$',...,
'$\sigma[1]$ =1, $\sigma[3199] = 3.865 \times 10^{-4}$',...,
'Location','southwest');
set(l, 'interpreter', 'latex')
saveas(gcf,fullfile(fig_dir, 'eigenvalues_gauss_hann.fig'));
saveas(gcf,fullfile(fig_dir, 'eigenvalues_gauss_hann.png'));
%% eigenvectors
figure;
eigs_gauss = evdn_gauss.Gauss256.U;
eigs_gauss = evdn_gauss.U;
figure;
set(gcf,'position',[1, 1 1000 800]);
subplot(221);
plot_spectrogram(eigs_gauss(:,147), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,1), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(222);
plot_spectrogram(eigs_gauss(:,86), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,100), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(223);
plot_spectrogram(eigs_gauss(:,3039), param_dict.('Gauss256').dgt_params,...,
subplot(223);
plot_spectrogram(eigs_gauss(:,3199), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
ylabel('Frequency (kHz)')
set(gca, 'FontSize', 20, 'fontName','Times');
subplot(224);
plot_spectrogram(eigs_gauss(:,3046), param_dict.('Gauss256').dgt_params,...,
plot_spectrogram(eigs_gauss(:,3072), param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
ylabel('Frequency (kHz)')
yticks([0,1000,2000,3000,4000]);
......@@ -76,7 +83,6 @@ saveas(gcf,fullfile(fig_dir, 'eigvectors_prop_illustration.png'));
saveas(gcf,fullfile(fig_dir, 'eigvectors_prop_illustration.fig'));
%% mask
mask_gauss = mask_dict.Gauss256;
figure;
plot_spectrogram(mask_gauss.mask, param_dict.('Gauss256').dgt_params,...,
param_dict.('Gauss256').signal_params, param_dict.('Gauss256').dgt);
......@@ -84,7 +90,7 @@ ylabel('Frequency (kHz)')
yticks([0,1000,2000,3000,4000]);
yticklabels([0,1,2,3,4]);
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(fig_dir, 'mask_gauss.pdf'));
saveas(gcf,fullfile(fig_dir, 'mask_gauss.pngcd sc'));
%%
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment