//
//   example program to read a BYU SIR image file
//
//   Written by DGL 16 Apr 2003
//
//   This simple program reads the header of a BYU sir-format input file 
//   using c++
//
//   should be linked with: libsir.cpp sir_ez.c sir_io.c
//   requires includes: libsir.h sir3.h
//

#include "libsir.h"  // cpp interface routines

#include <iostream.h>

void byte_file(char *fname,float smin,float smax, float *stval, 
	       struct sirheader *head);

int main(int argc, char **argv)
{
  cout << "BYU SIR c++ example program";
  if(argc < 2) {
    cout << "\nusage: %s file out1 out2 <dmin> <dmax>\n\n",argv[0];
    cout << " input parameters:\n";
    cout << "   file  =  input SIR file\n";
    cout << "   out1  =  output byte file\n";
    cout << "   out2  =  output SIR file\n";
    cout << "   dmin  =  min saturation value (optional)\n";
    cout << "   dmax  =  max saturation value (optional)\n"; 
    return 0;
  }

  // read SIR image file into memory

  float *stval;                 // pointer to image storage
  struct sirheader head;        // SIR file header structure

  stval=LoadSIR(argv[1], &head);
  if (stval == NULL) {
    cout << "*** Error reading SIR header from file" << endl;
    return -1;
  }

  // print out contents of sir file header
  print_sir_header(&head);

  // test transformations at image corners (may get errors for polar images)
  
  int nsx=head.nsx;
  int nsy=head.nsy;

  int   i, j, ix, iy, ix1, iy1;
  float x, y, alon, alat;

  cout << "\nSIR geometric transformation test results at corners of image\n";
  // note: x,y indexing is 1-based (rather than 0-based) to be
  // consistent with fortran and the documentation
  for (ix = 1; ix <= nsx; ix = ix + nsx - 1)
    for (iy = 1; iy <= nsy; iy = iy + nsy - 1) {
      i = SIRlex(ix, iy, &head);  // get image pixel index  0 ... nsx*nxy-1
      ipixtolatlon(&head, ix, iy, &alon, &alat);  // (ix,iy) -> (lat,lon)
      ilatlon2pix(&head, alon, alat, &ix1, &iy1); // (lat,lon) -> (x,y)
      printf(" Pixel (%d,%d) Lat,Lon (%3.3f,%7.3f) at (%d,%d) Value: %f\n",
	     ix,iy,alat,alon,ix1,iy1,*(stval+i)); 
     };

  if (argc < 3) return 0;

  float smin=1.e25, smax=-1.e25;
  for (i=0; i< nsx * nsy; i++)
    if (*(stval+i)> head.anodata) {
      smin=(smin > *(stval+i) ? *(stval+i) : smin);
      smax=(smax < *(stval+i) ? *(stval+i) : smax);
    }
  cout << "Data Min, Max: %f , %f\n",smin,smax;

  smin = head.v_min;
  smax = head.v_max;
  if (argc > 4) sscanf(argv[4],"%f",&smin); 
  if (argc > 5) sscanf(argv[5],"%f",&smax); 
  if (smax - smin == 0.0) smax = smin + 1.0;
  cout << "Byte Min, Max: %f , %f\n",smin,smax;

  byte_file(argv[2],smin,smax,stval,&head);  // create a byte format file
  
  if (argc > 3) {   // write SIR format file
    cout << "Writing output SIR file" << endl;

    int ierr = PutSIR(argv[3], &head, stval);
    if (ierr < 0) {
      cout << "*** ERROR writing output SIR file ***" << endl;
      return -1;
    }
  }

  return 0;
}



void byte_file(char *fname, float smin, float smax, float *stval, 
	       struct sirheader *head)
{
  FILE *omf;
  unsigned char *data;             // pointer to byte array storage */
  int i, n = head->nsx * head->nsy;  // pixels
  float am, scale = (smax - smin);

  omf = fopen(fname,"w"); 
  if (omf == NULL) {
    cout << "ERROR: cannot open output byte file\n"; 
    return;
  }
  cout << "Generating byte file...\n";

  // allocate space for byte array image and create it
  data = new unsigned char [n];
  if (data == NULL) {
    cout << "*** ERROR: byte memory allocation failure...\n";
    return;
  }

  if (scale > 0.) scale = 255./ scale;

  for (i = 0; i < n; i++) {
    am = scale * (*(stval+i) - smin);    // scale floating point to byte values
    if (am > 255.) am = 255.;		// check overflow
    if (am < 0.) am = 0.;	        // check underflow
    *(data+i) = (unsigned char)((int)(am));    // byte array
  }

  fwrite(data, sizeof(unsigned char), n, omf);
  fclose(omf);

  cout << "BYTE output file successfully written\n";
}