Home > wafo > spec > wnormspec.m

# wnormspec

## PURPOSE

Normalize a spectral density such that m0=m2=1

## SYNOPSIS

[S,mn4,m0,m2,m4,m1]=wnormspec(spectrum,m0,m2,plotflag)

## DESCRIPTION

``` WNORMSPEC Normalize a spectral density such that m0=m2=1

CALL:  [Sn,mn4,m0,m2,m4,m1]=wnormspec(S,m0,m2,plotflag);

Sn  = the normalized spectrum struct
mn4 = the fourth spectral moment of Sn.
mi  = the i'th spectral moment of S.
S   = the original spectrum struct
plotflag = 0, do not plot the normalized spectrum (default).
1, plot the normalized spectrum.

Normalization performed such that
INT S(freq) dfreq = 1       INT freq^2  S(freq) dfreq = 1
where integration limits are given by  freq  and  S(freq)  is the
spectral density; freq can be frequency or wave number.
The normalization is defined by
A=sqrt(m0/m2); B=1/A/m0; freq'=freq*A; S(freq')=S(freq)*B;

If S is a directional spectrum then a normalized gravity (.g) is added
to Sn, such that mxx normalizes to 1, as well as m0 and mtt.
(See spec2mom for notation of moments)

If S is complex-valued cross spectral density which has to be
normalized, then m0, m2 (suitable spectral moments) should be given.

Example:
S = jonswap;
[Sn,mn4] = wnormspec(S);
mts = spec2mom(S,2)     % Should be equal to one!

## CROSS-REFERENCE INFORMATION

This function calls:
 gravity returns the constant acceleration of gravity spec2spec Transforms between different types of spectra wspecplot Plot a spectral density datestr String representation of date. exist Check if variables or functions are defined. isfield True if field is in structure array. lower Convert string to lowercase. now Current date and time as date number. strcmp Compare strings. trapz Trapezoidal numerical integration.
This function is called by:
 dat2spec Estimate one-sided spectral density from data. dat2spec2 Estimate one-sided spectral density, version 2. spec2AcAt Evaluates survival function R(h1,h2)=P(Ac>h1,At>h2). spec2Acdf Evaluates cdf of crests P(Ac<=h) or troughs P(At<=h). spec2cmat Joint intensity matrix for cycles (max,min)-, rainflow- and (crest,trough) spec2mmtpdf Calculates joint density of Maximum, minimum and period. spec2sdat Simulates a Gaussian process and its derivative from spectrum spec2tccpdf Evaluates densities of wave period Tcc, wave lenght Lcc. spec2thpdf Joint density of amplitude and period/wave-length characteristics spec2tpdf Evaluates densities for crest-,trough-period, length. spec2tpdf2 Evaluates densities for various wave periods or wave lengths

## SOURCE CODE

```001 function [S,mn4,m0,m2,m4,m1]=wnormspec(spectrum,m0,m2,plotflag)
002 %WNORMSPEC Normalize a spectral density such that m0=m2=1
003 %
004 % CALL:  [Sn,mn4,m0,m2,m4,m1]=wnormspec(S,m0,m2,plotflag);
005 %
006 %        Sn  = the normalized spectrum struct
007 %        mn4 = the fourth spectral moment of Sn.
008 %        mi  = the i'th spectral moment of S.
009 %        S   = the original spectrum struct
010 %   plotflag = 0, do not plot the normalized spectrum (default).
011 %              1, plot the normalized spectrum.
012 %
013 % Normalization performed such that
014 %    INT S(freq) dfreq = 1       INT freq^2  S(freq) dfreq = 1
015 % where integration limits are given by  freq  and  S(freq)  is the
016 % spectral density; freq can be frequency or wave number.
017 % The normalization is defined by
018 % A=sqrt(m0/m2); B=1/A/m0; freq'=freq*A; S(freq')=S(freq)*B;
019 %
020 % If S is a directional spectrum then a normalized gravity (.g) is added
021 % to Sn, such that mxx normalizes to 1, as well as m0 and mtt.
022 % (See spec2mom for notation of moments)
023 %
024 % If S is complex-valued cross spectral density which has to be
025 % normalized, then m0, m2 (suitable spectral moments) should be given.
026 %
027 % Example:
028 %   S = jonswap;
029 %   [Sn,mn4] = wnormspec(S);
030 %   mts = spec2mom(S,2)     % Should be equal to one!
031 %
033
034 % Tested on: Matlab 5.3
035 %
036 % History:
037 % revised ir 31.08.01, introducing normalizing  g for w-spectrum
038 % Revised by jr 20.04.2001
039 % - Condition on nargout related to calculation of mn4; minor change
040 % - Updated information, added example
041 % By es 28.09.1999
042
043 if  nargin==2 & ~isempty(m0)
044   plotflag=m0;
045 end
046 if (nargin<4 & nargin ~=2)|isempty(plotflag)
047   if nargout==0,
048     plotflag=1;
049   else
050     plotflag=0;
051   end
052 end
053
054 if strcmp(lower(spectrum.type(end-2:end)),'k2d')
055   intype=spectrum.type;
056   spectrum=spec2spec(spectrum,'dir');
057 end
058
059 S=spectrum.S;   % size np x nf
060
061 if isfield(spectrum,'w')
062   f=spectrum.w(:); % length nf
063 elseif isfield(spectrum,'k')
064   f=spectrum.k(:);
065 else
066   f=2*pi*spectrum.f(:);
067   S=S/2/pi;
068 end
069
070 if (nargin<3)|isempty(m0)|isempty(m2)
071   S1=abs(S);
072   if strcmp(lower(spectrum.type(end-2:end)),'dir')
073     % integrate out theta:
074     S2=trapz(spectrum.theta(:),...
075          S1.*(cos(spectrum.theta(:)*ones(size(f'))).^2),1);
076
077     m20=trapz(f,f.^4/gravity^2.*S2.',1); % second order moment in x
078     m02=trapz(f,f.^4/gravity^2.*trapz(spectrum.theta(:),S1.*...
079     (sin(spectrum.theta(:)*ones(size(f'))).^2),1).',1); % second order moment in y
080     S1=trapz(spectrum.theta(:),S1,1).'; % integrate out theta
081     mom=trapz(f,[S1  f.*S1 f.^2.*S1 f.^4.*S1] ,1);
082     m0=mom(1); m1=mom(2); m2=mom(3);  m4=mom(4);
083   else
084     S1=S1(:);
085     mom=trapz(f,[S1  f.*S1 f.^2.*S1 f.^4.*S1] ,1);
086     m0=mom(1); m1=mom(2); m2=mom(3);  m4=mom(4);
087     if isfield(spectrum,'w')
088        m02=m4/gravity^2;
089         m20=m02;
090     end
091   end
092 end
093
094 SM0=sqrt(m0); SM2=sqrt(m2);
095 A=SM0/SM2;
096 B=SM2/(SM0*m0);
097
098 f=f*A;
099 S1=S*B;
100
101 if (nargout >= 2)&nargin<3
102   mn4 = m4*A^5*B;
103 elseif (nargout >= 2)& (nargin>3)
104   mn4=trapz(f,[f.^4.*S1] );
105 end
106
107 S=spectrum;
108 S.S=S1;
109 if isfield(spectrum,'w')
110   S.w=f;
111 elseif isfield(spectrum,'k')
112   S.k=f;
113 else
114   S.f=f/2/pi;
115   S.S=S.S*2*pi;
116 end
117 if strcmp(lower(spectrum.type(end-2:end)),'dir')
118   S.g=[gravity*sqrt(m0*m20)/m2 gravity*sqrt(m0*m02)/m2];
119   % normalization gravity in x and y, resp.
120 end
121 if isfield(spectrum,'w')
122   S.g=[gravity*sqrt(m0*m20)/m2];
123 end
124
125 S.norm=1;
126 S.date=datestr(now);
127 if exist('intype')
128   S=spec2spec(S,intype);
129 end
130 if plotflag
131   wspecplot(S)
132 end
133```

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