# HG changeset patch # User idamnjanovic # Date 1300285017 0 # Node ID c3eca463202dd7ca8cf7bd461c65539913962d81 # Parent 5f1f436057ca375ff53ecc7ebe79d6e5b1ad80f4 diff -r 5f1f436057ca -r c3eca463202d util/ksvd utils/add_dc.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/add_dc.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,33 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/addtocols.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/addtocols.c Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,85 @@ +/************************************************************************** + * + * 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 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] 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=m) { + mexErrMsgTxt("Row index in ROWS is out of range."); + } + } + + /* Do the actual computation */ + for (i=0; i=n) { + mexErrMsgTxt("Column index in COLS is out of range."); + } + for (j=0; j=n) { + mexErrMsgTxt("Column index in COLS is out of range."); + } + for (j=0; j + + +/* 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]=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 5f1f436057ca -r c3eca463202d util/ksvd utils/imnormalize.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/imnormalize.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,19 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/iswhole.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/iswhole.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,17 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/make.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/make.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,40 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/mexutils.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/mexutils.c Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,79 @@ +/************************************************************************** + * + * File name: mexutils.c + * + * Ron Rubinstein + * Computer Science Department + * Technion, Haifa 32000 Israel + * ronrubin@cs + * + * Last Updated: 15.8.2009 + * + *************************************************************************/ + +#include "mexutils.h" +#include + + + +/* 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 5f1f436057ca -r c3eca463202d util/ksvd utils/mexutils.h --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/mexutils.h Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,103 @@ +/************************************************************************** + * + * 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 5f1f436057ca -r c3eca463202d util/ksvd utils/normcols.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/normcols.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,17 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/printf.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/printf.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,26 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/reggrid.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/reggrid.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,136 @@ +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(active_dims(id))) + active_dims = find(n1); + 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(active_dims(id))) + active_dims = find(n 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=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=n) { + mexErrMsgTxt("Column index in COLS is out of range."); + } + } + + /* Do the actual computation */ + for (j=0; j 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 5f1f436057ca -r c3eca463202d util/ksvd utils/sprow.c --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/sprow.c Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,94 @@ +/************************************************************************** + * + * 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 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 5f1f436057ca -r c3eca463202d util/ksvd utils/timereta.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/timereta.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,98 @@ +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 5f1f436057ca -r c3eca463202d util/ksvd utils/timerinit.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/util/ksvd utils/timerinit.m Wed Mar 16 14:16:57 2011 +0000 @@ -0,0 +1,110 @@ +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