www.gusucode.com > wavelet工具箱matlab源码程序 > wavelet/wavelet/conofinf.m
function varargout = conofinf(wname,scales,LSig,x,ArgPLOT) %CONOFINF Cone of influence for CWT. % For a signal of length LSIG, analyzed with the continuous wavelet % transform (CWT) using the wavelet WNAME, and for the scales SCALES, % CONE = CONOFINF(WAME,SCALES,LSIG,X) returns the cone of influence % (COI) of the Xth value of the signal. The Xth value of the signal % may be outside the interval [1,LSIG]. % CONE is a NbScales by LSIG matrix such that CONE(i,j) = 1 % when the jth component of the CWT at the ith scale is influenced % by the Xth value of the signal, and CONE(i,j) = 0 otherwise. % % If X is a vector, CONE is a cell array. Each cell contains the COI % of the corresponding component of X. % % In addition, [CONE,PL,PR,PLmin,PRmax] = CONOFINF(...,X) returns several % polynomials of degree 1. PL and PR are the equations of the left and % right edges of the COI respectively. % If X is a vector of length LX, PL and PR are LX by 2 matrices. % (1-PL(:,2))./PL(:,1) and 1-PR(:,2))./PR(:,1) give respectively % the abscissa of left and right edges for the scale 1. % % PLmin and PRmax are polynomials of degree 1 which give the equation % of the minimal left vertex, and the equation of the maximal right % vertex respectively, completely included in the (NbScales by LSIG) % domain of the CWT. % % Note that [PLmin,PRmax] = CONOFINF(WAME,SCALES,LSIG) or % [PLmin,PRmax] = CONOFINF(WAME,SCALES,LSIG,[]) % returns only the PLmin and PRmax polynomials. % % In Addition, you can plot the cones of influence using: % [...] = CONOFINF(...,'plot') % M. Misiti, Y. Misiti, G. Oppenheim, J.M. Poggi 08-Feb-2010. % Last Revision: 04-Jan-2011. % Copyright 1995-2011 The MathWorks, Inc. % Check inputs. flagPLOT = false; if length(scales)<2 , scales = 1:scales; end if isempty(LSig) || ~isnumeric(LSig) || numel(LSig)>1 || ... (LSig<1) || LSig~=fix(LSig) error(message('Wavelet:FunctionInput:TypeOfbArg')) end if nargin<4 , x = []; end if nargin>4 , flagPLOT = ~isequal(ArgPLOT,false); end % Check outputs. if ~isempty(x) , maxOUT = 7; else maxOUT = 4; end error(nargoutchk(0,maxOUT,nargout, 'struct')) % Get the support of the wavelet. bounds = wavsupport(wname); % Compute the limits of cones of influence. xR = LSig-bounds(2)*scales(end); xL = -bounds(1)*scales(end)+1; Lmin = [0 scales(end);xL 0]; Rmax = [xR 0;LSig scales(end)]; % Compute the cone of influence for each x value. if ~isempty(x) x = x(:); scales = scales(:)'; nbS = length(scales); nbX = length(x); L = repmat(bounds(1)*scales,nbX,1) + repmat(x,1,nbS); R = repmat(bounds(2)*scales,nbX,1) + repmat(x,1,nbS); PL = zeros(nbX,2); PR = zeros(nbX,2); for k = 1:nbX PL(k,:) = getLine([L(k,[1 end])' , scales([1 end])']); PR(k,:) = getLine([R(k,[1 end])' , scales([1 end])']); end cone = cell(1,nbX); warning('Off','MATLAB:colon:nonIntegerIndex') for k = 1:nbX coneofx = zeros(nbS,LSig); L(k,:) = max(L(k,:),1); R(k,:) = min(R(k,:),LSig); for j=1:nbS coneofx(j,L(k,j):R(k,j)) = 1; end cone{k} = coneofx; end warning('On','MATLAB:colon:nonIntegerIndex') if nbX==1 , cone = cone{1}; end else L = []; R = []; cone = []; end Pmin = getLine(Lmin); Pmax = getLine(Rmax); xMin = (1-Pmin(2))/Pmin(1); xMax = (1-Pmax(2))/Pmax(1); if ~isempty(x) varargout = {cone,PL,PR,Pmin,Pmax,xMin,xMax}; else varargout = {Pmin,Pmax,xMin,xMax}; end % Plot the cones if requested. if ~flagPLOT , return; end LW = 2; hold on; plot(Lmin(:,1),Lmin(:,2),'b','LineWidth',LW); stem(xMin,scales(1),'filled',... 'MarkerSize',7,'MarkerEdgeColor','b','MarkerFaceColor','b'); plot(Rmax(:,1),Rmax(:,2),'r','LineWidth',LW); stem(xMax,scales(1),'filled',... 'MarkerSize',7,'MarkerEdgeColor','r','MarkerFaceColor','r'); if ~isempty(x) LW = 1; midPts = 0.5*(L(:,1)+R(:,1)); for k = 1:nbX plot(L(k,:),scales,'m:','LineWidth',LW); plot(R(k,:),scales,'m:','LineWidth',LW); stem(midPts(k,1),scales(1),'filled',... 'MarkerSize',7,'MarkerEdgeColor','k','MarkerFaceColor','k'); iL1 = find(L(k,:)==1,1,'first'); iL2 = find(L(k,:)==1,1,'last'); if isempty(iL1) , iL1 = nbS ; iL2 = nbS; end iR1 = find(R(k,:)==LSig,1,'first'); iR2 = find(R(k,:)==LSig,1,'last'); if isempty(iR1) , iR1 = nbS; iR2 = nbS; end t = [L(k,iL1) L(k,1) R(k,1) R(k,iR1) ... R(k,iR2) L(k,iL2) L(k,iL1)]; s = [scales(iL1) 1 1 scales(iR1) ... scales(iR2) scales(iL2) scales(iL1)]; hcoi = fill(t,s,[1 1 1]); set(hcoi,'alphadatamapping','direct','facealpha',.5) end end if isempty(x) t = [Lmin(1,1) , Lmin(2,1) , Rmax(1,1) , Rmax(2,1) , Lmin(1,1)]; s = [Lmin(1,2) , Lmin(2,2) , Rmax(1,2) , Rmax(2,2) , Lmin(1,2)]; hcoi = fill(t,s,[1 0.8 1]); set(hcoi,'alphadatamapping','direct','facealpha',.4) end if LSig<=20 xlim = [min([xMax,1]) max([xMin,LSig])]; plotXL = true; else xlim = [1,LSig]; plotXL = false; end if plotXL stem([1 LSig],[scales(1),scales(1)],'filled','MarkerSize',7,... 'MarkerEdgeColor',[0 0.8 0],'MarkerFaceColor',[0 0.8 0]); end set(gca,'XLim',xlim,'YLim',[1,scales(end)]); function P = getLine(L) dx = L(2,1)-L(1,1); dy = L(2,2)-L(1,2); P(1) = dy/dx; dxy = L(1,2)*L(2,1)-L(1,1)*L(2,2); P(2) = dxy/dx;