/*----------------------------MegaWave2 module----------------------------*/
/* mwcommand
name = {classif_kppv};
version = {"1.0"};
author = {"no author"};
function = {"Calcule la carte de classification k-plus proches voisins"};
usage = {
 'n':[n=500]->n       "taille de l'image (carrée), par défaut 500",
 k->k                 "nombre de voisins (k-PPV)",
 e1->e1               "échantillons de la classe 1 (2-Flist)",
 e2->e2               "échantillons de la classe 2 (2-Flist)",
 out<-classif_kppv    "résultat: carte de classification (Cimage)"
        };
*/

#include <stdio.h>
#include <math.h>
#include "mw.h"


void insert(dtab,ctab,d,c,k)
     double *dtab,d;
     int *ctab,c,k;
{
  int i,j;

  for (i=k;i>=1 && (ctab[i-1]==0 || dtab[i-1]>d);i--);
  if (i<k) {
    for (j=k-1;j>i;j--) {
      dtab[j]=dtab[j-1];
      ctab[j]=ctab[j-1];
    }
    dtab[i]=d;
    ctab[i]=c;
  }
}

Cimage classif_kppv(e1,e2,k,n)
     Flist e1,e2;
     int k,*n;
{
  Cimage out;
  int x,y,i,count,*class;
  float xx,yy,*p;
  double d,*dist;

  if (k>e1->size+e2->size) 
    mwerror(FATAL,1,"number of neighbors is larger than data size.");

  out = mw_change_cimage(NULL,*n,*n);
  dist = (double *)malloc(k*sizeof(double));
  class = (int *)malloc(k*sizeof(int));

  /* parcours des pixels de la carte de classification (out) */
  for (x=0;x<*n;x++)
    for (y=0;y<*n;y++) {
      xx = (float)x/(float)(*n);
      yy = (float)y/(float)(*n);
      
      /* initialisation */
      for (i=0;i<k;i++) class[i]=0;
      
      /* parcours des échantillons (classe 1) */
      for (i=0,p=e1->values;i<e1->size;i++,p+=2) {
	d = (double)( (*p-xx)*(*p-xx) + (*(p+1)-yy)*(*(p+1)-yy) );
	insert(dist,class,d,1,k);
      }
      /* parcours des échantillons (classe 2) */
      for (i=0,p=e2->values;i<e2->size;i++,p+=2) {
	d = (double)( (*p-xx)*(*p-xx) + (*(p+1)-yy)*(*(p+1)-yy) );
	insert(dist,class,d,2,k);
      }

      /* remplissage de la carte de classification */
      for (i=0,count=0;i<k;i++) count += (class[i]==1?1:0);
      out->gray[*n*y+x] = (count>k-count?80:170);
    }

  free(class); free(dist);

  return(out);
}