qt_aggregate

 PUPROSE:  Use quadtree decomposition on a matrix D and assign the
           aggregated value of the original matrix to a sparse matrix
           M, based on the sparse matrix S[0..1] generated using QTDECOMP
 -------------------------------------------------------------------
 USAGE: [M,S] = qt_aggregate(D,varargin)
 where: D is the input matrix    
        varargin are parameters options for the quadtreedecompostion
        M is the output sparse matrix with aggregated values of D
        S is the sparse (intensity) matrix generated using QTDECOMP

 OPTIONS: 
        'Sin'  is an input sparse matrix holding an existing quadtree
               structure created using QTDECOMP. If Sin is given, no
               additional quadtree decomposition is executed, but values
               of D are aggregated according to the structure given by Sin

        'crit' is the criterion for qtdecomp to split blocks of range [0..1],
               where the normalised maximum of a block minus the normalised
               minimum is greater than 'crit', a block is split further
        'minq' is the minimum block size created (default=1)
        'maxq' is the maximum block size created by qtdecomp (default=512)

 EXAMPLE:
       [M,S] = qt_aggregate(peaks(512),'crit',0.05,'minq',1,'maxq',16);
       - uses quadtree decomposition on a 512x512 example peaks surface,
         creating blocks betweeen 1 and 16 cell sizes.
       - approx. 3400 blocks are created, compared to 262144 cells of the
         original full input matrix
       - values from the input full matrix are aggregated (mean) for each
         block and written to sparse output matrix M
       - Sparse matrix S holds upper left corner cordinates and block size

 See also: QTDECOMP

 Note: Input matrix needs to be square, and size(D) needs to be multiple of
       min/max block size (minq/maxq). This is a requirement for qtdecomp.

 Felix Hebeler, Dept. of Geography, University Zurich, June 2008.

0001 function [M,S] = qt_aggregate(D,varargin)
0002 % PUPROSE:  Use quadtree decomposition on a matrix D and assign the
0003 %           aggregated value of the original matrix to a sparse matrix
0004 %           M, based on the sparse matrix S[0..1] generated using QTDECOMP
0005 % -------------------------------------------------------------------
0006 % USAGE: [M,S] = qt_aggregate(D,varargin)
0007 % where: D is the input matrix
0008 %        varargin are parameters options for the quadtreedecompostion
0009 %        M is the output sparse matrix with aggregated values of D
0010 %        S is the sparse (intensity) matrix generated using QTDECOMP
0011 %
0012 % OPTIONS:
0013 %        'Sin'  is an input sparse matrix holding an existing quadtree
0014 %               structure created using QTDECOMP. If Sin is given, no
0015 %               additional quadtree decomposition is executed, but values
0016 %               of D are aggregated according to the structure given by Sin
0017 %
0018 %        'crit' is the criterion for qtdecomp to split blocks of range [0..1],
0019 %               where the normalised maximum of a block minus the normalised
0020 %               minimum is greater than 'crit', a block is split further
0021 %        'minq' is the minimum block size created (default=1)
0022 %        'maxq' is the maximum block size created by qtdecomp (default=512)
0023 %
0024 % EXAMPLE:
0025 %       [M,S] = qt_aggregate(peaks(512),'crit',0.05,'minq',1,'maxq',16);
0026 %       - uses quadtree decomposition on a 512x512 example peaks surface,
0027 %         creating blocks betweeen 1 and 16 cell sizes.
0028 %       - approx. 3400 blocks are created, compared to 262144 cells of the
0029 %         original full input matrix
0030 %       - values from the input full matrix are aggregated (mean) for each
0031 %         block and written to sparse output matrix M
0032 %       - Sparse matrix S holds upper left corner cordinates and block size
0033 %
0034 % See also: QTDECOMP
0035 %
0036 % Note: Input matrix needs to be square, and size(D) needs to be multiple of
0037 %       min/max block size (minq/maxq). This is a requirement for qtdecomp.
0038 %
0039 % Felix Hebeler, Dept. of Geography, University Zurich, June 2008.
0040 
0041 %% configure inputs
0042 %default parameters
0043 crit=0.2;
0044 minq=1;
0045 maxq=512;
0046 
0047 %parse inputs
0048 if isempty(varargin)~=1     % check if any arguments are given
0049     [m1,n1]=size(varargin);
0050     opts={'crit';'minq';'maxq';'Sin'};
0051     for i=1:n1;             % check which parameters are given
0052         indi=strcmpi(varargin{i},opts);
0053         ind=find(indi==1);
0054         if isempty(ind)~=1
0055             switch ind
0056                 case 1
0057                     crit=varargin{i+1};
0058                 case 2
0059                     minq=varargin{i+1};
0060                 case 3
0061                     maxq=varargin{i+1};
0062                 case 4
0063                     S=varargin{i+1};
0064             end
0065         end
0066     end
0067 end
0068 %% do quadtree decomposition if Sin is not given
0069 if ~exist('S','var')
0070     % create intensity image
0071     I=(D-min(D(:)))./(max(D(:))-min(D(:)));
0072     % perfom standard quadtree decomposition on it
0073     S=qtdecomp(I,crit,[minq maxq]);
0074     clear I;
0075 end
0076 %% could use minq and maxq but it is possible not all block sizes are being used
0077 z_min=min(nonzeros(S));
0078 z_max=max(nonzeros(S));
0079 
0080 %% create matrix to hold result
0081 M=spalloc(size(S,1),size(S,2),nnz(S));
0082 
0083 % loop over all possible block sizes
0084 for dim=log2(z_min):log2(z_max)
0085     [vals,r,c] = qtgetblk(D,S,2^dim);
0086     for v=1:size(vals,3)
0087                M(r(v),c(v)) = aggregate(vals(:,:,v));
0088     end   
0089 end
0090 
0091 
0092 %% separated so different aggregation functions (eg weighted) can be used
0093 %% easily
0094 function a=aggregate(x)
0095 a=mean(mean(x));
0096