Home > wafo > spec > private > time2spa.m

time2spa

PURPOSE

Transform of spectrum from frequency to wave no. (used in spec2spec)

SYNOPSIS

Sk=time2spa(S,k,k2,g,rate)

DESCRIPTION

``` TIME2SPA Transform of spectrum from frequency to wave no. (used in spec2spec)

CALL:  Sk = time2spa(S,k,k2,g)

Sk = wave number spectrum, 1- or 2D depending on input
S  = spectrum struct
k,k2= wave number lag of output (default depending on input spectrum)
g  = constant of gravity (default see gravity)
rate = defining the ratio between the wave numbers and angular
frequencies. (default 2)

Transforms a frequency or directional spectrum to
a wave number spectrum.
Input 1D gives output 1D, input 2D gives output 2D
Rotation, transformation etc is preserved

CROSS-REFERENCE INFORMATION

This function calls:
 gravity returns the constant acceleration of gravity k2w Translates from wave number to frequency w2k Translates from frequency to wave number error Display message and abort function. interp1 1-D interpolation (table lookup) interp2 2-D interpolation (table lookup). isfield True if field is in structure array. linspace Linearly spaced vector. rmfield Remove fields from a structure array. strcmpi Compare strings ignoring case.
This function is called by:
 spec2spec Transforms between different types of spectra

SOURCE CODE

```001 function Sk=time2spa(S,k,k2,g,rate)
002 %TIME2SPA Transform of spectrum from frequency to wave no. (used in spec2spec)
003 %
004 % CALL:  Sk = time2spa(S,k,k2,g)
005 %
006 %   Sk = wave number spectrum, 1- or 2D depending on input
007 %   S  = spectrum struct
008 %  k,k2= wave number lag of output (default depending on input spectrum)
009 %   g  = constant of gravity (default see gravity)
010 % rate = defining the ratio between the wave numbers and angular
011 %        frequencies. (default 2)
012 %
013 % Transforms a frequency or directional spectrum to
014 % a wave number spectrum.
015 % Input 1D gives output 1D, input 2D gives output 2D
016 % Rotation, transformation etc is preserved
017 %
019
020 % Tested on Matlab 5.3
021 % History:  revised by så 03.10.2005: corrected bug when input k=dk
022 %           revised by es 29.11.1999: norm. in both x and y (length(g) 1 OR 2)
023 %           by es 99.08.17
024
025 error(nargchk(1,5,nargin))
026 if nargin<5|isempty(rate), rate = 2; end
027 rate = max(round(abs(rate)),1);
028
029 if nargin<3
030   if isfield(S,'g')
031     g=S.g;
032   else
033     g=gravity;
034   end
035 end
036 if strcmpi(S.type(end-2:end),'k1d')|strcmpi(S.type(end-2:end),'k2d')
037   error('Spectrum already in space domain')
038 end
039 Sk=S;
040 if isfield(S,'f')
041   w=2*pi*S.f(:);
042   Sf=S.S/2/pi;
043   Sk=rmfield(Sk,'f');
044 else
045   w=S.w(:);
046   Sf=S.S;
047   Sk=rmfield(Sk,'w');
048 end
049 if isfield(Sk,'v')
050   Sk=rmfield(Sk,{'v','phi'});
051   % enc.velocity  and direction make no sense for wave no spectrum
052 end
053 Sk.g=g;
054
055 if strcmpi(S.type,'freq')|strcmpi(S.type,'enc')
056   if nargin<2|isempty(k)
057     k=linspace(0,w2k(w(end),0,S.h,g(1)),rate*length(w))';
058   else
059     if length(k)==1 % then interpret k as dk (step length)
060       k=((0:length(w)-1)*k)';
061     end
062   end
063   length(k);
064   wk = k2w(k,0,Sk.h,g(1));
065   in = (wk>min(w))&(wk<=w(end));
066   Gk = zeros(size(w));
067   Gk(in)=interp1(w,Sf,wk(in));
068   Sk.S=zeros(size(k'));
069   Sk.S=(Gk.*dwdk(wk,k,Sk.h,g(1)))';
070   Sk.type='k1D';
071   Sk.k=k';
072 else % directional spectrum
073   if S.h<inf
074     disp('This transformation for finite depth is not available yet')
075     return
076   end
077   if nargin<2|(isempty(k)&isempty(k2)) % no arg-in for wave-numbers
078     nw=max(rate-1,1)*length(w);
079     k=linspace(0,w(end)^2/g(1),nw);
080     if g(end)>0
081       k2=linspace(-w(end)^2/g(end), w(end)^2/g(end), 2*nw-1)';
082     else
083        k2=linspace(-w(end)^2/g(1), w(end)^2/g(1), 2*nw-1)';
084     end
085   else % Some arg-in for k
086     if nargin<3|isempty(k2) % no arg-in for k2
087       k2=k;
088     end
089     if length(k)==1 % k scalar, ie dk
090       k=(0:(nw-1))*k;
091     end
092     if length(k2)==1 % ditto for k2
093       k2=(-(nw-1):(nw-1))'*k2;
094     end
095   end
096   % k,k2 vectors. Make matrices
097   k=k(:)'; % make k row vector
098   if (S.theta(1)<-pi/2)|(S.theta(end)>pi/2);
099     %thetas on the left half plane => k's<0
100     k=[-fliplr(k(:,2:end)) k];
101   end
102   k2=k2(:); % make k2 column vector
103   K1=ones(size(k2))*k;
104   K2=k2*ones(size(k));
105   Sk.k=k;
106   Sk.k2=k2;
107
108   % Matrices K1 and K2 are equidistant
109   W=sqrt(sqrt(g(1)^2*K1.^2+g(end)^2*K2.^2)); % Non-equidistant W derived from (K1,K2)
110   TH=atan2(g(end)*K2,g(1)*K1); % => -pi < TH <= pi, Non-eq.dist. TH ----"----
111   Gk=zeros(size(W));
112   in=(W>w(1))&(W<=w(end))&(TH>S.theta(1))&(TH<=S.theta(end));
113   Gk(in)=interp2(w,S.theta,Sf,W(in),TH(in),'*linear');
114   Sk.S=zeros(size(Gk));
115   Sk.S(W>0)=Gk(W>0)./abs(W(W>0)).^3/2*g(1)*g(end);  % Wave-number spectrum
116   Sk=rmfield(Sk,'theta');
117   Sk.type='k2D';
118 end
119 return
120
121 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
122
123 function dw=dwdk(w,k,h,g)
124
125 dw=zeros(size(w));
126 if h==inf
127   dw(w>0)=g/2./w(w>0);
128 else
129   dw(w>0)=(g*tanh(k(w>0)*h)+g*h*k(w>0).*(1-tanh(k(w>0)*h).^2))./w(w>0)/2;
130 end
131 return
132
133
134
135
136
137
138
139
140
141```

Mathematical Statistics
Centre for Mathematical Sciences
Lund University with Lund Institute of Technology

Comments or corrections to the WAFO group

Generated on Thu 06-Oct-2005 02:21:16 for WAFO by m2html © 2003