g2_addgrid.c

Go to the documentation of this file.

 
#include <stdio.h>
#include <stdlib.h>
#include "grib2.h"
 
g2int
g2_addgrid(unsigned char *cgrib, g2int *igds, g2int *igdstmpl, g2int *ideflist,
           g2int idefnum)
{
    g2int ierr = 0;
    static unsigned char G = 0x47;       /* 'G' */
    static unsigned char R = 0x52;       /* 'R' */
    static unsigned char I = 0x49;       /* 'I' */
    static unsigned char B = 0x42;       /* 'B' */
    static unsigned char seven = 0x37;   /* '7' */
 
    static g2int one = 1, three = 3, miss = 65535;
    g2int lensec3, iofst, ibeg, lencurr, len;
    g2int i, j, temp, ilen, isecnum, nbits;
    gtemplate *mapgrid = 0;
 
    /* Check to see if beginning of GRIB message exists. */
    if (cgrib[0] != G || cgrib[1] != R || cgrib[2] != I || cgrib[3] != B)
    {
        printf("g2_addgrid: GRIB not found in given message.\n");
        printf("g2_addgrid: Call to routine gribcreate required to initialize GRIB messge.\n");
        ierr = -1;
        return(ierr);
    }
 
    /* Get current length of GRIB message. */
    gbit(cgrib, &lencurr, 96, 32);
 
    /* Check to see if GRIB message is already complete. */
    if (cgrib[lencurr - 4] == seven && cgrib[lencurr - 3] == seven &&
        cgrib[lencurr - 2] == seven && cgrib[lencurr - 1] == seven)
    {
        printf("g2_addgrid: GRIB message already complete.  Cannot add new section.\n");
        ierr = -2;
        return(ierr);
    }
 
    /* Loop through all current sections of the GRIB message to find
     * the last section number. */
    len = 16;    /* length of Section 0 */
    for (;;)
    {
        /* Get section number and length of next section. */
        iofst = len * 8;
        gbit(cgrib, &ilen, iofst, 32);
        iofst = iofst + 32;
        gbit(cgrib, &isecnum, iofst, 8);
        len = len + ilen;
        /* Exit loop if last section reached. */
        if (len == lencurr) break;
        /* If byte count for each section doesn't match current total
         * length, then there is a problem. */
        if (len > lencurr)
        {
            printf("g2_addgrid: Section byte counts don''t add to total.\n");
            printf("g2_addgrid: Sum of section byte counts = %ld\n", len);
            printf("g2_addgrid: Total byte count in Section 0 = %ld\n", lencurr);
            ierr = -3;
            return(ierr);
        }
    }
 
    /* Section 3 can only be added after sections 1, 2 and 7. */
    if (isecnum != 1 && isecnum != 2 && isecnum != 7)
    {
        printf("g2_addgrid: Section 3 can only be added after Section 1, 2 or 7.\n");
        printf("g2_addgrid: Section ',isecnum,' was the last found in given GRIB message.\n");
        ierr = -4;
        return(ierr);
    }
 
    /* Add Section 3  - Grid Definition Section. */
    ibeg = lencurr * 8;        /*   Calculate offset for beginning of section 3 */
    iofst = ibeg + 32;         /*   leave space for length of section */
    sbit(cgrib, &three, iofst, 8);     /* Store section number (3) */
    iofst = iofst + 8;
    sbit(cgrib, igds+0, iofst, 8);     /* Store source of Grid def. */
    iofst = iofst + 8;
    sbit(cgrib, igds+1, iofst, 32);    /* Store number of data pts. */
    iofst = iofst + 32;
    sbit(cgrib, igds+2, iofst, 8);     /* Store number of extra octets. */
    iofst = iofst + 8;
    sbit(cgrib, igds+3, iofst, 8);     /* Store interp. of extra octets. */
    iofst = iofst + 8;
 
    /* if Octet 6 is not equal to zero, Grid Definition Template may
     * not be supplied. */
    if (igds[0] == 0)
        sbit(cgrib, igds+4, iofst, 16);  /* Store Grid Def Template num. */
    else
        sbit(cgrib, &miss, iofst, 16);   /* Store missing value as Grid Def Template num. */
    iofst = iofst + 16;
 
    /* Get Grid Definition Template. */
    if (igds[0] == 0)
    {
        if (!(mapgrid = getgridtemplate(igds[4])))
        {       /* undefined template */
            ierr = -5;
            return ierr;
        }
 
        /* Extend the Grid Definition Template, if necessary. The
         *  number of values in a specific template may vary depending
         *  on data specified in the "static" part of the template. */
        if (mapgrid->needext)
        {
            free(mapgrid);
            mapgrid = extgridtemplate(igds[4], igdstmpl);
        }
    }
 
    /* Pack up each input value in array igdstmpl into the
     * appropriate number of octets, which are specified in
     * corresponding entries in array mapgrid. */
    for (i = 0; i < mapgrid->maplen; i++)
    {
        nbits = abs(mapgrid->map[i]) * 8;
        if ((mapgrid->map[i] >= 0) || (igdstmpl[i] >= 0))
            sbit(cgrib, igdstmpl+i, iofst, nbits);
        else
        {
            sbit(cgrib, &one, iofst, 1);
            temp = abs(igdstmpl[i]);
            sbit(cgrib, &temp, iofst+1, nbits-1);
        }
        iofst = iofst + nbits;
    }
 
    /* Pack template extension, if appropriate. */
    j = mapgrid->maplen;
    if (mapgrid->needext && mapgrid->extlen > 0)
    {
        for (i = 0; i < mapgrid->extlen; i++)
        {
            nbits = abs(mapgrid->ext[i]) * 8;
            if (mapgrid->ext[i] >= 0 || igdstmpl[j] >= 0)
                sbit(cgrib, igdstmpl + j, iofst, nbits);
            else
            {
                sbit(cgrib, &one, iofst, 1);
                temp = abs(igdstmpl[j]);
                sbit(cgrib, &temp, iofst + 1, nbits - 1);
            }
            iofst = iofst + nbits;
            j++;
        }
    }
    free(mapgrid);
 
    /* If requested, insert optional list of numbers defining number
     * of points in each row or column. This is used for non regular
     * grids. */
    if (igds[2] != 0)
    {
        nbits = igds[2] * 8;
        sbits(cgrib, ideflist, iofst, nbits, 0, idefnum);
        iofst = iofst + (nbits * idefnum);
    }
 
    /* Calculate length of section 3 and store it in octets 1-4 of section 3. */
    lensec3 = (iofst - ibeg) / 8;
    sbit(cgrib, &lensec3, ibeg, 32);
 
    /* Update current byte total of message in Section 0. */
    lencurr += lensec3;
    sbit(cgrib, &lencurr, 96, 32);
 
    return lencurr;
}