#include <stdio.h>
#include <math.h>
#define n1 4096
#define n2 4032
FILE *fp;
char a[80];
char f1[60],f2[60],f3[60];
float x1[3000],x2[3000],x3[3000];
float u[1000],v[1000],w[1000];
int n=0,ip=0;
float b[n2][n1];
char  h[72][80];
ifm(float x0, float y0)
{
  int i;
  float z;
  for(i=0;i<=n;i++){
    z=x0-x2[i]; if(z<0.)z=-z;
    if(z>2.)continue;
    z=y0-x3[i]; if(z<0.)z=-z;
    if(z>2.)continue;
    return(0);
  }    
  return(1);
}

main(int ac, char **av)
{
  int i,j,k;
  float x,y,z,x0,y0;
  if(ac<2){
    printf("\n\t ***** auto produce Coordination_file for get seeing...***\n");
    printf("\n\t before run: d4 focus.0067.fits  (divide to 4 file)");
    printf("\n\t before run: e4 focus.0067.fit_? (-overscan,addftis,x-y)");
    printf("\n\t after  run: e4f eocus.0067_?\n");
    printf("\n\t Ex: e4dat eocus.0067_?  [!] (need jsex support)");
    printf("\n\t                [!] (produce match file, named a.fit)\n\n");
    exit(0);
  }
 
  if(av[1][0]!='e'){ printf("\n\tfits_file_name should be e*_?\n\n"); exit(0); }
  strcpy(f1,av[1]); k=strlen(f1); if(f1[k-3]=='f' && f1[k-2]=='i')f1[k-4]=0;
  strcpy(f2,f1); strcpy(f3,f1);
  strcat(f2,".fit");  strcat(f3,".dat");
  fp=fopen(f2,"rb"); if(fp==0){ printf("%s not found!\n",f2); exit(0); }
  fclose(fp); sprintf(a,"jsex %s",f2); system(a);
  fp=fopen("test.cat","r");
  for(i=0;i<3;i++)fgets(a,80,fp); 
l10:
  fgets(a,80,fp); if(feof(fp))goto l20;
  sscanf(a,"%f %f %f",&x1[n],&x2[n],&x3[n]);
  if(n<3000-1)n++;
  goto l10;
l20:
  fclose(fp);
  for(i=0;i<=n;i++){
    x=x2[i]; y=x3[i];
//  printf("%8.1f %8.1f\n",x,y);
    if(x<200 || x> n1-200)continue;
    if(x<2048)k=1; else k=-1;
    if(y<15  || y>n2-15)continue;
    for(j=0;j<=n;j++){
      z=(y-x3[j]); if(z<0.)z=-z;
      if(z>1.8)continue;
      z=(x-x2[j])*k-60; if(z<0.)z=-z;
      if(z>3.2)continue;
      x0=(x+x2[j])*0.5; y0=(y+x3[j])*0.5;
      k=ifm(x0,y0);
      if(k){ u[ip]=x; v[ip]=y; w[ip]=x1[i]; ip++; }
      break;
    }      
  }
  printf("total: %d\n",ip);
  if(ac>2){ mark(); printf("******* produced mark file named a.fit *****\n"); }
  fp=fopen(f3,"w");
  for(k=0;k<ip;k++){
    i=u[k]+.5; j=v[k]+.5;
    fprintf(fp,"%4d %4d %6.1fmag\n",i,j,w[k]);
  }
  fclose(fp);
}

mark()
{
  int i,j,k,x0,y0;
  fp=fopen(f2,"rb"); fread(h,72,80,fp); fread(b,n1*n2,4,fp); fclose(fp);
  for(k=0;k<ip;k++){
    x0=u[k]; y0=v[k];
    for(i=x0-10;i<x0+10;i++)b[y0][i]=0;
    for(j=y0-10;j<y0+10;j++)b[j][x0]=0;
  }
  fp=fopen("a.fit","wb"); fwrite(h,72,80,fp); fwrite(b,n1*n2,4,fp);
  fclose(fp);
}