# include "delstar.h"

void find_matches (Catalog *catalog, int photcode, int start, int end) {

  int i, j, k, n, m, N, M, drop;
  int *next, *next_miss, *ave_miss, last, last_miss;
  int Nave, Nmeas, Nmiss;
  int Nmeasfound, Nsecfilt;
  int this, prev;
  Measure *tmpmeasure;
  Missing *tmpmissing;

  /** allocate local arrays **/
  Nave = catalog[0].Naverage;

  Nmeas = catalog[0].Nmeasure;
  ALLOCATE (next, int, MAX(Nmeas,1));
  
  Nmiss = catalog[0].Nmissing;
  ALLOCATE (next_miss, int, MAX(Nmiss,1));
  ALLOCATE (ave_miss, int, MAX(Nmiss,1));
  
  if (VERBOSE) fprintf (stderr, "starting with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss);

  /* set up pointers for linked list of measure */
  for (i = 0; i < Nmeas - 1; i++) {
    next[i] = i+1;
  }
  next[i] = -1;
  last = i;
  /* set up pointers for linked list of missing */
  for (i = 0; i < Nmiss - 1; i++) {
    next_miss[i] = i+1;
  }
  next_miss[i] = -1;
  last_miss = i;
  /* set up references for missing to average */
  for (i = 0; i < Nave; i++) {
    for (j = 0; j < catalog[0].average[i].Nn; j++) {
      ave_miss[catalog[0].average[i].missing + j] = i;
    }
  }

  if (VERBOSE) fprintf (stderr, "deleting for range %d to %d, photcode %d\n", start, end, photcode);
  Nmeasfound = 0;
  Nsecfilt = catalog[0].Nsecfilt;

  if (VERBOSE) fprintf (stderr, "fixing measure...\n"); 

  /** find measure in time range **/
  this = prev = -1;
  for (i = 0; (i < Nmeas); i++) {
    if (VERBOSE && !(i % 10000)) fprintf (stderr, ". ");
    drop = TRUE;
    drop &= (catalog[0].measure[i].t >= start);
    drop &= (catalog[0].measure[i].t <= end);
    drop &= ((photcode == -1) || (photcode == catalog[0].measure[i].source));
    if (!drop) {
      prev = i;
      continue;
    }
    Nmeasfound ++;

    /* this star is in this image */
    this = next[i];
    next[i] = -2; /* we delete this one */
    if (prev != -1) { next[prev] = this; }

# if (0) 
    /* why is this section disabled? */
    /* fix the list links: connect the previous valid link to the next valid link */
    for (j = i; (j >= 0) && (next[j] == -2); j--); /* find previous entry to fix link */
    if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */
      if (next[j] != i) {
	fprintf (stderr, "error? (1)  this link seems to have been lost\n");
	fprintf (stderr, "j: %d, next[j]: %d, i: %d\n", j, next[j], i);
	exit (1);
      }
      /* find next valid entry to fix link */
      for (k = i; (k < Nmeas) && (next[k] == -2); k++);
      if (k < Nmeas)
	next[j] = k;
      else 
	next[j] = -1;  /* last link in list gets a -1 */
    }
# endif      

    /*** fix the corresponding average entry ***/
    n = catalog[0].measure[i].averef;
    if (catalog[0].average[n].Nm == 0) { /* this should never happen */
      fprintf (stderr, "error? we deleted one too many objects?\n");
      exit (1);
    }
    catalog[0].average[n].Nm --;
    /* this was only entry in list: will be deleted below.  meanwhile, delete all missing entries*/
    if ((catalog[0].average[n].Nm < 1) && (catalog[0].average[n].Nn > 0)) { 
      m = catalog[0].average[n].missing;
      for (j = 0; j < catalog[0].average[n].Nn; j++) {
	M = next_miss[m];
	next_miss[m] = -2;
	m = M;
      }
      m = catalog[0].average[n].missing;
      /* fix the list links: connect the previous valid link to the next valid link */
      for (j = m; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */
      if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */
	if (next_miss[j] != m) {
	  fprintf (stderr, "error? (2) this link seems to have been lost\n");
	  fprintf (stderr, "j: %d, next_miss[j]: %d, i: %d\n", j, next_miss[j], i);
	  exit (1);
	}
	/* find next valid entry to fix link */
	for (k = m; (k < Nmiss) && (next_miss[k] == -2); k++);
	if (k < Nmiss)
	  next_miss[j] = k;
	else 
	  next_miss[j] = -1;  /* last link in list gets a -1 */
      }
    }
    /* this was first entry in list */
    if ((catalog[0].average[n].offset == i) && (catalog[0].average[n].Nm > 0)) { 
      m = catalog[0].average[n].offset;
      /* find next valid entry -- notice lack of error checking... */
      for (j = 0; (j < Nmeas) && (next[m+j] == -2); j++);
      catalog[0].average[n].offset = m + j;
    }
  } 
  fprintf (stderr, "found %d meas to remove\n", Nmeasfound);

  if (VERBOSE) fprintf (stderr, "fixing missing..."); 
  /** find missing in time range of image **/
  for (i = 0; (i < Nmiss); i++) {
    if (next_miss[i] == -2) continue;
    if (catalog[0].missing[i].t < start) continue;
    if (catalog[0].missing[i].t > end) continue;

    next_miss[i] = -2; /* we delete this one */
    /* fix the list links: connect the previous valid link to the next valid link */
    for (j = i; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */
    if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */
      if (next_miss[j] != i) {
	fprintf (stderr, "error? (3) this link seems to have been lost\n");
	fprintf (stderr, "j: %d, next_miss[j]: %d, i: %d\n", j, next_miss[j], i);
	exit (1);
      }
      /* find next valid entry to fix link */
      for (k = i; (k < Nmiss) && (next_miss[k] == -2); k++);
      if (k < Nmiss)
	next_miss[j] = k;
      else 
	next_miss[j] = -1;  /* last link in list gets a -1 */
    }
      
    /*** fix the corresponding average entry ***/
    n = ave_miss[i];
    if (catalog[0].average[n].Nn == 0) { /* this should never happen */
      fprintf (stderr, "error? we deleted one too many missing?\n");
      exit (1);
    }
    catalog[0].average[n].Nn --;
    /* this was first entry in list */
    if ((catalog[0].average[n].missing == i) && (catalog[0].average[n].Nn > 0)) { 
      m = catalog[0].average[n].missing;
      for (j = 0; (j < Nmiss) && (next_miss[m+j] == -2); j++);
      catalog[0].average[n].missing = m + j;
    }
  }

  /* we should delete the references to blended image and cat stars ?? */
  /* or drop since we are changing this concept ?? */

  /* fix Average list: delete entries with Nm == 0 */
  for (i = j = 0; (i < Nave) && (j < Nave); i++, j++) {
    for (; (j < Nave) && (catalog[0].average[j].Nm == 0); j++);
    if ((i != j) && (j < Nave)) {
      catalog[0].average[i] = catalog[0].average[j];
      for (k = 0; k < catalog[0].Nsecfilt; k++) {
	catalog[0].secfilt[i*Nsecfilt + k] = catalog[0].secfilt[j*Nsecfilt + k];
      }
    }    
    if (j == Nave) i--;
  }
  Nave = i;
  REALLOCATE (catalog[0].average, Average, Nave);
  REALLOCATE (catalog[0].secfilt, SecFilt, MAX (1, Nave*Nsecfilt));

  /* fix order of Measure (memory intensive, but fast) */
  N = 0; 
  ALLOCATE (tmpmeasure, Measure, Nmeas);
  for (i = 0; i < Nave; i++) {
    n = catalog[0].average[i].offset;
    catalog[0].average[i].offset = N;
    for (k = 0; k < catalog[0].average[i].Nm; k++, N++) {
      if ((n == -1) || (n == -2)) {
	fprintf (stderr, "error: linked list is confused\n");
	exit (1);
      }
      tmpmeasure[N] = catalog[0].measure[n]; 
      tmpmeasure[N].averef = i;
      n = next[n];
    }
  }
  Nmeas = N;
  free (catalog[0].measure);
  catalog[0].measure = tmpmeasure;
  REALLOCATE (catalog[0].measure, Measure, Nmeas);
    
  /* fprintf (stderr, "fixing the mising order...\n"); */
  /* fix order of Missing (memory intensive, but fast) */
  N = 0; 
  ALLOCATE (tmpmissing, Missing, Nmiss);
  for (i = 0; i < Nave; i++) {
    if (catalog[0].average[i].Nn > 0) {
      n = catalog[0].average[i].missing;
      catalog[0].average[i].missing = N;
      for (k = 0; k < catalog[0].average[i].Nn; k++, N++) {
	if ((n == -1) || (n == -2)) {
	  fprintf (stderr, "error: linked list is confused\n");
	  exit (1);
	}
	tmpmissing[N] = catalog[0].missing[n]; 
	n = next_miss[n];
      }
    }
  }
  Nmiss = N;
  free (catalog[0].missing);
  catalog[0].missing = tmpmissing;
  REALLOCATE (catalog[0].missing, Missing, Nmiss);

  fprintf (stderr, "\n");

  catalog[0].Naverage = Nave;
  catalog[0].Nmeasure = Nmeas;
  catalog[0].Nmissing = Nmiss;
  if (VERBOSE) fprintf (stderr, "  ending with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss);

  free (next);
  free (next_miss);
  free (ave_miss);
  return;
}

