static char *IMREGIONS_C="\n/*\n    NB: MAKE SURE YOU EDIT THE TEMPLATE FILE!!!!\n*/\n\n#ifndef FILTER_PTYPE\n#include <regions.h>\n#endif\n\n#ifndef UNUSED\n#ifdef __GNUC__\n#  define UNUSED(x) UNUSED_ ## x __attribute__((__unused__))\n#else\n#  define UNUSED(x) UNUSED_ ## x\n#endif\n#endif\n\n\n#define USE_ASTRO_ANGLE 0\n\n\n#define USE_FPU_DOUBLE 0\n#if USE_FPU_DOUBLE\n#include <fpu_control.h>\nstatic fpu_control_t _cw;\n#define FPU_DOUBLE {fpu_control_t _cw2; _FPU_GETCW(_cw); _cw2 = _cw & ~_FPU_EXTENDED; _cw2 |= _FPU_DOUBLE; _FPU_SETCW(_cw2);}\n#define FPU_RESTORE {_FPU_SETCW(_cw);}\n#else\n#define FPU_DOUBLE\n#define FPU_RESTORE\n#endif\n\n\nstatic int imregno=0;\nvoid initimregions(void)\n{\n  imregno++;\n  return;\n}\n\nstatic void markx(GFilt g, int UNUSED(sno), int flag, int type, int x, int y)\n{\n  \n  if( type == TOK_EREG )\n    return;\n  \n  if( flag ){\n    if( x <= g->x0s[y] )\n      g->x0s[y] = max(x,g->x0);\n    if( x >= g->x1s[y] )\n      g->x1s[y] = min(x,g->x1);\n  }\n  \n  else{\n    g->x0s[y] = g->x0;\n    g->x1s[y] = g->x1;\n  }\n}\n\nstatic void marky(GFilt g, int sno, int flag, int type)\n{\n  int i;\n  \n  if( type == TOK_EREG )\n    return;\n  \n  if( flag ){\n    for(i=max(g->y0,g->shapes[sno].ystart);\n	i<=min(g->y1,g->shapes[sno].ystop);\n	i++)\n      g->ybuf[i] = 1;\n  }\n  \n  else{\n    for(i=g->y0; i<=g->shapes[sno].ystart-1; i++)\n      g->ybuf[i] = 1;\n    for(i=g->shapes[sno].ystop+1; i<=g->y1; i++)\n      g->ybuf[i] = 1;\n  }\n}\n\nstatic int\nimagemaskcmp(const void *s1, const void *s2)\n{\n  FilterMask f1 = (FilterMask)s1;\n  FilterMask f2 = (FilterMask)s2;\n  \n  if( f1->y < f2->y ){\n    return -1;\n  }\n  else if( f1->y > f2->y ){\n    return 1;\n  }\n  else{\n    if( f1->xstart < f2->xstart ){\n      return -1;\n    }\n    else{\n      return 1;\n    }\n  }\n}\n\n\n\nstatic void quadeq(double a, double b, double c,\n		   double *x1, double *x2, int *nr, int *nc)\n{\n  double dis, q;\n  if( feq(a,0.0) ){\n    *nc = 0;\n    if( feq(b,0.0) ){\n      *nr = 0; *x1 = 0.0;\n    }\n    else{\n      *nr = 1; *x1 = -c / b;\n    }\n    *x2 = *x1;\n  }\n  else{\n    dis = b*b - 4.0 * a * c;\n    if( dis > 0.0 ){\n      *nr = 2; *nc = 0;\n      dis = sqrt(dis);\n      if( b < 0.0 ) dis = -dis;\n      q = -0.5 * (b + dis);\n      *x1 = q/a; *x2 = c/q;\n      if(*x1 > *x2){\n	q = *x1; *x1 = *x2; *x2 = q;\n      }\n    } \n    else if( feq(dis,0.0) ){\n      *nr = 1; *nc = 0; *x1 = - 0.5 * b / a; *x2 = *x1;\n    }\n    else{\n      *nr = 0; *nc = 2; *x1 = - 0.5 * b / a; *x2 = 0.5 * sqrt(-dis) / a;\n    }\n  }\n}\n\nstatic void rgs_mark(GFilt g, Scan *scanlist,\n		     int sno, int flag, int type, int xval, int yval )\n{\n  Scan scanmark, mark;\n  \n  if(yval < g->y0) yval = g->y0;\n  if(yval > g->y1) yval = g->y1;\n  mark = (Scan)calloc(1, sizeof(ScanRec));\n  mark->x = xval;\n  \n  if( !scanlist ) return;\n  \n  if( !scanlist[yval] || ((scanlist[yval])->x > xval) ){\n    mark->next = scanlist[yval];\n    scanlist[yval] = mark;\n  } else {\n    scanmark = scanlist[yval];\n    while( scanmark->next && ((scanmark->next)->x < xval) )\n      scanmark = scanmark->next;\n    mark->next = scanmark->next;\n    scanmark->next = mark;\n  }\n  markx(g, sno, flag, type, xval, yval);\n}\n\nstatic void rgs_segment(GFilt g, Scan *scanlist, int width, int height,\n			int sno, int flag, int type,\n			double x1, double y1, double x2, double y2)\n{\n  int ystart, ystop, yval, xval;\n  double invslope, xoffset;\n\n  ystart = PIXINCL(y1);\n  if( ystart < 1 ) ystart = 1;\n  \n  ystop = PIXINCL(y2) - 1;\n  if( ystop > height ) ystop = height;\n  \n  if( (ystart > ystop) || (ystop < 1) )\n    return;\n  \n  invslope = (x1 - x2) / (y1 - y2);\n  xoffset = x1 + ((ystart - y1) * invslope);\n  for(yval=ystart; yval<=ystop; yval++){\n    xval = PIXINCL(xoffset);\n    \n    if(xval < 1)     xval = 1;\n    if(xval > width) xval = width + 1;\n    rgs_mark(g, scanlist, sno, flag, type, xval, yval);\n    xoffset = xoffset + invslope;\n  }\n}\n\nstatic void _polygoni(GFilt g, int qt, int UNUSED(rno), int sno, int flag,\n                      int type, double *vx, double *vy, int count)\n{\n  int i, j;\n  double xlo, xhi;\n  double ylo, yhi;\n\n#ifdef ALREADY_DONE\n  \n  for(i=0; i<count; i++){\n    vx[i] = (vx[i] - g->xmin)/g->block + 1.0;\n    vy[i] = (vy[i] - g->ymin)/g->block + 1.0;\n  }\n#endif\n  \n  xlo = vx[0];\n  xhi = xlo;\n  ylo = vy[0];\n  yhi = ylo;\n  for(i=0; i<count; i++){\n    if(vx[i] > xhi) xhi = vx[i];\n    if(vx[i] < xlo) xlo = vx[i];\n    if(vy[i] > yhi) yhi = vy[i];\n    if(vy[i] < ylo) ylo = vy[i];\n  }\n  if( qt && (sno > 1) ){\n    g->shapes[sno].ystart = g->shapes[sno-1].ystart;\n    g->shapes[sno].ystop = g->shapes[sno-1].ystop;\n  }\n  else{\n    g->shapes[sno].ystart = max(g->y0,PIXINCL(ylo));\n    g->shapes[sno].ystop = min(g->y1,PIXINCL(yhi) - 1);\n  }\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  marky(g, sno, flag, type);\n\n  \n  \n  j = count-1;\n  for(i=0; i<count; i++){\n    \n    if(vy[i] > vy[j]){\n      rgs_segment(g, g->shapes[sno].scanlist, g->x1, g->y1,\n		  sno, flag, type, vx[j], vy[j], vx[i], vy[i]);\n    }\n    else{\n      rgs_segment(g, g->shapes[sno].scanlist, g->x1, g->y1,\n		  sno, flag, type, vx[i], vy[i], vx[j], vy[j]);\n    }\n    j = i;\n  }\n}\n\nstatic int corner_vertex(int index, int width, int height,\n			 double *x, double *y)\n{\n  switch (index) {\n  case 1:\n    *x = 0.0;\n    *y = height + 1;\n    break;\n  case 2:\n    *x = 0.0;\n    *y = 0.0;\n    break;\n  case 3:\n    *x = width + 1;\n    *y = 0.0;\n    break;\n  case 4:\n    *x = width + 1;\n    *y = height + 1;\n  default:\n    break;\n  }\n  index = index + 1;\n  if(index > 4) index = 1;\n  return(index);\n}\n\nstatic int pie_intercept(int width, int height, double xcen, double ycen,\n			 double angle, double *xcept, double *ycept)\n{\n  double angl, slope;	\n  angl = angle;\n  \n  while (angl < 0.0)\n    angl = angl + 360.0;\n  while (angl >= 360.0)\n    angl = angl - 360.0;\n  \n#if USE_ASTRO_ANGLE\n  if(fabs(angl - 90.0) < SMALL_NUMBER) {\n#else\n  if(fabs(angl - 180.0) < SMALL_NUMBER) {\n#endif\n    *xcept = 0.0;\n    *ycept = ycen;\n    return(2);\n  }\n#if USE_ASTRO_ANGLE\n  if(fabs(angl - 270.0) < SMALL_NUMBER) {\n#else\n  if(fabs(angl - 0.0) < SMALL_NUMBER) {\n#endif\n    *xcept = width + 1;\n    *ycept = ycen;\n    return(4);\n  }\n#if USE_ASTRO_ANGLE\n  \n  angl = angl + 90.0;\n#endif\n  if(angl >= 360.0)\n    angl = angl - 360.0;\n  if(angl < 180.0) {\n    *ycept = height + 1;\n    \n    if(fabs(angl - 90.0) < SMALL_NUMBER) {\n      *xcept = xcen;\n      return(1);\n    }\n  } else {\n    *ycept = 0.0;\n    \n    if(fabs(angl - 270.0) < SMALL_NUMBER) {\n      *xcept = xcen;\n      return(3);\n    }\n  }\n  \n  angl = (angl / 180.0) * M_PI;\n  \n  slope = tan(angl);\n  \n  *xcept = xcen + ((*ycept - ycen) / slope);\n  if(*xcept < 0) {\n    *ycept = (ycen - (xcen * slope));\n    *xcept = 0.0;\n    return(2);\n  } else if(*xcept > (width + 1)) {\n    *ycept = (ycen + ((width + 1 - xcen) * slope));\n    *xcept = width + 1;\n    return(4);\n  } else {\n    if(*ycept < height)\n      return(3);\n    else\n      return(1);\n  }\n}\n\nvoid _impiei(GFilt g, int qt, int rno, int sno, int flag, int type,\n	     double UNUSED(x), double UNUSED(y),\n	     double xcen, double ycen, double angle1, double angle2)\n{\n  int width, height;		\n  double sweep;			\n  double vx[7], vy[7];		\n  int count;			\n  int intrcpt1, intrcpt2;	\n  double x2, y2;		\n\n  \n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  \n  width = g->x1;\n  height = g->y1;\n  \n  vx[0] = xcen;\n  vy[0] = ycen;\n  sweep = angle2 - angle1;\n  \n  if(fabs(sweep) < SMALL_NUMBER)\n    return;\n  if (sweep < 0.0) sweep = sweep + 360.0;\n  intrcpt1 = pie_intercept(width, height, xcen, ycen, angle1,\n			   &(vx[1]), &(vy[1]));\n  intrcpt2 = pie_intercept(width, height, xcen, ycen, angle2,\n			   &x2, &y2);\n  count = 2;\n  \n  \n  if((intrcpt1 != intrcpt2) || (sweep > 180.0)){\n    do{\n      intrcpt1 = corner_vertex(intrcpt1, width, height,  \n			       &(vx[count]), &(vy[count]));\n      count = count + 1;\n    }while(intrcpt1 != intrcpt2);\n  }\n  vx[count] = x2;\n  vy[count] = y2;\n  count++;\n  _polygoni(g, qt, rno, sno, flag, type, vx, vy, count);\n}\n\n\n\nvoid imannulusi(GFilt g, int rno, int sno, int flag, int type,\n		double x, double y,\n		double xcen, double ycen, double iradius, double oradius)\n{\n  int yy;\n  double dval;\n  double xoff, yoff;\n  Scan *scanlist;\n\n  \n\n  if( iradius == 0 ){\n    imcirclei(g, rno, sno, flag, type, x, y, xcen, ycen, oradius);\n    return;\n  }\n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  iradius /= (double)g->block;\n  oradius /= (double)g->block;\n  \n  if( PIXSTART(ycen - oradius) < PIXSTOP(ycen + oradius) ){\n    g->shapes[sno].ystart = max(g->y0,PIXSTART(ycen - oradius));\n    g->shapes[sno].ystop  = min(g->y1,PIXSTOP(ycen + oradius));\n  }\n  else{\n    g->shapes[sno].ystart = min(g->y1,PIXSTOP(ycen + oradius));\n    g->shapes[sno].ystop  = max(g->y0,PIXSTART(ycen - oradius));\n  }\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  scanlist = g->shapes[sno].scanlist;\n  marky(g, sno, flag, type);\n  \n  for(yy=g->shapes[sno].ystart; yy<=g->shapes[sno].ystop; yy++){\n    yoff = PIXCEN(yy) - ycen;\n    if( (dval=(oradius * oradius) - (yoff * yoff)) > 0.0 ){\n      xoff = sqrt(dval);\n      if( PIXSTART(xcen - xoff) <= PIXSTOP(xcen + xoff)  ){\n	rgs_mark(g, scanlist, sno, flag, type, PIXSTART(xcen - xoff), yy);\n	rgs_mark(g, scanlist, sno, flag, type, PIXSTOP(xcen + xoff), yy);\n      }\n      if( (dval=(iradius * iradius) - (yoff * yoff)) > 0.0 ){\n	xoff = sqrt(dval);\n	if( PIXSTART(xcen - xoff) <= PIXSTOP(xcen + xoff)  ){\n	  rgs_mark(g, scanlist, sno, flag, type, PIXSTART(xcen - xoff), yy);\n	  rgs_mark(g, scanlist, sno, flag, type, PIXSTOP(xcen + xoff), yy);\n	}\n      }\n    }\n  }\n}\n\nint imannulus(GFilt g, int rno, int sno, int flag, int type,\n	      double x, double y,\n	      double xcen, double ycen, double iradius, double oradius)\n{\n  Scan scan;\n\n  if( iradius == 0 ){\n    return(imcircle(g, rno, sno, flag, type, x, y, xcen, ycen, oradius));\n  }\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  scan = g->shapes[sno].scanlist[(int)y];\n  if( (scan							&& \n       ((y >= g->shapes[sno].ystart) && (y <= g->shapes[sno].ystop))	&&\n       (scan->next->next ?\n	(((x >= scan->x) && (x <= scan->next->next->next->x))	&&\n	 !((x >= scan->next->x) && (x <= scan->next->next->x))) :\n	((x >= scan->x) && (x <= scan->next->x))))\n      ==flag ){\n    if( rno && flag ) g->rid = rno;\n    return 1;\n  }\n  else\n    return 0;\n}\n\nvoid imboxi(GFilt g, int rno, int sno, int flag, int type,\n	    double UNUSED(x), double UNUSED(y),\n	    double xcen, double ycen, double xwidth, double yheight,\n	    double angle)\n{\n  double angl;			 \n  double half_width, half_height;\n  double cosangl, sinangl;	 \n  double hw_cos, hw_sin;	 \n  double hh_cos, hh_sin;	 \n  double cornerx[4], cornery[4]; \n\n  \n  if( (xwidth == 0) && (yheight==0) ){\n    return;\n  }\n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  xwidth /= (double)g->block;\n  yheight /= (double)g->block;\n#if USE_ASTRO_ANGLE\n  \n  angl = angle + 90.0;\n#else\n  angl = angle;\n#endif\n  while (angl >= 360.0) angl = angl - 360.0;\n  \n  angl = (angl / 180.0) * M_PI;\n  sinangl = sin (angl);\n  cosangl = cos (angl);\n#if USE_ASTRO_ANGLE\n  \n  \n  \n  half_width = yheight / 2.0;\n  half_height = xwidth / 2.0;\n#else\n  half_width = xwidth / 2.0;\n  half_height = yheight / 2.0;\n#endif\n  hw_cos = half_width * cosangl;\n  hw_sin = half_width * sinangl;\n  hh_cos = half_height * cosangl;\n  hh_sin = half_height * sinangl;\n#if USE_ASTRO_ANGLE\n  cornerx[0] = xcen - hw_cos - hh_sin;\n  cornery[0] = ycen - hw_sin + hh_cos;\n  cornerx[1] = xcen + hw_cos - hh_sin;\n  cornery[1] = ycen + hw_sin + hh_cos;\n  cornerx[2] = xcen + hw_cos + hh_sin;\n  cornery[2] = ycen + hw_sin - hh_cos;\n  cornerx[3] = xcen - hw_cos + hh_sin;\n  cornery[3] = ycen - hw_sin - hh_cos;\n#else\n  cornerx[0] = xcen - hw_cos + hh_sin;\n  cornery[0] = ycen - hh_cos - hw_sin;\n  cornerx[1] = xcen - hw_cos - hh_sin;\n  cornery[1] = ycen + hh_cos - hw_sin;\n  cornerx[2] = xcen + hw_cos - hh_sin;\n  cornery[2] = ycen + hh_cos + hw_sin;\n  cornerx[3] = xcen + hw_cos + hh_sin;\n  cornery[3] = ycen - hh_cos + hw_sin;\n#endif\n  _polygoni(g, 0, rno, sno, flag, type, cornerx, cornery, 4);\n}\n\nint imbox(GFilt g, int rno, int sno, int flag, int type,\n	  double x, double y,\n	  double UNUSED(xcen), double UNUSED(ycen),\n	  double xwidth, double yheight,\n	  double UNUSED(angle))\n{\n  if( (xwidth == 0) && (yheight==0) ){\n    return(!flag);\n  }\n  return impolygon(g, rno, sno, flag, type, x, y);\n}\n\nvoid imcirclei(GFilt g, int UNUSED(rno), int sno, int flag, int type,\n	       double UNUSED(x), double UNUSED(y),\n	       double xcen, double ycen, double radius)\n{\n  int yy;\n  double dval;\n  double xoff, yoff;\n  Scan *scanlist;\n\n  \n  if( radius == 0 ){\n    return;\n  }\n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  radius /= (double)g->block;\n  \n  if( PIXSTART(ycen - radius) < PIXSTOP(ycen + radius) ){\n    g->shapes[sno].ystart = max(g->y0,PIXSTART(ycen - radius));\n    g->shapes[sno].ystop  = min(g->y1,PIXSTOP(ycen + radius));\n  }\n  else{\n    g->shapes[sno].ystart = min(g->y1,PIXSTOP(ycen + radius));\n    g->shapes[sno].ystop  = max(g->y0,PIXSTART(ycen - radius));\n  }\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  scanlist = g->shapes[sno].scanlist;\n  marky(g, sno, flag, type);\n  \n  for(yy=g->shapes[sno].ystart; yy<=g->shapes[sno].ystop; yy++){\n    yoff = PIXCEN(yy) - ycen;\n    if( (dval=(radius * radius) - (yoff * yoff))>=0.0 ){\n      xoff = sqrt(dval);\n      if( PIXSTART(xcen - xoff) <= PIXSTOP(xcen + xoff)  ){\n	rgs_mark(g, scanlist, sno, flag, type, PIXSTART(xcen - xoff), yy);\n	rgs_mark(g, scanlist, sno, flag, type, PIXSTOP(xcen + xoff), yy);\n      }\n    }\n  }\n}\n\nint imcircle(GFilt g, int rno, int sno, int flag, int UNUSED(type),\n	     double x, double y,\n	     double UNUSED(xcen), double UNUSED(ycen), double radius)\n{\n  Scan scan;\n\n  if( radius == 0 ){\n    return(!flag);\n  }\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  scan = g->shapes[sno].scanlist[(int)y];\n  if( (scan 							&& \n      ((y>=g->shapes[sno].ystart) && (y<=g->shapes[sno].ystop))	&&\n      ((x >= scan->x) && (x <= (scan->next)->x))) == flag ){\n    if( rno && flag ) g->rid = rno;\n    return 1;\n  }\n  else{\n   return 0;\n  }\n}\n\nvoid imellipsei(GFilt g, int rno, int sno, int flag, int type,\n		double x, double y,\n		double xcen, double ycen, double xrad, double yrad,\n		double angle)\n{\n  int yy;\n  int nr, nc;\n  double yhi;\n  double yoff;\n  double xboff, xfoff;\n  double angl;\n  double sinangl, cosangl;\n  double cossq, sinsq, xradsq, yradsq;\n  double a, b_partial, c_partial;\n  double b, c;\n  Scan *scanlist;\n\n  \n\n  if( xrad == yrad ){\n    imcirclei(g, rno, sno, flag, type, x, y, xcen, ycen, xrad);\n    return;\n  }\n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  xrad /= (double)g->block;\n  yrad /= (double)g->block;\n  \n#if USE_ASTRO_ANGLE\n  \n  angl = angle + 90.0;\n#else\n  angl = angle;\n#endif\n  while( angl >= 360.0 )\n    angl = angl - 360.0;\n  \n  angl = (angl / 180.0) * M_PI;\n  sinangl = sin(angl);\n  cosangl = cos(angl);\n  \n  \n  FPU_DOUBLE\n#if USE_ASTRO_ANGLE\n  yhi = fabs(sinangl * yrad) + fabs(cosangl * xrad);\n#else\n  yhi = fabs(sinangl * xrad) + fabs(cosangl * yrad);\n#endif\n  yhi = min(yhi, max(yrad, xrad));\n  \n  if( PIXSTART(ycen - yhi) < PIXSTOP(ycen + yhi) ){\n    g->shapes[sno].ystart = max(g->y0,PIXSTART(ycen - yhi));\n    g->shapes[sno].ystop  = min(g->y1,PIXSTOP(ycen + yhi));\n  }\n  else{\n    g->shapes[sno].ystart = min(g->y1,PIXSTOP(ycen + yhi));\n    g->shapes[sno].ystop  = max(g->y0,PIXSTART(ycen - yhi));\n  }\n  FPU_RESTORE\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  scanlist = g->shapes[sno].scanlist;\n  marky(g, sno, flag, type);\n  \n  cossq = cosangl * cosangl;\n  sinsq = sinangl * sinangl;\n#if USE_ASTRO_ANGLE\n  \n  \n  \n  \n  xradsq = yrad * yrad;\n  yradsq = xrad * xrad;\n#else\n  xradsq = xrad * xrad;\n  yradsq = yrad * yrad;\n#endif\n  \n  a = (cossq / xradsq) + (sinsq / yradsq);\n  b_partial = (2.0 * sinangl) * ((cosangl / xradsq) - (cosangl / yradsq));\n  c_partial = (sinsq / xradsq) + (cossq / yradsq);\n  \n  for(yy=g->shapes[sno].ystart; yy<=g->shapes[sno].ystop; yy++){\n    yoff = yy - ycen;\n    b = b_partial * yoff;\n    c = (c_partial * yoff * yoff) - 1.0;\n    \n    quadeq (a, b, c, &xboff, &xfoff, &nr, &nc);\n    \n    if( nr != 0 ) {\n      \n      rgs_mark(g, scanlist, sno, flag, type, PIXSTART(xcen + xboff), yy);\n      rgs_mark(g, scanlist, sno, flag, type, PIXSTOP(xcen + xfoff), yy);\n    }\n  }\n}\n\nint imellipse(GFilt g, int rno, int sno, int flag, int type,\n	      double x, double y,\n	      double xcen, double ycen, double xrad, double yrad,\n	      double UNUSED(angle))\n{\n  Scan scan;\n\n  if( xrad == yrad ){\n    return(imcircle(g, rno, sno, flag, type, x, y, xcen, ycen, xrad));\n  }\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  scan = g->shapes[sno].scanlist[(int)y];\n  FPU_DOUBLE\n  if( (scan 							&& \n      ((y>=g->shapes[sno].ystart) && (y<=g->shapes[sno].ystop))	&&\n      ((x >= scan->x) && (x <= (scan->next)->x))) == flag ){\n    if( rno && flag ) g->rid = rno;\n    FPU_RESTORE\n    return 1;\n  }\n  else{\n    FPU_RESTORE\n    return 0;\n  }\n}\n\nvoid imfieldi(GFilt g, int UNUSED(rno), int sno, int flag, int type,\n	      double UNUSED(x), double UNUSED(y))\n{\n  int yy;\n  Scan *scanlist;\n\n  \n  g->shapes[sno].ystart = g->y0;\n  g->shapes[sno].ystop = g->y1;\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  scanlist = g->shapes[sno].scanlist;\n  marky(g, sno, flag, type);\n  \n  for(yy=g->shapes[sno].ystart; yy<=g->shapes[sno].ystop; yy++){\n    rgs_mark(g, scanlist, sno, flag, type, g->x0, yy);\n    rgs_mark(g, scanlist, sno, flag, type, g->x1, yy);\n  }\n}\n\nint imfield(GFilt g, int rno, int UNUSED(sno), int flag, int UNUSED(type),\n	    double UNUSED(x), double UNUSED(y))\n{\n  if( flag ){\n    if( rno && flag ) g->rid = rno;\n    return 1;\n  }\n  else{\n    return 0;\n  }\n}\n\nvoid imlinei(GFilt g, int UNUSED(rno), int sno, int flag, int type,\n	     double UNUSED(x), double UNUSED(y),\n	     double x1, double y1, double x2, double y2)\n{\n  double vx[2];\n  double vy[2];\n  int xval, yval;\n  double invslope, xoffset;\n\n  \n  \n  x1 = (x1 - g->xmin)/g->block + 1.0;\n  y1 = (y1 - g->ymin)/g->block + 1.0;\n  x2 = (x2 - g->xmin)/g->block + 1.0;\n  y2 = (y2 - g->ymin)/g->block + 1.0;\n  \n  if( y1 < y2 ){\n    vx[0] = x1;\n    vy[0] = y1;\n    vx[1] = x2;\n    vy[1] = y2;\n  }\n  else{\n    vx[0] = x2;\n    vy[0] = y2;\n    vx[1] = x1;\n    vy[1] = y1;\n  }\n  \n  g->shapes[sno].ystart = PIXNUM(vy[0]);\n  g->shapes[sno].ystop = PIXNUM(vy[1]);\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  marky(g, sno, flag, type);\n  if( feq(vy[0],vy[1]) ){\n    rgs_mark(g, g->shapes[sno].scanlist, sno, flag, type,\n	     (int)min(vx[0],vx[1]), (int)vy[0]);\n    rgs_mark(g, g->shapes[sno].scanlist, sno, flag, type,\n	     (int)max(vx[0],vx[1]), (int)vy[0]);\n  }\n  else{\n    \n    invslope = (vx[0] - vx[1]) / (vy[0] - vy[1]);\n    xoffset = vx[0];\n    for(yval=vy[0]; yval<=vy[1]; yval++){\n      xval = xoffset;\n      rgs_mark(g, g->shapes[sno].scanlist, sno, flag, type, xval, yval);\n      xoffset = xoffset + invslope;\n    }\n  }\n}\n\nint imline(GFilt g, int rno, int sno, int flag, int UNUSED(type),\n	   double x, double y,\n	   double UNUSED(x1), double UNUSED(y1),\n	   double UNUSED(x2), double UNUSED(y2))\n{\n  Scan scan;\n\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  scan = g->shapes[sno].scanlist[(int)y];\n  if( (scan                                     &&\n      ((x==(int)scan->x) 			||\n       (scan->next&&((x>=(int)scan->x)&&(x<=(int)scan->next->x)))))==flag ) {\n    if( rno && flag ) g->rid = rno;\n    return 1;\n  }\n  else\n    return 0;\n}\n\nvoid impiei(GFilt g, int rno, int sno, int flag, int type,\n	    double x, double y,\n	    double xcen, double ycen, double angle1, double angle2)\n{\n  _impiei(g, 0, rno, sno, flag, type, x, y, xcen, ycen, angle1, angle2);\n}\n\nint impie(GFilt g, int rno, int sno, int flag, int type,\n	  double x, double y,\n	  double UNUSED(xcen), double UNUSED(ycen),\n	  double UNUSED(angle1), double UNUSED(angle2))\n{\n  return impolygon(g, rno, sno, flag, type, x, y);\n}\n\nvoid imqtpiei(GFilt g, int rno, int sno, int flag, int type,\n	      double x, double y,\n	      double xcen, double ycen, double angle1, double angle2)\n{\n  _impiei(g, 1, rno, sno, flag, type, x, y, xcen, ycen, angle1, angle2);\n}\n\nint imqtpie(GFilt g, int rno, int sno, int flag, int type,\n	    double x, double y,\n	    double UNUSED(xcen), double UNUSED(ycen),\n	    double UNUSED(angle1), double UNUSED(angle2))\n{\n  return impolygon(g, rno, sno, flag, type, x, y);\n}\n\nvoid impointi(GFilt g, int UNUSED(rno), int sno, int flag, int type,\n	      double UNUSED(x), double UNUSED(y),\n	      double xcen, double ycen)\n{\n  \n  \n  xcen = (xcen - g->xmin)/g->block + 1.0;\n  ycen = (ycen - g->ymin)/g->block + 1.0;\n  \n  g->shapes[sno].ystart = PIXNUM(ycen);\n  g->shapes[sno].ystop = PIXNUM(ycen);\n  g->shapes[sno].scanlist = (Scan *)calloc(g->y1+1, sizeof(Scan));\n  marky(g, sno, flag, type);\n  rgs_mark(g, g->shapes[sno].scanlist, sno, flag, type,\n	   PIXNUM(xcen), PIXNUM(ycen));\n}\n\nint impoint(GFilt g, int rno, int sno, int flag, int UNUSED(type),\n	    double x, double y,\n	    double UNUSED(xcen), double UNUSED(ycen))\n{\n  Scan scan;\n\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  scan = g->shapes[sno].scanlist[(int)y];\n  if( (scan                                     &&\n      (y == (int)g->shapes[sno].ystart)		&& \n      (x == (int)scan->x))==flag )		{\n    if( rno && flag ) g->rid = rno;\n    return 1;\n  }\n  else\n    return 0;\n}\n\n#ifdef __STDC__\nvoid\nimpolygoni(GFilt g, int rno, int sno, int flag, int type,\n	   double UNUSED(x), double y, ...)\n{\n  double *vx=NULL, *vy=NULL;\n  int count, maxcount;\n  va_list args;\n  va_start(args, y);\n#else\nvoid impolygoni(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double *vx=NULL, *vy=NULL;\n  int count, maxcount;\n  va_list args;\n  va_start(args);\n  g  =  va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x = va_arg(args, double);\n  y = va_arg(args, double);\n#endif\n  \n  \n  maxcount = MASKINC;\n  vx = (double *)calloc(maxcount, sizeof(double));\n  vy = (double *)calloc(maxcount, sizeof(double));\n  \n  count = 0;\n  while( 1 ){\n    if( (count + 1) >= maxcount ){\n      maxcount += MASKINC;\n      vx = (double *)realloc(vx, maxcount*sizeof(double));\n      vy = (double *)realloc(vy, maxcount*sizeof(double));\n    }\n    vx[count] = va_arg(args, double);\n    vy[count] = va_arg(args, double);\n    if( feq(vx[count],PSTOP) && feq(vy[count],PSTOP) )\n      break;\n    vx[count] = (vx[count] - g->xmin)/g->block + 1.0;\n    vy[count] = (vy[count] - g->ymin)/g->block + 1.0;\n    count ++;\n  }\n  va_end(args);\n  \n  vx = (double *)realloc(vx, count*sizeof(double));\n  vy = (double *)realloc(vy, count*sizeof(double));\n  \n  _polygoni(g, 0, rno, sno, flag, type, vx, vy, count);\n  if( vx ) free(vx);\n  if( vy ) free(vy);\n}\n\n#ifdef __STDC__\nint\nimpolygon(GFilt g, int rno, int sno, int flag, int UNUSED(type),\n	  double x, double y, ...)\n{\n  int crossings;\n  Scan scan;\n  va_list args;\n  va_start(args, y);\n#else\nint impolygon(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  int crossings;\n  Scan scan;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x = va_arg(args, double);\n  y = va_arg(args, double);\n#endif\n  va_end(args);\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  if( flag ){\n    if(y < g->shapes[sno].ystart ) return 0;\n    if(y > g->shapes[sno].ystop ) return 0;\n  }\n  \n  if( (y>=g->shapes[sno].ystart)&&(y<=g->shapes[sno].ystop) ){\n    crossings = 0;\n    for(scan=g->shapes[sno].scanlist[(int)y]; scan; scan=scan->next){\n      if( x >= scan->x )\n	crossings++;\n      else\n	break;\n    }\n    if( (crossings%2) == flag ){\n      if( rno && flag ) g->rid = rno;\n      return 1;\n    }\n    else{\n      return 0;\n    }\n  }\n  else{\n    return !flag;\n  }\n}\n\nvoid imnannulusi(GFilt g, int rno, int sno, int flag, int type,\n		 double x, double y,\n		 double xcen, double ycen,\n		 double lo, double hi, int n)\n{\n  int i;\n  int xsno;\n  double dinc;\n\n  \n\n  \n  dinc = (hi - lo)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imannulusi(g, 0, xsno, flag, type, x, y, xcen, ycen, lo, hi);\n  for(i=0; i<n; i++){\n    imannulusi(g, rno+i, sno+i, flag, type, x, y,\n	       xcen, ycen, lo+(i*dinc), lo+((i+1)*dinc));\n  }\n}\n\nvoid imnboxi(GFilt g, int rno, int sno, int flag, int type,\n	     double x, double y,\n	     double xcen, double ycen,\n	     double lox, double loy, double hix, double hiy, int n,\n	     double ang)\n{\n  int i;\n  int xsno;\n  double dincx;\n  double dincy;\n\n  \n\n  \n  dincx = (hix - lox)/n;\n  dincy = (hiy - loy)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imboxi(g, 0, xsno, flag, type, x, y, xcen, ycen, hix, hiy, ang);\n  imboxi(g, 0, xsno+1, flag, type, x, y, xcen, ycen, lox, loy, ang);\n  for(i=0; i<n; i++){\n    imboxi(g, rno+i, sno+i, flag, type, x, y,\n	   xcen, ycen, lox+((i+1)*dincx), loy+((i+1)*dincy), ang);\n  }\n}\n\nvoid imnellipsei(GFilt g, int rno, int sno, int flag, int type,\n		 double x, double y,\n		 double xcen, double ycen,\n		 double lox, double loy, double hix, double hiy, int n,\n		 double ang)\n{\n  int i;\n  int xsno;\n  double dincx;\n  double dincy;\n\n  \n\n  \n  dincx = (hix - lox)/n;\n  dincy = (hiy - loy)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imellipsei(g, 0, xsno, flag, type, x, y, xcen, ycen, hix, hiy, ang);\n  imellipsei(g, 0, xsno+1, flag, type, x, y, xcen, ycen, lox, loy, ang);\n  for(i=0; i<n; i++){\n    imellipsei(g, rno+i, sno+i, flag, type, x, y,\n	       xcen, ycen, lox+((i+1)*dincx), loy+((i+1)*dincy), ang);\n  }\n}\n\nvoid imnpiei(GFilt g, int rno, int sno, int flag, int type,\n	     double x, double y,\n	     double xcen, double ycen,\n	     double lo, double hi, int n)\n{\n  int i;\n  int xsno;\n  double dinc;\n\n  \n\n  \n  while( lo > hi ) lo -= 360.0;\n  dinc = (hi - lo)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  impiei(g, 0, xsno, flag, type, x, y, xcen, ycen, lo, hi);\n  for(i=0; i<n; i++){\n    impiei(g, rno+i, sno+i, flag, type, x, y,\n	   xcen, ycen, lo+(i*dinc), lo+((i+1)*dinc));\n  }\n}\n\nvoid impandai(GFilt g, int rno, int sno, int flag, int type,\n	      double x, double y,\n	      double xcen, double ycen,\n	      double anglo, double anghi, double angn,\n	      double radlo, double radhi, double radn)\n{\n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, rinc;\n\n  \n\n  \n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  rinc = (radhi - radlo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imannulusi(g, 0, xsno, flag, type, x, y, xcen, ycen, radlo, radhi);\n  imqtpiei(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi);\n  for(a=0; a<ahi; a++){\n    for(r=0; r<rhi; r++){\n      imannulusi(g, rno+n, sno+(2*n), flag, type, x, y,\n		 xcen, ycen, radlo+(r*rinc), radlo+((r+1)*rinc));\n      imqtpiei(g, rno+n, sno+(2*n+1), flag, type, x, y,\n	       xcen, ycen, anglo+(a*ainc), anglo+((a+1)*ainc));\n      n++;\n    }\n  }\n}\n\nvoid imbpandai(GFilt g, int rno, int sno, int flag, int type,\n	       double x, double y,\n	       double xcen, double ycen,\n	       double anglo, double anghi, double angn,\n	       double xlo, double ylo, double xhi, double yhi, double radn,\n	       double ang)\n{\n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, xinc, yinc;\n\n  \n\n  \n  anglo += ang;\n  anghi += ang;\n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  xinc = (xhi - xlo)/radn;\n  yinc = (yhi - ylo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imboxi(g, 0, xsno, flag, type, x, y, xcen, ycen, xhi, yhi, ang);\n  imqtpiei(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi);\n  imboxi(g, 0, xsno+2, flag, type, x, y, xcen, ycen, xlo, ylo, ang);\n  for(a=0; a<ahi; a++){\n    for(r=1; r<=rhi; r++){\n      imboxi(g, rno+n, sno+(2*n), flag, type, x, y,\n	     xcen, ycen, xlo+(r*xinc), ylo+(r*yinc), ang);\n      imqtpiei(g, rno+n, sno+(2*n+1), flag, type, x, y,\n	       xcen, ycen, anglo+(a*ainc), anglo+((a+1)*ainc));\n      n++;\n    }\n  }\n}\n\nvoid imepandai(GFilt g, int rno, int sno, int flag, int type,\n	       double x, double y,\n	       double xcen, double ycen,\n	       double anglo, double anghi, double angn,\n	       double xlo, double ylo, double xhi, double yhi, double radn,\n	       double ang)\n{\n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, xinc, yinc;\n\n  \n\n  \n  anglo += ang;\n  anghi += ang;\n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  xinc = (xhi - xlo)/radn;\n  yinc = (yhi - ylo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  \n  imellipsei(g, 0, xsno, flag, type, x, y, xcen, ycen, xhi, yhi, ang);\n  imqtpiei(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi);\n  imellipsei(g, 0, xsno+2, flag, type, x, y, xcen, ycen, xlo, ylo, ang);\n  for(a=0; a<ahi; a++){\n    for(r=1; r<=rhi; r++){\n      imellipsei(g, rno+n, sno+(2*n), flag, type, x, y,\n		 xcen, ycen, xlo+(r*xinc), ylo+(r*yinc), ang);\n      imqtpiei(g, rno+n, sno+(2*n+1), flag, type, x, y,\n	       xcen, ycen, anglo+(a*ainc), anglo+((a+1)*ainc));\n      n++;\n    }\n  }\n}\n\n\n\nint imnannulus(GFilt g, int rno, int sno, int flag, int type,\n	       double x, double y, double xcen, double ycen,\n	       double lo, double hi, int n)\n{\n  int i;\n  int xsno;\n  double dinc;\n\n  \n  dinc = (hi - lo)/(double)n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imannulus(g, 0, xsno, flag, type, x, y, xcen, ycen, lo, hi) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( imannulus(g, rno+i, sno+i, flag, type, x, y,\n		    xcen, ycen, lo+(i*dinc), lo+((i+1)*dinc)) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  else{\n    \n    if( !imannulus(g, 0, xsno, 1, type, x, y, xcen, ycen, lo, hi) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\nint imnbox(GFilt g, int rno, int sno, int flag, int type,\n	   double x, double y, double xcen, double ycen,\n	   double lox, double loy, double hix, double hiy, int n,\n	   double ang)\n{\n  int i;\n  int xsno;\n  double dincx;\n  double dincy;\n\n  \n  dincx = (hix - lox)/n;\n  dincy = (hiy - loy)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imbox(g, 0, xsno, flag, type, x, y, xcen, ycen, hix, hiy, ang) ){\n      return(0);\n    }\n    \n    if( imbox(g, 0, xsno+1, flag, type, x, y, xcen, ycen, lox, loy, ang) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( imbox(g, rno+i, sno+i, flag, type, x, y,\n		xcen, ycen, lox+((i+1)*dincx), loy+((i+1)*dincy), ang) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !imbox(g, 0, xsno, 1, type, x, y, xcen, ycen, hix, hiy, ang) ){\n      return(1);\n    }\n    \n    if( imbox(g, 0, xsno+1, 1, type, x, y, xcen, ycen, lox, loy, ang) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\nint imnellipse(GFilt g, int rno, int sno, int flag, int type,\n	       double x, double y, double xcen, double ycen,\n	       double lox, double loy, double hix, double hiy, int n,\n	       double ang)\n{\n  int i;\n  int xsno;\n  double dincx;\n  double dincy;\n\n  \n  dincx = (hix - lox)/n;\n  dincy = (hiy - loy)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imellipse(g, 0, xsno, flag, type, x, y, xcen, ycen, hix, hiy, ang) ){\n      return(0);\n    }\n    \n    if( imellipse(g, 0, xsno+1, flag, type, x, y, xcen, ycen, lox, loy, ang) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( imellipse(g, rno+i, sno+i, flag, type, x, y,\n		    xcen, ycen, lox+((i+1)*dincx), loy+((i+1)*dincy), ang) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !imellipse(g, 0, xsno, 1, type, x, y, xcen, ycen, hix, hiy, ang) ){\n      return(1);\n    }\n    \n    if( imellipse(g, 0, xsno+1, 1, type, x, y, xcen, ycen, lox, loy, ang) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\nint imnpie(GFilt g, int rno, int sno, int flag, int type,\n	   double x, double y, double xcen, double ycen,\n	   double lo, double hi, int n)\n{\n  int i;\n  int xsno;\n  double dinc;\n\n  \n  while( lo > hi ) lo -= 360.0;\n  dinc = (hi - lo)/n;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !impie(g, 0, xsno, flag, type, x, y, xcen, ycen, lo, hi) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( impie(g, rno+i, sno+i, flag, type, x, y,\n		xcen, ycen, lo+(i*dinc), lo+((i+1)*dinc)) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  else{\n    \n    if( !impie(g, 0, xsno, 1, type, x, y, xcen, ycen, lo, hi) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\nint impanda(GFilt g, int rno, int sno, int flag, int type,\n	    double x, double y,\n	    double xcen, double ycen,\n	    double anglo, double anghi, double angn,\n	    double radlo, double radhi, double radn)\n{\n    \n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, rinc;\n\n  \n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  rinc = (radhi - radlo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imannulus(g, 0, xsno, flag, type, x, y, xcen, ycen, radlo, radhi) ||\n	!impie(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi)   ){\n      return(0);\n    }\n    \n    for(a=1; a<=ahi; a++){\n      for(r=1; r<=rhi; r++){\n	if( imannulus(g, rno+n, sno+(2*n), flag, type, x, y,\n		      xcen, ycen, radlo+((r-1)*rinc), radlo+(r*rinc)) &&\n	    impie(g, rno+n, sno+(2*n+1), flag, type, x, y,\n		  xcen, ycen, anglo+((a-1)*ainc), anglo+(a*ainc))     ){\n	  return(1);\n	}\n	n++;\n      }\n    }\n    return(0);\n  }\n  else{\n    \n    if( !imannulus(g, 0, xsno, 1, type, x, y, xcen, ycen, radlo, radhi) )\n      return(1);\n    else if( !impie(g, 0, xsno+1, 1, type, x, y, xcen, ycen, anglo, anghi) ){\n      return(1);\n    }\n    else{\n      return(0);\n    }\n  }\n}\n\nint imbpanda(GFilt g, int rno, int sno, int flag, int type,\n	     double x, double y,\n	     double xcen, double ycen,\n	     double anglo, double anghi, double angn,\n	     double xlo, double ylo, double xhi, double yhi, double radn,\n	     double ang)\n{\n    \n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, xinc, yinc;\n\n  \n  anglo += ang;\n  anghi += ang;\n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  xinc = (xhi - xlo)/radn;\n  yinc = (yhi - ylo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imbox(g, 0, xsno, flag, type, x, y, xcen, ycen, xhi, yhi,\n	       ang) ){\n      return(0);\n    }\n    \n    else if( imbox(g, 0, xsno+2, flag, type, x, y, xcen, ycen, xlo, ylo,\n		   ang) ){\n      return(0);\n    }\n    \n    else if( !impie(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi)){\n      return(0);\n    }\n    \n    for(a=0; a<ahi; a++){\n      for(r=1; r<=rhi; r++){\n	if( imbox(g, rno+n, sno+(2*n), flag, type, x, y,\n		  xcen, ycen, xlo+(r*xinc), ylo+(r*yinc), ang)   &&\n	    imqtpie(g, rno+n, sno+(2*n+1), flag, type, x, y,\n		    xcen, ycen, anglo+(a*ainc), anglo+((a+1)*ainc))  ){\n	  return(1);\n	}\n	n++;\n      }\n    }\n    return(0);\n  }\n  else{\n    \n    if( !imbox(g, 0, xsno, 1, type, x, y, xcen, ycen, xhi, yhi, ang) )\n      return(1);\n    \n    else if( !imbox(g, 0, xsno+2, 1, type, x, y, xcen, ycen, xlo, ylo,\n		    ang) )\n      return(1);\n    \n    else if( !impie(g, 0, xsno+1, 1, type, x, y, xcen, ycen, anglo, anghi) ){\n      return(1);\n    }\n    else{\n      return(0);\n    }\n  }\n}\n\nint imepanda(GFilt g, int rno, int sno, int flag, int type,\n	     double x, double y,\n	     double xcen, double ycen,\n	     double anglo, double anghi, double angn,\n	     double xlo, double ylo, double xhi, double yhi, double radn,\n	     double ang)\n{\n    \n  int a, r;\n  int ahi, rhi;\n  int xsno;\n  int n=0;\n  double ainc, xinc, yinc;\n\n  \n  anglo += ang;\n  anghi += ang;\n  while( anglo > anghi ) anglo -= 360.0;\n  ainc = (anghi - anglo)/angn;\n  ahi = (int)angn;\n  xinc = (xhi - xlo)/radn;\n  yinc = (yhi - ylo)/radn;\n  rhi = (int)radn;\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( flag ){\n    \n    if( !imellipse(g, 0, xsno, flag, type, x, y, xcen, ycen, xhi, yhi,\n		   ang) ){\n      return(0);\n    }\n    \n    else if( imellipse(g, 0, xsno+2, flag, type, x, y, xcen, ycen, xlo, ylo,\n		       ang) ){\n      return(0);\n    }\n    \n    else if( !impie(g, 0, xsno+1, flag, type, x, y, xcen, ycen, anglo, anghi)){\n      return(0);\n    }\n    \n    for(a=0; a<ahi; a++){\n      for(r=1; r<=rhi; r++){\n	if( imellipse(g, rno+n, sno+(2*n), flag, type, x, y,\n		      xcen, ycen, xlo+(r*xinc), ylo+(r*yinc), ang)   &&\n	    imqtpie(g, rno+n, sno+(2*n+1), flag, type, x, y,\n		    xcen, ycen, anglo+(a*ainc), anglo+((a+1)*ainc))  ){\n	  return(1);\n	}\n	n++;\n      }\n    }\n    return(0);\n  }\n  else{\n    \n    if( !imellipse(g, 0, xsno, 1, type, x, y, xcen, ycen, xhi, yhi, ang) )\n      return(1);\n    \n    else if( !imellipse(g, 0, xsno+2, 1, type, x, y, xcen, ycen, xlo, ylo,\n			ang) )\n      return(1);\n    \n    else if( !impie(g, 0, xsno+1, 1, type, x, y, xcen, ycen, anglo, anghi) ){\n      return(1);\n    }\n    else{\n      return(0);\n    }\n  }\n}\n\n#ifdef __STDC__\nvoid\nimvannulusi(GFilt g, int rno, int sno, int flag, int type,\n	    double x, double y, double xcen, double ycen, ...)\n{\n  int i, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvannulusi(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  \n  if( n == 2 ){\n    imannulusi(g, rno, sno, flag, type, x, y, xcen, ycen, xv[0], xv[1]);\n    return;\n  }\n  imannulusi(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[0], xv[n-1]);\n  for(i=0; i<(n-1); i++){\n    imannulusi(g, rno+i, sno+i, flag, type, x, y, xcen, ycen, xv[i], xv[i+1]);\n  }\n}\n\n#ifdef __STDC__\nvoid\nimvboxi(GFilt g, int rno, int sno, int flag, int type,\n	double x, double y, double xcen, double ycen, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double ang;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvboxi(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  ang = xv[--n];\n  \n  if( n == 2 ){\n    imboxi(g, rno, sno, flag, type, x, y, xcen, ycen, xv[0], xv[1], ang);\n    return;\n  }\n  imboxi(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[n-2], xv[n-1], ang);\n  imboxi(g, 0, xsno+1, flag, type, x, y, xcen, ycen, xv[0], xv[1], ang);\n  for(i=2, j=0; i<n; i+=2, j++){\n    imboxi(g, rno+j, sno+j, flag, type, x, y, xcen, ycen, xv[i], xv[i+1], ang);\n  }\n}\n\n#ifdef __STDC__\nvoid\nimvellipsei(GFilt g, int rno, int sno, int flag, int type,\n	    double x, double y, double xcen, double ycen, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double ang;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvellipsei(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  ang = xv[--n];\n  \n  if( n == 2 ){\n    imellipsei(g, rno, sno, flag, type, x, y, xcen, ycen, xv[0], xv[1], ang);\n    return;\n  }\n  imellipsei(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[n-2], xv[n-1], ang);\n  imellipsei(g, 0, xsno+1, flag, type, x, y, xcen, ycen, xv[0], xv[1], ang);\n  for(i=2, j=0; i<n; i+=2, j++){\n    imellipsei(g, rno+j, sno+j, flag, type, x, y, xcen, ycen,\n	       xv[i], xv[i+1], ang);\n  }\n}\n\n#ifdef __STDC__\nvoid\nimvpiei(GFilt g, int rno, int sno, int flag, int type,\n	double x, double y, double xcen, double ycen, ...)\n{\n  int i, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvpiei(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  \n  if( n == 2 ){\n    impiei(g, rno, sno, flag, type, x, y, xcen, ycen, xv[0], xv[1]);\n    return;\n  }\n  impiei(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[0], xv[n-1]);\n  for(i=0; i<(n-1); i++){\n    impiei(g, rno+i, sno+i, flag, type, x, y, xcen, ycen, xv[i], xv[i+1]);\n  }\n}\n\n#ifdef __STDC__\nvoid\nimvpointi(GFilt g, int rno, int sno, int flag, int type,\n	  double x, double y, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, y);\n#else\nint imvpointi(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  for(i=0, j=0; i<n; i+=2, j++){\n    impointi(g, rno+j, sno+j, flag, type, x, y, xv[i], xv[i+1]);\n  }\n}\n\n\n\n#ifdef __STDC__\nint\nimvannulus(GFilt g, int rno, int sno, int flag, int type,\n	   double x, double y, double xcen, double ycen, ...)\n{\n  int i, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvannulus(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double *xv;\n  int i, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  \n  if( n == 2 ){\n    return(imannulus(g, rno, sno, flag, type, x, y, xcen, ycen, xv[0], xv[1]));\n  }\n  if( flag ){\n    \n    if( !imannulus(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[0], xv[n-1]) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( imannulus(g, rno+i, sno+i, flag, type, x, y, xcen, ycen,\n		  xv[i], xv[i+1]) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !imannulus(g, 0, xsno, 1, type, x, y, xcen, ycen, xv[0], xv[n-1]) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\n#ifdef __STDC__\nint\nimvbox(GFilt g, int rno, int sno, int flag, int type,\n       double x, double y, double xcen, double ycen, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double ang;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvbox(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  ang = xv[--n];\n  \n  if( n == 2 ){\n    return(imbox(g, rno, sno, flag, type, x, y,\n		 xcen, ycen, xv[0], xv[1], ang));\n  }\n  if( flag ){\n    \n    if( !imbox(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[n-2], xv[n-1],\n	       ang) ){\n      return(0);\n    }\n    \n    if( imbox(g, 0, xsno+1, flag, type, x, y, xcen, ycen, xv[0], xv[1], ang) ){\n      return(0);\n    }\n    \n    for(i=2, j=0; i<n; i+=2, j++){\n      if( imbox(g, rno+j, sno+j, flag, type, x, y, xcen, ycen,\n	      xv[i], xv[i+1], ang) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !imbox(g, 0, xsno, 1, type, x, y, xcen, ycen, xv[n-2], xv[n-1], ang) ){\n      return(1);\n    }\n    \n    else if( imbox(g, 0, xsno+1, 1, type, x, y, xcen, ycen, xv[0], xv[1], ang) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\n#ifdef __STDC__\nint\nimvellipse(GFilt g, int rno, int sno, int flag, int type,\n	   double x, double y, double xcen, double ycen, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double ang;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvellipse(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double ang;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  ang = xv[--n];\n  \n  if( n == 2 ){\n    return(imellipse(g, rno, sno, flag, type, x, y,\n		     xcen, ycen, xv[0], xv[1], ang));\n  }\n  if( flag ){\n    \n    if( !imellipse(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[n-2], xv[n-1],\n		   ang) ){\n      return(0);\n    }\n    \n    if( imellipse(g, 0, xsno+1, flag, type, x, y, xcen, ycen, xv[0], xv[1],\n		  ang) ){\n      return(0);\n    }\n    \n    for(i=2, j=0; i<n; i+=2, j++){\n      if( imellipse(g, rno+j, sno+j, flag, type, x, y, xcen, ycen,\n		  xv[i], xv[i+1], ang) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !imellipse(g, 0, xsno, 1, type, x, y, xcen, ycen, xv[n-2], xv[n-1],\n		   ang) ){\n      return(1);\n    }\n    \n    if( imellipse(g, 0, xsno+1, 1, type, x, y, xcen, ycen, xv[0], xv[1], ang) ){\n      return(1);\n    }\n    return(0);\n  }\n}\n\n#ifdef __STDC__\nint\nimvpie(GFilt g, int rno, int sno, int flag, int type,\n       double x, double y, double xcen, double ycen, ...)\n{\n  int i, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, ycen);\n#else\nint imvpie(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double xcen, ycen;\n  double *xv;\n  int i, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n  xcen  = va_arg(args, double);\n  ycen  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  \n  if( n == 2 ){\n    return(impie(g, rno, sno, flag, type, x, y,\n		 xcen, ycen, xv[0], xv[1]));\n  }\n  if( flag ){\n    \n    if( !impie(g, 0, xsno, flag, type, x, y, xcen, ycen, xv[0], xv[n-1]) ){\n      return(0);\n    }\n    \n    for(i=0; i<n; i++){\n      if( impie(g, rno+i, sno+i, flag, type, x, y, xcen, ycen, xv[i], xv[i+1]) ){\n	return(1);\n      }\n    }\n    return(0);\n  }\n  \n  else{\n    \n    if( !impie(g, 0, xsno, 1, type, x, y, xcen, ycen, xv[0], xv[n-1]) ){\n      return(1);\n    }\n    return(1);\n  }\n}\n\n#ifdef __STDC__\nint\nimvpoint(GFilt g, int rno, int sno, int flag, int type,\n	 double x, double y, ...)\n{\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  double *xv;\n  va_list args;\n  va_start(args, y);\n#else\nint imvpoint(va_alist) va_dcl\n{\n  GFilt g;\n  int rno, sno, flag, type;\n  double x, y;\n  double *xv;\n  int i, j, n;\n  int maxpts;\n  int xsno;\n  va_list args;\n  va_start(args);\n  g  = va_arg(args, GFilt);\n  rno  = va_arg(args, int);\n  sno  = va_arg(args, int);\n  flag  = va_arg(args, int);\n  type  = va_arg(args, int);\n  x  = va_arg(args, double);\n  y  = va_arg(args, double);\n#endif\n  xsno = (g->nshapes+1)+((sno-1)*XSNO);\n  if( !g->shapes[xsno].xv ){\n    maxpts = MASKINC;\n    g->shapes[xsno].xv = (double *)calloc(maxpts, sizeof(double));\n    g->shapes[xsno].nv = 0;\n    while( 1 ){\n      if( g->shapes[xsno].nv >= maxpts ){\n        maxpts += MASKINC;\n        g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  					      maxpts*sizeof(double));\n      }\n      g->shapes[xsno].xv[g->shapes[xsno].nv] = va_arg(args, double);\n      if( feq(g->shapes[xsno].xv[g->shapes[xsno].nv],PSTOP)   &&\n  	feq(g->shapes[xsno].xv[g->shapes[xsno].nv-1],PSTOP) ){\n        g->shapes[xsno].nv--;\n        break;\n      }\n      g->shapes[xsno].nv++;\n    }\n    va_end(args);\n    g->shapes[xsno].xv = (double *)realloc(g->shapes[xsno].xv,\n  	  	        g->shapes[xsno].nv*sizeof(double));\n  }\n  n = g->shapes[xsno].nv;\n  xv = g->shapes[xsno].xv;\n  \n  for(i=0, j=0; i<n; i+=2, j++){\n    if( impoint(g, rno+j, sno+j, flag, type, x, y, xv[i], xv[i+1]) ){\n      return(1);\n    }\n  }\n  return(0);\n}\n\nvoid imimagemaski(GFilt g, int UNUSED(rno), int UNUSED(sno),\n                  int UNUSED(flag), int UNUSED(type),\n		  double UNUSED(x), double UNUSED(y))\n{\n  int i, j;\n  int got;\n  int fsize;\n  int mblock;\n  FilterMask xmasks;\n\n  \n  \n  mblock = (int)(((double)((g->xmax - g->xmin + 1)/g->block)/(double)g->maskdim) + 0.5);\n  if( mblock < 1 ){\n    g->masks = NULL;\n    g->nmask = 0;\n  }\n  \n  fsize = g->nmask * sizeof(FilterMaskRec) * mblock;\n  xmasks = (FilterMask)calloc(1, fsize);\n  for(got=0, i=0; i<g->nmask; i++){\n    xmasks[got].region = g->masks[i].region;\n    xmasks[got].y = (g->masks[i].y - 1.0) * mblock + 1.0;\n    xmasks[got].xstart = (g->masks[i].xstart - 1.0) * mblock + 1.0;\n    xmasks[got].xstop = (g->masks[i].xstop - 1.0) * mblock + 1.0;\n    \n    for(j=1; j<mblock; j++){\n      xmasks[got+j].region = xmasks[got].region;\n      xmasks[got+j].y = xmasks[got].y+j;\n      xmasks[got+j].xstart = xmasks[got].xstart;\n      xmasks[got+j].xstop = xmasks[got].xstop;\n    }\n    got += mblock;\n  }\n  \n  qsort(xmasks, got, sizeof(FilterMaskRec), imagemaskcmp);\n  /* we don't free previous masks, just overwrite, because they were not\n     allocated in the first place */\n  \n  \n  g->masks = xmasks;\n  g->nmask = got;\n\n  \n  for(i=0; i<got; i++){\n    if( !g->ybuf[g->masks[i].y] )\n      g->ybuf[g->masks[i].y] = i+1;\n  }\n}\n\nint imimagemask(GFilt g, int UNUSED(rno), int UNUSED(sno),\n                int UNUSED(flag), int UNUSED(type),\n		double x, double y)\n{\n  int i;\n  int ix, iy;\n\n  if( g->nmask == 0 )\n    return(0);\n  if( g->evsect ){\n    if( g->usebinsiz ){\n      x = (int)((x - g->tlminx)/g->binsizx + 1.0);\n      y = (int)((y - g->tlminy)/g->binsizy + 1.0);\n    }\n    else{\n      x = (int)((x - g->tlminx) + 1.0);\n      y = (int)((y - g->tlminy) + 1.0);\n    }\n    x = (int)((x - g->xmin)/g->block + 1.0);\n    y = (int)((y - g->ymin)/g->block + 1.0);\n    if(y < g->y0) return 0;\n    if(y > g->y1) return 0;\n    if(x < g->x0) return 0;\n    if(x > g->x1) return 0;\n  }\n  ix = (int)x;\n  iy = (int)y;\n  \n  i = g->ybuf[iy];\n  \n  if( i != 0 ){\n    i--;\n    for(; i<g->nmask; i++){\n      if( g->masks[i].y > iy ){\n	return(0);\n      }\n      if( (ix >= g->masks[i].xstart) && (ix <= g->masks[i].xstop) ){\n	g->rid = g->masks[i].region;\n	return(1);\n      }\n    }\n  }\n  return(0);\n}\n\n";