# include "psphot.h"

static psHash *timers = NULL;

// start/restart a named timer
bool psTimerStart (char *name) {

    psTime *start;

    if (timers == NULL) timers = psHashAlloc (16);

    start = psTimeGetTime (PS_TIME_UTC);

    psHashAdd (timers, name, start);

    return (TRUE);
}

// get current elapsed time on named timer
psF64 psTimerMark (char *name) {

    psTime *start;
    psTime *mark;
    psF64   delta;

    if (timers == NULL) return (0);

    start = psHashLookup (timers, name);
    if (start == NULL) return (0);

    mark = psTimeGetTime (PS_TIME_UTC);
    delta = psTimeDelta (mark, start);
    return (delta);
}

// we have log levels 1 (Error), 2 (Warning), 3 (Info), 4 (Details), 5 (Minutiae)
// 2 = default, -v = 3, -vv = 4, -vvv = 5 
psS32 psLogArguments (int *argc, char **argv) {
  
    int N, level;

    // default log level is 2
    level = 2;

    // set in order, so that -vvv overrides -vv overrides -v
    if ((N = get_argument (*argc, argv, "-v"))) {
	remove_argument (N, argc, argv);
	level = 3;
    }
    if ((N = get_argument (*argc, argv, "-vv"))) {
	remove_argument (N, argc, argv);
	level = 4;
    }
    if ((N = get_argument (*argc, argv, "-vvv"))) {
	remove_argument (N, argc, argv);
	level = 5;
    }

    if ((N = get_argument (*argc, argv, "-logfmt"))) {
	if (*argc < N + 2) {
	    psAbort ("psLogArguments", "USAGE: -logfmt (format)");
	}
	remove_argument (N, argc, argv);
	psLogSetFormat (argv[N]); // XXX EAM : this function should return an error if the log format is invalid
	remove_argument (N, argc, argv);
    }

    // set the level, return the level
    psLogSetLevel (level);
    return (level);
}

// set trace levels by facility
psS32 psTraceArguments (int *argc, char **argv) {
  
    int N;

    // argument format is: -trace (facil) (level) 
    while ((N = get_argument (*argc, argv, "-trace"))) {
	if (*argc < N + 3) {
	    psAbort ("psTraceArguments", "USAGE: -trace (facility) (level)");
	}
	remove_argument (N, argc, argv);
	psTraceSetLevel (argv[N], atoi(argv[N+1]));
	remove_argument (N, argc, argv);
	remove_argument (N, argc, argv);
    }
    if ((N = get_argument (*argc, argv, "-trace-levels"))) {
	psTracePrintLevels ();
	exit (2);
    }
    return (TRUE);
}

// find the location of the specified argument
int psArgumentGet (int argc, char **argv, char *arg) {

    int i;

    for (i = 0; i < argc; i++) {
	if (!strcmp(argv[i], arg))
	    return (i);
    }
  
    return ((int)NULL);
}

// remove the specified argument (by location)
int psArgumentRemove (int N, int *argc, char **argv) {

    int i;

    if ((N != (int)NULL) && (N != 0)) {
	(*argc)--;
	for (i = N; i < *argc; i++) {
	    argv[i] = argv[i+1];
	}
    }

    return (N);
}

// find the location of the specified argument
int get_argument (int argc, char **argv, char *arg) {

    int i;

    for (i = 0; i < argc; i++) {
	if (!strcmp(argv[i], arg))
	    return (i);
    }
  
    return ((int)NULL);
}

// remove the specified argument (by location)
int remove_argument (int N, int *argc, char **argv) {

    int i;

    if ((N != (int)NULL) && (N != 0)) {
	(*argc)--;
	for (i = N; i < *argc; i++) {
	    argv[i] = argv[i+1];
	}
    }

    return (N);
}

// alternate implementation of this function from pmObjects.c
psVector *psGetRowVectorFromImage(psImage *image,
				  psU32 row)
{
    PS_IMAGE_CHECK_NULL(image, NULL);
    PS_IMAGE_CHECK_TYPE(image, PS_TYPE_F32, NULL);

    psVector *tmpVector = psVectorAlloc(image->numCols, PS_TYPE_F32);
    memcpy (tmpVector->data.F32, image->data.F32[row], image->numCols*sizeof(psF32));
    return(tmpVector);
}

// extract config informatin from config data or from image header, if specified
psF32 pmConfigLookupF32 (bool *status, psMetadata *config, psMetadata *header, char *name) {

    char *source;
    char *keyword;
    psF32 value;
    psMetadataItem *item;

    // look for the entry in the config collection
    item = psMetadataLookup (config, name);
    if (item == NULL) {
	psLogMsg ("pmConfigLookupF32", 2, "no key %s in config data, trying as header keyword", name);
	value = psMetadataLookupF32 (status, header, name);
	if (*status == false) {
	    psAbort ("pmConfigLookupF32", "no key %s in header", name);
	}
	*status = true;
	return (value);
    }	

    // I'm either expecting a string, with the name "HD:keyword"...
    if (item->type == PS_META_STR) {
	source = item->data.V;
	if (!strncasecmp (source, "HD:", 3)) {
	    keyword = &source[3]; 
	    value = psMetadataLookupF32 (status, header, keyword);
	    if (*status == false) {
		psAbort ("pmConfigLookupF32", "no key %s in config", name);
	    }
	    *status = true;
	    return (value);
	}	
    }

    //  ... or a value (F32?)
    if (item->type == PS_META_F32) {
	value = item->data.F32;
	return (value);
    }

    *status = false;
    psError(PS_ERR_UNKNOWN, true, "invalid item");
    return (0);
}

void psImageSmooth (psImage *image, float sigma, float Nsigma) {

    int Nx, Ny, Npixel, Nrange;
    float factor, g, s;
    psVector *temp;

    // relevant terms
    Nrange = sigma*Nsigma + 0.5;
    Npixel = 2*Nrange + 1;
    factor = -0.5/(sigma*sigma);

    Nx = image->numCols;
    Ny = image->numRows;

    // generate gaussian
    psVector *gaussnorm = psVectorAlloc (Npixel, PS_TYPE_F32);
    for (int i = -Nrange; i < Nrange + 1; i++) {
	gaussnorm->data.F32[i+Nrange] = exp (factor*i*i);
    }
    psF32 *gauss = &gaussnorm->data.F32[Nrange];

    // smooth in X direction
    temp = psVectorAlloc (Nx, PS_TYPE_F32);
    for (int j = 0; j < Ny; j++) {
	psF32 *vi = image->data.F32[j];
	psF32 *vo = temp->data.F32;
	for (int i = 0; i < Nx; i++) {
	    g = s = 0;
	    for (int n = -Nrange; n < Nrange + 1; n++) {
		if (i+n < 0) continue;
		if (i+n >= Nx) continue;
		s += gauss[n]*vi[i+n];
		g += gauss[n];
	    }
	    vo[i] = s / g;
	}
	memcpy (image->data.F32[j], temp->data.F32, Nx*sizeof(psF32));
    }
    psFree (temp);

    // smooth in Y direction
    temp = psVectorAlloc (image->numRows, PS_TYPE_F32);
    for (int i = 0; i < Nx; i++) {
	psF32  *vo = temp->data.F32;
	psF32 **vi = image->data.F32;
	for (int j = 0; j < Ny; j++) {
	    g = s = 0;
	    for (int n = -Nrange; n < Nrange + 1; n++) {
		if (j+n < 0) continue;
		if (j+n >= Ny) continue;
		s += gauss[n]*vi[j+n][i];
		g += gauss[n];
	    }
	    vo[j] = s / g;
	}
	// replace temp in image
	for (int j = 0; j < Ny; j++) {
	    vi[j][i] = vo[j];
	}
    }
    psFree (temp);
    psFree (gaussnorm);
}
    
bool psImageInit (psImage *image,...) {

  va_list argp;
  psU8  vU8;
  psF32 vF32;
  psF64 vF64;

  if (image == NULL) return (false);

  va_start (argp, image);

  switch (image->type.type) {
    case PS_TYPE_U8:
      vU8 = va_arg (argp, psU32);

      for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	  image->data.U8[iy][ix] = vU8;
	}
      }
      break;

    case PS_TYPE_F32:
      vF32 = va_arg (argp, psF64);
      
      for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	  image->data.F32[iy][ix] = vF32;
	}
      }
      return (true);

    case PS_TYPE_F64:
      vF64 = va_arg (argp, psF64);

      for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	  image->data.F64[iy][ix] = vF64;
	}
      }
      return (true);

    default:
      psAbort ("psphot.psUtils", "datatype %d not defined in psImageInit\n", image->type);
      return (false);
  }
  return (false);
}	    

// define a square region centered on the given coordinate
psRegion *psRegionSquare (psF32 x, psF32 y, psF32 radius) {
    psRegion *region;
    region = psRegionAlloc (x - radius, x + radius + 1,    
			    y - radius, y + radius + 1);
    return (region);
}

// set actual region based on image parameters:
// compensate for negative upper limits
// XXX this is inconsistent: the coordindates should always be in the parent
//     frame, which means the negative values should subtract from Nx,Ny of
//     the parent, not the child.  but, we don't carry the dimensions of the 
//     parent in the psImage container.  for now, us the child Nx,Ny
// force range to be on this subimage 
// XXX EAM : this needs to be changes to use psRegion rather than psRegion*
psRegion *psRegionForImage (psRegion *out, psImage *image, psRegion *in) {
    
    // x0,y0, x1,y1 are in *parent* units

    if (out == NULL) {
	out = psRegionAlloc(in->x0, in->x1, in->y0, in->y1);
    } else {
	*out = *in;
    }
    // XXX these are probably wrong (see above)
    if (out->x1 <= 0) {
	out->x1 = image->col0 + image->numCols + out->x1;
    }
    if (out->y1 <= 0) {
	out->y1 = image->row0 + image->numRows + out->y1;
    }

    // force the lower-limits to be on the child
    out->x0 = PS_MAX(image->col0, out->x0);
    out->y0 = PS_MAX(image->row0, out->y0);

    // force the upper-limits to be on the child
    out->x1 = PS_MIN(image->col0 + image->numCols, out->x1);
    out->y1 = PS_MIN(image->row0 + image->numRows, out->y1);
    return (out);
}

// mask the area contained by the region
// the region is defined wrt the parent image
void psImageMaskRegion (psImage *image, psRegion *region, bool logical_and, int maskValue) {

    for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	    if (ix + image->col0 <  region->x0) continue;
	    if (ix + image->col0 >= region->x1) continue;
	    if (iy + image->row0 <  region->y0) continue;
	    if (iy + image->row0 >= region->y1) continue;
	    if (logical_and) {
		image->data.U8[iy][ix] &= maskValue;
	    } else {
		image->data.U8[iy][ix] |= maskValue;
	    }
	}
    }
}

// mask the area not contained by the region
// the region is defined wrt the parent image
void psImageKeepRegion (psImage *image, psRegion *region, bool logical_and, int maskValue) {

    for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	    if (ix + image->col0 <  region->x0) goto maskit;
	    if (ix + image->col0 >= region->x1) goto maskit;
	    if (iy + image->row0 <  region->y0) goto maskit;
	    if (iy + image->row0 >= region->y1) goto maskit;
	    continue;
	maskit:
	    if (logical_and) {
		image->data.U8[iy][ix] &= maskValue;
	    } else {
		image->data.U8[iy][ix] |= maskValue;
	    }
	}
    }
}

// mask the area contained by the region
// the region is defined wrt the parent image
void psImageMaskCircle (psImage *image, double x, double y, double radius, bool logical_and, int maskValue) {

    double dx, dy, r2, R2;

    R2 = PS_SQR(radius);

    for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	    dx = ix + image->col0 - x;
	    dy = iy + image->row0 - y;
	    r2 = PS_SQR(dx) + PS_SQR(dy);
	    if (r2 > R2) continue;
	    if (logical_and) {
		image->data.U8[iy][ix] &= maskValue;
	    } else {
		image->data.U8[iy][ix] |= maskValue;
	    }
	}
    }
}

// mask the area contained by the region
// the region is defined wrt the parent image
void psImageKeepCircle (psImage *image, double x, double y, double radius, bool logical_and, int maskValue) {

    double dx, dy, r2, R2;

    R2 = PS_SQR(radius);

    for (int iy = 0; iy < image->numRows; iy++) {
	for (int ix = 0; ix < image->numCols; ix++) {
	    dx = ix + image->col0 - x;
	    dy = iy + image->row0 - y;
	    r2 = PS_SQR(dx) + PS_SQR(dy);
	    if (r2 < R2) continue;
	    if (logical_and) {
		image->data.U8[iy][ix] &= maskValue;
	    } else {
		image->data.U8[iy][ix] |= maskValue;
	    }
	}
    }
}

psVector *psVectorCreate (double lower, double upper, double delta, psElemType type) {

  int nBin = (upper - lower) / delta;

  psVector *out = psVectorAlloc (nBin, type);
  
  for (int i = 0; i < nBin; i++) {
    out->data.F64[i] = lower + i * delta;
  }

  return (out);
}