/* 
   ERIK S. WHEELER 
   cross_corr.cc
   Cross Correlation using linked lists read 
   from 2 binary data file.
   16 BIT DATA
*/

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

#define FALSE 0
#define MAXLENGTH_DEFAULT 32768
#define MAXHD 10000
typedef double datatype;
 
typedef struct listelement_type {
  datatype dataitem;
  struct listelement_type *next;
  struct listelement_type *prev;
} listelement;

listelement_type *AddItemOldList (listelement_type *finalptr, datatype data) {
  listelement_type *tempptr;

  if ((finalptr -> next = (struct listelement_type *) malloc (sizeof (listelement_type))) == NULL) {
    printf("\nMemory Allocation Error -- Execution Aborted\n");
    exit(1);
  }
  tempptr = finalptr;
  finalptr = finalptr -> next;
  finalptr -> next = NULL;
  finalptr -> prev = tempptr;
  finalptr -> dataitem = data;
  return finalptr;
}

listelement_type * AddItemNewList (listelement_type * listptr, datatype data) {

  if ((listptr = (struct listelement_type *) malloc (sizeof (listelement_type))) == NULL) {
    printf("\nMemory Allocation Error -- Execution Aborted\n");
    exit(1);
  }
  listptr -> next = NULL;
  listptr -> prev = NULL;
  listptr -> dataitem = data;
  return listptr;
}

datatype MultItem (listelement_type *ptr1, listelement_type *ptr2) {
  datatype returndata = 0;
  
  returndata = ptr1 -> dataitem * ptr2 -> dataitem;
  while (((ptr1 = ptr1-> next) != NULL) && ((ptr2 = ptr2 -> next) != NULL)) {
    returndata = returndata + ((ptr1 -> dataitem) * (ptr2 -> dataitem));
  }
  return returndata;
}

listelement_type * LargestItem(listelement_type * listptr) {
  listelement_type * tempptr = listptr;
  while (listptr != NULL) {
    if (abs((int)listptr -> dataitem) > abs((int)tempptr -> dataitem)) {
      tempptr = listptr;
    }
    listptr =  listptr -> next;
  }
  return tempptr;
}

void NormalizeToFirst (listelement_type * listptr) {
  datatype normalfactor = (listptr -> dataitem) / MAXHD;
  while (listptr != NULL) {
    listptr -> dataitem = (listptr -> dataitem) / normalfactor ;
    listptr = listptr -> next;
  }
}

void NormalizeToLargest (listelement_type * listptr) {
  listelement_type *tempptr;
  datatype normalfactor;
  tempptr = LargestItem(listptr);
  normalfactor = (tempptr -> dataitem) / MAXHD;
  while (listptr != NULL) {
    listptr -> dataitem = (listptr -> dataitem) / normalfactor ;
    listptr = listptr -> next;
  }
}

void PrintQueue (listelement_type * listptr) {
  if (listptr == NULL)
    printf ("Queue is empty!\n");
  else
    while (listptr != NULL) {
      printf ("%f\t", listptr -> dataitem);
      listptr = listptr -> next;
    }
  printf ("\n");
}

listelement_type * RemoveFirstItem (listelement_type * listptr) {
  listelement_type * tempptr;

  tempptr = listptr -> next;
  if (tempptr != NULL)  tempptr -> prev = NULL;
  free (listptr);
  return tempptr;
}

listelement_type * RemoveLastItem (listelement_type * finalptr) {
  listelement_type * tempptr;

  tempptr = finalptr -> prev;
  tempptr -> next = NULL;
  free (finalptr);
  return tempptr;
}

/*
void ClearQueueLtoF (listelement_type * finalptr) {
  while (finalptr != NULL) 
    finalptr = RemoveLastItem (finalptr);
}
*/

void ClearQueue (listelement_type * listptr) {
  while (listptr != NULL)
    listptr = RemoveFirstItem (listptr);
}

void Usage(char *argv[]) {
  printf("\nUsage: %s X-FILE H-FILE Y-FILE [-lxX-WINDOW] [-lhH-WINDOW] [-n]\n", argv[0]);
  printf("****************************************************\n");
  printf("X-FILE, H-FILE = binary files to be cross correlated\n");
  printf("Y-FILE         = normalized cross correlation of X-FILE & H-FILE\n");
  printf("[-xlX-WINDOW] = maximum X-FILE window size (integer)\n");
  printf("[-hlH-WINDOW] = maximum H-FILE window size (integer)\n");
  printf("Default window size = %d\n",MAXLENGTH_DEFAULT);
  printf("[-n] = Normalize (if NOT present, output normalizes to largest value\n"); 
  printf("Brackets [] denote optional\n");
  printf("Input and output files are SHORT INT (16 bit) (max value = %d)\n\n",MAXHD);
  exit(1); 
}

int main (int argc, char *argv[]) {
  listelement_type *listptrx = NULL, *currptrx = NULL, *finalptrx = NULL,
                   *listptrh = NULL, *currptrh = NULL, *finalptrh = NULL,
                   *listptry = NULL, *currptry = NULL, *finalptry = NULL;
  int i, fileflag, xwindow=MAXLENGTH_DEFAULT, hwindow=MAXLENGTH_DEFAULT;
  short int hd_data=0, norm=1;
  char *xfile = "xtemp.dat";
  char *hfile = "htemp.dat";
  char *yfile = "ytemp.dat";
  FILE *fp;

  if((argc<8)&&(argc>3)) {
    fileflag=0;
    for(i=1;i<argc;i++) {
      if (argv[i][0]=='-') {
	if ((argv[i][1]=='l')&&(argv[i][2]=='x')) {
	  printf("X-WINDOW = %d\n",(xwindow = atoi(argv[i] + 3) - 1));
	}
	else {
	  if ((argv[i][1]=='l')&&(argv[i][2]=='h')) {
	    printf("H-WINDOW = %d\n",(hwindow = atoi(argv[i] + 3) - 1));
	  }
	  else {
	    if ((argv[i][1]=='n')||(argv[i][1]=='N')) {
	      printf("Normalization OFF -- Warning Overflow may occur\n");
	      norm = FALSE;
	    }
	    else {
	      Usage(argv);
	    }
	  }
	}
      }
      else {
	fileflag++;
	switch (fileflag) {
	case 1:
	  printf("\nFile 1: %s\n",xfile = argv[i]);
	  break;
	case 2:
	  printf("File 2: %s\n",hfile = argv[i]);
	  break;
	case 3:
	  printf("File 3: %s\n",yfile = argv[i]);
	  break;
	default:
	  Usage(argv);
	}
      }
    }
    if (fileflag != 3) Usage(argv);
  }
  else Usage(argv);

  listptrx = NULL;
  currptrx = NULL;
  finalptrx = NULL;
  listptrh = NULL;
  currptrh = NULL;
  finalptrh = NULL;


/*  XFILE OPEN BINARY FILE FOR READING */
  if ((fp=fopen(xfile,"rb"))==NULL){
    printf("Cannot open file %s\n",xfile);
    exit(1);
  }

/*  XFILE READ DATA INTO LINKED LIST OF xwindow <= FILELENGTH */
  if (fread(&hd_data, sizeof(short int), 1, fp) != 1) {
    printf("Read Error\n");
    exit(1);
  }
  listptrx = AddItemNewList (listptrx, (datatype)hd_data);
  finalptrx = listptrx;

  i = 0;
  while (i++ < xwindow) {
    if (fread(&hd_data, sizeof(short int), 1, fp) != 1) {
      if (!feof(fp)) {
	printf("Read Error\n");
	exit(1);
      }
      else i = xwindow;
    }
    else finalptrx = AddItemOldList (finalptrx, (datatype)hd_data);
  }
  fclose(fp);

/*  HFILE OPEN BINARY FILE FOR READING */
  if ((fp=fopen(hfile,"rb"))==NULL){
    printf("Cannot open file %s\n",hfile);
    exit(1);
  }

/*  HFILE READ DATA INTO LINKED LIST OF hwindow <= FILELENGTH */
  if (fread(&hd_data, sizeof(short int), 1, fp) != 1) {
    printf("Read Error\n");
    exit(1);
  }
  listptrh = AddItemNewList (listptrh, (datatype)hd_data);
  finalptrh = listptrh;

  i = 0;
  while (i++ < hwindow) {
    if (fread(&hd_data, sizeof(short int), 1, fp) != 1) {
      if (!feof(fp)) {
	printf("Read Error\n");
	exit(1);
      }
      else i = hwindow;
    }
    else finalptrh = AddItemOldList (finalptrh, (datatype)hd_data);
  }
  fclose(fp);

/* CROSS CORRELATION FUNCTION */
  listptry = NULL;
  currptry = NULL;
  finalptry = NULL;
  currptrx = listptrx;
  currptrh = finalptrh;

  listptry = AddItemNewList (listptry, MultItem(currptrx,currptrh));
  finalptry = listptry;
  while (currptrh -> prev != NULL) {
    currptrh = currptrh -> prev;
    finalptry = AddItemOldList (finalptry, MultItem (currptrx, currptrh));
  }
  while (currptrx -> next != NULL) {
    currptrx = currptrx -> next;
    finalptry = AddItemOldList (finalptry, MultItem (currptrx,currptrh));
  }

  if (norm != FALSE) NormalizeToLargest (listptry);

/* Store Cross correlation in file fp as SHORT INT */
  if ((fp = fopen(yfile, "wb")) == NULL) {
    printf("Cannot Open File %s \n",yfile);
    exit(1);
  }
  currptry = listptry;
  while (currptry != NULL) {
    hd_data = (short int)(currptry -> dataitem);
    if (fwrite(&hd_data, sizeof(short int), 1, fp) != 1) {
      printf ("Read Error\n");
      exit(1);
    }
    currptry = currptry -> next;
  }
  fclose(fp);

  ClearQueue (listptrx); 
  ClearQueue (listptrh);
  ClearQueue (listptry);
  printf("\nCROSS CORRELATION of %s and %s to %s successful\n\n",xfile,hfile,yfile);
}