/***************************************************************************
 *
 * excmd - Extract a CMD from two PLANET photometry archives.
 *
 * Usage:
 *    excmd Varchive Iarchive CMDfile [-axv]
 *
 * where:
 *    Varchive = name of the V-band photometry archive
 *    Iarchive = name of the I-band photometry archive
 *    CMDfile  = name of the output ascii CMD file to create
 *
 * options:
 *    -a = extract all stars, not just types 11-13, and extra info
 *    -x = use internal cross-references (misaligned archives)
 *    -v = send verbose output to tty
 *
 * Note that V and I above are merely placeholders to illustrate
 * the command, in principle any two filter archives may be used.
 *
 * R. Pogge, OSU Astronomy
 * pogge@astronomy.ohio-state.edu
 * 2000 May 18
 *
 * Modification History:
 *   2000 May 22: enabled -x flag to handle the case of mis-aligned archives
 *                with cross-reference IDs set in the StarInfo data structure.
 *    2000 Aug 2: Fixed the TypMean problem, whereby JP thought he
 *                would save space by coding the observation number and
 *                the mean DoPhot type as one 2-byte integer.  [rwp/osu]
 *
 ***************************************************************************/

#include "fsui.h"

void Usage();

main(int narg, char *argv[]) 
{ 
  SUIVIFIP *sfp1,*sfp2;   /* archive file pointers */
  FILE *fpCMD;            /* output file pointer */

  char ArchRoot1[80];   /* archive and output filenames, with and w/o paths */
  char Archive1[256];
  char ArchRoot2[80];
  char Archive2[256]; 
  char CMDRoot[80];       
  char CMDFile[256];  

  STARINFO *star_i1;       /* StarInfo, TimeInfo, and Mesure structs */
  STARINFO *star_i2;
  TIMEINFO timeinfo;
  MESURE mesure;

  /* Photometric data */

  int NStars1, NObs1;
  float Mag1;
  float MagErr1, RMS1; 
  int MeanTyp1, NObsl1;
  int NRefs1; 
  int RefList1[10];

  int NStars2, NObs2;
  float Mag2;
  float MagErr2, RMS2; 
  int MeanTyp2, NObsl2;
  int NRefs2;
  int RefList2[10];

  float Color, ColErr;

  /* Working Variables */

  int i, j, iarg;
  int nstar;

  int allstars = 0;    /* default: only type 11-13 stars   */
  int docross = 0;     /* default: assume aligned archives */
  int verbose = 0;     /* default: silent operation        */

  int a;
  int NMax;
  int LensID;

  /* Some initializations */

  star_i1 = &StarInfo;
  star_i2 = &StarInfo2; 
  NStars1 = 0;

  /* The command line must have at least the command and 3 arguments */

  if (narg < 4) { 
    Usage();
    exit(1);
  }

  /* Get the default input and output directories from the respective 
   * environment variables if defined, otherwise assume the cwd */

  if ((InDir=getenv("IN_FILE")) == NULL)
    InDir  = "";
  if ((OutDir=getenv("OUT_FILE")) == NULL)
    OutDir = "";

  /* Parse the command line.  The first 3 arguments, in order, are
   *   name of the first photometry archive
   *   name of the second photometry archive
   *   name of the output merged file to create
   */

  iarg = 1;
  strcpy(ArchRoot1,argv[iarg++]);          
  sprintf(Archive1 ,"%s%s",InDir,ArchRoot1);

  strcpy(ArchRoot2,argv[iarg++]);          
  sprintf(Archive2 ,"%s%s",InDir,ArchRoot2);

  strcpy(CMDRoot,argv[iarg++]);         
  sprintf(CMDFile ,"%s%s",OutDir,CMDRoot);

  /* the last arguments are optional flags to control other functions */

  for (a=iarg; a<narg; a++) 
    switch(tolower(argv[a][1])){             
    case 'v' :  /* verbose mode */
      verbose = 1;
      break;
    case 'a' :  /* select all stars irrespective of DoPhot type + extras */
      allstars = 1;
      break;
    case 'x' :  /* archives are misaligned, use cross-reference lists */
      docross = 1;
      break;
    default : 
      fprintf(stderr,"Error: unknown option %s\n", argv[a]); 
      Usage();
      exit(1);
      break;
    }

  if (verbose) {
    fprintf(stdout,"Archive file 1: %s\n", Archive1);
    fprintf(stdout,"Archive file 2: %s\n", Archive2);
    fprintf(stdout,"Output CMD file: %s\n", CMDFile);
  }

  /* Open the archives, and test to see if they are byteswapped.  If
   * so, close them and gripe.
   */

  if ((sfp1 = SuiviOpen(Archive1,SUOF_RO)) == NULL) {
    fprintf(stderr,"Error opening archive file %s\n",ArchRoot1); 
    exit(1); 
  }
    
  if ((sfp2 = SuiviOpen(Archive2,SUOF_RO)) == NULL) {
    fprintf(stderr,"Error opening archive file %s\n",ArchRoot2); 
    exit(1); 
  }
  
  NStars1 = SuiviGetNbStars(sfp1);   /* is archive data byte swapped? */
  if ((NStars1 < 0) || (NStars1 > 100000)) { 
    fprintf(stderr,"Error: one or both archives are byteswapped.\n");
    fprintf(stderr,"       Please convert to unswapped and try again.\n");
    exit(1);
  }

  /* open an ascii file as an "append" text stream */

  if ((fpCMD =  fopen(CMDFile,"a")) == NULL) {
    fprintf(stderr,"Error opening output file %s\n",CMDRoot);
    exit(1);   
  }  

  NObs1 = SuiviGetNbMesures(sfp1);      
  NStars1 = SuiviGetNbStars(sfp1);
  NMax = NStars1;

  /* Extract the IDs of the microlensing event and reference stars */
  
  GetMicrostarNumPL(sfp1,&LensID);

  NObs2  = SuiviGetNbMesures(sfp2);     
  NStars2 = SuiviGetNbStars(sfp2);

  /* Get the number of reference stars in each archive */
 
  GetRefstarNumPL(sfp1,RefList1,(int *) &NRefs1); 
  GetRefstarNumPL(sfp2,RefList2,(int *) &NRefs2); 

  if (verbose) {
    fprintf(stdout,"%s: NObs %d, NStars %d, NRefs %d\n ",
	    ArchRoot1,NObs1,NStars1,NRefs1);
    fprintf(stdout,"%s: NObs %d, NStars %d, NRefs %d\n ",
	    ArchRoot2,NObs2,NStars2,NRefs2);
  }

  /* write a header into the output file, use # for SM compatibility */

  if (allstars) {
    fprintf(fpCMD,"#\n");
    fprintf(fpCMD,"# Archive 1: %s\n",ArchRoot1);
    fprintf(fpCMD,"# Archive 2: %s\n",ArchRoot2);
    fprintf(fpCMD,"# Lens ID: %d\n",LensID);
    fprintf(fpCMD,"#\n");
    fprintf(fpCMD,
    "#Star   X1     Y1   Mag1   Err1  rms1   N1  T1 Mag2   Err2  rms2   N2  T2  Color  Err\n"); 
  } else {
    fprintf(fpCMD,
	    "#Star Mag1   Err1  rms1  Mag2   Err2  rms2   Color  Err\n"); 
  }
  
  /**********************************************************************
   * Read the photometry data from the archives, and output the 
   * merged 2-filter data.  Read only up to the max number of stars in
   * the first archive.  If "docross" is set, use the cross-ref ID's
   * in the starinfo block instead of assuming the archives are aligned
   */

  for (nstar=1; nstar<=NMax; nstar++) {
    NObs2 = SuiviGetNbMesures(sfp2);
    if (SuiviReadStarInfo(sfp1,nstar,(char *)(star_i1)) != 0) {
      fprintf(stderr,"Error: cannot read StarInfo in archive %s\n",ArchRoot1); 
      exit(1);
    }  
  
    /* 
     * The default assumption is that that the archives are aligned, which
     * means that they use a common template and so star NNN in one is the
     * same as star NNN in the other.
     *
     * If docross is set, the archives are mis-aligned but cross-references
     * are stored in the StarInfo data block.  The cross-ref ID is in the
     * first element of the XRef[] array in the starinfo struct of the first
     * archive.
     */

    if (docross) {
      if ((star_i1->XRef[0]) > 0)
	if (SuiviReadStarInfo(sfp2,(star_i1->XRef[0]),(char *)(star_i2)) != 0){
	  printf("Error: cannot read StarInfo for archive %s\n",ArchRoot2); 
	  exit(1);
	}
    } else {
      if (SuiviReadStarInfo(sfp2,nstar,(char *)(star_i2)) != 0) {
	printf("Error: cannot read StarInfo for archive %s\n",ArchRoot2); 
	exit(1);
      }
    }

    /* Store the StarInfo data from each archive in working variables */

    NObsl1   = (int) (star_i1->TypMean/100.0) ;
    Mag1     = (float) star_i1->FlxMean;
    MagErr1  = (float) star_i1->SigRef/1000.0;
    RMS1  = (float) star_i1->SigMean/1000.;
    MeanTyp1 = star_i1->TypMean-(NObsl1*100.0);

    /* sigh... fix JPs TypMean bug */

    if (NObsl1 < 0) NObsl1 = NObsl1 + 655;
    if (MeanTyp1 < 0) MeanTyp1 = MeanTyp1 + 36;

    NObsl2   = (int) (star_i2->TypMean/100.0) ;
    Mag2     = (float) star_i2->FlxMean;
    MagErr2  = (float) star_i2->SigRef/1000.0;
    RMS2  = (float) star_i2->SigMean/1000.;
    MeanTyp2 = star_i2->TypMean-(NObsl2*100.0);

    /* and again... */

    if (NObsl2 < 0) NObsl2 = NObsl2 + 655;
    if (MeanTyp2 < 0) MeanTyp2 = MeanTyp2 + 36;

    /* compute the color and its formal uncertainty */

    Color    = Mag1 - Mag2;
    ColErr   = (float)(sqrt((float)(MagErr1*MagErr1 + MagErr2*MagErr2))) ;

    /* 
     * if a mis-aligned archive, handle the case that the star in
     * question has no measurement in the other archive 
     */

    if (docross && (star_i1->XRef[0] < 1)) { 
      NObs2 = 0;
      MagErr2 = 0;
      RMS2 = 0;
      NObsl2 = 0; 
      MeanTyp2 = 0; 
      Mag2 = 0; 
      star_i1->NumStar = -1;
      star_i2->NumStar = -1; 
      star_i2->FlxMean = 0; 
      star_i2->X = 0.0;
      star_i2->Y = 0.0;
    }

    if (verbose && (star_i1->NumStar != -1)) {
      if (docross) {
	printf("Star %d (%d):\n",star_i1->NumStar, star_i1->XRef[0]);
      } else {
	printf("Star %d:\n",star_i1->NumStar);
      }
      printf("V: %f %f %f %d %d \n",Mag1,MagErr1,RMS1,NObsl1,MeanTyp1); 
      printf("I: %f %f %f %d %d \n",Mag2,MagErr2,RMS2,NObsl2,MeanTyp2); 
    }

    /* Write the data into the CMD file */

    if (allstars) {
      if ((star_i1->NumStar > 0) && ((MeanTyp1 <= 17) && (MeanTyp2 <= 17)))
	fprintf(fpCMD,
	"%5d %6.1f %6.1f %6.3f %5.3f %5.3f %3d %3d %6.3f %5.3f %5.3f %3d %3d %6.3f %5.3f\n",
		star_i1->NumStar, star_i1->X, star_i1->Y, 
		Mag1, MagErr1, RMS1, NObsl1, MeanTyp1,
		Mag2, MagErr2, RMS2, NObsl2, MeanTyp2, Color, ColErr);
    } else {
      if ((MeanTyp1 <= 13) && (MeanTyp2 <= 13)) {
	if (star_i1->NumStar > 0) 
	  fprintf(fpCMD,
		  "%5d %6.3f %5.3f %5.3f %6.3f %5.3f %5.3f %6.3f %5.3f\n",
		  star_i1->NumStar, Mag1, MagErr1, RMS1,
		  Mag2, MagErr2, RMS2, Color, ColErr);
      }
    }
    
  } /* loop over all stars in archive 1 */

  /* Close the photometry archives and the CMD file */

  SuiviClose(sfp1);
  SuiviClose(sfp2);
  fclose(fpCMD);

  fprintf(stdout,"Done.\n");

}

/***************************************************************************/

void Usage()
{

  fprintf(stdout,"\nUsage: excmd Varchive Iarchive CMDfile [-axv]\n");
  fprintf(stdout,"\nwhere:\n");
  fprintf(stdout,"   Varchive = name of the V-band photometry archive\n");
  fprintf(stdout,"   Iarchive = name of the I-band photometry archive\n");
  fprintf(stdout,"   CMDfile  = name of the output ascii file to create\n");
  fprintf(stdout,"   -a = extract all stars, not just types 11-13, and\n");
  fprintf(stdout,"        add other info to the CMD file\n");
  fprintf(stdout,"   -x = use internal cross-references (misaligned archives)\n");
  fprintf(stdout,"   -v = verbose output to tty\n");
  fprintf(stdout,"\nExtract a CMD by merging two PLANET photometry archives.\n");
  fprintf(stdout,"Creates an ascii data file with the mean photometry in\n");
  fprintf(stdout,"each filter for each star and the instrumental color.\n");
  fprintf(stdout,"V and I are just placeholds, any 2 filters will work, provided\n");
  fprintf(stdout,"that the archives have a common template.\n");
  fprintf(stdout,"If the two archives do not share a common template but have\n");
  fprintf(stdout,"have been cross-referenced, use the -x option.\n");

  return;

}