% ttide_adcp_comb.m An mfile to run T_TIDE on South Carolina ADCP % velocity and pressure data. % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Use of this program is described in: % % Sullivan, C.M., Warner, J.C., Martini, M.A., Voulgaris, G., % Work, P.A., Haas, K.A., and Hanes, D.H. (2006) % South Carolina Coastal Erosion Study Data Report for Observations % October 2003 - April 2004., USGS Open-File Report 2005-1429. % % Program written in Matlab v7.1.0 SP3 % Program ran on PC with Windows XP Professional OS. % % "Although this program has been used by the USGS, no warranty, % expressed or implied, is made by the USGS or the United States % Government as to the accuracy and functioning of the program % and related program material nor shall the fact of distribution % constitute any such warranty, and no responsibility is assumed % by the USGS in connection therewith." %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % C. Sullivan 05/15/06, Water flow data in bins likely exposed to air has % been set to the NetCDF FillValue_. Editing code to % consider the surface bin the upward-most good bin. % Also ttide velocities and pressure at Site 4. Also % found and fixed typos that would have produced % weird results % C. Sullivan 02/10/05, Perform harmonic analysis on entire timeseries at % all sites with good ADCP data. Sites with good % data include Site 1 (velocities only), Site 2, % Site 3 (velocities only), and Site 7 % (velocities only). % C. Sullivan 01/11/05, This version use velocity data from surface, % bottom, and user-selected mid-depths. The pressure % used is the bottom pressure (P_4). The harmonic % analysis is run using the entire timeseries % thus combining the data from deployments 1 and 2. % We only run this version at site 2 b/c there was % good quality data during both deply. 1 and 2. more off % Define data and output directories dataDirectory = 'E:\CSULLIVAN\SOUTH_CAROLINA\DATA_DVD\DATAFILES'; plotDirectory = 'E:\CSULLIVAN\SOUTH_CAROLINA\PLOTS\TTIDE'; % Define ADCP data filenames. dataFiles={'7201adc-a.nc','7401adc-a.nc',... '7221adc-a.nc','7421adc-a.nc',... '7241adc-a.nc','7441adc-a.nc',... '7331adc-a.nc','7531adc-a.nc',... '7321adc-a.nc','7521wh-a.nc'}; % The following loop loads the entire timeseries of % pressure and velocities at site 2 and runs harmonic % analysis on the data. for f = 2:2:length(dataFiles) nc1=netcdf(fullfile(dataDirectory,dataFiles{f-1}), 'nowrite'); nc2=netcdf(fullfile(dataDirectory,dataFiles{f}), 'nowrite'); site=nc1.DESCRIPT(6); latitude=nc1.latitude(:); binDepths1=nc1{'depth'}(:); binDepths2=nc2{'depth'}(:); %time in julian days jd1=nc1{'time'}(:) + (nc1{'time2'}(:)/3600/1000/24); jd2=nc2{'time'}(:) + (nc2{'time2'}(:)/3600/1000/24); rawTime=[jd1;jd2]; %time intervals for deployments 1 and 2 should %be 15 minutes each interval1=mean(diff(jd1)); interval2=mean(diff(jd2)); if abs(interval1-interval2)>15/60/24 disp('Time intervals don''t match!') end interval=interval1; %will interpolate data onto array that is time in julian days interpTime=[rawTime(1):interval:rawTime(end)]; %bottom pressure converted to height using H=P/(rho*g) %where rho=1025 km/m^3 and g=9.8 m/s^2 P_41=nc1{'P_4'}(:); P_41(P_41>=1e35)=nan; Height1=P_41/(1025*9.8); P_42=nc2{'P_4'}(:); P_42(P_42>=1e35)=nan; Height2=P_42/(1025*9.8); rawHeight=[Height1;Height2]; %%%%%%run t_tide on combined bottom pressure data %%%%%%only at Site 2 and Site 4 if strcmp(site,'2') || strcmp(site,'4') %interpolate data [interpHeight]=smart_interp(rawTime,rawHeight,interpTime,3); startInd = find(~isnan(interpHeight),1,'first'); startTime = julian2datenum(interpTime(startInd)); interval = mean(diff(julian2datenum(interpTime)))*24; ttideOutFile = fullfile(plotDirectory,['ttide_pres_site',site,'.txt']); [nameu,fu,tidecon,xout]=t_tide(interpHeight,'start time',startTime,... 'interval',interval,... 'latitude',latitude,... 'output',ttideOutFile); end %%%%%%run t_tide on velocity data at the bottom ubot1 = [nc1{'u_1205'}(:,1)]; vbot1 = [nc1{'v_1206'}(:,1)]; bads=find(ubot1+vbot1 >= 1e35); ubot1(bads) = nan; vbot1(bads) = nan; ubot2 = [nc2{'u_1205'}(:,1)]; vbot2 = [nc2{'v_1206'}(:,1)]; bads=find(ubot2+vbot2 >= 1e35); ubot2(bads) = nan; vbot2(bads) = nan; raw_uvbot=[ubot1 + i*vbot1; ubot2 + i*vbot2]; %interpolate data [interp_uvbot]=smart_interp(rawTime,raw_uvbot,interpTime,3); startInd = find(~isnan(interp_uvbot),1,'first'); startTime = julian2datenum(interpTime(startInd)); interval = mean(diff(julian2datenum(interpTime)))*24; ttideOutFile = fullfile(plotDirectory,['ttide_curr_site',site,'_bottom.txt']); [nameu,fu,tidecon,xout]=t_tide(interp_uvbot,'start time',startTime,... 'interval',interval,... 'latitude',latitude,... 'output',ttideOutFile); %%%%%%run t_tide on velocity data at the surface %%%%%%C. Sullivan, 05/15/06: the surface is the upward-most good bin questBins1 = nc1{'u_1205'}.bins_questioned(:); lastGoodBin1 = questBins1(1)-1; depthLastGoodBin1 = binDepths1(lastGoodBin1); questBins2 = nc2{'u_1205'}.bins_questioned(:); lastGoodBin2 = questBins2(1)-1; depthLastGoodBin2 = binDepths2(lastGoodBin2); %%%%%%we know the last good bins between deployments at each site are %%%%%%less than 0.5 m (a bin width) apart usfc1 = [nc1{'u_1205'}(:,lastGoodBin1)]; vsfc1 = [nc1{'v_1206'}(:,lastGoodBin1)]; bads=find(usfc1+vsfc1 >= 1e35); usfc1(bads) = nan; vsfc1(bads) = nan; usfc2 = [nc2{'u_1205'}(:,lastGoodBin2)]; vsfc2 = [nc2{'v_1206'}(:,lastGoodBin2)]; bads=find(usfc2+vsfc2 >= 1e35); usfc2(bads) = nan; vsfc2(bads) = nan; raw_uvsfc=[usfc1 + i*vsfc1; usfc2 + i*vsfc2]; %interpolate data [interp_uvsfc]=smart_interp(rawTime,raw_uvsfc,interpTime,3); startInd = find(~isnan(interp_uvsfc),1,'first'); startTime = julian2datenum(interpTime(startInd)); interval = mean(diff(julian2datenum(interpTime)))*24; ttideOutFile = fullfile(plotDirectory,['ttide_curr_site',site,'_surface.txt']); [nameu,fu,tidecon,xout]=t_tide(interp_uvsfc,'start time',startTime,... 'interval',interval,... 'latitude',latitude,... 'output',ttideOutFile); %%%%%%run t_tide on velocity data at 4 mbs and 6 mbs targetDepths=[4,6]; %mbs for b=1:size(targetDepths,2) [index_num,dist_to_hit] = value2Index(binDepths1,targetDepths(b),0.5); u1 = [nc1{'u_1205'}(:,index_num)]; v1 = [nc1{'v_1206'}(:,index_num)]; bads=find(u1+v1 >= 1e35); u1(bads) = nan; v1(bads) = nan; [index_num,dist_to_hit] = value2Index(binDepths2,targetDepths(b),0.5); u2 = [nc2{'u_1205'}(:,index_num)]; v2 = [nc2{'v_1206'}(:,index_num)]; bads=find(u2+v2 >= 1e35); u2(bads) = nan; v2(bads) = nan; raw_uv=[u1 + i*v1; u2 + i*v2]; %interpolate data [interp_uv]=smart_interp(rawTime,raw_uv,interpTime,3); startInd = find(~isnan(interp_uv),1,'first'); startTime = julian2datenum(interpTime(startInd)); interval = mean(diff(julian2datenum(interpTime)))*24; ttideOutFile = fullfile(plotDirectory,['ttide_curr_site',site,'_',num2str(targetDepths(b)),'mbs.txt']); [nameu,fu,tidecon,xout]=t_tide(interp_uv,'start time',startTime,... 'interval',interval,... 'latitude',latitude,... 'output',ttideOutFile); end nc1=close(nc1); nc2=close(nc2); end