Home > wafo > wstats > dist2dsmfun.m

dist2dsmfun

PURPOSE ^

Smooths the conditional DIST2D distribution parameters.

SYNOPSIS ^

[Av , Bv , Cv ]=dist2dsmfun(phat,h2,CSMA,lin_extrap)

DESCRIPTION ^

 DIST2DSMFUN Smooths the conditional DIST2D distribution parameters.  
  
  CALL:  [A, B, C] = dist2dsmfun(phat,x2,csm,lin) 
  
    A, B, C = scale shape and location parameter, respectively 
              evaluated at x2 
    phat    = parameter structure array (see dist2dfit) 
    x2      = evaluation points (default phat.x{1}(:,1))  
    csm     = [csma csmb csmc], smoothing parameter vector which defines 
              the smoothing of parameter A,B and C, respectively  
              0 -> LS straight line 
              1 -> cubic spline interpolant (default [1 1 1]) 
    lin     = [lina linb linc], vector defining the extrapolation of 
              the parameters A,B and C,  respectively (default [1 1 1]) 
              0 No linear extrapolation outside the range of data 
              1 Linear extrapolation outside the range of data  
  
    The size of A B and C is the size of the argument  X2. 
     
  See also   dist2dfit, dist2dsmfun2, smooth

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

001 function [Av , Bv , Cv ]=dist2dsmfun(phat,h2,CSMA,lin_extrap) 
002 %DIST2DSMFUN Smooths the conditional DIST2D distribution parameters.  
003 % 
004 % CALL:  [A, B, C] = dist2dsmfun(phat,x2,csm,lin) 
005 % 
006 %   A, B, C = scale shape and location parameter, respectively 
007 %             evaluated at x2 
008 %   phat    = parameter structure array (see dist2dfit) 
009 %   x2      = evaluation points (default phat.x{1}(:,1))  
010 %   csm     = [csma csmb csmc], smoothing parameter vector which defines 
011 %             the smoothing of parameter A,B and C, respectively  
012 %             0 -> LS straight line 
013 %             1 -> cubic spline interpolant (default [1 1 1]) 
014 %   lin     = [lina linb linc], vector defining the extrapolation of 
015 %             the parameters A,B and C,  respectively (default [1 1 1]) 
016 %             0 No linear extrapolation outside the range of data 
017 %             1 Linear extrapolation outside the range of data  
018 % 
019 %   The size of A B and C is the size of the argument  X2. 
020 %    
021 % See also   dist2dfit, dist2dsmfun2, smooth 
022  
023  
024 % Tested on matlab 5.x 
025 % history 
026 % revised pab July2004 
027 %  -made useSig work correctl 
028 % revised pab 03.12.2000 
029 % - added truncated weibull and truncated rayleigh 
030 % revised pab 12.11.2000 
031 %  - added wggampdf option 
032 %  Per A. Brodtkorb 28.10.98 
033  
034 PVH   = phat.x{1}; 
035 ind   = find(~isnan(sum(PVH(:,:),2)));% & PVH(:,1)<2); 
036 PVH   = PVH(ind,:); 
037 CDIST = phat.dist{1};   
038 if nargin<2 |isempty(h2), 
039   h2=PVH(:,1); 
040 end 
041  
042 CSMB=1; 
043 CSMC=1; 
044 if nargin<3|isempty(CSMA) 
045   CSMA = 1; 
046 elseif length(CSMA)==2 
047   CSMB=CSMA(2); 
048   CSMA=CSMA(1); 
049 elseif length(CSMA)>=2 
050   CSMB=CSMA(2); 
051   CSMC=CSMA(3); 
052   CSMA=CSMA(1); 
053 end 
054 linB=1; 
055 linC=1; 
056 if (nargin<4) | isempty(lin_extrap), 
057     linA=1; 
058 elseif length(lin_extrap)==2, 
059   linA=lin_extrap(1); 
060   linB=lin_extrap(2); 
061 elseif length(lin_extrap)>=2, 
062   linA=lin_extrap(1); 
063   linB=lin_extrap(2); 
064   linC=lin_extrap(3); 
065 end 
066  
067 %smooth extrapolates linearly outside the range of PVH 
068 useSIG = (1 & isfield(phat,'CI')) ; % 
069 useSTATS = (useSIG & isfield(phat,'stats1') & ~isempty(phat.stats1)); 
070 if useSTATS 
071   Npts =  phat.stats1{end}(ind); 
072 end 
073 if useSIG 
074   useSIG = any(~isnan(phat.CI{1}(:))); 
075 end 
076 K5=1; 
077 if all(PVH(:,1)>=0) 
078   da = (PVH(:,1)).^K5.*ones(size(PVH(:,1))); 
079 else 
080   da = ones(size(PVH(:,1))); 
081 end 
082 if useSIG,  
083   CIa = phat.CI{1}(ind,1:2); 
084   if useSTATS 
085     NptsA = Npts; 
086   end 
087   k1 =  find(PVH(:,2)<CIa(:,1)| CIa(:,2)< PVH(:,2)); 
088   if any(k1) 
089     CIa(k1,:) = []; 
090     PVH(k1,:)=[]; 
091     ind(k1) = []; 
092     da(k1)  = []; 
093     if useSTATS 
094       NptsA(k1) = []; 
095     end 
096   end 
097   tmp = (diff(CIa,1,2)/4).^2; % approximate 2 standard 
098                                            % deviations 
099    if ~isempty(tmp) 
100       k0 = find(~isnan(tmp)); 
101       maxStd = 10; 
102       if any(k0) 
103     da(k0) = tmp(k0); 
104     maxStd = max(maxStd,10*max(tmp(k0))); 
105       end 
106       k1 = find(isnan(tmp)); 
107       if any(k1) 
108     da(k1) = maxStd; 
109       end 
110       %size(PVH),size(da) 
111        
112        da = da+[0; da(1:end-1)]+... 
113         [da(2:end);da(end)]+... 
114        PVH(:,1).*[0; abs(diff(PVH(:,2)))];; 
115        if useSTATS 
116      da = da./sqrt(NptsA); 
117        end 
118       da = PVH(:,1).*da./(PVH(:,2)+eps); 
119    end                      
120 end 
121  
122 Av=smooth(PVH(:,1),PVH(:,2),CSMA,h2,linA,da); 
123 switch lower(CDIST(1:2)), 
124   case {'ra','tr','tg','gu','ga','we','tw','gg'},  
125     if (any(Av<0)), 
126       da2 = abs(log(PVH(:,2)) - log(PVH(:,2)+ da)); 
127       Av=exp(smooth(PVH(:,1),log(PVH(:,2)),CSMA,h2,linA,da2)); 
128     end 
129 end 
130 Bv=[];       
131 if ~strcmp(lower(CDIST(1:2)),'ra') 
132   if all(PVH(:,1)>=0) 
133     db = PVH(:,1).^K5.*ones(size(PVH(:,1))); 
134   else 
135     db = ones(size(PVH(:,1))); 
136   end 
137   if useSTATS 
138     NptsB = Npts; 
139   end 
140   if useSIG  
141     CIb = phat.CI{1}(ind,3:4); 
142     k1 =  find(PVH(:,3)<CIb(:,1)| CIb(:,2)< PVH(:,3)); 
143     if any(k1) 
144       CIb(k1,:) = []; 
145       PVH(k1,:)=[]; 
146       ind(k1) = []; 
147       db(k1) = []; 
148       if useSTATS 
149     NptsB(k1) = []; 
150       end 
151     end 
152      
153     % approximate 2 standard deviations 
154     tmp=(diff(CIb,1,2)/4).^2;  
155     if ~isempty(tmp) 
156       k0 = find(~isnan(tmp)); 
157       maxStd = 10; 
158       if any(k0) 
159     db(k0) = tmp(k0); 
160     maxStd = max(maxStd,10*max(tmp(k0))); 
161       end 
162       k1 = find(isnan(tmp)); 
163       if any(k1) 
164     db(k1) = maxStd; 
165       end 
166       db = (db+[0; db(1:end-1)]+... 
167        [db(2:end);db(end)])/3+... 
168        PVH(:,1).*[0; abs(diff(PVH(:,3)))]; 
169       if useSTATS 
170      db = db./sqrt(NptsB); 
171        end 
172        db = PVH(:,1).*db./(PVH(:,3)+eps); 
173        %db = db./(PVH(:,3)+eps); 
174     end 
175   end 
176    
177   Bv=smooth(PVH(:,1),PVH(:,3),CSMB,h2,linB,db); 
178   switch CDIST(1:2) , 
179     case {'lo','ga','we','tw','gg'},   
180       if  (any(Bv<0)), 
181     db2 = abs(log(PVH(:,3)) - log(PVH(:,3)+ db)); 
182     Bv=exp(smooth(PVH(:,1),log(PVH(:,3)),CSMB,h2,linB,db2)); 
183       end 
184   end  
185 end 
186  
187 Cv=0; 
188 if (((size(PVH,2)>=4) & (~strcmpi(CDIST(1:2),'ra'))) | ... 
189     ((size(PVH,2)>=3) & (strcmpi(CDIST(1:2),'ra')))), 
190   Cv=smooth(PVH(:,1),PVH(:,4),CSMC,h2,linC); 
191    
192   if all(PVH(:,1)>=0) 
193     dc = PVH(:,1).^K5.*ones(size(PVH(:,1))); 
194   else 
195     dc = ones(size(PVH(:,1))); 
196   end 
197  
198   if useSIG  
199     CIc = phat.CI{1}(ind,5:6); 
200     if useSTATS 
201       NptsC = Npts; 
202     end 
203     k1 =  find(PVH(:,4)<CIc(:,1)| CIc(:,2)< PVH(:,4)); 
204     if any(k1) 
205       CIc(k1,:) = []; 
206       PVH(k1,:)=[]; 
207       ind(k1) = []; 
208       dc(k1) = []; 
209       NptsC(k1) = []; 
210     end 
211     tmp = (diff(CIc,1,2)/4).^2;  
212     % approximate 2 standard deviations 
213     if ~isempty(tmp) 
214       k0 = find(~isnan(tmp)); 
215       maxStd = 10; 
216       if any(k0) 
217     dc(k0) = tmp(k0); 
218     maxStd = max(maxStd,10*max(tmp(k0))); 
219       end 
220       k1 = find(isnan(tmp)); 
221       if any(k1) 
222     dc(k1) = maxStd; 
223       end 
224        dc = dc+[0; dc(1:end-1)]+... 
225         [dc(2:end);dc(end)];%+... 
226        %PVH(:,1).*[0; abs(diff(PVH(:,4)))];; 
227        if useSTATS 
228      dc = dc./sqrt(NptsC); 
229        end 
230       dc = PVH(:,1).*dc./(PVH(:,4)+eps); 
231    end                      
232 end 
233   switch lower(CDIST(1:2)) , 
234     case {'ra','tg','lo','ga','we'},   
235       if  (any(Cv<0)), 
236     Cv(Cv<0)=0; 
237       end 
238     case 'gg', 
239        if  (any(Cv<0)), 
240     Cv=exp(smooth(PVH(:,1),log(PVH(:,4)),CSMC,h2,linC,dc)); 
241       end  
242   end  
243 end 
244 return 
245  
246  
247

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