Changeset 334
- Timestamp:
- 06/14/10 12:19:31
- Files:
-
- DPWP/trunk/DPWP/ADCP_splitter/adcpwaves.py (modified) (1 diff)
- DPWP/trunk/DPWP/adcp_matlab/radialtouvw.m (modified) (3 diffs)
- DPWP/trunk/DPWP/adcp_matlab/specmultiplot.m (modified) (8 diffs)
- DPWP/trunk/DPWP/adcp_matlab/waveplot.m (modified) (5 diffs)
- DPWP/trunk/DPWP/diwasp_1_1GD/dirspec.m (modified) (5 diffs)
- DPWP/trunk/DPWP/diwasp_1_1GD/private/IMLM.m (modified) (5 diffs)
- DPWP/trunk/DPWP/diwasp_1_1GD/private/radial.m (modified) (3 diffs)
- DPWP/trunk/DPWP/diwasp_1_1GD/private/velx.m (modified) (1 diff)
- DPWP/trunk/DPWP/diwasp_1_1GD/private/vely.m (modified) (1 diff)
- DPWP/trunk/DPWP/diwasp_1_1GD/private/velz.m (modified) (1 diff)
- DPWP/trunk/DPWP/nortek/nortek_radialtouvw.m (modified) (1 diff)
- DPWP/trunk/DPWP/nortek/nortek_specmultiplot.m (modified) (5 diffs)
- DPWP/trunk/DPWP/nortek/nortek_waveplot.m (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
- Modified
- Copied
- Moved
DPWP/trunk/DPWP/ADCP_splitter/adcpwaves.py
r315 r334 4 4 5 5 NOTE: Transducer height above bottom is hardwired into this code at .4m. 6 Need to update code in the future to be a user selected option""" 6 Need to update code in the future to be a user selected option. It is also 7 hardwired into the Matlab code radialtouvw.m and radial.m""" 7 8 8 9 DPWP/trunk/DPWP/adcp_matlab/radialtouvw.m
r314 r334 18 18 %orbit=load('orbital data'),sysinfo=load('sysinfo data') and range 19 19 20 %what is the transducer face height off the bottom? 20 %what is the transducer face height off the bottom 21 %NOTE: also hardwired into radial.m 21 22 adcpheight=0.4; 22 23 … … 251 252 ID.fs=2; 252 253 %depth 253 ID.depth= (0.4+mean(press)); %note changed for work with the FRF ADCP**********************254 ID.depth=avgdepth; 254 255 %the spectral matrix structure 255 256 SM.freqs=[0.01:0.01:0.4]; 256 SM.dirs=[0:2:360]; 257 SM.dirs=[0:2:360]; %note that SM.dirs is hardwired into specmultiplot (line 123) 257 258 SM.xaxisdir= 90; 258 259 %the estimation parameter 259 EP.method= 'EMEP'; 260 EP.iter=100; 260 EP.method= 'IMLM'; 261 EP.iter=50; 262 EP.nfft=256; 261 263 262 264 263 265 if data_type == 1 264 %For uvw and pressure for bins 2,3,4 266 %For uvw and pressure for the highest three bins (3,4,5) 267 265 268 % the datatypes 266 269 ID.datatypes={'pres' 'velx' 'vely' 'velz' 'velx' 'vely' 'velz' 'velx' 'vely' 'velz'}; … … 286 289 % the datatypes 287 290 ID.datatypes={'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial' 'radial'}; 288 % the layout using bins 2,3,4 289 %ID.layout = [xyzpositions(1,1,2) xyzpositions(1,2,2) xyzpositions(1,3,2) xyzpositions(1,4,2) xyzpositions(1,1,3) xyzpositions(1,2,3) xyzpositions(1,3,3) xyzpositions(1,4,3) xyzpositions(1,1,4) xyzpositions(1,2,4) xyzpositions(1,3,4) xyzpositions(1,4,4); 290 %xyzpositions(2,1,2) xyzpositions(2,2,2) xyzpositions(2,3,2) xyzpositions(2,4,2) xyzpositions(2,1,3) xyzpositions(2,2,3) xyzpositions(2,3,3) xyzpositions(2,4,3) xyzpositions(2,1,4) xyzpositions(2,2,4) xyzpositions(2,3,4) xyzpositions(2,4,4); 291 %xyzpositions(3,1,2) xyzpositions(3,2,2) xyzpositions(3,3,2) xyzpositions(3,4,2) xyzpositions(3,1,3) xyzpositions(3,2,3) xyzpositions(3,3,3) xyzpositions(3,4,3) xyzpositions(3,1,4) xyzpositions(3,2,4) xyzpositions(3,3,4) xyzpositions(3,4,4)]; 291 292 % the layout using highest 3 bins depending on the number of bins 293 % recorded 294 295 %highest bin 296 bin3=5; 297 %second highest 298 bin2=4; 299 %third highest 300 bin1=3; 301 302 ID.layout = [xyzpositions(1,1,bin3) xyzpositions(1,2,bin3) xyzpositions(1,3,bin3) xyzpositions(1,4,bin3) xyzpositions(1,1,bin2) xyzpositions(1,2,bin2) xyzpositions(1,3,bin2) xyzpositions(1,4,bin2) xyzpositions(1,1,bin1) xyzpositions(1,2,bin1) xyzpositions(1,3,bin1) xyzpositions(1,4,bin1); 303 xyzpositions(2,1,bin3) xyzpositions(2,2,bin3) xyzpositions(2,3,bin3) xyzpositions(2,4,bin3) xyzpositions(2,1,bin2) xyzpositions(2,2,bin2) xyzpositions(2,3,bin2) xyzpositions(2,4,bin2) xyzpositions(2,1,bin1) xyzpositions(2,2,bin1) xyzpositions(2,3,bin1) xyzpositions(2,4,bin1); 304 xyzpositions(3,1,bin3) xyzpositions(3,2,bin3) xyzpositions(3,3,bin3) xyzpositions(3,4,bin3) xyzpositions(3,1,bin2) xyzpositions(3,2,bin2) xyzpositions(3,3,bin2) xyzpositions(3,4,bin2) xyzpositions(3,1,bin1) xyzpositions(3,2,bin1) xyzpositions(3,3,bin1) xyzpositions(3,4,bin1)]; 292 305 % the data 293 %ID.data = horzcat(orbitnew(:,5),orbitnew(:,6),orbitnew(:,7),orbitnew(:,8),orbitnew(:,9),orbitnew(:,10),orbitnew(:,11),orbitnew(:,12),orbitnew(:,13),orbitnew(:,14),orbitnew(:,15),orbitnew(:,16)); 294 295 % the layout using bins 3,4,5 296 ID.layout = [xyzpositions(1,1,5) xyzpositions(1,2,5) xyzpositions(1,3,5) xyzpositions(1,4,5) xyzpositions(1,1,4) xyzpositions(1,2,4) xyzpositions(1,3,4) xyzpositions(1,4,4) xyzpositions(1,1,3) xyzpositions(1,2,3) xyzpositions(1,3,3) xyzpositions(1,4,3); 297 xyzpositions(2,1,5) xyzpositions(2,2,5) xyzpositions(2,3,5) xyzpositions(2,4,5) xyzpositions(2,1,4) xyzpositions(2,2,4) xyzpositions(2,3,4) xyzpositions(2,4,4) xyzpositions(2,1,3) xyzpositions(2,2,3) xyzpositions(2,3,3) xyzpositions(2,4,3); 298 xyzpositions(3,1,5) xyzpositions(3,2,5) xyzpositions(3,3,5) xyzpositions(3,4,5) xyzpositions(3,1,4) xyzpositions(3,2,4) xyzpositions(3,3,4) xyzpositions(3,4,4) xyzpositions(3,1,3) xyzpositions(3,2,3) xyzpositions(3,3,3) xyzpositions(3,4,3)]; 299 % the data 300 ID.data = horzcat(orbitnew(:,17),orbitnew(:,18),orbitnew(:,19),orbitnew(:,20),orbitnew(:,13),orbitnew(:,14),orbitnew(:,15),orbitnew(:,16),orbitnew(:,9),orbitnew(:,10),orbitnew(:,11),orbitnew(:,12)); 306 %set up values based on what bins are available 307 orb=orbOut; 308 ID.data = horzcat(orbitnew(:,orb-3),orbitnew(:,orb-2),orbitnew(:,orb-1),orbitnew(:,orb),orbitnew(:,orb-7),orbitnew(:,orb-6),orbitnew(:,orb-5),orbitnew(:,orb-4),orbitnew(:,orb-11),orbitnew(:,orb-10),orbitnew(:,orb-9),orbitnew(:,orb-8)); 301 309 302 310 DPWP/trunk/DPWP/adcp_matlab/specmultiplot.m
r314 r334 1 function [Winfo,spec2d,fd,thetad,td]=specmultiplot(dateinfo1,dateinfo2,timestep) ;1 function [Winfo,spec2d,fd,thetad,td]=specmultiplot(dateinfo1,dateinfo2,timestep) 2 2 %where dateinfo1 and 2 are the strings of the first and last date extension 3 3 %of the pressure,range,etc files and timestep is the time in hours between … … 8 8 %peak period, direction of peak period and dominant dir. It also generates 9 9 %polar plots and 2d freq and directional plots of a variety of different 10 %estimation methods, resolutions and the wavesmon output. 10 %estimation methods. 11 12 %%% IMPT NOTE!!! The ADCP height above bottom is manually set in the 13 %%% radialtouvw.m and the radial.m code. In this case the default is set 14 %%% at .4m This will be changed in a future release!! 11 15 12 16 … … 35 39 36 40 %set up data structure for wave info 37 Winfo.setup={' EMEPradial'};41 Winfo.setup={'IMLM radial'}; 38 42 Winfo.hsig=[]; 39 43 Winfo.hconf=[]; … … 64 68 65 69 % run diwasp to generate this spectrum with EMEP 66 [radialE, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0});70 [radialE,~]=dirspec(ID,SM,EP,{'MESSAGE',0,'PLOTTYPE',0}); 67 71 68 72 %Use waveplot to plot the polar, freq and dir plots of the spectra 69 73 70 [E_radial_WI ,fig1,fig2]=waveplot(radialE,sysinfo);74 [E_radial_WI]=waveplot(radialE); 71 75 72 saveas(fig1,['polar1_' dateinfo '.fig']); 73 saveas(fig2,['dirfreq1_' dateinfo '.fig']); 74 75 movefile(['polar1_' dateinfo '.fig'],'tempdir'); 76 movefile(['dirfreq1_' dateinfo '.fig'],'tempdir'); 77 close all; 78 79 76 %rename the files 77 movefile('tempdir/fig1.fig',['tempdir/polar1_' dateinfo '.fig']); 78 movefile('tempdir/fig2.fig',['tempdir/dirfreq1_' dateinfo '.fig']); 79 80 80 81 %Fill up the Winfo data structure 81 82 … … 93 94 Winfo.Ddir=[horzcat(Winfo.Ddir,Ddir)]; 94 95 Winfo.Spectrum.EMEPradial=[cat(3,Winfo.Spectrum.EMEPradial,radialE.S)]; 96 97 %save a continuing backup of Winfo in case of an error 98 save('backup_Winfo','Winfo'); 99 movefile('backup_Winfo.mat','tempdir'); 95 100 end 96 101 … … 159 164 %************************************************************************** 160 165 166 %first delete the backup_Winfo file 167 %delete('tempdir/backup_Winfo.mat'); 168 161 169 %now create the time series plots of the wave info 162 170 [fig_sigh,fig_pp,fig_dp,fig_dd]=Winfo_plot(Winfo); … … 167 175 saveas(fig_dd,['domdir_' dateinfo '.fig']); 168 176 177 close all; 178 169 179 movefile(['sigh_' dateinfo '.fig'],'tempdir'); 170 180 movefile(['peakp_' dateinfo '.fig'],'tempdir'); … … 174 184 movefile('tempdir',strcat('procdata_',dateinfo)); 175 185 176 close all;177 186 178 187 DPWP/trunk/DPWP/adcp_matlab/waveplot.m
r231 r334 1 function [E_radial_WI ,fig1,fig2]=waveplot(radialE,sysinfo)1 function [E_radial_WI]=waveplot(radialE) 2 2 3 3 %make sure to load the EMEP, IMLM and wavesmon samples … … 17 17 scrsz = get(0,'ScreenSize'); 18 18 fig1=figure('Position',[scrsz]); 19 subplot(1,1,1);20 19 subplotspec(radialE,4); 21 20 title('radial velocity data'); 21 22 saveas(fig1,'tempdir/fig1.fig'); 23 close all force; 22 24 23 25 … … 71 73 end 72 74 73 legend([h1],' EMEPradial','location','best');75 legend([h1],'IMLM radial','location','best'); 74 76 title('directional wave spectrum integrated over frequency'); 75 77 xlabel('axis angle (degrees true)'); … … 82 84 h2=plot(freqs,EMEPradialfreqUP,'b--'); 83 85 h2a=plot(freqs,EMEPradialfreqLOW,'b--'); 84 legend([h1,h2],' EMEPradial','95% confidence limits','location','best');86 legend([h1,h2],'IMLM radial','95% confidence limits','location','best'); 85 87 title('directional wave spectrum integrated over direction'); 86 88 xlabel('frequency in Hz'); 87 89 ylabel('m^2 / hz'); 90 91 saveas(fig2,'tempdir/fig2.fig'); 92 close all force; 88 93 89 94 % ______Calculate and display the wave parameters SigH, Tp, Dp,DTp_______ … … 110 115 111 116 %Display on the screen the SigH,Tp,Dp,DTp 112 disp([' EMEPradial']);117 disp(['IMLM radial']); 113 118 disp(['SigH (meters): ' num2str(EMEP_radial_Hsig)]); 114 119 disp(['SigH 95% confidence limits: ' num2str(EMEP_radial_HsConf)]); DPWP/trunk/DPWP/diwasp_1_1GD/dirspec.m
r314 r334 159 159 %have Ss = the spectral average if the radial velocities are used without a 160 160 %pressure or vertical range datatype first 161 lags=50; 161 162 compare=strcmp('radial',ID.datatypes(1)); 162 163 if compare == 1 … … 165 166 %the total number of sensors or szd 166 167 dataUsed=szd; 167 end 168 168 169 %also need to now calculate the effective number of windows to 170 %accurately calculate the degrees of freedom below. 171 Tau=ones(dataUsed); %Tau Matrix 172 Nstar=Tau; %Nstar Matrix 173 for i=1:dataUsed 174 for k=1:dataUsed 175 covi=xcov(data(:,i),lags,'coeff'); 176 covk=xcov(data(:,k),lags,'coeff'); 177 covik=xcov(data(:,i),data(:,k),lags,'coeff'); 178 covki=xcov(data(:,k),data(:,i),lags,'coeff'); 179 covicovk=covi.*covk; 180 covikcovki=covik.*covki; 181 %sign1=sign(covicovk); 182 %sign2=sign(covikcovki); 183 %part1=(abs(covicovk)).^.5; 184 %part2=(abs(covikcovki)).^.5; 185 Tau(i,k)=sum(covicovk+covikcovki); 186 Nstar=dataLength./Tau; 187 end 188 end 189 %calculate the average value of Nstar for the off-diagnal values 190 totalN=sum(sum(Nstar)); 191 for i=1:dataUsed 192 out(i)=Nstar(i,i); 193 end 194 autoOut=sum(out); %don't count the auto values 195 totalN=totalN-autoOut; 196 avgN=totalN/((dataUsed^2)-dataUsed);%find the average N value for all co values 197 Nstar=avgN; 198 else 199 %if there is only ONE time series used for the frequency spec 200 covi=xcov(data(:,1),lags,'coeff'); 201 Tau=sum(covi.*covi); 202 Nstar=dataLength./Tau; 203 end 204 205 %now find the effective number windows for each time series given Nstar 206 numWindows=Nstar/nfft; 169 207 170 208 … … 172 210 %See Emery and Thomson, 2004, added on 9/17/08 173 211 windowConstant=2.5164; %for Hamming Window (Matlab default for cpsd) 174 degF= dataUsed*2*numWindows*windowConstant;212 degF=2*numWindows*dataUsed*windowConstant; 175 213 SM.degF=degF; 176 214 … … 180 218 % call appropriate estimation function 181 219 disp(['directional spectra using' blanks(1) EP.method ' method']);disp(' ') 182 Specout=feval(EP.method,xps(:,:,1:ffs),trm(:,1:ffs,:),kx(:,:,1:ffs,:),Ss(:,1:ffs),pidirs,EP.iter,displ); 220 221 [Specout]=feval(EP.method,xps(:,:,1:ffs),trm(:,1:ffs,:),kx(:,:,1:ffs,:),Ss(:,1:ffs),pidirs,EP.iter,displ); 183 222 184 223 Specout(find(isnan(Specout)))=0.0; … … 196 235 197 236 %check Hsig of mapped spectrum and check sufficiently close to original 198 Hs2=Hsig(S,SM.freqs(2)-SM.freqs(1),SM.dirs(2)-SM.dirs(1)); 199 if (Hs2-Hs)/Hs >0.01 200 warning('User defined grid may be too coarse; try increasing resolution of ''SM.freqs'' or ''SM.dirs'''); 201 end 237 %*****As far as I can tell this warning is useless**** Edited out on 3/2/10 238 239 %Hs2=Hsig(S,SM.freqs(2)-SM.freqs(1),SM.dirs(2)-SM.dirs(1)); 240 %if (Hs2-Hs)/Hs >0.01 241 % warning('User defined grid may be too coarse; try increasing resolution of ''SM.freqs'' or ''SM.dirs'''); 242 %end 202 243 203 244 %smooth spectrum DPWP/trunk/DPWP/diwasp_1_1GD/private/IMLM.m
r168 r334 1 1 function [S]=IMLM(xps,trm,kx,Ss,pidirs,miter,displ) 2 2 3 gamma=0.1; 3 %this code was edited significantly on 5/2010 --- see comments for details 4 5 gamma=.1; 4 6 beta=1.0; 5 7 alpha=0.1; … … 32 34 end 33 35 34 35 invcps= inv(xps(:,:,ff));36 I=eye(m)*sqrt(eps); 37 invcps=pinv(xps(:,:,ff)+I); %changed on 5/5/10 to pinv and +I 36 38 Sftmp=zeros(ddirs,1); 37 39 for m=1:szd … … 41 43 end 42 44 end 43 Eo=(1./Sftmp(:)); 44 kappa=1./(ddir*sum(Eo)); 45 Eo=kappa*Eo; 45 Eo=(1./real(Sftmp(:))); %added 5/5/10 to take just the real part 46 [Eo]=normalize(Eo,ddir); %added 5/11/10 to remove negative values before normalizing 46 47 E=Eo; 47 48 T=Eo; … … 54 55 end 55 56 end 56 invcps= inv(ixps);57 invcps=pinv(ixps+I); % changed on 5/5/10 to pinv and +I 57 58 Sftmp=zeros(ddirs,1); 58 59 for m=1:szd … … 63 64 end 64 65 Told=T; 65 T=(1./Sftmp(:)); 66 kappa=1./(ddir*sum(T)); 67 T=kappa*T; 66 T=(1./real(Sftmp(:))); %added on 5/5/10 to take just real part 67 [T]=normalize(T,ddir); %added 5/11/10 to remove negative values before normalizing 68 68 69 69 %lambda=ones(size(T))-(T./Eo) 70 70 %ei=gamma*lambda.*E; 71 71 ei=gamma*((Eo-T)+alpha*(T-Told)); 72 %lambda=Eo-T; %edited on 5/5/10 to calculate ei slightly differently 73 %ei=(abs(lambda).^(beta+1))./(gamma.*lambda); 72 74 E=E+ei; 73 kappa=1./(ddir*sum(E)); 74 E=kappa*E; 75 [E]=normalize(E,ddir);%added 5/11/10 to remove negative values before normalizing 75 76 76 77 DPWP/trunk/DPWP/diwasp_1_1GD/private/radial.m
r231 r334 1 function [trm]=radial(ffreqs,dirs,wns,z,depth,xpos,ypos) ;1 function [trm]=radial(ffreqs,dirs,wns,z,depth,xpos,ypos) 2 2 3 3 %This transfer function takes the along beam radial velocities from an ADCP … … 21 21 22 22 % compute the axis angle 23 if xpos > 0 & ypos > 0 % 1st quadrant23 if xpos > 0 && ypos > 0 % 1st quadrant 24 24 B= acos(xpos/dist); 25 elseif xpos < 0 & ypos < 0 % 3 quadrant25 elseif xpos < 0 && ypos < 0 % 3 quadrant 26 26 B= pi + acos(abs(xpos/dist)); 27 elseif xpos < 0 & ypos > 0 % 2nd quadrant27 elseif xpos < 0 && ypos > 0 % 2nd quadrant 28 28 B= pi/2 + asin(abs(xpos/dist)); 29 elseif xpos > 0 & ypos < 0 % 4th quadrant29 elseif xpos > 0 && ypos < 0 % 4th quadrant 30 30 B = 3*pi/2 + asin(abs(xpos/dist)); 31 31 end … … 66 66 67 67 %now put them together to create the transfer matrix 68 trm=(FLKZ*(sin(a)*cos(dirs-B))-( i*FLKZi)*A1);68 trm=(FLKZ*(sin(a)*cos(dirs-B))-(1i*FLKZi)*A1); 69 69 70 70 DPWP/trunk/DPWP/diwasp_1_1GD/private/velx.m
r168 r334 1 1 function trm=velx(ffreqs,dirs,wns,z,depth,xpos,ypos) 2 2 3 Kz=cosh(z*wns)./sinh(depth*wns); 4 %include a maximum cuttoff for the velocity response function 5 Kz(find(Kz<0.1))=0.1; 6 Kz(find(isnan(Kz)))=1; 7 trm=(ffreqs.*Kz)*cos(dirs); 3 4 %Edits were made on 4/20/10 to fix a bug in the cutoff limits in the KZ 5 %parts of the transfer function. 6 7 %set the cutoff value (for high frequencies) 8 cutoff=0.1; 9 10 11 KZ=cosh(z*wns)./sinh(depth*wns); 12 FLKZ=ffreqs.*KZ; 13 14 %find the cutoff and replace the values below it 15 cutindex=find(FLKZ < cutoff); 16 FLKZ(cutindex)=cutoff; 17 18 trm=(FLKZ)*cos(dirs); DPWP/trunk/DPWP/diwasp_1_1GD/private/vely.m
r168 r334 1 1 function trm=vely(ffreqs,dirs,wns,z,depth,xpos,ypos) 2 2 3 Kz=cosh(z*wns)./sinh(depth*wns); 4 %include a maximum cuttoff for the velocity response function 5 Kz(find(Kz<0.1))=0.1; 6 Kz(find(isnan(Kz)))=1; 7 trm=(ffreqs.*Kz)*sin(dirs); 3 %Edits were made on 4/20/10 to fix a bug in the cutoff limits in the KZ 4 %parts of the transfer function. 5 6 %set the cutoff value (for high frequencies) 7 cutoff=0.1; 8 9 10 KZ=cosh(z*wns)./sinh(depth*wns); 11 FLKZ=ffreqs.*KZ; 12 13 %find the cutoff and replace the values below it 14 cutindex=find(FLKZ < cutoff); 15 FLKZ(cutindex)=cutoff; 16 17 18 trm=(FLKZ)*sin(dirs); DPWP/trunk/DPWP/diwasp_1_1GD/private/velz.m
r168 r334 1 1 function trm=velz(ffreqs,dirs,wns,z,depth,xpos,ypos) 2 2 3 Kz=sinh(z*wns)./sinh(depth*wns); 4 %include a maximum cuttoff for the velocity response function 5 Kz(find(Kz<0.1))=0.1; 6 Kz(find(isnan(Kz)))=1; 7 trm=-i*(ffreqs.*Kz)*ones(size(dirs)); 3 %Edits were made on 4/20/10 to fix a bug in the cutoff limits in the KZ 4 %parts of the transfer function. 5 6 %set the cutoff value (for high frequencies) 7 cutoff=0.1; 8 9 10 KZ=sinh(z*wns)./sinh(depth*wns); 11 FLKZ=ffreqs.*KZ; 12 13 %find the cutoff and replace the values below it 14 cutindex=find(FLKZ < cutoff); 15 FLKZ(cutindex)=cutoff; 16 17 18 trm=-1i*(FLKZ)*ones(size(dirs)); DPWP/trunk/DPWP/nortek/nortek_radialtouvw.m
r314 r334 171 171 %the spectral matrix structure 172 172 SM.freqs=[0.01:0.01:0.4]; 173 SM.dirs=[ -180:2:180];173 SM.dirs=[0:2:360]; 174 174 SM.xaxisdir= 90; 175 175 %the estimation parameter 176 176 EP.method= 'IMLM'; 177 EP.iter=3; 177 EP.iter=50; 178 EP.nfft=256; 178 179 179 180 if data_type == 1 DPWP/trunk/DPWP/nortek/nortek_specmultiplot.m
r314 r334 38 38 39 39 %set up data structure for wave info 40 Winfo.setup={' EMEP'}';40 Winfo.setup={'IMLM'}'; 41 41 Winfo.hsig=[]; 42 42 Winfo.hconf=[]; … … 44 44 Winfo.dirP=[]; 45 45 Winfo.Ddir=[]; 46 Winfo.Spectrum .IMLM=[];46 Winfo.Spectrum=[]; 47 47 48 48 %load the wave data from the nortek quick wave software … … 70 70 %For IMLM 71 71 EP.method='IMLM'; 72 EP.iter= 3;73 [I_F01_D2, EPout]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0});72 EP.iter=50; 73 [I_F01_D2, ~]=dirspec(ID,SM,EP,{'PLOTTYPE',0,'MESSAGE',0}); 74 74 75 75 … … 77 77 78 78 %[E_WI,I_WI,nortek_WI,fig1,fig2]=nortek_waveplot(E_F01_D2,I_F01_D2,sysinfo,nortekspec); 79 [I_WI,fig1,fig2]=nortek_waveplot 1(I_F01_D2);79 [I_WI,fig1,fig2]=nortek_waveplot(I_F01_D2); 80 80 81 81 saveas(fig1,['polar1_' dateinfo '.fig']); … … 100 100 Winfo.dirP=[horzcat(Winfo.dirP,dirP)]; 101 101 Winfo.Ddir=[horzcat(Winfo.Ddir,Ddir)]; 102 Winfo.Spectrum .IMLM=[cat(3,Winfo.Spectrum.IMLM,I_F01_D2.S)];102 Winfo.Spectrum=[cat(3,Winfo.Spectrum,I_F01_D2.S)]; 103 103 end 104 104 DPWP/trunk/DPWP/nortek/nortek_waveplot.m
r314 r334 14 14 freqs=[0.01:0.01:.4]; 15 15 dirres=2; 16 dirs=[ -180:2:180];16 dirs=[0:2:360]; 17 17 18 18 % plot the spectrum generated through DIWASP