precision.c

precision.c How to control precision when coding a grib field.

00001 
00010 /*
00011  * C Implementation: precision
00012  *
00013  * Description: how to control decimal precision when packing fields.
00014  *
00015  *
00016  * Author: Enrico Fucile
00017  *
00018  *
00019  */
00020 #include <stdio.h>
00021 #include <stdlib.h>
00022 #include <math.h>
00023 
00024 #include "grib_api.h"
00025 
00026 int main(int argc, char** argv) {
00027   int err = 0;
00028   size_t size=0;
00029 
00030   FILE* in = NULL;
00031   char* infile = "../../data/regular_latlon_surface.grib1";
00032   FILE* out = NULL;
00033   char* outfile = "out.grib1";
00034   grib_handle *h = NULL;
00035   const void* buffer = NULL;
00036   double* values1=NULL;
00037   double* values2=NULL;
00038   double maxa=0,a=0;
00039   double maxv=0,minv=0;
00040   double maxr=0,r=0;
00041   long decimalPrecision;
00042   long bitsPerValue1=0, bitsPerValue2=0;
00043   int i=0;
00044 
00045   in = fopen(infile,"r");
00046   if(!in) {
00047     printf("ERROR: unable to open file %s\n",infile);
00048     return 1;
00049   }
00050 
00051   out = fopen(outfile,"w");
00052   if(!in) {
00053     printf("ERROR: unable to open file %s\n",outfile);
00054     return 1;
00055   }
00056 
00057   /* create a new handle from a message in a file */
00058   h = grib_handle_new_from_file(0,in,&err);
00059   if (h == NULL) {
00060     printf("Error: unable to create handle from file %s\n",infile);
00061   }
00062 
00063   /* bitsPerValue before changing the packing parameters */
00064   GRIB_CHECK(grib_get_long(h,"bitsPerValue",&bitsPerValue1),0);
00065 
00066   /* get the size of the values array*/
00067   GRIB_CHECK(grib_get_size(h,"values",&size),0);
00068 
00069   values1 = malloc(size*sizeof(double));
00070   /* get data values before changing the packing parameters*/
00071   GRIB_CHECK(grib_get_double_array(h,"values",values1,&size),0);
00072 
00073   /* changing decimal precition to 2 means that 2 decimal digits
00074      are preserved when packing.  */
00075   decimalPrecision=2;
00076   GRIB_CHECK(grib_set_long(h,"changeDecimalPrecision",decimalPrecision),0);
00077    
00078   /* bitsPerValue after changing the packing parameters */
00079   GRIB_CHECK(grib_get_long(h,"bitsPerValue",&bitsPerValue2),0);
00080 
00081   values2 = malloc(size*sizeof(double));
00082   /* get data values after changing the packing parameters*/
00083   GRIB_CHECK(grib_get_double_array(h,"values",values2,&size),0);
00084 
00085   /* computing error */
00086   maxa=0;
00087   maxr=0;
00088   maxv=values2[0];
00089   minv=maxv;
00090   for (i=0;i<size;i++) {
00091      a=fabs(values2[i]-values1[i]);
00092      if ( values2[i] > maxv ) maxv=values2[i];
00093      if ( values2[i] < maxv ) minv=values2[i];
00094      if ( values2[i] !=0 ) r=fabs((values2[i]-values1[i])/values2[i]);
00095      if ( a > maxa ) maxa=a;
00096      if ( r > maxr ) maxr=r;
00097   }
00098   printf("max absolute error = %g\n",maxa);
00099   printf("max relative error = %g\n",maxr);
00100   printf("min value = %g\n",minv);
00101   printf("max value = %g\n",maxv);
00102 
00103   printf("old number of bits per value=%ld\n",(long)bitsPerValue1);
00104   printf("new number of bits per value=%ld\n",(long)bitsPerValue2);
00105 
00106   /* get the coded message in a buffer */
00107   GRIB_CHECK(grib_get_message(h,&buffer,&size),0);
00108 
00109   /* write the buffer in a file*/
00110   if(fwrite(buffer,1,size,out) != size) 
00111   {
00112      perror(argv[1]);
00113      exit(1);
00114   }
00115 
00116   /* delete handle */
00117   grib_handle_delete(h);
00118 
00119   fclose(in);
00120   fclose(out);
00121 
00122   return 0;
00123 }
00124 

Generated on Tue Sep 22 15:18:21 2009 for grib_api by  doxygen 1.5.3