NCCOOS Trac Projects: Top | Web | Platforms | Processing | Viz | Sprints | Sandbox | (Wind)

Changeset 215

Show
Ignore:
Timestamp:
11/02/08 15:12:09
Author:
gdusek
Message:

--

Files:

Legend:

Unmodified
Added
Removed
Modified
Copied
Moved

  • DPWavesProc/trunk/DPWavesProc/RELEASENOTES.txt

    r40 r215  
     1This is the first BETA release of DPWavesProc. 

  • DPWavesProc/trunk/DPWavesProc/adcp_matlab/radialtouvw.m

    r168 r215  
    182182%transform the original x,y,z positions to the new positions accounting for 
    183183%heading, pitch and roll... we also add adcpheight back in 
    184 newxyzpos=ones(3,4,5); 
     184new_xyzpos=ones(3,4,5); 
    185185new_xyzpos(1,:,:)=xyzpos(1,:,:)*a+xyzpos(2,:,:)*b+xyzpos(3,:,:)*c; 
    186186new_xyzpos(2,:,:)=xyzpos(1,:,:)*d+xyzpos(2,:,:)*e+xyzpos(3,:,:)*f; 

  • DPWavesProc/trunk/DPWavesProc/adcp_matlab/radialwaveplot.m

    r168 r215  
    5858IMLMfreq=sum(real(radialI.S)')*dirres; 
    5959wmonfreq=sum(wmon.S')*wmondirres; 
     60 
     61%Compute the coefficient for the upper and lower error bounds for the 
     62%frequency 
     63%spectrum assuming 95% confidence.  Added on 9/17/08 
     64degF=radialE.degF; 
     65chiUp=chi2inv(.975,degF); 
     66chiLow=chi2inv(.025,degF); 
     67coeffUp=degF/chiLow; 
     68coeffLow=degF/chiUp; 
     69 
     70%calculate the conf limits throughout the frequency spectrum 
     71EMEPfreqUP=EMEPfreq*coeffUp; 
     72EMEPfreqLOW=EMEPfreq*coeffLow; 
    6073 
    6174%Find the maximum for the directional spectrum so we can set up the proper 
    …… 
    115128%plot the frequency energy spectrum 
    116129subplot(1,2,2); 
    117 plot(freqs,EMEPfreq,'b'); 
     130h1 = plot(freqs,EMEPfreq,'b'); 
    118131hold on 
    119 plot(wmon.freqs,wmonfreq,'k'); 
     132h2 = plot(freqs,EMEPfreqLOW,'b--'); 
     133h3 = plot(freqs,EMEPfreqUP,'b--'); 
     134h4 = plot(wmon.freqs,wmonfreq,'k'); 
    120135axis(axis); 
    121 plot(freqs,IMLMfreq,'c'); 
    122 legend('EMEP radial vel','wavesmon','IMLM radial vel','location','best'); 
     136%h5 = plot(freqs,IMLMfreq,'c'); 
     137legend([h1,h2,h4],'EMEP radial vel','EMEP radial 95% conf limits', 'wavesmon','location','best'); 
    123138title('directional wave spectrum integrated over direction'); 
    124139xlabel('frequency in Hz'); 
    …… 
    135150% Calculate the Sig wave height 
    136151EMEP_radial_Hsig=4*sqrt(m0); 
     152%Use the function HsigConf.m to calculate the sigH confidence limits 
     153EMEP_HsConf=HsigConf(radialE); 
    137154 
    138155% Calculate the peak period Tp 
    …… 
    151168disp(['EMEP radial vel']); 
    152169disp(['SigH (meters): ' num2str(EMEP_radial_Hsig)]); 
     170disp(['SigH 95% confidence limits: ' num2str(EMEP_HsConf)]); 
    153171disp(['peak period (seconds): ' num2str(EMEP_radial_Tp)]); 
    154172disp(['Dir of peak period: ' num2str(compangle(EMEP_radial_DTp, radialE.xaxisdir))]); 
    …… 
    157175 
    158176E_radial_WI.hsig=EMEP_radial_Hsig; 
     177E_radial_WI.hconf=EMEP_HsConf; 
    159178E_radial_WI.tp=EMEP_radial_Tp; 
    160179E_radial_WI.dtp=compangle(EMEP_radial_DTp, radialE.xaxisdir); 
    …… 
    169188% Calculate the Sig wave height 
    170189IMLM_radial_Hsig=4*sqrt(m0); 
     190%Use the function HsigConf.m to calculate the sigH confidence limits 
     191IMLM_HsConf=HsigConf(radialI); 
    171192 
    172193% Calculate the peak period Tp 
    …… 
    185206disp(['IMLM radial vel']); 
    186207disp(['SigH (meters): ' num2str(IMLM_radial_Hsig)]); 
     208disp(['SigH 95% confidence limits: ' num2str(IMLM_HsConf)]); 
    187209disp(['peak period (seconds): ' num2str(IMLM_radial_Tp)]); 
    188210disp(['Dir of peak period: ' num2str(compangle(IMLM_radial_DTp, radialI.xaxisdir))]); 
    …… 
    191213 
    192214I_radial_WI.hsig=IMLM_radial_Hsig; 
     215I_radial_WI.hconf=IMLM_HsConf; 
    193216I_radial_WI.tp=IMLM_radial_Tp; 
    194217I_radial_WI.dtp=compangle(IMLM_radial_DTp, radialI.xaxisdir); 

  • DPWavesProc/trunk/DPWavesProc/adcp_matlab/specmultiplot.m

    r168 r215  
    3535Winfo.setup={'EMEP UVW','IMLM UVW','EMEP range','IMLM range','EMEP radial','IMLM radial','wavesmon'}'; 
    3636Winfo.hsig=[]; 
     37Winfo.hconf=[]; 
    3738Winfo.peakP=[]; 
    3839Winfo.dirP=[]; 
    …… 
    6667 
    6768    %run diwasp to generate this spectrum 
    68     [E_uvw_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     69    [E_uvw_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    6970 
    7071    %now for IMLM 
    7172    EP.method='IMLM'; 
    7273    EP.iter=3; 
    73     [I_uvw_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     74    [I_uvw_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    7475 
    7576    %set up data with ranges, freq and dir at default 
    …… 
    7778 
    7879    % run diwasp to generate this spectrum with EMEP 
    79     [E_range_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     80    [E_range_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    8081 
    8182    %now with IMLM 
    8283    EP.method='IMLM'; 
    8384    EP.iter=3; 
    84     [I_range_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     85    [I_range_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    8586 
    8687    %set up data with radial velocities, freq and dir at default 
    …… 
    8889 
    8990    % run diwasp to generate this spectrum with EMEP 
    90     [E_radial_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     91    [E_radial_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    9192 
    9293    %now with IMLM 
    9394    EP.method='IMLM'; 
    9495    EP.iter=3; 
    95     [I_radial_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     96    [I_radial_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    9697 
    9798    %Use waveplot to plot the polar, freq and dir plots of the spectra 
    …… 
    115116 
    116117    hsig=[vertcat(E_UVW_WI.hsig,I_UVW_WI.hsig,E_range_WI.hsig,I_range_WI.hsig,E_radial_WI.hsig,I_radial_WI.hsig,wmon_WI.hsig)]; 
     118    hconf=[vertcat(E_UVW_WI.hconf,I_UVW_WI.hconf,E_range_WI.hconf,I_range_WI.hconf,E_radial_WI.hconf,I_radial_WI.hconf,[NaN NaN])]; 
    117119    peakP=[vertcat(E_UVW_WI.tp,I_UVW_WI.tp,E_range_WI.tp,I_range_WI.tp,E_radial_WI.tp,I_radial_WI.tp,wmon_WI.tp)]; 
    118120    dirP=[vertcat(E_UVW_WI.dtp,I_UVW_WI.dtp,E_range_WI.dtp,I_range_WI.dtp,E_radial_WI.dtp,I_radial_WI.dtp,wmon_WI.dtp)]; 
    …… 
    121123     
    122124    Winfo.hsig=[horzcat(Winfo.hsig,hsig)]; 
     125    Winfo.hconf=[horzcat(Winfo.hconf,hconf)]; 
    123126    Winfo.peakP=[horzcat(Winfo.peakP,peakP)]; 
    124127    Winfo.dirP=[horzcat(Winfo.dirP,dirP)]; 

  • DPWavesProc/trunk/DPWavesProc/adcp_matlab/waveplot.m

    r168 r215  
    7070IMLMrangefreq=sum(real(rangeI.S)')*dirres; 
    7171wmonfreq=sum(wmon.S')*wmondirres; 
     72 
     73%Compute the coefficient for the upper and lower error bounds for the power 
     74%spectrum assuming 95% confidence.  Added on 9/17/08 
     75degF=EMEP.degF; 
     76chiUp=chi2inv(.975,degF); 
     77chiLow=chi2inv(.025,degF); 
     78coeffUp=degF/chiLow; 
     79coeffLow=degF/chiUp; 
    7280 
    7381%Find the maximum for the directional spectrum so we can set up the proper 
    …… 
    148156 
    149157%For EMEP uvw 
    150  
    151 %calculate the 0,1,2 moments 
    152 m0=sum(EMEPfreq*freqres);                               
    153 m1=sum(freqs.*EMEPfreq*freqres); 
    154 m2=sum((freqs.^2).*EMEPfreq*freqres); 
    155 % Calculate the Sig wave height 
    156 EMEP_UVW_Hsig=4*sqrt(m0); 
     158% Calculate the Sig wave height 
     159EMEP_UVW_Hsig=Hsig(EMEP); 
     160%Use the function HsigConf.m to calculate the sigH confidence limits 
     161EMEP_HsConf=HsigConf(EMEP); 
    157162 
    158163% Calculate the peak period Tp 
    …… 
    171176disp(['EMEP uvw']); 
    172177disp(['SigH (meters): ' num2str(EMEP_UVW_Hsig)]); 
     178disp(['SigH 95% confidence limits: ' num2str(EMEP_HsConf)]); 
    173179disp(['peak period (seconds): ' num2str(EMEP_UVW_Tp)]); 
    174180disp(['Dir of peak period: ' num2str(compangle(EMEP_UVW_DTp, EMEP.xaxisdir))]); 
    …… 
    177183 
    178184E_UVW_WI.hsig=EMEP_UVW_Hsig; 
     185E_UVW_WI.hconf=EMEP_HsConf; 
    179186E_UVW_WI.tp=EMEP_UVW_Tp; 
    180187E_UVW_WI.dtp=compangle(EMEP_UVW_DTp, EMEP.xaxisdir); 
    …… 
    184191% For IMLM uvw 
    185192 
    186 %calculate the 0,1,2 moments 
    187 m0=sum(IMLMfreq*freqres);                               
    188 m1=sum(freqs.*IMLMfreq*freqres); 
    189 m2=sum((freqs.^2).*IMLMfreq*freqres); 
    190 % Calculate the Sig wave height 
    191 IMLM_UVW_Hsig=4*sqrt(m0); 
     193% Calculate the Sig wave height 
     194IMLM_UVW_Hsig=Hsig(IMLM); 
     195%Use the function HsigConf.m to calculate the sigH confidence limits 
     196IMLM_HsConf=HsigConf(IMLM); 
    192197 
    193198% Calculate the peak period Tp 
    …… 
    206211disp(['IMLM uvw']); 
    207212disp(['SigH (meters): ' num2str(IMLM_UVW_Hsig)]); 
     213disp(['SigH 95% confidence limits: ' num2str(IMLM_HsConf)]); 
    208214disp(['peak period (seconds): ' num2str(IMLM_UVW_Tp)]); 
    209215disp(['Dir of peak period: ' num2str(compangle(IMLM_UVW_DTp, IMLM.xaxisdir))]); 
    …… 
    212218 
    213219I_UVW_WI.hsig=IMLM_UVW_Hsig; 
     220I_UVW_WI.hconf=IMLM_HsConf; 
    214221I_UVW_WI.tp=IMLM_UVW_Tp; 
    215222I_UVW_WI.dtp=compangle(IMLM_UVW_DTp, IMLM.xaxisdir); 
    …… 
    218225% for Range of EMEP 
    219226 
    220 %calculate the 0,1,2 moments 
    221 m0=sum(EMEPrangefreq*freqres);                               
    222 m1=sum(freqs.*EMEPrangefreq*freqres); 
    223 m2=sum((freqs.^2).*EMEPrangefreq*freqres); 
    224 % Calculate the Sig wave height 
    225 EMEP_range_Hsig=4*sqrt(m0); 
     227% Calculate the Sig wave height 
     228EMEP_range_Hsig=Hsig(rangeE); 
     229%Use the function HsigConf.m to calculate the sigH confidence limits 
     230EMEP_range_HsConf=HsigConf(rangeE); 
    226231 
    227232% Calculate the peak period Tp 
    …… 
    240245disp(['EMEP range']); 
    241246disp(['SigH (meters): ' num2str(EMEP_range_Hsig)]); 
     247disp(['SigH 95% confidence limits: ' num2str(EMEP_range_HsConf)]); 
    242248disp(['peak period (seconds): ' num2str(EMEP_range_Tp)]); 
    243249disp(['Dir of peak period: ' num2str(compangle(EMEP_range_DTp, rangeE.xaxisdir))]); 
    …… 
    246252 
    247253E_range_WI.hsig=EMEP_range_Hsig; 
     254E_range_WI.hconf=EMEP_range_HsConf; 
    248255E_range_WI.tp=EMEP_range_Tp; 
    249256E_range_WI.dtp=compangle(EMEP_range_DTp, rangeE.xaxisdir); 
    …… 
    253260%for Range of IMLM 
    254261 
    255 %calculate the 0,1,2 moments 
    256 m0=sum(IMLMrangefreq*freqres);                               
    257 m1=sum(freqs.*IMLMrangefreq*freqres); 
    258 m2=sum((freqs.^2).*IMLMrangefreq*freqres); 
    259 % Calculate the Sig wave height 
    260 IMLM_range_Hsig=4*sqrt(m0); 
     262% Calculate the Sig wave height 
     263IMLM_range_Hsig=Hsig(rangeI); 
     264%Use the function HsigConf.m to calculate the sigH confidence limits 
     265IMLM_range_HsConf=HsigConf(rangeI); 
    261266 
    262267% Calculate the peak period Tp 
    …… 
    275280disp(['IMLM range']); 
    276281disp(['SigH (meters): ' num2str(IMLM_range_Hsig)]); 
     282disp(['SigH 95% confidence limits: ' num2str(IMLM_range_HsConf)]); 
    277283disp(['peak period (seconds): ' num2str(IMLM_range_Tp)]); 
    278284disp(['Dir of peak period: ' num2str(compangle(IMLM_range_DTp, rangeI.xaxisdir))]); 
    …… 
    281287 
    282288I_range_WI.hsig=IMLM_range_Hsig; 
     289I_range_WI.hconf=IMLM_range_HsConf; 
    283290I_range_WI.tp=IMLM_range_Tp; 
    284291I_range_WI.dtp=compangle(IMLM_range_DTp, rangeI.xaxisdir); 

  • DPWavesProc/trunk/DPWavesProc/diwasp_1_1GD/dirspec.m

    r168 r215  
    1 function [SM,EP]=dirspec(ID,SM,EP,varargin); 
     1function [SM,EP]=dirspec(ID,SM,EP,varargin) 
    22 
    33%DIWASP V1.1 function with GD edit 
    …… 
    8484end 
    8585 
    86 % NOTE: changed from csd to cpsd 6/15/07 
    87 %calculate the cross-power spectra 
     86%NOTE: changed on 9/5/08 window length for cpsd set to be a window equal to 
     87%or slightly larger than nfft so padding is limited 
     88dataLength=size(ID.data,1); 
     89numWindows=ceil(dataLength/nfft); 
     90window=ceil(dataLength/numWindows); 
     91 
     92%Do a version check for the upcoming cpsd function. If the version is 2006b 
     93%or later (date Aug,3,2006 or datenum 732892 cpsd will be handled 
     94%differently.  added on 9/17/08 
     95[v d]=version; 
     96ver_date=datenum(d); 
     97 
     98 
     99% NOTE: changed cross power spectral function from csd to cpsd 6/15/07 
    88100for m=1:szd 
    89101for n=1:szd 
    90    [xpstmp,Ftmp]=cpsd(data(:,m),data(:,n),[],[],nfft,ID.fs); 
     102    % Need to have 2 different cpsd functions depending on Matlab version. 
     103    % added of 9/17/08 
     104    if ver_date < 732892 
     105        [xpstmp,Ftmp]=cpsd(data(:,m),data(:,n),window,[],nfft,ID.fs); 
     106    else 
     107        [xpstmp,Ftmp]=cpsd(data(:,n),data(:,m),window,[],nfft,ID.fs); 
     108    end; 
    91109   xps(m,n,:)=xpstmp(2:(nfft/2)+1); 
    92  end 
     110end 
    93111end 
    94112 
    …… 
    119137end 
    120138 
     139%What are the number of sensors (data series) used to compute the power 
     140%spectrum?  If we aren't using radial velocities it is just 1. 
     141dataUsed=1; 
     142 
     143 
     144%NOTE: edited on 8/07 
    121145%have Ss = the spectral average if the radial velocities are used without a 
    122146%pressure or vertical range datatype first 
    123147compare=strcmp('radial',ID.datatypes(1)); 
    124  
    125148if compare == 1 
    126149    Ss=mean(Ss); 
    127 end 
    128      
     150    %since this is using radial velocities the number of data series used 
     151    %the total number of sensors or szd 
     152    dataUsed=szd; 
     153end 
     154 
     155 
     156 
     157%Compute the degrees of freedom depending on cpsd options and the data used 
     158%See Emery and Thomson, 2004, added on 9/17/08 
     159windowConstant=2.5164; %for Hamming Window (Matlab default for cpsd) 
     160degF=dataUsed*2*numWindows*windowConstant; 
     161SM.degF=degF; 
     162 
     163 
    129164ffs=sum(F<=max(SM.freqs)); 
    130165 

  • DPWavesProc/trunk/DPWavesProc/diwasp_1_1GD/infospec.m

    r168 r215  
    1 function [H,Tp,DTp,Dp]=infospec(SM) 
     1function [H,HsConf,Tp,DTp,Dp]=infospec(SM) 
     2 
     3%Edited on 9/17/08 to include confidence limits for the sig wave height 
     4 
    25 
    36%DIWASP V1.1 function 
    …… 
    811%Outputs: 
    912%Hsig           Signficant wave height (Hmo) 
     13%HsConf     95% conf limits for Hsig 
    1014%Tp                     Peak period 
    1115%DTp            Direction of spectral peak 
    …… 
    1418%Inputs: 
    1519%SM             A spectral matrix structure containing the file data 
     20% 
    1621% 
    1722%Hsig is the significant wave height. Tp is the peak frequency, the highest point in the one dimensional spectrum.  
    …… 
    2530SM=check_data(SM,2);if isempty(SM) return;end; 
    2631 
    27 H=HSig(SM); 
     32H=Hsig(SM); 
     33HsConf=HsigConf(SM); 
    2834 
    2935S=sum(real(SM.S),2); % S is freq spectrum 
    …… 
    3844disp(['Infospec::']); 
    3945disp(['Significant wave height (Hmo): ' num2str(H)]); 
     46disp(['Sig wave height 95% confidence limits: ' num2str(HsConf)]); 
    4047disp(['Peak period: ' num2str(Tp)]); 
    4148disp(['Direction of peak period: ' num2str(DTp) ' axis angle / ' num2str(compangle(DTp,SM.xaxisdir)) ' compass bearing']); 

  • DPWavesProc/trunk/DPWavesProc/diwasp_1_1GD/private/Hsig.m

    r168 r215  
    1 function Hs=Hmo(varargin) 
     1function Hs=Hsig(varargin) 
    22%DIWASP V1.1 function to calculate significant wave height 
    33% 

  • DPWavesProc/trunk/DPWavesProc/nortek/nortek_Winfo_plot.m

    r168 r215  
    1010h1=plot(Winfo.time,Winfo.hsig(1,:),'b.-'); 
    1111hold on 
     12h1a=plot(Winfo.time,Winfo.hconf(1,1:2:end),'b--'); 
     13plot(Winfo.time,Winfo.hconf(1,2:2:end),'b--'); 
    1214h2=plot(Winfo.time,Winfo.hsig(2,:),'r.-'); 
     15h2a=plot(Winfo.time,Winfo.hconf(2,1:2:end),'r--'); 
     16plot(Winfo.time,Winfo.hconf(2,2:2:end),'r--'); 
    1317h3=plot(Winfo.time,Winfo.hsig(3,:),'k.-'); 
    1418 
    1519datetick('x',15,'keepticks'); 
    16 legend([h1, h2, h3],'EMEP','IMLM','nortek','location','best'); 
     20legend([h1,h1a,h2,h2a,h3],'EMEP','EMEP 95% confidence','IMLM','IMLM 95% confidence','nortek','location','best'); 
    1721title(['Significant Wave Height starting on ',datestr(Winfo.time(:,1),0)]); 
    1822xlabel('time'); 

  • DPWavesProc/trunk/DPWavesProc/nortek/nortek_specmultiplot.m

    r168 r215  
    1111%peak period, direction of peak period and dominant dir.  It also generates 
    1212%polar plots and 2d freq and directional plots of a variety of different 
    13 %estimation methods, resolutions and the wavesmon output.  
     13%estimation methods, resolutions and the nortek quickwave output.  
    1414 
    1515 
    …… 
    3838Winfo.setup={'EMEP','IMLM' 'nortek'}'; 
    3939Winfo.hsig=[]; 
     40Winfo.hconf=[]; 
    4041Winfo.peakP=[]; 
    4142Winfo.dirP=[]; 
    …… 
    6869 
    6970    %run diwasp to generate this spectrum 
    70     [E_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     71    [E_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    7172 
    7273    %now for IMLM 
    7374    EP.method='IMLM'; 
    7475    EP.iter=3; 
    75     [I_F01_D2, EPout]=dirspectest(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
     76    [I_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 
    7677 
    7778 
    …… 
    8283    saveas(fig1,['polar1_' dateinfo '.fig']); 
    8384    saveas(fig2,['dirfreq1_' dateinfo '.fig']); 
     85     
     86    close all; 
    8487    
    85   
    8688    %Fill up the Winfo data structure 
    8789 
    8890    hsig=[vertcat(E_WI.hsig,I_WI.hsig,nortek_WI.hsig)]; 
     91    hconf=[vertcat(E_WI.hconf,I_WI.hconf,nortek_WI.hconf)]; 
    8992    peakP=[vertcat(E_WI.tp,I_WI.tp,nortek_WI.tp)]; 
    9093    dirP=[vertcat(E_WI.dtp,I_WI.dtp,nortek_WI.dtp)]; 
    …… 
    9295     
    9396    Winfo.hsig=[horzcat(Winfo.hsig,hsig)]; 
     97    Winfo.hconf=[horzcat(Winfo.hconf,hconf)]; 
    9498    Winfo.peakP=[horzcat(Winfo.peakP,peakP)]; 
    9599    Winfo.dirP=[horzcat(Winfo.dirP,dirP)]; 

  • DPWavesProc/trunk/DPWavesProc/nortek/nortek_waveplot.m

    r168 r215  
    5454IMLMfreq=sum(real(IMLM.S)')*dirres; 
    5555nortekfreq=sum(nortek.S')*nortekdirres; 
     56 
     57%Compute the coefficient for the upper and lower error bounds for the 
     58%frequency 
     59%spectrum assuming 95% confidence.  Added on 9/17/08 
     60degF=2*EMEP.degF; % 2* since for the Nortek the range beam sampled at 4hz (4096) is decimated or averaged to 2hz (2048) 
     61chiUp=chi2inv(.975,degF); 
     62chiLow=chi2inv(.025,degF); 
     63coeffUp=degF/chiLow; 
     64coeffLow=degF/chiUp; 
     65 
     66%calculate the conf limits throughout the frequency spectrum 
     67EMEPfreqUP=EMEPfreq*coeffUp; 
     68EMEPfreqLOW=EMEPfreq*coeffLow; 
     69IMLMfreqUP=IMLMfreq*coeffUp; 
     70IMLMfreqLOW=IMLMfreq*coeffLow; 
    5671 
    5772%Find the maximum for the directional spectrum so we can set up the proper 
    …… 
    110125%plot the frequency energy spectrum 
    111126subplot(1,2,2); 
    112 plot(freqs,EMEPfreq,'b'); 
     127h1=plot(freqs,EMEPfreq,'b'); 
    113128hold on 
    114 plot(freqs,IMLMfreq,'r'); 
    115 plot(nortekfreqs,nortekfreq,'k'); 
    116 legend('EMEP','IMLM','nortek quickwave','location','best'); 
     129h2=plot(freqs,EMEPfreqUP,'b--'); 
     130plot(freqs,EMEPfreqLOW,'b--'); 
     131h3=plot(freqs,IMLMfreq,'r'); 
     132h4=plot(freqs,IMLMfreqUP,'r--'); 
     133plot(freqs,IMLMfreqLOW,'r--'); 
     134h5=plot(nortekfreqs,nortekfreq,'k'); 
     135legend('EMEP','EMEP 95% confidence','IMLM','IMLM 95% confidence','nortek quickwave','location','best'); 
    117136title('directional wave spectrum integrated over direction'); 
    118137xlabel('frequency in Hz'); 
    …… 
    129148% Calculate the Sig wave height 
    130149EMEP_Hsig=4*sqrt(m0); 
     150%Use the function HsigConf.m to calculate the sigH confidence limits 
     151EMEP_HsConf=nortek_HsigConf(EMEP); 
    131152 
    132153% Calculate the peak period Tp 
    …… 
    145166disp(['EMEP']); 
    146167disp(['SigH (meters): ' num2str(EMEP_Hsig)]); 
     168disp(['SigH 95% conf. (meters): ' num2str(EMEP_HsConf)]); 
    147169disp(['peak period (seconds): ' num2str(EMEP_Tp)]); 
    148170disp(['Dir of peak period: ' num2str(compangle(EMEP_DTp, EMEP.xaxisdir))]); 
    …… 
    151173 
    152174E_WI.hsig=EMEP_Hsig; 
     175E_WI.hconf=EMEP_HsConf; 
    153176E_WI.tp=EMEP_Tp; 
    154177E_WI.dtp=compangle(EMEP_DTp, EMEP.xaxisdir); 
    …… 
    164187% Calculate the Sig wave height 
    165188IMLM_Hsig=4*sqrt(m0); 
     189%Use the function HsigConf.m to calculate the sigH confidence limits 
     190IMLM_HsConf=nortek_HsigConf(IMLM); 
    166191 
    167192% Calculate the peak period Tp 
    …… 
    180205disp(['IMLM']); 
    181206disp(['SigH (meters): ' num2str(IMLM_Hsig)]); 
     207disp(['SigH 95% conf. (meters): ' num2str(IMLM_HsConf)]); 
    182208disp(['peak period (seconds): ' num2str(IMLM_Tp)]); 
    183209disp(['Dir of peak period: ' num2str(compangle(IMLM_DTp, IMLM.xaxisdir))]); 
    …… 
    186212 
    187213I_WI.hsig=IMLM_Hsig; 
     214I_WI.hconf=IMLM_HsConf; 
    188215I_WI.tp=IMLM_Tp; 
    189216I_WI.dtp=compangle(IMLM_DTp, IMLM.xaxisdir); 
    …… 
    222249 
    223250nortek_WI.hsig=nortek_Hsig; 
     251nortek_WI.hconf=[NaN NaN]; 
    224252nortek_WI.tp=nortek_Tp; 
    225253nortek_WI.dtp=compangle(nortek_DTp, nortek.xaxisdir); 

  • DPWavesProc/trunk/DPWavesProc/radialtouvw1.m

    r168 r215  
    1 function [ID,SM,EP,orbitnew]=radialtouvw(pressure,range,orbit,sysinfo,data_type) 
     1function [ID,SM,EP,orbitnew]=radialtouvw1(pressure,range,orbit,sysinfo,data_type) 
    22 
    33