/* [+GenLib General Purpose Library+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved. */
# include "genlib.h"

void mexFunction(int nlhs,       mxArray *plhs[],
                 int nrhs, const mxArray *prhs[])
{

# define MEXNAME ellipfastmex

 /* f = ellipfastmex(m,kind); */

# define F	pout[0]
	
# define M	prhs[0] /* double */
# define KIND	prhs[1] /* numeric, scalar */

 CHECK_NARGIN_EQ(2);

 CHECK_DOUBLE(M);
 CHECK_NUMERIC(KIND);
 CHECK_SCALAR(KIND);

 CHECK_NARGOUT_LE(1);

 mxArray *pout[1];

 int kind;
 register int k;
 register double *m, *f;

 F = mxCreateDoubleMatrix(mxGetM(M), mxGetN(M), mxREAL);
 kind = mxGetScalar(KIND);
 m = mxGetPr(M);
 f = mxGetPr(F);
 for (k = mxGetM(M)*mxGetN(M); k--;)
  *f++ = ellipfast(*m++,kind);

 ASSIGN_PLHS;
}
