function[ws, wd, date_time, station_id]= PK_WND(full_filename) if ~exist('full_filename') || isempty(full_filename) [filename, pathname]= uigetfile({'*.op','OP (ASOS records)(*.op)'},'Select ASOS file'); if any(pathname~=0) full_filename = fullfile( pathname, filename ); else error('File selection cancelled.'); end end fid = fopen( full_filename ); if fid==-1 err = errordlg({'The following file could not be opened for reading:',full_filename},'Error reading file'); uiwait(err); return; end h = waitbar(0,['File name: ' strrep(strrep(full_filename,'\','\\'),'_','\_')],... 'Name','Extracting peak wind data...', 'WindowStyle','modal'); drawnow; file_contents = textscan( fid, '%s', 'whitespace', '', 'delimiter', '\n'); file_contents = file_contents{1}; fclose(fid); ws = zeros(length(file_contents),1); wd = zeros(length(file_contents),1); date_time = zeros(length(file_contents),1); station_id = zeros(length(file_contents),1); last_ws = 0; last_timestr = 'zero'; last_time = 0; pkwindow = 1/12; % time interval (days) within which equivalent peak wind reports (same speed and time, or missing time) are considered repeated count = 0; warn_count = 0; warn_cell{4} = 'Errors were encountered in the "PK WND" reports on the following lines, and these reports were rejected:'; for i=1:length(file_contents) textline = file_contents{i}; if length(textline)<22 continue end ac = str2num(textline(5:10)); datetimestr=textline(11:22); yr=str2num(datetimestr(1:4)); mo=str2num(datetimestr(5:6)); dy=str2num(datetimestr(7:8)); hr=str2num(datetimestr(9:10)); mn=str2num(datetimestr(11:12)); linetime=datenum(yr,mo,dy,hr,mn,0); k = min(strfind(textline,'PK')); if ~isempty(k) windstr_ex = textline(k+2:length(textline)); j = min(strfind(windstr_ex,'/')); if isempty(j) || j>30 warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': No slash found within 30 characters of the "PK" code']; continue end jj = min(find(ismember(windstr_ex(1:j-1),'0123456789'))); if isempty(jj) warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': No numeric characters found preceding slash' ]; continue end windstr = windstr_ex(jj:j-1); if any(~ismember(windstr_ex(jj:j-1),'0123456789')) warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': Non-consecutive numeric characters found preceding slash: ' windstr]; continue end n_windstr=length(windstr); y = min(find(~ismember(windstr_ex(j+1:length(windstr_ex)),'0123456789'))); timestr = windstr_ex(j+1:(j+(y-1))); n_timestr=length(timestr); switch n_windstr case 3 wdir=str2num(windstr(1))*10; wspd=str2num(windstr(2:3)); case 4 if str2num(windstr(1))>3 | str2num(windstr(3))<=1 wdir = str2num(windstr(1))*10; wspd = str2num(windstr(2:4)); else wdir=str2num(windstr(1:2))*10; % wind directions were in tens of degrees until July 1,1996 for all NYC stations wspd=str2num(windstr(3:4)); end case 5 if str2num(windstr(3)) == 0 wdir=str2num(windstr(1:3)); wspd=str2num(windstr(4:5)); elseif str2num(windstr(3)) == 1 wdir = str2num(windstr(1:2))*10; wspd = str2num(windstr(3:5)); else warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': Ambiguous 5-digit code found before slash: ' windstr]; continue end case 6 wdir=str2num(windstr(1:3)); wspd=str2num(windstr(4:6)); otherwise warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': ' num2str(n_windstr) ... ' numeric characters found before slash (3 to 6 expected): ' windstr]; continue end date_num_i = case_time(timestr, datetimestr); if wspd == last_ws && (linetime-last_time) skip to next line end if wdir<=360 & wspd>=25 & wspd<199 count = count+1; wd(count)=wdir; ws(count)=wspd; date_time(count)=date_num_i; station_id(count) = ac; if ishandle(h), waitbar(i/length(file_contents),h); drawnow; end % store temporary values to check for repeated reports: last_ws = wspd; last_timestr = timestr; last_time = date_num_i; else warn_count = warn_count+1; warn_cell{warn_count+4}= ['Line ' num2str(i) ': Reported wind speed (' num2str(wspd) ' knot) ' ... 'or wind direction (' num2str(wdir) ' degrees) exceeds allowable range: ' windstr]; continue end end end wd=wd(1:count); ws=ws(1:count); date_time=date_time(1:count); station_id = station_id(1:count); if ishandle(h), waitbar(1,h); drawnow; close(h); end warn_cell{1} = ['Filename: ' full_filename]; warn_cell{2} = ['Number of "PK WND" reports extracted: ' num2str(length(ws))]; % updated below warn_cell{3} = ' '; if warn_count>0 wrn = helpdlg(warn_cell,['Extraction of peak wind data completed']); uiwait(wrn); end