#include <stdio.h>
FILE *fp,*fp1;

char  head[144][80],h[72][80],f1[30];
short dd[512];
short a[2048][2036],b[4096][4072],c[4096][40];
int   n1,n2,kh;
main(int ac, char **av)
{
  int i,j,k,id,ip;
  float x;
  if(ac<2){
     printf("\n\t*** change orgion_image to 4_(4032+20*2)*4096  **\n");
     printf("\n\tUsage: d4 *.in_fits\n\n"); exit(0);
  }
  for(ip=1;ip<ac;ip++){
    fp=fopen(av[ip],"rb");
    if(fp==0){ printf("\n\tinput file not found ! %s\n\n",av[ip]); continue; }
    x=2036*2048/1440.; i=x; x-=i; if(x!=0.)x=1.-x; id=x*1440.+0.5;
//    printf("%d\n",id);   512;
    n1=4072;   n2=4096;
    fread(head,108,80,fp);
    sscanf(&head[2][28],"%d",&k);
    if(k){ printf("\n\tinput_file not a raw_data\n\n"); continue; }
    k=indexpos(head,"CCDBIN1 ",108);
      if(k!=108){
        sscanf(&head[k][20],"%d",&k);
        if(k!=1){ printf("\n\t ****** COMBINE CCD ******* \n"); exit(0); }
      }
   
    kh=108; k=indexpos(head,"END     ",108);
    if(k==108){ kh+=36; fread(&head[k][0],36,80,fp); }
    for(i=0;i<7;i++)for(j=0;j<80;j++)head[i][j]=32;
    strcpy(&head[ 0][0], "SIMPLE  =                    T ");
    strcpy(&head[ 1][0], "BITPIX  =                   16 ");
    strcpy(&head[ 2][0], "NAXIS   =                    2 ");
    strcpy(&head[ 3][0], "NAXIS1  =                 4072 / 4032+20*2");
    strcpy(&head[ 4][0], "NAXIS2  =                 4096 ");
    strcpy(&head[ 5][0], "BZERO   =                32768 ");
    strcpy(&head[ 6][0], "BSCALE  =                    1 ");
    strcpy(&head[ 7][0], "CCD_NO: =                    ? ");
    for(i=0;i<8;i++)head[i][31]='/';  head[3][42]=32;
    for(i=33;i<60;i++)head[7][i]=32;
    for(k=1;k<=4;k++){ head[7][29]=k+48;
      fread(h,72,80,fp); for(j=0;j<78;j++)printf("%c",h[24][j]); printf("\n");
      fread(a,2036*2048,2,fp); fread(dd,2,id,fp);
      if(k<3)for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[n2-1-i][n1-1-j]=a[i][j];
      else   for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[i][j]=a[i][j];

      fread(h,72,80,fp); for(j=0;j<78;j++)printf("%c",h[24][j]); printf("\n");
      fread(a,2036*2048,2,fp); fread(dd,2,id,fp);
      if(k<3)for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[n2-1-i][j]=a[i][j];
      else   for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[i][n1-1-j]=a[i][j];

      fread(h,72,80,fp); for(j=0;j<78;j++)printf("%c",h[24][j]); printf("\n");
      fread(a,2036*2048,2,fp); fread(dd,2,id,fp);
      if(k<3)for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[i][n1-1-j]=a[i][j];
      else   for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[n2-1-i][j]=a[i][j];

      fread(h,72,80,fp); for(j=0;j<78;j++)printf("%c",h[24][j]); printf("\n");
      fread(a,2036*2048,2,fp); fread(dd,2,id,fp);
      if(k<3)for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[i][j]=a[i][j];
      else   for(i=0;i<2048;i++)for(j=0;j<2036;j++)b[n2-1-i][n1-1-j]=a[i][j];

      for(i=0;i<4096;i++)for(j=0;j<40;j++)c[i][j]=b[i][2016+j];
      for(i=0;i<4096;i++)for(j=2016;j<4032;j++)b[i][j]=b[i][j+40];
      for(i=0;i<4096;i++)for(j=0;j<40;j++)b[i][4032+j]=c[i][j];

      strcpy(f1,av[ip]); j=strlen(f1)-5;
      sprintf(&f1[j],"_%d.fit",k);
      fp1=fopen(f1,"wb");
      fwrite(head,80,kh,fp1);
      fwrite(b,4072*4096,2,fp1);
      fclose(fp1);
      printf("%s produced!\n",f1);
    }
    fclose(fp);
  }
}
