Home > wafo > trgauss > mc2rfc.m

# mc2rfc

## PURPOSE

Calculates a rainflow matrix given a Markov chain with kernel f_xy;

## SYNOPSIS

f_rfc = mc2rfc(f_xy,paramv,paramu)

## DESCRIPTION

``` MC2RFC  Calculates a rainflow matrix given a Markov chain with kernel f_xy;
f_rfc = f_xy + F_mc(f_xy).

CALL: f_rfc = mc2rfc(f_xy);

where

f_rfc = the rainflow matrix,
f_xy =  the frequency matrix of Markov chain (X0,X1)
but only the triangular part for X1>X0.

Further optional input arguments;

CALL:  f_rfc = mc2rfc(f_xy,paramx,paramy);

paramx = the parameter matrix defining discretization of x-values,
paramy = the parameter matrix defining discretization of y-values,```

## CROSS-REFERENCE INFORMATION

This function calls:
 levels Calculates discrete levels given the parameter matrix. axis Control axis scaling and appearance. clf Clear current figure. pcolor Pseudocolor (checkerboard) plot. rot90 Rotate matrix 90 degrees. subplot Create axes in tiled positions. title Graph title. xlabel X-axis label. ylabel Y-axis label.
This function is called by:
 iter_mc Calculates a kernel of a MC given a rainflow matrix

## SOURCE CODE

```001 function f_rfc = mc2rfc(f_xy,paramv,paramu)
002 %MC2RFC  Calculates a rainflow matrix given a Markov chain with kernel f_xy;
003 %       f_rfc = f_xy + F_mc(f_xy).
004 %
005 %  CALL: f_rfc = mc2rfc(f_xy);
006 %
007 %  where
008 %
009 %        f_rfc = the rainflow matrix,
010 %        f_xy =  the frequency matrix of Markov chain (X0,X1)
011 %                but only the triangular part for X1>X0.
012 %
013 %  Further optional input arguments;
014 %
015 %  CALL:  f_rfc = mc2rfc(f_xy,paramx,paramy);
016 %
017 %       paramx = the parameter matrix defining discretization of x-values,
018 %       paramy = the parameter matrix defining discretization of y-values,
019 %
020 if nargin<2
021 paramv=[-1, 1, length(f_xy)];
022 paramu=paramv;
023 end
024
025 if nargin<3
026 paramu=paramv;
027 end
028 dd=diag(rot90(f_xy));
029 N=length(f_xy);
030 Splus=sum(f_xy');
031 Sminus=fliplr(sum(f_xy));
032 Max_rfc=zeros(N,1);
033 Min_rfc=zeros(N,1);
034 norm=zeros(N,1);
035 for i=1:N
036   Spm=Sminus(i)+Splus(i)-dd(i);
037   if Spm>0.
038     Max_rfc(i)=(Splus(i)-dd(i))*(Splus(i)-dd(i))/(1-dd(i)/Spm)/Spm;
039     Min_rfc(i)=(Sminus(i)-dd(i))*(Sminus(i)-dd(i))/(1-dd(i)/Spm)/Spm;
040     norm(i)=Spm;
041   end
042 end
043
044 %cross=zeros(N,1);
045 %for i=2:N
046 %   cross(N-i+1)=cross(N-i+2)+Sminus(N-i+2)-Splus(N-i+2);
047 %end
048
049 f_rfc=zeros(N,N);
050 f_rfc(N-1,1)=Max_rfc(N-1);
051 f_rfc(1,N-1)=Min_rfc(2);
052
053 for k=3:N-1
054      for i=2:k-1,
055
056 %       AAe= f_xy(1:N-k,1:k-i);
057 %       SAAe=sum(sum(AAe));
058        AA = f_xy(N-k+1:N-k+i,k-i+1:k);
059        RAA=f_rfc(N-k+1:N-k+i,k-i+1:k);
060        nA=length(AA);
061        MA= Splus(N-k+1:N-k+i);
062        mA=Sminus(N-k+1:N-k+i);
063        normA=norm(N-k+1:N-k+i);
064        MA_rfc=Max_rfc(N-k+1:N-k+i);
065        mA_rfc=Min_rfc(k-i+1:k);
066        SA=sum(sum(AA));
067        SRA=sum(sum(RAA));
068        SMA_rfc=sum(MA_rfc);
069        SMA=sum(MA);
070        DRFC=SA-SMA-SRA+SMA_rfc;
071
072        NT=MA_rfc(1)-sum(RAA(1,:));
073
074 %       if k==35
075 %          check=[MA_rfc(1) sum(RAA(1,:))]
076 %          pause
077 %       end
078
079     NT=max(NT,0);
080
081     if NT>1e-6*MA_rfc(1)
082
083      NN=MA-sum(AA');
084      e=(fliplr(mA)-sum(AA))';
085      e=flipud(e);
086      AA=AA+flipud(rot90(AA,-1));
087      AA=rot90(AA);
088          AA=AA-0.5*diag(diag(AA));
089
090
091      for j=1:nA,
092           if normA(j)~=0
093             AA(j,:)=AA(j,:)/normA(j);
094             e(j)=e(j)/normA(j);
095           end
096        end
097        fx=0.;
098
099        if max(abs(e))>1e-7 & max(abs(NN))>1e-7*MA_rfc(1)
100
101
102        I=eye(size(AA));
103
104          if nA==1
105         fx=NN/(1-AA)*e;
106          else
107         fx=NN*((I-AA)\e);
108          end
109        end
110
111     f_rfc(N-k+1,k-i+1)=DRFC+fx;
112     end
113    end
114   m0=max(0,Min_rfc(N)-sum(f_rfc(N-k+2:N,1)));
115   M0=max(0,Max_rfc(N-k+1)-sum(f_rfc(N-k+1,2:k)));
116   f_rfc(N-k+1,1)=min(m0,M0);
117 %  n_loops_left=N-k+1
118  end
119
120 for k=2:N
121   M0=max(0,Max_rfc(1)-sum(f_rfc(1,N-k+2:N)));
122   m0=max(0,Min_rfc(k)-sum(f_rfc(2:k,N-k+1)));
123   f_rfc(1,N-k+1)=min(m0,M0);
124 end
125 f_rfc=f_rfc+rot90(diag(dd),-1);
126 clf
127 subplot(1,2,2)
128 pcolor(levels(paramv),levels(paramu),flipud(f_xy+flipud(rot90(f_xy,-1))))
129   axis([paramv(1), paramv(2), paramu(1), paramu(2)])
130   title('MC-kernel  f(x,y)')
131   ylabel('y'), xlabel('x')
132 axis('square')
133
134 subplot(1,2,1)
135 pcolor(levels(paramv),levels(paramu),flipud(f_rfc))
136   axis([paramv(1), paramv(2), paramu(1), paramu(2)])
137   title('Rainflow matrix')
138   ylabel('max'), xlabel('rfc-min')
139 axis('square')
140
141
142
143
144
145
146
147
148```

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