#include <stdio.h>
#include <math.h>
char c1[14],c2[14],head[72][80];
double a8[8],b8[8],pi,cx,cy;
FILE *fp,*fp1;
int  n1,n2;
#define bok 48.

main(ac,av)
int ac; char *av[];
{
  char a[90],f1[50],f2[50];
  double x,y,xa,xd,epoch,alpha,delta,z,z1;
  int i,j,k,iz;
  pi=4.*atan(1.0); cx=pi/12.; cy=pi/180.;
  if(ac<2){
    printf("\n\t usage: u_xy6ad filename");
    printf("\n\t                       2010.10\n");
    printf("\n\t    file is come from pip231");
    printf("\n\t input:   *.als|lst, *.mes|fit");
    printf("\n\t output:  *.cat6");
    printf("\n\t  epoch = 2000.0\n\n");
    exit(0);
  }
  epoch=2000.0;
  fp1=fopen(av[1],"r"); if(fp1==0){printf("\ninputfile not found!\n"),exit(0);}
  strcpy(f1,av[1]);  j=nindex(f1,"."); if(j)f1[j]=0; 
  strcpy(f2,f1); strcat(f2,".mes");
  fp=fopen(f2,"r");
  if(fp){
    for(i=0;i<19;i++){
      fgets(a,90,fp);
      if(a[0]!='A' || a[1]!='8' || a[2]!='1')continue;
      sscanf(&a[10],"%le",&a8[0]);
      for(j=1;j<8;j++){
        fgets(a,90,fp); sscanf(&a[10],"%le",&a8[j]);
      }
      break;
    }
    fclose(fp);
  } else {
    strcpy(f2,f1); strcat(f2,".fit");
    fp=fopen(f2,"rb"); if(fp==0){ printf("cannot find %s or .mes\n\n");exit(0);}
    fread(head,72,80,fp); fclose(fp);
    k=indexpos(head,"A81     ",72);
    for(i=0;i<8;i++)sscanf(&head[k+i][10],"%le",&a8[i]);
  }
  xytoad(a8,b8);

  strcpy(f2,f1); strcat(f2,".cat6");
  fp=fopen(f2,"w");
  toms2(a8[6]/cx,c1,1); toms2(a8[7]/cy,c2,0);
  printf("#  center: %s %s (2000.0)\n",c1,c2);
  fprintf(fp,"#  center: %s %s (2000.0)\n",c1,c2); 
//  printf("#    R.A.          Dec.      Mag.    Err_mag. Number      X        Y\n");
  fprintf(fp,"#    R.A.          Dec.      Mag.    Err_mag. Number      X        Y\n");
l10:
  fgets(a,90,fp1);  if(feof(fp1))goto l20;
  sscanf(a,"%d %lf %lf %lf %lf",&iz,&x,&y,&z,&z1);
  xy_rad(a8[6],a8[7],x,y,&alpha,&delta);
  toms2(alpha,c1,1); toms2(delta,c2,0);
//  printf("%s %s %8.3lf %8.3lf %6d %9.2lf %9.2lf\n",c1,c2,z,z1,iz,x,y);  
  fprintf(fp,"%s %s %8.3lf %8.3lf %6d %9.2lf %9.2lf\n",c1,c2,z,z1,iz,x,y);  
  goto l10;
l20:
  fclose(fp1); fclose(fp);
}

xytoad(x,a)
double *x,*a;
{
  double z;
  z=x[0]*x[3]-x[2]*x[1];
  a[0]= x[3]/z;
  a[1]=-x[1]/z;
  a[2]=-x[2]/z;
  a[3]= x[0]/z;
  a[4]= (x[2]*x[5]-x[4]*x[3])/z;
  a[5]= (x[4]*x[1]-x[0]*x[5])/z;
}

toms2(aa,c,k)
double aa;           // input
char c[14];          // output
int k;               // input
/*
 hour or degree to char_line
 if k=1 (hour case) in **:**:**.***
    k=0 (degree       -**:**:**.**
*/
{
  int  i,j;
  double a,b,x;
  a=aa;
  if(k==1){ if(a>24.)a-=24.; if(a<0.)a+=24.; }
  c[0]=32; if(a<0.){ a=-a; c[0]='-'; }
  i=a; b=(a-i)*60.;
  j=b; x=(b-j)*60.;
  if(j==60){ j=0; i++; }
  if(k==1)sprintf(&c[1],"%2.2d:%2.2d:%6.3f",i,j,x);
  if(k==0)sprintf(&c[1],"%2.2d:%2.2d:%5.2f",i,j,x);
  if(c[7]==32)c[7]='0';   if(c[8]==32)c[8]='0';
  if(c[7]=='6'){
    c[7]='0'; j++; sprintf(&c[4],"%2.2d",j);
    if(c[4]=='6'){
      c[4]='0'; i++; sprintf(&c[1],"%2.2d",i);
    }
    c[3]=':'; c[6]=':';
  }
}

rad_xy(ac,dc,ra,de,x,y)
double ac,dc,ra,de,*x,*y;
{
  double xi,xn, tmp,rar,der;
  rar=ra*cx; der=de*cy;
  tmp=atan(tan(der)/cos(rar-ac));
  xi=cos(tmp)*tan(rar-ac)/cos(tmp-dc);
  xn=tan(tmp-dc);
  tmp=1.+bok*(xi*xi+xn*xn);
  xi*=tmp; xn*=tmp;
  *x=b8[0]*xi+b8[2]*xn+b8[4];
  *y=b8[1]*xi+b8[3]*xn+b8[5];
}

xy_rad(ac,dc,x,y,ra,de)
double ac,dc,x,y,*ra,*de;
{
  double xi,xn,tmp,tmp1;
  double x1,y1,r1,d1;
  int i;
  xi=a8[0]*x+a8[2]*y+a8[4];
  xn=a8[1]*x+a8[3]*y+a8[5];
  tmp=1.+bok*(xi*xi+xn*xn);
  for(i=0;i<4;i++){
    x1=xi/tmp; y1=xn/tmp;
    tmp=tan(dc);  tmp1=1.-y1*tmp;
    r1=atan(x1/cos(dc)/tmp1);
    d1=atan((y1+tmp)*cos(r1)/tmp1);
    r1+=ac;
    tmp=atan(tan(d1)/cos(r1-ac));
    x1=cos(tmp)*tan(r1-ac)/cos(tmp-dc);
    y1=tan(tmp-dc);
    tmp=1.+bok*(x1*x1+y1*y1);
  }
  xi/=tmp; xn/=tmp;
  tmp=tan(dc);  tmp1=1.-xn*tmp;
  *ra=atan(xi/cos(dc)/tmp1);
  *de=atan((xn+tmp)*cos(*ra)/tmp1);
  *ra+=ac;
  tmp=(*de)-dc; if(tmp<0.)tmp=-tmp;
  if(tmp>2.){  *de=-(*de);  *ra+=pi; }
  if(*ra<0.)*ra+=(pi+pi);
  *ra/=cx; *de/=cy;
}