/* Convert JCAMP-DX spectrum to LISE
*/

#include <stdio.h>
#include <spec.h>

#define TLN10 0.434294482
#define LN10   2.30258509

float x,y,*spc,*err,*tim;
float yfactor, xfactor;
char EOL[4];

help()
{
  printf("jcampio spectrum [-l2j] [-j2l] [-nocr]\n");
  printf("  converts spectra from LISE to JCAMP and vice versa\n");
  printf("  options:\n");
  printf("     -l2j     convert from LISE to JCAMP\n");
  printf("     -j2l     convert from JCAMP to LISE\n");
  printf("     -yf fpnr divide Y-values by factor on jcamp output\n");
  printf("     -xf fpnr divide X-values by factor on jcamp output\n");
  printf("     -nocr    suppress CR before LF\n");
  return(0);
}

main(argc,argv)
int argc;
char *argv[];
{
int n,i,max,memsiz;
char z[80],comment[80];

   yfactor = 0.0;
   xfactor = 0.0;

   memsiz = (_MAXSPCLEN+2) * sizeof(float);
   spc= (float *)malloc(memsiz);
   err= (float *)malloc(memsiz);
   tim= (float *)malloc(memsiz);
   if(tim==NULL) {
      printf("sorry, not enough memory\n");
      exit(-1);
   }
   EOL[0]=13; EOL[1]=10; EOL[2]=0;
   if(checkopt(argc,argv,"-nocr",z)) {
      EOL[0]=10; EOL[1]=0;
   }
   if(checkopt(argc,argv,"-yf",z)) {
      yfactor = atosf(z);
   }
   if(checkopt(argc,argv,"-xf",z)) {
      xfactor = atosf(z);
   }
   if(checkopt(argc,argv,"-l2j",z)) {
      max=readspec(argv[1],spc,err,tim,comment);
      strcpy(z,argv[1]); n = strlen(z);
      for(i = 0; i < n; i++) if(z[i] == '.') z[i] = 0;
      writejcamp(z,spc,err,_utim,max,comment);  /* we use "_utim" as a Bugfix for the SAS-C */
   }

   if(checkopt(argc,argv,"-j2l",z)) {
      max = readjcamp(argv[1],spc,err,tim,comment);
      strcpy(z,argv[1]); n = strlen(z);
      for(i = 0; i < n; i++) if(z[i] == '.') z[i] = 0;
      writespec(z,spc,err,max,2,comment);
      if(_spc_onam[0] != 0) strcpy(z,_spc_onam);
      strcat(z,".tim");
      writespec(z,tim,err,max,2,comment);
   }
   free(err); free(spc); free(tim);
   exit(0);
}

float getnumber(s)
char *s;
{
char z[80];
int i,n,l;
float erg;

   n = 0; l = strlen(s);
   erg = 0.0;
   while(s[n++] != '=') if(n >= l) return(0.0);
   i = 0; while(s[n] > 31) z[i++] = s[n++];
   z[i] = 0;
   erg = atosf(z);
   return(erg);
}

float calcfactor(max,min)
float max,min;
{
float y;
int imax,imin;
int i;

   if(max < 0.0) max = -1.0 * max;
   if(min < 0.0) min = -1.0 * min;

   imax = (int) (log(max) * TLN10);
   imin = (int) (log(min) * TLN10);
   i = imax;
   if((i > 60) || (i < -60)) i = imin;
   if((i > 60) || (i < -60)) return(1.0);
   y = exp(LN10 * i);
   y = 1.0 / (y * 10000.0);
   return(1.0); /* makes some trouble */
}

writejcamp(name,spc,err,tim,max,comment)
char *name, *comment;
float *spc, *err, *tim;
int max;
{
FILE *fp;
int i,n,m;
char *s, *z;
float maxy, miny;
float maxx, minx;

   s = (char *) malloc(256);
   z = (char *) malloc(256);

   sprintf(s,"%s.jcamp",name);
   fp = fopen(s,"w");
   if(fp == NULL) {
      fprintf(stderr,"jcampio: could not open >%s< for write\n",s);
      free(s);
      return(0);
   }

   maxy = spc[0]; miny = spc[0];
   maxx = tim[0]; minx = tim[0];
   for(n = 0; n < max; n++) {
      if(spc[n] > maxy) maxy = spc[n];
      if(spc[n] < miny) miny = spc[n];
      if(tim[n] > maxx) maxx = tim[n];
      if(tim[n] < minx) minx = tim[n];
   }

   if(xfactor == 0.0) xfactor = calcfactor(maxx,minx);
   if(yfactor == 0.0) yfactor = calcfactor(maxy,miny);

   fprintf(fp,"##TITLE=%s%s",comment,EOL);
   fprintf(fp,"##JCAMP-DX=4.24%s",EOL);
   fprintf(fp,"##DATA TYPE= LISE%s",EOL);
   fprintf(fp,"##ORIGIN= (C) 1993 by Rainer Kowallik%s",EOL);
   fprintf(fp,"##OWNER= Public Domain%s",EOL);
   fprintf(fp,"##XUNITS= arbitrary%s",EOL);
   fprintf(fp,"##YUNITS= arbitrary%s",EOL);
   fpfmt(s,xfactor);
   fprintf(fp,"##XFACTOR=%s%s",s,EOL);
   fpfmt(s,yfactor);
   fprintf(fp,"##YFACTOR=%s%s",s,EOL);
   fpfmt(s,tim[0]);
   fprintf(fp,"##FIRSTX=%s%s",s,EOL);
   fpfmt(s,tim[max-1]);
   fprintf(fp,"##LASTX=%s%s",s,EOL);
   fprintf(fp,"##NPOINTS=%d%s",max,EOL);
   fpfmt(s,spc[0]);
   fprintf(fp,"##FIRSTY=%s%s",s,EOL);
   fpfmt(s,maxy);
   fprintf(fp,"##MAXY=%s%s",s,EOL);
   fpfmt(s,miny);
   fprintf(fp,"##MINY=%s%s",s,EOL);
   fpfmt(s,_tica);
   fprintf(fp,"##DELTAX=%s%s",s,EOL);
   fprintf(fp,"##XYDATA=(X++(Y..Y))%s",EOL);
   n = 0;
   while(n < max) {
      sprintf(z,"%d",(int)(tim[n] / xfactor));
      while(strlen(z) < 70) {
         sprintf(s,"%+d",(int)(spc[n] / yfactor));
         strcat(z,s);
         if((n++) >= max) break;
      }
      fprintf(fp,"%s%s",z,EOL);
   }
   fprintf(fp,"##END=%s",EOL);
   fclose(fp);

   free(z); free(s);
}

fpfmt(s,x)
char *s;
float x;
{
double y;
int i;

   y = (double) x;
   i = (int) x;
   sprintf(s,"%g",y);
   return(0);
}

readjcamp(name,spc,err,tim,comment)
char *name, *comment;
float *spc, *err, *tim;
{
FILE *fp;
int i,n,l,max;
char *s, *z;
float firstx, deltax, yfactor;

   s = (char *) malloc(256); z = (char *) malloc(256);

   firstx = 0.0; deltax = 1.0;
   yfactor = 1.0;
   strcpy(comment,"no comment");

   strcpy(s,name);
   fp = fopen(s,"r");
   if(fp == NULL) {
      strcat(s,".jcamp");
      fp = fopen(s,"r");
      if(fp == NULL) {
         fprintf(stderr,"jcampio: could not open >%s< for read\n",name);
         free(s); free(z);
         return(0);
      }
   }
   while(!feof(fp)) {
      fgets(s,200,fp);
      if(instr("TITLE",s) > 0) {
         n = 0; i = 0;
         while(s[i++] != '=') if(i > 100) break;
         while(s[i] > 30) comment[n++] = s[i++];
         comment[n] = 0;
      }
      if(instr("FIRSTX",s) > 0) firstx = getnumber(s);
      if(instr("DELTAX",s) > 0) {
         deltax = getnumber(s);
         _tica = deltax;
      }
      if(instr("YFACTOR",s) > 0) yfactor = getnumber(s);
      if(instr("XYDATA",s) > 0) break;
   }
   max = 0;

   while(!feof(fp)) {
      fgets(s,200,fp);
      n = 0; while(s[n] > 20) n++;
      s[n] = 0;
      if(s[0] == '#') continue;
      if(strlen(s) < 3) break;
      if(instr("END",s) > 0) break;
      n = 0; l = strlen(s);
      while(s[n] < 42) n++;           /* skip blanks and sign of X*/
      while(s[n++] > 45) if(n >= l) break; /* first skip X */
      while( 1 == 1) {
         i = 0; if(n >= l) break;
         n = n - 1;
         while(s[n] == 32) n++;           /* skip blanks */
         z[i++] = s[n++];                 /* read sign */
         while(s[n] > 45) z[i++] = s[n++];
         z[i] = 0; n++;
         spc[max] = ((float) atoi(z)) * yfactor;
         err[max] = 0.001 * spc[max]; /* you'll need an error for the fit */
         tim[max] = firstx + deltax * ((float) max);
         max = max + 1;
         if(max >= _MAXSPCLEN) {
            printf("Warning, spectrum is not fully read, INCREASE MAXSPCLEN\n");
            fclose(fp); free(s); free(z); return(max);
         }
      }
   }
   fclose(fp);
   free(s); free(z);
   return(max);
}