Index: /issm/trunk/src/c/Makefile.am
===================================================================
--- /issm/trunk/src/c/Makefile.am	(revision 7859)
+++ /issm/trunk/src/c/Makefile.am	(revision 7860)
@@ -392,4 +392,5 @@
 					./io/FetchNodeSets.cpp\
 					./io/OptionParse.cpp\
+					./io/MatlabNArrayToNArray.cpp\
 					./io/pfopen.cpp\
 					./io/pfclose.cpp\
Index: /issm/trunk/src/c/io/FetchData.cpp
===================================================================
--- /issm/trunk/src/c/io/FetchData.cpp	(revision 7859)
+++ /issm/trunk/src/c/io/FetchData.cpp	(revision 7860)
@@ -87,4 +87,47 @@
 }
 /*}}}*/
+/*FUNCTION FetchData(double** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){{{1*/
+void FetchData(double** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){
+
+	double*  outmatrix=NULL;
+	int      outmatrix_numel,outmatrix_ndims;
+	int*     outmatrix_size=NULL;
+
+	if(mxIsEmpty(dataref) ){
+		/*Nothing to pick up. Just initialize matrix pointer to NULL: */
+		outmatrix_numel=0;
+		outmatrix_ndims=0;
+		outmatrix_size =NULL;
+		outmatrix=NULL;
+	}
+	else if (mxIsClass(dataref,"double") ){
+
+		/*Check dataref is not pointing to NaN: */
+		if ( mxIsNaN(*(mxGetPr(dataref))) && (mxGetNumberOfElements(dataref)==1) ){
+			outmatrix_numel=0;
+			outmatrix_ndims=0;
+			outmatrix_size =NULL;
+			outmatrix=NULL;
+		}
+		else{
+
+			/*Convert matlab n-dim array to double* matrix: */
+			MatlabNArrayToNArray(&outmatrix,&outmatrix_numel,&outmatrix_ndims,&outmatrix_size,dataref);
+		}
+	}
+	else{
+		/*This is an error: we don't have the correct input!: */
+		_error_("Input parameter of class %s not supported yet",mxGetClassName(dataref));
+	}
+			
+	/*Assign output pointers:*/
+	*pmatrix=outmatrix;
+	if (pnumel)*pnumel=outmatrix_numel;
+	if (pndims)*pndims=outmatrix_ndims;
+	if (psize )*psize =outmatrix_size;
+	else xfree((void**)&outmatrix_size);
+
+}
+/*}}}*/
 /*FUNCTION FetchData(int** pmatrix,int* pM,int *pN,const mxArray* dataref){{{1*/
 void FetchData(int** pmatrix,int* pM,int *pN,const mxArray* dataref){
@@ -169,4 +212,69 @@
 	if (pM)*pM=outmatrix_rows;
 	if (pN)*pN=outmatrix_cols;
+}
+/*}}}*/
+/*FUNCTION FetchData(bool** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){{{1*/
+void FetchData(bool** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){
+
+	int      i;
+	int      outmatrix_numel,outmatrix_ndims;
+	int*     outmatrix_size=NULL;
+	double*  doublematrix=NULL;
+	bool*    outmatrix=NULL;
+
+	if(mxIsEmpty(dataref) ){
+		/*Nothing to pick up. Just initialize matrix pointer to NULL: */
+		outmatrix_numel=0;
+		outmatrix_ndims=0;
+		outmatrix_size =NULL;
+		outmatrix=NULL;
+	}
+	else if (mxIsClass(dataref,"logical") ){
+
+		/*Check dataref is not pointing to NaN: */
+		if ( mxIsNaN(*((bool*)mxGetData(dataref))) && (mxGetNumberOfElements(dataref)==1) ){
+			outmatrix_numel=0;
+			outmatrix_ndims=0;
+			outmatrix_size =NULL;
+			outmatrix=NULL;
+		}
+		else{
+
+			/*Convert matlab n-dim array to bool* matrix: */
+			MatlabNArrayToNArray(&outmatrix,&outmatrix_numel,&outmatrix_ndims,&outmatrix_size,dataref);
+		}
+	}
+	else if (mxIsClass(dataref,"double") ){
+
+		/*Check dataref is not pointing to NaN: */
+		if ( mxIsNaN(*(mxGetPr(dataref))) && (mxGetNumberOfElements(dataref)==1) ){
+			outmatrix_numel=0;
+			outmatrix_ndims=0;
+			outmatrix_size =NULL;
+			outmatrix=NULL;
+		}
+		else{
+
+			/*Convert matlab n-dim array to double* matrix: */
+			MatlabNArrayToNArray(&doublematrix,&outmatrix_numel,&outmatrix_ndims,&outmatrix_size,dataref);
+
+			/*Convert double matrix into bool matrix: */
+			outmatrix=(bool*)xmalloc(outmatrix_numel*sizeof(bool));
+			for(i=0;i<outmatrix_numel;i++)outmatrix[i]=(bool)doublematrix[i];
+			xfree((void**)&doublematrix);
+		}
+	}
+	else{
+		/*This is an error: we don't have the correct input!: */
+		_error_("Input parameter of class %s not supported yet",mxGetClassName(dataref));
+	}
+			
+	/*Assign output pointers:*/
+	*pmatrix=outmatrix;
+	if (pnumel)*pnumel=outmatrix_numel;
+	if (pndims)*pndims=outmatrix_ndims;
+	if (psize )*psize =outmatrix_size;
+	else xfree((void**)&outmatrix_size);
+
 }
 /*}}}*/
@@ -366,4 +474,46 @@
 	/*Assign output pointers:*/
 	*pstring=outstring;
+}
+/*FUNCTION FetchData(char** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){{{1*/
+void FetchData(char** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref){
+
+	int      outmatrix_numel,outmatrix_ndims;
+	int*     outmatrix_size=NULL;
+	char*    outmatrix=NULL;
+
+	if(mxIsEmpty(dataref) ){
+		/*Nothing to pick up. Just initialize matrix pointer to NULL: */
+		outmatrix_numel=0;
+		outmatrix_ndims=0;
+		outmatrix_size =NULL;
+		outmatrix=NULL;
+	}
+	else if (mxIsClass(dataref,"char") ){
+
+		/*Check dataref is not pointing to NaN: */
+		if ( mxIsNaN(*(mxGetPr(dataref))) && (mxGetNumberOfElements(dataref)==1) ){
+			outmatrix_numel=0;
+			outmatrix_ndims=0;
+			outmatrix_size =NULL;
+			outmatrix=NULL;
+		}
+		else{
+
+			/*Convert matlab n-dim array to char* matrix: */
+			MatlabNArrayToNArray(&outmatrix,&outmatrix_numel,&outmatrix_ndims,&outmatrix_size,dataref);
+		}
+	}
+	else{
+		/*This is an error: we don't have the correct input!: */
+		_error_("Input parameter of class %s not supported yet",mxGetClassName(dataref));
+	}
+			
+	/*Assign output pointers:*/
+	*pmatrix=outmatrix;
+	if (pnumel)*pnumel=outmatrix_numel;
+	if (pndims)*pndims=outmatrix_ndims;
+	if (psize )*psize =outmatrix_size;
+	else xfree((void**)&outmatrix_size);
+
 }
 /*}}}*/
Index: /issm/trunk/src/c/io/MatlabNArrayToNArray.cpp
===================================================================
--- /issm/trunk/src/c/io/MatlabNArrayToNArray.cpp	(revision 7860)
+++ /issm/trunk/src/c/io/MatlabNArrayToNArray.cpp	(revision 7860)
@@ -0,0 +1,272 @@
+/* \file MatlabConversions.cpp
+ * \brief: convert a sparse or dense matlab array to cpp array
+ */
+
+
+#ifdef HAVE_CONFIG_H
+	#include "config.h"
+#else
+#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
+#endif
+
+#include "./io.h"
+#include "../shared/shared.h"
+#include "../include/include.h"
+
+#ifdef _SERIAL_
+#include <mex.h>
+
+int MatlabNArrayToNArray(double** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix){
+
+	int  i,j,rows,cols;
+	int  numel,ndims;
+	int *size,*dims;
+	double* mxmatrix_ptr=NULL;
+	const mwSize* ipt=NULL;
+
+	/*output: */
+	double* matrix=NULL;
+
+	/*matlab indices: */
+	mwIndex*    ir=NULL;
+	mwIndex*    jc=NULL;
+	double* pr=NULL;
+	int     count;
+	int     nnz;
+	int     nz;
+
+	/*get Matlab matrix information: */
+	numel=mxGetNumberOfElements(mxmatrix);
+	ndims=mxGetNumberOfDimensions(mxmatrix);
+	ipt  =mxGetDimensions(mxmatrix);
+	size =(int *) xcalloc(ndims,sizeof(int));
+	for (i=0; i<ndims; i++) size[i]=(int)ipt[i];
+
+	/*Ok, first check if we are dealing with a sparse or full matrix: */
+	if (mxIsSparse(mxmatrix)){
+
+		/*Dealing with sparse matrix: recover size first: */
+		rows=mxGetM(mxmatrix);
+		cols=mxGetN(mxmatrix);
+		nnz=mxGetNzmax(mxmatrix);
+		nz=(int)((double)nnz/(double)rows);
+
+		matrix=(double*)xcalloc(rows*cols,sizeof(double));
+
+		/*Now, get ir,jc and pr: */
+		ir=mxGetIr(mxmatrix);
+		jc=mxGetJc(mxmatrix);
+		pr=mxGetPr(mxmatrix);
+
+		/*Now, start inserting data into double* matrix: */
+		count=0;
+		for(i=0;i<cols;i++){
+			for(j=0;j<(jc[i+1]-jc[i]);j++){
+				*(matrix+rows*ir[count]+i)=pr[count];
+				count++;
+			}
+		}
+
+	}
+	else{
+
+		/*Dealing with dense matrix: recover pointer and size: */
+		mxmatrix_ptr=(double*)mxGetPr(mxmatrix);
+		
+		/*Create serial matrix: */
+		matrix=(double*)xcalloc(numel,sizeof(double));
+
+		dims=(int *) xcalloc(ndims,sizeof(int));
+		for(i=0;i<numel;i++){
+			ColumnWiseDimsFromIndex(dims,i,size,ndims);
+			j=IndexFromRowWiseDims(dims,size,ndims);
+			*(matrix+j)=*(mxmatrix_ptr+i);
+		}
+		xfree((void**)&dims);
+		
+	}
+
+	/*Assign output pointer: */
+	*pmatrix=matrix;
+	*pmatrix_numel=numel;
+	*pmatrix_ndims=ndims;
+	*pmatrix_size=size;
+
+	return 1;
+}
+
+int MatlabNArrayToNArray(bool** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix){
+
+	int  i,j,rows,cols;
+	int  numel,ndims;
+	int *size,*dims;
+	bool* mxmatrix_ptr=NULL;
+	const mwSize* ipt=NULL;
+
+	/*output: */
+	bool* matrix=NULL;
+
+	/*matlab indices: */
+	mwIndex*    ir=NULL;
+	mwIndex*    jc=NULL;
+	bool*   pm=NULL;
+	int     count;
+	int     nnz;
+	int     nz;
+
+	/*get Matlab matrix information: */
+	numel=mxGetNumberOfElements(mxmatrix);
+	ndims=mxGetNumberOfDimensions(mxmatrix);
+	ipt  =mxGetDimensions(mxmatrix);
+	size =(int *) xcalloc(ndims,sizeof(int));
+	for (i=0; i<ndims; i++) size[i]=(int)ipt[i];
+
+	/*Ok, first check if we are dealing with a sparse or full matrix: */
+	if (mxIsSparse(mxmatrix)){
+
+		/*Dealing with sparse matrix: recover size first: */
+		rows=mxGetM(mxmatrix);
+		cols=mxGetN(mxmatrix);
+		nnz=mxGetNzmax(mxmatrix);
+		nz=(int)((double)nnz/(double)rows);
+
+		matrix=(bool*)xcalloc(rows*cols,sizeof(bool));
+
+		/*Now, get ir,jc and pm: */
+		ir=mxGetIr(mxmatrix);
+		jc=mxGetJc(mxmatrix);
+		pm=(bool*)mxGetData(mxmatrix);
+
+		/*Now, start inserting data into bool* matrix: */
+		count=0;
+		for(i=0;i<cols;i++){
+			for(j=0;j<(jc[i+1]-jc[i]);j++){
+				*(matrix+rows*ir[count]+i)=pm[count];
+				count++;
+			}
+		}
+
+	}
+	else{
+
+		/*Dealing with dense matrix: recover pointer and size: */
+		mxmatrix_ptr=(bool*)mxGetData(mxmatrix);
+		
+		/*Create serial matrix: */
+		matrix=(bool*)xcalloc(numel,sizeof(bool));
+
+		dims=(int *) xcalloc(ndims,sizeof(int));
+		for(i=0;i<numel;i++){
+			ColumnWiseDimsFromIndex(dims,i,size,ndims);
+			j=IndexFromRowWiseDims(dims,size,ndims);
+			*(matrix+j)=(bool)*(mxmatrix_ptr+i);
+		}
+		xfree((void**)&dims);
+		
+	}
+
+	/*Assign output pointer: */
+	*pmatrix=matrix;
+	*pmatrix_numel=numel;
+	*pmatrix_ndims=ndims;
+	*pmatrix_size=size;
+
+	return 1;
+}
+
+int MatlabNArrayToNArray(char** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix){
+
+	int  i,j,rows,cols;
+	int  numel,ndims;
+	int *size,*dims;
+	mxChar* mxmatrix_ptr=NULL;
+	const mwSize* ipt=NULL;
+
+	/*output: */
+	char* matrix=NULL;
+
+	/*matlab indices: */
+	mwIndex*    ir=NULL;
+	mwIndex*    jc=NULL;
+	char*   pm=NULL;
+	int     count;
+	int     nnz;
+	int     nz;
+
+	/*get Matlab matrix information: */
+	numel=mxGetNumberOfElements(mxmatrix);
+	ndims=mxGetNumberOfDimensions(mxmatrix);
+	ipt  =mxGetDimensions(mxmatrix);
+	size =(int *) xcalloc(ndims,sizeof(int));
+	for (i=0; i<ndims; i++) size[i]=(int)ipt[i];
+
+	/*Ok, first check if we are dealing with a sparse or full matrix: */
+	if (mxIsSparse(mxmatrix)){
+
+		/*Dealing with sparse matrix: recover size first: */
+		rows=mxGetM(mxmatrix);
+		cols=mxGetN(mxmatrix);
+		nnz=mxGetNzmax(mxmatrix);
+		nz=(int)((double)nnz/(double)rows);
+
+		matrix=(char*)xcalloc(rows*cols,sizeof(double));
+
+		/*Now, get ir,jc and pm: */
+		ir=mxGetIr(mxmatrix);
+		jc=mxGetJc(mxmatrix);
+		pm=(char*)mxGetData(mxmatrix);
+
+		/*Now, start inserting data into char* matrix: */
+		count=0;
+		for(i=0;i<cols;i++){
+			for(j=0;j<(jc[i+1]-jc[i]);j++){
+				*(matrix+rows*ir[count]+i)=(char)pm[count];
+				count++;
+			}
+		}
+
+	}
+	else{
+
+		/*Dealing with dense matrix: recover pointer and size: */
+		mxmatrix_ptr=mxGetChars(mxmatrix);
+		
+		/*Create serial matrix: */
+		matrix=(char*)xcalloc(numel+1,sizeof(mxChar));
+
+		/*looping code adapted from Matlab example explore.c: */
+		int elements_per_page = size[0] * size[1];
+		/* total_number_of_pages = size[2] x size[3] x ... x size[N-1] */
+		int total_number_of_pages = 1;
+		for (i=2; i<ndims; i++) {
+			total_number_of_pages *= size[i];
+		}
+
+		i=0;
+		for (int page=0; page < total_number_of_pages; page++) {
+			int row;
+			/* On each page, walk through each row. */
+			for (row=0; row<size[0]; row++)  {
+				int column;
+				j = (page * elements_per_page) + row;
+
+				/* Walk along each column in the current row. */
+				for (column=0; column<size[1]; column++) {
+					*(matrix+i++)=(char)*(mxmatrix_ptr+j);
+					j += size[0];
+				}
+			}
+		}
+
+	}
+
+	/*Assign output pointer: */
+	*pmatrix=matrix;
+	*pmatrix_numel=numel;
+	*pmatrix_ndims=ndims;
+	*pmatrix_size=size;
+
+	return 1;
+}
+
+#endif
Index: /issm/trunk/src/c/io/OptionParse.cpp
===================================================================
--- /issm/trunk/src/c/io/OptionParse.cpp	(revision 7859)
+++ /issm/trunk/src/c/io/OptionParse.cpp	(revision 7860)
@@ -42,5 +42,4 @@
 
 	OptionDouble *odouble = NULL;
-	const mwSize  *ipt     = NULL;
 
 	/*check and parse the name  */
@@ -54,16 +53,5 @@
 	}
 
-	odouble->numel=mxGetNumberOfElements(prhs[0]);
-	odouble->ndims=mxGetNumberOfDimensions(prhs[0]);
-	ipt           =mxGetDimensions(prhs[0]);
-	odouble->size =(int *) xmalloc(odouble->ndims*sizeof(int));
-	for (int i=0; i<odouble->ndims; i++) odouble->size[i]=(int)ipt[i];
-
-	//  note that FetchData does not correctly handle ndims >= 3
-	if (odouble->ndims > 2) {
-		_printf_(true,"WARNING -- option \"%s\" of class \"%s\" has ndims=%d and will be skipped.\n",odouble->name,mxGetClassName(prhs[0]),odouble->ndims);
-		odouble->numel=0;
-	}
-	else FetchData(&odouble->values,NULL,NULL,prhs[0]);
+	FetchData(&odouble->values,&odouble->numel,&odouble->ndims,&odouble->size,prhs[0]);
 
 	return(odouble);
@@ -73,6 +61,4 @@
 
 	OptionLogical *ological = NULL;
-	const mwSize   *ipt      = NULL;
-	bool            btemp;
 
 	/*check and parse the name  */
@@ -86,21 +72,5 @@
 	}
 
-	ological->numel=mxGetNumberOfElements(prhs[0]);
-	ological->ndims=mxGetNumberOfDimensions(prhs[0]);
-	ipt            =mxGetDimensions(prhs[0]);
-	ological->size =(int *) xmalloc(ological->ndims*sizeof(int));
-	for (int i=0; i<ological->ndims; i++) ological->size[i]=(int)ipt[i];
-
-	//  note that FetchData does not correctly handle non-scalar logicals
-	if (ological->ndims > 2 || ological->size[0] > 1 || ological->size[1] > 1) {
-		_printf_(true,"WARNING -- option \"%s\" of class \"%s\" is more than [1x1] and will be skipped.\n",ological->name,mxGetClassName(prhs[0]));
-		ological->numel=0;
-	}
-	else {
-		//FetchData(&ological->values,prhs[0]);
-		//could be memory leak until FetchData handles logical arrays
-		ological->values=(bool *) xmalloc(sizeof(bool));
-		FetchData(ological->values,prhs[0]);
-	}
+	FetchData(&ological->values,&ological->numel,&ological->ndims,&ological->size,prhs[0]);
 
 	return(ological);
@@ -110,5 +80,4 @@
 
 	OptionChar  *ochar = NULL;
-	const mwSize *ipt   = NULL;
 
 	/*check and parse the name  */
@@ -122,16 +91,5 @@
 	}
 
-	ochar->numel=mxGetNumberOfElements(prhs[0]);
-	ochar->ndims=mxGetNumberOfDimensions(prhs[0]);
-	ipt         =mxGetDimensions(prhs[0]);
-	ochar->size =(int *) xmalloc(ochar->ndims*sizeof(int));
-	for(int i=0; i<ochar->ndims; i++) ochar->size[i]=(int)ipt[i];
-
-	//note that FetchData does not correctly handle ndims >= 2 or multiple rows
-	if (ochar->ndims > 2 || ochar->size[0] > 1) {
-		_printf_(true,"WARNING -- option \"%s\" of class \"%s\" is more than [1xn] and will be skipped.\n",ochar->name,mxGetClassName(prhs[0]));
-		ochar->numel=0;
-	}
-	else FetchData(&ochar->values,prhs[0]);
+	FetchData(&ochar->values,&ochar->numel,&ochar->ndims,&ochar->size,prhs[0]);
 
 	return(ochar);
Index: /issm/trunk/src/c/io/io.h
===================================================================
--- /issm/trunk/src/c/io/io.h	(revision 7859)
+++ /issm/trunk/src/c/io/io.h	(revision 7860)
@@ -33,6 +33,8 @@
 void FetchData(DataSet** pdataset,const mxArray* dataref);
 void FetchData(double** pmatrix,int* pM,int *pN,const mxArray* dataref);
+void FetchData(double** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref);
 void FetchData(int** pmatrix,int* pM,int *pN,const mxArray* dataref);
 void FetchData(bool** pmatrix,int* pM,int *pN,const mxArray* dataref);
+void FetchData(bool** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref);
 void FetchData(Mat* pmatrix,const mxArray* dataref);
 void FetchData(int** pvector,int* pM,const mxArray* dataref);
@@ -42,4 +44,5 @@
 void FetchData(Vec* pvector,const mxArray* dataref);
 void FetchData(char** pstring,const mxArray* dataref);
+void FetchData(char** pmatrix,int* pnumel,int* pndims,int** psize,const mxArray* dataref);
 void FetchData(double* pscalar,const mxArray* dataref);
 void FetchData(int* pinteger,const mxArray* dataref);
@@ -64,4 +67,8 @@
 OptionStruct*   OptionStructParse( char* name, const mxArray* prhs[]);
 OptionCell*     OptionCellParse( char* name, const mxArray* prhs[]);
+
+int MatlabNArrayToNArray(double** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix);
+int MatlabNArrayToNArray(bool** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix);
+int MatlabNArrayToNArray(char** pmatrix,int* pmatrix_numel,int* pmatrix_ndims,int** pmatrix_size,const mxArray* mxmatrix);
 
 #endif