added corrected NDBC_netcdf

NDBC_hdf5/angle_onlns.m

@@ -0,0 +1,9 @@
+function [sint cost sin2t cos2t] = angle_onlns(dtheta,nang)
+
+anguse = 270 - dtheta.*[0:nang-1];
+theta = 270 - anguse;
+tmp1 = pi/180 * theta;
+sint = sin(tmp1);
+cost = cos(tmp1);
+sin2t = sin(2.*tmp1);
+cos2t = cos(2.*tmp1);

NDBC_hdf5/create_onlns.m

@@ -0,0 +1,60 @@
+function create_onlns(stat,year,mon,aa,sdir)
+
+onlform = ['%5s%6i%3i%3i%3i%3i%5i%5i%5i%5i%5i%5i%7.1f%7.2f%7.2f%7.2f', ...
+    '%7.1f%7.2f%7.2f%7.1f%8.2f%8.2f%8.2f%8.2f%7.1f%7.1f%7.1f%7.1f', ...
+    '%8.2f%8.2f%8.2f%8.2f%8.2f%10.4f%7.1f%7.1f%8.2f% 14.5E% 14.5E', ...
+    '% 14.5E% 14.5E\n'];
+sp1form1 = ['%5s%6i%3i%3i%3i%3i'];
+
+if ~exist(sdir,'dir')
+    mkdir(sdir);
+end
+fout = [sdir,'\n',stat,'_',year,'_',mon];
+foutone = [fout,'.onlns'];
+fid = fopen(foutone,'w');
+for qq = 1:size(aa,1)
+   fprintf(fid,onlform,stat,aa(qq,:));
+end
+fclose(fid);
+        % print spec1d file
+%     if p == 1
+%         fspec1 = [fout,'.spe1d'];
+%         if nn == 1
+%             fid = fopen(fspec1,'w');
+%         else
+%             fid = fopen(fspec1,'a+');
+%         end
+%         for qq = 1:size(cc.c11,2)
+%             fprintf(fid,sp1form1,stat,aa(qq,1:5));
+%             for jj = 1:length(cc.freq)
+%                 fprintf(fid,'%8.4f%12.6f',cc.freq(jj),cc.c11(jj,qq));
+%             end
+%             fprintf(fid,'\n');
+%         end
+%         fclose(fid);
+%     else
+%     % print spec2d file
+%         fspec1 = [fout,'.spe1d'];
+%         fspec2 = [fout,'.spe2e'];
+%         if nn == 1
+%             fid  = fopen(fspec1,'w');
+%             fid2 = fopen(fspec2,'w');
+%         else
+%             fid = fopen(fspec1,'a+');
+%             fid2 = fopen(fspec2,'a+');
+%         end
+%         for qq = 1:size(cc.c11,2)
+%             fprintf(fid,sp1form1,stat,aa(qq,1:5));
+%             fprintf(fid2,sp1form1,stat,aa(qq,1:5));
+%             for jj = 1:length(cc.freq)
+%                 fprintf(fid,'%8.4f%12.6f',cc.freq(jj),cc.c11(jj,qq));
+%                 fprintf(fid2,'%8.4f%8.4f%8.4f%6.1f%6.1f%12.6f',cc.freq(jj), ...
+%                     cc.r1(jj,qq),cc.r2(jj,qq),cc.alpha1(jj,qq), ...
+%                     cc.alpha2(jj,qq),cc.c11(jj,qq));
+%             end
+%             fprintf(fid,'\n');
+%             fprintf(fid2,'\n');
+%         end
+%         fclose(fid);
+%         fclose(fid2);
+%     end

NDBC_hdf5/delftsdrs.m

@@ -7,23 +7,23 @@
 idt = -1;
 
 dire = 270 - dtheta.*[0:ide-1];
-dfsq = repmat(delf',[1 nang]);
+dfsq = repmat(delf,[1 nang]);
 etf = sum(dfsq(1:ife,:).*(ef2d(1:ife,:) + ef2d(2:ife+1,:)));
 
 ii = find(etf == max(etf));
 etop = etf(ii);
 idt = ii;
 
 if idt > 0
-    delftpdr = dire(idt);
+    delftpdr = dire(max(idt));
     if delftpdr < 0
         delftpdr = delftpdr + 360;
     else
         delftpdr = -999.;
     end
 end
 
-[sint cost sin2t cos2t] = angle(dtheta,nang);
+[sint cost sin2t cos2t] = angle_onlns(dtheta,nang);
 edt = etf*ainc;
 sum0 = sum(edt);
 eex = sum(cost.*edt);

NDBC_hdf5/delftsums.m

@@ -21,7 +21,7 @@
     delfthmo = 4.0*sqrt(sum0);
     ii = find(ef == max(ef));
     etop = ef(ii);
-    ift = ii;
+    ift = ii(1);
 
     if ift <= nfre
         delfttp = 1./fr(ift);

NDBC_hdf5/directns.m

@@ -20,11 +20,11 @@
 tkh(ii) = 20;
 cgg = 0.5.*cphz.*(1.0 + tkh./sinh(tkh));
 
-[sint cost sin2t cos2t] = angle(dtheta,nang);
+[sint cost sin2t cos2t] = angle_onlns(dtheta,nang);
 sa = length(fr);
 sint = repmat(sint,[sa 1]);cost = repmat(cost,[sa 1]);
 sin2t = repmat(sin2t,[sa 1]);cos2t = repmat(cos2t,[sa 1]);
-delsq = repmat(delf',[1 nang]);
+delsq = repmat(delf,[1 nang]);
 
 xsum = sum(cost.*ef2d.*delsq,2);
 ysum = sum(sint.*ef2d.*delsq,2);
@@ -68,20 +68,24 @@
 
 wavpd = peakdr(ef2d,mp,nfre,nang,dtheta);
 wavpd = 180/pi*wavpd;
-
+% 
 wavdvt = 270. -  wavdvt;
+%wavdvt = wavdvt - 180;
 if wavdvt < 0
     wavdvt = wavdvt + 360;
 end
 wavdvfm = 270 - wavdvfm;
+%wavdvfm = wavdvfm - 180;
 if wavdvfm < 0
     wavdvfm = wavdvfm + 360;
 end
 wavdmx = 270 - wavdmx;
+%wavdmx = wavdmx - 180;
 if wavdmx < 0
     wavdmx = wavdmx + 360;
 end
 wavpd = 270 - wavpd;
+%wavpd = wavpd - 180;
 if wavpd < 0
     wavpd = wavpd + 360;
 end

NDBC_hdf5/get_NDBC_air.m

@@ -1,23 +1,26 @@
 function [aa payload] = get_NDBC_air(fname)
-
+%
+%   read NODC netCDF4 file to get met information
+%   created 07/12 by TJ. Hesser
+% 
 fprintf(1,'Analyizing Station %5s\n',fname(6:10));
 %Get attributes distributed
-fileatt = h5info(fname);
+fileatt = ncinfo(fname);
 fatt = {fileatt.Attributes.Name};
 
 % Get time out of the file and convert to year mon day hour min sec
-times = double(h5read(fname,'/time'));
+times = double(ncread(fname,'/time'));
 [year,mon,day,hour,min,sec] = datevec( ...
         datenum(1970,1,1)+times./86400);
 
 % Get lat and lon out of file and convert to deg min sec  
 latm = strcmp(fatt,'nominal_latittude');
-lat = h5readatt(fname,'/','nominal_latitude');
+lat = ncreadatt(fname,'/','nominal_latitude');
 [lad,lam,las] = deg2dms(abs(lat));
 lam = abs(lam);las = round(las);las = abs(las);
 lad = repmat(lad,[length(year) 1]);lam = repmat(lam,[length(year) 1]);
 las = repmat(las,[length(year) 1]);
-lon = h5readatt(fname,'/','nominal_longitude');
+lon = ncreadatt(fname,'/','nominal_longitude');
 [lod,lom,los] = deg2dms(abs(lon));
 lom = abs(lom);los = abs(los);los = round(los);
 lod = repmat(lod,[length(year) 1]);lom = repmat(lom,[length(year) 1]);
@@ -32,44 +35,44 @@
 
 % Get depth out of file
 try
-    dep = h5readatt(fname,'/','sea_floor_depth_below_sea_level');
+    dep = ncreadatt(fname,'/','sea_floor_depth_below_sea_level');
 catch
     dep = -999.00;
 end
 dep = repmat(dep,[length(year) 1]);
 
 % Get windspeed, windgust, wind direction
 try
-    wspd = h5read(fname,['/payload_1/anemometer_1/wind_speed']);
+    wspd = ncread(fname,['/payload_1/anemometer_1/wind_speed']);
     payload = 'payload_1';
-    payatt = h5readatt(fname,'/payload_1','description');
+    payatt = ncreadatt(fname,'/payload_1','description');
     np = length(payatt);sp = ['%',num2str(np),'s\n'];
     fprintf(1,'Using payload_1\n');
     fprintf(1,sp,payatt);
 catch
-    wspd = h5read(fname,['/payload_2/anemometer_1/wind_speed']);
+    wspd = ncread(fname,['/payload_2/anemometer_1/wind_speed']);
     payload = 'payload_2';
-    payatt = h5readatt(fname,'/payload_2','description');
+    payatt = ncreadatt(fname,'/payload_2','description');
     np = length(payatt);sp = ['%',num2str(np),'s\n'];
     fprintf(1,'Using payload_2\n');
     fprintf(1,sp,payatt);
 end
-wspd_q = h5read(fname,['/',payload,'/anemometer_1/wind_speed_qc']);
+wspd_q = ncread(fname,['/',payload,'/anemometer_1/wind_speed_qc']);
 ii = wspd_q ~= 0;
 wspd(ii) = -99.99;
 if sum(ii) > 0
-  wspd2 = h5read(fname,['/',payload,'/anemometer_2/wind_speed']);
-  wspd_q2 = h5read(fname,['/',payload,'/anemometer_2/wind_speed_qc']);
+  wspd2 = ncread(fname,['/',payload,'/anemometer_2/wind_speed']);
+  wspd_q2 = ncread(fname,['/',payload,'/anemometer_2/wind_speed_qc']);
   jj = wspd_q2(ii) == 0;
   wspd(ii(jj)) = wspd2(ii(jj));
 end
-wspd_info = h5info(fname,['/',payload,'/anemometer_1/wind_speed']);
+wspd_info = ncinfo(fname,['/',payload,'/anemometer_1/wind_speed']);
 ii = wspd == wspd_info.FillValue;
 wspd(ii) = -99.99;
 
 % Get anemometer height
 try
-    anemh = h5readatt(fname,['/',payload,'/anemometer_1','height_of_instrument']);
+    anemh = ncreadatt(fname,['/',payload,'/anemometer_1','height_of_instrument']);
 catch
     anemh = 'NaN';
     re = 0;
@@ -88,58 +91,58 @@
 end
 anemh = repmat(anemh,[length(year) 1]);
 
-wgus = h5read(fname,['/',payload,'/anemometer_1/wind_gust']);
-wgus_q = h5read(fname,['/',payload,'/anemometer_1/wind_gust_qc']);
+wgus = ncread(fname,['/',payload,'/anemometer_1/wind_gust']);
+wgus_q = ncread(fname,['/',payload,'/anemometer_1/wind_gust_qc']);
 ii = wgus_q ~= 0;
 wgus(ii) = -99.99;
 if sum(ii) > 0
-  wgus2 = h5read(fname,['/',payload,'/anemometer_2/wind_gust']);
-  wgus_q2 = h5read(fname,['/',payload,'/anemometer_2/wind_gust_qc']);
+  wgus2 = ncread(fname,['/',payload,'/anemometer_2/wind_gust']);
+  wgus_q2 = ncread(fname,['/',payload,'/anemometer_2/wind_gust_qc']);
   jj = wgus_q2(ii) == 0;
   wgus(ii(jj)) = wgus2(ii(jj));
 end
-wgus_info = h5info(fname,['/',payload,'/anemometer_1/wind_gust']);
+wgus_info = ncinfo(fname,['/',payload,'/anemometer_1/wind_gust']);
 ii = wgus == wgus_info.FillValue;
 wgus(ii) = -99.99;
 
-wdir = h5read(fname,['/',payload,'/anemometer_1/wind_direction']);
-wdir_q = h5read(fname,['/',payload,'/anemometer_1/wind_direction_qc']);
+wdir = ncread(fname,['/',payload,'/anemometer_1/wind_direction']);
+wdir_q = ncread(fname,['/',payload,'/anemometer_1/wind_direction_qc']);
 ii = wdir_q ~= 0;
 wdir(ii) = -999.0;
 if sum(ii) > 0
-  wdir2 = h5read(fname,['/',payload,'/anemometer_2/wind_direction']);
-  wdir_q2 = h5read(fname,['/',payload,'/anemometer_2/wind_direction_qc']);
+  wdir2 = ncread(fname,['/',payload,'/anemometer_2/wind_direction']);
+  wdir_q2 = ncread(fname,['/',payload,'/anemometer_2/wind_direction_qc']);
   jj = wdir_q2(ii) == 0;
   wdir(ii(jj)) = wdir2(ii(jj));
 end
-wdir_info = h5info(fname,['/',payload,'/anemometer_1/wind_direction']);
+wdir_info = ncinfo(fname,['/',payload,'/anemometer_1/wind_direction']);
 ii = wdir == wdir_info.FillValue;
 wdir(ii) = -999.0;
 wdir = single(wdir);
 
 % Get Air temperature and Sea tempurature
-atemp = h5read(fname,['/',payload,'/air_temperature_sensor_1/air_temperature']);
-atemp_q = h5read(fname,['/',payload,'/air_temperature_sensor_1/air_temperature_qc']);
-atemp_info = h5info(fname,['/',payload,'/air_temperature_sensor_1/air_temperature']);
+atemp = ncread(fname,['/',payload,'/air_temperature_sensor_1/air_temperature']);
+atemp_q = ncread(fname,['/',payload,'/air_temperature_sensor_1/air_temperature_qc']);
+atemp_info = ncinfo(fname,['/',payload,'/air_temperature_sensor_1/air_temperature']);
 ii = atemp == atemp_info.FillValue;
 atemp = atemp - 272.15;
 atemp(ii) = -99.99;
 ii = atemp_q ~= 0;
 atemp(ii) = -99.99;
 try
     if sum(ii) > 0
-    atemp2 = h5read(fname,['/',payload,'/air_temperature_sensor_2/air_temperature']);
-    atemp_q2 = h5read(fname,['/',payload,'/air_temperature_sensor_2/air_temperature']);
+    atemp2 = ncread(fname,['/',payload,'/air_temperature_sensor_2/air_temperature']);
+    atemp_q2 = ncread(fname,['/',payload,'/air_temperature_sensor_2/air_temperature']);
     jj = atemp_q2(ii) == 0;
     atemp(ii(jj)) = atemp2(ii(jj));
     end
 catch
 end
 
 try
-    stemp = h5read(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature']);
-    stemp_q = h5read(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature_qc']);
-    stemp_info = h5info(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature']);
+    stemp = ncread(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature']);
+    stemp_q = ncread(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature_qc']);
+    stemp_info = ncinfo(fname,['/',payload,'/ocean_temperature_sensor_1/sea_surface_temperature']);
     ii = stemp == stemp_info.FillValue;
     stemp = stemp - 272.15;
     stemp(ii) = -99.99;
@@ -149,17 +152,17 @@
     stemp = repmat(-99.99,[length(times) 1]);
 end
 % Get barometric pressure from file
-barp = h5read(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level']);
-barp_q = h5read(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level_qc']);
-barp_info = h5info(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level']);
+barp = ncread(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level']);
+barp_q = ncread(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level_qc']);
+barp_info = ncinfo(fname,['/',payload,'/barometer_1/air_pressure_at_sea_level']);
 ii = barp == barp_info.FillValue;
 barp = double(barp)./100;
 barp(ii) = -999.;
 ii = barp_q ~= 0;
 barp(ii) = -999.;
 if sum(ii) ~= 0
-    barp2 = h5read(fname,['/',payload,'/barometer_2/air_pressure_at_sea_level']);
-    barp_q2 = h5read(fname,['/',payload,'/barometer_2/air_pressure_at_sea_level_qc']);
+    barp2 = ncread(fname,['/',payload,'/barometer_2/air_pressure_at_sea_level']);
+    barp_q2 = ncread(fname,['/',payload,'/barometer_2/air_pressure_at_sea_level_qc']);
     jj = barp_q2(ii) == 0;
     barp(ii(jj)) = barp2(ii(jj));
 end