time   calls  line
                  1 function [res1,res2] = spdiags(arg1,arg2,arg3,arg4)
                  2 %SPDIAGS Sparse matrix formed from diagonals.
                  3 %   SPDIAGS, which generalizes the function "diag", deals with three
                  4 %   matrices, in various combinations, as both input and output.
                  5 %
                  6 %   [B,d] = SPDIAGS(A) extracts all nonzero diagonals from the m-by-n
                  7 %   matrix A.  B is a min(m,n)-by-p matrix whose columns are the p
                  8 %   nonzero diagonals of A.  d is a vector of length p whose integer
                  9 %   components specify the diagonals in A.
                 10 %
                 11 %   B = SPDIAGS(A,d) extracts the diagonals specified by d.
                 12 %   A = SPDIAGS(B,d,A) replaces the diagonals of A specified by d with
                 13 %       the columns of B.  The output is sparse.
                 14 %   A = SPDIAGS(B,d,m,n) creates an m-by-n sparse matrix from the
                 15 %       columns of B and places them along the diagonals specified by d.
                 16 %
                 17 %   Roughly, A, B and d are related by
                 18 %       for k = 1:p
                 19 %           B(:,k) = diag(A,d(k))
                 20 %       end
                 21 %
                 22 %   Example: These commands generate a sparse tridiagonal representation
                 23 %   of the classic second difference operator on n points.
                 24 %       e = ones(n,1);
                 25 %       A = spdiags([e -2*e e], -1:1, n, n)
                 26 %
                 27 %   Some elements of B, corresponding to positions "outside" of A, are
                 28 %   not actually used.  They are not referenced when B is an input and
                 29 %   are set to zero when B is an output.  If a column of B is longer than
                 30 %   the diagonal it's representing, elements of super-diagonals of A
                 31 %   correspond to the lower part of the column of B, while elements of
                 32 %   sub-diagonals of A correspond to the upper part of the column of B.
                 33 %
                 34 %   Example: This uses the top of the first column of B for the second
                 35 %   sub-diagonal and the bottom of the third column of B for the first
                 36 %   super-diagonal.
                 37 %       B = repmat((1:n)',1,3);
                 38 %       S = spdiags(B,[-2 0 1],n,n);
                 39 %
                 40 %   See also DIAG, SPEYE.
                 41 
                 42 %   Rob Schreiber
                 43 %   Copyright 1984-2007 The MathWorks, Inc.
                 44 %   $Revision: 5.12.4.9 $  $Date: 2010/09/02 13:37:02 $
                 45 
                 46 
  0.96  119397   47 if nargin <= 2 
                 48    if size(arg1,1) ~= size(arg1,2) && min(size(arg1))==1,
                 49       error(message('MATLAB:spdiags:VectorInput'));
                 50    end
                 51    % Extract diagonals
                 52    A = arg1;
                 53    if nargin == 1
                 54       % Find all nonzero diagonals
                 55       [i,j] = find(A);
                 56       d = sort(j-i);
                 57       d = d(find(diff([-inf; d])));
                 58    else
                 59       % Diagonals are specified
                 60       d = arg2(:);
                 61    end
                 62    [m,n] = size(A);
                 63    p = length(d);
                 64    B = zeros(min(m,n),p);
                 65    for k = 1:p
                 66       if m >= n
                 67          i = max(1,1+d(k)):min(n,m+d(k));
                 68       else
                 69          i = max(1,1-d(k)):min(m,n-d(k));
                 70       end
                 71       if ~isempty(i), B(i,k) = diag(A,d(k)); end
                 72    end
                 73    res1 = B;
                 74    res2 = d;
                 75 end
                 76 
  0.94  119397   77 if nargin >= 3 
                 78    % Create a sparse matrix from its diagonals
  0.22  119397   79    B = arg1; 
  0.19  119397   80    d = arg2(:); 
  0.14  119397   81    p = length(d); 
  0.05  119397   82    moda = (nargin == 3); 
  0.07  119397   83    if moda 
                 84       % Modify a given matrix
                 85       A = arg3;
  0.09  119397   86    else 
                 87       % Start from scratch
  1.36  119397   88       A = sparse(arg3,arg4); 
  0.21  119397   89    end 
                 90    
                 91    % Process A in compact form
  0.98  119397   92    [i,j,a] = find(A); 
  0.89  119397   93    a = [i j a]; 
  0.80  119397   94    [m,n] = size(A); 
                 95    
  0.40  119397   96    if moda 
                 97       for k = 1:p     
                 98          % Delete current d(k)-th diagonal
                 99          i = a(:,2)-a(:,1) == d(k);
                100          a(i,:) = [];
                101          % Append new d(k)-th diagonal to compact form
                102          i = (max(1,1-d(k)):min(m,n-d(k)))';
                103          a = [a; i i+d(k) B(i+(m>=n)*d(k),k)];
                104       end
  0.15  119397  105    else 
  0.59  119397  106       len=zeros(p+1,1); 
  0.15  119397  107       for k = 1:p 
  0.33  119411  108           len(k+1) = len(k)+length(max(1,1-d(k)):min(m,n-d(k))); 
  0.11  119411  109       end 
  0.40  119397  110       a = zeros(len(p+1),3); 
  0.09  119397  111       for k = 1:p 
                112          % Append new d(k)-th diagonal to compact form
  0.68  119411  113          i = (max(1,1-d(k)):min(m,n-d(k)))'; 
  1.79  119411  114          a((len(k)+1):len(k+1),:) = [i i+d(k) B(i+(m>=n)*d(k),k)]; 
  0.51  119411  115       end 
  0.13  119397  116    end 
  2.27  119397  117    res1 = sparse(a(:,1),a(:,2),full(a(:,3)),m,n); 
  0.15  119397  118    if (moda && islogical(A)) || islogical(B) 
                119        res1 = (res1~=0);
                120    end
  0.43  119397  121 end