
function [theResult, duplicate_O] = average_and_remove_duplicates_oxy (output, outfile)
%
% Find duplicate samples in output array
%
% output = seven column array of... 
%   [output_sta output_niskin output_bottle output_twiceratio output_salinity output_quality output_cast];
%
% this routine looks for duplicate output_niskin for each output_sta
% Averages these duplicates together provided they are less than 2*stddev
% of all duplicates away from the mean duplicate difference (i.e. tries to remove outliers
% before averaging).
%
% Output:
% outfile = text file of seven columns as in output array, but with
% duplicates removed.
% Averaged duplicates qC flags are replaced with a "6"
%
% version 1 called average_and_remove_duplicates for clivar with columns of
% output different
% version 2 called average_and_remove_duplicates_aoml for the format as
% described above.
%
Icolumn_O = 7;
Icolumn_Omml = 6;
Icolumn_thio = 5;
Icolumn_quality = 8;

[pathstr,this_name,ext] = fileparts(outfile);

duplicate_O = [];
duplicate_Omml = [];
index_dups = find_repeats2 (output);

if isempty(index_dups)
    fprintf(1,['\nNo Duplicates Found.  Writing ' outfile '\n\n'])
else
    duplicate_O     = zeros (size(index_dups,1), size(index_dups,2));
    duplicate_Omml  = zeros (size(index_dups,1), size(index_dups,2));
    duplicate_ratio = zeros (size(index_dups,1), size(index_dups,2));
    for II = 1:size(index_dups,1)
        duplicate_O (II,:) = output(index_dups(II,:),Icolumn_O);
        duplicate_Omml (II,:) = output(index_dups(II,:),Icolumn_Omml);
        duplicate_ratio (II,:) = output(index_dups(II,:),Icolumn_thio);
    end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                CREATES A THE OXYGEN DUPLICATE TABLE                     %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if 1,    
    k=1;
    staa=min(output(:,1));
    stab=max(output(:,1));
    for i=staa:1:stab;
        aux_sta=find(output(:,1)==i & output(:,8)==2);
        for j=1:24
            aux=find(output(aux_sta,2)==j);
            if length(aux)>1;
                sta(k,1)=i;
                niskin(k,1)=j;
                botoxy1(k,1)=output(aux_sta(aux(1)),7);
                botoxy2(k,1)=output(aux_sta(aux(2)),7);
                k=k+1;
            end
        end  
    end
    fid=fopen(['duplicates_oxy.tab'],'w');
    fprintf(fid,'%s','\begin{longtable}{ccccc}');
    fprintf(fid,'\n');    
    fprintf(fid,'%s','\caption{Duplicate dissolved oxygen samples collected during the ABACO cruise (values in ?mol/kg).\label{long} }');
    fprintf(fid,'\n');
    fprintf(fid,'%s','\hline');
    fprintf(fid,'\n');
    bp=[sprintf('%s\t','Station &','Niskin &','Oxygen1 &','Oxygen2 &','Differences \\')];
    fprintf(fid,bp);
    fprintf(fid,'\n');    
    fprintf(fid,'%s','\hline');
    fprintf(fid,'\n');
    for ii=1:length(sta)
    fprintf(fid,'%2i\t%s%2i\t%s%3.1f\t%s%3.1f\t%s%3.3f%s\n',sta(ii), '&',...
        niskin(ii), '&', botoxy1(ii), '&', botoxy2(ii), '&', botoxy1(ii)-botoxy2(ii),'\\');
    end;
    fprintf(fid,'\n');    
    fprintf(fid,'%s','\hline');
    fprintf(fid,'\n');          
    fprintf(fid,'%s','\end{longtable}');
    fclose(fid);
end;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %
    % Average duplicates and plot differences
    %
%    plot(mean(duplicate_O'),diff(duplicate_O'),'+')
    diff_dup_O = diff(duplicate_O');
    Index_good_dups = find(abs(diff_dup_O) <= median(diff_dup_O)                 + 2*std(diff_dup_O));
    Index_good_dups = find(abs(diff_dup_O) <= median(diff_dup_O(Index_good_dups))+ 2*std(diff_dup_O(Index_good_dups)));
    Index_good_dups = find(abs(diff_dup_O) <= median(diff_dup_O(Index_good_dups))+ 2*std(diff_dup_O(Index_good_dups)));

    figure;
        plot(diff_dup_O,'+')
        hold on;
        [xlimits]=get(gca,'xlim');
        plot(xlimits, [median(diff_dup_O(Index_good_dups)) median(diff_dup_O(Index_good_dups))],'k--')
        plot(xlimits, [median(diff_dup_O(Index_good_dups))+2*std(diff_dup_O(Index_good_dups)) median(diff_dup_O(Index_good_dups))+2*std(diff_dup_O(Index_good_dups))],'r:')
        plot(xlimits, [median(diff_dup_O(Index_good_dups))-2*std(diff_dup_O(Index_good_dups)) median(diff_dup_O(Index_good_dups))-2*std(diff_dup_O(Index_good_dups))],'r:')
        [ylimits] = get(gca,'ylim');
        text([xlimits(1)+0.1*diff(xlimits)], [ylimits(2)-0.1*diff(ylimits)], ['Median = ' num2str(median(diff_dup_O(Index_good_dups)),5) ' +/- ' num2str(std(diff_dup_O(Index_good_dups)),5)])
        title (['Duplicate Oxygen'])% from ' upper(strrep(this_name,'_',' '))])
        set(gca,'xlim',[0.5 length(diff_dup_O)+0.5])
        xlabel('Station')
        ylabel('Dissolved Oxygen Differences (umol/kg)')
        set(gcf, 'PaperPositionMode', 'auto');
        print -depsc2 oxy_duplicates.eps
    %
    % Average good duplicates
    %
    avg_O = meanmiss(duplicate_O')';
    avg_Omml = meanmiss(duplicate_Omml')';
    avg_ratio = meanmiss(duplicate_ratio')';
    for II = 1: length(Index_good_dups);
        index_dups(Index_good_dups(II),1);
        output (index_dups(Index_good_dups(II),1),Icolumn_thio)    = avg_ratio(Index_good_dups(II));
        output (index_dups(Index_good_dups(II),1),Icolumn_O)       = avg_O(Index_good_dups(II));
        output (index_dups(Index_good_dups(II),1),Icolumn_Omml)    = avg_Omml(Index_good_dups(II));
        output (index_dups(Index_good_dups(II),1),Icolumn_quality) = 6;
    end;
    %
    % Remove the duplicates (by convention remove the second duplicate).
    %
    output(index_dups(:,2),:) = [];
end

fid=fopen(outfile,'w+');
fprintf(fid,'%s','%Sta   Niskin  Sample_Btl   Draw T  Thio_vol      Oxygen_mml       Oxygen        Flag    Cast\n');
fprintf(fid,'\n');
%fprintf(fid,'%d\t%d\t%d\t%d\t%1.5f\t%2f5\t%2f5\t%d\n',output');
fprintf(fid,'%d\t %d\t %d\t %2.1f\t %1.1d\t %1.5f\t %1.5f\t %d %d\n',output');
fclose(fid);

theResult = output;
return;


