/************************************************************************/
/*                                                                      */
/* LECTURE DES ENTREES                                                  */
/*                                                                      */
/* Fonction qui verifie la validite des arguments passes a la fonction  */
/* EMD et qui stocke les entrees dans des vecteurs C.                   */
/* Si le troisieme argument est superieur a 1, c'est le nombre d'IMFs a */
/* extraire du signal, s'il appartient a ]0,0.1], c'est le parametre    */
/* epsilon. De meme pour le quatrieme.                                  */
/*                                                                      */
/************************************************************************/

entree_t lire_entrees(int nlhs,int nrhs,const mxArray *prhs[]) {
  entree_t entree;
  int n,i,choix_eps=0;
  double *x,*y,*y_temp,troisieme,quatrieme;

  entree.choix=0;
  entree.eps=EPS_DEFAUT;

  /* verification des arguments */
  if (nrhs<2 || nrhs>4)
    mexErrMsgTxt("EMD requiert deux, trois ou quatre arguments en entree.");
  if (nlhs>1)
    mexErrMsgTxt("EMD requiert un seul argument en sortie.");
  if (!mxIsNumeric(prhs[0]) || mxIsComplex(prhs[0]) ||
        mxIsSparse(prhs[0]) || !mxIsDouble(prhs[0]) ||
      !mxIsNumeric(prhs[1]) || mxIsComplex(prhs[1]) ||
        mxIsSparse(prhs[1]) || !mxIsDouble(prhs[1]))
    mexErrMsgTxt("x et y doivent contenir des reels double precision.");

  /* lecture des entrees */
  n=PLUSGRAND(mxGetN(prhs[0]),mxGetM(prhs[0])); /* longueur du vecteur x */
  x=mxGetPr(prhs[0]);
  y_temp=mxGetPr(prhs[1]);

  /* troisieme arguments */
  if (nrhs==3 || nrhs==4) {
    if (!mxIsNumeric(prhs[2]) || mxIsComplex(prhs[2]) || mxIsSparse(prhs[2])
        || !mxIsDouble(prhs[2]) || mxGetN(prhs[2])!=1 || mxGetM(prhs[2])!=1)
      mexErrMsgTxt("epsilon et nb_imfs doivent etre des scalaires reels");
    troisieme=*mxGetPr(prhs[2]);
    if (troisieme<=0 || (troisieme>0.1 && troisieme<1))
      mexErrMsgTxt("le troisieme parametre a une valeur incorrecte");
    if (troisieme>=1) { /* le troisieme argument est nb_imfs */
      entree.choix=1;
      entree.nb_imfs=(int)troisieme;
    } else {            /* le troisieme argument est epsilon */
      choix_eps=1;
      entree.eps=troisieme;
    }
  }

  /* quatrieme argument */
  if (nrhs==4) {
    if (!mxIsNumeric(prhs[3]) || mxIsComplex(prhs[3]) || mxIsSparse(prhs[3])
        || !mxIsDouble(prhs[3]) || mxGetN(prhs[3])!=1 || mxGetM(prhs[3])!=1)
      mexErrMsgTxt("epsilon et nb_imfs doivent etre des scalaires reels");
    quatrieme=*mxGetPr(prhs[3]);
    if (quatrieme<=0 || (quatrieme>0.1 && quatrieme<1))
      mexErrMsgTxt("le quatrieme parametre a une valeur incorrecte");
    if (quatrieme>=1) { /* le quatrieme argument est nb_imfs */
      if (entree.choix)
        mexErrMsgTxt("le quatrieme argument doit etre inferieur a 0.1");
      entree.choix=1;
      entree.nb_imfs=(int)quatrieme;
    } else {            /* le quatrieme argument est epsilon */
      if (choix_eps)
        mexErrMsgTxt("le quatrieme argument doit etre superieur a 1");
      entree.eps=quatrieme;
    }
  }

  /* verification plus precise de x et de y */
  if (PLUSPETIT(mxGetN(prhs[0]),mxGetM(prhs[0]))!=1 ||
      PLUSPETIT(mxGetN(prhs[1]),mxGetM(prhs[1]))!=1)
    mexErrMsgTxt("x et y doivent etre des vecteurs");
  if (PLUSGRAND(mxGetN(prhs[1]),mxGetM(prhs[1]))!=n)
    mexErrMsgTxt("x et y doivent avoir la meme longueur.");
  i=1;
  while (i<n && x[i]>x[i-1]) i++;
  if (i<n) mexErrMsgTxt("Le vecteur x doit etre ordonne.");

  /* lecture et recopie de y (puisque l'on va modifier le */
  /* vecteur y, on ne travaille pas directement sur le    */
  /* vecteur Matlab passe en argument de la fonction EMD) */
  y=(double *)malloc(n*sizeof(double));
  for (i=0;i<n;i++) y[i]=y_temp[i];

  entree.n=n;
  entree.x=x;
  entree.y=y;
  return entree;
}


/************************************************************************/
/*                                                                      */
/* INITIALISATION DE LA LISTE                                           */
/*                                                                      */
/* Fonction qui initialise la liste chainee dans laquelle on stocke les */
/* imfs extraites du signal.                                            */
/*                                                                      */
/************************************************************************/

liste_t initialiser_liste(int n) {
  liste_t liste;
  liste.debut=NULL;
  liste.fin=NULL;
  liste.n=n;
  liste.m=0;
  return liste;
}


/************************************************************************/
/*                                                                      */
/* AJOUT D'UN MODE DANS LA LISTE                                        */
/*                                                                      */
/* Fonction qui ajoute une imf extraite du signal dans une liste        */
/* chainee.                                                             */
/*                                                                      */
/************************************************************************/

void ajouter(liste_t *liste,double *p) {
  double *v=(double *)malloc(liste->n*sizeof(double));
  int i;
  imf_t *mode=(imf_t *)malloc(sizeof(imf_t));
  for (i=0;i<liste->n;i++) v[i]=p[i];
  mode->pointeur=v;
  mode->suivant=NULL;
  if (!liste->debut) {
    liste->debut=mode;
  } else {
    (liste->fin)->suivant=mode;
  }
  liste->fin=mode;
  liste->m++;
}


/************************************************************************/
/*                                                                      */
/* LIBERATION DE LA LISTE                                               */
/*                                                                      */
/* Fonction qui libere la memoire allouee pour la liste chainee.        */
/*                                                                      */
/************************************************************************/

void liberer_liste(liste_t liste) {
  imf_t *courant=liste.debut, *precedent;  
  while (courant) {
    precedent=courant;
    courant=courant->suivant;
    free(precedent->pointeur);
    free(precedent);
  }
}


/************************************************************************/
/*                                                                      */
/* ECRITURE DES SORTIES                                                 */
/*                                                                      */
/* Fonction qui cree un matrice Matlab dans laquelle on recopie le      */
/* contenu de la liste chainee. Les matrices Matlab sont stockees       */
/* colonne par colonne, or on veut stocker la liste de telle sorte que  */
/* l'on ait une imf par ligne de la matrice Matlab. Voila pourquoi on   */
/* ne stocke pas la liste imf par imf.                                  */
/*                                                                      */
/************************************************************************/

void ecrire_sorties(liste_t liste,mxArray *plhs[]) {
  double *sortie;
  imf_t *parcours;
  int i=0,j,m=liste.m,n=liste.n;
  plhs[0]=mxCreateDoubleMatrix(m,n,mxREAL);
  sortie=mxGetPr(plhs[0]);
  for (parcours=liste.debut;parcours;parcours=parcours->suivant) {

    /* on stocke en memoire les elements de l'imf en    */
    /* les espacant de m cases memoires, ou "m" est le  */
    /* nombre d'imfs extraits du signal (residu inclus) */
    for (j=0;j<n;j++) *(sortie+j*m+i)=parcours->pointeur[j];

    i++;
  }
}
