sgwt_cheby_op : Chebyshev polynomial of Laplacian applied to vector function r=sgwt_cheby_op(f,L,c,arange) Compute (possibly multiple) polynomials of laplacian (in Chebyshev basis) applied to input. Coefficients for multiple polynomials may be passed as a cell array. This is equivalent to setting r{1}=sgwt_cheby_op(f,L,c{1},arange); r{2}=sgwt_cheby_op(f,L,c{2},arange); ... but is more efficient as the Chebyshev polynomials of L applied to f can be computed once and shared. Inputs: f- input vector L - graph laplacian (should be sparse) c - Chebyshev coefficients. If c is a plain array, then they are coefficients for a single polynomial. If c is a cell array, then it contains coefficients for multiple polynomials, such that c{j}(1+k) is k'th Chebyshev coefficient the j'th polynomial. arange - interval of approximation Outputs: r - result. If c is cell array, r will be cell array of vectors size of f. If c is a plain array, r will be a vector the size of f.
0001 % sgwt_cheby_op : Chebyshev polynomial of Laplacian applied to vector 0002 % 0003 % function r=sgwt_cheby_op(f,L,c,arange) 0004 % 0005 % Compute (possibly multiple) polynomials of laplacian (in Chebyshev 0006 % basis) applied to input. 0007 % 0008 % Coefficients for multiple polynomials may be passed as a cell array. This is 0009 % equivalent to setting 0010 % r{1}=sgwt_cheby_op(f,L,c{1},arange); 0011 % r{2}=sgwt_cheby_op(f,L,c{2},arange); 0012 % ... 0013 % 0014 % but is more efficient as the Chebyshev polynomials of L applied 0015 % to f can be computed once and shared. 0016 % 0017 % Inputs: 0018 % f- input vector 0019 % L - graph laplacian (should be sparse) 0020 % c - Chebyshev coefficients. If c is a plain array, then they are 0021 % coefficients for a single polynomial. If c is a cell array, 0022 % then it contains coefficients for multiple polynomials, such 0023 % that c{j}(1+k) is k'th Chebyshev coefficient the j'th polynomial. 0024 % arange - interval of approximation 0025 % 0026 % Outputs: 0027 % r - result. If c is cell array, r will be cell array of vectors 0028 % size of f. If c is a plain array, r will be a vector the size 0029 % of f. 0030 0031 % This file is part of the SGWT toolbox (Spectral Graph Wavelet Transform toolbox) 0032 % Copyright (C) 2010, David K. Hammond. 0033 % 0034 % The SGWT toolbox is free software: you can redistribute it and/or modify 0035 % it under the terms of the GNU General Public License as published by 0036 % the Free Software Foundation, either version 3 of the License, or 0037 % (at your option) any later version. 0038 % 0039 % The SGWT toolbox is distributed in the hope that it will be useful, 0040 % but WITHOUT ANY WARRANTY; without even the implied warranty of 0041 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 0042 % GNU General Public License for more details. 0043 % 0044 % You should have received a copy of the GNU General Public License 0045 % along with the SGWT toolbox. If not, see <http://www.gnu.org/licenses/>. 0046 0047 function r=sgwt_cheby_op(f,L,c,arange) 0048 0049 if ~iscell(c) 0050 r=sgwt_cheby_op(f,L,{c},arange); 0051 r=r{1}; 0052 return; 0053 end 0054 0055 Nscales=numel(c); 0056 M=zeros(size(Nscales)); 0057 for j=1:Nscales 0058 M(j)=numel(c{j}); 0059 end 0060 assert(all(M>=2)); 0061 0062 maxM=max(M); 0063 %Twf_new = T_j(L) f 0064 %Twf_cur T_{j-1}(L) f 0065 %TWf_old T_{j-2}(L) f 0066 0067 a1=(arange(2)-arange(1))/2; 0068 a2=(arange(2)+arange(1))/2; 0069 0070 Twf_old=f; %j=0; 0071 Twf_cur=(L*f-a2*f)/a1; % j=1; 0072 for j=1:Nscales 0073 r{j}=.5*c{j}(1)*Twf_old + c{j}(2)*Twf_cur; 0074 end 0075 0076 for k=2:maxM 0077 Twf_new = (2/a1)*(L*Twf_cur-a2*Twf_cur)-Twf_old; 0078 for j=1:Nscales 0079 if 1+k<=M(j) 0080 r{j}=r{j}+c{j}(k+1)*Twf_new; 0081 end 0082 end 0083 Twf_old=Twf_cur; 0084 Twf_cur=Twf_new; 0085 end