/* void.c                                                       */
           /* to find the voids in the random number coverage of the plane */
           /*                                                              */
           /* algorithm: find the circles in the plane whitin which        */
           /*            no random points is found                         */
           /*                                                              */

/* Syntax:                                                                 */
/* random_demo_2_series 2 1000 | void > data_file                          */

                                                     /* Michel Vallieres */


#include <stdio.h>
#include <math.h>
#include <stdlib.h>

#define NMAX 100000
#define dR 0.01
#define Nx 100
#define Ny 100
#define Npoint 55

int main ( int argc, char *argv[] )
{
  int    N;
  double x[NMAX], y[NMAX];
  double dxc, dyc, d;
  double R, R2, xc, yc, dist_x, dist_y;
  int    count, icx, icy, ir;
  double radius, angle, xp, yp;
  int    ic;

                                                 /* read in the random points coordinates */
  N = 0;
  while ( scanf( "%lf %lf", &x[N], &y[N] ) != EOF )
      N++;
  fprintf( stderr,  " Number of random numbers %d \n", N );
                                               
  dxc = 1.0/(double)(Nx-1);                     /* scan parameters for circle center coordinates */
  dyc = 1.0/(double)(Ny-1);
  R = 0.1;                                      /* scan circle radius */
  while ( R > 0.05 )
    {
      R2 = R*R;
      for ( icx=0 ; icx<Nx ; icx++ )            /* scan x direction (circle center) */
        {
          xc = icx*dxc;
          for ( icy=0 ; icy<Ny ; icy++ )        /* scan y direction (circle center) */
            {
              yc = icy*dyc;
              count = 0;
              for ( ir=0 ; ir<N ; ir++ )        /* scan over random points */
                {
                  dist_x = xc - x[ir];          /* check if within circle */
                  dist_y = yc - y[ir];
                  if ( dist_x*dist_x + dist_y*dist_y < R2 )
                    {
                      count++;
                      break;                    /* found a point - quit scanning random points */
                    }
                }
              if ( count == 0 )                 /* empty circle */
		{
                  fprintf(stderr,  " Circle of radius %f at ( %f, %f )\n", R, xc, yc );
		                                /* output dots within circle for plotting filled circle */
                  for ( ic=0 ; ic<Npoint ; ic++ )
		    {
                      angle = ic*2.0*M_PI/(Npoint-1);
                      radius = R;
                      while ( radius > dR )
			{
                           xp = xc + radius*cos(angle);
                           yp = yc + radius*sin(angle);
                           if ( xp>0.  & xp<1.0 && yp>0.0 && yp<1.0 )
                               printf( "%f %f\n", xp, yp );
                           radius = radius - dR;
			}
		    }
		} 
            }
        }
      R = R - dR;                               /* next smaller circle size */
    }
}