whp_sum

 WHP_SUM reads in standard WOCE *.sum files 
         (as designated in table 3.5 in the
         'Requirements for WOCE Hydrograhic Programme Data Reporting'
 
 Usage: [lat,lon,station,time,depth] = whp_sum(flname)

 Inputs:
   flname: full name and path of *.sum file
   
 Outputs: (vectors with matched entries for each station)
   lat: latitude 
   lon: longitude
   station: station number
   time: decimal year time (e.g. 1992.34982)
   
   If the Sum file contains multiple entries for each station the values
   correspnding to the bottom of the cast are returned.  
   
   This file is a bit finicky about the format of the sum file.  Specifically
   if there are any letters (a-z;A-Z) in the body of the data it will crash.
   This usually happens when bottle numbers are specified with letter
   designations

 Paul E. Robbins  copyright 1995
 last update Dec 4 1996

0001 function [lat,lon,station,time,depth] = whp_sum(flname)
0002  
0003 % WHP_SUM reads in standard WOCE *.sum files
0004 %         (as designated in table 3.5 in the
0005 %         'Requirements for WOCE Hydrograhic Programme Data Reporting'
0006 %
0007 % Usage: [lat,lon,station,time,depth] = whp_sum(flname)
0008 %
0009 % Inputs:
0010 %   flname: full name and path of *.sum file
0011 %
0012 % Outputs: (vectors with matched entries for each station)
0013 %   lat: latitude
0014 %   lon: longitude
0015 %   station: station number
0016 %   time: decimal year time (e.g. 1992.34982)
0017 %
0018 %   If the Sum file contains multiple entries for each station the values
0019 %   correspnding to the bottom of the cast are returned.
0020 %
0021 %   This file is a bit finicky about the format of the sum file.  Specifically
0022 %   if there are any letters (a-z;A-Z) in the body of the data it will crash.
0023 %   This usually happens when bottle numbers are specified with letter
0024 %   designations
0025 %
0026 % Paul E. Robbins  copyright 1995
0027 % last update Dec 4 1996
0028 
0029 fid = fopen(flname,'r');
0030 if fid == -1
0031    %try to find file in matlab_data directory
0032   fid = fopen([getenv('MATLAB_DATA'),flname],'r');
0033   if fid == -1
0034     fid = fopen([getenv('MATLAB_DATA'),'bottle/woce/',flname],'r');    
0035     if fid == -1;
0036      disp(['File ',flname,' not found'])
0037      break
0038     end
0039   end    
0040 end
0041 
0042 %read in four lines of header data
0043 hdr1 = fgetl(fid);
0044 hdr2 = fgetl(fid);
0045 hdr3 = fgetl(fid);
0046 hdr4 = fgetl(fid);
0047 
0048 %use line 3 to figure out positions of desired columns
0049 iexpo    = findstr(hdr3,'EXPOCODE');
0050 isect    = findstr(hdr3,'SECT');
0051 istation = findstr(hdr3,'STNNBR');
0052 if length(istation)==0
0053   istation = findstr(hdr3,'STA');
0054 end
0055 
0056 icast    = findstr(hdr3,'CASTNO');
0057 itype    = findstr(hdr3,'TYPE');
0058 idate    = findstr(hdr3,'DATE');
0059 iutc     = findstr(hdr3,'TIME');
0060 icode    = findstr(hdr3,' CODE')+1;
0061 ilat     = findstr(hdr3,'LATITUDE');
0062 ilon     = findstr(hdr3,'LONGITUDE');
0063 inav     = findstr(hdr3,'NAV');
0064 idepth    = findstr(hdr3,'DEPTH');
0065 if length(idepth)==0
0066   idepth = min(findstr(hdr3,'MAX'));
0067 end
0068 
0069 i = 0;
0070 
0071 %  while 1
0072 for j = 1:1000;
0073   ln = [];ln = fgetl(fid);
0074   if ~isstr(ln) | all(ln == ' ') | length(ln) < 10 break; break; end % check if EOF or blank line
0075   %pad the line to full length if not long enough
0076   if length(ln) < inav;  
0077     ln = [ln blanks(inav-length(ln))];
0078   end;
0079   
0080   sect(j,:) = ln(isect:isect+4);
0081   station(j) = str2num(ln(isect+5:istation+6));
0082   cast(j) = str2num(ln(istation+7: icast+6));
0083   type(j,:) = ln(itype:itype+3);    
0084   
0085   % convert data and utc to a decimal year
0086   month= str2num(ln(idate:idate+1));
0087   day = str2num(ln(idate+2:idate+3));
0088   year = 1900+str2num(ln(idate+4:idate+5));
0089   hour = str2num(ln(iutc:iutc+1));
0090   if hour == []; hour = 0; end
0091   minutes  = str2num(ln(iutc+2:iutc+3));
0092   if minutes == []; minutes = 0; end
0093   time(j) =  year + cal2dec(month,day,hour,minutes)/365.25;
0094   code(j,1:2) =ln(icode+3:icode+4);
0095   if nargout > 4
0096     depth(j) = str2num(ln(idepth:idepth+4));
0097   end
0098    
0099 %  if ~strcmp(junk,10)
0100    %  read in latitude
0101    latstr = ln(ilat:ilon-1);
0102    %find blanks and retain only index to first blank if multiple blanks
0103    blnks = find(latstr == ' '); 
0104    blnks([0 diff(blnks) == 1]) = []; blnks(blnks == 1) = [];
0105    if length(blnks) < 3     %assume empty or courrpted
0106      lat(j) = nan;
0107    else
0108      latint = str2num(latstr(1:blnks(1)));
0109      latmin = str2num(latstr(blnks(1):blnks(2)));
0110      hemisphere = latstr(blnks(3)-1);
0111      if strcmp(hemisphere,'N')
0112        lat(j) = latint + latmin/60;
0113      elseif  strcmp(hemisphere,'S')
0114        lat(j) = -latint - latmin/60;
0115      else
0116        lat(j) = nan;     %if can't determine hemisphere    consider bogus
0117      end
0118    end
0119     %  read in longitude
0120    lonstr = ln(ilon:inav-1);
0121    blnks = find(lonstr == ' ');
0122    blnks([0 diff(blnks) == 1]) = [];blnks(blnks == 1) = [];
0123    if length(blnks) < 3     %assume empty or courrpted
0124      lon(j)= nan;
0125    else
0126      lonint = str2num(lonstr(1:blnks(1)));
0127      lonmin = str2num(lonstr(blnks(1):blnks(2)));
0128      hemisphere = lonstr(blnks(3)-1);   
0129      if strcmp(hemisphere,'E')
0130        lon(j) = lonint + lonmin/60;
0131      elseif strcmp(hemisphere,'W')
0132        lon(j) = -lonint - lonmin/60;    
0133      else
0134        lon(j)= nan;
0135      end
0136    end
0137 end
0138 
0139 %chuck the non-rosette entries (e.g. drifter deployments, mooring work,etc)
0140 ok =0*ones(length(lon),1);
0141 for j = 1:length(lon)
0142   if any(findstr(type(j,:),'RO')) 
0143     ok(j) = 1;
0144   end
0145 end
0146 
0147 lon = lon(ok); lat = lat(ok); station = station(ok); time = time(ok);
0148 
0149 fclose(fid);
0150 % check for multiple entries per station
0151 if any(diff(station) == 0) ;
0152   newstation = station;
0153   newstation(find(diff(newstation) == 0)) = [];
0154   newlat  = nan*newstation; newlon = newlat; newtime  = newlat;
0155   for i = 1:length(newstation);
0156     fs = find(newstation(i) == station);
0157     %check to see if any station entries are for bottom and use that time if
0158     %it exists
0159     
0160     found = 0;
0161     for j = 1:length(fs)
0162       if strcmp(code(fs(j),:),'BO')
0163     found = 1;
0164     newtime(i) = time(fs(j));
0165     newlat(i) = lat(fs(j));
0166     newlon(i) = lon(fs(j));
0167       end
0168     end
0169     if found == 0;  % if no bottom time found use start time
0170       for j = 1:length(fs)
0171     if strcmp(code(fs(j),:),'BE')
0172       found = 1;
0173       newtime(i) = time(fs(j));
0174       newlat(i) = lat(fs(j));
0175       newlon(i) = lon(fs(j));
0176       if nargout > 4
0177         newdepth(i) = depth(fs(j));
0178       end
0179     end
0180       end
0181     end
0182     %check to see if a lat and lon value is found
0183     if isnan(newlat(i)) | isnan(newlon(i));
0184       %use mean value of all non nan entries
0185       alllon  =lon(fs);
0186       alllat = lat(fs);
0187       newlon(i)  = mean(alllon(~isnan(alllon)));
0188       newlat(i)  = mean(alllat(~isnan(alllat)));
0189       if nargout > 4
0190     alldepth = depth(fs);
0191     newdepth(i) = max(alldepth(~isnan(alldepth)));
0192       end
0193     end
0194     if isnan(newtime(i))
0195       alltime  = time(fs);      
0196       newtime(i) =  mean(alltime(~isnan(alltime)));
0197     end
0198     
0199  end   
0200  station = newstation;  
0201  lat = newlat;
0202  lon = newlon;
0203  time = newtime;
0204  depth = newdepth;
0205 end
0206 
0207 
0208 
0209 
0210 
0211 
0212 
0213