/* [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved. */
# include "meq.h"

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

 /* [FI,BRI,BZI,BRRI,BRZI,BZZI] = qint(RX,ZX,FX,RI,ZI,INM) */

# define MEXNAME qintmex

# define FI   pout[0]
# define BRI  pout[1]
# define BZI  pout[2]
# define BRRI pout[3]
# define BRZI pout[4]
# define BZRI pout[5]
# define BZZI pout[6]

# define RX   prhs[0]  /* double, def size=nr */
# define ZX   prhs[1]  /* double, def size=nz */
# define FX   prhs[2]  /* double, size=nz*nr */
# define RI   prhs[3]  /* double, def size=ni */
# define ZI   prhs[4]  /* double, size=ni */
# define INM  prhs[5]  /* double, def size=[n,*], size=[n,6] */

 int ni, nr, nz, n;
 int mri, nri;
  
 CHECK_NARGIN_EQ(6);

 ni = mxGetNumberOfElements(RI);
 nr = mxGetNumberOfElements(RX);
 nz = mxGetNumberOfElements(ZX);
 n  = mxGetM(INM);

 CHECK_DOUBLE(RX);
 CHECK_DOUBLE(ZX);
 CHECK_DOUBLE(FX);
 CHECK_NUMEL(FX,nr*nz);
 CHECK_DOUBLE(RI);
 CHECK_DOUBLE(ZI);
 CHECK_NUMEL(ZI,ni);
 CHECK_DOUBLE(INM);
 CHECK_NCOLS(INM,6);

 CHECK_NARGOUT_LE(7);

 mxArray *pout[7];
 
 mri = mxGetM(RI);
 nri = mxGetN(RI);
 
 FI   = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BRI  = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BZI  = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BRRI = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BRZI = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BZRI = mxCreateDoubleMatrix(mri,nri,mxREAL);
 BZZI = mxCreateDoubleMatrix(mri,nri,mxREAL);

 dqint(mxGetPr(FI), mxGetPr(BRI), mxGetPr(BZI), mxGetPr(BRRI), mxGetPr(BRZI), mxGetPr(BZRI), mxGetPr(BZZI),
       mxGetPr(RX), mxGetPr(ZX), mxGetPr(FX), mxGetPr(RI), mxGetPr(ZI), mxGetPr(INM),
       ni, nr, nz, n);


 ASSIGN_PLHS;
}
