From 54baef77335b98c92977d28a6a626cd31a3bec01 Mon Sep 17 00:00:00 2001
From: Marina Kreme <amamarinak@gmail.com>
Date: Tue, 1 Dec 2020 17:04:26 +0100
Subject: [PATCH] load mask
---
.../tfgm/scripts/exp_gabmul_eigs_properties.m | 120 ++++++++----------
.../rank_estimation_halko_vs_eigs_gausswin.m | 69 ++++++----
.../scripts/script_gabmul_eigs_properties.m | 48 ++++---
3 files changed, 124 insertions(+), 113 deletions(-)
diff --git a/matlab/tfgm/scripts/exp_gabmul_eigs_properties.m b/matlab/tfgm/scripts/exp_gabmul_eigs_properties.m
index 0a3a1d6..ebae14c 100644
--- a/matlab/tfgm/scripts/exp_gabmul_eigs_properties.m
+++ b/matlab/tfgm/scripts/exp_gabmul_eigs_properties.m
@@ -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');
diff --git a/matlab/tfgm/scripts/rank_estimation_halko_vs_eigs_gausswin.m b/matlab/tfgm/scripts/rank_estimation_halko_vs_eigs_gausswin.m
index ae7c2a8..42a925b 100644
--- a/matlab/tfgm/scripts/rank_estimation_halko_vs_eigs_gausswin.m
+++ b/matlab/tfgm/scripts/rank_estimation_halko_vs_eigs_gausswin.m
@@ -1,38 +1,50 @@
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');
diff --git a/matlab/tfgm/scripts/script_gabmul_eigs_properties.m b/matlab/tfgm/scripts/script_gabmul_eigs_properties.m
index 078d8ab..0590e0f 100644
--- a/matlab/tfgm/scripts/script_gabmul_eigs_properties.m
+++ b/matlab/tfgm/scripts/script_gabmul_eigs_properties.m
@@ -1,5 +1,5 @@
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
--
GitLab