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

remove folder

parent 1a5662ba
No related branches found
No related tags found
No related merge requests found
clc; clear; close all;
%%
% experiments with
% - wideband sound: engine(Genesis)
% - localized sound: birdsong (Genesis)
%%
make_wav_pairs() % make sure alls pairs are ok
%% Figures Directory
pwd;
pathname ='figures_JSTSP';
if ~exist('figures_JSTSP','dir')
mkdir('figures_JSTSP');
end
addpath('figures_JSTSP')
%%
resampling_fs=8000;
sig_len = 16384;
win_type_list={'hann', 'gauss'};
win_len_list = [256, 128];
mask_params_hann =[4, 0.001, 1];
mask_params_gauss =[3, 0.002, 3];
mask_params = {mask_params_hann, mask_params_gauss};
eigsvect_hann = [45, 18, 1945, 2051];
eigsvect_gauss = [45, 8, 1650, 1438];
eigsvect = {eigsvect_hann, eigsvect_gauss };
e_target_list = [0.04, 0.049];
%%
for k=1:length(win_type_list)
k=2;
%% DGT parameters - gabor frames operators - analysis window
win_type = win_type_list{k};
win_len = win_len_list(k);
fprintf('analysis window : %s\n', win_type);
fprintf('win_len : %.f \n', win_len);
signal_params = generate_signal_parameters(resampling_fs, sig_len);
dgt_params = generate_dgt_parameters(win_type, win_len);
[dgt, idgt] = get_stft_operators(dgt_params, signal_params);
% plot associated window
figure;
plot_win(dgt_params.win, signal_params.fs, signal_params.sig_len, win_type)
title([num2str(win_type), ' - window']);
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, [num2str(win_type), num2str(win_len),'_window.png']));
%% Load signals - get mixtures - their spectrograms
ind_engine = 3;
ind_bird = 5;
deb = 0;
[x_engine, x_bird] = load_pairs(ind_engine, ind_bird, resampling_fs, signal_params.sig_len, deb);
signals = generate_mix_signal(x_engine, x_bird);
% compute dgt
tf_mat_engine = compute_dgt(signals.target, dgt );
tf_mat_bird = compute_dgt(signals.noise, dgt );
tf_mat_mix = compute_dgt(signals.mix, dgt );
% mix -sdr
sdr_mix = sdr(signals.target, signals.mix);
fprintf('The SDR of the mixture is : %e\n', sdr_mix)
%% Generate mask
% mask
alpha= mask_params{k}(1) ;
seuil = mask_params{k}(2);
radius = mask_params{k}(3);
mask = generate_mask(tf_mat_engine, tf_mat_bird, alpha, seuil, radius);
[mask_area, mask_area_ratio] = get_mask_area(mask);
%% Baselines reconstruction
%zero value method
x_zero = solver_tfgm_zero(tf_mat_mix, mask, idgt);
fprintf('Zeros filling SDR is : %e\n', sdr(signals.target, x_zero));
%interpolation + random phases method
x_interp= solver_tfgm_interp(tf_mat_mix, mask, idgt);
fprintf('interp + random phases filling SDR is : %e\n', sdr(signals.target, x_interp));
%% generate Gabor mutliplier
gab_mul = gen_gabmul_operator(dgt, idgt, mask);
%% evd via halko
% halko parameters
tolerance_arrf = 1e-6;
proba_arrf = 0.9999;
r = 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, r);
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;
%% tf filtering reconstruction
u_mat = evdn.U;
s_vec = diag(evdn.D);
ut_x = U_transpose_observation( signals.mix, u_mat);
%%
e_target = e_target_list(k);
x_rec= solver_tfgm( signals.mix, u_mat,s_vec, ut_x);
obj_fun = @(lambda_coef) abs(e_target - norm(mask.*dgt(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 Results
%% Eigenvalue plot - suplots
l1 = eigsvect{k}(1);
l2 = eigsvect{k}(2);
l3 = eigsvect{k}(3);
l4 = eigsvect{k}(4);
figure;
%set(gcf,'position',[1, 1 1100 400]);
%subplot(121);
plot_spectrogram(mask, dgt_params,signal_params, dgt);
set(gca, 'FontSize', 20, 'fontName','Times');
axis square
saveas(gcf,fullfile(pathname, ['mask_cuicui_', num2str(win_type), '.png']));
%subplot(122);
figure;
semilogy(diag(evdn.D), 'Linewidth',3);
hold on;
plot(l1,s_vec(l1),'k-*','Linewidth',3);
plot(l2,s_vec(l2),'m-*','Linewidth',3);
plot(l3,s_vec(l3),'g-*','Linewidth',3);
plot(l4,s_vec(l4),'c-*','Linewidth',3);
grid on;
xlabel('$k$','Interpreter','latex')
legend({'$\sigma_k$',['$\lambda$ =',num2str(l1)],['$\lambda$ =', num2str(l2)],....,
['$\lambda$ =', num2str(l3)], ['$\lambda$ =', num2str(l4)]},...,
'Interpreter','latex','Location','southwest')
axis square
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['eigval_cuicui_', num2str(win_type), '.png']));
%%
figure;
%set(gcf,'position',[1, 1 900 400]);
%subplot(221);
plot_spectrogram(evdn.U(:,l1), dgt_params,signal_params, dgt);
axis square
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['eigvect_', num2str(l1), '.png']));
%subplot(222);
figure;
plot_spectrogram(evdn.U(:,l2), dgt_params,signal_params, dgt);
axis square
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['eigvect_', num2str(l2), '.png']));
%subplot(223);
figure;
plot_spectrogram(evdn.U(:,l3), dgt_params,signal_params, dgt);
axis square
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['eigvect_', num2str(l3), '.png']));
%subplot(224);
figure;
plot_spectrogram(evdn.U(:,l4), dgt_params,signal_params, dgt);
axis square
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['eigvect_', num2str(l4), '.png']));
%%
%figure;
%set(gcf,'position',[1, 1 950 400]);
figure;
plot_spectrogram(signals.noise, dgt_params, signal_params, dgt)
title('perturbation: birdsong' );
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['birdsong_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(signals.mix, dgt_params, signal_params, dgt)
title(['car+birdsong : SDR= ',num2str(sdr_mix),'dB']);
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['mix_birdsong_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(signals.target, dgt_params, signal_params,dgt)
title('true source: car')
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['car_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(mask, dgt_params,signal_params, dgt);
title(['mask car +birdsong : area = ',num2str(mask_area)]);
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['mask_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(x_zero, dgt_params, signal_params,dgt)
title(['Zero fill SDR= ', num2str(sdr_zero),'dB'])
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['x_zero_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(x_interp, dgt_params, signal_params,dgt)
title(['interp SDR= ', num2str(sdr_interp),'dB'])
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['interp_' , num2str(win_type),'_' , num2str(win_len),'.png']));
figure;
plot_spectrogram(x_est, dgt_params, signal_params,dgt)
title(['\lambda\_op = ' ,num2str(lambda_opt,2), ' - SDR=' , num2str(sdr_opt), 'dB'])
set(gca, 'FontSize', 20, 'fontName','Times');
saveas(gcf,fullfile(pathname, ['x_est_' , num2str(win_type),'_' , num2str(win_len),'.png']));
%%
file_name = ['exp_engine_cuicui_1area_', num2str(win_type),'.mat'];
save(file_name,'signal_params','dgt_params', 'signals','mask','mask_area',....,
'gab_mul', 'q_mat', 'evdn','sdr_mix', 'sdr_opt','sdr_zero','sdr_interp','x_zero','x_interp','x_est','lambda_opt','e_target')
end
clc; clear; close all;
%
%%
pwd;
pathname ='figures_JSTSP';
if ~exist('figures_JSTSP','dir')
mkdir('figures_JSTSP');
end
addpath('figures_JSTSP')
%%
ind_loc = 5;
ind_wd = 3;
deb_ind_loc = 0;
deb_ind_wd=0;
resampling_fs = 8000;
sig_len = 16384;
%% DGT params - signals - mask
clc;
param_gauss = get_win_gauss_param();
win_len = param_gauss.win_len;
win_type = param_gauss.win_type;
alpha = param_gauss.alpha;
seuil = param_gauss.seuil;
radius = param_gauss.radius;
approx_win_len = 128;
hop =32;
nbins=512;
[signals, dgt_params, signal_params, mask, dgt,idgt] = get_mix(ind_loc, ...,
ind_wd, deb_ind_loc, deb_ind_wd, resampling_fs, sig_len,...,
approx_win_len,hop, nbins, win_type, alpha, seuil, radius);
[mask_area, mask_area_ratio] = get_mask_area(mask);
fprintf("We work with %s window of length %.f\n", win_type, win_len);
fprintf("Gabor transform parameters are: \n")
fprintf('hop :%2.f\n', dgt_params.hop);
fprintf('n_bins: %2.f\n', dgt_params.nbins);
fprintf("The parameters for smoothing the mask are: \n")
fprintf("alpha = %f\n", alpha);
fprintf("seuil = %f\n", seuil);
fprintf("radius = %f\n", radius);
% plot associated window
figure;
plot_win(dgt_params.win, signal_params.fs, signal_params.sig_len, win_type)
title([num2str(win_type), ' - window']);
set(gca, 'FontSize', 20, 'fontName','Times');
%saveas(gcf,fullfile(pathname, [num2str(win_type),'_window.png']));
%% Spectrogrammes - Mask
tf_mat_engine = compute_dgt(signals.target, dgt);
tf_mat_bird = compute_dgt(signals.noise, dgt);
tf_mat_mix = compute_dgt(signals.mix, dgt);
%plot_spectrogram(x, dgt_params, signal_params, dgt, dynrange, clim)
dynrange=90;
figure('name','engine'); plot_spectrogram(tf_mat_engine, dgt_params, signal_params, dgt);
title('Source')
set(gca, 'FontSize', 20, 'fontName','Times');
figure('name','bird'); plot_spectrogram(tf_mat_bird, dgt_params, signal_params, dgt);
title('perturbation')
set(gca, 'FontSize', 20, 'fontName','Times');
figure('name','mix');
plot_spectrogram(tf_mat_mix, dgt_params,signal_params, dgt);
set(gca, 'FontSize', 20, 'fontName','Times');
%title('mix')
%saveas(gcf,fullfile(pathname, 'engine_bird_mixture.png'));
%plot mask
figure('name','mask'); plot_spectrogram(mask, dgt_params,signal_params, dgt );
%title(['mask : mask-area = ',num2str(mask_area)]);
set(gca, 'FontSize', 20, 'fontName','Times');
%saveas(gcf,fullfile(pathname, 'mask_cuicui_gauss.png'));
% figure;
% plot_spectrogram((1-mask).*tf_mat_bird, dgt_params,signal_params, dgt);
fprintf("mask area is: %f\n",mask_area);
%% Mix - SDR
sdr_mix = sdr(signals.target, signals.mix);
fprintf('The SDR of the mixture is : %e\n', sdr_mix)
%%
%[mask_limast,mask_labeled, mask_area_list] = make_subregions(mask, dgt_params, signal_params);
%pq_norms_val = pq_norms(sig_len,dgt,idgt,mask_list);
%pq_norms_val1 = get_pq_norms(sig_len, dgt, idgt, mask_labeled);
%pq_norms = get_pq_norms(sig_len, dgt, idgt, mask_list);
%%
tol=10^(-3);
[pq_norms_val, mask_labeled] = create_subregions(mask, dgt, idgt, ...,
dgt_params, signal_params, tol);
%%
final_mask_labeled = mask_labeled;
[gabmul_list, mask_list] = get_P_gabmul(final_mask_labeled, dgt, idgt);
%%
x_mix = signals.mix;
tolerance_arrf = 10^(-3);
proba_arrf = 0.999;
[t_arrf,t_evdn, t_ut_x, rank_q, s_vec_list, u_mat_list,...,
ut_x_list,r] = compute_decomposition(x_mix, mask_list, gabmul_list,...,
tolerance_arrf, proba_arrf);
%%
lambda_coef=0.01;
x=compute_estimate(lambda_coef, x_mix, s_vec_list, u_mat_list, ut_x_list);
%%
[lambda_est, t_est] = compute_lambda(x_mix, mask, dgt_params,...,
signal_params, dgt, s_vec_list, u_mat_list, ut_x_list,...,
gabmul_list);
%% figure
figure;
imagesc(real(log10(pq_norms_val)))
ylabel('$p$','Interpreter','latex')
xlabel('$q$', 'Interpreter', 'latex')
colorbar()
set(gca, 'FontSize', 20, 'fontName','Times');
%title('Final norms of Gabor multiplier composition')
%saveas(gcf,fullfile(pathname, 'norm_mulpq.png'));
%% fixed
gabmul_list = get_P_gabmul(mask_labels, dgt, idgt);
%%
rank = 10;
x_mix = signals.mix;
[t_rrf,t_evdn, t_ut_x, s_vec_list, u_mat_list,...,
ut_x_list]=compute_decomposition_fixed_rank(x_mix, mask_labels, gabmul_list, rank);
%%
tolerance_arrf = 0.1;
proba_arrf = 0.9;
[mask_labels, mask_area_list,n_labels] = make_subregions(mask, dgt_params, signal_params);
[t_arrf,t_evdn, t_ut_x, rank_q, s_vec_list, u_mat_list,...,
ut_x_list,r] = compute_decomposition(x_mix, mask_labels, gabmul_list,...,
tolerance_arrf, proba_arrf);
%%
n_areas = length(mask_area);
lambda_coef = ones(n_areas,1);
x=compute_estimate(lambda_coef, x_mix, n_areas,s_vec_list, u_mat_list, ut_x_list);
%%
x_target = signals.target;
[lambda_oracle, t_oracle] = compute_lambda_oracle_sdr(x_target, x_mix,...,
n_areas,s_vec_list, u_mat_list, ut_x_list);
%%
[lambda_est, t_est] = compute_lambda(x_mix, mask, dgt_params,...,
signal_params, dgt, mask_labels, s_vec_list, u_mat_list, ut_x_list,...,
gabmul_list);
%% parametres de la EVD halko
tolerance_arrf = 1e-6;
proba_arrf = 1 - 1e-4;
x_mix = signals.mix;
compute_decomposition(x_mix,mask_labels, gabmul_list,tolerance_arrf, proba_arrf);
%% create subregions
[mask_labels, mask_area_list,n_labels] = make_subregions(mask, dgt_params, signal_params);
%%
e_target = zeros(n_labels, 1);
for k=1:n_labels
e_target_k = norm(mask_labels{k}.*tf_mat_engine) ;
e_target(k) = e_target_k;
end
%%
x_mix = signals.mix;
x_target = signals.target;
[x_rec, t_arrf, t_evdn, t_ut_x, rank_q, s_vec_list, u_mat_list,...,
ut_x_list, lambda_vec_opt] = filtering_out_Pareas(x_mix, mask_labels, dgt, idgt, x_target, tolerance_arrf, proba_arrf);
%%
i_p=1;
pq_norms_val = update_pq_norms(mask_labels, pq_norms, i_p, signal_params, dgt, idgt);
%%
i_p = 2;
i_q=1;
[ mask, pq_norms_val] = merge_subregions(mask, pq_norms_val, i_p, i_q);
%%
figure;
plot_spectrogram(x_rec, dgt_params, signal_params,dgt);
%%
%%
% %% Generate list of Gabor multipliers
%
% gab_mul_list = get_P_gabmul(mask_labels, dgt, idgt);
%
% %% EVD decomposition
%
% [t_arrf, t_evdn, t_uh_x, s_vec_list, u_mat_list, t_uh_x_list, ...,
% rank_q] = compute_decomposition(mask_labels, gab_mul_list, dgt_params, signal_params, signals.mix, tolerance_arrf,r);
%
% %% Tuning Lambda
% n_areas = n_labels;
% uh_x_list = t_uh_x_list;
%
%
%%
%%
% %%
% %obj_fun = @(lambda_vec) norm(signals.target - compute_estimate(lambda_vec, s_vec_list, u_mat_list, uh_x_list,n_areas, signals.mix));
% obj_fun = @(lambda_vec) abs(e_target -norm( mask.*dgt(compute_estimate(lambda_vec, s_vec_list, u_mat_list, uh_x_list,n_areas, signals.mix))));
%
% %% Generate and save msk for each regions
% x0= ones(n_labels,1);
% tic;
% sol = fmincon(obj_fun,x0);
% t1 =toc;
%
% %%
% lambda_opt = obj_fun(x0);
% x_est = compute_estimate(lambda_opt, s_vec_list, u_mat_list, uh_x_list, n_areas,signals.mix);
%
% figure;
% plot_spectrogram(x_est, dgt_params, signal_params,dgt);
% %%
% %
% % all_mask = zeros(size(mask,1)*size(mask,2),n_labels);
% %
% %
% % for k =1:n_labels
% % [mask_label,~] = bwlabel(mask);
% %
% % % on construit chaque mask
% % mask_label(mask_label~=k)=0;
% % mask_ =mask_label;
% % all_mask(:,k) = mask_(:);
% % 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_);
% %
% % end
%
%
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment