TRANGOOD Makes a transformation that is suitable for efficient transforms. CALL: f = trangood(ff,nmin,mini,maxi,nmax); f = the good transform function, [X f(X)]. ff = the input transform function, [X f(X)]. nmin = the minimum number of points in the good transform. (Default size(ff,1)) mini = the minimum data value to transform. (Default min(ff(:,1))) maxi = the maximum data value to transform. (Default max(ff(:,1))) nmax = then maximum number of points in the good transform (default inf) TRANGOOD interpolates ff linearly and optionally extrapolate ff linearly outside the range of ff(:,1) with X uniformly spaced. See also tranproc, interp1q
Clear variables and functions from memory. | |
Difference and approximate derivative. | |
Display message and abort function. | |
Quick 1-D linear interpolation. |
reconstruct the spurious/missing points of timeseries | |
Transforms process X and up to four derivatives |
001 function f = trangood(ff,nmin,mini,maxi,nmax) 002 %TRANGOOD Makes a transformation that is suitable for efficient transforms. 003 % 004 % CALL: f = trangood(ff,nmin,mini,maxi,nmax); 005 % 006 % f = the good transform function, [X f(X)]. 007 % ff = the input transform function, [X f(X)]. 008 % nmin = the minimum number of points in the good transform. 009 % (Default size(ff,1)) 010 % mini = the minimum data value to transform. 011 % (Default min(ff(:,1))) 012 % maxi = the maximum data value to transform. 013 % (Default max(ff(:,1))) 014 % nmax = then maximum number of points in the good transform 015 % (default inf) 016 % 017 % TRANGOOD interpolates ff linearly and optionally 018 % extrapolate ff linearly outside the range of ff(:,1) 019 % with X uniformly spaced. 020 % 021 % See also tranproc, interp1q 022 023 024 % History: 025 % revised pab 07.02.2001 026 % - added nmax 027 % - replaced x = linspace(f(1,1),f(nf,1),nmin)' with x = (f(1,1):df:f(nf,1))' 028 % revised pab 07.01.2001 029 % -fixed a bug: x = linspace(f(1,1),f(nf,1),nmin)' is safer than using 030 % x=f(1,1)+(0:nmin-1)'/(nmin-1)*(f(nf,1)-f(1,1)); with interp1q 031 % revised pab 12.11.2000 032 % - updated header info: A more detailed description of what TRANGOOD does. 033 % - changed interpolation with a call to interp1q which is much faster 034 % added nargchk and isempty(maxi),....isempty(nmin) 035 % by ??? 036 037 error(nargchk(1,5,nargin)) 038 if (size(ff,2)~=2) 039 error('ff must be a two column matrix.') 040 end 041 if (size(ff,1)<2) 042 error('ff must have at least two rows.') 043 end 044 045 [f,i] = sort(ff(:,1)); 046 f = [f ff(i,2)]; 047 clear i; 048 df = diff(f(:,1)); 049 if ( any(df<=0)), %eps 050 error('Duplicate x-values in ff not allowed.') 051 end 052 053 nf = size(f,1); 054 if (nargin<5)|isempty(nmax), nmax = inf; end 055 if (nargin<4)|isempty(maxi), maxi = f(nf,1); end 056 if (nargin<3)|isempty(mini), mini = f(1,1); end 057 if (nargin<2)|isempty(nmin), nmin = nf; end 058 if (nmin<2), nmin = 2; end 059 if (nmax<2), nmax = 2; end 060 061 if ( (nf<nmin) |(nmax<nf) | any(abs(diff(df))>10*eps*(f(nf,1)-f(1,1))) ) 062 % pab 07.01.2001: Always choose the stepsize df so that 063 % it is an exactly representable number. 064 % This is important when calculating numerical derivatives and is 065 % accomplished by the following. 066 df = (f(nf,1)-f(1,1))/(min(nmin,nmax)-1); 067 df = donothing(df+2)-2; 068 x = (f(1,1):df:f(nf,1)).'; 069 % New call pab 11.11.2000: This is much quicker 070 f = [ x interp1q(f(:,1),f(:,2),x)]; 071 %f = [ x interp1(f(:,1),f(:,2),x,'linear')]; 072 end 073 % f(:,1) is now uniformly spaced 074 df = f(2,1)-f(1,1); 075 076 % Extrapolate linearly outside the range of ff 077 %---------------------------------------------- 078 if (mini<f(1,1)), 079 f1 = df*(floor((mini-f(1,1))/df):1:-1)'; 080 f2 = f(1,2)+f1*(f(2,2)-f(1,2))/(f(2,1)-f(1,1)); 081 f = [f1+f(1,1) f2;f]; 082 end 083 n = size(f,1); 084 if (maxi>f(n,1)) 085 f1 = df*(1:1:ceil((maxi-f(n,1))/df))'; 086 f2 = f(n,2)+f1*(f(n,2)-f(n-1,2))/(f(n,1)-f(n-1,1)); 087 f = [f;f1+f(n,1) f2]; 088 end 089 090 return 091 function y=donothing(x) 092 y=x; 093 return 094 095 % Old call: (Saved it just in case...) 096 x=f(1,1)+(0:nmin-1)'/(nmin-1)*(f(nf,1)-f(1,1)); 097 % "y=interp1(f(:,1),f(:,2),x,'linear');" is slow: 098 % Use the fact that transforms often are "continuous" to 099 % find y_k=f(x_k) incrementally. 100 y=zeros(nmin,1); 101 i = 1; 102 for k=1:nmin 103 xx=x(k); 104 if (xx>f(i+1,1)) 105 while (xx>f(i+1,1)) 106 i=i+1; 107 if (i>=nf), i=nf-1; break; end; 108 end 109 else % xx<=f(i+1,1) 110 while (xx<=f(i,1)) 111 i=i-1; 112 if (i<1), i=1; break; end; 113 end 114 end 115 x0 = f(i,1); x1 = f(i+1,1); 116 y0 = f(i,2); y1 = f(i+1,2); 117 y(k) = (xx-x0)*(y1-y0)/(x1-x0)+y0; 118 end 119 f=[x y]; 120 clear x y; 121 122 123
Comments or corrections to the WAFO group