/* Programme Pyramides Morphologiques */

#include <stdio.h>
#include <stdlib.h>
#include "cimage.h"

#define ENTIERE 0
#define QUINCONCE 1
#define NMAX 120

/* type Image : comporte un champ echelle
                         (seuls les points de coordonnees multiples de echelle appartiennent a l'image)
                          un champ type_grille
                          (QUINCONCE ou ENTIERE)
*/

typedef struct {
   int nrow;
   int ncol;
   int nb_coul;
   int type_grille;
   int echelle;
   unsigned char* gray;
} Cimageprov;

/*****************************************/
/* Fonction de lecture d'une image de type PGM */
/*****************************************/

FILE* lire_image(char* name, int* width, int* height, int* nb_coul,
               Cimageprov* image, int* pos) {
long nb_lu;
char ce_char;
FILE* fid;
char   start_of_line;
char   magic [3];
int    in_header = 1;

fid = fopen (name, "r");

   if (((int) fid) == -1) {
      printf ("Pb fopen\n");
      return (0);
   }

nb_lu = fread (magic, 1, 3, fid);

while (in_header) {
nb_lu = fread (&start_of_line, 1, 1, fid);

      if (start_of_line != '#') {
         printf ("in header no more #\n");
         fseek (fid, -1, SEEK_CUR);
         break;
      }
      else {
         printf ("in header lire comment\n");
         nb_lu = fread (&ce_char, 1, 1, fid);
         while ((ce_char != 10) && (ce_char !=13))
            nb_lu = fread (&ce_char, 1, 1, fid);
         nb_lu = fread (&ce_char, 1, 1, fid);
         if ((ce_char != 10) && (ce_char !=13))
            fseek (fid, -1, SEEK_CUR);
         printf ("in header fin lire comment\n");
      }
   }

printf ("FIN comment\n");

   fscanf (fid, "%d", width);
   printf ("Lecture width ok\n");
   fscanf (fid, "%d", height);
   printf ("Lecture h ok\n");
   fscanf (fid, "%d", nb_coul);
   printf ("Lecture coul ok\n");

printf ("%d %d %d \n", *height, *width, *nb_coul);

fread (&ce_char, 1, 1, fid);
*pos = ftell (fid);

image->gray = malloc ((*height) * (*width));
nb_lu = fread (image->gray, *width, *height, fid);

   image->nrow = *height;
   image->ncol = *width;
   image->nb_coul = *nb_coul;

return (fid);
}

/************************************************/
/* Fonction de copie d'une variable image dans une autre */
/* (formatte de plus une image type QUINCONCE pour */
/* affichage en mettant a blanc les points hors treillis) */
/************************************************/

void sauver_image_niveau_n (Cimageprov* image, Cimageprov* image_niveau_n) {
int i,j;
   image_niveau_n->ncol = image->ncol/image->echelle;
   image_niveau_n->nrow = image->nrow/image->echelle;
   image_niveau_n->gray = malloc (image_niveau_n->ncol * image_niveau_n->nrow);

   if (image->type_grille == QUINCONCE)
      for (i = 0; i < image->ncol; i+=image->echelle)
         for (j = 0; j < image->nrow; j+=image->echelle)
            if (((i+j) / image->echelle) % 2 != 0)
              image->gray [i+j*image->ncol] = 255;

   for (i = 0; i < image->ncol/image->echelle; i++) {
      for (j = 0; j < image->nrow/image->echelle; j++) {
         image_niveau_n->gray [i+ j*image_niveau_n->ncol] =
              image->gray [image->echelle*i+ j*image->ncol*image->echelle];
         }
      }
}

/**************************************************/
/* Passage de Xn a Xn+1 */
/* le resultat est mis dans la meme variable image */
/**************************************************/

void calcul_niveau_suivant (Cimageprov* image) {
long k;
int i,j;
   if (image->type_grille == ENTIERE) {
      for (i = 0; i < image->ncol; i+=image->echelle) {
         for (j = 0; j < image->nrow; j+=image->echelle) {
            if (((i+j) / image->echelle) % 2 == 0) {
               k = i+j*image->ncol;
               if ((i>=image->echelle) && (image->gray [k-image->echelle]<image->gray [k]))
                  image->gray [k] = image->gray [k-image->echelle];
               if ((i<image->ncol-image->echelle) && (image->gray [k+image->echelle]<image->gray [k]))
                  image->gray [k] = image->gray [k+image->echelle];
               if ((j>= image->echelle) && (image->gray [k-image->echelle*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [k-image->echelle*image->ncol];
               if ((j<image->nrow-image->echelle) && (image->gray [k+image->echelle*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [k+image->echelle*image->ncol];
            }
         }
      }
      image->type_grille = QUINCONCE;
   } else {
      for (i = 0; i < image->ncol; i+=2*image->echelle) {
         for (j = 0; j < image->nrow; j+=2*image->echelle) {
               k = i+j*image->ncol;
               if ((i>=image->echelle) && (j>=image->echelle) &&
                    (image->gray [i-image->echelle+ (j-image->echelle)*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [i-image->echelle+ (j-image->echelle)*image->ncol];
               if ((i>=image->echelle) && (j<image->nrow-image->echelle) &&
                    (image->gray [i-image->echelle+ (j+image->echelle)*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [i-image->echelle+ (j+image->echelle)*image->ncol];
               if ((i<image->ncol-image->echelle) && (j>=image->echelle) &&
                   (image->gray [i+image->echelle+ (j-image->echelle)*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [i+image->echelle+ (j-image->echelle)*image->ncol];
               if ((i<image->ncol-image->echelle) && (j<image->nrow-image->echelle) &&
                   (image->gray [i+image->echelle+ (j+image->echelle)*image->ncol]<image->gray [k]))
                  image->gray [k] = image->gray [i+image->echelle+ (j+image->echelle)*image->ncol];
         }
      }
      image->echelle = 2*image->echelle;
      image->type_grille = ENTIERE;
   }
}

/*************************************************/
/* Operateur de synthese : passage Xn a Xn-1 */
/*************************************************/

void calcul_niveau_prec (Cimageprov* image) {
int max,i,j;
long k;
if (image->type_grille == QUINCONCE) {

      for (i = 0; i < image->ncol; i+=image->echelle) {
         for (j = 0; j < image->nrow; j+=image->echelle) {
            if (((i+j)/image->echelle) % 2 != 0) {
               max = 0;
               k = i+j*image->ncol;
               if ((i>=image->echelle) && (image->gray [k-image->echelle]>max))
                  max = image->gray [k-image->echelle];
               if ((i<image->ncol-image->echelle) && (image->gray [k+image->echelle]>max))
                  max = image->gray [k+image->echelle];
               if ((j>= image->echelle) && (image->gray [k-image->echelle*image->ncol]>max))
                  max = image->gray [k-image->echelle*image->ncol];
               if ((j<image->nrow-image->echelle) && (image->gray [k+image->echelle*image->ncol]>max))
                  max = image->gray [k+image->echelle*image->ncol];
               image->gray [k] = max;
            }
         }
      }

image->type_grille = ENTIERE;
} else {

image->echelle = image->echelle/2;

      for (i = 0; i < image->ncol; i+=image->echelle) {
         for (j = 0; j < image->nrow; j+=image->echelle) {
             if ((((i/image->echelle) % 2) != 0) && (((j/image->echelle) % 2) != 0)) {
                  max = 0;
                  k = i+j*image->ncol;
                  if ((i>=image->echelle) && (j>=image->echelle) &&
                       (image->gray [i-image->echelle+ (j-image->echelle)*image->ncol]>max))
                     max = image->gray [i-image->echelle+ (j-image->echelle)*image->ncol];
                  if ((i>=image->echelle) && (j<image->nrow-image->echelle) &&
                       (image->gray [i-image->echelle+ (j+image->echelle)*image->ncol]>max))
                     max = image->gray [i-image->echelle+ (j+image->echelle)*image->ncol];
                  if ((i<image->ncol-image->echelle) && (j>=image->echelle) &&
                       (image->gray [i+image->echelle+ (j-image->echelle)*image->ncol]>max))
                     max = image->gray [i+image->echelle+ (j-image->echelle)*image->ncol];
                  if ((i<image->ncol-image->echelle) && (j<image->nrow-image->echelle) &&
                      (image->gray [i+image->echelle+ (j+image->echelle)*image->ncol]>max))
                     max = image->gray [i+image->echelle+ (j+image->echelle)*image->ncol];
                  image->gray [k] = max;
            }
         }
      }

image->type_grille = QUINCONCE;
}
}

int main () {
   int height;
   int width;
   int nb_coul;
   int pos;
   int delta;
   char   ce_char;
   Cimageprov image_in;
   Cimageprov image_courante;
   Cimageprov image_niveau_n;

   int k,n,i,j;
   long nb_lu, nb_ec;
   char name[30];
   char name_out[3];
   FILE* fid_in;      /* image originale */
   FILE* fid_out;     /* approximation */
   FILE* fid_diff;    /* erreur de l'approximation */
   FILE* fid_xn;      /* image niveau n */

   printf ("Entrez le nom de l'image : \n");
   scanf("%s", name);
   printf ("Entrez un nb pair de niveaux de la pyramide :\n");
   scanf("%d", &n);

   if (n > NMAX) {
      printf ("Nb trop grand !\n");
      return 0;
   }

   fid_in = lire_image (name, &width, &height, &nb_coul, &image_in, &pos);
   fid_out = fopen ("approx.pgm", "w");
   if (((int) fid_out) == -1) {
      printf ("Pb fopen\n");
      return (0);
   }
   fid_diff = fopen ("diff.pgm", "w");
   if (((int) fid_diff) == -1) {
      printf ("Pb fopen\n");
      return (0);
   }
   fid_xn = fopen ("xn.pgm", "w");
   if (((int) fid_xn) == -1) {
      printf ("Pb fopen\n");
      return (0);
   }

   image_courante.gray = malloc (image_in.ncol*image_in.nrow);
   image_courante.ncol = image_in.ncol;
   image_courante.nrow = image_in.nrow;
   image_courante.echelle = 1;
   image_courante.type_grille = ENTIERE;
   for (k = 0; k < image_in.ncol*image_in.nrow; k++)
      image_courante.gray [k] = image_in.gray [k];

   for (i=0; i<n; i++) {
     calcul_niveau_suivant (&image_courante);
   }

sauver_image_niveau_n (&image_courante, &image_niveau_n);

   fprintf (fid_xn, "P5\n# CREATOR : FIFI\n%d %d\n255\n", image_niveau_n.ncol,
            image_niveau_n.nrow, fid_xn);
   nb_ec = fwrite (image_niveau_n.gray, image_niveau_n.ncol, image_niveau_n.nrow, fid_xn);

   for (i=0; i<n; i++) {
     calcul_niveau_prec (&image_courante);
   }

   fseek (fid_in, 0, SEEK_SET);
   for (k = 1; k <= pos; k++) {
      nb_lu = fread (&ce_char, 1, 1, fid_in);
      nb_ec = fwrite (&ce_char, 1, 1, fid_out);
   }

nb_ec = fwrite (image_courante.gray, width, height, fid_out);

for (k = 0; k < width*height; k++)
image_courante.gray [k] = image_in.gray [k] - image_courante.gray [k];

   fseek (fid_in, 0, SEEK_SET);
   for (k = 1; k <= pos; k++) {
      nb_lu = fread (&ce_char, 1, 1, fid_in);
      nb_ec = fwrite (&ce_char, 1, 1, fid_diff);
   }

nb_ec = fwrite (image_courante.gray, width, height, fid_diff);

   fclose (fid_in);
   fclose (fid_out);
   fclose (fid_diff);

free (image_courante.gray);
return (0);
}