function oxy_cal_stats(cruiseid)
% oxy_cal_stats -  Statistics about the condcutivity calibration to be used on the Data Report.
%
% CTD Calibration toolbox.
%
% USAGE:
% oxy_cal_stats(cruiseid);
%
%INPUT: cruiseid: Cruise identification name.
%   
% OUTPUT: Statistics about the oxygen calibration to be used on the Data Report..
%          
% DESCRIPTION: This function uses the summary conductivity calibration file created by  
% make_cond_sumfile....
%
% AUTHORS:Carlos Fonseca 
%         UM/CIMAS
%         Derrick Snowden
%         NOAA/AOML/PhOD
%         Tue May 11 11:23:37 EDT 2004
%
% CHANGELOG: 
%   23-Jul-2004, Version 1.0
%		* Initial version.

% Retrieve the necessary preferences for this cruise

%%%

%%%
group = group_name(cruiseid);

% cruise_dir is the main calibration data directory 
if ispref(group,'cal_data_dir')
	cruise_dir=getpref(group,'cal_data_dir');
else
	error(['The cal_data_dir preference was not set for cruise: ' cruiseid '.  Run register_cruise.m']);
end

% Load data from the cruise database
db_file = fullfile(cruise_dir,[cruiseid '_db.mat']);
if exist(db_file)==2
ld_file=(['load ',db_file,' oxy_match_cal oxy_match']);
eval(ld_file)
else
	error(['The cruise data base file was not found on the path.  ' db_file]);
end
n=who('oxy_match_cal');
 if isempty(n)==1
 error(['oxy_match_cal not found. Run make_cal_sumfile.m ']);
 return;
 end
disp('****************************************** ')
disp('PRELIMINARY COMPARISONS')
disp('****************************************** ')
disp('The data found in cond_match goes from') 
disp(['the station ',num2str(min(oxy_match_cal.stnnbr)),' to the station ',num2str(max(oxy_match_cal.stnnbr)),'.'] ) 
disp('****************************************** ')
%
staa=input('Input first station: ');
stab=input('Input last station: ');

%%%%%%%%
%
disp(' 0 = use the primary sensor (DEFAULT)' )
disp(' 1 = use the secondary sensor' )
pri_sec=input('Do the calibration was applied to the primary or to the secondary sensor: ');
while (isempty(pri_sec)==1 | pri_sec<0 | pri_sec>1)
	disp('Please answer 0 (zero) for primary or  1(one) for secondary ')
	pri_sec=input('Do the calibration was applied to the primary or to the secondary sensor ');
end
disp('***************************************************** ')
YN=input('Do you want to remove any niskin bottle? y/n ','s');
bad_niskin=1;
if strcmpi(YN,'y'),
    while bad_niskin~=0,
    bad_niskin=input(['Enter the niskin number to remove or ', ...
	'0 to skip: ']);  
    bad_index=find(oxy_match_cal.niskinbtlpos==bad_niskin);
    oxy_match_cal.quality1(bad_index)=1;
    end;
end
%s_index=find( cond_match_cal.stnnbr>= staa & cond_match_cal.stnnbr  <= stab   & cond_match_cal.quality1==2 & cond_match_cal.quality2==2);
s_index=find( oxy_match.stnnbr>= staa & oxy_match.stnnbr  <= stab);
disp('**************************************************************************************')
disp('First the Raw Values...')
disp(['Number of good data points found: ',num2str(length(s_index))])
aux_gd=find( oxy_match.stnnbr>= staa & oxy_match.stnnbr  <= stab);
per_gd=(length(s_index)/length(aux_gd))*100;
disp(['This represents ', num2str(per_gd),'% of the data']); 
%
 station_num_raw=oxy_match.stnnbr(s_index);
  if pri_sec==1
 bottle_oxy_raw=ox_units(oxy_match.btloxy(s_index),oxy_match.ctdsal2(s_index),sw_ptmp(oxy_match.ctdsal2(s_index),oxy_match.ctdtmp2(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');   
 ctd_oxy_raw=ox_units(oxy_match.ctdoxy2(s_index),oxy_match.ctdsal2(s_index),sw_ptmp(oxy_match.ctdsal2(s_index),oxy_match.ctdtmp2(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');
  else
 bottle_oxy_raw=ox_units(oxy_match.btloxy(s_index),oxy_match.ctdsal1(s_index),sw_ptmp(oxy_match.ctdsal1(s_index),oxy_match.ctdtmp1(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');       
 ctd_oxy_raw=ox_units(oxy_match.ctdoxy1(s_index),oxy_match.ctdsal1(s_index),sw_ptmp(oxy_match.ctdsal1(s_index),oxy_match.ctdtmp1(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');   
 end
 ctd_pressure_raw=oxy_match.ctdprs(s_index);
aux_500=find(ctd_pressure_raw>=500);
aux_3000=find(ctd_pressure_raw>=3000);
disp(['average residual ',num2str(mean(bottle_oxy_raw-ctd_oxy_raw)),' sd ',num2str(std(bottle_oxy_raw-ctd_oxy_raw))])
disp(['average residual (below 500) ',num2str(mean(bottle_oxy_raw(aux_500)-ctd_oxy_raw(aux_500))),' sd ',num2str(std(bottle_oxy_raw(aux_500)-ctd_oxy_raw(aux_500)))])
disp(['average residual (below 3000)',num2str(mean(bottle_oxy_raw(aux_3000)-ctd_oxy_raw(aux_3000))),' sd ',num2str(std(bottle_oxy_raw(aux_3000)-ctd_oxy_raw(aux_3000)))])
 
s_index=find( oxy_match.stnnbr>= staa & oxy_match.stnnbr  <= stab   & oxy_match.quality1==2 & oxy_match.quality2==2);
disp('**************************************************************************************')
disp('Now the QC Values...')
disp(['Number of good data points found: ',num2str(length(s_index))])
aux_gd=find( oxy_match.stnnbr>= staa & oxy_match.stnnbr  <= stab);
per_gd=(length(s_index)/length(aux_gd))*100;
disp(['This represents ', num2str(per_gd),'% of the data']); 
%
 station_num_raw=oxy_match.stnnbr(s_index);
if pri_sec==1
 bottle_oxy_raw=ox_units(oxy_match.btloxy(s_index),oxy_match.ctdsal2(s_index),sw_ptmp(oxy_match.ctdsal2(s_index),oxy_match.ctdtmp2(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');   
 ctd_oxy_raw=ox_units(oxy_match.ctdoxy2(s_index),oxy_match.ctdsal2(s_index),sw_ptmp(oxy_match.ctdsal2(s_index),oxy_match.ctdtmp2(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');
  else
 bottle_oxy_raw=ox_units(oxy_match.btloxy(s_index),oxy_match.ctdsal1(s_index),sw_ptmp(oxy_match.ctdsal1(s_index),oxy_match.ctdtmp1(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');       
 ctd_oxy_raw=ox_units(oxy_match.ctdoxy1(s_index),oxy_match.ctdsal1(s_index),sw_ptmp(oxy_match.ctdsal1(s_index),oxy_match.ctdtmp1(s_index),oxy_match.ctdprs(s_index),0),'mll2mmk');   
 end
 ctd_pressure_raw=oxy_match.ctdprs(s_index);
aux_500=find(ctd_pressure_raw>=500);
aux_3000=find(ctd_pressure_raw>=3000);
disp(['average residual ',num2str(mean(bottle_oxy_raw-ctd_oxy_raw)),' sd ',num2str(std(bottle_oxy_raw-ctd_oxy_raw))])
disp(['average residual (below 500) ',num2str(mean(bottle_oxy_raw(aux_500)-ctd_oxy_raw(aux_500))),' sd ',num2str(std(bottle_oxy_raw(aux_500)-ctd_oxy_raw(aux_500)))])
disp(['average residual (below 3000)',num2str(mean(bottle_oxy_raw(aux_3000)-ctd_oxy_raw(aux_3000))),' sd ',num2str(std(bottle_oxy_raw(aux_3000)-ctd_oxy_raw(aux_3000)))])

s_index=find( oxy_match_cal.stnnbr>= staa & oxy_match_cal.stnnbr  <= stab   & oxy_match_cal.quality1==2 & oxy_match_cal.quality2==2);
disp('**************************************************************************************')
disp('And  finally the Calibrated Values...')
disp(['Number of good data points found: ',num2str(length(s_index))])
aux_gd=find( oxy_match_cal.stnnbr>= staa & oxy_match_cal.stnnbr  <= stab);
per_gd=(length(s_index)/length(aux_gd))*100;
disp(['This represents ', num2str(per_gd),'% of the data']); 
%
 station_num=oxy_match_cal.stnnbr(s_index);
 bottle_oxy=ox_units(oxy_match_cal.btloxy(s_index),oxy_match_cal.ctdsal(s_index),sw_ptmp(oxy_match_cal.ctdsal(s_index),oxy_match_cal.ctdtmp(s_index),oxy_match_cal.ctdprs(s_index),0),'mll2mmk');
 ctd_oxy_1=oxy_match_cal.ctdoxyumolkg(s_index);
 niskin_num=oxy_match_cal.niskinbtlpos(s_index);
 sample_bottle_num=oxy_match_cal.sampnbr(s_index);
 ctd_pressure=oxy_match_cal.ctdprs(s_index);

 aux_500=find(ctd_pressure>=500);
aux_3000=find(ctd_pressure>=3000);
disp(['average residual ',num2str(mean(bottle_oxy-ctd_oxy_1)),' sd ',num2str(std(bottle_oxy-ctd_oxy_1))])
disp(['average residual (below 500) ',num2str(mean(bottle_oxy(aux_500)-ctd_oxy_1(aux_500))),' sd ',num2str(std(bottle_oxy(aux_500)-ctd_oxy_1(aux_500)))])
disp(['average residual (below 3000)',num2str(mean(bottle_oxy(aux_3000)-ctd_oxy_1(aux_3000))),' sd ',num2str(std(bottle_oxy(aux_3000)-ctd_oxy_1(aux_3000)))])


%%%%%%%%%%%%%%%
%Start Ploting%
%%%%%%%%%%%%%%%
return