external/io_grib1/MEL_grib1/grib_dec.c

   1 #include <stdio.h>              /* standard I/O header file          */
   2 #include <stdlib.h>
   3 #include <string.h>
   4 #include "dprints.h"            /* for dprints & func prototype*/
   5 #include "gribfuncs.h"          /* prototypes */
   6
   7 /* PROGRAMMER : Steve Lowe and Todd Kienitz, SAIC Monterey
   8    DATE       : February 7, 1996
   9                 Oct. 1996 by Alice Nakajima, SAIC Monterey
  10 *
  11 *********************************************************************
  12 * A.  FUNCTION:  grib_dec
  13 *     decode a Gridded Binary (GRIB edition 1) format message
  14 *
  15 *    INTERFACE:
  16 *      int grib_dec (curr_ptr, pds, gds, bds_head, bms, ppgrib_data, errmsg)
  17 *
  18 *    ARGUMENTS (I=input, O=output, I&O=input and output):
  19 *      (I)  char *curr_ptr;
  20 *           pointer to block containing GRIB message to decode;
  21 *      (O)  PDS_INPUT  *pds ;
  22 *           to be filled with decoded Product Defn Section info;
  23 *      (O)  grid_desc_sec  *gds;
  24 *           to be filled with decoded Binary Data Section info;
  25 *      (O)  BDS_HEAD_INPUT *bds_head;
  26 *           to be filled with decoded Binary Data Section info;
  27 *      (O)  BMS_INPUT *bms;
  28 *           to be filled with decoded Bitmap Section info;
  29 *      (O)  float **ppgrib_data;
  30 *           points to NULL upon entry; upon successful exit, points to newly
  31 *           malloced Float array filled with unpacked and restored data;
  32 *      (O)  char *errmsg;
  33 *           Empty array, Returned filled if error occurred;
  34 *
  35 *     RETURN CODE:
  36 *        0> Success, **ppgrib_data now points to a block containing
  37 *           the unpacked & restored data (float);
  38 *        1> Fail: first 4 bytes of curr_ptr is not 'GRIB'
  39 *        2> Fail: last 4 bytes of curr_ptr is not '7777'
  40 *        3> Fail: not Grib Edition 1
  41 *        4> Fail: unknown projection type;
  42 ***********************************************************************
  43 */
  44
  45 #if PROTOTYPE_NEEDED
  46 int   grib_dec (char *curr_ptr, PDS_INPUT  *pds, grid_desc_sec  *gds,
  47                 BDS_HEAD_INPUT *bds_head, BMS_INPUT *bms, float **ppgrib_data,
  48                 char *errmsg)
  49 #else
  50 int grib_dec (curr_ptr, pds, gds, bds_head, bms, ppgrib_data, errmsg)
  51   char          *curr_ptr; /*input= ptr to 1st byte of GRIB message block*/
  52   PDS_INPUT     *pds;      /* output=ptr to Internal PDS struct*/
  53   grid_desc_sec *gds;      /* output=ptr to Internal GDS struct*/
  54   BDS_HEAD_INPUT*bds_head; /*out=ptr to Internal BDS header struct*/
  55   BMS_INPUT     *bms;      /*output=ptr to Internal bitmap section struct*/
  56   float         **ppgrib_data; /*outp=ptr to nothing upon entry; upon exit, */
  57                            /* points to a newly malloced array of floats;  */
  58   char          *errmsg;   /* output= empty unless Error happens */
  59 #endif
  60 {
  61   char *func="grib_dec";
  62   unsigned long lMessageSize;       /* message and section size */
  63   long edition;                     /* GRIB edition number */
  64   int flag;                         /* tests if a condition has happened */
  65   int gds_flag;                     /* set if Gds present */
  66   int nReturn = 0;
  67   unsigned long skip;
  68   float *outdata;
  69   int xsize;
  70   int j;
  71
  72 /*
  73 *
  74 * A.0     DEBUG printing
  75 */
  76  DPRINT1 ("Entering %s\n", func);
  77  DPRINT6 (
  78   "curr_ptr=%ld, pds=%ld, gds=%ld\nbds_head=%ld, bms=%ld, ppgrib_data=%ld\n",
  79   curr_ptr, pds, gds, bds_head, bms, ppgrib_data);
  80 /*
  81 *
  82 * A.1     IF (incoming pointer is not at 'GRIB')
  83 *            RETURN 1  !errmsg filled
  84 *         ENDIF
  85 */
  86 if(strncmp(curr_ptr,"GRIB",4) != 0) {
  87   sprintf (errmsg,"%s:  no 'GRIB' at beg. of this msg\n", func);
  88   nReturn= (1);   /* GRIB not found */
  89  }
  90
  91 /*
  92 *
  93 * A.2     FUNCTION gbyte   !get total message length from IDS
  94 */
  95 skip=32;
  96 gbyte(curr_ptr,&lMessageSize,&skip,24);
  97 DPRINT0 ("lMessageSize\n");
  98
  99 /*
 100 *
 101 * A.3     IF (Message does not end with '7777')
 102 *            RETURN 2  !errmsg filled
 103 *         ENDIF
 104 */
 105 if(strncmp((curr_ptr + lMessageSize - 4),"7777",4)!=0) {
 106   DPRINT1 ("%s:  no '7777' at end of this msg\n", func);
 107   sprintf (errmsg,"%s:  no '7777' at end of this msg\n", func);
 108   nReturn= 2; goto BYE;
 109   }
 110
 111 /*
 112 *
 113 * A.4     EXTRACT the  GRIB edition out of Section 0
 114 *         IF (not GRIB edition 1)
 115 *            RETURN 3  !errmsg filled
 116 *         ENDIF
 117 */
 118 edition = (long) curr_ptr[7];        /* get edition */
 119 pds->usEd_num = (unsigned short) edition;
 120 if(edition != 1) {
 121   DPRINT1 ("%s:  error, not Grib Edition 1 \n", func);
 122   sprintf (errmsg,"%s:  not Grib Edition 1 \n", func);
 123   nReturn=(3);   goto BYE;
 124   }
 125
 126 /*
 127 *
 128 * A.5     MOVE pointer to the Product Definition section
 129 */
 130 curr_ptr = curr_ptr + 8;
 131
 132 /*
 133 *
 134 * A.6     FUNCTION gribgetpds  !decode the PDS
 135 *         RETURN error code if fails  !errmsg filled
 136 */
 137 if( nReturn= gribgetpds(curr_ptr, pds, errmsg)) {
 138    DPRINT2 ("%s:  error=%d  in grib get pds;\n", func, nReturn);
 139    upd_child_errmsg (func, errmsg);
 140    goto BYE;
 141   }
 142
 143 /*
 144 *
 145 * A.7     MOVE pointer to the end of PDS
 146 */
 147 curr_ptr += pds->uslength;
 148
 149 /*
 150 *
 151 * A.8     IF (GDS is present)
 152 */
 153 gds_flag = pds->usGds_bms_id >> 7 & 1;
 154 if(gds_flag)  /* grid description section present */
 155   {
 156 /*
 157 * A.8.1      FUNCTION gribgetgds   !decode GDS
 158 *            RETURN error code if fails  !errmsg filled
 159 */
 160    if( nReturn=gribgetgds(curr_ptr, gds, errmsg)) {
 161       DPRINT2 ("%s:  error=%d  in grib get GDS;\n", func, nReturn);
 162       upd_child_errmsg (func, errmsg);
 163       goto BYE;
 164        }
 165
 166 /*
 167 * A.8.2      MOVE the cursor to the next section (either BMS/BDS)
 168 */
 169    curr_ptr += gds->head.uslength;
 170 /*
 171 * A.8.3      SET the number of data points depending on Projection
 172 */
 173    switch(gds->head.usData_type)
 174   {
 175      case LATLON_PRJ:           /* Lat/Lon Grid */
 176      case GAUSS_PRJ:            /* Gaussian Latitude/Longitude grid */
 177      case ROT_LATLON_PRJ:       /* Rotated Lat/Lon */
 178      case ROT_GAUSS_PRJ:        /* Rotated Gaussian */
 179      case STR_LATLON_PRJ:       /* Stretched Lat/Lon */
 180      case STR_GAUSS_PRJ :       /* Stretched Gaussian */
 181      case STR_ROT_LATLON_PRJ :  /* Stretched and Rotated Lat/Lon */
 182      case STR_ROT_GAUSS_PRJ :   /* Stretched and Rotated Gaussian */
 183         bds_head->ulGrid_size = gds->llg.usNi * gds->llg.usNj;
 184         break;
 185
 186      case MERC_PRJ:             /* Mercator Grid */
 187         bds_head->ulGrid_size = gds->merc.cols * gds->merc.rows;
 188         break;
 189
 190      case LAMB_PRJ:             /* Lambert Conformal */
 191      case ALBERS_PRJ:           /* Albers equal-area */
 192      case OBLIQ_LAMB_PRJ:       /* Oblique Lambert Conformal */
 193         bds_head->ulGrid_size = gds->lam.iNx * gds->lam.iNy;
 194         break;
 195
 196      case POLAR_PRJ:            /* Polar Stereographic */
 197         bds_head->ulGrid_size = gds->pol.usNx * gds->pol.usNy;
 198         break;
 199
 200      default: /* unknown */
 201        DPRINT2 ("%s: unknown usData_type=%d\n",func,gds->head.usData_type);
 202        sprintf(errmsg,"%s: unknown usData_type=%d\n",
 203         func, gds->head.usData_type);
 204        nReturn= (4); goto BYE;
 205      }
 206   }
 207 /*
 208 * A.8     ENDIF (GDS is present)
 209 */
 210
 211 /*
 212 *
 213 * A.9     IF (bitmap Section is present)
 214 */
 215 flag = pds->usGds_bms_id >> 6 & 1;
 216 if(flag)  /* bit map section present */
 217   {
 218 /*
 219 * A.9.1      FUNCTION gribgetbms   !decode BMS
 220 *            RETURN error code if fails  !errmsg filled
 221 */
 222    if( nReturn=gribgetbms(curr_ptr,bms,gds_flag, bds_head->ulGrid_size,errmsg))
 223    {
 224      DPRINT2 ("%s:  error=%d  in grib get BMS;\n",func,nReturn);
 225      upd_child_errmsg (func, errmsg);
 226      goto BYE;
 227    }
 228
 229 /*
 230 * A.9.2      MOVE the cursor to beginning of Binary Data Section
 231 */
 232      curr_ptr += bms->uslength;
 233
 234   } /* Bms present */
 235 /*
 236 * A.9     ENDIF  !bms present
 237 */
 238
 239
 240 /*
 241 *
 242 * A.10    FUNCTION  gribgetbds()
 243 *         RETURN error code if failed  !errmsg filled
 244 */
 245   if(nReturn=gribgetbds(curr_ptr, pds->sDec_sc_fctr, bms, gds, ppgrib_data,
 246                 bds_head, errmsg))
 247    {
 248      DPRINT2 ("%s:  error=%d  in grib get BDS;\n",func,nReturn);
 249      upd_child_errmsg (func, errmsg);
 250      goto BYE;
 251     }
 252
 253 /*
 254 *
 255 * A.11    SET return code to 0  !no errors
 256 */
 257   nReturn = 0;
 258
 259 /*
 260 *
 261 * A.12    RETURN return code;
 262 */
 263 BYE:
 264   DPRINT2  ("Exit %s, Stat=%d\n", func, nReturn);
 265   return(nReturn);
 266 /*
 267 *
 268 * END OF FUNCTION
 269 *
 270 *
 271 */
 272 }