Overloaded double method where efficient.

mpf · Jul 30, 2009 · 5ae5e18 · 5ae5e18
1 parent 63a686b
commit 5ae5e18
Show file tree

Hide file tree

Showing 15 changed files with 288 additions and 203 deletions.
diff --git a/TODO b/TODO
@@ -1,19 +1,10 @@
-TODO ($Id$)
+TODO
 
 [] opSpot class: move funhandle and matrix props, which are def'd in
    many (most?!) of the subclasses, into opSpot. (These are used too
    often.) Then can have opSpot/double.m do a simple check if
    ~isempty(op.matrix).
 
-[] isposintmat lives in both private and @opSpot/private -- is this
-really neces sary!
-
-[] Read the Rice Vector Library manual.
-
-[] Method for generating random elements (ie, vectors) from domain and
-   range of operator. Note: special cases such as real FFT. Useful for
-   dot-tests.
-
 [] opBlockDiag.m: Avoid replicating operator explicitly?  How to
                   display repeated blocks:
                   opBlockDiag(ones(10,1),opDCT(4))?  make sure help
@@ -30,7 +21,6 @@ really neces sary!
 [] opWavelet.m: Help comment
 [] opToeplitz.m, opToepGauss.m, opToepSign.m: Help comment
 [] Comments in the @opSpot files
-[] Set ID property in SVN
 [] opSurfacelet.m: Haven't got a clue about how to use? Still have to test.
 [] For each operator, check if complexity, linear flag needs to be set or not
 [] opSubsRef.m, opSubsAsgn.m, opExcise.m: Add proper help
@@ -55,96 +45,3 @@ randn('state'). These are not compatible with the latest stream-based
 calls. Need to document this.
 
 
-DONE
-- opBernoulli
-- opBinary
-- opBlockDiag
-- opClass
-- opConj
-- opConvolve
-- opCurvelet
-- opDiag
-- opDictionary
-- opDirac
-- opEmpty
-- opExcise
-- opEye
-- opFFT
-- opFoG
-- opFunction
-- opGaussian
-- opHaar
-- opHadamard
-- opHeaviside
-- opImag
-- opInverse
-- opKron
-- opMatrix
-- opMask
-- opMinus
-- opOnes
-- opPower
-- opReal
-- opRestriction
-- opSparseBinary
-- opStack
-- opSubsAsgn
-- opSubsRef
-- opSum
-- opSurfacelet
-- opToepGauss
-- opToeplitz
-- opToepSign
-- opTranspose
-- opUnaryMinus
-- opWavelet
-- opWindow
-- opZeros
-
-DONE
-- bicg
-- bicgstab
-- blkdiag
-- cgs
-- char
-- conj
-- ctranspose
-- display
-- disp
-- double
-- end
-- full
-- gmres
-- horzcat
-- imag
-- inv
-- isscalar
-- isreal
-- isempty
-- kron
-- lsqr
-- minres
-- minus
-- mldivide
-- mtimes
-- mpower
-- mrdivide
-- ndims
-- pcg
-- plus
-- qmr
-- real
-- size
-- subsasgn
-- subsref
-- symmlq
-- transpose
-- uminus
-- uplus
-- vertcat
-
-[] opFFT: should this be unnormalized so that it mimics Matlab's
-          built-in fft? In other words, should we expect fft(y) ==
-          opFFT*y?
-
-   NO! Normalizing makes it an orthogonal matrix.
diff --git a/opDiag.m b/opDiag.m
@@ -4,41 +4,48 @@
 %   diagonal matrix with D on its diagonal.
 
 %   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opDiag < opSpot
 
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     % Properties
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     properties (SetAccess = private)
-       diag = []; % Diagonal entries
+       diag     % Diagonal entries
     end % Properties
 
-
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     % Methods - Public
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     methods
 
-        % Constructor
-        function op = opDiag(d)
-           if nargin ~= 1
-              error('Invalid number of arguments.');
-           end
-
-           % Vectorize d and get size
-           d = d(:);
-           n = length(d);
-
-           % Construct operator
-           op = op@opSpot('Diag',n,n);
-           op.cflag      = ~isreal(d);
-           op.linear     = 1;
-           op.diag       = d;
-           op.precedence = 1;
-        end % Constructor
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       % Constructor
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       function op = opDiag(d)
+          if nargin ~= 1
+             error('Invalid number of arguments.');
+          end
+
+          % Vectorize d and get size
+          d = d(:);
+          n = length(d);
+
+          % Construct operator
+          op = op@opSpot('Diag',n,n);
+          op.cflag      = ~isreal(d);
+          op.linear     = 1;
+          op.diag       = d;
+          op.precedence = 1;
+       end % Constructor
+
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       % double
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       function A = double(op)
+          A = diag(op.diag);
+       end % double
 
     end % Methods
 

diff --git a/opDictionary.m b/opDictionary.m
@@ -1,20 +1,17 @@
 %opDictionary   Dictionary of concatenated operators
 %
-%   opDictionary(WEIGHTS,OP1,OP2,...OPn) creates a dictionary
-%   operator consisting of the concatenation of all operators;
-%   [WEIGHT1*OP1 | WEIGHT2*OP2 | ... WEIGHTn*OPn]. Vector WEIGHTS
-%   can be of size nx1 or 1xn. If the same weight is to be applied
-%   to each operator, set WEIGHTS to a scalar. When WEIGHTS is
-%   empty [], it is set to  one. The WEIGHT parameter can be
-%   omitted as long as OP1 is not a vector of length (n-1); in
-%   which case there is no way to decide whether it is a weight
-%   vector or operator.
-% 
-%   See also opFoG, opStack, opSum.
+%   D = opDictionary(OP1,OP2,...OPn) creates a dictionary
+%   operator consisting of the concatenation of all operators, i.e.,
+%   
+%       D = [ OP1, OP2, ..., OPn ].
+%
+%   In general, it's best to use Matlab's horizonal concatenation
+%   operations instead of calling opDictionary. (The two are equivalent.)
+%
+%   See also opFoG, opStack, opSum, @opSpot/horzcat.
 
-%   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   Copyright 2008-2009, Ewout van den Berg and Michael P. Friedlander
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opDictionary < opSpot
 
@@ -27,33 +24,16 @@
        % Constructor
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        function op = opDictionary(varargin)
-          % Checks weights parameter
-          if ~isnumeric(varargin{1})
-             weights = ones(nargin,1);
-             opList = varargin;
-          else
-             weights = varargin{1};
-             if isempty(weights), weights = 1; end;
-             [m,n] = size(weights);
-             if (((m == 1) && (n == nargin-1)) || ...
-                 ((n == 1) && (m == nargin-1)) || ...
-                 ((m == 1) && (n == 1)))
-               weights = ones(nargin-1,1).*weights(:); 
-               opList = varargin(2:end);
-             else
-               weights = ones(nargin,1);
-               opList = varargin;
-             end
-          end
 
           % Check number of operators
-          if (length(opList) < 1)
-            error('At least one operator must be specified.');
+          opList = varargin;
+          if isempty(opList)
+             error('At least one operator must be specified.');
           end
 
           % Convert all arguments to operators
           for i=1:length(opList)
-             if ~isa(opList{i},'opSpot') 
+             if ~isa(opList{i},'opSpot')
                 opList{i} = opMatrix(opList{i});
              end
           end
@@ -110,6 +90,20 @@
 
           str = [str, ']'];
        end % Display
+
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       % Double
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       function A = double(op)
+          A = zeros(size(op));
+          k = 0;
+          for i=1:length(op.children)
+             child = op.children{i};
+             n = size(child,2);
+             A(:,k+1:k+n) = double(child);
+             k = k + n;
+          end
+       end % double
 
     end % Methods
 

diff --git a/opDirac.m b/opDirac.m
@@ -1,12 +1,11 @@
 %opDirac  Dirac basis
 %
-%   opDirac(N) creates the square N by N identity operator. Without
+%   opDirac(N) creates the square N-by-N identity operator. Without
 %   any arguments an operator corresponding to the scalar 1 is
 %   created.
 
 %   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opDirac < opSpot
 
@@ -22,6 +21,10 @@
             op = op@opSpot('Dirac',n,n);
         end
 
+        function A = double(op)
+           A = eye(size(op));
+        end
+
     end % Methods
 
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

diff --git a/opEmpty.m b/opEmpty.m
@@ -4,8 +4,7 @@
 %   matrix. At least one of M and N must be zero.
 
 %   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opEmpty < opSpot
 
@@ -28,7 +27,14 @@
 
           op = op@opSpot('Empty',m,n);
        end % Constructor
-
+
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       % Double
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       function A = double(op)
+          A = [];          
+       end       
+
     end % Methods
 
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

diff --git a/opEye.m b/opEye.m
@@ -1,12 +1,11 @@
 %opEye  Identity operator
 %
-%   opEye(M,N) creates the identity operator of size M by N. If N is
+%   opEye(M,N) creates the identity operator of size M-by-N. If N is
 %   omitted it is set to M by default. Without any arguments an
 %   operator corresponding to the scalar 1 is created.
 
 %   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opEye < opSpot
 
@@ -21,6 +20,10 @@
            if nargin < 2, n = m; end
            op = op@opSpot('Eye',m,n);
         end % Constructor
+
+        function A = double(op)
+           A = eye(size(op));
+        end % double
 
     end % Methods
 

diff --git a/opFoG.m b/opFoG.m
@@ -10,8 +10,7 @@
 %   See also opDictionary, opStack, opSum.
 
 %   Copyright 2009, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id$
+%   http://www.cs.ubc.ca/labs/scl/spot
 
 classdef opFoG < opSpot
 
@@ -71,10 +70,19 @@
 
           % Preprocess operators
           op.operators = {A,B};
-          if isscalar(A), op.operators{1} = opMatrix(double(A)); end;
-          if isscalar(B), op.operators{2} = opMatrix(double(B)); end;
+          if isscalar(A), op.operators{1} = opMatrix(double(A)); end
+          if isscalar(B), op.operators{2} = opMatrix(double(B)); end
        end % Constructor
-
+
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       % double
+       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+       function A = double(op)
+          C1 = op.children{1};
+          C2 = op.children{2};
+          A  = double(C1)*double(C2);
+       end
+
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        % Display
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%