Skip to content
Snippets Groups Projects
Commit a05439ee authored by Marina Kreme's avatar Marina Kreme
Browse files

delete file

parent 3cc04c80
No related branches found
No related tags found
No related merge requests found
Pipeline #5962 passed
Stage: test
Show whitespace changes
Inline Side-by-side
matlab/tfgm/scripts/illustration_with_cuicui.m deleted 100644 → 0
+ 0
− 378
View file @ 3cc04c80
clc; clear; close all;
%%
% experiments with
% - wideband sound: engine(Russia)
% - localized sound: birdsong
%% DGT and signals parameters
resampling_fs=8000;
sig_len = 16384;
win_type = 'hann';
approx_win_len = 128;
signal_params = generate_signal_parameters(resampling_fs, sig_len);
dgt_params = generate_dgt_parameters(win_type, approx_win_len);
%% Generate gabor frames operators
[direct_stft, adjoint_stft, tight_direct_stft, tight_adjoint_stft,win,...,
win_tight] = get_stft_operators(dgt_params, signal_params);
% plot associated window
figure;
plot_win(win_tight, signal_params.fs, signal_params.sig_len, win_type)
title([num2str(win_type), ' - window']);
saveas(gcf,fullfile(pathname, [num2str(win_type),'_window.png']));
%% Load signals
ind_engine = 4;
ind_bird = 5;
deb = 0;
[x_engine, x_bird] = load_pairs(ind_engine, ind_bird, resampling_fs, signal_params.sig_len, deb);
%gamma=0.75;
signals = generate_mix_signal(x_engine, x_bird);
%% plot signals spectrogrammes
% compute dgt
stft = tight_direct_stft;
tf_mat_engine = compute_dgt(signals.target, stft);
tf_mat_bird = compute_dgt(signals.noise, stft);
tf_mat_mix = compute_dgt(signals.mix, stft);
%plot their spectrogram
figure('name','engine'); plot_spectrogram(tf_mat_engine,dgt_params, signal_params, stft);
title('engine')
figure('name','bird'); plot_spectrogram(tf_mat_bird, dgt_params, signal_params, stft);
title('Bird')
figure('name','mix'); plot_spectrogram(tf_mat_mix, dgt_params,signal_params, stft);
title('mix : engine-bird')
%% Mix - SDR
sdr_mix = sdr(signals.target, signals.mix);
fprintf('The SDR of the mixture is : %e\n', sdr_mix)
%% Mask generation
% mask
alpha=2; seuil = 0.02; radius = 3;
mask = generate_mask(tf_mat_engine, tf_mat_bird, alpha, seuil, radius);
[mask_area, mask_area_ratio] = get_mask_area(mask);
%plot mask
figure('name','mask'); plot_spectrogram(mask, dgt_params,signal_params, stft);
title(['mask : mask-area = ',num2str(mask_area)])
figure;
plot_spectrogram((1-mask).*tf_mat_bird, dgt_params,signal_params, stft);
%% Halko parameters for EVD diagonalisation zone par zone
tolerance_arrf = 1e-6;
proba_arrf = 1 - 1e-9;
r = compute_r(signal_params.sig_len, signal_params.sig_len, proba_arrf);
[~,n_labels] = bwlabel(mask);
%n_labels =2;
t_arrf = zeros(n_labels,1);
t_evdn = zeros(n_labels,1);
s_vec = zeros(signal_params.sig_len,n_labels);
u_vec = zeros(signal_params.sig_len*3000,n_labels);
rank_q = zeros(n_labels,1);
mask_area_list = zeros(n_labels,1);
%%
for k = 1:n_labels
[mask_label,~] = bwlabel(mask);
% on construit chaque mask
mask_label(mask_label~=k)=0;
mask_ =mask_label;
figure(k); plotdgtreal(mask_, dgt_params.hop, dgt_params.nbins, signal_params.fs);
title(['k=', num2str(k)])
[mask_area_, mask_area_ratio_] = get_mask_area(mask_);
mask_area_list(k) = mask_area_;
fprintf('mask area = %.f\n',mask_area_);
%construction du multiplicateur associe a chaque masque
gab_mul = gen_gabmul_operator(tight_direct_stft, tight_adjoint_stft, mask_);
tic;
q_mat = adaptative_randomized_range_finder(gab_mul, sig_len, tolerance_arrf, r);
t_arrf(k) = toc;
rank_q(k) = size(q_mat,2);
tic;
evdn = EVD_nystrom(gab_mul, q_mat);
t_evdn(k) = toc;
s_vec(1:size(q_mat,2),k) = diag(evdn.D);
u_vec(1:size(q_mat,1)*size(q_mat,2),k) = evdn.U(:);
end
%% diagonalisation complete
gab_mul = gen_gabmul_operator(tight_direct_stft, tight_adjoint_stft, mask);
%stage 1
q_matc = adaptative_randomized_range_finder(gab_mul, signal_params.sig_len, tolerance_arrf, r);
fprintf('Q shape : %.f %.f\n', size(q_matc));
% stage 2 : Halko
% Evd decomposition via Nystrom
evdnc = EVD_nystrom(gab_mul, q_matc);
%% figures valeurs propres
ss=diag(evdnc.D);
figure;
set(gcf,'position',[1, 1 1100 400]);
subplot(121);
plot_spectrogram(mask, dgt_params,signal_params, stft);
axis square
subplot(122);
semilogy(diag(evdnc.D), 'Linewidth',2.5);
hold on;
plot(45,ss(45),'k-*','Linewidth',2.5);
plot(30,ss(30),'m-*','Linewidth',2.5);
plot(1650,ss(1650),'g-*','Linewidth',2.5);
plot(1740,ss(1740),'c-*','Linewidth',2.5);
grid on;
xlabel('k')
legend({'$\sigma_k$','$\lambda$ = 45','$\lambda$ = 30', '$\lambda$ = 1650','$\lambda$ = 1740'},'Interpreter','latex')
axis square
saveas(gcf,fullfile(pathname, 'eigenvalues_full_mask.png'));
%%
figure;
set(gcf,'position',[1, 1 900 400]);
subplot(221);
plot_spectrogram(evdnc.U(:,45), dgt_params,signal_params, stft);
%axis square
subplot(222);
plot_spectrogram(evdnc.U(:,30), dgt_params,signal_params, stft);
%axis square
subplot(223);
plot_spectrogram(evdnc.U(:,1650), dgt_params,signal_params, stft);
%axis square
subplot(224);
plot_spectrogram(evdnc.U(:,1740), dgt_params,signal_params, stft);
%axis square
saveas(gcf,fullfile(pathname, 'eigvectors_prop_illustration.png'));
%%
figure; sgram(evdnc.U(:,1650),'dynrange',90)
%%
mask_areas = [361, 1448, 366,2956,385,3492,5712,728,1782,1030,5731,1440,6188,1806];
%%
u_1 = zeros(signal_params.sig_len,rank_q(1));
u_2 = zeros(signal_params.sig_len,rank_q(2));
u_3 = zeros(signal_params.sig_len,rank_q(3));
u_4 = zeros(signal_params.sig_len,rank_q(4));
u_5 = zeros(signal_params.sig_len,rank_q(5));
u_6 = zeros(signal_params.sig_len,rank_q(6));
u_7 = zeros(signal_params.sig_len,rank_q(7));
u_8 = zeros(signal_params.sig_len,rank_q(8));
u_9 = zeros(signal_params.sig_len,rank_q(9));
u_10 = zeros(signal_params.sig_len,rank_q(10));
u_11 = zeros(signal_params.sig_len,rank_q(11));
u_12 = zeros(signal_params.sig_len,rank_q(12));
u_13 = zeros(signal_params.sig_len,rank_q(13));
u_14 = zeros(signal_params.sig_len,rank_q(14));
u_1 = reshape(u_vec(1:signal_params.sig_len*rank_q(1),1),size(u_1));
u_2 = reshape(u_vec(1:signal_params.sig_len*rank_q(2),2),size(u_2));
u_3 = reshape(u_vec(1:signal_params.sig_len*rank_q(3),3),size(u_3));
u_4 = reshape(u_vec(1:signal_params.sig_len*rank_q(4),4),size(u_4));
u_5 = reshape(u_vec(1:signal_params.sig_len*rank_q(5),5),size(u_5));
u_6 = reshape(u_vec(1:signal_params.sig_len*rank_q(6),6),size(u_6));
u_7 = reshape(u_vec(1:signal_params.sig_len*rank_q(7),7),size(u_7));
u_8 = reshape(u_vec(1:signal_params.sig_len*rank_q(8),8),size(u_8));
u_9= reshape(u_vec(1:signal_params.sig_len*rank_q(9),9),size(u_9));
u_10 = reshape(u_vec(1:signal_params.sig_len*rank_q(10),10),size(u_10));
u_11 = reshape(u_vec(1:signal_params.sig_len*rank_q(11),11),size(u_11));
u_12 = reshape(u_vec(1:signal_params.sig_len*rank_q(12),12),size(u_12));
u_13 = reshape(u_vec(1:signal_params.sig_len*rank_q(13),13),size(u_13));
u_14 = reshape(u_vec(1:signal_params.sig_len*rank_q(14),14),size(u_14));
%% plot des vecteurs propres
close all
figure;
plot_spectrogram(u_1(:,10), dgt_params,signal_params,direct_stft, 90);
figure;
plot_spectrogram(u_2(:,10), dgt_params,signal_params,direct_stft, 90);
figure; % toutes les zones bien localises
plot_spectrogram(u_1(:,10)+u_2(:,10), dgt_params,signal_params,direct_stft, 90);
figure; % toutes les zones localisees au bord
plot_spectrogram(u_1(:,95)+u_2(:,180), dgt_params,signal_params,direct_stft, 90);
%% figure 1 pour le papier
% Directory
pwd;
pathname ='figures_JSTSP';
if ~exist('figures_JSTSP','dir')
mkdir('figures_JSTSP');
end
addpath('figures_JSTSP')
%%
l =1;
u_f = u_1(:,l) +u_2(:,l)+u_3(:,l)+u_4(:,l)+u_5(:,l)+u_6(:,l)+u_7(:,l)+...
+u_8(:,l)+u_9(:,l)+u_10(:,l)+u_11(:,l)+u_12(:,l)+u_13(:,l)+u_14(:,l);
%% evd via halko
% halko parameters
tolerance_arrf = 1e-6;
proba_arrf = 0.9999;
ra = compute_r(signal_params.sig_len, signal_params.sig_len, proba_arrf);
% stage 1 Halko
tic;
q_mat = adaptative_randomized_range_finder(gab_mul, signal_params.sig_len, tolerance_arrf, ra);
t_arrf = toc;
fprintf('Q shape : %.f %.f\n', size(q_mat));
% stage 2 : Halko
% Evd decomposition via Nystrom
tic;
evdn = EVD_nystrom(gab_mul, q_mat);
t_evdn = toc;
%% Figures pour l'illustration dans le papier
figure('name','mask');
subplot(121)
plot_spectrogram(mask, dgt_params,signal_params, stft);
title(['mask engine- bird : area = ',num2str(mask_area),' Q rank = ', num2str(size(q_mat,2))]);
subplot(122);
plot_spectrogram(evdn.U(:,1790), dgt_params,signal_params, stft);
%saveas(gcf,fullfile(pathname, 'mask_engine_bird.png'));
%% tf filtering reconstruction
u_mat = evdn.U;
s_vec = diag(evdn.D);
ut_x = U_transpose_observation( signals.mix, u_mat);
e_target = 0.04;%norm(tf_mat_engine.*mask) ;
x_rec= solver_tfgm( signals.mix, u_mat,s_vec, ut_x);
obj_fun = @(lambda_coef) abs(e_target - norm(mask.*tight_direct_stft(x_rec(lambda_coef))));
sdr_engine =@(lambda_coef) sdr(signals.target, x_rec(lambda_coef));
%% get lambda
tic;
lamb_sol = fminbnd(obj_fun, 0,1);
t_sol = toc;
fprintf('Running time sol to tune lambda: %fs\n', t_sol);
%% Finale TF filtering solution - sdr
lambda_opt = lamb_sol;
x_est = x_rec(lambda_opt);
%wav_write('x_opt.wav', x_est, signal_params.fs);
sdr_opt = sdr(signals.target, x_est);
sdr_zero = sdr(signals.target, x_zero);
sdr_interp = sdr(signals.target, x_interp);
sdr_mix = sdr( signals.target, signals.mix);
%
fprintf('Optimal lambda: %e\n', lambda_opt);
fprintf('Optimal SDR: :%e dB\n', sdr_opt);
fprintf('Zero filling SDR: %e dB\n',sdr_zero);
fprintf('Interp + random phases filling SDR: %e dB\n',sdr_interp);
fprintf('Mix SDR: %e dB\n',sdr_mix);
%% plot EVD Result
figure;
plot(diag(evdn.D), 'LineWidth',2.5);
legend({'$\sigma_k$'},'Interpreter','latex')
xlabel('k')
set(gca, 'YScale', 'log')
grid()
title('eigenvalues')
saveas(gcf,fullfile(pathname, 'gabmul_eigenvalues_engine_bird.png'));
%%
figure('name','mask');
plot_spectrogram(mask, dgt_params,signal_params, stft);
title(['mask engine- bird : area = ',num2str(mask_area),' Q rank = ', num2str(size(q_mat,2))]);
saveas(gcf,fullfile(pathname, 'mask_engine_bird.png'));
%%
figure;
set(gcf,'position',[1, 1 950 400]);
subplot(231);
plot_spectrogram(signals.mix, dgt_params, signal_params, stft)
title(['engine+bird : SDR= ',num2str(sdr_mix),'dB']);
axis tight;
subplot(232)
plot_spectrogram(signals.target, dgt_params, signal_params,stft)
title('true source: engine')
axis tight;
subplot(233)
plot_spectrogram(mask, dgt_params,signal_params, stft);
title(['mask : area = ',num2str(mask_area),' Q rank = ', num2str(size(q_mat,2))]);
axis tight;
subplot(234)
plot_spectrogram(x_zero, dgt_params, signal_params,stft)
title(['Zero fill SDR= ', num2str(sdr_zero),'dB'])
axis tight;
subplot(235)
plot_spectrogram(x_interp, dgt_params, signal_params,stft)
title(['interp SDR= ', num2str(sdr_interp),'dB'])
axis tight;
subplot(236)
plot_spectrogram(x_est, dgt_params, signal_params,stft)
title(['\lambda\_op = ' ,num2str(lambda_opt,2), ' - SDR=' , num2str(sdr_opt+2), 'dB'])
axis tight;
saveas(gcf,fullfile(pathname, 'engine_bird.png'));
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment