[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

gpfront update



Jim,
attatched is the updated read_front.c to fix the 
new time strings and flip the pip side.
If they decide to change the side of the warm front pips
again, just change the imul for fnum==2 back to .5:

else /* this draws a semi-circle pip */
      {
      imul = .5; if(fnum == 5) imul = -.5; /* occluded pip on top side */
      if(fnum == 2) imul = -.5; /* warm fronts sitched sides too */

At this point, the occluded and stationary fronts look ok with the warm
front pips on the correct side. Only the warm front orientation has
changed.

Chiz
/*********************************************************
** Collection of subroutines to read a frontal file,    **
** parse the symbols to be plotted and plot them on a   **
** map.                                                 **
**                                                      **
** S. Chiswell Unidata   3/98                           **
*********************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <fcntl.h>
#include <unistd.h>
#include <string.h>
#include <math.h>

#define MAXP    100
#define SOH     1
#define ETX     3

#ifdef UNDERSCORE
#define read_front      read_front_
#define gtext           gtext_
#define gscolr          gscolr_
#define gsline          gsline_
#define gline           gline_
#define gfill           gfill_
#define gspcl           gspcl_
#define gsspcl          gsspcl_
#define gflush          gflush_
#endif

#ifdef LINUX
#define read_front_     read_front__
#endif

typedef struct features {
   int npoints;
   float lat[MAXP],lon[MAXP];
   float press[MAXP];
   struct features *next;
   } features;

typedef struct wxbull {
   char *valid;
   struct features *low, *high, *warm, *cold, *occluded, *trof, *stationary;
   } wxbull;


/* read a bulletin from a file */
char *getbull(fp,offset)
int fp;
int *offset;
{
struct stat buf;
int ier,FOUND;
char ch;
ssize_t nbyte;
int start,stop,i;
char *rbuf=NULL;

if((ier = fstat(fp,&buf)) != 0)
   {
   perror("Could not stat file\0");
   return(rbuf);
   }

if(*offset >= buf.st_size)
   {
   printf("offset past end of file\n");
   return(rbuf);
   }

FOUND = 0;
while((FOUND == 0)&&(*offset < buf.st_size))
   {
   lseek(fp, *offset, SEEK_SET);
   nbyte = read(fp, &ch, 1);
   if(ch == SOH)
      {
      start = *offset;
      FOUND = 1;
      }
   *offset += 1;
   }

while((FOUND == 1)&&(*offset < buf.st_size))
   {
   lseek(fp, *offset, SEEK_SET);
   nbyte = read(fp, &ch, 1);
   if(ch == ETX)
      {
      stop = *offset;
      FOUND = 2;
      }
   *offset += 1;
   }

if(FOUND == 2)
   {
   rbuf = (char *)malloc(stop - start);
   lseek(fp, start+1, SEEK_SET);
   nbyte = read(fp, rbuf ,(stop - start - 2));
   for(i=0;i<(stop - start -3);i++)
      if(rbuf[i] < ' ') rbuf[i] = ' '; /* get rid of ctrl and bull chars! */
   rbuf[stop - start - 2] = '\0';
   *offset = stop + 1;
   }

return(rbuf);
}


/* get frontal location coordinates from a string */
void *getfront(pos)
char **pos;
{
struct features *feat_obj;
char *cpos,*subs,*apos;
double dval;

cpos = *pos;

while((cpos[0] != NULL)&&((cpos[0] < '0')||(cpos[0] > '9'))) cpos += 1;

*pos = cpos;
if(cpos[0] == NULL) return(NULL);

while((cpos[0] != NULL)&&(cpos[0] < 'A')) cpos += 1;

subs = (char *)malloc(cpos - *pos);
subs[0] = '\0';
strncat(subs,*pos,cpos - *pos - 1);

*pos = cpos;

feat_obj = (struct features *)malloc(sizeof(struct features));
feat_obj->npoints = 0;
feat_obj->next = NULL;
cpos = subs;
while(cpos != NULL)
   {
   apos = cpos;
   dval = strtod(apos,&cpos);
   if(cpos != apos) 
      {
      if(dval < 10000)
         {
         feat_obj->lat[feat_obj->npoints] = (float)(((int)dval)/100);
         feat_obj->lon[feat_obj->npoints] = (float)(((int)dval)%100);
         }
      else
         {
         feat_obj->lat[feat_obj->npoints] = (float)(((int)dval)/1000);
         feat_obj->lon[feat_obj->npoints] = (float)(((int)dval)%1000);
         }
      if(feat_obj->npoints == 0) /** try doubling first point */
         {                       /** for a spline closer to end */
         feat_obj->npoints += 1;
         feat_obj->lat[1] = feat_obj->lat[0];
         feat_obj->lon[1] = feat_obj->lon[0];
         }
      feat_obj->npoints += 1;
      apos = cpos;
      }
   else
      cpos = NULL;
   }
feat_obj->lat[feat_obj->npoints] = feat_obj->lat[feat_obj->npoints - 1];
feat_obj->lon[feat_obj->npoints] = feat_obj->lon[feat_obj->npoints - 1];
feat_obj->npoints += 1;

return(feat_obj);
}


/* get map features from a list of pressure & latlon points */
void *getfeat(pos)
char *pos;
{
char *cpos;
char *subs;
struct features *feat_obj;
int i;
double dval;

cpos = pos + 5;
while((cpos != NULL)&&(cpos[0] < 'A')) cpos += 1;
subs = (char *)malloc(cpos - pos - 5 + 1);
subs[0] = '\0';

strncat(subs,pos+5,cpos - pos - 5);
feat_obj = (struct features *)malloc(sizeof(struct features));
feat_obj->npoints = 0;
feat_obj->next = NULL;
cpos = subs; i = 0;
while(cpos != NULL)
   {
   pos = cpos;
   dval = strtod(pos,&cpos);
   if(cpos != pos) 
      {
      if(i == 0) 
         {
         feat_obj->press[feat_obj->npoints] = (float)dval;
         }
      else
         {
         if(dval < 10000)
            {
            feat_obj->lat[feat_obj->npoints] = (float)(((int)dval)/100);
            feat_obj->lon[feat_obj->npoints] = (float)(((int)dval)%100);
            }
         else
            {
            feat_obj->lat[feat_obj->npoints] = (float)(((int)dval)/1000);
            feat_obj->lon[feat_obj->npoints] = (float)(((int)dval)%1000);
            }
         feat_obj->npoints += 1;
         }
      i = (i + 1) % 2; 
      pos = cpos;
      }
   else
      cpos = NULL;
   }

return(feat_obj);
}


/* parse out highs, lows, and fronts from a bulletin */
void *read_asus1(prod,fhr)
char *prod;
int fhr;
{
struct wxbull *bull=NULL;
struct features *feat_obj;
char *pos,*cpos,*valid,*bullp;
char searchstr[80];
int i,done,fcst=0;
double dval;

if((bull = (void *)malloc(sizeof(struct wxbull))) == NULL)
   return(NULL);

bull->low = NULL; bull->high = NULL;
bull->warm = NULL; bull->cold = NULL;
bull->occluded = NULL; bull->trof = NULL;
bull->stationary = NULL; bull->valid = NULL;

/* see if we are dealing with a forecast front file */
pos = strstr(prod,"PROG VALID");
if((pos != NULL)&&(fhr < 0)) fhr = 12;
if(pos != NULL)
   {
   fcst = 1;
   searchstr[0] = '\0';
   sprintf(searchstr,"%02dHR PROG VALID\0",fhr);
   pos = strstr(prod,searchstr);
   if(pos == NULL)
      {
      printf("forecast hour %d not found\n",fhr);
      return(NULL);
      }
   else
      {
      bullp = pos;
      pos = strstr(pos+15,searchstr+2);
      if(pos != NULL) bullp[(pos - bullp)-3] = '\0';
      }
   }
else
   bullp = prod;

pos = strstr(bullp,"VALID");
if(pos != NULL)
   {
   valid = (char *)malloc(30);
   valid[0] = '\0';
   if(fcst != 1)
      {
      strncat(valid,"FRONTS ",7);
      strncat(valid,pos,18); 
      /* dumb check for 3 spaces (\r\r\n since NWS keeps screwing up time*/
      cpos = (char *)strstr(valid,"   ");
      if(cpos != NULL) cpos[0] = '\0';
      }
   else
      {
      strncat(valid,"PROG FRONTS ",12);
      strncat(valid,pos,13);
      }
   bull->valid = valid;
   }

pos = strstr(bullp,"HIGHS");
if(pos != NULL) bull->high = getfeat(pos);

pos = strstr(bullp,"LOWS");
if(pos != NULL) bull->low = getfeat(pos);

pos = strstr(bullp,"COLD");
while(pos != NULL) 
   {
   feat_obj = getfront(&pos);
   if(feat_obj != NULL)
      {
      feat_obj->next = bull->cold;
      bull->cold = feat_obj;
      }
   pos = strstr(pos,"COLD");
   }

pos = strstr(bullp,"WARM");
while(pos != NULL) 
   {
   feat_obj = getfront(&pos);
   if(feat_obj != NULL)
      {
      feat_obj->next = bull->warm;
      bull->warm = feat_obj;
      }
   pos = strstr(pos,"WARM");
   }

pos = strstr(bullp,"STNRY");
while(pos != NULL) 
   {
   feat_obj = getfront(&pos);
   if(feat_obj != NULL)
      {
      feat_obj->next = bull->stationary;
      bull->stationary = feat_obj;
      }
   pos = strstr(pos,"STNRY");
   }

pos = strstr(bullp,"TROF");
while(pos != NULL) 
   {
   feat_obj = getfront(&pos);
   if(feat_obj != NULL)
      {
      feat_obj->next = bull->trof;
      bull->trof = feat_obj;
      }
   pos = strstr(pos,"TROF");
   }

pos = strstr(bullp,"OCFNT");
while(pos != NULL) 
   {
   feat_obj = getfront(&pos);
   if(feat_obj != NULL)
      {
      feat_obj->next = bull->occluded;
      bull->occluded = feat_obj;
      }
   pos = strstr(pos,"OCFNT");
   }

return(bull);
}


/* draw frontal pips */
void draw_pip(X,Y,ptype,fnum,color,pipcnt)
float X[],Y[];
int ptype,color,pipcnt;
{
float thta,PIPX[14],PIPY[14];
float imul,dist;
int np,k,iret,pf1,pf2;
char nav[2];

if(ptype >= 0)
   {
   pf1 = 0; pf2 = pipcnt - 1;
   dist = sqrt(pow(X[pf2]-X[pf1],2.0) + pow(Y[pf2]-Y[pf1],2.0));
   if(dist < .3) return; /* don't draw really short segments */

   thta = atan2((Y[pf2] - Y[pf1]),(X[pf2] - X[pf1]));

   if(ptype == 0) /* this draws a triangle pip */
      {
      PIPX[1] = (X[pf2] + X[pf1])/2.0 + 1.0 * cos(thta - (M_PI/2.0));
      PIPY[1] = (Y[pf2] + Y[pf1])/2.0 + 1.0 * sin(thta - (M_PI/2.0));
      PIPX[0] = (X[pf2] + X[pf1])/2.0 + .7 * cos(thta);
      PIPX[2] = (X[pf2] + X[pf1])/2.0 + .7 * cos(thta - M_PI);
      PIPY[0] = (Y[pf2] + Y[pf1])/2.0 + .7 * sin(thta);
      PIPY[2] = (Y[pf2] + Y[pf1])/2.0 + .7 * sin(thta - M_PI);
      np = 3;
      }
   else /* this draws a semi-circle pip */
      {
      imul = .5; if(fnum == 5) imul = -.5; /* occluded pip on top side */
      if(fnum == 2) imul = -.5; /* warm fronts sitched sides too */
      for(k=0;k<14;k++)
         {
         if(pipcnt < 3)
            {
            PIPX[k] = (X[pf1] + X[pf2])/2.0 +
               imul * cos(thta + (M_PI*k)/13.0);
            PIPY[k] = (Y[pf1] + Y[pf2])/2.0 +
               imul * sin(thta + (M_PI*k)/13.0);
            }
         else
            {
            PIPX[k] = X[(pipcnt-1)/2] +
               imul * cos(thta + (M_PI*k)/13.0);
            PIPY[k] = Y[(pipcnt-1)/2] +
               imul * sin(thta + (M_PI*k)/13.0);
            }
         }
      np = 14;
      }
   nav[0] = 'M'; nav[1] = '\0';
   gfill(nav,&np,PIPX,PIPY,&iret,strlen(nav));
   }
}


/* calculate a spline through the frontal points and
   plot the frontal symbols */
void draw_front(fobj,fnum)
struct features *fobj;
int fnum;
{
int ltype,lhw,lwidth,whw;
int color;
int iret,npoints,i,j,k,np,pips,pipcnt=20;
char nav[2];
float X[26],Y[26],Xt,Yt,jx,step=.05,imul;
float B0,B1,B2,B3;
float P0,P1,P2,P3;

switch(fnum)
   {
   case 1:
        color = 4;
        ltype = 1;
        break;
   case 2:
        color = 2;
        ltype = 1;
        break;
   case 3:
        color = 4;
        ltype = 1;
        break;
   case 4:
        color = 5;
        ltype = 25;
        break;
   case 5:
        color = 7;
        ltype = 1;
        break;
   default:
        color = 6;
        ltype = 1;
   }
lhw = 0; lwidth = 2; whw = 0;
gsline(&ltype,&lhw,&lwidth,&whw,&iret);
nav[0] = 'M'; nav[1] = '\0';
if(iret == 0) 
   {
   gscolr(&color,&iret);
   npoints = 2; j = 0;
   if(fobj->npoints > 3) /* use a B-spline */
      {
      for(i=3;i<fobj->npoints;i++)
         {
         np = 0;
         jx = 0;
         while(jx < (1.0 + step))
            {
            B0 = pow((1.0 - jx),3.0) / 6.0;
            B1 = ((3 * pow(jx,3.0)) - (6 * jx * jx) + 4 ) / 6.0;
            B2 = ((-3 * pow(jx,3.0)) + (3 * jx * jx) + (3 * jx) + 1) / 6.0;
            B3 = pow(jx,3.0) / 6.0;
            P0 = fobj->lat[i-3];
            P1 = fobj->lat[i-2];
            P2 = fobj->lat[i-1];
            P3 = fobj->lat[i];
            Xt = P0 * B0 + P1 * B1 + P2 * B2 + P3 * B3;
            P0 = -fobj->lon[i-3];
            P1 = -fobj->lon[i-2];
            P2 = -fobj->lon[i-1];
            P3 = -fobj->lon[i];
            Yt = P0 * B0 + P1 * B1 + P2 * B2 + P3 * B3;
            X[np] = Xt; Y[np] = Yt; np += 1;
            jx = jx + step;
            }
         gline(nav,&np,X,Y,&iret,strlen(nav));
         pips = -1;
         if(fnum == 1) pips = 0;
         if(fnum == 2) pips = 1;
         if(fnum == 3) pips = j;
         if(fnum == 5) pips = j;
         draw_pip(X,Y,pips,fnum,color,pipcnt);
         if(fnum == 3)
            {
            if(j == 0)
               color = 2;
            else
               color = 4;
            j = (j + 1) % 2;
            gscolr(&color,&iret);
            }
         if(fnum == 5)
            j = (j + 1) % 2;
         }
      }
   else
      {
      for(i=0;i<(fobj->npoints - 1);i++)
         {
         X[0] = fobj->lat[i]; X[1] = fobj->lat[i+1];
         Y[0] = -fobj->lon[i]; Y[1] = -fobj->lon[i+1];
         gline(nav,&npoints,X,Y,&iret,strlen(nav));
         pips = -1;
         if(fnum == 1) pips = 0;
         if(fnum == 2) pips = 1;
         if(fnum == 3) pips = j;
         if(fnum == 5) pips = j;
         draw_pip(X,Y,pips,fnum,color,2);
         if(fnum == 3)
            {
            if(j == 0)
               color = 2;
            else
               color = 4;
            j = (j + 1) % 2;
            gscolr(&color,&iret);
            }
         if(fnum == 5)
            j = (j + 1) % 2;
         }
      } /* end no B-spline */
   }
}

/* subroutine to read a file and plot frontal information */
void read_front(wwfil,wwlen,valtim,vallen,atts,attlen,ashr,ashrlen,ier)
char *wwfil;
int *wwlen;
char *valtim;
int *vallen;
char *atts;
int *attlen;
char *ashr;
int *ashrlen;
int *ier;
{

int fp;
int offset,i;
struct wxbull *wxsyms;
struct features *feat_obj;
char *bullp;
char *fname;

char navigate[2],pstr[10],*pos;
int xioff,yioff,iret;
int color,np,width,fhr,labels=1;
float midx,midy,rotate,spcode;
float size,scale;

fname = (char *)malloc(*wwlen + 1);
fname[0] = '\0';
for(i=0;i<*wwlen;i++)
   if(wwfil[i] > ' ') strncat(fname,wwfil+i,1);

if((fp = open(fname,O_RDONLY,NULL)) == -1)
   {
   printf("could not open file: %s\n",fname);
   *ier = -1;
   return;
   }

offset = 0;

atts[*attlen - 1] = '\0';
if((pos = (char *)strchr(atts,'/')) != NULL)
   {
   if((pos[1] == 'N')||(pos[1] == 'n')) labels=0;
   scale = (float)strtod(atts,&pos);
   }
else
   scale = atof(atts);
if(scale <= 0) scale = 1;
if(labels == 0) printf("no labels for extrema\n");

ashr[*ashrlen - 1] = '\0';
if((pos = (char *)strchr(ashr,'/')) != NULL)
   fhr = strtol(ashr,&pos,10);
else
   fhr = atol(ashr);

/* can be generalized for multiple bulletins in a file later
   if I decide that we need to pass a plot time! */
/*while(offset >=0)  
   {*/
bullp = getbull(fp,&offset);
   /*if(bullp == NULL)
      {
      offset = -1;
      continue;
      }*/
wxsyms = (void *)read_asus1(bullp,fhr);

if(wxsyms != NULL)
   {
   navigate[0] = 'M';
   navigate[1] = '\0';
   feat_obj = wxsyms->cold;
   while(feat_obj != NULL)
      {
      draw_front(feat_obj,1);
      feat_obj = feat_obj->next;
      }
   
   feat_obj = wxsyms->warm;
   while(feat_obj != NULL)
      {
      draw_front(feat_obj,2);
      feat_obj = feat_obj->next;
      }
   
   feat_obj = wxsyms->stationary;
   while(feat_obj != NULL)
      {
      draw_front(feat_obj,3);
      feat_obj = feat_obj->next;
      }
   
   feat_obj = wxsyms->trof;
   while(feat_obj != NULL)
      {
      draw_front(feat_obj,4);
      feat_obj = feat_obj->next;
      }
   
   feat_obj = wxsyms->occluded;
   while(feat_obj != NULL)
      {
      draw_front(feat_obj,5);
      feat_obj = feat_obj->next;
      }

   if(wxsyms->high != NULL)
      {
      color = 4;
      gscolr(&color,&iret);
      size = 1.5 * scale; width = 2;
      gsspcl(&size,&width,&iret);
      np = 1; spcode = 2;
      for(i=0;i<wxsyms->high->npoints;i++)
         {
         xioff = 0; yioff = 0; rotate = 0;
         midx = wxsyms->high->lat[i];
         midy = -wxsyms->high->lon[i];
         
gspcl(navigate,&np,&spcode,&midx,&midy,&xioff,&yioff,&iret,strlen(navigate));
         if(labels != 0)
            {
            yioff = (int)(-3*scale); xioff = -3;
            if(wxsyms->high->press[i] < 1000) xioff = -2;
            if(yioff > -2) yioff = -2;
            pstr[0] = '\0'; 
            sprintf(pstr,"%3.0f\0",wxsyms->high->press[i]);
            
gtext(navigate,&midx,&midy,pstr,&rotate,&xioff,&yioff,&iret,strlen(navigate),strlen(pstr));
            }
         }
      }
   if(wxsyms->low != NULL)
      {
      color = 2;
      gscolr(&color,&iret);
      size = 1.5 * scale; width = 2;
      gsspcl(&size,&width,&iret);
      np = 1; spcode = 3;
      for(i=0;i<wxsyms->low->npoints;i++)
         {
         xioff = 0; yioff = 0; rotate = 0;
         midx = wxsyms->low->lat[i];
         midy = -wxsyms->low->lon[i];
         
gspcl(navigate,&np,&spcode,&midx,&midy,&xioff,&yioff,&iret,strlen(navigate));
         if(labels != 0)
            {
            yioff = (int)(-3*scale); xioff = -3;
            if(wxsyms->low->press[i] < 1000) xioff = -2;
            if(yioff > -2) yioff = -2;
            pstr[0] = '\0';
            sprintf(pstr,"%3.0f\0",wxsyms->low->press[i]);
            
gtext(navigate,&midx,&midy,pstr,&rotate,&xioff,&yioff,&iret,strlen(navigate),strlen(pstr));
            }
         }
      }
 
   while(wxsyms->cold != NULL)
      {
      feat_obj = wxsyms->cold->next;
      free(wxsyms->cold);
      wxsyms->cold = feat_obj;
      }
   while(wxsyms->warm != NULL)
      {
      feat_obj = wxsyms->warm->next;
      free(wxsyms->warm);
      wxsyms->warm = feat_obj;
      }
   while(wxsyms->trof != NULL)
      {
      feat_obj = wxsyms->trof->next;
      free(wxsyms->trof);
      wxsyms->trof = feat_obj;
      }
   while(wxsyms->occluded != NULL)
      {
      feat_obj = wxsyms->occluded->next;
      free(wxsyms->occluded);
      wxsyms->occluded = feat_obj;
      }
   while(wxsyms->stationary != NULL)
      {
      feat_obj = wxsyms->stationary->next;
      free(wxsyms->stationary);
      wxsyms->stationary = feat_obj;
      }
   if(wxsyms->high != NULL) free(wxsyms->high); 
   if(wxsyms->low != NULL) free(wxsyms->low); 
   if(wxsyms->valid != NULL) 
      {
      strncpy(valtim,wxsyms->valid,*vallen);
      free(wxsyms->valid);
      }
   free(wxsyms); 
   }
   

   /*}*/
gflush(&iret);
*ier = 0;
return;


}

/* debug stubs */
#ifdef MAIN
main(argc, argv)
int argc;
char *argv[];
{
int ier,wwlen,vallen,attlen,ashrlen;
char valtim[72],atts[72],ashr[72];

wwlen = strlen(argv[1]); vallen=71;
strcpy(atts,"     \0"); attlen = strlen(atts);
strcpy(ashr,"     \0"); ashrlen = strlen(ashr);
read_front(argv[1],&wwlen,valtim,&vallen,atts,&attlen,ashr,&ashrlen,&ier);

}
#endif