# include "dvo.h"
void copy_region (GSCRegion *in, GSCRegion *out);

int imdata (int argc, char **argv) {
  
  int i, j, k;
  int Nimage, N, NPTS, skip, mode, TimeSelect;
  int n, Nregion, Ntregion, TimeFormat;
  int *subset, Nsubset;
  char filename[256];
  double trange;
  unsigned long tzero, start, stop, TimeReference;
  Image *image;
  Catalog catalog;
  GSCRegion *region, *tregion;
  Vector *vec;

  TimeSelect = FALSE;
  if (N = get_argument (argc, argv, "-time")) {
    remove_argument (N, &argc, argv);
    if (!str_to_time (argv[N], &tzero)) { 
      fprintf (stderr, "syntax error\n");
      return (FALSE);
    }
    remove_argument (N, &argc, argv);
    if (!str_to_dtime (argv[N], &trange)) { 
      fprintf (stderr, "syntax error\n");
      return (FALSE);
    }
    remove_argument (N, &argc, argv);
    TimeSelect = TRUE;
    fprintf (stderr, "searching in range %ds - %ds (%f seconds)\n", (int)tzero, (int)(tzero + trange), trange);
    if (trange > 0) {
      start = tzero;
      stop  = tzero + trange;
    } else {
      stop = tzero;
      start  = tzero + trange;
    }      
  }

  if (argc != 4) {
    fprintf (stderr, "USAGE: imdata (value) [-time t dt]\n");
    return (FALSE);
  }

  /* identify selection */
  mode = 0;
  if (!strcasecmp (argv[1], "ra")) 
    mode = 1;
  if (!strcasecmp (argv[1], "dec")) 
    mode = 2;
  if (!strcasecmp (argv[1], "mag")) 
    mode = 3;
  if (!strcasecmp (argv[1], "dmag")) 
    mode = 4;
  if (!strcasecmp (argv[1], "Mcal")) 
    mode = 5;
  if (!strcasecmp (argv[1], "Mrel")) 
    mode = 6;
  if (!strcasecmp (argv[1], "source")) 
    mode = 7;
  if (!strcasecmp (argv[1], "x")) 
    mode = 8;
  if (!strcasecmp (argv[1], "y")) 
    mode = 9;
  if (!strcasecmp (argv[1], "time")) 
    mode = 10;
  if (mode == 0) {
    fprintf (stderr, "value may be one of the following:\n");
    fprintf (stderr, " ra dR dec dD mag dmag Mrel Mcal source time\n");
    return (FALSE);
  }
  if ((vec = SelectVector (argv[1], ANYVECTOR, TRUE)) == NULL) return (FALSE);

  image = LoadImages (&Nimage);
  if (Nimage == 0) {
    return (FALSE);
  }
  image_subset (image, Nimage, &subset, &Nsubset, (Graphdata *) NULL, FALSE, tzero, trange, TimeSelect);
  GetTimeFormat (&TimeReference, &TimeFormat);

  /* for each image of interest, find the appropriate region files */
  Nregion = 0;
  ALLOCATE (region, GSCRegion, 1);
  for (i = 0; i < Nsubset; i++) {
    tregion = get_regions (&image[subset[i]], &Ntregion);
    REALLOCATE (region, GSCRegion, MAX (1, Nregion + Ntregion));
    for (j = 0; j < Ntregion; j++) {
      skip = FALSE;
      for (k = 0; (k < Nregion) && !skip; k++) {
	if (!strcmp (region[k].filename, tregion[j].filename)) skip = TRUE;
      }
      if (!skip) {
	copy_region (&tregion[j], &region[Nregion]);
	Nregion ++;
      }
    }
    free (tregion);
  }	
  free (subset);
  for (i = 0; i < Nregion; i++) {
    fprintf (stderr, "try %s\n", region[i].filename);
  } 
  
  /* create output vector */
  NPTS = 1000;
  REALLOCATE (vec[0].elements, float, NPTS);
  vec[0].Nelements = N = 0;

  /* for each region file, extract the data of interest in the right time range */
  for (j = 0; j < Nregion; j++) {

    /* get file name and open */
    strcpy (filename, region[j].filename);
    catalog.filename = filename;
    switch (lock_catalog (&catalog, LCK_SOFT)) {
    case 2:
      unlock_catalog (&catalog);
    case 0:
      continue;
    }
    if (!load_catalog (&catalog, LOAD_AVES, TRUE)) {
      unlock_catalog (&catalog);
      continue;
    }
    unlock_catalog (&catalog);

    /* assign vector values */
    switch (mode) {
    case (1):  /* ra */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	n = catalog.measure[i].averef;
	vec[0].elements[N] = catalog.average[n].R - catalog.measure[i].dR / 360000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (2):  /* dec */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	n = catalog.measure[i].averef;
	vec[0].elements[N] = catalog.average[n].D - catalog.measure[i].dD / 360000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (3):  /* mag */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	vec[0].elements[N] = catalog.measure[i].M / 1000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (4):  /* dmag */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	vec[0].elements[N] = catalog.measure[i].dM / 1000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (5):  /* Mcal */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	vec[0].elements[N] = catalog.measure[i].Mcal / 1000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (6):  /* Mrel */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	n = catalog.measure[i].averef;
	vec[0].elements[N] = catalog.average[n].M / 1000.0;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (7):  /* source */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	vec[0].elements[N] = catalog.measure[i].source;
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    case (10):  /* time */
      for (i = 0; i < catalog.Nmeasure; i++) {
	if ((catalog.measure[i].t < start) || (catalog.measure[i].t > stop)) continue;
	vec[0].elements[N] = TimeValue (catalog.measure[i].t, TimeReference, TimeFormat);
	N++;
	if (N == NPTS - 1) {
	  NPTS += 1000;
	  REALLOCATE (vec[0].elements, float, NPTS);
	}
      }
      break;
    }
    if (catalog.average != 0) {
      free (catalog.average);
    }
    if (catalog.measure != 0) {
      free (catalog.measure);
    }
  }
  
  vec[0].Nelements = N;
  REALLOCATE (vec[0].elements, float, MAX(1,N));
  return (TRUE);
}


void copy_region (GSCRegion *in, GSCRegion *out) {
    
  out[0].RA[0]  = in[0].RA[0];
  out[0].RA[1]  = in[0].RA[1];
  out[0].DEC[0] = in[0].DEC[0];
  out[0].DEC[1] = in[0].DEC[1];
  strcpy (out[0].filename, in[0].filename);
}

