#include "stdio.h"
#include "string.h"
#include "stdlib.h"
#include <math.h>

unsigned char quantify(unsigned char greyLevel, unsigned char *histo, int n);
void calcIntervalles( int *histo, unsigned char **bornes, int *n);
void calcLogScale(unsigned char **bornes, int *n);
void calcUnifScale(unsigned char **bornes, int *n);


void main(int argc, char **argv)
{
  char     nameIma[255], nameDim[255], nameIma2[255], nameDim2[255], *val;
  FILE     *image, *dim, *image2, *dim2, *fhisto;
  int      width, heigth, i, j, histo[256], max, N;
  unsigned char   **array, **array2, greyLevel, greyLevel2, *bornes;
  double    stdev;



  sprintf(nameIma, "%s.ima", argv[1]);
  sprintf(nameDim, "%s.dim", argv[1]);
  sprintf(nameIma2, "%sQ.ima", argv[1]);
  sprintf(nameDim2, "%sQ.dim", argv[1]);


  sprintf(nameIma, "%s.ima", argv[1]);
  sprintf(nameDim, "%s.dim", argv[1]);
  sprintf(nameIma2, "%sQ.ima", argv[1]);
  sprintf(nameDim2, "%sQ.dim", argv[1]);

  image = fopen(nameIma, "rb");
  dim = fopen(nameDim, "r");
  image2 = fopen(nameIma2, "wb");
  dim2 = fopen(nameDim2, "w");
  fhisto = fopen("histo", "w");
  fscanf(dim, "%d %d", &width, &heigth);
  fprintf(dim2, "%d %d", width, heigth);
  
  for (i=0; i<256; i++) histo[i] = 0;


  printf("Width = %d, Heigth = %d\n", width, heigth);

  array = (unsigned char **) malloc(width * sizeof(unsigned char *));
  array2 = (unsigned char **) malloc(width * sizeof(unsigned char *));

  for(i=0; i<width; i++)
  {
    array[i] = (unsigned char *) malloc(heigth * sizeof(unsigned char));
    array2[i] = (unsigned char *) malloc(heigth * sizeof(unsigned char));
    for(j=0; j<heigth; j++)
      { 
        fread(&array[i][j], sizeof(unsigned char), 1, image);
      	greyLevel = array[i][j];
        histo[greyLevel] ++;

      }
  }
  max = 0;
  for(i=0; i< 256; i++)
    if (histo[i] > max) max = histo[i];
 
  
  for (i=0; i<256; i++)
    {
      val = (char *) malloc(82 * sizeof(char));

      for (j=0; j<(histo[i]  * 80 / max); j++)
	{
	  val[j] = '*';
	}
      val[j] = '\0';
      fprintf(fhisto, "%s %d \n", val, i);
    }
        calcIntervalles(histo, &bornes, &N); 
      
  /*   calcUnifScale(&bornes, &N);  */
/*   calcLogScale(&bornes, &N);  */

  for(i=0; i<width; i++)
    {
    for(j=0; j<heigth; j++)
      {       
 	greyLevel2 = quantify(array[i][j], bornes, N);
	array2[i][j] = greyLevel2;
      	fwrite(&greyLevel2, sizeof(unsigned char), 1, image2);
      }
    } 
  stdev = 0;
  for(i=0;i<width; i++)
    {
    for(j=0; j<heigth; j++)
      {       
	stdev += (array[i][j] - array2[i][j]) * (array[i][j] - array2[i][j]);
      }
    }
  stdev /= width * heigth;
  stdev = sqrt(stdev);

  printf("Erreur de quantification moyenne =  %lf ( pourc/erreur = %lf)\n", stdev, 100 *stdev / 256);

  fclose(fhisto);
  fclose(image);
  fclose(dim);
  fclose(image2);
  fclose(dim2);

}


unsigned char quantify(unsigned char greyLevel, unsigned char *bornes, int n)
{
  unsigned char greyLevel2;
  int     i;

 
  i=1;
  while(bornes[i] < greyLevel)
    {
      i++;

    }
  greyLevel2 = (unsigned char) ((bornes[i-1]+bornes[i])/2);


  return(greyLevel2);
}
  
  
void calcIntervalles( int *histo, unsigned char **bornes, int *n)
{
  int         Nmediane = 8, i, j, *NparInter, *Nordre, temp;
  int         indice;

  printf("Methode de la mediane , nintervalles = %d\n", ((Nmediane)*(Nmediane + 1))/2);
  NparInter = (int *) malloc(Nmediane * sizeof(int));
  Nordre = (int *) malloc(Nmediane * sizeof(int));

  for (i=0; i<Nmediane; i++)
    {
      NparInter[i] = 0;
    }


  for (i=0; i<256; i++)
    {
    NparInter[i / Nmediane] += histo[i];
    }

   for (i=0; i<Nmediane; i++)
    {
      Nordre[i] = i;
    } 
 
   for (i=0; i<Nmediane; i++)
    {
      for (j=1; j <Nmediane - i; j++)
	{
	  if (NparInter[j] < NparInter[j-1]) 
	    {
	      temp = NparInter[j];
	      NparInter[j] = NparInter[j-1];
	      NparInter[j-1] = temp;
	      temp = Nordre[j];
	      Nordre[j] = Nordre[j-1];
	      Nordre[j-1] = temp;
	    }
	}
    } 

   for (i=0; i<Nmediane; i++)
     {
       printf("%d rang= %d\n", NparInter[i], Nordre[i]);
     }

 *bornes = (unsigned char *) malloc(( (Nmediane * (Nmediane +1)) / 2  + 1) * sizeof(unsigned char));

 indice = 0;


for (i=0; i<Nmediane; i++)
  {
    for (j=0; j<=Nordre[i]; j++)
      {
	(*bornes)[indice] =  (unsigned char) ( i * 255 / Nmediane + j * 255 / Nmediane / (Nordre[i]+1) );
	printf("%d\n", (*bornes)[indice]);
	indice ++;
      }
}

    (*bornes)[indice] = 255;
    indice ++;

    *n = indice;
}
  
  

 

 
void calcLogScale(unsigned char **bornes, int *n)
{
  int          i, LOG = 20;
  double   PARAM = 1;

 (*bornes) = (unsigned char *) malloc((LOG + 1) * sizeof(unsigned char));

  printf("Methode log , nintervalles = %d\n", LOG);


for (i=0; i<= LOG; i++)
  {
   (* bornes)[i] = (unsigned char) floor( (double) ( 255 * (exp (PARAM * i/(double)(LOG)) -1) 
				                 / (exp(PARAM )- 1)));
    printf ("bornes[%d]=%d\n", i , (*bornes)[i]);

  }
    *n = LOG + 1;
}
  
  
 
void calcUnifScale(unsigned char **bornes, int *n)
{
  int          i, NUNIF = 20;

 (*bornes) = (unsigned char *) malloc((NUNIF + 1) * sizeof(unsigned char));

  printf("Methode unif , nintervalles = %d\n", NUNIF);




for (i=0; i<= NUNIF; i++)
  {
   (* bornes)[i] = (unsigned char) floor( (double) ( 255 * i / ((double)(NUNIF))));
    printf ("bornes[%d]=%d\n", i , (*bornes)[i]);

  }
    *n = NUNIF + 1;
}