7 This is the TRMM Office Radar Software Library.
8 Copyright (C) 1996, 1997
10 Space Applications Corporation
13 This library is free software; you can redistribute it and/or
14 modify it under the terms of the GNU Library General Public
15 License as published by the Free Software Foundation; either
16 version 2 of the License, or (at your option) any later version.
18 This library is distributed in the hope that it will be useful,
19 but WITHOUT ANY WARRANTY; without even the implied warranty of
20 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
21 Library General Public License for more details.
23 You should have received a copy of the GNU Library General Public
24 License along with this library; if not, write to the Free
25 Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
27 /*----------------------------------------------------------------------
30 * Ingest a Lassen file and fill the Lassen_volume data structure.
31 * It can handle version 1.3 and 1.4 files (MCTEX data).
33 * The algorthm is pretty simple, call an 'xdr' routine each time for
34 * each and every member of Lassen_volume and all substructures.
36 * The order and type of xdr routine called defines the file organization.
37 *----------------------------------------------------------------------
43 #include <netinet/in.h>
46 /* xdr_destroy is broken on HPUX, SGI, and SUN; Linux, the only working one? */
47 /* This broken behavior could be from old xdr versions, too. */
48 #undef BROKEN_XDR_DESTROY
49 #if defined(__hpux) || defined(__sgi) || defined(__sun)
50 #define BROKEN_XDR_DESTROY
52 /*************************************************************/
54 /* read_lassen_head */
56 /*************************************************************/
57 int read_lassen_head(FILE *f, Lassen_head *head)
63 memset(head->magic, 0, sizeof(head->magic));
65 xdrstdio_create(&xdrs, f, XDR_DECODE);
68 * Can I read the first string? If not, then this is probably
69 * not a valid header. If this happens, then return.
72 if((rc &= xdr_string(&xdrs, &tmp, 8))==0) {
73 #ifndef BROKEN_XDR_DESTROY
78 rc &= xdr_u_char(&xdrs, &head->mdate.year);
79 rc &= xdr_u_char(&xdrs, &head->mdate.month);
80 rc &= xdr_u_char(&xdrs, &head->mdate.day);
81 rc &= xdr_u_char(&xdrs, &head->mdate.hour);
82 rc &= xdr_u_char(&xdrs, &head->mdate.minute);
83 rc &= xdr_u_char(&xdrs, &head->mdate.second);
84 rc &= xdr_u_char(&xdrs, &head->cdate.year);
85 rc &= xdr_u_char(&xdrs, &head->cdate.month);
86 rc &= xdr_u_char(&xdrs, &head->cdate.day);
87 rc &= xdr_u_char(&xdrs, &head->cdate.hour);
88 rc &= xdr_u_char(&xdrs, &head->cdate.minute);
89 rc &= xdr_u_char(&xdrs, &head->cdate.second);
90 rc &= xdr_int(&xdrs, &head->type);
92 rc &= xdr_string(&xdrs, &tmp, 16);
94 rc &= xdr_string(&xdrs, &tmp, 16);
95 rc &= xdr_int(&xdrs, &head->protection);
96 rc &= xdr_int(&xdrs, &head->checksum);
97 tmp=head->description;
98 rc &= xdr_string(&xdrs, &tmp, 40);
99 rc &= xdr_int(&xdrs, &head->id);
101 rc &= xdr_int(&xdrs, &head->spare[i]);
103 #ifndef BROKEN_XDR_DESTROY
110 /*************************************************************/
112 /* free_lassen_volume */
114 /*************************************************************/
116 void free_lassen_volume(Lassen_volume *vol)
120 if (vol == NULL) return;
121 for(i=0; i<vol->numsweeps; i++ ) {
122 if(vol->index[i] == NULL) continue;
123 for(j=0; j<vol->index[i]->numrays; j++ )
124 if (vol->index[i]->ray[j] != NULL)
125 free( vol->index[i]->ray[j] );
126 free( vol->index[i] );
132 static void short_shift(unsigned short *n)
137 /*************************************************************/
139 /* read_lassen_volume */
141 /*************************************************************/
142 int read_lassen_volume( XDR *xdr, Lassen_volume *vol)
147 unsigned short ushort_tmp;
149 /* Read the version number. */
150 rc &= xdr_u_short(xdr, &vol->version);
153 * Currently, Version 1.3 and 1.4 seem to be the same on disk.
154 * I have no doc's. The only thing that I know of to be different
155 * is how the ray azimuth angles are calculated. For MCTEX data,
156 * this is time dependant. The correction will happen in the
157 * application; after RSL_lassen_to_radar is called, for example.
159 * Version 1.3 is coded as 13 and version 1.4 is coded as 14.
161 if( (int)vol->version != 13 && (int)vol->version != 14 ) {
162 (void)fprintf(stderr, "Incompatible version of LASSEN file.\n" );
163 (void)fprintf(stderr, "File version: %d.%d\n",
164 (int)vol->version/10, (int)vol->version%10);
168 /* Read in the Lassen_volume. */
169 rc &= xdr_short(xdr, &vol->filled);
170 rc &= xdr_u_int(xdr, &vol->volume);
171 rc &= xdr_u_short(xdr, &vol->sweep);
172 rc &= xdr_u_short(xdr, &vol->sweep_type);
173 rc &= xdr_u_short(xdr, &vol->max_height);
174 rc &= xdr_u_short(xdr, &vol->status);
175 rc &= xdr_u_short(xdr, &vol->min_fangle);
176 rc &= xdr_u_short(xdr, &vol->max_fangle);
177 rc &= xdr_u_short(xdr, &vol->min_var);
178 rc &= xdr_u_short(xdr, &vol->max_var);
179 rc &= xdr_u_short(xdr, &vol->a_start);
180 rc &= xdr_u_short(xdr, &vol->a_stop);
181 rc &= xdr_u_short(xdr, &vol->numsweeps);
182 for(i=0;i<LASSEN_MAX_SWEEPS;i++)
183 rc &= xdr_u_short(xdr, &vol->fangles[i]);
184 rc &= xdr_u_short(xdr, &vol->gatewid);
185 rc &= xdr_u_short(xdr, &vol->rangeg1);
186 for(i=0;i<LASSEN_MAX_SWEEPS;i++)
187 rc &= xdr_u_short(xdr, &vol->numgates[i]);
188 rc &= xdr_u_short(xdr, &vol->maxgates);
189 rc &= xdr_u_short(xdr, &ushort_tmp);
190 rc &= xdr_u_short(xdr, &ushort_tmp);
191 rc &= xdr_u_short(xdr, &ushort_tmp);
192 rc &= xdr_u_short(xdr, &ushort_tmp);
193 rc &= xdr_u_short(xdr, &vol->prf);
194 rc &= xdr_u_short(xdr, &vol->prflow);
195 rc &= xdr_u_int(xdr, &vol->freq);
196 rc &= xdr_u_short(xdr, &vol->n_pulses);
198 for(i=0;i<LASSEN_MAX_SWEEPS;i++)
199 for(j=0;j<NUMOFFSETS;j++)
200 rc &= xdr_u_short(xdr, &vol->offset[i][j]);
202 rc &= xdr_u_char(xdr, &vol->year);
203 rc &= xdr_u_char(xdr, &vol->month);
204 rc &= xdr_u_char(xdr, &vol->day);
205 rc &= xdr_u_char(xdr, &vol->shour);
206 rc &= xdr_u_char(xdr, &vol->sminute);
207 rc &= xdr_u_char(xdr, &vol->ssecond);
208 rc &= xdr_u_char(xdr, &vol->ehour);
209 rc &= xdr_u_char(xdr, &vol->eminute);
210 rc &= xdr_u_char(xdr, &vol->esecond);
211 /* The flags use bits. The data type length is short; I expect. */
212 rc &= xdr_u_short(xdr, (unsigned short *)&vol->volflags);
213 short_shift((unsigned short *)&vol->volflags);
214 /* This is radar info. */
215 tmp=vol->radinfo.radar_name;
216 rc &= xdr_string(xdr, &tmp, 8);
217 tmp=vol->radinfo.site_name;
218 rc &= xdr_string(xdr, &tmp, 8);
219 rc &= xdr_u_short(xdr, &vol->radinfo.antenna_height);
220 rc &= xdr_short(xdr, &vol->radinfo.latitude.degree);
221 rc &= xdr_short(xdr, &vol->radinfo.latitude.minute);
222 rc &= xdr_short(xdr, &vol->radinfo.latitude.second);
223 rc &= xdr_short(xdr, &vol->radinfo.longitude.degree);
224 rc &= xdr_short(xdr, &vol->radinfo.longitude.minute);
225 rc &= xdr_short(xdr, &vol->radinfo.longitude.second);
228 * The return value is the logical and of all the xdr reads.
233 /*************************************************************/
235 /* read_lassen_ray */
237 /*************************************************************/
238 int read_lassen_ray(XDR *xdr, Lassen_ray *ray )
242 /* Read in Lassen_ray */
243 rc &= xdr_u_short( xdr, &ray->vangle);
244 rc &= xdr_u_short( xdr, &ray->fanglet);
245 rc &= xdr_u_short( xdr, &ray->fanglea);
246 rc &= xdr_u_short( xdr, &ray->a_start);
247 rc &= xdr_u_short( xdr, &ray->a_stop);
248 rc &= xdr_u_char( xdr, &ray->max_height);
249 rc &= xdr_u_char( xdr, &ray->volume);
250 rc &= xdr_u_char( xdr, &ray->sweep);
251 rc &= xdr_u_char( xdr, &ray->sweep_type);
252 rc &= xdr_u_short( xdr, &ray->gatewid);
253 rc &= xdr_u_short( xdr, &ray->rangeg1);
254 rc &= xdr_u_short( xdr, &ray->numgates);
255 rc &= xdr_u_short( xdr, &ray->prf);
256 rc &= xdr_u_short( xdr, &ray->prflow);
257 rc &= xdr_u_short( xdr, &ray->n_pulses);
258 rc &= xdr_u_char( xdr, &ray->p_width);
259 rc &= xdr_u_char( xdr, &ray->cfilter);
260 rc &= xdr_u_short( xdr, &ray->status);
261 memset(&ray->flags, 0, 4);
262 /* I expect the length of the bit field to be a 'short', never 'long' */
263 rc &= xdr_u_short( xdr, (unsigned short *)&ray->flags); /* a bit field */
264 short_shift((unsigned short *)&ray->flags);
266 fprintf(stderr, "short_shift: flags = %x\n", ray->flags);
267 fprintf(stderr, "packed = %d", (int)ray->flags.packed);
268 fprintf(stderr, "good_data = %d", (int)ray->flags.good_data);
269 fprintf(stderr, "uz = %d", (int)ray->flags.uz);
270 fprintf(stderr, "cz = %d", (int)ray->flags.cz);
271 fprintf(stderr, "vel= %d", (int)ray->flags.vel);
272 fprintf(stderr, "wid= %d", (int)ray->flags.wid);
273 fprintf(stderr, "zdr= %d", (int)ray->flags.zdr);
274 fprintf(stderr, "phi= %d", (int)ray->flags.phi);
275 fprintf(stderr, "rho= %d", (int)ray->flags.rho);
276 fprintf(stderr, "ldr= %d", (int)ray->flags.ldr);
277 fprintf(stderr, "kdp= %d", (int)ray->flags.kdp);
278 fprintf(stderr, "tim= %d", (int)ray->flags.time);
279 fprintf(stderr, "spa= %x", (int)ray->flags.spares);
281 for(i=0;i<NUMOFFSETS;i++)
282 rc &= xdr_u_short( xdr, &ray->offset[i]);
284 rc &= xdr_u_char( xdr, &ray->year);
285 rc &= xdr_u_char( xdr, &ray->month);
286 rc &= xdr_u_char( xdr, &ray->day);
287 rc &= xdr_u_char( xdr, &ray->hour);
288 rc &= xdr_u_char( xdr, &ray->minute);
289 rc &= xdr_u_char( xdr, &ray->second);
295 /*************************************************************/
297 /* read_entire_lassen_file */
299 /*************************************************************/
300 int read_entire_lassen_file(FILE *f, Lassen_volume *vol)
308 int lastray=0,numrays=0, rc;
314 /* Skip the header. */
315 read_lassen_head(f, &head);
317 /* Connect to the XDR stream. */
318 xdrstdio_create(&xdr, f, XDR_DECODE);
320 /* Check the volume header. Is the version is correct? */
321 if( read_lassen_volume( &xdr, vol ) == 0) {
322 #ifndef BROKEN_XDR_DESTROY
330 /* Read until XDR error. */
334 * Read in the ray header. If there is an error in reading
335 * in the ray header, then make sure the rc value is
338 if( read_lassen_ray( &xdr, &ray ) == 0 )
341 if( (short)ANGLE_CONVERT(ray.vangle) < 0 ||
342 (short)ANGLE_CONVERT(ray.vangle) > 360 ) {
343 fprintf(stderr, "read_entire_lassen_file: angle out of range (%d)\n",
344 ANGLE_CONVERT(ray.vangle));
349 if(ray.volume != (vol->volume&0x00ff)) {
350 fprintf(stderr, "read_entire_lassen_file: Volume serial number out of sync:\n");
351 fprintf(stderr, " Ray->volume = %u, Should be %u\n",
352 ray.volume, (unsigned char)vol->volume&0x00ff );
353 fprintf(stderr, " (Ray->volume == Vol->volume(%d) & 0x00ff)\n",
355 if(isweep == -99) rc = 0;
360 * If we have a new sweep, then create
361 * the new sweep data structure. The 'isweep'
362 * contains the offset number in the volume. 'ray.sweep'
363 * has the sweep number.
365 if((ray.sweep-1) != isweep) {
367 fprintf(stderr,"Sweep number in a ray is 0!\n");
368 fprintf(stderr,"Aborting read_entire_lassen_file()\n");
369 #ifndef BROKEN_XDR_DESTROY
377 isweep = ray.sweep-1;
380 /* Reset the last ray when passing 2 pi. */
382 if( lastray < 0 ) lastray = 359;
384 sweep = vol->index[isweep-1];
387 fprintf(stderr, "vol->index[%d] is NULL.\n", isweep-1);
390 sweep->ehour = sweep->ray[lastray]->hour;
391 sweep->eminute = sweep->ray[lastray]->minute;
392 sweep->esecond = sweep->ray[lastray]->second;
393 sweep->numrays = numrays;
394 sweep->max_var = sweep->ray[lastray]->vangle;
398 /* We have a new sweep. */
399 if((vol->index[isweep] = (Lassen_sweep *)calloc
400 (1,sizeof(Lassen_sweep))) == 0) {
401 perror("read_entire_lassen_file");
408 /* Copy to the Sweep. */
409 sweep = vol->index[isweep];
410 sweep->sweep = isweep;
411 sweep->shour = ray.hour;
412 sweep->sminute = ray.minute;
413 sweep->ssecond = ray.second;
414 sweep->volume = ray.volume;
415 sweep->sweep = ray.sweep;
416 sweep->sweep_type = ray.sweep_type;
417 sweep->max_height = ray.max_height;
418 sweep->fangle = ray.fanglet;
419 sweep->a_start = ray.a_start;
420 sweep->a_stop = ray.a_stop;
421 sweep->gatewid = ray.gatewid;
422 sweep->rangeg1 = ray.rangeg1;
423 sweep->numgates = ray.numgates;
424 sweep->prf = ray.prf;
425 sweep->prflow = ray.prflow;
426 sweep->n_pulses = ray.n_pulses;
427 sweep->p_width = ray.p_width;
428 sweep->cfilter = ray.cfilter;
429 for(i=0; i<NUMOFFSETS; i++)
430 sweep->offset[i] = ray.offset[i];
431 sweep->year = ray.year;
432 sweep->month = ray.month;
433 sweep->day = ray.day;
434 sweep->min_var = ray.vangle;
440 * If everything is ok, then add the ray to the
444 vol->index[isweep]->status |= ray.status;
446 /* Determine the 'size' of the ray and allocate memory. */
448 /* fprintf(stderr,"\nOFFSETS (%d):", NUMOFFSETS); */
449 for(i=0; i<NUMOFFSETS; i++) {
450 size += (ray.numgates*(ray.offset[i] !=0 ));
451 /* fprintf(stderr, " %d", ray.numgates*(ray.offset[i])); */
454 if((vol->index[isweep]->ray[lastray] =
455 (Lassen_ray *) calloc(1,sizeof(Lassen_ray)+size))==0) {
456 perror("read_entire_lassen_file");
458 #ifndef BROKEN_XDR_DESTROY
465 p = (unsigned char *)vol->index[isweep]->ray[lastray];
467 /* Copy the ray header. */
468 memcpy(p, &ray, sizeof(Lassen_ray));
469 p = (unsigned char *)p + sizeof(Lassen_ray);
472 /* printf("Previous rc=%d ... ", rc); */
473 rc &= xdr_bytes(&xdr, (char **)&p, &size,size);
475 /* printf("after rc=%d. xdr_bytes %d bytes, tbytes %d.\n", rc, size, tbytes);
478 if(lastray > 359) lastray = 0;
482 /* Clean up and get out. */
483 free_lassen_volume(vol);
484 #ifndef BROKEN_XDR_DESTROY
487 fprintf(stderr,"\tAborting read_entire_lassen_file.\n\n");
492 * Check to see if we read a sweep!
498 if( lastray < 0 ) lastray = 359;
500 /* Update the last sweep. */
501 sweep = vol->index[isweep];
502 sweep->ehour = sweep->ray[lastray]->hour;
503 sweep->eminute = sweep->ray[lastray]->minute;
504 sweep->esecond = sweep->ray[lastray]->second;
505 sweep->max_var = ray.vangle;
506 sweep->numrays = numrays;
510 /* Check the number of tilts read vs. the number in the file. */
511 if(vol->numsweeps != isweep+1) {
512 fprintf(stderr, "read_entire_lassen_file: Warning, read error. The number of sweeps ");
513 fprintf(stderr, "detected is different from the number written.");
516 vol->numsweeps = isweep+1;
518 #ifndef BROKEN_XDR_DESTROY