#include <math.h>
#include <stdlib.h>
#include <stdio.h>
  FILE *fd, *fb;
void write_data(void)
{
  int dummy;
  int i, k, pc;
  double x,y,z,v;
  unsigned long N;

  N = 256;
  x =0.;
  y =128.;

  fprintf(fd,"%12g %12g\n",x,y);

  fwrite(&N, sizeof (unsigned int),1,fb);
  fwrite(&x, sizeof (double),1,fb);
  fwrite(&y, sizeof (double),1,fb);
  fwrite(&N, sizeof (unsigned int),1,fb); 
}
void read_header(void)
{
  int dummy;
  int i, k, pc;
  double x,y,z,v;
  unsigned long Nn;
  unsigned long Npart[6];
  double        mass[6];
  double        time;
  double        redshift;
  long          flag_sfr;
  long          flag_feedback;
  unsigned long npartTotal[6];
  long          flag_cooling;
  long          num_files;
  double        BoxSize;
  double        Omega0;
  double        OmegaLambda;
  double        HubbleParam;
  long          flag_stellarage;
  long          flag_metals;
  long          hashtabsize;
  char          fill[84];
  int          NpartN[6];
 
  printf(" Start reading \n");

  fread(&Nn,         sizeof (unsigned long),1,fb);
  printf("N     = %lu \n ",Nn); 
  
  fread(&Npart,      sizeof (unsigned long),6,fb);
  fread(&mass,       sizeof (double),6,fb);
  fread(&time,       sizeof (double),1,fb);
  fread(&redshift,   sizeof (double),1,fb);
  fread(&flag_sfr,   sizeof (long),1,fb);
  fread(&flag_feedback, sizeof (long),1,fb);
  fread(&npartTotal,   sizeof (unsigned long),6,fb);
  fread(&flag_cooling, sizeof (long),1,fb);
  fread(&num_files,    sizeof (long),1,fb);
  fread(&BoxSize,      sizeof (double),1,fb);
  fread(&Omega0,       sizeof (double),1,fb);
  fread(&OmegaLambda,  sizeof (double),1,fb);
  fread(&HubbleParam,  sizeof (double),1,fb);
  fread(&flag_stellarage, sizeof (long),1,fb);
  fread(&flag_metals,  sizeof (long),1,fb);
  fread(&hashtabsize,  sizeof (long),1,fb);
  fread(&fill,         sizeof (char),84,fb);
  fread(&Nn,           sizeof (unsigned long),1,fb);

  for (i =0; i<6; ++i) {
    NpartN[i] = Npart[i];
  }
  
  printf("N     = %lu \n ",Nn); 
  printf("Npart 0 = %d %lu\n",NpartN[0],npartTotal[0]);
  printf("      1 = %d %lu\n",NpartN[1],npartTotal[1]);
  printf("      2 = %d %lu\n",NpartN[2],npartTotal[2]);
  printf("      3 = %d %lu\n",NpartN[3],npartTotal[3]);
  printf("      4 = %d %lu\n",NpartN[4],npartTotal[4]);
  printf("      5 = %d %lu\n",NpartN[5],npartTotal[5]);

  printf("Mass 0  = %lf %lf %lf\n",mass[0],mass[1],mass[2]);
  printf("N files = %lu\n", num_files);
  printf("N files = %lf\n", HubbleParam );
    }

int main(int argc, char **argv)
{
  fd = fopen("test.dat","r");
  fb = fopen("ics_1000mpc_2048_2lpt.0","rb");
  /* fb = fopen("ics_1000mpc_2048_2lpt.0 ","wb"); */

   read_header();  
    /*rwrite_data(); */
  
  fclose(fd);
  fclose(fb);
    exit(0);
}