% Read DMS files for HiWinGS 2013 data

disp(['read_DMS: ']);
fclose('all');
jd=sprintf('%03i',doy);

%declare arrays
jax=1;                              %line counter
nr2=13;
dmsdata=ones(nr2+1,36000)*NaN;


for hhh=0:23,					    %cycle thru 24 hourly stats files
    hr=sprintf('%02i',hhh);
    e=['G:\Data\CRUISES\HiWinGS_2013\Knorr\UH_DMS\' jd '\DMS_2013_' jd '_' hr '00.txt'];
    %e=[way_raw_data_flux, 'day',jd,'\met0' year(3:4) jd hr '_raw.txt'];
    disp(['Reading DMS file from hour ',e]);
    flist=fopen(e,'r');
    temp=[];stx=[];
    
    if flist>0,					    %if the file exists,
        while ~feof(flist)
            try
                temp =textscan(flist,'%f %*f %*f %*f %*f %*f %f %f %f %f %f %f %f %f %f %f %f %*f %f %*f %*f %*f %*f %*f %*[^\n]', 'delimiter', ',', 'headerlines', 1,  'emptyvalue', NaN);
                stx=[stx cell2mat(temp)'];
            catch
                [lnr2,cnr2]=size(temp{1,nr2});
                for ll=1:nr2
                     if length(temp{1,ll})~=lnr2
                         temp{1,ll}(lnr2+1)=[];
                     end;
                end;
                if isempty(temp{1,1}); continue; end
                if isempty(temp{1,1});
                    disp('Should never do this!')
                    temp{1,nr2}=temp{1,nr2}';
                end;
             stx=[stx cell2mat(temp)'];
            end;
        end;  %end while
    
        [nr1,nl1]=size(stx);
        hhrs=ones(1,nl1)*hhh;
        dmsdata(:,jax:jax+nl1-1)=[stx' hhrs']';
        jax=jax+nl1;
    
    else
        disp(sprintf('\b<--File not found.'))
        stx=ones(nr2,3600)*NaN;
%         stx(1,:)=ones(1,3600)*jam; stx(2,:)=ones(1,3600)*jam;
        [nr1,nl1]=size(stx);
        hhrs=ones(1,nl1)*jam;
        dmsdata(:,jax:jax+nl1-1)=[stx' hhrs']';
        jax=jax+nl1;
    end;      %     end if list
  
end;                                %end for jam

if doy==282
    dmsdata(1,1)= 24278400.0239;
end
jd_UH= doy+(dmsdata(1,:)-dmsdata(1,1))/3600/24; % time in decimal yearday
U_UH=(dmsdata(2,:)-2.5).*12;     % scale wind to physical units
V_UH=(dmsdata(3,:)-2.5).*12;  
W_UH=(dmsdata(4,:)-2.5).*12; 

uv = [U_UH;V_UH];
theta = (-30*(pi/180));
Tm = [(cos(theta)) (-sin(theta));(sin(theta)) (cos(theta))];
UVrot = Tm*uv;
U_UH = UVrot(1,:);
V_UH = UVrot(2,:);
% Reverse sign of U to agree with motion coordinates (forward positive)
U_UH = -U_UH;
   
Spd_Snd=dmsdata(5,:); 
%compute sonic temperature, NOAA R2 formula
SpdSnd = Spd_Snd.*16 + 290;  % speed of sound in m/s
Tsonic_UH = (SpdSnd.*SpdSnd + (U_UH.*U_UH + W_UH.*W_UH)/2. + V_UH.*V_UH)/403-273.15;
