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

int main(int argc, char* argv[])
{

  char ch;
  int i;

  double zero_B=-0.12, zero_V=0.0, zero_R=.18;

  double zero;


//determining which filter the data is in so the zero point can be correctly set
  for(i=0;i<13;i++){
    ch=argv[1][i];
    switch(ch){
    case 'B': zero=zero_B;
      break;
    case 'V': zero=zero_V;
      break;
    case 'R': zero=zero_R;
      break;}
  }  

  FILE *fake_mag, *real_mag, *out, *txt_out_m, *txt_out_f;
  FILE *txt_out_mg, *txt_out_fg;

  char julian[40];
  double JD,MJD;

  int m=100;
  int n;

  struct star{
    double mag;
    double err;};

  struct source{
    double mag;
    double err;
    double gerr;};

  star ref[m];
  source mag_source;
  source flux_source;
  double err_tmp;

  star ref_actual[m];
  
  double offset[m];
  double offset_err[m];
  
  double offset_sum;
  double offset_mean;
  double err_sumsq;
  double mean_err;
  double sumsq_res;
  double mean_err_gauss;
  
  char in_raw[500],in_cal[500],out_all[500],out_txt_m[500],out_txt_f[500];
  char out_txt_mg[500],out_txt_fg[500];

  //setting all of the proper input/output files

  sprintf(in_raw,"%s.dat",argv[1]);
  sprintf(in_cal,"%scal.dat",argv[1]);
  sprintf(out_all,"../corrected_magnitudes/%s.dat",argv[1]);
  sprintf(out_txt_m,"../corrected_magnitudes/%s.merr",argv[1]);
  sprintf(out_txt_f,"../corrected_magnitudes/%s.ferr",argv[1]);
  sprintf(out_txt_mg,"../corrected_magnitudes/%s.mgerr",argv[1]);
  sprintf(out_txt_fg,"../corrected_magnitudes/%s.fgerr",argv[1]);

  fake_mag=fopen(in_raw,"r");
  real_mag=fopen(in_cal,"r");
  out=fopen(out_all,"w");

  txt_out_m=fopen(out_txt_m,"w");
  txt_out_f=fopen(out_txt_f,"w");
  txt_out_mg=fopen(out_txt_mg,"w");
  txt_out_fg=fopen(out_txt_fg,"w");

  int stars=0;


  //getting the photmetry from IRAF and the calibration data
  do{
    fscanf(real_mag,"%lf %lf",&ref_actual[stars].mag,&ref_actual[stars].err);
    stars++;}while(!feof(real_mag));

  while(1){

    fscanf(fake_mag,"%s = %lf",julian,&JD);
    MJD=JD-2400000.0;
    for(i=0;i<stars;i++){ 
      fscanf(fake_mag,"%lf %lf",&ref[i].mag,&ref[i].err);}
    fscanf(fake_mag,"%lf %lf",&mag_source.mag,&err_tmp);
    if(feof(fake_mag)){break;}
    
    fprintf(out,"%s = %lf\n",julian,MJD);
  
    offset_sum=0;
    err_sumsq=0;
    sumsq_res=0;
    n=0;

    //calculating the offset between reality and the IRAF photometry
    for(i=0;i<stars;i++){
      if(ref[i].mag!=00.000){
	offset[i]=ref[i].mag-ref_actual[i].mag;
	fprintf(out,"%lf ",offset[i]);
	offset_sum+=offset[i];
	offset_err[i]=sqrt(ref[i].err*ref[i].err+ref_actual[i].err*ref_actual[i].err);
	fprintf(out,"%lf\n",offset_err[i]);
	err_sumsq+=offset_err[i]*offset_err[i];
	n++;}
      //dealing with case of no photmetry from a given star
      else{
      fprintf(out,"%lf %lf\n",ref[i].mag,ref[i].err);}
    }
    
    offset_mean=offset_sum/(double)n;
    fprintf(out,"%lf ",offset_mean);
    mean_err=sqrt(err_sumsq)/(double)n;
    fprintf(out,"%lf ",mean_err);
    
    for(i=0;i<stars;i++){
      if(ref[i].mag!=00.000){
	sumsq_res+=(offset_mean-offset[i])*(offset_mean-offset[i]);}
    }
    
    if(n>1){
      mean_err_gauss=sqrt(sumsq_res/((double)n-1.0));
      fprintf(out,"%lf\n",mean_err_gauss);}

    else{
      fprintf(out,"Not enough reference stars!\n");
      mean_err_gauss=mean_err;}

    mag_source.mag-=offset_mean;
    mag_source.err=sqrt(err_tmp*err_tmp+mean_err*mean_err);
    mag_source.gerr=sqrt(err_tmp*err_tmp+mean_err_gauss*mean_err_gauss);
    flux_source.mag=pow(10.0,(mag_source.mag-8.90+zero)/-2.5);
    flux_source.err=-1.0*flux_source.mag*(pow(10.0,-0.4*mag_source.err)-1.0);
    flux_source.gerr=-1.0*flux_source.mag*(pow(10.0,-0.4*mag_source.gerr)-1.0);


    fprintf(out,"%lf %lf %lf\n",mag_source.mag,mag_source.err,
	    mag_source.gerr);
    fprintf(out,"%lf %lf %lf\n",flux_source.mag,flux_source.err,
	    flux_source.gerr);
    fprintf(txt_out_m,"%lf %lf %lf\n",MJD,mag_source.mag,mag_source.err);
    fprintf(txt_out_f,"%lf %lf %lf\n",MJD,flux_source.mag,flux_source.err);
    fprintf(txt_out_mg,"%lf %lf %lf\n",MJD,mag_source.mag,mag_source.gerr);
    fprintf(txt_out_fg,"%lf %lf %lf\n",MJD,flux_source.mag,flux_source.gerr);
    
  }

}