clear all;
close all;
clc;

addpath(genpath(pwd));

%% Load the data paths
fid = fopen('_file_list.txt', 'r');
pth = textscan(fid, '%s');
pth = pth{1, 1};
fclose(fid);

clear fid;

%% Load UPDRS III
[~, ~, clin] = xlsread('info.xlsx', 'paraclinical_info');

%% Perform a parametrization
% Each features that has a form o vector or a matrix is post-processed:
% - 1. quartile
% - 2. quartile (a.k.a. median)
% - 3. quartile
% - 4. mean
% - 5. std
data = struct();

i = 0;
for dat = 1:1:length(pth)
    
    % Load the recording
    wav_pth = pth{dat, 1};
    [x, fs] = audioread(wav_pth);
    
    [path, name, ext] = fileparts(wav_pth);
    disp([num2str(dat) '. Processing: ' name])
    
    % Check the recording's validity
    % <recording is only valid if it is recorded from a PD patient,
    %  and the patient has an UPDRS III value atached>
    is_valid = 0;
    
    for j = 1:size(clin, 1)
        if ~isempty(strfind(name, clin{j, 1}))
            is_valid = 1;
            break;
        end
    end
    
    % Continue only for a valid recording
    if is_valid
        i = i + 1;
        
        % Save the paraclinical data
        data(i).pid = name;
        data(i).gender = find_gender(name);
        data(i).age = clin{j, 2};
        data(i).pd_duration = clin{j, 3};
        data(i).updrs_3 = clin{j, 4};
        
        % Compute the speech features
        % - Short-time energy
        % - Zero-crossing rate
        % - Formants and bandidths (praat!
        % - F0 (praat)
        % - jitter variants
        % - shimmer variants
        % - hnr/nhr
        % - LPC coeff.
        % - MFCC coeff.
        % - PLP coeff.
        praat_str = praat_info(x, fs); 

        zcr = zero_crossing_rate(x, hamming(0.02*fs), 0.01*fs); 
        data(i).zcr_q1 = prctile(zcr, 25);
        data(i).zcr_median = prctile(zcr, 50);
        data(i).zcr_q3 = prctile(zcr, 75);
        data(i).zcr_mean = mean(zcr);
        data(i).zcr_std = std(zcr);

        ste = energy(x, hamming(0.02*fs), 0.01*fs);
        data(i).ste_q1 = prctile(ste, 25);
        data(i).ste_median = prctile(ste, 50);
        data(i).ste_q3 = prctile(ste, 75);
        data(i).ste_mean = mean(ste);
        data(i).ste_std = std(ste);

        formants = praat_formants(x, fs);
        data(i).F1_q1 = prctile(formants{1, 1}, 25);
        data(i).F1_median = prctile(formants{1, 1}, 50);
        data(i).F1_q3 = prctile(formants{1, 1}, 75);
        data(i).F1_mean = mean(formants{1, 1});
        data(i).F1_std = std(formants{1, 1});
        data(i).BW1_q1 = prctile(formants{1, 2}, 25);
        data(i).BW1_median = prctile(formants{1, 2}, 50);
        data(i).BW1_q3 = prctile(formants{1, 2}, 75);
        data(i).BW1_mean = mean(formants{1, 2});
        data(i).BW1_std = std(formants{1, 2});

        data(i).F2_q1 = prctile(formants{1, 3}, 25);
        data(i).F2_median = prctile(formants{1, 3}, 50);
        data(i).F2_q3 = prctile(formants{1, 3}, 75);
        data(i).F2_mean = mean(formants{1, 3});
        data(i).F2_std = std(formants{1, 3});
        data(i).BW2_q1 = prctile(formants{1, 4}, 25);
        data(i).BW2_median = prctile(formants{1, 4}, 50);
        data(i).BW2_q3 = prctile(formants{1, 4}, 75);
        data(i).BW2_mean = mean(formants{1, 4});
        data(i).BW2_std = std(formants{1, 4});

        data(i).F3_q1 = prctile(formants{1, 5}, 25);
        data(i).F3_median = prctile(formants{1, 5}, 50);
        data(i).F3_q3 = prctile(formants{1, 5}, 75);
        data(i).F3_mean = mean(formants{1, 5});
        data(i).F3_std = std(formants{1, 5});
        data(i).BW3_q1 = prctile(formants{1, 6}, 25);
        data(i).BW3_median = prctile(formants{1, 6}, 50);
        data(i).BW3_q3 = prctile(formants{1, 6}, 75);
        data(i).BW3_mean = mean(formants{1, 6});
        data(i).BW3_std = std(formants{1, 6});

        F0 = praat_F0(x, fs);
        data(i).F0_q1 = prctile(F0, 25);
        data(i).F0_median = prctile(F0, 50);
        data(i).F0_q3 = prctile(F0, 75);
        data(i).F0_mean = mean(F0);
        data(i).F0_std = std(F0); 

        data(i).jitter_local = praat_str.jitter_local;
        data(i).jitter_localabs = praat_str.jitter_localabs;
        data(i).jitter_rap = praat_str.jitter_rap;
        data(i).jitter_ppq5 = praat_str.jitter_ppq5;
        data(i).jitter_ddp = praat_str.jitter_ddp;
        data(i).shimm_local = praat_str.shimm_local;
        data(i).shimm_localdb = praat_str.shimm_localdb;
        data(i).shimm_apq3 = praat_str.shimm_apq3;
        data(i).shimm_apq5 = praat_str.shimm_apq5;
        data(i).shimm_apq11 = praat_str.shimm_apq11;
        data(i).shimm_dda = praat_str.shimm_dda;

        data(i).hnr_q1 = prctile(praat_str.hnr_hn, 25);
        data(i).hnr_median = prctile(praat_str.hnr_hn, 50);
        data(i).hnr_q3 = prctile(praat_str.hnr_hn, 75);
        data(i).hnr_mean = mean(praat_str.hnr_hn);
        data(i).hnr_std = std(praat_str.hnr_hn); 

        data(i).nhr_q1 = prctile(praat_str.hnr_nh, 25);
        data(i).nhr_median = prctile(praat_str.hnr_nh, 50);
        data(i).nhr_q3 = prctile(praat_str.hnr_nh, 75);
        data(i).nhr_mean = mean(praat_str.hnr_nh);
        data(i).nhr_std = std(praat_str.hnr_nh);
    end
end

%% Save the results
save('database_info.mat', 'data');
