# include "dvo.h"

int pmeasure (int argc, char **argv) {
  
  FILE *f;
  int i, j, k, m, N;
  int InRegion, Nregions, Ngraph;
  char filename[128], catdir[256];
  double Mz, Mr, mag;
  double Radius, Rmin, Rmax;
  unsigned IDclip, IDchoice, LimExclude;
  unsigned FlagChoice, FlagClip, flags;

  Catalog catalog;
  RegionFile *regions;
  Graphdata graphmode;
  Vector Xvec, Yvec, Zvec;

  Ngraph = -1;
  if (!GetGraph (&graphmode, NULL, &Ngraph)) return (FALSE);
  if (VarConfig ("CATDIR", "%s", catdir) == NULL) return (FALSE);
  if (!InitPhotcodes ()) return (FALSE);

  regions = (RegionFile *) NULL;
  f = (FILE *) NULL;
  Mz = 17.0;
  Mr = -5.0;
  Rmin = graphmode.coords.crval1 - 182.0;
  Rmax = graphmode.coords.crval1 + 182.0;

  IDclip = FALSE;
  if (N = get_argument (argc, argv, "-ID")) {
    remove_argument (N, &argc, argv);
    IDchoice  = atoi(argv[N]);
    remove_argument (N, &argc, argv);
    IDclip = TRUE;
  }
  FlagClip = FALSE;
  if (N = get_argument (argc, argv, "-flag")) {
    remove_argument (N, &argc, argv);
    FlagChoice  = atoi(argv[N]);
    remove_argument (N, &argc, argv);
    FlagClip = TRUE;
  }

  LimExclude = FALSE;
  if (N = get_argument (argc, argv, "-x")) {
    remove_argument (N, &argc, argv);
    LimExclude = TRUE;
  }

  if (N = get_argument (argc, argv, "-m")) {
    remove_argument (N, &argc, argv);
    Mr  = atof(argv[N]);
    remove_argument (N, &argc, argv);
    Mz = atof(argv[N]);
    Mr = Mr - Mz;
    remove_argument (N, &argc, argv);
  }

  InRegion = FALSE;
  if (N = get_argument (argc, argv, "-all")) {
    remove_argument (N, &argc, argv);
    InRegion = TRUE;
  }

  if ((InRegion || (argc != 2)) && (!InRegion || (argc != 1))) {
    fprintf (stderr, "USAGE: catalog (filename / -all) [-m M M] [-n N N] [-g] [-a RA DEC MAG] \n");
    return (FALSE);
  }
  
  if (InRegion) {
    Radius = MAX (fabs(graphmode.xmax), fabs(graphmode.ymax));
    regions = find_regions (graphmode.coords.crval1, graphmode.coords.crval2, Radius, &Nregions);
  } else {
    Nregions = 1;
  }
  
  for (j = 0; j < Nregions; j++) {
    catalog.average = 0;
    
    if (InRegion) {
      sprintf (filename, "%s/%s", catdir, regions[j].name);
    } else {
      sprintf (filename, "%s/%s", catdir, argv[1]);
    }
    
    /* lock, load, unlock catalog */
    catalog.filename = filename;
    switch (lock_catalog (&catalog, LCK_SOFT)) {
    case 2:
      unlock_catalog (&catalog);
    case 0:
      continue;
    }
    if (!load_catalog (&catalog, LOAD_AVES | LOAD_MEAS, TRUE)) {
      unlock_catalog (&catalog);
      continue;
    }
    unlock_catalog (&catalog);

    /* data has been loaded, get ready to plot it */
    Xvec.Nelements = catalog.Nmeasure;
    Yvec.Nelements = catalog.Nmeasure;
    Zvec.Nelements = catalog.Nmeasure;
    ALLOCATE (Xvec.elements, float, Xvec.Nelements);
    ALLOCATE (Yvec.elements, float, Yvec.Nelements);
    ALLOCATE (Zvec.elements, float, Zvec.Nelements);

    /* project stars to screen display coords */
    Xvec.Nelements = 0;
    for (N = i = 0; i < catalog.Naverage; i++) {
      if (IDclip && (catalog.average[i].code != IDchoice)) continue;
      while (catalog.average[i].R < Rmin) catalog.average[i].R += 360.0;
      while (catalog.average[i].R > Rmax) catalog.average[i].R -= 360.0;
      m = catalog.average[i].offset;
      for (k = 0; k < catalog.average[i].Nm; k++) {
	if (FlagClip) {
	  flags = catalog.measure[m+k].flags;
	  if (!(flags & FlagChoice)) continue;
	}
	mag = PhotCat (&catalog.measure[m+k]);
	Zvec.elements[N] = MIN (1.0, MAX (0.01, (mag - Mz) / Mr));
	if (LimExclude && (Zvec.elements[N] > 0.99)) continue;
	if (Zvec.elements[N] < 0.011) continue;
	fRD_to_XY (&Xvec.elements[N], &Yvec.elements[N], catalog.average[i].R - 2.777e-6*catalog.measure[m+k].dR, catalog.average[i].D - 2.777e-6*catalog.measure[m+k].dD, &graphmode.coords);
	N ++;
      }
    }
    Zvec.Nelements = Yvec.Nelements = Xvec.Nelements = N;

    graphmode.style = 2; /* set style to points */
    graphmode.size = -1; /* point size determined by Zvec */
    graphmode.etype = 0; /* no errorbars */
    PrepPlotting (N, &graphmode);
    
    PlotVector (N, Xvec.elements);
    PlotVector (N, Yvec.elements);
    PlotVector (N, Zvec.elements);

    free (Xvec.elements);
    free (Yvec.elements);
    free (Zvec.elements);

    if (catalog.average != 0) free (catalog.average);
    if (catalog.measure != 0) free (catalog.measure);
  }
  return (TRUE);

}
  
