Pileus Git - ~andy/rsl/blob - rapic_routines.c

   1 #include <stdio.h>
   2 #include "rapic_routines.h"
   3 #include <string.h>
   4
   5 #include <netinet/in.h>
   6 void rapic_decode(unsigned char *inbuf, int inbytes, unsigned char *outbuf, int *outbytes,
   7                                   float *azim, float *elev, int *delta_time)
   8 {
   9   /* Decode RLE rapic buffer.
  10    *
  11    * Output to 'outbuf' and indicate the size w/ 'nout'.
  12    */
  13
  14   /* There is some risidual parsing to do:
  15    *
  16  *      @               <data>                  *
  17  *      <data>          AAA.A,EEE.E,TTT=LEN16D1D1D1NLD1D1D1NLNL
  18  *
  19  *      AAA.A           Azimuth in degrees
  20  *      EEE.E           Elevation in degress
  21  *      TTT             Delta time in seconds since the start of this scan
  22  *      LEN16           2 byte long length of radial
  23  *
  24  *      D1              Run length coded Data.
  25  *      NL              Null The Next Byte is count of NULL data.
  26  *      NLNL            End of this Radial
  27  *
  28  *         eg.          There will be no white space in the actual radial data.
  29  *
  30  *                      @066.1,010.6,004=002082B2817F8C84830048D72D0038
  31  *                                       999C0036202D35FD2C00238A99008AFE920000
  32  *                      Azimuth = 66.1
  33  *                      Elev    = 10.6
  34  *                      Dt      = 4 sec since the start
  35  *                      0020    = Bytes to follow
  36  *                      Data    = 82,B2,81,7F,8C,84,83
  37  *                      0048    = 48H null bytes
  38  *                      Data    = D7,2D
  39  *                      0038    = 38H null bytes
  40  *                      Data    = 99,9C
  41  *                      0036    = 36H Null bytes
  42  *                      ........................
  43  *                      0000    = End of Data.
  44  *      In  versions before 10.1
  45  *      @               <data>
  46  *      <data>          AAALEN16D1D1D1NLD1D1D1NLNL
  47  *
  48  *      AAA             Azimuth in degrees
  49  *      LEN16           2 byte long length of radial
  50  *
  51  *      D1              Run length coded Data.
  52  *      NL              Null The Next Byte is count of NULL data.
  53  *      NLNL            End of this Radial
  54  *
  55  *         eg.          There will be no white space in the actual radial data.
  56  *
  57  *                      @066002082B2817F8C84830048D72D0038
  58  *                          999C0036202D35FD2C00238A99008AFE920000
  59  *                      Azimuth = 66
  60  *                      0020    = Bytes to follow
  61  *                      Data    = 82,B2,81,7F,8C,84,83
  62  *                      0048    = 48H null bytes
  63  *                      Data    = D7,2D
  64  *                      0038    = 38H null bytes
  65  *                      Data    = 99,9C
  66  *                      0036    = 36H Null bytes
  67  *                      ........................
  68  *                      0000    = End of Data.
  69  *
  70  */
  71
  72
  73   /* The parser won't give us a line w/o '@' at the begining nor '\0\0'
  74    * at the end.  So we can be sure of that.
  75    */
  76
  77   int i;
  78   char prefix[16];
  79   unsigned short nnulls;
  80   short i16;
  81
  82   /* Find the '=' and start RLE decode from there. */
  83   *outbytes = 0;
  84   memset(prefix, '\0', sizeof(prefix));
  85   memcpy(prefix, &inbuf[1], 15);
  86
  87   sscanf(prefix, "%f,%f,%d", azim, elev, delta_time);
  88   /*  fprintf(stderr, "AZIM=%f, ELEV=%f, TTT=%d\n", *azim, *elev, *delta_time); */
  89
  90   /* Now, decode RLE. Don't care about 17,18 (they are the RLE buf size) */
  91   memcpy(&i16, &inbuf[17], 2);
  92   i16 = ntohs(i16);
  93   /*  fprintf(stderr, "Expecting %d bins\n", (int)i16); */
  94   i = 19;
  95   while (i<inbytes-2) {
  96         /*      fprintf(stderr, "i=%d byte=%2.2x(next %2.2x%2.2x) outbytes=%d\n", i, (unsigned char)inbuf[i], (unsigned char)inbuf[i+1], (unsigned char)inbuf[i+2], *outbytes); */
  97         if (inbuf[i] == '\0') { /* Next byte is a count of NULL's */
  98           i++;
  99           nnulls = (int)inbuf[i];
 100           /* fprintf(stderr, "NULL .. number of nulls=%4.4x\n", nnulls); */
 101           memset(&outbuf[*outbytes], '\0', (int)nnulls);
 102           *outbytes += nnulls;
 103         } else {
 104           /* fprintf(stderr, "Data\n"); */
 105           /* Data. */
 106           outbuf[*outbytes] = inbuf[i];
 107           *outbytes += 1;
 108         }
 109         i++;
 110   }
 111
 112   /* fprintf(stderr, "Decoded RLE: len=%d buf=<", *outbytes); */
 113   /* binprint(outbuf, *outbytes); */
 114   /* fprintf(stderr, ">\n"); */
 115 }
 116
 117
 118
 119 /*---------------------------------------------------------*/
 120 /*                                                         */
 121 /*                    binprint                             */
 122 /*                                                         */
 123 /*---------------------------------------------------------*/
 124 void binprint(char *s, int n)
 125 {
 126 int i;
 127
 128 for (i=0; i<n; i++)
 129          fprintf(stderr,"%c", s[i]);
 130
 131 }
 132
 133 /* Shut the damn linker up for librsl.so on Linux. This reference is
 134  * in the HDF library; don't need it, don't care.
 135  */
 136 void i_len(void)
 137 {
 138 }
 139
 140
 141 #include <time.h>
 142 void rapic_fix_time (Ray *ray)
 143 {
 144   struct tm the_time;
 145   float fsec;
 146   /* Fixes possible overflow values in month, day, year, hh, mm, ss */
 147   /* Normally, ss should be the overflow.  This code ensures end of
 148    * month, year and century are handled correctly by using the Unix
 149    * time functions.
 150    */
 151   if (ray == NULL) return;
 152   memset(&the_time, 0, sizeof(struct tm));
 153   the_time.tm_sec = ray->h.sec;
 154   fsec = ray->h.sec - the_time.tm_sec;
 155   the_time.tm_min  = ray->h.minute;
 156   the_time.tm_hour = ray->h.hour;
 157   the_time.tm_mon  = ray->h.month - 1;
 158   the_time.tm_year = ray->h.year - 1900;
 159   the_time.tm_mday = ray->h.day;
 160   the_time.tm_isdst = -1;
 161   (void) mktime(&the_time);
 162   /* The time is fixed. */
 163   ray->h.sec    = the_time.tm_sec;
 164   ray->h.sec   += fsec;
 165   ray->h.minute = the_time.tm_min;
 166   ray->h.hour   = the_time.tm_hour;
 167   ray->h.month  = the_time.tm_mon + 1;
 168   ray->h.year   = the_time.tm_year + 1900;
 169   ray->h.day    = the_time.tm_mday;
 170   return;
 171 }
 172
 173 void rapic_load_ray_header(Rapic_sweep_header rh, int iray, int isweep, float elev, float azim, Ray_header *h)
 174 {
 175   sscanf(rh.yyyymoddhhmmss,"%4d%2d%2d%2d%2d%2f",
 176                  &h->year,&h->month,&h->day,
 177                  &h->hour, &h->minute, &h->sec);
 178
 179   h->unam_rng = rh.end_range/1000.0;  /* Unambiguous range. (KM). */
 180   h->azimuth  = azim;   /* Azimuth angle. (degrees). Must be positive
 181                                         * 0=North, 90=east, -90/270=west.
 182                     * This angle is the mean azimuth for the whole ray.
 183                                         * Eg. for NSIG the beginning and end azimuths are
 184                                         *     averaged.
 185                                         */
 186   h->ray_num  = iray;     /* Ray no. within elevation scan. */
 187   h->elev     = rh.elev;  /* Elevation angle. (degrees). */
 188   h->elev_num = isweep;   /* Elevation no. within volume scan. */
 189
 190   h->range_bin1 = rh.start_range;       /* Range to first gate.(meters) */
 191   h->gate_size  = rh.range_resolution;  /* Data gate size (meters)*/
 192
 193   h->vel_res    = rh.range_resolution;   /* Doppler velocity resolution */
 194   h->sweep_rate = rh.anglerate/6.0;      /* Sweep rate. Full sweeps/min. */
 195
 196   h->prf       = rh.prf;          /* Pulse repetition frequency, in Hz. */
 197   h->azim_rate = rh.anglerate;    /* degrees/sec */
 198   h->fix_angle = elev;
 199   h->pitch     = 0;       /* Pitch angle. */
 200   h->roll      = 0;       /* Roll  angle. */
 201   h->heading   = 0;       /* Heading. */
 202   h->pitch_rate   = 0;            /* (angle/sec) */
 203   h->roll_rate    = 0;            /* (angle/sec) */
 204   h->heading_rate = 0;            /* (angle/sec) */
 205   h->lat          = rh.lat;       /* Latitude (degrees) */
 206   h->lon          = rh.lon;       /* Longitude (degrees) */
 207   h->alt          = rh.height;    /* Altitude (m) */
 208   h->rvc          = 0;            /* Radial velocity correction (m/sec) */
 209   h->vel_east     = 0;            /* Platform velocity to the east  (m/sec) */
 210   h->vel_north    = 0;            /* Platform velocity to the north (m/sec) */
 211   h->vel_up       = 0;            /* Platform velocity toward up    (m/sec) */
 212   h->pulse_count  = 0;            /* Pulses used in a single dwell time. */
 213   h->pulse_width  = rh.pulselen;  /* Pulse width (micro-sec). */
 214   h->beam_width   = rh.angle_resolution; /* Beamwidth in degrees. */
 215   h->frequency    = rh.freq/1000.0; /* Carrier freq. GHz. */
 216   h->wavelength   = 0;            /* Wavelength. Meters. */
 217   h->nyq_vel      = rh.nyquist;   /* Nyquist velocity. m/s */
 218
 219   return;
 220 }
 221
 222 extern float rapic_nyquist;
 223
 224 float RAPIC_DZ_F(unsigned char x) {
 225   if (x == 0) return NOECHO;
 226   return (((float)x-64)/2.0);  /* rapic -> float */
 227 }
 228
 229 float RAPIC_VR_F(unsigned char x) {
 230   if (x == 0) return BADVAL;
 231   return (((float)((int)x-128))/128.0*rapic_nyquist);  /* rapic -> float */
 232 }
 233
 234 float RAPIC_SW_F(unsigned char x) {
 235   if (x == 0) return NOECHO;
 236   return (((float)x)/256.0*rapic_nyquist);  /* rapic -> float */
 237 }
 238
 239 float RAPIC_ZT_F(unsigned char x) {
 240   return RAPIC_DZ_F(x);
 241 }
 242
 243 float RAPIC_ZD_F(unsigned char x) {
 244   if (x == 0) return NOECHO;
 245   return (((float)x-128)/16.0);  /* rapic -> float */
 246 }
 247
 248 /* USE RSL INDEXING! */
 249 static  float (*RAPIC_f_list[])(unsigned char x) = {RAPIC_DZ_F,
 250                                                                                                         RAPIC_VR_F,
 251                                                                                                         RAPIC_SW_F,
 252                                                                                                         NULL,
 253                                                                                                         RAPIC_ZT_F,
 254                                                                                                         NULL,
 255                                                                                                         NULL,
 256                                                                                                         RAPIC_ZD_F};
 257
 258
 259 void rapic_load_ray_data(unsigned char *buf, int bufsize, int ifield, Ray *ray)
 260 {
 261   /* ifield is the RSL numbering scheme for field types.  The conversion
 262    * is done in rapic.y.  In other words, we've already converted rapic
 263    * index to RSL indexes.
 264    */
 265   int i;
 266   for (i=0; i<bufsize; i++) {
 267         ray->range[i] = ray->h.invf(RAPIC_f_list[ifield](buf[i]));
 268         /*      fprintf(stderr,"i=%d ifield=%d, buf[%d]=%2.2x, ray->range[%d]=%4.4x value=%f\n", i,ifield,i,buf[i],i,ray->range[i], RAPIC_f_list[ifield](buf[i]) ); */
 269   }
 270   ray->h.nbins = bufsize;
 271 }
 272
 273 Radar *fill_header(Radar *radar)
 274 {
 275   /* Learn as much as possible from the Ray headers.  Place this
 276    * information into radar->h.xxxxxx
 277    */
 278   Ray *ray;
 279   Volume *volume;
 280   int i;
 281   float tmp;
 282
 283   volume = NULL;
 284   if (radar == NULL) return NULL;
 285   for (i=0; i<radar->h.nvolumes && !(volume = radar->v[i]); i++)
 286         ;
 287   if (volume == NULL) return NULL;
 288
 289   ray = RSL_get_first_ray_of_volume(volume);
 290   if (ray  == NULL) return NULL;
 291
 292   radar->h.month = ray->h.month;
 293   radar->h.day   = ray->h.day;
 294   radar->h.year  = ray->h.year;
 295   radar->h.hour  = ray->h.hour;
 296   radar->h.minute= ray->h.minute;
 297   radar->h.sec   = ray->h.sec; /* Second plus fractional part. */
 298   sprintf(radar->h.radar_type, "rapic"); /* Type of radar. */
 299                   /* nvolumes   is already filled in YACC. */
 300                   /* number     is already filled in YACC. */
 301                   /* name       is already filled in YACC. */
 302                   /* radar_name is already filled in YACC. */
 303   /*  radar->h.city[15];     */ /* Not available from RAPIC. */
 304   /*  radar->h.state[3];     */ /* Not available from RAPIC. */
 305   /*  radar->h.country[15];  */ /* Not available from RAPIC. */
 306
 307    /** Latitude deg, min, sec **/
 308   radar->h.latd = (int)ray->h.lat;
 309   tmp = (ray->h.lat - radar->h.latd) * 60.0;
 310   radar->h.latm = (int)tmp;
 311   radar->h.lats = (int)((tmp - radar->h.latm) * 60.0);
 312    /** Longitude deg, min, sec **/
 313   radar->h.lond = (int)ray->h.lon;
 314   tmp = (ray->h.lon - radar->h.lond) * 60.0;
 315   radar->h.lonm = (int)tmp;
 316   radar->h.lons = (int)((tmp - radar->h.lonm) * 60.0);
 317
 318   radar->h.height = ray->h.alt; /* height of site in meters above sea level*/
 319   radar->h.spulse = 0; /* length of short pulse (ns)*/
 320   radar->h.lpulse = 0; /* length of long pulse (ns) */
 321
 322   return radar;
 323 }