%%% loadhurrwind.m %%%%%%%%%%%%%%%%%

if ~exist('stormfile')
    stormfile=uigetfile('*.txt', 'Select storm forecast file');
end;


if ~exist('Onestorm')
    Onestorm=0;
end;

% read file
fprintf(1,'Read storm file %s\n',stormfile);
a=textread(stormfile,'%s');


ind=22;
day=str2num(char(a(1:ind:length(a))));
hr=char(a(2:ind:length(a)));
lats=char(a(3:ind:length(a)));
lons=char(a(4:ind:length(a)));
r64_ne=char(a(7:ind:length(a)));
r64_se=char(a(8:ind:length(a)));
r64_sw=char(a(9:ind:length(a)));
r64_nw=char(a(10:ind:length(a)));
r50_ne=char(a(13:ind:length(a)));
r50_se=char(a(14:ind:length(a)));
r50_sw=char(a(15:ind:length(a)));
r50_nw=char(a(16:ind:length(a)));
r34_ne=char(a(19:ind:length(a)));
r34_se=char(a(20:ind:length(a)));
r34_sw=char(a(21:ind:length(a)));
r34_nw=char(a(22:ind:length(a)));
clear ind a
  
% determine month of earliest value in stormfile 
if ~exist('dep')
    month=input(['month (8=Aug, 9=Sep, etc.) of earliest value in ',stormfile,'? ']);
    % assume that the year is the current year
    year=str2num(datestr(now,10));
else
    % assume that the day of deployment, dep(3), is 
    % in the forecast period 
    q=find(day==dep(3)); 
    if (day(1)<=day(q))  % no jump across month
        month=dep(2);
    else
        month=dep(2)-1;
    end; 
    year=dep(1);
end; 
 
lat=NaN*ones(length(lats),1);
lon=lat; t=lat;
r64=NaN*ones(length(lat),4);
r50=r64; r34=r64; r=r64;
 
for i=1:length(lat);
    lat(i)=str2num(lats(i,1:size(lats,2)-1));
    ind=1:size(lons,2);
    for j=1:max(ind);
        if (strcmp('W',lons(i,j)))
            ind=1:j-1;
        end;
    end;
    lon(i)=-str2num(lons(i,ind));

    for radfield=1:3;
        r=NaN*ones(size(r));
        if (radfield==1)
            rne=r64_ne; rse=r64_se;
            rsw=r64_sw; rnw=r64_nw;
        elseif (radfield==2)
            rne=r50_ne; rse=r50_se;
            rsw=r50_sw; rnw=r50_nw;
        else
            rne=r34_ne; rse=r34_se;
            rsw=r34_sw; rnw=r34_nw;
        end;

        for j=1:size(rne,2)-1
            if strcmp(rne(i,j:j+1),'NE')
                r(i,1)=str2num(rne(i,1:j-1));
            end;
        end;
        for j=1:size(rse,2)-1
            if strcmp(rse(i,j:j+1),'SE')
                r(i,2)=str2num(rse(i,1:j-1));
            end;
        end;
        for j=1:size(rsw,2)-1
            if strcmp(rsw(i,j:j+1),'SW')
                r(i,3)=str2num(rsw(i,1:j-1));
            end;
        end;
        for j=1:size(rnw,2)-1
            if strcmp(rnw(i,j:j+1),'NW')
                r(i,4)=str2num(rnw(i,1:j-1));
            end;
        end;
        if (radfield==1) r64(i,:)=r(i,:);
        elseif (radfield==2) r50(i,:)=r(i,:);
        else r34(i,:)=r(i,:);
        end;
    end;    % loop over radfield
    clear radfield
end;
clear r i lats lons *_ne *_nw *_se *_sw rne rnw rse rsw 

% calculate time t
hour=str2num(hr(:,1:2));
minute=str2num(hr(:,3:4));
clear hr;
year=year*ones(length(hour),1);
month=month*ones(length(hour),1);
% if the day suddenly decreases to 1, assume
% that you've jumped forward a month
q=find(diff(day)<0);
if (length(q))
    if (day(q+1)~=1)
        disp(['warning: time not increasing monotonically in ',stormfile,'\n'])
    end;
    month(q(1)+1:length(month))=month(q(1)+1:length(month))+1;
end;
clear q;

    
t=juliand([year month day hour minute zeros(size(minute))]);
clear year month day hour minute 

% visualize results
if Onestorm~=1
    plot(lon,lat);
    hold on;
    windpatch(t,lon,lat,r64,r50,r34)
    hold on;
end

%%%%%%%%%  END load hurrwind.m  %%%%%%%

