# include "relphot.h"
double BigDecBounds[] = {0.0, 7.5, 15.0, 22.5, 30.0, 37.5, 45.0, 
			 52.5, 60.0, 67.5, 75.0, 82.5, 90.0,
			 0.0, -7.5, -15.0, -22.5, -30.0, -37.5, -45.0, 
			 -52.5, -60.0, -67.5, -75.0, -82.5, -90.0};

double Areas[] = {2700.0, 2665.4505, 2563.567, 2399.5344, 2181.6343, 1920.7299, 1629.5942, 1322.1251, 1012.5, 714.32555, 439.83856, 199.21043, 0.0,
		  2700.0, 2665.4505, 2563.567, 2399.5344, 2181.6343, 1920.7299, 1629.5942, 1322.1251, 1012.5, 714.32555, 439.83856, 199.21043, 0.0};

char *DecSections[] = {"N0000", "N0730", "N1500", "N2230", "N3000", "N3730", "N4500", 
		       "N5230", "N6000", "N6730", "N7500", "N8230", "weirdness", 
		       "S0000", "S0730", "S1500", "S2230", "S3000", "S3730", "S4500", 
		       "S5230", "S6000", "S6730", "S7500", "S8230", "weirdness"};

char *Dec2Sections[] = {"n0000", "n0730", "n1500", "n2230", "n3000", "n3730", "n4500", 
			"n5230", "n6000", "n6730", "n7500", "n8230", "weirdness", 
			"s0000", "s0730", "s1500", "s2230", "s3000", "s3730", "s4500", 
			"s5230", "s6000", "s6730", "s7500", "s8230", "weirdness"};

char *disk[] = {"disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "disk 1", 
		"disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "weirdness", 
		"disk 1", "disk 2", "disk 2", "disk 2", "disk 2", "disk 2", "disk 2", 
		"disk 2", "disk 2", "disk 2", "disk 2", "disk 2", "weirdness"};

int NBigRASections [] = {48, 47, 45, 43, 40, 36, 32, 28, 21, 15, 9, 3, 3, 48, 47, 45, 43, 40, 36, 32, 28, 21, 15, 9, 3, 3};

int NDecLines[] = {593, 584, 551, 530, 522, 465, 406, 362, 280, 198, 123, 24, 
                   0, 597, 578, 574, 577, 534, 499, 442, 376, 294, 212, 144, 48};

/* returns a list of region files within the desired RA, DEC region */
GSCRegion *get_regions (double minRa, double maxRa, double minDec, double maxDec, int *Nregions) {
  
  char buffer[28800], temp[50];
  GSCRegion *regions;
  FILE *f;
  double RA0, RA1, DEC0, DEC1;
  int i, j, NBigDec;
  int NLINES, done, NREGIONS, nregions;
  
  f = fopen (GSCFILE, "r");
  if (f == NULL) {
    fprintf (stderr, "ERROR: can't find regions file %s\n", GSCFILE);
    exit (0);
  }
  
  NREGIONS = 50;
  ALLOCATE (regions, GSCRegion, NREGIONS);
  nregions = 0;

  if ((minDec <= -90) || (maxDec >= 90)) {
    minRa = 0;
    maxRa = 360;
  }
  
  /* use the pole regions, if near pole */
  if (maxDec > 86.25) {
    sprintf (regions[nregions].filename, "n8230/pole.cpt\0");
    regions[nregions].RA[0] = 0;
    regions[nregions].RA[1] = 360;
    regions[nregions].DEC[0] = 86.25;
    regions[nregions].DEC[1] = 90.0;
    nregions ++;
    if (nregions == NREGIONS) {
      NREGIONS += 50;
      REALLOCATE (regions, GSCRegion, NREGIONS);
    }
  }

  if (minDec > 86.25) {
    *Nregions = nregions;
    return (regions);
  }
    
  if ((minDec < 0) && (maxDec > 0)) {
    /* Search Both Sides */
    NBigDec = 0;
  } else {
    /* find large DEC region (directory) */
    NBigDec = -1;
    for (i = 0; i < 12; i++) {
      if ((minDec >= BigDecBounds[i]) && (minDec < BigDecBounds[i+1])) {
	NBigDec = i;
	break;
      }
    }
    if (NBigDec < 0) {
      for (i = 13; i < 24; i++) {
	if ((maxDec < BigDecBounds[i]) && (maxDec >= BigDecBounds[i+1])) {
	  NBigDec = i;
	  break;
	}
      }
    }
  }
  if (NBigDec < 0) {
    fprintf (stderr, "ERROR: Dec out of range: %f\n", minDec);
    *Nregions = 0;
    return ((GSCRegion *) NULL);
  }
  
  /* count lines before section */
  NLINES = 0;
  for (i = 0; i < NBigDec; i++) {
    NLINES += NDecLines[i];
  }
  fseek (f, 5*2880 + 48*NLINES, SEEK_SET);
  
  /* should be in this section.  if not, there is a problem counting... */
  /* careful with the 0,360.0 boundary **/
  done = FALSE;
  for (j = 0; !done && (NBigDec + j < 25); j++) {
    fread (buffer, 48*NDecLines[NBigDec + j], 1, f);
    for (i = 0; (i < NDecLines[NBigDec + j]); i++) {
      strncpy (temp, &buffer[i*48], 48);
      temp[49] = 0;
      hms_to_deg (&RA0, &RA1, &DEC0, &DEC1, &temp[7]);
      if (RA1 < RA0) RA1 += 360.0;
      if ((DEC1 > 0) && (minDec < DEC1) && (maxDec > DEC0) && (minRa < RA1) && (maxRa > RA0)) {
	temp[5] = 0;
	sprintf (regions[nregions].filename, "%s/%s.cpt\0", Dec2Sections[NBigDec + j], &temp[1]);
	regions[nregions].RA[0] = RA0;
	regions[nregions].RA[1] = RA1;
	regions[nregions].DEC[0] = DEC0;
	regions[nregions].DEC[1] = DEC1;
	nregions ++;
	if (nregions == NREGIONS) {
	  NREGIONS += 50;
	  REALLOCATE (regions, GSCRegion, NREGIONS);
	}
      }
      if ((DEC1 < 0) && (minDec < DEC0) && (maxDec > DEC1) && (minRa < RA1) && (maxRa > RA0)) {
	temp[5] = 0;
	sprintf (regions[nregions].filename, "%s/%s.cpt\0", Dec2Sections[NBigDec + j], &temp[1]);
	regions[nregions].RA[0] = RA0;
	regions[nregions].RA[1] = RA1;
	regions[nregions].DEC[0] = DEC0;
	regions[nregions].DEC[1] = DEC1;
	nregions ++;
	if (nregions == NREGIONS) {
	  NREGIONS += 50;
	  REALLOCATE (regions, GSCRegion, NREGIONS);
	}
      }
      if (((DEC1 > 0) && (maxDec <= DEC1)) || ((DEC1 < 0) && (minDec >= DEC1))) {
	done = TRUE;
      }
    }
    if (done && (minDec < 0) && (maxDec > 0) && (BigDecBounds[NBigDec + j + 1] > 0)) {
      /* skip remaining north sections, try south sections */
      /* count lines before section */
      NLINES = 0;
      for (i = 0; i < 13; i++) { 
	NLINES += NDecLines[i];
      }
      fseek (f, 5*2880 + 48*NLINES, SEEK_SET);
      done = FALSE;
      j = 12;
    }
  }

  REALLOCATE (regions, GSCRegion, MAX (1,nregions));
  *Nregions = nregions;
  
  fclose (f);
  return (regions);
  
}

/**********/
int hms_to_deg (double *h0, double *h1, double *d0, double *d1, char *string) {
  
  int flag_d0, flag_d1, flag_h0, flag_h1;
  double tmp;
  
  *d0 = *h0 = *d1 = *h1 = 0;

  flag_h0 = dparse (h0, 1, string);
  flag_h1 = dparse (h1, 4, string);
  flag_d0 = dparse (d0, 7, string);
  flag_d1 = dparse (d1, 9, string);
  *h0 *= flag_h0;
  *h1 *= flag_h1;
  *d0 *= flag_d0;
  *d1 *= flag_d1;

  dparse (&tmp, 2, string);
  *h0 += tmp/60.0;
  dparse (&tmp, 3, string);
  *h0 += tmp/3600.0;
  
  dparse (&tmp, 5, string);
  *h1 += tmp/60.0;
  dparse (&tmp, 6, string);
  *h1 += tmp/3600.0;
  
  dparse (&tmp, 8, string);
  *d0 += tmp/60.0;

  dparse (&tmp, 10, string);
  *d1 += tmp/60.0;

  *h0 *= 15*flag_h0;
  *h1 *= 15*flag_h1;
  *d0 *= flag_d0;
  *d1 *= flag_d1;

  return (TRUE);
}
