p% JBA June 2011, after code by others (Roger, Woojae, Riya?)
%
% Stacked error barplots (as a funscion of SNRs) and confusion pattern for HI ears
% CM is recreated after taking off the bad utterances (24/96*100 = 25%) from Exps. I and III.

% There are 4 data conditions, which are
% Exp.I: 45 ears; Retest of Exp.I: 6 ears; Exp.II: 17 ears; Exp.III(NAL-R): 43 ears
% help on: string = input('hi there','s') which returns a string

%clear all %this will kill the debugging code
close all

% 0 for all SNRS: 1,...,SNR_COUNT, for all plots and entropy calculations
snr_i = 0; %do all   
%snr_i =2;  %e.g.: only do 12 dB SNR
%snr_i=1:2; %do Q and 12

addpath('bin');

SUBJ_I ={'AS-R-I','AW-R-I','BD-L-I','BD-R-I','BG-L-I','BG-R-I','CL-L-I','CL-R-I','CP-L-I','CP-R-I'...
'DG-L-I','DG-R-I','DN-L-I','DN-R-I','EG-L-I','EG-R-I','EM-L-I','HV-L-I','HV-R-I','JG-L-I','JG-R-I'...
'JS-L-I','JS-R-I','JT-L-I','KW-R-I','LB-L-I','LB-R-I','LW-L-I','LW-R-I','MA-L-I','MA-R-I','MB-L-I','MB-R-I'...
'MC-L-I','MC-R-I','MJ-L-I','MJ-R-I','PB-L-I','PB-R-I','SN-L-I','SN-R-I','TB-L-I','TB-R-I','VS-L-I','VS-R-I'};

SUBJ_IR1={'HV-L-retest','HV-R-retest','LW-L-retest','LW-R-retest','TB-L-retest','TB-R-retest'};

SUBJ_IR2={'HV-L-retest2','LW-L-retest2','TB-L-retest2'};

%Exp II in sorted order of overall-all confusions: Best to worst
SUBJ_II={ ...   
 '44L/BG-L-II',...   % 4
 '36L/JG-L-II',...   % 4
 '44R/BG-R-II',...   % 4
 '36R/JG-R-II',...   % 4
 '40L/PB-L-II',...   % 4
 '40R/PB-R-II',...	% 4
 '32R/DN-R-II',...	% 4
 '32L/DN-L-II',...	% 4
 '34L/VS-L-II',...	% 4
 '34R/VS-R-II',...	% 4
 '46L/EG-L-II',...	% 4
 '46R/EG-R-II',...	% 4
 '14R/AW-R-II',...	% 4
 '30L/DG-L-II',...	% 4
 '30R/DG-R-II',...	% 4
 '01L/AS-L-II',...	% 4
 '01R/AS-R-II'...	% 4 %This one has an error: size is 12x10
};

SNR_LIST = {'Q', '12', '6', '0', '-6', '-12'};
%define an snr_i for each of the 17 subjects
% in the plot order (the above list)
%map from DB order to presented order
%SUBJ_SNR = {1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1};
%SUBJ_SNR = {2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2};
SUBJ_SNR = {3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3};
%SUBJ_SNR = {4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4};

%SUBJ_CV = {'pa'};

SUBJ_III={'AW-R-III','BD-L-III','BD-R-III','BG-L-III','BG-R-III','CL-L-III','CL-R-III','CP-L-III','CP-R-III'...
'DG-L-III','DG-R-III','DN-L-III','DN-R-III','EG-L-III','EG-R-III','EM-L-III','HV-L-III','HV-R-III','JG-L-III','JG-R-III'...
'JS-L-III','JS-R-III','JT-L-III','KW-R-III','LW-L-III','LW-R-III','MA-L-III','MA-R-III','MB-L-III','MB-R-III'...
'MC-L-III','MC-R-III','MJ-L-III','MJ-R-III','PB-L-III','PB-R-III','SN-L-III','SN-R-III','TB-L-III','TB-R-III','VS-L-III','VS-R-III'};

STIM_LIST={'pa','ta','ka','fa','xta','sa','xsa','ba','da','ga','va','xda','za','xza','ma','na'};
% order     1    2    3    4     5    6     7    8    9    10   11    12   13   14    15   16    

RESP_LIST = {'pa','ta','ka','fa','xta','sa','xsa','ba','da','ga','va','xda','za','xza','ma','na', '?'};

RESP_LIST_D = {'pa','ta','ka','fa','Ta','sa','Sa','ba','da','ga','va','Da','za','Za','ma','na', '?'};
STIM_LIST_D = {'pa','ta','ka','fa','Ta','sa','Sa','ba','da','ga','va','Da','za','Za','ma','na'};

%You can uncomment the following two lines to manually specify the file
%subject = 'MJ-R'; 

SNR_COUNT = 4; %Good for Exp II

for subj=1:length(SUBJ_II)   %loop over subjects
%k=strfind(F(subj).name,'.')-1 %find location of "."-1
close all;


%% start ExpII
%subject = 'AS-R-II';
subject = SUBJ_II{subj};
%disp([subject()]);
disp([subject]);

%subject = 'AS-R-I' for Exp.I (non NAL-R)
%subject = 'HV-L-retest' for retest for Exp.I (Same stimuli)
%subject = 'AS-R-II' for Exp.II
%subject = 'LW-R-III' for Exp.III (NAL-R)

%JBA: Load 6 sets of CMs and the corresponding SNRs
[CM, SNR] = r_loadcm(subject(5:end)); %remove the "ID/"

%Remove rows and cols given zero presentations
rowsum=sum(CM{1},2);
for i=2:SNR_COUNT  %JBA: i=1 was done on the previous line
  rowsum=rowsum+sum(CM{i},2);
end %for i=2:SNR_COUNT
KEEP=find(rowsum>0)'; %dont keep when all the rowsums are zero

%JBA: We need a loop over SNRs, but there doesnt seem to be one
%JBA: This sets the SNR for the entropy + stack'ed bar plot
for snr_i = 1:SNR_COUNT; %%%SUBJ_SNR{subj}; %maping from subj to snr_i

%Stacked-Errors Bar Plots

%For the next few lines we will compute the total entropy
CM_total = zeros(16,17); % collapse 5 CM's into one

%JBA: snr_i is a flag to indicate an average over all SNRs
if snr_i == 0
    for i=1:SNR_COUNT;
        CM_total = CM_total + CM{i};
    end
else
    CM_total = CM{snr_i};
end %if snr_i == 0

if (0)
    % hack to distribute ? responses uniformly
    for i=1:16
        q = CM_total(i, 17) / 16;
        CM_total(i,17) = 0;
        CM_total(i,1:16) = CM_total(i,1:16) + q;
    end
end

% Finally, compute row entropy
N = sum(CM_total, 2);
Nrep = repmat(N, 1, 17);   % create a normalization matrix for ./

% only keep those consonants that had presentations
%This is now done above, based on row-sums
% tmp = find(N == 0);
% KEEP = [1:16];  % keep all initially
% for i=1:length(tmp);
%     KEEP(find(KEEP == tmp(i))) = '';   % delete if not presented
% end

%Note that the entropy is  computed including the last 17th row, of '?'
Prow = CM_total ./ Nrep; 
Hrow = sum(log2(Prow.^-Prow), 2);  % contains entropy for each row

%Hrow(find(isnan(Hrow))) = 0;

%
fprintf('2^{Entropy} =');
fprintf(' %2.1f', 2.^Hrow); fprintf('\n');

% compute overall response entropy
%Pall = CM_total / sum(sum(CM_total));   % P 
%Hall = sum(sum(log2(Pall.^-Pall)));


tmp = Hrow(KEEP);
tmp = tmp(find(~isnan(tmp)));
Havg = mean(tmp);  % only compute entropy average if presentation > 0.

% make a bar graph
Pe_row = zeros(size(Prow));  %last row for question mark response

pe_order = zeros(16,2);
for i=1:16
    pe = 1 - Prow(i,i);
    pe_order(i,:) = [pe, i];
    Pe_row(i,:) = Prow(i,:);
    Pe_row(i,i) = 0;
end

a = find(isnan(pe_order(:,1)));
pe_order(a, 1) = zeros(length(a), 1);   % remove nan's (xta, xda)

pe_order = pe_order(KEEP, :);

pe_order = sortrows(pe_order);  % comment out to disable sorting of responses
s_order = pe_order(:,2);

fig(1);
bar(Pe_row(s_order, :)*100, 'stack')

x = 1:length(Pe_row(s_order));
p = sum(Pe_row(s_order,:), 2);
e = sqrt(p.*(1-p)./N(s_order));
hold on
errorbar(x-0.45,p*100,e*100, '.', 'Color', [0,0,0])
hold off


set(gca,'XTick',1:1:length(s_order),'XTickLabel',STIM_LIST_D(s_order),'Fontsize',16);
%set(gca,'XDir','reverse');
xlim([0, length(s_order)] + 0.5);
ylim([0, 120])
set(gca,'YTick',0:20:120,'YTickLabel',{0, 20, 40, 60, 80, 100, ''},'Fontsize',16);
xlabel('Ranked-Order Consonants as Stimuli','Fontsize',18);
ylabel('Sum of Confusions in P_e [%]','Fontsize',18);

legend(STIM_LIST_D,'Location','NorthEastoutside','Fontsize',14);

%Fix the title to remove the '/'
%subject(3)='-';

if snr_i == 0
    title(subject);
else
    title(sprintf('%s at %s dB-SNR', subject, SNR_LIST{snr_i}),'Fontsize',18);
end

for x=1:16;
    x2 = find(KEEP == x); % this translates the index in the original data set to the index in the reduced data set
    if isempty(x2)
        continue;
    end
    tok = s_order(x2);  % index translation, x along axis corresonds to a different token due to sorting
    c = cumsum(Pe_row(tok,:));
    large_error = find(Pe_row(tok,:) > 0.10);
    ypos = c(large_error);
    for j=1:length(large_error)
        y = ypos(j);
        txt = sprintf('%2.0f%%', 100*Pe_row(tok, large_error(j)));
        col = [0,0,0];
        if any([1,2,3,4, 15,16 17] - large_error(j) == 0)
            col = [1,1,1];
        end
        
        text(x2,y*100 - 10, txt, 'rotation', 90, 'color', col); 
       
    end
    if ~isnan(Hrow(tok))
        txt = sprintf('%1.3f', Hrow(tok));
    else
        txt = 'NO PRESENT';
    end
    if isnan(c(end))
        c(end) = 0;
    end
    
    text(x2, c(end)*100 + 4, txt, 'rotation', 90);
 end


txt = sprintf('H_{avg} = %2.3f', Havg);
text(1.5, 110, txt, 'FontSize', 16);

%title('BD-R: Exp.I', 'FontSize',22);
if snr_i > 0
    snr_txt = sprintf('_SNR_%s', SNR_LIST{snr_i});
else
    snr_txt = '';
end
%fname = sprintf('stackedbar_%s%s.eps', subject, snr_txt); %put the SNR in the filename
subject(4)='_';
fname = sprintf('stackedbar_%s_SNR_%s.eps', subject, SNR_LIST{snr_i}); %name with no SNR
print(gcf, fullfile('graphics', fname), '-depsc2');
end % for snr_i = 1:SNR_COUNT;


%%
% Plotting Individual CVs Confusion Pattern
% the comparison @ between flat gain(even#) and nalr(odd#)

for token=1:16
%token = 2;  % Adjust this value between 1 and 16 for each stimulus
%STIM_LIST = {'pa','ta','ka','fa','xta','sa','xsa','ba','da','ga','va','xda','za','xza','ma','na'};
% order        1    2    3    4     5    6     7    8    9    10   11    12   13   14    15   16   

resp = zeros(SNR_COUNT, 17);
for i=1:SNR_COUNT
    resp(i,:) = CM{i}(token, :);
end

% normalize responses so that row sums equal 1
resp_n = resp ./ repmat(sum(resp, 2), 1, 17);

%plotting 1 by 2
%subplot(221) %**1=left; **2=right
snr = [18, 12, 6, 0, -6, -12]';
snr = snr(1:SNR_COUNT);

eb = zeros(SNR_COUNT,1);

fig(2);
C = get(gca(), 'ColorOrder');
for i=1:17
    ci = mod((i-1), length(C)) + 1;
    if (i ~= token)
        semilogy(snr, 100*resp_n(:,i), '-*', 'color', C(ci,:))
    else
        errorbar(snr, 100*resp_n(:,i), eb, '-*', 'color', C(ci,:))
    end
    hold on
end

% display N's
resp_sum = sum(resp, 2);
for i=1:4
    text(snr(i)+1, 130, sprintf('%i', resp_sum(i)));
end

set(gca(), 'YScale', 'log');

ylim([5,140]);
set(gca, 'XTick',[-18:6:18],'XTickLabel',{'-18';'-12';'-6';'0';'6';'12';'Q'},'YTick',[0:5:10:20:100],'Fontsize',16);
% make x axis pretty with SNR values
xlim([-15, 20]); 
grid on
%title(sprintf('Confusions for /%s/ in %s', STIM_LIST{token},subject),'Fontsize',18)
%legend(sprintf('stim: /%s/', STIM_LIST{token}),'Location','SouthWest');
txt = sprintf('stim: /%s/', STIM_LIST_D{token});
text(-14, 107, txt, 'Fontsize', 16)
xlabel('Signal-to-Noise Ratio [dB]','Fontsize',18)
ylabel('Percent Correct [%]','Fontsize',18)


%x1 = snr;
%y1 = 100-100*wonalr;
%e1 = std(y1)*ones(size(x1));
%errorbar(x1,y1,e1,'xb','LineWidth',2);


% label each point on the CM
for i=1:SNR_COUNT
    for j=1:17
        ci = mod((j-1), length(C)) + 1;
        p = 100*resp_n(i,j);
        text(snr(i), p * (1.05^(ci/2)), RESP_LIST_D{j}, 'color', C(ci, :), 'FontSize', 16);
    end
end
subject(4)='_';
fname = sprintf('cp%s_%s.eps', subject, RESP_LIST{token});
lab=subject;lab(4)='/';
title(lab,'FontSize',18);%clear lab;
print(gcf, fullfile('graphics', fname), '-depsc2');
close(gcf);
end   %for token=1:17   %Loop over tokens to build CP(snr)
%
%%

% We already have CM{}: Reloading CM is not necessary (& thus not desireable). 
%CM = r_loadcm(subject, 'cm_data_filter_rds');
% CM = r_loadcm(subject); 'cm_data_filter_rds' is the default

fid = 1;  % to output to console
%fid = fopen(sprintf('CSV-%s.CSV', fname), 'w');

%clc

fprintf(fid, 'LISTENER: %s, \n', subject);

for snr_i=1:SNR_COUNT  %for production code
%for snr_i=3 %for debugging code
  
PCM = CM{snr_i}; %make a working copy of CM

Sum=sum(PCM,2);  %Keep true row sum
L=Sum~=0;  %Logical variable: false for zero rowsum
Lc=[L;true]; %Need col = 17 than row to hold the '?' responses
PCM=PCM(L,Lc); %remove rows and columns for sounds never presented
SL=STIM_LIST_D(L); %prune STIM LIST removing, e.g. T,D
RL=RESP_LIST_D(Lc); %prune RESP LIST
Sum=Sum(L);  %pruned true rowsums 
clear L Lc; %we are done with these

%JBA: fixed the problem of removing rows and col if diag <3 (June 3, 2011)
D=diag(PCM); %Prepare to remove PCM<3, but save the Diagional!
L=PCM<3;
%dont remove anything on the diagional:
I=eye(size(L)); %1 on the diagional
L(I==1)=0; %L=1 where we will zero elements <3
PCM(L)=0; %zero out all elements <3 (e.g.: 0, 1, 2) except the DIAG
clear L; %done with this

% The goal is to display the pruned CM: 
% 1) Remove all rows with no (less than 3) errors
% 2) sort the columns that appear as rows (and only these columns)
% 3) Not now: we may wish to sort the error along the diagional

% Step 1) remove rows with no errors:
% Method 1: compare col sum to row sum: If equal remove that row
% Method 2: remove rows where Diag == row sum (D==S)
S=sum(PCM,2); %recompute the RowSum of modified CM (this will differ from Sum due to PCM<3
%Find rows where the Diag~=Rowsum => WITH errors
Lr= S~=D;
%if ~isempty(L) %I removed find() so its now logical => never empty
  % Remove rows; We keep the cols due to possible misheard sounds
PCM=PCM(Lr,:); %Prune perfect rows (no errors) (keep columns)
%if ~isempty(PCM); %dont continue if PCM is empty
Sum=Sum(Lr,:); %Prune sums of perfect rows (keep columns)
SL=SL(Lr);  %Remove zero error Stimuli (i.e., row) labels
clear Lr; %done pruning the rows (arrays and labels)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%We are done pruning rows. Time to move on to the columns
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% This is only required because of a bug in the order of doing things.
% It should be unnecessary for working code
if length(SL)>length(RL);  %Look for "impossible" case #rows > #cols
  disp('Error with dimensions of matrix');
  disp(['STIM:' sprintf('%s ',SL{:} )]);
  disp(['RESP:' sprintf('%s ',RL{:} )]);
  disp('Trapped on line 399 due to logic flow in code.');
  pause
 %k=length(SL);i=length(RL);
 %disp(sprintf('\n|SL|=%i, |RL|=%i',k,i))
 %disp(sprintf('RL{k} does not exist:  %i,%s,%i,%s,%s',k,SL{k},i,RL{i},subject));
 %XCM = PCM(:,L);  %This is an error condition
end %if length(SL)>...
  
  
if ~isempty(PCM); %dont continue if PCM is empty
  
% Step 2) Sort the columns that appear as rows
% I.E.: REORDER: force CM to be diagonal at the front

%Lc=1:length(PCM); %Lc is the col index of PCM

K=length(SL); I=length(RL); Lc=1:I;
for k=1:K; %Loop over rows
   for i=1:I %Loop over cols
  if RL{i}==SL{k}; %swap cols
  TMP = RL{k};  RL{k}=RL{i};  RL{i}=TMP;
  TMP = Lc(k);   Lc(k) =Lc(i);   Lc(i) =TMP;
  end %if RL{i}... 
 end %i=1:length(RL)
 end %for k=1:K
 
 %From here refer to PCM(:,Lc) to have the right order
 PCM=PCM(:,Lc);
 clear Lc;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Sort done: CM is fat with K rows and columns aligned as a sqare CM
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%This is the correct place to sort the rows in terms of decreasing error

%Step 3) remove columns >K having zero errors
S=sum(PCM,1); %compute column sums, to find and remove remove unwanted columns

%Remove cols(K+1:end) with no errors
% This is where we go wrong!!! Restrict this operation to cols > K
Lc=S>0; 
Lc(1:K)=true; %Never !!! remove the first K columns
%L=find(L); %Use logical (NOT array indexing)
RL=RL(Lc);  %Not sure this is correct (as with SL above)

%dont need to do this, as the Sum is passed in the argument of cm2tex
%PCM(:,end+1)=Sum; %append the "true" row sums, after '?'

%format the CM for LaTex output and write it in './tex/...dB.txt' file
name = ['tex/',subject,['@'],SNR_LIST{snr_i},'dB.txt'];
cm2tex([PCM(:,Lc),Sum],name,SL,[RL,'\Sigma'], SNR_LIST(snr_i))
clear Lc; %done with this

end %if is empty(L)
end %for snr_i =
end %for subj=


% if (fid > 2)
%     fclose(fid);
%end %fid>2

%STIM_LIST = {'pa','ta','ka','fa','xta','sa','xsa','ba','da','ga','va','xda','za','xza','ma','na'};
% order        1    2    3    4     5    6     7    8    9    10   11   12    13   14    15   16    

%print(gcf,'-depsc2',['CP_AS-L_pa']);


%%
disp(['DONE (normal end)'])
