/*
 * (C) 1997 René Scholz <mrz@informatik.uni-jena.de>
 *
 *
 * Bestimmt die Bild-Differenz zweier PGM5-Bilder
 * und gibt diesen als RMSE-Wert und PSNR-Wert aus.
 *
 * Außerdem kann das Differenzbild (ABS der Grauwerte der entsprechenden
 * Pixel) in eine Datei geschrieben werden.
 *
 * Es wird zusätzlich noch die Entropie der Bilder berechnet, allerdings
 * nur 1dim (nach der Def. der Entropie = -\sum p_i ld(p_i) )
 *
 * Kommentare nach "P5" sind möglich.
 * Beruht auf einem Code von C. Rahn
 *
*/

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

#define S_LEN 1024     /* Größe des Strings, der zum Lesen des Headers
                          benutzt wird */

#define NO_PGM5_FORMAT          1  /* Fehlercodes vereinbaren */
#define CANT_READ_FILE          2
#define CANT_WRITE_FILE         3
#define WRONG_SIZE_OF_PICTURES  4
#define FAILED_MALLOC           5

int    DEBUG=0;

unsigned int X=0,Y=0;     /* Bildgröße: XxY */
int    DiffPic=0;         /* if 1: schreibe Differenz-Bild */
int    Named_Pic2=0;      /* if 1: lese nicht von STDIO, sondern von -p2 */
int    Raw=0;             /* if 1: lese RAW PGM's */
long   Color_Position;    /* Stelle im Diff-File, an der der Color-Wert
                             des PGM-Files geschrieben werden muß */
char   *dp;               /* Filename des Differenz-Bilds */
char   *p1,*p2;           /* Filename der Bilder 1 und 2 */
FILE   *FDp1,*FDp2,*FDdp; /* die Filedescriptoren der Bilder */


void show_options(int return_code) {
  printf("Usage: pgmdiff -p1 picture1 [ -p2 picture2]  [options]\n\n");
  printf("Options:\n\n");
  printf("  -h | --help        Show this help.\n");
  printf("  -p1 pic1           read picture 1 from file <pic1>.\n");
  printf("  -p2 pic2 | < pic2  read picture 2 from file <pic2> or from STDIO.\n");
  printf("  [-dp diffpic]      write the difference picture <diffpic> (as PGM5).\n");
  printf("  -raw XxY           read RAW PGM5 pictures with size = XxY. (0<X,Y<=2048)\n\n");
  exit(return_code);
}



void show_help() {
  system("clear");
  printf("\n\t Calculates the RMSE and PSNR of two PGM5 pictures\n");
  printf("\t and (optionally) writes out the difference picture. \n\n");
  printf("\t (C) 1997   René Scholz    <http://www.thur.de/~Voland/>\n");
  printf("\t          & Christoph Rahn <rahn@informatik.uni-jena.de>\n\n");
  show_options(0);
}



void scan_options(int argc, char *argv[]) {
  int i; char *s,*s1;
  printf("\n");
  if(argc==1) show_help();

  for(i=1; i<argc; i++) {
    if(!strcmp("--help", argv[i]) || !strcmp ("-h", argv[i]) )
      show_help();

    if(!strcmp("-debug", argv[i]) ) { DEBUG=1; continue;         }

    if(!strcmp("-raw", argv[i]) ) {
      if(++i<argc) {
        s=argv[i];
        if ( (s1=strtok(s,"x")) !=NULL)    X=atoi(s1);  else goto M1;
        if ( (s1=strtok(NULL,"x")) !=NULL) Y=atoi(s1);  else goto M1;
        if ( (X*Y>0) && (X<=2048) && (Y<=2048) )   Raw=1;
        else  {
        M1: fprintf(stderr,"Error: wrong picture size!\n");
        show_options(-1);
        }
      }
      else {
        fprintf(stderr,"Error: picture size missing!\n");
        show_options(-1);
      }
      continue;
    }

    if(!strcmp("-dp", argv[i])) {
      if (++i<argc) { dp=argv[i]; DiffPic=1; }
      else {
        printf("\nError: file name of difference picture missing!\n");
        show_options(-1);
      }
      continue;
    }

    if(!strcmp("-p1", argv[i])) {
      if (++i<argc)  p1=argv[i];
      else {
        printf("\nError: file name of picture 1 missing!\n");
        show_options(-1);
      }
      continue;
    }

    if(!strcmp("-p2", argv[i])) {
      if (++i<argc) { p2=argv[i]; Named_Pic2=1; }
      else {
        printf("\nError: file name of picture 2 missing!\n");
        show_options(-1);
      }
      continue;
    }

    /* Wenn alles nicht passt: */
    printf("\nWrong Option [%s] !\n\n",argv[i]);
    show_options(-1);
  }
}



/* schließt bei einem Fehler alle nötigen Filedescriptoren: */
void close_FDs(int sub_error) {
  if ( (sub_error & 1) == 1) fclose(FDp1);
  if ( (sub_error & 2) == 2) if (Named_Pic2==1) fclose(FDp2);
  if ( (sub_error & 4) == 4) fclose(FDdp);
  exit(1);
}



/* sub_error gibt an, welche Filedescriptoren geschlossen werden sollen */
void print_error(int error_code,char *name,int sub_error) {
  switch (error_code) {
    case NO_PGM5_FORMAT: {
      fprintf(stderr,"%s is not a valid PGM5 file!\n",name);
      fprintf(stderr,"\nCorrect Header:\n\nP5\n[# Comments]\n");
      fprintf(stderr,"Width Height\nMax-Colors\n\n");
      fprintf(stderr,"Example:\n\nP5\n# bird.pgm from Picasso\n800 600\n256\n");
      fprintf(stderr,"\nis a picture with size 800x600 and 256 gray values.\n\n");
      close_FDs(sub_error);
    }
    case CANT_READ_FILE: {
      fprintf(stderr,"Cannot read from file %s !\n",name);
      close_FDs(sub_error);
    }
    case CANT_WRITE_FILE: {
      fprintf(stderr,"Cannot write to file %s !\n",name);
      close_FDs(sub_error);
    }
    case WRONG_SIZE_OF_PICTURES: {
      fprintf(stderr,"%s is not a valid PGM5 file,\n",name);
      fprintf(stderr,"or the pictures 1 and 2 have different sizes.!\n");
      fprintf(stderr,"\nCorrect Header:\n\nP5\n[# Comments]\n");
      fprintf(stderr,"Width Height\nMax-Colors\n\n");
      fprintf(stderr,"Example:\n\nP5\n# bird.pgm from Picasso\n800 600\n256\n");
      fprintf(stderr,"\nis a picture with size 800x600 and 256 gray values.\n\n");
      close_FDs(sub_error);
    }
    case FAILED_MALLOC: {
      fprintf(stderr,"Error: failed to allocate %u bytes!\n",S_LEN);
      exit(1);
    }
  }
}



/* versucht, die Bilddateien zu öffnen und überprüft diese, ob sie
 * korrekte PGM5-Bilder sind.
 * Anschließend werden die Header ausgelesen und verglichen.
 */
void open_pictures(void) {
  unsigned int X1,Y1,X2,Y2,max_gray;   /* die Bildgrößen */
  char *s;                             /* zum Lesen des Headers */

  if((FDp1=fopen(p1, "r"))==NULL)       print_error(CANT_READ_FILE,p1,0);
  if(Named_Pic2==1)
    { if ((FDp2=fopen(p2, "r"))==NULL)  print_error(CANT_READ_FILE,p2,1); }
  else { FDp2=stdin; p2="<stdin>"; }

  if (Raw==0) {  /* falls kein RAW-PGM, sondern PGM mit Header */
    if( (s=(char *)malloc(S_LEN*sizeof(char))) == NULL)
      print_error(FAILED_MALLOC,"",0);
    /* Öffne Bild1: */
    fscanf(FDp1,"%s\n",s);
    if(strcmp(s,"P5")) print_error(NO_PGM5_FORMAT,p1,1);
    /* Auslesen der Kommentare: */
    do fgets(s,S_LEN,FDp1); while ((char)s[0] == '#');
    
    sscanf(s,"%u %u\n",&X1,&Y1);
    fscanf(FDp1,"%u",&max_gray); fgetc(FDp1);

    if( (X1<=0) || (Y1<=0) || (max_gray<0) || (max_gray>255) )
      print_error(NO_PGM5_FORMAT,p1,1);
    
    /* Öffne Bild2: */
    fscanf(FDp2,"%s\n",s);
    if(strcmp(s,"P5")) print_error(NO_PGM5_FORMAT,p2,3);
    /* Auslesen der Kommentare: */
    do fgets(s,S_LEN,FDp2); while ((char)s[0] == '#');

    sscanf(s,"%u %u\n",&X2,&Y2);
    fscanf(FDp2,"%u",&max_gray); fgetc(FDp2);

    if( (X2!=X1) || (Y2!=Y1) || (max_gray<0) || (max_gray>255) )
      print_error(WRONG_SIZE_OF_PICTURES,p2,3);

    X=X2; Y=Y2;
  }

  if(DiffPic==1) {  /* falls diffpic geschrieben werden soll */
    if((FDdp=fopen(dp, "w+"))==NULL) print_error(CANT_WRITE_FILE,dp,3);
    fprintf(FDdp,"P5\n%s\n#\t%s\n#\t%s\n%u %u\n",
            "# difference picture of the pictures:",p1,p2,X,Y);
    Color_Position=ftell(FDdp);
    printf("Color_Position: %u\n", Color_Position);
    fprintf(FDdp,"%u\n",255);
    /* am Programmende 255 durch den hellsten Pixel ersetzen ! */
  }
}



/* Bildpunkte vergleichen und differenzen aufsummieren, PSNR berechnen */
void calculate() {
  register unsigned char A,B,t;
  register unsigned long int i,Pixels;
  register unsigned long int sum=0,p;
  unsigned long int histogramA[256]={0};
  unsigned long int histogramB[256]={0};
  double RMSE, PSNR;
  register double eA=0,eB=0,f, l=log(2);
  register unsigned int Max_Gray_Value=0, Hellster_Pixel=0;

  Pixels=(unsigned long)X*(unsigned long)Y;
  printf("calculate ...\n");

  if(DiffPic==0)
    for(i=0; i<Pixels; i++){
      A=fgetc(FDp1); B=fgetc(FDp2);
      if(A>Max_Gray_Value) Max_Gray_Value=A;
      if(B>Max_Gray_Value) Max_Gray_Value=B;

      t=(A>=B) ? A-B : B-A;
      sum+=t*t;
      histogramA[A]++; histogramB[B]++;
    }
  else /* Differenzbild erzeugen */
    for(i=0; i<Pixels; i++) {
      A=fgetc(FDp1); B=fgetc(FDp2);
      if(A>Max_Gray_Value) Max_Gray_Value=A;
      if(B>Max_Gray_Value) Max_Gray_Value=B;

      t=(A>=B) ? A-B : B-A;
      sum+=t*t;
      histogramA[A]++; histogramB[B]++;
      if(t>Hellster_Pixel) Hellster_Pixel=t;
      fprintf(FDdp,"%c",(char)t);
    }

  /* Ergebnisse ausrechnen: */
  RMSE=sqrt( (double)sum / (double)Pixels );
  /* PSNR=20*log10( (double)Max_Gray_Value / RMSE ); */
  if(RMSE>0) PSNR=20*log10( (double)255.0 / RMSE );
  else PSNR=0;

  for(i=0; i<256; i++)
    if(p=histogramA[i])
      { f=(double)p/(double)Pixels; eA-=f*log(f)/l; }

  for(i=0; i<256; i++)
    if(p=histogramB[i])
      { f=(double)p/(double)Pixels; eB-=f*log(f)/l; }


  fclose(FDp1);
  if(Named_Pic2==1) fclose(FDp2);

  if(DiffPic==1) {
    if(Hellster_Pixel<255) {
      if(Hellster_Pixel==0) Hellster_Pixel=1;
      fseek(FDdp,Color_Position,SEEK_SET);
      (Hellster_Pixel<100) ?
        (Hellster_Pixel<10) ? fprintf(FDdp,"00%u",Hellster_Pixel)
        : fprintf(FDdp,"0%u",Hellster_Pixel)
        : fprintf(FDdp,"%u",Hellster_Pixel);
      fseek(FDdp,0,SEEK_END);
    }
    fclose(FDdp);
  }


  
  if(PSNR>0) printf("\nRMSE=%f\nPSNR=%f\n\n",RMSE,PSNR);
  else printf("\nRMSE=%f\nPSNR=Inf\n\n",RMSE);
  printf("Entropie (1dim) von Bild %s = %f bit.\n",p1,eA);
  printf("Entropie (1dim) von Bild %s = %f bit.\n",p2,eB);
  if(DEBUG) printf("\npixels=%li sum=%li\n",Pixels,sum);
}



int main(int argc, char *argv[])
{
  scan_options(argc,argv);
  open_pictures();
  calculate();
  return(0);
}