%This file reads the good.list file and loads the good stations
%from a folder containing a section of the xbt AX8 line.
%This script is created by Marlos
%Last change: 09/28/2010

function [data_full,lat,lon,variables,cruise,month,TIME] = load_xbt_file(cruise,month)
warning('off')

addpath /phodnet/share/mgoes/TS_PACKAGE_THACKER_noyear_globe   %ORIG
%addpath /data/mgoes/SA_gyre/Argo/2015/Atlantic                  %FOR NEWSAL
%addpath /data/mgoes/SA_gyre/Argo/2015/Atlantic2                  %FOR NEWSAL2

%Choose the cruise you want to use (stored in file bad.list)
pout = [5:10:750 800:50:2000 2100:100:6000]; %OUTPUT depth

%Define data directory
%suffix = ['/WSmounts/phodnet/mobster/m1/data/xbt/aoml/shallow/',cruise,'_2db/'];
suffix = ['/WSmounts/phodnet/mobster/m1/data/xbt/aoml/shallow/',cruise,'/']; %WOA05

preffix = month2date_temp(month);
%preffix = datestr(datenum(month,'yymm'),'yyyymm');

preffix = preffix([5 6 3 4]);
direc   = strcat(suffix,month,'/');
filepreffix = [direc,cruise,preffix];

eval(['!basename -a ',filepreffix,'.* > fiels']);
files = importdata('fiels');

%Get good data info
lista  = strcat(suffix,month,'/good.list');

 % retrieve only good data
  if exist(lista,'file') 
     filen = load(lista);
     filen = filen(filen~=0);
   station = num2str(filen*1e-3,'%1.3f');
   station = station(:,2:end);
     files = strcat(cruise,preffix,station);
  end
% Declare and/or initialize a few parameters
lat = [];
lon = [];
data_full = [];       %matrix where all data is stored
TIME = [];

for i = 1:size(files,1);
% Open Data File
  filename = [direc,files(i,:)];
 
   dat = importdata(filename,' ',12);
   header = dat.textdata;
   x = dat.data;
   maxp = size(x,1);

   variables = header{11,:}; 
   [ca,latn,lonn,cb,dat,cc]=strread(header{2,:},'%d%f%f%d%s%d'); 
    datx{1} = char(dat);
    lat(i) = latn;
    lon(i) = lonn;
   TIME(i) = str2num(datestr(datx,'yyyymmdd'));
  data_full(1:maxp,:,i) =  x;

 end

 %fix direction if needed
  idir = nanmedian(diff(lon));
 if idir < 0
    dirf = i:-1:1;
     lon = lon(dirf);
     lat = lat(dirf)
    TIME = TIME(dirf);
    data_full = data_full(:,:,dirf);
 end 

 s = squeeze(data_full(:,4,:));
 t = squeeze(data_full(:,2,:));
 d = squeeze(data_full(:,6,:));
 p = squeeze(data_full(:,1,:));

 anao = t==0;
 t(anao)=nan;
 d(anao)=nan;
 s(anao)=nan;
 p(anao)=nan;
clear data_full

%CALCULATE SALINITY WITH PADDING AND VERTICAL INTERPOLATION
pad = 900; %PAD WOA13 CLIMATOLOGY AT 900 meters

[S0,~,t0,pp]=Calc_sal_Thacker_Goes_EmDr_Stom_svd_globe(t,d,lat,lon,TIME,900,pout,'goes');
 nlat = length(lat);
   nz = length(pout);
 pres = sw_pres(pp(ones(nlat,1),:)',lat(ones(1,nz),:));
 anao = isnan(t0)|isnan(S0);
 t0(anao) = nan;
 S0(anao) = nan;
 ptmp = sw_ptmp(S0,t0,pp',0);

%CALCULATE VELOCITY
 ref_lev = 1000;                                               %REF LEV FOR GEOSTROPHIC VEL
 [DH,vg] = calc_dynH_fixlev(lon,lat,nanmean(pres,2),t0,S0,ref_lev);
 data_full(1:length(pout),1,:) = pres;                         %PRESSURE
 data_full(1:length(pout),2,:) = t0;                           %TEMP
 data_full(1:length(pout),3,:) = ptmp;                         %POTENTIAL TEMP
 data_full(1:length(pout),4,:) = S0;                           %SALINITY
 data_full(1:length(pout),5,:) = DH;                           %DYNAMIC HEIGHT (CM)
 data_full(1:length(pout),6,:) = pp(ones(1,nlat),:)';          %DEPTH
 data_full(1:length(pout),7,:) = vg;                          %GEOSTROPHIC VELOCITY
 
 pat = '\s+';
     s = regexp(variables, pat, 'split');
  ind1 = find(strcmp(s,'pr'));
variables = s(ind1:end);variables{7} = 'vg'

return

pcolor(lon,-pout,vg);shading flat
%caxis([-1 1])
colorbar
drawnow
end




