# HG changeset patch
# User idamnjanovic
# Date 1300284897 0
# Node ID 5f1f436057ca375ff53ecc7ebe79d6e5b1ad80f4
# Parent  f2ca04532eac0551c39be4d58efd0ad940adab81


diff -r f2ca04532eac -r 5f1f436057ca Problems/private/add_dc.m
--- a/Problems/private/add_dc.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,33 +0,0 @@
-function x = add_dc(y,dc,columns)
-%ADD_DC Add DC channel to signals.
-%   X = ADD_DC(Y,DC) adds the specified DC value to the (possibly
-%   multi-dimensional) signal Y, returning the result as X. DC should be a
-%   scalar value.
-%
-%   X = ADD_DC(Y,DC,'columns') where Y is a 2D matrix and DC is an array of
-%   length size(Y,2), treats the columns of Y as individual 1D signals, 
-%   adding to each one the corresponding DC value from the DC array. X is
-%   the same size as Y and contains the resulting signals.
-%
-%   See also REMOVE_DC.
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2009
-
-
-if (nargin==3 && strcmpi(columns,'columns')), columns = 1;
-else columns = 0;
-end
-
-if (columns)
-  x = addtocols(y,dc);
-else
-  x = y + dc;
-end
-
-
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/addtocols.c
--- a/Problems/private/addtocols.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,85 +0,0 @@
-/**************************************************************************
- *
- * File name: addtocols.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 19.4.2009
- *
- *************************************************************************/
-
-
-#include "mex.h"
-
-
-/* Input Arguments */
-
-#define	X_IN	prhs[0]
-#define V_IN  prhs[1]
-
-
-/* Output Arguments */
-
-#define	Y_OUT	plhs[0]
-
-
-void mexFunction(int nlhs, mxArray *plhs[], 
-		             int nrhs, const mxArray*prhs[])
-     
-{ 
-    double *x, *y, *v, *xend;  
-    mwSize m,n,m1,n1;
-    mwIndex counter; 
-    
-    
-    /* Check for proper number of arguments */
-    
-    if (nrhs != 2) {
-      mexErrMsgTxt("Two input arguments required."); 
-    } else if (nlhs > 1) {
-      mexErrMsgTxt("Too many output arguments."); 
-    } 
-    
-    
-    /* Check the the input dimensions */ 
-    
-    m = mxGetM(X_IN);
-    n = mxGetN(X_IN);
-    if (!mxIsDouble(X_IN) || mxIsComplex(X_IN) || mxGetNumberOfDimensions(X_IN)>2) {
-      mexErrMsgTxt("ADDTOCOLS requires that X be a double matrix.");
-    }
-    m1 = mxGetM(V_IN);
-    n1 = mxGetN(V_IN);
-    if (!mxIsDouble(V_IN) || mxIsComplex(V_IN) || (m1!=1 && n1!=1)) {
-      mexErrMsgTxt("ADDTOCOLS requires that V be a double vector.");
-    }
-    if ((m1==1 && n1!=n) || (n1==1 && m1!=n)) {
-      mexErrMsgTxt("Error in ADDTOCOLS: dimensions of V and X must agree.");
-    }
-    
-    
-    /* Create a matrix for the return argument */ 
-    Y_OUT = mxCreateDoubleMatrix(m, n, mxREAL);
-         
-    
-    /* Assign pointers to the various parameters */ 
-    x = mxGetPr(X_IN); 
-    v = mxGetPr(V_IN);
-    y = mxGetPr(Y_OUT);
-            
-    
-    /* Do the actual computation */
-    
-    xend = x+(m*n);
-    counter = 0;
-    while (x<xend) {
-      (*y) = (*x) + (*v);
-      y++; x++; counter++;
-      if (counter==m) {v++; counter=0;}
-    }
-    
-    return;
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/addtocols.m
--- a/Problems/private/addtocols.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,13 +0,0 @@
-%ADDTOCOLS Add values to the columns of a matrix.
-%  Y=ADDTOCOLS(X,V) adds to each column of the MxN matrix X the
-%  corresponding value from the N-element vector V.
-%
-%  See also NORMCOLS.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  June 2005
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/addtocols.mexa64
Binary file Problems/private/addtocols.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/col2imstep.c
--- a/Problems/private/col2imstep.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,146 +0,0 @@
-/**************************************************************************
- *
- * File name: col2imstep.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 31.8.2009
- *
- *************************************************************************/
-
-
-#include "mex.h"
-
-
-/* Input Arguments */
-
-#define	B_IN	 prhs[0]
-#define N_IN   prhs[1]
-#define SZ_IN  prhs[2]
-#define S_IN   prhs[3]
-
-
-/* Output Arguments */
-
-#define	X_OUT	plhs[0]
-
-
-void mexFunction(int nlhs, mxArray *plhs[], 
-		             int nrhs, const mxArray*prhs[])
-     
-{ 
-    double *x, *b, *s;
-    mwSize sz[3], stepsize[3], n[3], ndims;
-    mwIndex i, j, k, l, m, t, blocknum;
-    
-    
-    /* Check for proper number of arguments */
-    
-    if (nrhs < 3 || nrhs > 4) {
-      mexErrMsgTxt("Invalid number of input arguments."); 
-    } else if (nlhs > 1) {
-      mexErrMsgTxt("Too many output arguments."); 
-    } 
-    
-    
-    /* Check the the input dimensions */ 
-    
-    if (!mxIsDouble(B_IN) || mxIsComplex(B_IN) || mxGetNumberOfDimensions(B_IN)>2) {
-      mexErrMsgTxt("B should be a double matrix.");
-    }
-    if (!mxIsDouble(N_IN) || mxIsComplex(N_IN) || mxGetNumberOfDimensions(N_IN)>2) {
-      mexErrMsgTxt("Invalid output matrix size.");
-    }
-    ndims = mxGetM(N_IN)*mxGetN(N_IN);
-    if (ndims<2 || ndims>3) {
-      mexErrMsgTxt("Output matrix can only be 2-D or 3-D.");
-    }
-    if (!mxIsDouble(SZ_IN) || mxIsComplex(SZ_IN) || mxGetNumberOfDimensions(SZ_IN)>2 || mxGetM(SZ_IN)*mxGetN(SZ_IN)!=ndims) {
-      mexErrMsgTxt("Invalid block size.");
-    }
-    if (nrhs == 4) {
-      if (!mxIsDouble(S_IN) || mxIsComplex(S_IN) || mxGetNumberOfDimensions(S_IN)>2 || mxGetM(S_IN)*mxGetN(S_IN)!=ndims) {
-        mexErrMsgTxt("Invalid step size.");
-      }
-    }
-    
-    
-    /* Get parameters */
-    
-    s = mxGetPr(N_IN);
-    if (s[0]<1 || s[1]<1 || (ndims==3 && s[2]<1)) {
-      mexErrMsgTxt("Invalid output matrix size.");
-    }
-    n[0] = (mwSize)(s[0] + 0.01);
-    n[1] = (mwSize)(s[1] + 0.01);
-    n[2] = ndims==3 ? (mwSize)(s[2] + 0.01) : 1;
-    
-    s = mxGetPr(SZ_IN);
-    if (s[0]<1 || s[1]<1 || (ndims==3 && s[2]<1)) {
-      mexErrMsgTxt("Invalid block size.");
-    }
-    sz[0] = (mwSize)(s[0] + 0.01);
-    sz[1] = (mwSize)(s[1] + 0.01);
-    sz[2] = ndims==3 ? (mwSize)(s[2] + 0.01) : 1;
-    
-    if (nrhs == 4) {
-      s = mxGetPr(S_IN);
-      if (s[0]<1 || s[1]<1 || (ndims==3 && s[2]<1)) {
-        mexErrMsgTxt("Invalid step size.");
-      }
-      stepsize[0] = (mwSize)(s[0] + 0.01);
-      stepsize[1] = (mwSize)(s[1] + 0.01);
-      stepsize[2] = ndims==3 ? (mwSize)(s[2] + 0.01) : 1;
-    }
-    else {
-      stepsize[0] = stepsize[1] = stepsize[2] = 1;
-    }
-    
-    if (n[0]<sz[0] || n[1]<sz[1] || (ndims==3 && n[2]<sz[2])) {
-      mexErrMsgTxt("Block size too large.");
-    }
-
-    if (mxGetN(B_IN) != ((n[0]-sz[0])/stepsize[0]+1)*((n[1]-sz[1])/stepsize[1]+1)*((n[2]-sz[2])/stepsize[2]+1)) {
-      mexErrMsgTxt("Invalid number of columns in B. Please use IM2COLSTEP to compute B.");
-    }
-    
-    
-    /* Create a matrix for the return argument */
-    
-    X_OUT = mxCreateNumericArray(ndims, n, mxDOUBLE_CLASS, mxREAL);
-    
-    
-    /* Assign pointers */
-    
-    b = mxGetPr(B_IN);
-    x = mxGetPr(X_OUT);
-            
-    
-    /* Do the actual computation */
-    
-    blocknum = 0;
-    
-    /* iterate over all blocks */
-    for (k=0; k<=n[2]-sz[2]; k+=stepsize[2]) {
-      for (j=0; j<=n[1]-sz[1]; j+=stepsize[1]) {
-        for (i=0; i<=n[0]-sz[0]; i+=stepsize[0]) {
-          
-          /* add single block */
-          for (m=0; m<sz[2]; m++) {
-            for (l=0; l<sz[1]; l++) {
-              for (t=0; t<sz[0]; t++) {
-                (x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i)[t] += (b + blocknum*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0])[t];
-              }
-            }
-          }
-          blocknum++;
-          
-        }
-      }
-    }
-    
-    return;
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/col2imstep.m
--- a/Problems/private/col2imstep.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,25 +0,0 @@
-%COL2IMSTEP Rearrange matrix columns into blocks.
-%  A = COL2IMSTEP(B,[MM NN],[N1 N2]) rearranges the columns of B into
-%  sliding N1-by-N2 blocks producing the matrix A of size MM-by-NN. B is
-%  usually the result of calling IM2COLSTEP(...) with a stepsize of 1, or
-%  using Matlab's IM2COL(..,'sliding'). Overlapping blocks are summed in A.
-%
-%  A = COL2IMSTEP(B,[MM NN],[N1 N2],[S1 S2]) arranges the blocks in A with
-%  a step size of (S1,S2) between them. The first block is at A(1:N1,1:N2),
-%  and the rest are at A((1:N1)+i*S1,(1:N2)+j*S2). Overlapping blocks are
-%  summed in A. Note that B is usually the result of calling
-%  IM2COLSTEP(...) with a stepsize of [S1 S2].
-%
-%  A = IM2COLSTEP(B,[MM NN KK],[N1 N2 N3],[S1 S2 S3]) generates a 3-D
-%  output matrix A. The step size [S1 S2 S3] may be ommitted, and defaults
-%  to [1 1 1].
-%
-%  See also IM2COLSTEP, IM2COL, COUNTCOVER.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  August 2009
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/col2imstep.mexa64
Binary file Problems/private/col2imstep.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/collincomb.c
--- a/Problems/private/collincomb.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,165 +0,0 @@
-/**************************************************************************
- *
- * File name: collincomb.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 21.5.2009
- *
- *************************************************************************/
-
-
-#include "mex.h"
-
-
-/* Input Arguments */
-
-#define	A_IN	   prhs[0]
-#define ROWS_IN  prhs[1]
-#define COLS_IN1 prhs[1]
-#define COLS_IN2 prhs[2]
-#define X_IN1    prhs[2]
-#define X_IN2    prhs[3]
-
-
-/* Output Arguments */
-
-#define	Y_OUT	plhs[0]
-
-
-void mexFunction(int nlhs, mxArray *plhs[],
-int nrhs, const mxArray*prhs[])
-
-{
-  double *A, *x, *y, *rows, *cols;
-  mwSize m,n,m1,n1,m2,n2,rownum,colnum;
-  mwIndex col_id,*row_ids,i,j;
-  int rownumspecified=0;
-  
-  
-  /* Check for proper number of arguments */
-  
-  if (nrhs!=3 && nrhs!=4) {
-    mexErrMsgTxt("Invalid number of arguments.");
-  } else if (nlhs > 1) {
-    mexErrMsgTxt("Too many output arguments.");
-  }
-  
-  
-  /* Check the input dimensions */
-  
-  m = mxGetM(A_IN);
-  n = mxGetN(A_IN);
-  if (!mxIsDouble(A_IN) || mxIsComplex(A_IN) || mxGetNumberOfDimensions(A_IN)>2) {
-    mexErrMsgTxt("COLLINCOMB requires that A be a double matrix.");
-  }
-  
-  if (nrhs==3) {
-    
-    m1 = mxGetM(COLS_IN1);
-    n1 = mxGetN(COLS_IN1);
-    if (!mxIsDouble(COLS_IN1) || mxIsComplex(COLS_IN1) || (m1!=1 && n1!=1)) {
-      mexErrMsgTxt("COLLINCOMB requires that COLS be an index vector of type double.");
-    }
-    colnum = (m1 > n1) ? m1 : n1;   /* the number of columns in the linear combination */
-    
-    m2 = mxGetM(X_IN1);
-    n2 = mxGetN(X_IN1);
-    if (!mxIsDouble(X_IN1) || mxIsComplex(X_IN1) || (m2!=1 && n2!=1)) {
-      mexErrMsgTxt("COLLINCOMB requires that X be a double vector.");
-    }
-    
-    if (m2!=colnum && n2!=colnum) {
-      mexErrMsgTxt("The length of X does not match the number of columns in COLS.");
-    }
-    
-    rows = 0;
-    Y_OUT = mxCreateDoubleMatrix(m, 1, mxREAL);
-    cols = mxGetPr(COLS_IN1);
-    x = mxGetPr(X_IN1);
-  }
-  else {
-    
-    m1 = mxGetM(ROWS_IN);
-    n1 = mxGetN(ROWS_IN);
-    if (!mxIsDouble(ROWS_IN) || mxIsComplex(ROWS_IN) || (m1!=1 && n1!=1)) {
-      mexErrMsgTxt("COLLINCOMB requires that ROWS be an index vector of type double.");
-    }
-    rownum = (m1 > n1) ? m1 : n1;   /* the number of rows in the linear combination */
-    rownumspecified = 1;
-    rows = mxGetPr(ROWS_IN);
-    
-    m1 = mxGetM(COLS_IN2);
-    n1 = mxGetN(COLS_IN2);
-    if (!mxIsDouble(COLS_IN2) || mxIsComplex(COLS_IN2) || (m1!=1 && n1!=1)) {
-      mexErrMsgTxt("COLLINCOMB requires that COLS be an index vector of type double.");
-    }
-    colnum = (m1 > n1) ? m1 : n1;   /* the number of columns in the linear combination */
-    
-    m2 = mxGetM(X_IN2);
-    n2 = mxGetN(X_IN2);
-    if (!mxIsDouble(X_IN2) || mxIsComplex(X_IN2) || (m2!=1 && n2!=1)) {
-      mexErrMsgTxt("COLLINCOMB requires that X be a double vector.");
-    }
-    
-    if (m2!=colnum && n2!=colnum) {
-      mexErrMsgTxt("The length of X does not match the number of columns in COLS.");
-    }
-    
-    Y_OUT = mxCreateDoubleMatrix(rownum, 1, mxREAL);
-    cols = mxGetPr(COLS_IN2);
-    x = mxGetPr(X_IN2);
-  }
-  
-  
-  /* Assign pointers to the various parameters */
-  A = mxGetPr(A_IN);
-  y = mxGetPr(Y_OUT);
-  
-  
-  if (rownumspecified) {
-    
-     /* check row indices */
-    
-    row_ids = (mwIndex*)mxMalloc(rownum*sizeof(mwIndex));
-    
-    for (i=0; i<rownum; ++i) {
-      row_ids[i] = (mwIndex)(rows[i]+0.1)-1;
-      if (row_ids[i]<0 || row_ids[i]>=m) {
-        mexErrMsgTxt("Row index in ROWS is out of range.");
-      }
-    }
-    
-    /* Do the actual computation */
-    for (i=0; i<colnum; ++i) {
-      col_id = (mwIndex)(cols[i]+0.1)-1;
-      if (col_id<0 || col_id>=n) {
-        mexErrMsgTxt("Column index in COLS is out of range.");
-      }
-      for (j=0; j<rownum; ++j) {
-        y[j] += A[m*col_id+row_ids[j]]*x[i];
-      }
-    }
-    
-    mxFree(row_ids);
-  }
-  
-  else {
-    
-    /* Do the actual computation */
-    for (i=0; i<colnum; ++i) {
-      col_id = (mwIndex)(cols[i]+0.1)-1;
-      if (col_id<0 || col_id>=n) {
-        mexErrMsgTxt("Column index in COLS is out of range.");
-      }
-      for (j=0; j<m; ++j) {
-        y[j] += A[m*col_id+j]*x[i];
-      }
-    }
-  }
-  
-  return;
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/collincomb.m
--- a/Problems/private/collincomb.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,27 +0,0 @@
-%COLLINCOMB Linear combination of matrix columns.
-%  Y = COLLINCOMB(A,COLS,X) computes a linear combination of the columns of
-%  the matrix A. The column indices are specified in the vector COLS, and
-%  the correspoinding coefficients are specified in the vector X. The
-%  vectors COLS and X must be of the same length. 
-%  
-%  The call Y = COLLINCOMB(A,COLS,X) is essentially equivalent to
-%
-%         Y = A(:,COLS)*X .
-%
-%  However, it is implemented much more efficiently.
-%
-%  Y = COLLINCOMB(A,ROWS,COLS,X) only works on the rows of A specified
-%  in ROWS, returning a vector of length equal to ROWS. This call is
-%  essentially equivalent to the command
-%
-%         Y = A(ROWS,COLS)*X .
-%
-%  See also ROWLINCOMB.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2009
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/collincomb.mexa64
Binary file Problems/private/collincomb.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/completeOps.m
--- a/Problems/private/completeOps.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,78 +0,0 @@
-function data = completeOps(data)
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: completeOps.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-operators = {};
-flagM = 0; if isfield(data,'M'), flagM = 1; end;
-flagB = 0; if isfield(data,'B'), flagB = 1; end;
-
-if (~flagM) && (~flagB)
-   error('At least one of the operators M or B has be to given.');
-end
-
-% Define measurement matrix if needed
-if ~flagM
-  info   = data.B([],0);
-  data.M = opDirac(info{1});
-else
-  operators{end+1} = data.M;
-end
-
-% Define sparsity basis if needed
-if ~flagB
-  info   = data.M([],0);
-  data.B = opDirac(info{2});
-else
-  operators{end+1} = data.B;
-end
-
-% Define operator A if needed
-if ~isfield(data,'A')
-  if (length(operators) > 1)
-     data.A = opFoG(operators{:});
-  else
-     data.A = operators{1};
-  end
-end
-
-% Define empty solution if needed
-if ~isfield(data,'x0')
-  data.x0 = [];
-end
-
-% Define the operator size and string
-opInfo       = data.A([],0);
-data.sizeA   = [opInfo{1},opInfo{2}];
-opInfo       = data.B([],0);
-data.sizeB   = [opInfo{1},opInfo{2}];
-opInfo       = data.M([],0);
-data.sizeM   = [opInfo{1},opInfo{2}];
-data.op.strA = opToString(data.A);
-data.op.strB = opToString(data.B);
-data.op.strM = opToString(data.M);
-
-% Get the size of the desired signal
-if ~isfield(data,'signalSize')
-   if ~isfield(data,'signal')
-     error(['At least one of the fields signal ', ...
-            'or signalSize must be given.']);
-   end
-   data.signalSize = size(data.signal);
-end
-
-% Reconstruct signal from sparse coefficients
-if ~isfield(data,'reconstruct')
-   data.reconstruct = @(x) reshape(data.B(x,1),data.signalSize);
-end
-
-% Reorder the fields (sort alphabetically)
-m = fieldnames(data);
-m = sort(m);
-dataReorder = struct();
-for i=1:length(m)
-  eval(sprintf('dataReorder.%s = data.%s;',m{i},m{i}));
-end
-
-data = dataReorder;
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/countcover.m
--- a/Problems/private/countcover.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,33 +0,0 @@
-function cnt = countcover(sz,blocksize,stepsize)
-%COUNTCOVER Covering of signal samples by blocks
-%  CNT = COUNTCOVER(SZ,BLOCKSIZE,STEPSIZE) assumes a p-dimensional signal
-%  of size SZ=[N1 N2 ... Np] covered by (possibly overlapping) blocks of
-%  size BLOCKSIZE=[M1 M2 ... Mp]. The blocks start at position (1,1,..,1)
-%  and are shifted between them by steps of size STEPSIZE=[S1 S2 ... Sp].
-%  COUNTCOVER returns a matrix the same size as the signal, containing in
-%  each entry the number of blocks covering that sample.
-%
-%  See also IM2COLSTEP, COL2IMSTEP, IM2COL.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  August 2008
-
-
-cnt = ones(sz);
-for k = 1:length(sz)
-  
-  % this code is modified from function NDGRID, so it computes one
-  % output argument of NDGRID at a time (to conserve memory)
-  ids = (1:sz(k))';
-  s = sz; s(k) = [];
-  ids = reshape(ids(:,ones(1,prod(s))),[length(ids) s]);
-  ids = permute(ids,[2:k 1 k+1:length(sz)]);
-  
-  cnt = cnt .* max( min(floor((ids-1)/stepsize(k)),floor((sz(k)-blocksize(k))/stepsize(k))) - ...
-                    max(ceil((ids-blocksize(k))/stepsize(k)),0) + 1 , 0 );
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/diag_ids.m
--- a/Problems/private/diag_ids.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,30 +0,0 @@
-function id = diag_ids(x,k)
-%DIAG_IDS Indices of matrix diagonal elements.
-%  ID = DIAG_IDS(X) returns the indices of the main diagonal of X.
-%
-%  ID = DIAG_IDS(X,K) returns the indices of the K-th diagonal. K=0
-%  represents the main diagonal, positive values are above the main
-%  diagonal and negative values are below the main diagonal.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  September 2006
-
-
-if (nargin==1), k=0; end
-
-[m,n] = size(x);
-l = max(m,n);
-
-if (k<=0)
-  id = (0:l-1)*m + (1:l) - k;
-else
-  id = (0:l-1)*m + (1:l) + k*m;
-end
-
-if (l-k>m), id = id(1:end-(l-k-m)); end
-if (l+k>n), id = id(1:end-(l+k-n)); end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/dictdist.m
--- a/Problems/private/dictdist.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,61 +0,0 @@
-function [dist,ratio] = dictdist(approxD,D,epsilon)
-%DICTDIST Distance between dictionaries.
-%  [DIST,RATIO] = DICTDIST(APPROXD,D) computes the distance between the
-%  approximate dictionary APPROXD and the true dictionary D, where APPROXD
-%  is NxK and D is NxM.
-%
-%  The distance between the dictionary APPROXD and a single atom A of D is
-%  defined as:
-%
-%      DIST(APPROXD,A) = min  { 1-abs(APPROXD(:,i)' * A) }
-%                         i
-%
-%  The distance between the dictionaries APPROXD and D is defined as:
-%
-%      DIST(APPROXD,D) = sum { dist(APPROXD, D(:,k)) } / M
-%                         k
-%
-%  Note that 0 <= DIST(APPROXD,D) <= 1, where 0 implies that all atoms in D
-%  appear in APPROXD, and 1 implies that the atoms of D are orthogonal to
-%  APPROXD.
-%
-%  The similarity ratio between APPROXD and D is defined as:
-%
-%      RATIO(APPROXD,D) = #(atoms in D that appear in APPROXD) / M
-%
-%  where two atoms are considered identical when DIST(A1,A2) < EPSILON with
-%  EPSILON=0.01 by default. Note that 0 <= RATIO(APPROXD,D) <= 1, where 0
-%  means APPROXD and D have no identical atoms, and 1 means that all atoms
-%  of D appear in APPROXD.
-%
-%  [DIST,RATIO] = DICTDIST(DICT1,DICT2,EPSILON) specifies a different value
-%  for EPSILON.
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  October 2007
-
-
-if (nargin < 3), epsilon = 0.01; end
-
-[n,m] = size(D);
-
-approxD = normcols(approxD*spdiag(sign(approxD(1,:))));
-D = normcols(D*spdiag(sign(D(1,:))));
-
-identical_atoms = 0;
-dist = 0;
-
-for i = 1:m
-  atom = D(:,i);
-  distances = 1-abs(atom'*approxD);
-  mindist = min(distances);
-  dist = dist + mindist;
-  identical_atoms = identical_atoms + (mindist < epsilon);
-end
-
-dist = dist / m;
-ratio = identical_atoms / m;
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/genPDF.m
--- a/Problems/private/genPDF.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,113 +0,0 @@
-function [pdf,val] = genPDF(imSize,p,pctg,distType,radius,disp)
-
-%[pdf,val] = genPDF(imSize,p,pctg [,distType,radius,disp])
-%
-%	generates a pdf for a 1d or 2d random sampling pattern
-%	with polynomial variable density sampling
-%
-%	Input:
-%		imSize - size of matrix or vector
-%		p - power of polynomial
-%		pctg - partial sampling factor e.g. 0.5 for half
-%		distType - 1 or 2 for L1 or L2 distance measure
-%		radius - radius of fully sampled center
-%		disp - display output
-%
-%	Output:
-%		pdf - the pdf
-%		val - min sampling density
-%
-% 
-%	Example:
-%	[pdf,val] = genPDF([128,128],2,0.5,2,0,1);
-%
-%	(c) Michael Lustig 2007
-
-%   This file is used with the kind permission of Michael Lustig
-%   (mlustig@stanford.edu), and originally appeared in the
-%   SparseMRI toolbox, http://www.stanford.edu/~mlustig/ .
-%
-% $Id: genPDF.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-
-if nargin < 4
-	distType = 2;
-end
-
-if nargin < 5
-	radius = 0;
-end
-
-if nargin < 6
-	disp = 0;
-end
-
-
-minval=0;
-maxval=1;
-val = 0.5;
-
-if length(imSize)==1
-	imSize = [imSize,1];
-end
-
-sx = imSize(1);
-sy = imSize(2);
-PCTG = floor(pctg*sx*sy);
-
-
-if sum(imSize==1)==0  % 2D
-	[x,y] = meshgrid(linspace(-1,1,sy),linspace(-1,1,sx));
-	switch distType
-		case 1
-			r = max(abs(x),abs(y));
-		otherwise
-			r = sqrt(x.^2+y.^2);
-			r = r/max(abs(r(:)));			
-	end
-
-else %1d
-	r = abs(linspace(-1,1,max(sx,sy)));
-end
-
-
-
-
-idx = find(r<radius);
-
-pdf = (1-r).^p; pdf(idx) = 1;
-if floor(sum(pdf(:))) > PCTG
-	error('infeasible without undersampling dc, increase p');
-end
-
-% begin bisection
-while(1)
-	val = minval/2 + maxval/2;
-	pdf = (1-r).^p + val; pdf(find(pdf>1)) = 1; pdf(idx)=1;
-	N = floor(sum(pdf(:)));
-	if N > PCTG % infeasible
-		maxval=val;
-	end
-	if N < PCTG % feasible, but not optimal
-		minval=val;
-	end
-	if N==PCTG % optimal
-		break;
-	end
-end
-
-if disp
-	figure, 
-	subplot(211), imshow(pdf)
-	if sum(imSize==1)==0
-		subplot(212), plot(pdf(end/2+1,:));
-	else
-		subplot(212), plot(pdf);
-	end
-end
-
-
-
-
-
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/genSampling.m
--- a/Problems/private/genSampling.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,53 +0,0 @@
-function [minIntrVec,stat,actpctg] = genSampling(pdf,iter,tol)
-
-%[mask,stat,N] = genSampling(pdf,iter,tol)
-%
-% a monte-carlo algorithm to generate a sampling pattern with 
-% minimum peak interference. The number of samples will be
-% sum(pdf) +- tol
-%
-%	pdf - probability density function to choose samples from
-%	iter - number of tries
-%	tol  - the deviation from the desired number of samples in samples
-%
-% returns:
-%	mask - sampling pattern
-%	stat - vector of min interferences measured each try
-%	actpctg    - actual undersampling factor
-%
-%	(c) Michael Lustig 2007
-
-%   This file is used with the kind permission of Michael Lustig
-%   (mlustig@stanford.edu), and originally appeared in the
-%   SparseMRI toolbox, http://www.stanford.edu/~mlustig/ .
-%
-% $Id: genSampling.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-% h = waitbar(0);
-
-pdf(find(pdf>1)) = 1;
-K = sum(pdf(:));
-
-minIntr = 1e99;
-minIntrVec = zeros(size(pdf));
-
-for n=1:iter
-	tmp = zeros(size(pdf));
-	while abs(sum(tmp(:)) - K) > tol
-		tmp = rand(size(pdf))<pdf;
-	end
-	
-	TMP = ifft2(tmp./pdf);
-	if max(abs(TMP(2:end))) < minIntr
-		minIntr = max(abs(TMP(2:end)));
-		minIntrVec = tmp;
-	end
-	stat(n) = max(abs(TMP(2:end)));
-      % waitbar(n/iter,h);
-end
-
-actpctg = sum(minIntrVec(:))/prod(size(minIntrVec));
-
-% close(h);
-
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/im2colstep.c
--- a/Problems/private/im2colstep.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,131 +0,0 @@
-/**************************************************************************
- *
- * File name: im2colstep.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 31.8.2009
- *
- *************************************************************************/
-
-
-#include "mex.h"
-#include <string.h>
-
-
-/* Input Arguments */
-
-#define	X_IN	 prhs[0]
-#define SZ_IN  prhs[1]
-#define S_IN   prhs[2]
-
-
-/* Output Arguments */
-
-#define	B_OUT	plhs[0]
-
-
-void mexFunction(int nlhs, mxArray *plhs[], 
-		             int nrhs, const mxArray*prhs[])
-     
-{ 
-    double *x, *b, *s;
-    mwSize sz[3], stepsize[3], n[3], ndims;
-    mwIndex i, j, k, l, m, blocknum;
-    
-    
-    /* Check for proper number of arguments */
-    
-    if (nrhs < 2 || nrhs > 3) {
-      mexErrMsgTxt("Invalid number of input arguments."); 
-    } else if (nlhs > 1) {
-      mexErrMsgTxt("Too many output arguments."); 
-    } 
-    
-    
-    /* Check the the input dimensions */ 
-    
-    ndims = mxGetNumberOfDimensions(X_IN);
-    
-    if (!mxIsDouble(X_IN) || mxIsComplex(X_IN) || ndims>3) {
-      mexErrMsgTxt("X should be a 2-D or 3-D double matrix.");
-    }
-    if (!mxIsDouble(SZ_IN) || mxIsComplex(SZ_IN) || mxGetNumberOfDimensions(SZ_IN)>2 || mxGetM(SZ_IN)*mxGetN(SZ_IN)!=ndims) {
-      mexErrMsgTxt("Invalid block size.");
-    }
-    if (nrhs == 3) {
-      if (!mxIsDouble(S_IN) || mxIsComplex(S_IN) || mxGetNumberOfDimensions(S_IN)>2 || mxGetM(S_IN)*mxGetN(S_IN)!=ndims) {
-        mexErrMsgTxt("Invalid step size.");
-      }
-    }
-    
-    
-    /* Get parameters */
-    
-    s = mxGetPr(SZ_IN);
-    if (s[0]<1 || s[1]<1 || (ndims==3 && s[2]<1)) {
-      mexErrMsgTxt("Invalid block size.");
-    }
-    sz[0] = (mwSize)(s[0] + 0.01);
-    sz[1] = (mwSize)(s[1] + 0.01);
-    sz[2] = ndims==3 ? (mwSize)(s[2] + 0.01) : 1;
-    
-    if (nrhs == 3) {
-      s = mxGetPr(S_IN);
-      if (s[0]<1 || s[1]<1 || (ndims==3 && s[2]<1)) {
-        mexErrMsgTxt("Invalid step size.");
-      }
-      stepsize[0] = (mwSize)(s[0] + 0.01);
-      stepsize[1] = (mwSize)(s[1] + 0.01);
-      stepsize[2] = ndims==3 ? (mwSize)(s[2] + 0.01) : 1;
-    }
-    else {
-      stepsize[0] = stepsize[1] = stepsize[2] = 1;
-    }
-    
-    n[0] = (mxGetDimensions(X_IN))[0];
-    n[1] = (mxGetDimensions(X_IN))[1];
-    n[2] = ndims==3 ? (mxGetDimensions(X_IN))[2] : 1;
-    
-    if (n[0]<sz[0] || n[1]<sz[1] || (ndims==3 && n[2]<sz[2])) {
-      mexErrMsgTxt("Block size too large.");
-    }
-    
-    
-    /* Create a matrix for the return argument */
-    
-    B_OUT = mxCreateDoubleMatrix(sz[0]*sz[1]*sz[2], ((n[0]-sz[0])/stepsize[0]+1)*((n[1]-sz[1])/stepsize[1]+1)*((n[2]-sz[2])/stepsize[2]+1), mxREAL);
-    
-    
-    /* Assign pointers */
-    
-    x = mxGetPr(X_IN);
-    b = mxGetPr(B_OUT);
-            
-    
-    /* Do the actual computation */
-    
-    blocknum = 0;
-    
-    /* iterate over all blocks */
-    for (k=0; k<=n[2]-sz[2]; k+=stepsize[2]) {
-      for (j=0; j<=n[1]-sz[1]; j+=stepsize[1]) {
-        for (i=0; i<=n[0]-sz[0]; i+=stepsize[0]) {
-          
-          /* copy single block */
-          for (m=0; m<sz[2]; m++) {
-            for (l=0; l<sz[1]; l++) {
-              memcpy(b + blocknum*sz[0]*sz[1]*sz[2] + m*sz[0]*sz[1] + l*sz[0], x+(k+m)*n[0]*n[1]+(j+l)*n[0]+i, sz[0]*sizeof(double));
-            }
-          }
-          blocknum++;
-          
-        }
-      }
-    }
-    
-    return;
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/im2colstep.m
--- a/Problems/private/im2colstep.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,31 +0,0 @@
-%IM2COLSTEP Rearrange matrix blocks into columns.
-%  B = IM2COLSTEP(A,[N1 N2]) converts each sliding N1-by-N2 block of the
-%  2-D matrix A into a column of B, with no zero padding. B has N1*N2 rows
-%  and will contain as many columns as there are N1-by-N2 neighborhoods in
-%  A. Each column of B contains a neighborhood of A reshaped as NHOOD(:),
-%  where NHOOD is a matrix containing an N1-by-N2 neighborhood of A.
-%
-%  B = IM2COLSTEP(A,[N1 N2],[S1 S2]) extracts neighborhoods of A with a
-%  step size of (S1,S2) between them. The first extracted neighborhood is
-%  A(1:N1,1:N2), and the rest are of the form A((1:N1)+i*S1,(1:N2)+j*S2).
-%  Note that to ensure coverage of all A by neighborhoods,
-%  (size(A,i)-Ni)/Si must be whole for i=1,2. The default function behavior
-%  corresponds to [S1 S2] = [1 1]. Setting S1>=N1 and S2>=N2 results in no
-%  overlap between the neighborhoods.
-%
-%  B = IM2COLSTEP(A,[N1 N2 N3],[S1 S2 S3]) operates on a 3-D matrix A. The
-%  step size [S1 S2 S3] may be ommitted, and defaults to [1 1 1].
-%
-%  Note: the call IM2COLSTEP(A,[N1 N2]) produces the same output as
-%  Matlab's IM2COL(A,[N1 N2],'sliding'). However, it is significantly
-%  faster.
-%
-%  See also COL2IMSTEP, IM2COL, COUNTCOVER.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  August 2009
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/im2colstep.mexa64
Binary file Problems/private/im2colstep.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/imnormalize.m
--- a/Problems/private/imnormalize.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,19 +0,0 @@
-function y = imnormalize(x)
-%IMNORMALIZE Normalize image values.
-%  Y = IMNORMALIZE(X) linearly transforms the intensity values of the image
-%  X to tightly cover the range [0,1]. If X has more than one channel, the
-%  channels are handled as one and normalized using the same transform.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  May 2004
-
-
-maxval = max(x(:));
-minval = min(x(:));
-
-y = (x-minval) / (maxval-minval);
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/iswhole.m
--- a/Problems/private/iswhole.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,17 +0,0 @@
-function z = iswhole(x,epsilon)
-%ISWHOLE Determine whole numbers with finite precision.
-%  Z = ISWHOLE(X,EPSILON) returns a matrix the same size as X, containing
-%  1's where X is whole up to precision EPSILON, and 0's elsewhere. 
-%
-%  Z = ISWHOLE(X) uses the default value of EPSILON=1e-6.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  May 2005
-
-if (nargin<2), epsilon = 1e-6; end
-z = abs(round(x)-x) < epsilon;
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/make.m
--- a/Problems/private/make.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,40 +0,0 @@
-function make
-%MAKE Build the KSVDBox MEX support files.
-%  MAKE compiles the KSVDBox supporting MEX functions, using Matlab's
-%  default MEX compiler. If the MEX compiler has not been set-up before,
-%  please run 
-%
-%    mex -setup
-%
-%  before using this MAKE file.
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  August 2009
-
-
-% detect platform 
-
-compstr = computer;
-is64bit = strcmp(compstr(end-1:end),'64');
-
-
-% compilation parameters
-
-compile_params = cell(0);
-if (is64bit)
-  compile_params{1} = '-largeArrayDims';
-end
-
-
-% Compile files %
-
-sourcefiles = {{'addtocols.c'}, {'collincomb.c'}, {'rowlincomb.c'}, {'sprow.c','mexutils.c'}, {'im2colstep.c'}, {'col2imstep.c'}};
-
-for i = 1:length(sourcefiles)
-  printf('Compiling %s...', sourcefiles{i}{1});
-  mex(sourcefiles{i}{:},compile_params{:});
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/mexutils.c
--- a/Problems/private/mexutils.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,79 +0,0 @@
-/**************************************************************************
- *
- * File name: mexutils.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 15.8.2009
- *
- *************************************************************************/
-
-#include "mexutils.h"
-#include <math.h>
-
-
-
-/* verify that the mxArray contains a double matrix */
-
-void checkmatrix(const mxArray *param, char *fname, char *pname)
-{
-  char errmsg[100];
-  sprintf(errmsg, "%.15s requires that %.25s be a double matrix.", fname, pname);
-  if (!mxIsDouble(param) || mxIsComplex(param) || mxGetNumberOfDimensions(param)>2) {
-    mexErrMsgTxt(errmsg);
-  }
-}
-
-
-/* verify that the mxArray contains a 1-D double vector */
-
-void checkvector(const mxArray *param, char *fname, char *pname)
-{
-  char errmsg[100];
-  sprintf(errmsg, "%.15s requires that %.25s be a double vector.", fname, pname);
-  if (!mxIsDouble(param) || mxIsComplex(param) || mxGetNumberOfDimensions(param)>2 || (mxGetM(param)!=1 && mxGetN(param)!=1)) {
-    mexErrMsgTxt(errmsg);
-  }
-}
-
-
-/* verify that the mxArray contains a double scalar */
-
-void checkscalar(const mxArray *param, char *fname, char *pname)
-{
-  char errmsg[100];
-  sprintf(errmsg, "%.15s requires that %.25s be a double scalar.", fname, pname);
-  if (!mxIsDouble(param) || mxIsComplex(param) || mxGetNumberOfDimensions(param)>2 || 
-      mxGetM(param)!=1 || mxGetN(param)!=1) 
-  {
-    mexErrMsgTxt(errmsg);
-  }
-}
-
-
-/* verify that the mxArray contains a sparse matrix */
-
-void checksparse(const mxArray *param, char *fname, char *pname)
-{
-  char errmsg[100];
-  sprintf(errmsg, "%.15s requires that %.25s be sparse.", fname, pname);
-  if (!mxIsSparse(param)) {
-    mexErrMsgTxt(errmsg);
-  }
-}
-
-
-/* verify that the mxArray contains a 1-dimensional cell array */
-
-void checkcell_1d(const mxArray *param, char *fname, char *pname)
-{
-  char errmsg[100];
-  sprintf(errmsg, "%.15s requires that %.25s be a 1-D cell array.", fname, pname);
-  if (!mxIsCell(param) || mxGetNumberOfDimensions(param)>2 || (mxGetM(param)!=1 && mxGetN(param)!=1)) {
-    mexErrMsgTxt(errmsg);
-  }
-}
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/mexutils.h
--- a/Problems/private/mexutils.h	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,103 +0,0 @@
-/**************************************************************************
- *
- * File name: mexutils.h
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 18.8.2009
- *
- * Utility functions for MEX files.
- *
- *************************************************************************/
-
-
-#ifndef __MEX_UTILS_H__
-#define __MEX_UTILS_H__
-
-#include "mex.h"
-
-
-
-/**************************************************************************
- * Function checkmatrix:
- *
- * Verify that the specified mxArray is real, of type double, and has 
- * no more than two dimensions. If not, an error message is printed
- * and the mex file terminates.
- * 
- * Parameters:
- *   param - the mxArray to be checked
- *   fname - the name of the function where the error occured (15 characters or less)
- *   pname - the name of the parameter (25 characters or less)
- *
- **************************************************************************/
-void checkmatrix(const mxArray *param, char *fname, char *pname);
-
-
-/**************************************************************************
- * Function checkvector:
- *
- * Verify that the specified mxArray is 1-D, real, and of type double. The
- * vector may be a column or row vector. Otherwise, an error message is
- * printed and the mex file terminates.
- * 
- * Parameters:
- *   param - the mxArray to be checked
- *   fname - the name of the function where the error occured (15 characters or less)
- *   pname - the name of the parameter (25 characters or less)
- *
- **************************************************************************/
-void checkvector(const mxArray *param, char *fname, char *pname);
-
-
-/**************************************************************************
- * Function checkscalar:
- *
- * Verify that the specified mxArray represents a real double scalar value. 
- * If not, an error message is printed and the mex file terminates.
- * 
- * Parameters:
- *   param - the mxArray to be checked
- *   fname - the name of the function where the error occured (15 characters or less)
- *   pname - the name of the parameter (25 characters or less)
- *
- **************************************************************************/
-void checkscalar(const mxArray *param, char *fname, char *pname);
-
-
-/**************************************************************************
- * Function checksparse:
- *
- * Verify that the specified mxArray contains a sparse matrix. If not,
- * an error message is printed and the mex file terminates.
- * 
- * Parameters:
- *   param - the mxArray to be checked
- *   fname - the name of the function where the error occured (15 characters or less)
- *   pname - the name of the parameter (25 characters or less)
- *
- **************************************************************************/
-void checksparse(const mxArray *param, char *fname, char *pname);
-
-
-/**************************************************************************
- * Function checkcell_1d:
- *
- * Verify that the specified mxArray is a 1-D cell array. The cell array 
- * may be arranged as either a column or a row. If not, an error message 
- * is printed and the mex file terminates.
- * 
- * Parameters:
- *   param - the mxArray to be checked
- *   fname - the name of the function where the error occured (15 characters or less)
- *   pname - the name of the parameter (25 characters or less)
- *
- **************************************************************************/
-void checkcell_1d(const mxArray *param, char *fname, char *pname);
-
-
-#endif
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/normcols.m
--- a/Problems/private/normcols.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,17 +0,0 @@
-function y = normcols(x)
-%NORMCOLS Normalize matrix columns.
-%  Y = NORMCOLS(X) normalizes the columns of X to unit length, returning
-%  the result as Y.
-%
-%  See also ADDTOCOLS.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2009
-
-
-y = x*spdiag(1./sqrt(sum(x.*x)));
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/printf.m
--- a/Problems/private/printf.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,26 +0,0 @@
-function str = printf(varargin)
-%PRINTF Print formatted text to screen.
-%  PRINTF(FMT,VAL1,VAL2,...) formats the data in VAL1,VAL2,... according to
-%  the format string FMT, and prints the result to the screen.
-%
-%  The call to PRINTF(FMT,VAL1,VAL2,...) simply invokes the call
-%  DISP(SPRINTF(FMT,VAL1,VAL2,...)). For a complete description of the
-%  format string options see function SPRINTF.
-%
-%  STR = PRINTF(...) also returns the formatted string.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2008
-
-
-if (nargout>0)
-  str = sprintf(varargin{:});
-  disp(str);
-else
-  disp(sprintf(varargin{:}));
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/pwsmoothfield.m
--- a/Problems/private/pwsmoothfield.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,32 +0,0 @@
-function f = pwsmoothfield(n,var,alpha)
-% PWSMOOTHFIELD(N,VAR,ALPHA)
-%   Generate an image of piecewise smooth filtered white noise
-%
-%   N = sidelength of the field in pixels
-%   VAR = variance of original white nose
-%   ALPHA = fraction of FFT coefficents to keep
-%
-%   Returns an N-by-N array
-%
-%   This file is used with the kind permission of Stephen J. Wright
-%   (swright@cs.wisc.edu), and was originally included in the GPSR
-%   v1.0 distribution: http://www.lx.it.pt/~mtf/GPSR .
-
-% $Id: pwsmoothfield.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-f = sqrt(var)*randn(n);
-F = fft2(f);
-a = ceil(n*alpha/2);
-b = fix(n*(1-alpha));
-F(a+1:a+b,:) = 0;
-F(:,a+1:a+b) = 0;
-f = real(ifft2(F));
-
-for i = 1:n
-    for j = 1:n
-        if (j/n >= 15*(i/n - 0.5)^3 + 0.4)
-            f(i,j) = f(i,j) + sqrt(var);
-        end
-    end
-end
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/reggrid.m
--- a/Problems/private/reggrid.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,136 +0,0 @@
-function [varargout] = reggrid(sz,num,mode)
-%REGGRID Regular sampling grid.
-%  [I1,I2,...,Ip] = REGGRID([N1 N2 ... Np], NUM) returns the indices
-%  of a regular uniform sampling grid over a p-dimensional matrix with
-%  dimensions N1xN2x...xNp. NUM is the minimal number of required samples,
-%  and it is ensured that the actual number of samples, given by
-%  length(I1)xlength(I2)x...xlength(Ip), is at least as large as NUM.
-%
-%  [I1,I2,...,Ip] = REGGRID([N1 N2 ... Np], NUM,'MODE') specifies the
-%  method for distributing the samples along each dimension. Valid modes
-%  include 'eqdist' (the default mode) and 'eqnum'. 'eqdist' indicates an
-%  equal distance between the samples in each dimension, while 'eqnum'
-%  indicates an equal number of samples in each dimension.
-%
-%  Notes about MODE:
-%
-%    1. The 'eqnum' mode will generally fail when the p-th root of NUM
-%    (i.e. NUM^(1/p)) is larger than min([N1 N2 ... Np]). Thus 'eqdist' is
-%    the more useful choice for sampling an arbitrary number of samples
-%    from the matrix (up to the total number of matrix entries).
-%  
-%    2. In both modes, the equality (of the distance between samples, or
-%    the number of samples in each dimension) is only approximate. This is
-%    because REGGRID attempts to maintain the appropriate equality while at
-%    the same time find a sampling pattern where the total number of
-%    samples is as close as possible to NUM. In general, the larger {Ni}
-%    and NUM are, the tighter the equality.
-%
-%  Example: Sample a set of blocks uniformly from a 2D image.
-%
-%    n = 512; blocknum = 20000; blocksize = [8 8];
-%    im = rand(n,n);
-%    [i1,i2] = reggrid(size(im)-blocksize+1, blocknum);
-%    blocks = sampgrid(im, blocksize, i1, i2);
-%
-%  See also SAMPGRID.
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  November 2007
-
-dim = length(sz);
-
-if (nargin<3)
-  mode = 'eqdist';
-end
-
-if (any(sz<1))
-  error(['Invalid matrix size : [' num2str(sz) ']']);
-end
-
-if (num > prod(sz))
-  warning(['Invalid number of samples, returning maximum number of samples.']);
-elseif (num <= 0)
-  if (num < 0)
-    warning('Invalid number of samples, assuming 0 samples.');
-  end
-  for i = 1:length(sz)
-    varargout{i} = [];
-  end
-  return;
-end
-
-
-if (strcmp(mode,'eqdist'))
-  
-  % approximate distance between samples: total volume divided by number of
-  % samples gives the average volume per sample. then, taking the p-th root
-  % gives the average distance between samples
-  d = (prod(sz)/num)^(1/dim);
-  
-  % compute the initial guess for number of samples in each dimension.
-  % then, while total number of samples is too large, decrese the number of
-  % samples by one in the dimension where the samples are the most crowded.
-  % finally, do the opposite process until just passing num, so the final
-  % number of samples is the closest to num from above.
-  
-  n = min(max(round(sz/d),1),sz);   % set n so that it saturates at 1 and sz
-  
-  active_dims = find(n>1);    % dimensions where the sample num can be reduced
-  while(prod(n)>num && ~isempty(active_dims))
-    [y,id] = min((sz(active_dims)-1)./n(active_dims));
-    n(active_dims(id)) = n(active_dims(id))-1;
-    if (n(active_dims(id)) < 2)
-      active_dims = find(n>1);
-    end
-  end
-
-  active_dims = find(n<sz);    % dimensions where the sample num can be increased
-  while(prod(n)<num && ~isempty(active_dims))
-    [y,id] = max((sz(active_dims)-1)./n(active_dims));
-    n(active_dims(id)) = n(active_dims(id))+1;
-    if (n(active_dims(id)) >= sz(active_dims(id)))
-      active_dims = find(n<sz);
-    end
-  end
-
-  for i = 1:dim
-    varargout{i} = round((1:n(i))/n(i)*sz(i));
-    varargout{i} = varargout{i} - floor((varargout{i}(1)-1)/2);
-  end
-  
-elseif (strcmp(mode,'eqnum'))
-  
-  % same idea as above
-  n = min(max( ones(size(sz)) * round(num^(1/dim)) ,1),sz);
-
-  active_dims = find(n>1);
-  while(prod(n)>num && ~isempty(active_dims))
-    [y,id] = min((sz(active_dims)-1)./n(active_dims));
-    n(active_dims(id)) = n(active_dims(id))-1;
-    if (n(active_dims(id)) < 2)
-      active_dims = find(n>1);
-    end
-  end
-  
-  active_dims = find(n<sz);
-  while(prod(n)<num && ~isempty(active_dims))
-    [y,id] = max((sz(active_dims)-1)./n(active_dims));
-    n(active_dims(id)) = n(active_dims(id))+1;
-    if (n(active_dims(id)) >= sz(active_dims(id)))
-      active_dims = find(n<sz);
-    end
-  end
-  
-  for i = 1:dim
-    varargout{i} = round((1:n(i))/n(i)*sz(i));
-    varargout{i} = varargout{i} - floor((varargout{i}(1)-1)/2);
-  end
-else
-  error('Invalid sampling mode');
-end
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/remove_dc.m
--- a/Problems/private/remove_dc.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,31 +0,0 @@
-function [y,dc] = remove_dc(x,columns)
-%REMOVE_DC Remove DC channel from signals.
-%   [Y,DC] = REMOVE_DC(X) removes the DC channel (i.e. the mean) from the
-%   specified (possibly multi-dimensional) signal X. Y is the DC-free
-%   signal and is the same size as X. DC is a scalar containing the mean of
-%   the signal X.
-%
-%   [Y,DC] = REMOVE_DC(X,'columns') where X is a 2D matrix, treats the
-%   columns of X as a set of 1D signals, removing the DC channel from each
-%   one individually. Y is the same size as X and contains the DC-free
-%   signals. DC is a row vector of length size(X,2) containing the means of
-%   the signals in X.
-%
-%   See also ADD_DC.
-
-
-if (nargin==2 && strcmpi(columns,'columns')), columns = 1;
-else columns = 0;
-end
-
-if (columns)
-  dc = mean(x);
-  y = addtocols(x,-dc);
-else
-  if (ndims(x)==2)  % temporary, will remove in future
-    warning('Treating 2D matrix X as a single signal and not each column individually');
-  end
-  dc = mean(x(:));
-  y = x-dc;
-end
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/rowlincomb.c
--- a/Problems/private/rowlincomb.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,148 +0,0 @@
-/**************************************************************************
- *
- * File name: rowlincomb.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 21.5.2009
- *
- *************************************************************************/
-
-#include "mex.h"
-
-
-/* Input Arguments */
-
-#define	X_IN	   prhs[0]
-#define A_IN     prhs[1]
-#define ROWS_IN  prhs[2]
-#define COLS_IN  prhs[3]
-
-
-/* Output Arguments */
-
-#define	Y_OUT	plhs[0]
-
-
-void mexFunction(int nlhs, mxArray *plhs[],
-int nrhs, const mxArray*prhs[])
-
-{
-  double *A, *x, *y, *rows, *cols;
-  mwSize m,n,m1,n1,m2,n2,rownum,colnum;
-  mwIndex *row_ids,*col_ids,i,j;
-  int colnumspecified=0;
-  
-  
-  /* Check for proper number of arguments */
-  
-  if (nrhs!=3 && nrhs!=4) {
-    mexErrMsgTxt("Invalid number of input arguments.");
-  } else if (nlhs > 1) {
-    mexErrMsgTxt("Too many output arguments.");
-  }
-  
-  
-  /* Check the input dimensions */
-  
-  m = mxGetM(A_IN);
-  n = mxGetN(A_IN);
-  if (!mxIsDouble(A_IN) || mxIsComplex(A_IN) || mxGetNumberOfDimensions(A_IN)>2) {
-    mexErrMsgTxt("ROWLINCOMB requires that A be a double matrix.");
-  }
-  
-  m1 = mxGetM(ROWS_IN);
-  n1 = mxGetN(ROWS_IN);
-  if (!mxIsDouble(ROWS_IN) || mxIsComplex(ROWS_IN) || (m1!=1 && n1!=1)) {
-    mexErrMsgTxt("ROWLINCOMB requires that ROWS be an index vector of type double.");
-  }
-  rownum = (m1 > n1) ? m1 : n1;   /* the number of rows in the linear combination */
-  
-  m2 = mxGetM(X_IN);
-  n2 = mxGetN(X_IN);
-  if (!mxIsDouble(X_IN) || mxIsComplex(X_IN) || ((m2!=1) && (n2!=1))) {
-    mexErrMsgTxt("ROWLINCOMB requires that X be a double vector.");
-  }
-  
-  if (m2 != rownum && n2 != rownum) {
-    mexErrMsgTxt("The length of X does not match the number of rows in ROWS.");
-  }
-  
-  if (nrhs==4) {
-    m1 = mxGetM(COLS_IN);
-    n1 = mxGetN(COLS_IN);
-    if (!mxIsDouble(COLS_IN) || mxIsComplex(COLS_IN) || (m1!=1 && n1!=1)) {
-      mexErrMsgTxt("ROWLINCOMB requires that COLS be an index vector of type double.");
-    }
-    colnum = (m1 > n1) ? m1 : n1;   /* the number of columns */
-    colnumspecified = 1;
-    cols = mxGetPr(COLS_IN);
-    
-    Y_OUT = mxCreateDoubleMatrix(1, colnum, mxREAL);
-  }
-  else {
-    cols = 0;
-    Y_OUT = mxCreateDoubleMatrix(1, n, mxREAL);
-  }
-  
-  
-  /* Assign pointers to the various parameters */
-  A = mxGetPr(A_IN);
-  rows = mxGetPr(ROWS_IN);
-  x = mxGetPr(X_IN);
-  y = mxGetPr(Y_OUT);
-  
-  
-  /* check row indices */
-  
-  row_ids = (mwIndex*)mxMalloc(rownum*sizeof(mwIndex));
-  
-  for (i=0; i<rownum; ++i) {
-    row_ids[i] = (mwIndex)(rows[i]+0.1)-1;
-    if (row_ids[i]<0 || row_ids[i]>=m) {
-      mexErrMsgTxt("Row index in ROWS is out of range.");
-    }
-  }
-  
-  
-  
-  if (colnumspecified) {
-    
-    /* check column indices */
-    col_ids = (mwIndex*)mxMalloc(colnum*sizeof(mwIndex));
-    
-    for (i=0; i<colnum; ++i) {
-      col_ids[i] = (mwIndex)(cols[i]+0.1)-1;
-      if (col_ids[i]<0 || col_ids[i]>=n) {
-        mexErrMsgTxt("Column index in COLS is out of range.");
-      }
-    }
-    
-    /* Do the actual computation */
-    for (j=0; j<colnum; ++j) {
-      for (i=0; i<rownum; ++i) {
-        y[j] += A[m*col_ids[j]+row_ids[i]]*x[i];
-      }
-    }
-    
-    mxFree(col_ids);
-  }
-  
-  else {
-    
-    /* Do the actual computation */
-    for (j=0; j<n; ++j) {
-      for (i=0; i<rownum; ++i) {
-        y[j] += A[m*j+row_ids[i]]*x[i];
-      }
-    }
-  }
-  
-  
-  mxFree(row_ids);
-  
-  return;
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/rowlincomb.m
--- a/Problems/private/rowlincomb.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,26 +0,0 @@
-%ROWLINCOMB Linear combination of matrix rows.
-%  Y = ROWLINCOMB(X,A,ROWS) computes a linear combination of the rows of
-%  the matrix A. The row indices are specified in the vector ROWS, and the
-%  correspoinding coefficients are specified in the vector X. The vectors
-%  ROWS and X must be of the same length. The call Y = ROWLINCOMB(X,A,ROWS)
-%  is essentially equivalent to the command
-%
-%         Y = X'*A(ROWS,:) .
-%
-%  However, it is implemented much more efficiently.
-%
-%  Y = ROWLINCOMB(X,A,ROWS,COLS) only works on the columns of A specified
-%  in COLS, returning a vector of length equal to COLS. This call is
-%  essentially equivalent to the command
-%
-%         Y = X'*A(ROWS,COLS) .
-%
-%  See also COLLINCOMB.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2009
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/rowlincomb.mexa64
Binary file Problems/private/rowlincomb.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/sampgrid.m
--- a/Problems/private/sampgrid.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,72 +0,0 @@
-function y = sampgrid(x,blocksize,varargin)
-%SAMPGRID Sample a multi-dimensional matrix on a regular grid.
-%  Y = SAMPGRID(X,BLOCKSIZE,I1,I2,...,Ip) extracts block samples of size
-%  BLOCKSIZE from the p-dimensional matrix X, arranging the samples as the
-%  column vectors of the matrix Y. The locations of the (1,1,..,1)-th
-%  elements of each block are given in the index vectors I1,I2,..Ip. The
-%  total number of samples taken is length(I1)xlength(I2)x...xlength(Ip).
-%  BLOCKSIZE should either be a p-element vector of the form [N1,N2,...Np],
-%  or a scalar N which is shorthand for the square block size [N N ... N].
-%
-%  Example: Sample a set of blocks uniformly from a 2D image.
-%
-%    n = 512; blocknum = 20000; blocksize = [8 8];
-%    im = rand(n,n);
-%    [i1,i2] = reggrid(size(im)-blocksize+1, blocknum);
-%    blocks = sampgrid(im, blocksize, i1, i2);
-%
-%  See also REGGRID.
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  November 2007
-
-
-p = ndims(x);
-if (p==2 && any(size(x)==1) && length(blocksize)==1)
-  p = 1;
-end
-
-if (numel(blocksize)==1)
-  blocksize = ones(1,p)*blocksize;
-end
-
-n = zeros(1,p);
-for i = 1:p
-  n(i) = length(varargin{i});
-end
-
-nsamps = prod(n);
-
-% create y of the same class as x
-y = zeros(prod(blocksize),nsamps,class(x));
-
-% ids() contains the index of the current block in I1..Ip
-ids = ones(p,1);
-
-% block_ids contains the indices of the current block in X
-block_ids = cell(p,1);
-for j = 1:p
-  block_ids{j} = varargin{j}(1) : varargin{j}(1)+blocksize(j)-1;
-end
-
-for k = 1:nsamps
-  block = x(block_ids{:});
-  y(:,k) = block(:);
-  
-  % increment ids() and block_ids{}
-  if (k<nsamps)
-    j = 1;
-    while (ids(j) == n(j))
-      ids(j) = 1;
-      block_ids{j} = varargin{j}(1) : varargin{j}(1)+blocksize(j)-1;
-      j = j+1;
-    end
-    ids(j) = ids(j)+1;
-    block_ids{j} = varargin{j}(ids(j)) : varargin{j}(ids(j))+blocksize(j)-1;
-  end
-end
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/scalarToRGB.m
--- a/Problems/private/scalarToRGB.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,18 +0,0 @@
-function s = scalarToRGB(x,colors)
-% input values are assumed to lie between 0 and 1
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: scalarToRGB.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-l  = size(colors,1);
-m  = size(x,1);
-n  = size(x,2);
-s  = zeros(m,n,3);
-
-for i=1:m
-  for j=1:n
-     idx = max(1,min(l,1+floor((l-1) * x(i,j))));
-     s(i,j,:) = colors(idx,:);
-  end
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/secs2hms.m
--- a/Problems/private/secs2hms.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,28 +0,0 @@
-function [h,m,s] = secs2hms(t)
-%SECS2HMS Convert seconds to hours, minutes and seconds.
-%  [H,M,S] = SECS2HMS(T) converts the specified number of seconds T to
-%  hours, minutes and seconds. H and M are whole numbers, and S is real.
-%
-%  Example: Estimate the remaining time of a loop
-%
-%    n = 10; tic;
-%    for i = 1:n
-%      pause(1);
-%      [h,m,s] = secs2hms( (n-i)*toc/i );
-%      printf('estimated remaining time: %02d:%02d:%05.2f',h,m,s);
-%    end
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  April 2008
-
-
-s = t;
-h = fix(s/3600);
-s = rem(s,3600);
-m = fix(s/60);
-s = rem(s,60);
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/spdiag.m
--- a/Problems/private/spdiag.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,38 +0,0 @@
-function Y = spdiag(V,K)
-%SPDIAG Sparse diagonal matrices.
-%   SPDIAG(V,K) when V is a vector with N components is a sparse square
-%   matrix of order N+ABS(K) with the elements of V on the K-th diagonal. 
-%   K = 0 is the main diagonal, K > 0 is above the main diagonal and K < 0
-%   is below the main diagonal. 
-%
-%   SPDIAG(V) is the same as SPDIAG(V,0) and puts V on the main diagonal.
-%
-%   See also DIAG, SPDIAGS.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  June 2008
-
-
-if (nargin<2)
-  K = 0;
-end
-
-n = length(V) + abs(K);
-
-if (K>0)
-  i = 1:length(V);
-  j = K+1:n;
-elseif (K<0)
-  i = -K+1:n;
-  j = 1:length(V);
-else
-  i = 1:n;
-  j = 1:n;
-end
-
-Y = sparse(i,j,V(:),n,n);
\ No newline at end of file
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/sprow.c
--- a/Problems/private/sprow.c	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,94 +0,0 @@
-/**************************************************************************
- *
- * File name: sprow.c
- *
- * Ron Rubinstein
- * Computer Science Department
- * Technion, Haifa 32000 Israel
- * ronrubin@cs
- *
- * Last Updated: 24.8.2009
- *
- *************************************************************************/
-
-
-#include "mex.h"
-#include "mexutils.h"
-
-
-/* Input Arguments */
-
-#define	A_IN	prhs[0]
-#define J_IN  prhs[1]
-
-
-/* Output Arguments */
-
-#define	X_OUT	  plhs[0]
-#define	ID_OUT	plhs[1]
-
-
-void mexFunction(int nlhs, mxArray *plhs[], 
-		             int nrhs, const mxArray*prhs[])
-     
-{ 
-  double *pr, *x, *id, rowid;
-  mwIndex *ir, *jc;
-  mwSize m, n;
-  mwIndex i, j, k, l, rowlen;
-  
-  if (nrhs != 2) {
-    mexErrMsgTxt("GETSPROW requires two input arguments.");
-  } else if (nlhs > 2) {
-    mexErrMsgTxt("Too many output arguments.");
-  }
-  
-  checkmatrix(A_IN, "GETSPROW", "A");
-  checksparse(A_IN, "GETSPROW", "A");
-  checkscalar(J_IN, "GETSPROW", "J");
-  
-  m = mxGetM(A_IN);
-  n = mxGetN(A_IN);
-  
-  rowid = mxGetScalar(J_IN);
-  if (rowid < 0) {
-    mexErrMsgTxt("Invalid row index.");
-  }
-  j = (mwIndex)(rowid + 1e-2);
-  if (j<1 || j>m) {
-    mexErrMsgTxt("Row index is out of range.");
-  }
-  j--;
-  
-  pr = mxGetPr(A_IN);
-  ir = mxGetIr(A_IN);
-  jc = mxGetJc(A_IN);
-  
-  /* Determine length of row */
-  rowlen = 0;
-  for (i=0; i<jc[n]; ++i) {
-    rowlen += (ir[i]==j);
-  }
-  
-  /* Allocate output parameters */
-  X_OUT = mxCreateDoubleMatrix(1, rowlen, mxREAL);
-  ID_OUT = mxCreateDoubleMatrix(1, rowlen, mxREAL);
-  
-  x = mxGetPr(X_OUT);
-  id = mxGetPr(ID_OUT);
-  
-  /* Compute j-th row */
-  k=0;
-  for (l=1; l<=n; ++l) {
-    i = jc[l-1];
-    while (i<jc[l] && ir[i]<j) {
-      i++;
-    }
-    if (i<jc[l] && ir[i]==j) {
-      x[k] = pr[i];
-      id[k] = l;
-      k++;
-    }
-  }
-  
-}
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/sprow.m
--- a/Problems/private/sprow.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,19 +0,0 @@
-%SPROW Extract row of sparse matrix.
-%  X = SPROW(A,J) where A is a sparse matrix, returns the nonzero values in
-%  the row A(J,:).
-%
-%  [X,ID] = SPROW(A,J) also returns the column indices of the nonzeros.
-%
-%  Note that the call [X,ID] = SPROW(A,J) is equivalent (but more efficient
-%  than) the Matlab code
-%
-%    IDS = find(A(J,:)); 
-%    X = A(J,IDS);
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  August 2009
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/sprow.mexa64
Binary file Problems/private/sprow.mexa64 has changed
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/thumbFromOp.m
--- a/Problems/private/thumbFromOp.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,42 +0,0 @@
-function P = thumbFromOp(op,m,n,sm,sn,grayscale)
-% Output matrix P, of size m x n, sampled at
-% sm x sn upper top
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: thumbFromOp.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-if nargin < 6, grayscale = 0; end
-    
-info = op([],0);
-
-sm = min(info{1},sm);
-sn = min(info{2},sn);
-M  = zeros(sm,sn);
-
-for i=1:sn
-   v = zeros(info{2},1); v(i) = 1;
-   w = real(op(v,1));
-
-   M(:,i) = w(1:sm);
-end
-
-mn = min(min(M));
-mx = max(max(M));
-M  = (M - mn) / (mx-mn);
-
-idxm = floor(linspace(1,sm+1,m+1)); idxm = idxm(1:end-1);
-idxn = floor(linspace(1,sn+1,n+1)); idxn = idxn(1:end-1);
-
-if grayscale
-   P = 1-M(idxm,idxn);
-else
-   clrmap = hsv;
-   M = 1 + round(M * (length(clrmap)-1));
-   P = zeros(m,n,3);
-   for j1=1:m
-       for j2=1:n
-           P(j1,j2,:) = clrmap(M(idxm(j1),idxn(j2)),:);
-       end
-   end
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/thumbPlot.m
--- a/Problems/private/thumbPlot.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,38 +0,0 @@
-function P = thumbPlot(P,x,y,color)
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: thumbPlot.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-m = size(P,1);
-n = size(P,2);
-if (size(P,3) == 0) & (length(color) == 3)
-  % Convert to gray-scale
-  color = 0.30*color(1) + 0.59*color(2) + 0.11*color(3);
-end
-
-mnx = min(x);  % Minimum x
-mxx = max(x);  % Maximum x
-mny = min(y);  % Minimum y
-mxy = max(y);  % Maximum y
-dy  = (mxy - mny) * 0.1;   % Offset on vertical axis
-sx  = (mxx - mnx) * 1.0;   % Scale of horizontal axis
-sy  = (mxy - mny) * 1.2;   % Scale of vertical axis
-
-if (sx < 1e-6), sx = 1; end
-if (sy < 1e-6), sy = 1; end
-
-for i=1:length(x)-1
-   x0 = floor(1 + (n-1) * (x(i  ) - mnx) / sx);
-   x1 = floor(1 + (n-1) * (x(i+1) - mnx) / sx);
-   y0 = floor(    (n-1) * (y(i  ) - mny + dy) / sy);
-   y1 = floor(    (n-1) * (y(i+1) - mny + dy) / sy);
-   
-   samples = 1+2*max(abs(x1-x0)+1,abs(y1-y0)+1);
-   c       = linspace(0,1,samples);
-   idx     = round((1-c)*x0 + c*x1);
-   idy     = n - round((1-c)*y0 + c*y1);
-   for j=1:samples
-      P(idy(j),idx(j),:) = color;
-   end
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/thumbwrite.m
--- a/Problems/private/thumbwrite.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,18 +0,0 @@
-function thumbwrite(data,name,opts)
-
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: thumbwrite.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-%
-% data       Thumbnail data in range 0-1
-% name       Name of file (no extension)
-% opts
-%   .thumbtype  Type of image (png, eps, ps, ...)
-%   .thumbdir   Output directory
-%
-
-[type,ext] = getFigureExt(opts.thumbtype);
-data = round(data * 255) / 255;
-imwrite(data,[opts.thumbpath,name,'.',ext],type);
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/timerclear.m
--- a/Problems/private/timerclear.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,37 +0,0 @@
-function timerclear()
-%TIMERCLEAR Clear all timers.
-%   TIMERCLEAR clears all currenly registered timers, invalidating all
-%   timer ids.
-%
-%   Note: since registered timers do not consume CPU power except for when
-%   the TIMER<*> functions are called, this function is only useful in
-%   situations where a large number of timers have been initialized, and
-%   there is a need to reclaim memory.
-%
-%   See also TIMERINIT, TIMERETA.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  June 2008
-
-
-global utiltbx_timer_start_times    % start times
-global utiltbx_time_lastdisp        % last display times
-global utiltbx_timer_iternums       % iteration numbers
-global utiltbx_timer_lastiter       % last queried iteration numbers
-global utiltbx_timer_name           % timer names
-global utiltbx_timer_callfun        % timer calling functions
-
-
-% clear all timers %
-
-utiltbx_timer_start_times = [];
-utiltbx_time_lastdisp = [];
-utiltbx_timer_iternums = [];
-utiltbx_timer_lastiter = [];
-utiltbx_timer_name = [];
-utiltbx_timer_callfun = [];
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/timereta.m
--- a/Problems/private/timereta.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,98 +0,0 @@
-function varargout = timereta(tid,iter,delay)
-%TIMERETA Estimated remaining time.
-%   S = TIMERETA(TID,ITER) returns the estimated remaining time (in
-%   seconds) for the process associated with timer TID, assuming the
-%   process has completed ITER iterations. Note: the function will exit
-%   with an error if the timer TID does not have an associated number of
-%   iterations (see function TIMERINIT).
-%
-%   [H,M,S] = TIMERETA(TID,ITER) returns the estimated remaining time in
-%   hours, minutes and seconds.
-%
-%   TIMERETA(TID,ITER), with no output arguments, prints the estimated
-%   remaining time to the screen. The time is displayed in the format
-%
-%     TIMERNAME: iteration ITER / ITERNUM, estimated remaining time: HH:MM:SS.SS
-%
-%   If the timer has no assigned name, the display format changes to
-%
-%     Iteration ITER / ITERNUM, estimated remaining time: HH:MM:SS.SS
-%
-%   TIMERETA(TID,ITER,DELAY) only displays the remaining time if the
-%   time elapsed since the previous printout is at least DELAY seconds.
-%
-%   See also TIMERINIT, TIMERCLEAR.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  June 2008
-
-
-global utiltbx_timer_start_times 
-global utiltbx_timer_iternums
-global utiltbx_timer_lastiter
-global utiltbx_time_lastdisp
-global utiltbx_timer_name
-
-
-if (tid<1 || tid>length(utiltbx_timer_iternums))
-  error('Unknown timer id');
-end
-
-if (utiltbx_timer_iternums(tid) < 0)
-  error('Specified timer does not have an associated number of iterations');
-end
-
-% update last reported iteration number
-utiltbx_timer_lastiter(tid) = iter;
-
-% compute elapsed time
-starttime = utiltbx_timer_start_times(tid,:);
-currtime = clock;
-timediff = etime(currtime, starttime);
-
-% total iteration number
-iternum = utiltbx_timer_iternums(tid);
-
-% compute eta
-timeremain = (iternum-iter)*timediff/iter;
-
-% return eta in seconds
-if (nargout==1)
-  varargout{1} = timeremain;
-  
-% return eta in hms
-elseif (nargout==3)
-  [varargout{1}, varargout{2}, varargout{3}] = secs2hms(timeremain);
-  
-  
-% print eta
-elseif (nargout==0)
-  
-  % check last display time
-  lastdisptime = utiltbx_time_lastdisp(tid,:);
-  if (nargin>2 && etime(currtime,lastdisptime) < delay)
-    return;
-  end
-  
-  % update last display time
-  utiltbx_time_lastdisp(tid,:) = currtime;
-
-  % display timer
-  [hrs,mins,secs] = secs2hms(timeremain);
-  if (isempty(utiltbx_timer_name{tid}))
-    printf('Iteration %d / %d, estimated remaining time: %02d:%02d:%05.2f', iter, iternum, hrs, mins, secs);
-  else
-    timername = utiltbx_timer_name{tid};
-    printf('%s: iteration %d / %d, estimated remaining time: %02d:%02d:%05.2f', timername, iter, iternum, hrs, mins, secs);
-  end
-   
-% invalid number of outputs
-else
-  error('Invalid number of output arguments');
-end
-
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/timerinit.m
--- a/Problems/private/timerinit.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,110 +0,0 @@
-function tid = timerinit(par1,par2)
-%TIMERINIT Initialize a new timer.
-%   TID = TIMERINIT() initializes a new timer for counting elapsed time,
-%   and returns its id.
-%
-%   TID = TIMERINIT('TIMERNAME') sets the timer name to the specified
-%   string for display purposes.
-%
-%   TID = TIMERINIT(ITERNUM) initializes a new ETA timer for a process with
-%   ITERNUM iterations. An ETA timer can be used for both counting elapsed
-%   time and estimating remaining time.
-%
-%   TID = TIMERINIT('TIMERNAME',ITERNUM) sets the ETA timer name to the
-%   specified string for display purposes.
-%
-%   Example:
-%
-%     tid = timerinit(100); 
-%     for i = 1:100
-%       pause(0.07);
-%       timereta(tid,i,1);
-%     end
-%     timereta(tid,i);
-%
-%   See also TIMERETA, TIMERCLEAR.
-
-
-%  Ron Rubinstein
-%  Computer Science Department
-%  Technion, Haifa 32000 Israel
-%  ronrubin@cs
-%
-%  June 2008
-
-
-global utiltbx_timer_start_times    % start times
-global utiltbx_time_lastdisp        % last display times
-global utiltbx_timer_iternums       % iteration numbers
-global utiltbx_timer_lastiter       % last queried iteration numbers
-global utiltbx_timer_name           % timer names
-global utiltbx_timer_callfun        % timer calling functions
-
-
-% parse function arguments %
-
-if (nargin==0)
-
-  iternum = -1;
-  timername = '';
-
-elseif (nargin==1)
-
-  if (ischar(par1))
-    iternum = -1;
-    timername = par1;
-
-  elseif (isnumeric(par1) && numel(par1)==1 && par1>0)
-    iternum = par1;
-    timername = '';
-
-  else
-    error('Invalid number of iterations');
-  end
-
-elseif (nargin==2)
-
-  if (ischar(par1) && isnumeric(par2))
-    if (numel(par2)==1 && par2>0)
-      timername = par1;
-      iternum = par2;
-    else
-      error('Invalid number of iterations');
-    end
-  else
-    error('Invalid function syntax');
-  end
-
-else
-  error('Too many function parameters');
-end
-
-
-% register the timer %
-
-if (isempty(utiltbx_timer_start_times))
-  utiltbx_timer_start_times = clock;
-  utiltbx_time_lastdisp = utiltbx_timer_start_times;
-  utiltbx_timer_iternums = double(iternum);
-  utiltbx_timer_lastiter = 0;
-  utiltbx_timer_name = { timername };
-  utiltbx_timer_callfun = {};
-  tid = 1;
-else
-  utiltbx_timer_start_times(end+1,:) = clock;
-  utiltbx_time_lastdisp(end+1,:) = utiltbx_timer_start_times(end,:);
-  utiltbx_timer_iternums(end+1) = double(iternum);
-  utiltbx_timer_lastiter(end+1) = 0;
-  utiltbx_timer_name{end+1} = timername;
-  tid = size(utiltbx_timer_start_times,1);
-end
-
-
-% detect timer calling function %
-
-st = dbstack;
-if (length(dbstack) >= 2)
-  utiltbx_timer_callfun{end+1} = st(2).name;
-else
-  utiltbx_timer_callfun{end+1} = '';
-end
diff -r f2ca04532eac -r 5f1f436057ca Problems/private/updateFigure.m
--- a/Problems/private/updateFigure.m	Wed Mar 16 14:13:48 2011 +0000
+++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
@@ -1,93 +0,0 @@
-function updateFigure(opts, figTitle, figFilename)
-
-%   Copyright 2008, Ewout van den Berg and Michael P. Friedlander
-%   http://www.cs.ubc.ca/labs/scl/sparco
-%   $Id: updateFigure.m 1040 2008-06-26 20:29:02Z ewout78 $
-
-% Ensure default values are available
-opts.linewidth  = getOption(opts,'linewidth', []);
-opts.fontsize   = getOption(opts,'fontsize',  []);
-opts.markersize = getOption(opts,'markersize',[]);
-
-% Output the plots
-if opts.update
-  % Set the line width, font size and marker size
-  chld = [gca; get(gca,'Children')];
-  lnwd = ones(length(chld),1) * NaN;
-  fnts = ones(length(chld),1) * NaN;
-  mrks = ones(length(chld),1) * NaN;
-  for i=1:length(chld)
-    conf = get(chld(i));
-    if ~isempty(opts.linewidth) && isfield(conf,'LineWidth')
-      lnwd(i) = get(chld(i),'LineWidth');
-      if (lnwd(i) == 0.5) % Default
-        set(chld(i),'Linewidth',opts.linewidth);
-      end
-    end
-    if ~isempty(opts.fontsize) && isfield(conf,'FontSize')
-      fnts(i) = get(chld(i),'FontSize');
-      if (fnts(i) == 10) % Default
-        set(chld(i),'FontSize',opts.fontsize);
-      end
-    end
-    if ~isempty(opts.markersize) && isfield(conf,'MarkerSize')
-      mrks(i) = get(chld(i),'MarkerSize');
-      if (mrks(i) == 6) % Default
-        set(chld(i),'MarkerSize',opts.markersize);
-      end
-    end
-  end
-  
-  for i=1:length(opts.figtype)
-     updateFigureType(opts.update, 0, opts.figtype{i}, ...
-                      opts.figpath, figTitle, figFilename);
-  end
-
-  % Restore the line-widths, font size
-  for i=1:length(chld)
-    if ~isnan(lnwd(i))
-      set(chld(i),'LineWidth',lnwd(i));
-    end
-    if ~isnan(fnts(i))
-      set(chld(i),'FontSize',fnts(i));
-    end
-    if ~isnan(mrks(i))
-      set(chld(i),'MarkerSize',mrks(i));
-    end
-  end
-    
-end
-
-% Show the plot
-if opts.show
-  updateFigureType(0,opts.show,'','',figTitle,'');
-end
-
-
-
-function updateFigureType(update,show,figtype,figpath,figTitle,figFilename)
-filename = [figpath,figFilename];
-
-switch lower(figtype)
- case {'pdf'}
-   cmdPostprocess = sprintf('!pdfcrop %s.pdf %s.pdf >& /dev/null', ...
-                            filename, filename);
- otherwise
-   cmdPostprocess = [];
-end
-
-[figtype,figext] = getFigureExt(figtype);
-
-% Print the figure for output (without title)
-if update
-  evalc(sprintf('print -d%s %s.%s;', figtype, filename, figext));
-  
-  if ~isempty(cmdPostprocess)
-    eval(cmdPostprocess);
-  end
-end
-
-% Add title if needed
-if show
-  title(figTitle);
-end