/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *  
*  nmcraobs.c
*
*  sample read program for nmc rawinsonde observation data
*
*  only categories 51 (land surface data) and 8 (additional data)
*  are expected in the KuDA land rawinsonde data set.
*
*  Copyright 1993, University Corporation for
*  Atmospheric Research, All Rights Reserved
*
*  ADM (NCAR/ATD/RAF) 93 FEB 26
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */ 

/* ---- compiler #includes ---- */
#include <fcntl.h>
#include <unistd.h>
#include <stdio.h> 


#define MXLV8 16
#define IDD 50
#define MX /*128*/ 32

#define TRUE 1
#define FALSE 0
#define OK 1
#define NOTOK 0

/* forward function declarations */
int stkadp(), trp_aft_adp(), wzzadp();
int is_quit();

/* globals */
float
    p[MX],z[MX],h[MX],d[MX],
    f[MX],t[MX];

float pm[/*20*/] = {
      1000.0, 850.0, 700.0, 500.0, 400.0, 300.0, 250.0,
      200.0, 150.0, 100.0, 70.0, 50.0, 30.0, 20.0, 10.0,
      7.0, 5.0, 3.0, 2.0, 1.0};

char
    nmsg[/*40*/]="9999999999999999999999999999999999999999",
    nbf[10000],
    q[MX][5],
    asta[6];

int
  istat=0,
  nu=99999,
  na=0,
  lthr=0;

/* ****** local global data ****** */
  static int
      infh,                      /* input file handle */
      outfh;                     /* output file handle */

/*
*-------------------------------------------------------------------
*
* main program
*
*   program adp77
*
*  asks user for filename and date, searches file for record
*  of that date, displays data, and asks for next record
*-------------------------------------------------------------------
*/
main()
{

    char
        ic25, ic26, ic27,
        flnm[128],
        ymd[20];

    int nro, nr, idsrch, idstyp, irtyp, iintyp, lth, nlvw,
        cc,   /* condition code */
        icc,
        i,    /* general purpose index */
        year,
        month,
        day,
        notfnd;
int iyr, imo, idy, ihr;
int nlv;
    float
        xlat, xlon, xhr, xelev;



/* --- set no matching records found flag */
    notfnd = TRUE;

/* --- set condition code to all's well */
    icc=0;


    do { /* infinite loop */
      bzero(q, sizeof(q));

      /* - - - - - - - - - - - - - - - - - - -*/
      do { /* while the input file is not open or does not exist */
         /* 
          *  ask user for desired file
          */
         bzero(flnm, sizeof(flnm));
         printf("\n");

         printf("  Enter file name to read or 'quit' -> ");
         gets(flnm);

         if ((i=is_quit(flnm))==TRUE) {
             printf("exiting program\n");
             exit();
         }

         /*
          *  -- open input file
          */

        if ((infh = open (flnm,O_RDONLY,0)) == -1) {
            printf(" cannot find file %s\n", flnm);
            printf(" please select another file.\n");
        } /* end if fd */
      } while (infh == -1);
      /* - - - - - - - - - - - - - - - - - - -*/



      nro=0;

      /* --- now enter yymmdd */
      icc=0;
      do {
        cc=OK;
        printf("\n");
        printf("  Enter YYMMDD or 'quit' -> ");
        gets(ymd);

        if ((i=is_quit(ymd))==TRUE) {
            printf("exiting program\n");
            exit();
        }
        i=
  sscanf(ymd, "%2d%2d%2d", &year, &month, &day);
        if (i==EOF || i!=3) {
            continue;
        }
        if (year < 90) {
            printf(" error:  year must be greater than 90\n");
            cc=NOTOK;
        }
        if (month < 1  ||  month > 12) {
            printf(" error:  month not between 1 and 12\n");
            cc=NOTOK;
        }
        if (day < 1  ||  day > 31) {
            printf(" error:  day not between 1 and 31\n");
            cc=NOTOK;
        }
      } while(cc==NOTOK);
      printf("\n  SEARCHING INPUT FILE\n\n");
          

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/*   read a record (end of record marked by carraige return)     */

      nr=0;
/*   70 continue */
   do { /* to read all records for a given day */
   do {
       istat = rdadp(infh,&idsrch,&iyr,&imo,&idy,&ihr,&idstyp,&xlat,&xlon,
           &xhr,&ic25,&ic26,&ic27,&irtyp,&xelev,&iintyp,&lth);
       icc = istat;
        if (icc == 0) {
          nr=nr+1;
          nro=nro+1;
    
          if (iyr > year) icc = 2;
          else if (imo > month  &&  iyr == year) icc = 2;
          else if (idy > day  &&  imo == month  &&  iyr == year)
              icc = 2;
       } /* end if(icc==0) */
   } while (icc==0 && (iyr != year  ||  imo != month  ||  idy != day));

/* ---------------------------------------------------------- */
   if (icc == 0) {
      printf(" ID - STATION %6s%10d\n",asta,1);
      printf(" %3d%3d%3d%5d%6.2f %7.2f%7.2f%6.0f %4d\n",
          iyr,imo,idy,ihr,xhr,xlat,xlon,xelev,irtyp);

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      nlv = stkadp(MX);
      if (istat == 0 && nlv  >=  1) {
          printf(" MAN AND SIG LEVELS %10d\n", nlv);
          for (i=0; i<nlv; i++) {
            printf(" %4d%8.1f%8.1f%8.1f%8.1f%8.1f%8.1f %4s\n", 
                i+1,p[i],z[i],t[i],h[i],d[i],f[i],q[i]);
          }
      } /* end if (nlv >= 1) */

/*
 - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     wind by height
*/
      if (istat == 0) nlvw = wzzadp(MX);
      if (istat == 0 && nlvw >= 1) {
          printf(" WIND BY HEIGHT     %10d\n", nlv);
          for (i=0; i<nlvw; i++) {
            printf(" %4d%8.1f%8.1f%8.1f %4s\n", 
                i+1,z[i],d[i],f[i],q[i]);
          }
      } /* end if (nlvwz >= 1) */

/*
 - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     tropopause
*/
      if (istat == 0) nlv=trp_aft_adp(MX,5);
      if (istat == 0 && nlv >= 1) {
          printf(" TROPOPAUSE         %10d\n", nlv);
          for (i=0; i<nlv; i++) {
            printf(" %4d%8.1f%8.1f%8.1f%8.1f%8.1f %4s\n", 
                i+1,p[i],t[i],h[i],d[i],f[i],q[i]);
          }
      } /* end if (nlv >= 1) */

/*
 - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     aircraft reports
*/
      if (istat == 0) nlv=trp_aft_adp(MX,0);
      if (istat == 0 && nlv  >=  1) {
          printf(" AIRCRAFT REPORT    %10d\n", nlv);
          for(i=0; i<nlv; i++) {
              printf(" %4d%8.1f%8.1f%8.1f%8.1f%8.1f %4s\n", 
                  i+1,p[i],t[i],h[i],d[i],f[i],q[i]);
          }
      } /* end if (nlv >= 1) */

/* ---------------------------------------------------------- */
    /* get user to hit a key for next reading */
         printf("\n");
         if (istat != 5) {
         printf("  Press return to read next record or enter 'quit' -> ");
           gets(ymd);

           if ((i=is_quit(ymd))==TRUE) icc=4;
           else printf("\n\n");
         } /* end if (istat!=50) */
    } /* end if (icc == 0) */
    } /* end do while */ while (icc == 0);

/* close the input file - a different one may be opened next */
      close(infh);
      icc=istat;
      switch(icc) {
        case 0:
            break;
        case 2:
            if (notfnd) {
                printf("no records found for date entered\n");
            }
            break;
        case 4:
            /* user entered quit so stop program */
            printf("exiting program\n");
            exit();
            break;
        case 5:
            /* end of file */
            printf(" end of file reached - # of records=%d\n",nro);
            break;
        case 6:
            /* record error */
            printf(" record %d error - reading next record\n",nro);
            break;
        default:
            /* unknown system error on open or read - stop program */
            printf("unrecoverable error\n");
            exit();
            break;
      } /* end switch(icc) */
    } while(TRUE); /* end infinite loop */

} /* end main() */

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


/*
c--------------------------------------------------------------------
c
c   subroutine rdadp
c
c   reads 60 character header info:
c       10 character yymmddhh block
c       40 character report identification block (lat,lon)
c--------------------------------------------------------------------
*/
int
rdadp(fh,idsrch,iyr,imo,idy,ihr,idstyp,xlat,xlon,
    xhr,ic25,ic26,ic27,irtyp,xelev,iintyp,lth)
int fh, *idsrch, *iyr, *imo, *idy, *ihr, *idstyp, *irtyp,*lth;
float *xlat, *xlon, *xhr, *xelev;
char  *ic25, *ic26, *ic27;
{
  int ist, i, 
    nbytes,
    lat, lon, hr,elev;


  nu=nu+lthr;

  do {
      ist=1;
      i=0;
      /* read a record */
      do {
          nbytes = read(fh,&nbf[i],1);
          if (nbytes == 0) {
              istat=5;
              return(istat);
          }
          else if (nbytes < 0) {
              printf(" err in rdadp- nu=%d\n", nu);
              istat=5;
              return(istat);
          }
      } while(nbf[i++] != '\n');

      ist=0;
      nu=1;

/* --- the first 10 characters of each record contain the year, */
/* --- month, day, and hour (2x,4i2) of each observation        */
/* --- get date from leader */

      i=sscanf(nbf, "  %2d%2d%2d%2d", iyr,imo,idy,ihr);
      if (i==EOF) {
          istat=5;
          return(istat);
      }
      if (i!=4) {
          istat=6;
          return(istat);
      }

/* --- get id for current report */
      nu=11;

/* --- the next 40 characters of each record contain the latitude, */
/* --- longitude, etc. */

      i=sscanf(&nbf[nu-1], "%5d%5d%6s%4d%*4s%c%c%c%3d%5d%2d%3d",
          &lat,&lon,asta,&hr,
          ic25,ic26,ic27,irtyp,&elev,iintyp,lth );
      if (i==EOF) {
          istat=5;
          return(istat);
      }
      if (i!=11) {
          istat=6;
          return(istat);
      }
      *xlat=lat*0.01;
      *xlon=lon*0.01;
      *xhr=hr*0.01;
      *xelev=elev*1.0;

/**/
      lthr= (*lth) * 10;
      istat=0;
      if(lthr > 40) return(istat);

/* --- error encountered */

      printf(" err in rdadp- nu=%d\n", nu);
  } while(TRUE);
    /*  go to 5 */

/* -----end of file reached */
  istat=5;
  return(istat);
} /* end rdadp() */

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

/*
c--------------------------------------------------------------------
c  mandatory? pressure levels
c--------------------------------------------------------------------
*/
int
manadp(maxlev)
int maxlev;
{
   int i,j;
   int ncnt, ntu, nlv, ntn, categ, nent, categc;
   int zz, tt, hh, dd, ff;

      ntu=nu+40;
      ncnt=0;
      nlv=0;
/*
 * look for category  1
 */
    do {
      i=strncmp(&nbf[ntu-1],"END REPORT",10);
      if(i == 0) {
          ncnt=nlv;
          return(ncnt);
      }
      i = sscanf(&nbf[ntu-1], "%2d%3d%2d%3d", &categ,&ntn,&nent,&categc);
      if (i==EOF) {
          istat=5;
          printf(" decode error in cat01\n");
          return(ncnt);
      }
      if (i!=4) {
          istat=6;
          return(ncnt);
      }
      if(ntn <= 0) {
          printf(" decode error in cat01\n");
          istat=5;
          return(ncnt);
      }
      if(categ !=  1) ntu=nu+10*(ntn-1);
    } while (categ != 1);

    nlv=nent;
    if (nlv  >  maxlev) {
        printf(" too many levels in manadp\n",nlv);
        nlv=maxlev;
    }
    ntu=ntu+10;

    for(i=0;i<nlv;i++) {
          p[i]=pm[i];
          j=sscanf(&nbf[ntu-1], "%5d%4d%3d%3d%3d%c%c%c%c",
              &zz,&tt,&hh,&dd,&ff,&q[i][0],&q[i][1],&q[i][2],&q[i][3]);
          if (j==EOF) {
              istat=5;
              printf(" decode error in cat01\n");
              return(ncnt);
          }
          if (j!=9) {
              istat=6;
              printf(" decode error in cat01\n");
              return(ncnt);
          }
          z[i] = zz * 1.0;
          t[i] = tt * 0.1;
          h[i] = hh * 0.1;
          d[i] = dd * 1.0;
          f[i] = ff * 1.0;

          ntu=ntu+22;
    }

   ncnt=nlv;
   return(ncnt);
} /* end manadp() */

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

/*
c--------------------------------------------------------------------
c  significant pressure levels
c--------------------------------------------------------------------
*/      
int
sigadp(kk, maxlev)
int kk, maxlev;
{

  int i,j;

  int ncnt, ntu, nlv, categ, ntn, nent, categc;
  int pp, tt, hh;

  ntu=nu+40;
  ncnt=0;
  nlv=0;
/*
 look for category  2
*/
    do {
      if(i=(strncmp(&nbf[ntu-1], "END REPORT",9)) == 0) {
          ncnt=nlv;
          return(ncnt);
      }
      j=sscanf(&nbf[ntu-1], "%2d%3d%2d%3d", &categ, &ntn, &nent, &categc);
      if (j==EOF) {
          printf(" decode error in cat02\n");
          ncnt=-1;
          istat=5;
          return(ncnt);
      }
      if ( j != 4) {
          printf(" decode error in cat02\n");
          ncnt=-1;
          istat=6;
          return(ncnt);
      }
      if(ntn <= 0) { /*go to 95*/
          printf(" decode error in cat02\n");
          ncnt=-1;
          return(ncnt);
      }
      if(categ !=  2) ntu=nu+10*(ntn-1);
    } while(categ != 2);

      nlv=nent;
      if (nlv  >  maxlev) {
          printf(" too many levels in sigadp%10d\n",nlv);
          nlv=maxlev;
      }
      ntu=ntu+10;

      for (i=kk-1; i<(nlv+kk-1); i++) {
          j=sscanf(&nbf[ntu-1], "%5d%4d%3d%c%c%c",
            &pp, &tt, &hh, &q[i][0], &q[i][1], &q[i][2]);
          q[i][3]=' ';
          if (j==EOF) {
              printf(" decode error in cat02\n");
              ncnt=-1;
              istat=5;
              return(ncnt);
          }
          if (j!=6) {
              printf(" decode error in cat02\n");
              ncnt=-1;
              istat=6;
              return(ncnt);
          }
          p[i] = pp * 0.1;
          t[i] = tt * 0.1;
          h[i] = hh * 0.1;
          ntu=ntu+15;
      }

    ncnt=nlv;
    return(ncnt);
} /* end sigadp() */


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

/*
c--------------------------------------------------------------------
c
c   subroutine stkadp
c
c--------------------------------------------------------------------
*/
int
stkadp(maxlev)
int maxlev;
{
  int i,j,k;

  int nlv, nmlv, mlev, nslv, il, is, flag, ib, ie, l;
  char qn[5];
  float xp, xz, xt, xh,xd, xf, xmul;
  nlv=0;
  nmlv=manadp(maxlev);
  if (istat!=0) return(-1);
  j=-1;

  if (nmlv  >  0) {
    for (i=0; i<nmlv; i++) {
      if (t[i] > 90.0  &&  h[i] > 90.0 && z[i] > 90000.0  && 
          d[i] > 900.0) continue;
      j=j+1;
      if (j  !=  i) {
          p[j]=p[i];
          z[j]=z[i];
          t[j]=t[i];
          h[j]=h[i];
          d[j]=d[i];
          f[j]=f[i];
          strncpy(q[j],q[i],5);
      }
    }
    nmlv=j+1;
  }


  mlev=maxlev-nmlv;
  i=nmlv+1;
  nslv=sigadp(i, mlev);
  if (istat!=0) return(-1);

  nlv=nmlv+nslv;
  if (nslv  ==  0) return(nlv);
  for (j=i-1; j<nlv; j++) {
      d[j]=999.0;
      f[j]=999.0;
      z[j]=99999.0;
  }

  if (nmlv == 0) return(nlv);

  il=nmlv;

  do {
      il=il+1;

      flag=TRUE;
      do {
          if (il  >  nlv) return(nlv);
          is=il;

          do {
              is=is-1;
              if (is  <  1) {
                  flag=FALSE;
                  break;
              }
          } while ( p[il-1]  >  p[is-1] );
          if (flag==FALSE) break;

          if (p[il-1]  !=  p[is-1]) {
              flag=FALSE;
              break;
          }

          nlv=nlv-1;
          if (il  >  nlv) return(nlv);

          for (k=il-1; k<nlv; k++) {
              j=k+1;
              p[k]=p[j];
              z[k]=z[j];
              t[k]=t[j];
              h[k]=h[j];
              d[k]=d[j];
              f[k]=f[j];
              strncpy(q[k],q[j],4); /*q[k]=q[j]; */
          }
      } while(flag==TRUE);

      is=is+1;
      xp=p[il-1];
      xz=p[il-1];
      xz=z[il-1];
      xt=t[il-1];
      xh=h[il-1];
      xd=d[il-1];
      xf=f[il-1];
      strncpy(qn,q[il-1],4); /*qn=q[il];*/
      ib=1;
      ie=il-is;

      for (j=ib; j<=ie; j++) {
          k=il-j-1;
          l=k+1;
          p[l]=p[k];
          z[l]=z[k];
          t[l]=t[k];
          h[l]=h[k];
          d[l]=d[k];
          f[l]=f[k];
          strncpy(q[l],q[k],4); /*q[l]=q[k];*/
      }

      p[is-1]=xp;
      z[is-1]=xz;
      t[is-1]=xt;
      h[is-1]=xh;
      d[is-1]=xd;
      f[is-1]=xf;
      strncpy(q[is-1],qn,4); /*q[is]=qn;*/
  } while (TRUE);
} /* end stkadp() */

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

/*
c--------------------------------------------------------------------
c
c--------------------------------------------------------------------
*/
int
trp_aft_adp(maxlev, lcat)
int maxlev, lcat;
{
  int i,j;
  int  ncnt, ntu, nlv, categ, ntn, nent, categc;
  int pp,tt,hh,dd,ff;
  float xmul;

  if (lcat==5) /*trpadp*/ xmul=0.1;
  else /*aftadp*/ xmul=1.0;

  ntu=nu+40;
  ncnt=0;
  nlv=0;
/*
  look for category  lcat
*/
  do {
      if(strncmp(&nbf[ntu-1], "END REPORT", 10) == 0) {
          ncnt=nlv;
          return(ncnt);
      }
      j=sscanf(&nbf[ntu-1], "%2d%3d%2d%3d", &categ, &ntn, &nent, &categc);
      if (j==EOF) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=5;
          return(ncnt);
      }
      if (j != 4) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=6;
          return(ncnt);
      }
      if(ntn <= 0) { /*go to 95*/
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          return(ncnt);
      }
      if(categ !=  lcat) ntu=nu+10*(ntn-1);
  } while (categ != lcat);

  nlv=nent;
  if (nlv  >  maxlev) {
      printf(" too many levels in txxadp%10d",nlv);
      nlv=maxlev;
  }

      ntu=ntu+10;
  for(i=0; i<nlv; i++) {
      j = sscanf(&nbf[ntu-1], "%5d%4d%3d%3d%3d%c%c%c%c",
          &pp,&tt,&hh,&dd,&ff,&q[i][0],&q[i][1],&q[i][2],&q[i][3]);

      if (j==EOF) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=5;
          return(ncnt);
      }
      if (j != 9) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=6;
          return(ncnt);
      }
      p[i] = pp * 1.0;
      t[i] = tt * 0.1;
      h[i] = hh * 0.1;
      d[i] = dd * 1.0;
      f[i] = ff * 1.0;
      p[i]=p[i]*xmul;
      ntu=ntu+22;
  }
  ncnt=nlv;
  return(ncnt);
} /* end trp_aft_adp() */
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */

/* 
c--------------------------------------------------------------------
c
c   subroutine wzzadp
c
c   reads category 4 data
c--------------------------------------------------------------------
*/
int
wzzadp(maxlev)
int maxlev;
{
  int i,j;

  int
     ncnt, lcat, ntu, nlv, categ, ntn, nent, categc;
  int zz,dd,ff;

  lcat=4;
  ntu=nu+40;
  ncnt=0;
  nlv=0;

/*
 * look for category  lcat
 */
  do {
      if (strncmp(&nbf[ntu-1], "END REPORT", 10) == 0) {
          ncnt=nlv;
          return(ncnt);
      }
      j=sscanf(&nbf[ntu-1], "%2d%3d%2d%3d", &categ,&ntn,&nent,&categc);
      if (j==EOF) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=5;
          return(ncnt);
      }
      if (j!=4) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=6;
          return(ncnt);
      }
      if(ntn <= 0) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          return(ncnt);
      }
      if(categ !=  lcat) ntu=nu+10*(ntn-1);
  } while(categ != lcat);
/*
*/
  printf("%2d ",lcat);
  nlv=nent;
  if (nlv  >  maxlev) {
      printf(" too many levels in wxxadp%10d\n",nlv);
      nlv=maxlev;
  }

  ntu=ntu+10;

  for(i=0; i<nlv; i++) {
      j=sscanf(&nbf[ntu-1], "%5d%3d%3d%c%c", 
          &zz,&dd,&ff,&q[i][0],&q[i][1]);
      q[i][2]=q[i][3]=' ';

      if (j==EOF) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=5;
          return(ncnt);
      }
      if (j!=5) {
          printf(" decode error in cat%2d\n",lcat);
          ncnt=-1;
          istat=6;
          return(ncnt);
      }
      z[i] = zz * 1.0;
      d[i] = dd * 1.0;
      f[i] = ff * 1.0;
      ntu=ntu+13;
  }

/* ---- 'END REPORT found' */
  ncnt=nlv;
  return(ncnt);
} /* end wzzadp() */

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 * is_quit()
 *
 * returns FALSE if string is not "quit"
 * else returns TRUE
 */
int
is_quit(ptr)
char *ptr;
{
   int i;
   char *c;
   c=ptr;

   if (isalpha(*c)) {

       for(i=0;i<4;i++,c++) *c=toupper(*c);

       if (strncmp(ptr,"QUIT",4) == 0) return(TRUE);
   }
   return(FALSE);
} /* end is_quit() */
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 
 *        end of file nmcraobs.c                                   *
 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */