* v1.41 In progress
*
* 1. wsr88d_m31.c: Simplified the WSR-88D ray structure and supporting code.
- * Removed function read_data_moment.
+ * Removed function read_data_moment. Added support for dual-polarization
+ * data fields.
* Added function get_wsr88d_unamb_and_nyq_vel.
- * 2. Added support for Sigmet 2-byte data types, as well as HydroClass 1 and 2 byte
- * types.
- * Files involved: nsig_to_radar.c, nsig.c, nsig.h, volume.c.
- * [TODO: dual PRF mods]
- * 3. Completed DORADE ingest. Files involved: dorade_to_radar.c, dorade.c, volume.c.
- * 4. prune.c: Added global variable prune_radar and two functions to set or unset:
- * RSL_prune_radar_on()
- * RSL_prune_radar_off()
+ * 2. Added support for Sigmet 2-byte data types, as well as HydroClass 1 and 2
+ * byte types. Files involved: nsig_to_radar.c, nsig.c, nsig.h, volume.c.
+ * Modified nsig_to_radar.c and rsl.h to handle Sigmet dual PRF. Much
+ * thanks to Fabio Sato and Cesar Beneti for the dual PRF code.
+ * 3. Completed DORADE ingest. Files involved: dorade_to_radar.c, dorade.c,
+ * volume.c.
*---------------------------------------------------------------------
* v1.40 Released 10/10/2008
*
lib_LTLIBRARIES = librsl.la
-librsl_la_LDFLAGS = -version-info 1:40
+librsl_la_LDFLAGS = -version-info 1:41
librsl_la_SOURCES = \
$(rapic_c) $(radtec_c)\
dorade.c dorade_print.c dorade_to_radar.c\
dnl Process this file with autoconf to produce a configure script.
-AC_INIT(rsl, v1.40.1)
+AC_INIT(rsl, v1.40.3)
AC_CONFIG_SRCDIR(volume.c)
AM_INIT_AUTOMAKE
<hr>
<h3>See also</h3>
-<a href=../examples/Example_RSL_clear.c>Example</a><br>
<a href=RSL_free.html>RSL_free_volume</a>, <a href=RSL_free.html>RSL_free_sweep</a>, <a href=RSL_free.html>RSL_free_ray</a>
<hr>
-<p>Author: <a href=john.merritt.html>John H. Merritt</a>
+<p>Author: John H. Merritt
</body>
<hr>
<h3>See also</h3>
-<a href=../examples/Example_RSL_copy.c>Example</a><br>
<a href=RSL_free.html>RSL_free_volume</a>, <a href=RSL_free.html>RSL_free_sweep</a>, <a href=RSL_free.html>RSL_free_ray</a>
<hr>
-<p>Author: <a href=john.merritt.html>John H. Merritt</a>
+<p>Author: John H. Merritt
</body>
<h3>
<hr>Description</h3>
-<b>RSL_get_window_from_radar </b>calls<b> RSL_get_window_from_volume </b>for the number of volumes. <b>RSL_get_window_from_volume</b> calls <b>RSL_get_window_from_sweep</b> for the number of sweeps. <b>RSL_get_window_from_sweep</b> calls <b>RSL_get_window_from_ray</b> for the number of rays. And <b>RSL_get_window_from_ray</b> gets window, between low_azim and hi_azim, from min_range to max_range.
+<b>RSL_get_window_from_radar </b>calls<b> RSL_get_window_from_volume </b>for the number of volumes.<br>
+<b>RSL_get_window_from_volume</b> calls <b>RSL_get_window_from_sweep</b> for the number of sweeps.<br>
+<b>RSL_get_window_from_sweep</b> calls <b>RSL_get_window_from_ray</b> for the number of rays.<br>
+<b>RSL_get_window_from_ray</b> gets window, between low_azim and hi_azim (degrees), from min_range to max_range (km).
<hr>
<h3>Return value</h3>
<hr>
<h3>See also</h3>
-<a href=../examples/test_get_win.c>Example</a>
+rsl/examples/test_get_win.c
<hr>
-<p>Author: <a href=john.merritt.html>John H. Merritt</a>
+<p>Author: John H. Merritt
</body>
<hr>
<h3>See also</h3>
-<a href=../examples/Example_RSL_new.c>Example</a><br>
<a href=RSL_free.html>RSL_free_volume</a>, <a href=RSL_free.html>RSL_free_sweep</a>, <a href=RSL_free.html>RSL_free_ray</a> <br>
<a href=RSL_copy.html>RSL_copy_volume</a>, <a href=RSL_copy.html>RSL_copy_sweep</a>, <a href=RSL_copy.html>RSL_copy_ray</a>
<hr>
-<p>Author: <a href=john.merritt.html>John H. Merritt</a>
+<p>Author: John H. Merritt
</body>
<hr>
<h3>See also</h3>
-<a href=../examples/Example_RSL_read.c>Example</a><br>
<a href=RSL_new.html>RSL_new_volume</a>, <a href=RSL_new.html>RSL_new_sweep</a>, <a href=RSL_new.html>RSL_new_ray</a>,<br>
<a href=RSL_write.html>RSL_write_volume</a>, <a href=RSL_write.html>RSL_write_sweep</a>, <a href=RSL_write.html>RSL_write_ray</a>,<br>
<a href=RSL_read_radar.html>RSL_read_radar</a>, <a href=RSL_write_radar.html>RSL_write_radar</a>,<br>
<a href=RSL_radar_file_format.html>File format</a>.
<hr>
-<p>Author: <a href=john.merritt.html>John H. Merritt</a>
+<p>Author: John H. Merritt
</body>
#include <stdlib.h>
#include <unistd.h>
#include <strings.h>
+#include <string.h>
#define USE_RSL_VARS
#include "rsl.h"
#include "dorade.h"
int find_rsl_field_index(char *dorade_field_name)
{
/*
- * Dorade: VE, DM, SW, DZ, ZDR, PHI, RHO, LDR, DX, CH, AH, CV, AV
- * RSL: VR, DM, SW, DZ, ZD, PH, RH, LR, *DX, *CH, *AH, *CV, *AV.
+ * Dorade: VE, DM, SW, DBZ, ZDR, PHI, RHO, LDR, DX, CH, AH, CV, AV
+ * RSL: VR, DM, SW, DZ, ZD, PH, RH, LR, *DX, *CH, *AH, *CV, *AV.
*/
- if (strncasecmp(dorade_field_name, "ve", 2) == 0) return VR_INDEX;
- if (strncasecmp(dorade_field_name, "dm", 2) == 0) return DM_INDEX;
- if (strncasecmp(dorade_field_name, "sw", 2) == 0) return SW_INDEX;
- if (strncasecmp(dorade_field_name, "dz", 2) == 0) return DZ_INDEX;
+ if (strncasecmp(dorade_field_name, "ve", 2) == 0) return VR_INDEX;
+ if (strncasecmp(dorade_field_name, "vr", 2) == 0) return VR_INDEX;
+ if (strncasecmp(dorade_field_name, "dm", 2) == 0) return DM_INDEX;
+ if (strncasecmp(dorade_field_name, "sw", 2) == 0) return SW_INDEX;
+ if (strncasecmp(dorade_field_name, "dbz", 3) == 0) return DZ_INDEX;
if (strncasecmp(dorade_field_name, "zdr", 3) == 0) return ZD_INDEX;
if (strncasecmp(dorade_field_name, "phi", 3) == 0) return PH_INDEX;
if (strncasecmp(dorade_field_name, "rho", 3) == 0) return RH_INDEX;
if (strncasecmp(dorade_field_name, "ldr", 3) == 0) return LR_INDEX;
- if (strncasecmp(dorade_field_name, "dx", 2) == 0) return DX_INDEX;
- if (strncasecmp(dorade_field_name, "ch", 2) == 0) return CH_INDEX;
- if (strncasecmp(dorade_field_name, "ah", 2) == 0) return AH_INDEX;
- if (strncasecmp(dorade_field_name, "cv", 2) == 0) return CV_INDEX;
- if (strncasecmp(dorade_field_name, "av", 2) == 0) return AV_INDEX;
+ if (strncasecmp(dorade_field_name, "dx", 2) == 0) return DX_INDEX;
+ if (strncasecmp(dorade_field_name, "ch", 2) == 0) return CH_INDEX;
+ if (strncasecmp(dorade_field_name, "ah", 2) == 0) return AH_INDEX;
+ if (strncasecmp(dorade_field_name, "cv", 2) == 0) return CV_INDEX;
+ if (strncasecmp(dorade_field_name, "av", 2) == 0) return AV_INDEX;
+ if (strncasecmp(dorade_field_name, "vs", 2) == 0) return VS_INDEX;
+ if (strncasecmp(dorade_field_name, "vl", 2) == 0) return VL_INDEX;
+ if (strncasecmp(dorade_field_name, "vg", 2) == 0) return VG_INDEX;
+ if (strncasecmp(dorade_field_name, "vt", 2) == 0) return VT_INDEX;
+ if (strncasecmp(dorade_field_name, "ncp", 2) == 0) return NP_INDEX;
- fprintf(stderr, "Unknown DORADE type <%s>\n", dorade_field_name);
return -1;
}
+void prt_skipped_field_msg(char *dorade_field_name)
+{
+ char prtname[9];
+ int i, already_printed;
+#define MAXFIELDS 20
+ static int nskipped = 0;
+ static char skipped_list[MAXFIELDS][9];
+
+ /* Make sure name is a properly formed string. */
+ strncpy(prtname, dorade_field_name, 8);
+ prtname[8] = '\0';
+
+ /* We don't want to repeat messages for the same field, so keep a list of
+ * fields already printed.
+ */
+ already_printed = 0;
+ i = 0;
+ while (!already_printed && i < nskipped) {
+ if (strncmp(prtname, skipped_list[i], 2) == 0) already_printed = 1;
+ i++;
+ }
+ if (!already_printed) {
+ fprintf(stderr, "Unknown DORADE parameter type <%s> -- skipping.\n", prtname);
+ strcpy(skipped_list[nskipped], prtname);
+ nskipped++;
+ }
+}
+
+/* For date conversion routines. */
+static int daytab[2][13] = {
+ {0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334, 365},
+ {0, 31, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335, 366}
+};
+
+/*************************************************************/
+/* */
+/* julian */
+/* */
+/*************************************************************/
+static int julian(int year, int mo, int day)
+{
+/* Converts a calendar date (month, day, year) to a Julian date.
+ Returns:
+ Julian day.
+*/
+ int leap;
+
+ leap = (year%4 == 0 && year%100 != 0) || year%400 == 0;
+ return(day + daytab[leap][mo-1]);
+}
+
+/*************************************************************/
+/* */
+/* ymd */
+/* */
+/*************************************************************/
+static void ymd(int jday, int year, int *mm, int *dd)
+{
+ /* Input: jday, yyyy */
+ /* Output: mm, dd */
+ /* Copied from hdf_to_radar.c, written by Mike Kolander. */
+
+ int leap;
+ int i;
+
+ leap = (year%4 == 0 && year%100 != 0) || year%400 == 0;
+ for (i=0; daytab[leap][i]<jday; i++) continue;
+ *mm = i;
+ i--;
+ *dd = jday - daytab[leap][i];
+}
+
/* Secretly defined in uf_to_radar.c */
Volume *copy_sweeps_into_volume(Volume *new_volume, Volume *old_volume);
+/**********************************************************************/
+/* */
+/* RSL_dorade_to_radar */
+/* */
+/**********************************************************************/
Radar *RSL_dorade_to_radar(char *infile)
{
Radar *radar;
Sweep *sweep;
Ray *ray;
int iv, iray, iparam;
+ int range_bin1, gate_size;
+
+ int nbins, data_len, word_size;
+ short *ptr2bytes;
+ int *ptr4bytes;
+ float scale, offset, value;
+ int datum, missing_data_flag;
+ float prf;
+
+ float (*f)(Range x);
+ Range (*invf)(float x);
FILE *fp;
+ Comment_block *cb;
Volume_desc *vd;
Sensor_desc **sd;
Sweep_record *sr;
Radar_desc *rd;
Data_ray *dray;
Parameter_data *pd;
+ Ray_info *ray_info;
+ Platform_info *platform_info;
int nsweep;
int i;
int degree, minute;
float second;
+ int year, month, day, jday, jday_vol;
+
radar = NULL;
if (infile == NULL) {
- int save_fd;
- save_fd = dup(0);
- fp = fdopen(save_fd, "r");
+ int save_fd;
+ save_fd = dup(0);
+ fp = fdopen(save_fd, "r");
} else
if((fp=fopen(infile, "r"))==(FILE *)NULL) {
- perror(infile);
- return radar;
+ perror(infile);
+ return radar;
}
fp = uncompress_pipe(fp); /* Transparently, use gunzip. */
+ cb = dorade_read_comment_block(fp);
+
/**********************************************************************/
vd = dorade_read_volume_desc(fp); /* R E A D */
- if (radar_verbose_flag) dorade_print_volume_desc(vd); /* P R I N T */
+ if (radar_verbose_flag) dorade_print_volume_desc(vd); /* P R I N T */
/* R E A D */
sd = (Sensor_desc **) calloc(vd->nsensors, sizeof(Sensor_desc *));
for (i=0; i<vd->nsensors; i++) {
- sd[i] = dorade_read_sensor(fp);
+ sd[i] = dorade_read_sensor(fp);
}
/* P R I N T */
if (radar_verbose_flag) {
- for (i=0; i<vd->nsensors; i++) {
- fprintf(stderr, "============ S E N S O R # %d =====================\n", i);
- dorade_print_sensor(sd[i]);
- }
+ for (i=0; i<vd->nsensors; i++) {
+ fprintf(stderr, "============ S E N S O R # %d =====================\n", i);
+ dorade_print_sensor(sd[i]);
+ }
}
/* R E A D sweeps. */
if (vd->nsensors > 1) {
- fprintf(stderr, "RSL_dorade_to_radar: Unable to process for more than 1 sensor.\n");
- fprintf(stderr, "RSL_dorade_to_radar: Number of sensors is %d\n", vd->nsensors);
- return NULL;
+ fprintf(stderr, "RSL_dorade_to_radar: Unable to process for more than 1 sensor.\n");
+ fprintf(stderr, "RSL_dorade_to_radar: Number of sensors is %d\n", vd->nsensors);
+ return NULL;
}
/* Use sensor 0 for vitals. */
rd = sd[0]->radar_desc;
+ range_bin1 = sd[0]->cell_range_vector->range_cell[0];
+ gate_size = sd[0]->cell_range_vector->range_cell[1] - range_bin1;
radar = RSL_new_radar(MAX_RADAR_VOLUMES);
radar->h.month = vd->month;
radar->h.height = rd->altitude * 1000.0;
radar->h.spulse = 0; /* FIXME */
radar->h.lpulse = 0; /* FIXME */
-
+
+ year = vd->year;
+ month = vd->month;
+ day = vd->day;
+ jday_vol = julian(year, month, day);
+
+ /* TODO:
+ Get any threshold values present in parameter descriptors.
+
+ Note: The parameter descriptor contains an 8-character string which may
+ contain the, "Name of parameter upon which this parameter is thresholded".
+ According to documentation, if there is no thresholding parameter, the
+ string will be "NONE". In practice, this string has occasionally been
+ seen to contain only spaces for this condition. The corresponding
+ missing-data-flag for threshold may also vary, the nominal value being
+ -999, but has also been -32768.
+
+ Pseudo code:
+ nparam = rd->nparam_desc; [<--useable code, i.e., not pseudo]
+ create string array for threshold parameter names, size nparams.
+ for each param:
+ strcpy thresh param name from param desc to threshold param array.
+ if string is all blanks, replace with 'NONE'.
+ endfor
+ */
+
/* Begin volume code. */
/* We don't know how many sweeps per volume exist, until we read
* the file. So allocate a large number of pointers, hope we don't
*/
if (radar_verbose_flag)
- fprintf(stderr, "Number of parameters: %d\n", rd->nparam_desc);
+ fprintf(stderr, "Number of parameters: %d\n", rd->nparam_desc);
/* All the parameters are together, however, their order within
* the ray is not guarenteed. For instance, VE could appear after
* DM. For this we'll keep a list of parameter names and perform
* a (linear) search. The result will be an index into the RSL
- * volume array (radar->v[i]). It is likely that the order will
+ * volume array (radar->v[i]). It is likely that the order will be
* consistant within a file, therefore, we'll keep track of the index of
* our previous parameter type and begin the search from there; the next
* index should be a match.
*
* The dorade parameter names and the rsl mapping is:
*
- * Dorade: VE, DM, SW, DZ, ZDR, PHI, RHO, LDR, DX, CH, AH, CV, AV
- * RSL: VR, DM, SW, DZ, ZD, PH, RH, LR, *DX, *CH, *AH, *CV, *AV.
+ * Dorade: VE, DM, SW, DBZ, ZDR, PHI, RHO, LDR, DX, CH, AH, CV, AV
+ * RSL: VR, DM, SW, DZ, ZD, PH, RH, LR, *DX, *CH, *AH, *CV, *AV.
*
* * means this is a new RSL name.
*/
#define DORADE_MAX_SWEEP 20
nsweep = 0;
while((sr = dorade_read_sweep(fp, sd))) {
- for(iray = 0; iray < sr->nrays; iray++) {
- dray = sr->data_ray[iray];
-
- /* Now, loop through the parameters and fill the rsl structures. */
- for (iparam = 0; iparam < dray->nparam; iparam++) {
- pd = dray->parameter_data[iparam];
- iv = find_rsl_field_index(pd->name);
- if (radar->v[iv] == NULL) {
- radar->v[iv] = RSL_new_volume(DORADE_MAX_SWEEP); /* Expandable */
- } else if (nsweep >= radar->v[iv]->h.nsweeps) {
- /* Must expand the number of sweeps. */
- /* Expand by another DORADE_MAX_SWEEP. */
- if (radar_verbose_flag) {
- fprintf(stderr, "nsweeps (%d) exceeds radar->v[%d]->h.nsweeps (%d).\n", nsweep, iv, radar->v[iv]->h.nsweeps);
- fprintf(stderr, "Increading it to %d sweeps\n", radar->v[iv]->h.nsweeps+DORADE_MAX_SWEEP);
- }
- new_volume = RSL_new_volume(radar->v[iv]->h.nsweeps+DORADE_MAX_SWEEP);
- /* Look in uf_to_radar.c for 'copy_sweeps_into_volume' */
- new_volume = copy_sweeps_into_volume(new_volume, radar->v[iv]);
- radar->v[iv] = new_volume;
- }
-
- /* Allocate the ray and load the parameter data. */
- if ((sweep = radar->v[iv]->sweep[nsweep]) == NULL) {
- sweep = radar->v[iv]->sweep[nsweep] = RSL_new_sweep(sr->s_info->nrays);
- sweep->h.sweep_num = sr->s_info->sweep_num;
- sweep->h.elev = sr->s_info->fixed_angle;
+ for(iray = 0; iray < sr->nrays; iray++) {
+ dray = sr->data_ray[iray];
+
+ /* Now, loop through the parameters and fill the rsl structures. */
+ for (iparam = 0; iparam < dray->nparam; iparam++) {
+ pd = dray->parameter_data[iparam];
+ iv = find_rsl_field_index(pd->name);
+ if (iv < 0) {
+ prt_skipped_field_msg(pd->name);
+ continue;
+ }
+ if (radar->v[iv] == NULL) {
+ radar->v[iv] = RSL_new_volume(DORADE_MAX_SWEEP); /* Expandable */
+ } else if (nsweep >= radar->v[iv]->h.nsweeps) {
+ /* Must expand the number of sweeps. */
+ /* Expand by another DORADE_MAX_SWEEP. */
+ if (radar_verbose_flag) {
+ fprintf(stderr, "nsweeps (%d) exceeds radar->v[%d]->h.nsweeps (%d)."
+ "\n", nsweep, iv, radar->v[iv]->h.nsweeps);
+ fprintf(stderr, "Increasing it to %d sweeps\n",
+ radar->v[iv]->h.nsweeps+DORADE_MAX_SWEEP);
+ }
+ new_volume = RSL_new_volume(radar->v[iv]->h.nsweeps+DORADE_MAX_SWEEP);
+ /* Look in uf_to_radar.c for 'copy_sweeps_into_volume' */
+ new_volume = copy_sweeps_into_volume(new_volume, radar->v[iv]);
+ radar->v[iv] = new_volume;
+ }
+ /* Assign f and invf
+ switch (iv) . . .
+ * Or just use RSL_f_list[iv] and RSL_invf_list[iv] as in sweep.h below.
+ */
+
+ /* Allocate the ray and load the parameter data. */
+ if ((sweep = radar->v[iv]->sweep[nsweep]) == NULL) {
+ sweep = radar->v[iv]->sweep[nsweep] = RSL_new_sweep(sr->s_info->nrays);
+ sweep->h.sweep_num = sr->s_info->sweep_num;
+ sweep->h.elev = sr->s_info->fixed_angle;
sweep->h.beam_width = rd->horizontal_beam_width;
- sweep->h.vert_half_bw = radar->v[iv]->sweep[nsweep]->h.beam_width / 2.0;
- sweep->h.horz_half_bw = rd->horizontal_beam_width / 2.0;
- sweep->h.f = RSL_f_list[iv];
- sweep->h.invf = RSL_invf_list[iv];
- }
-
-
- if ((ray = sweep->ray[iray]) == NULL) {
- if (radar_verbose_flag)
- fprintf(stderr, "Allocating %d bins for ray %d\n", dray->data_len[iparam], iray);
- ray = sweep->ray[iray] = RSL_new_ray(dray->data_len[iparam]);
- }
-
- /* Copy the ray data into the RSL ray. */
+ sweep->h.vert_half_bw = rd->vertical_beam_width / 2.0;
+ /*
+ sweep->h.vert_half_bw = radar->v[iv]->sweep[nsweep]->h.beam_width / 2.0;
+ */
+ sweep->h.horz_half_bw = rd->horizontal_beam_width / 2.0;
+ sweep->h.f = RSL_f_list[iv];
+ sweep->h.invf = RSL_invf_list[iv];
+ }
+
+
+ data_len = dray->data_len[iparam];
+ word_size = dray->word_size[iparam];
+ if (word_size != 2 && word_size != 4) {
+ fprintf(stderr,"RSL_dorade_to_radar: dray->word_size[%d] = %d\n",
+ iparam, dray->word_size[iparam]);
+ fprintf(stderr,"Expected value is 2 or 4.\n");
+ return NULL;
+ }
+ nbins = data_len / word_size;
+
+ if ((ray = sweep->ray[iray]) == NULL) {
+ if (radar_verbose_flag)
+ fprintf(stderr, "Allocating %d bins for ray %d\n",
+ dray->data_len[iparam], iray);
+ ray = sweep->ray[iray] = RSL_new_ray(nbins);
+ }
+
+ /* Load ray header. */
+
+ ray_info = dray->ray_info;
+ jday = ray_info->jday;
+ if (jday != jday_vol) {
+ if (jday > 0 && jday < 367) {
+ /* Recompute year month day */
+ if (jday < jday_vol) year++;
+ ymd(jday, year, &month, &day);
+ jday_vol = jday;
+ }
+ else { /* Invalid jday */
+ }
+ }
+ ray->h.year = year;
+ ray->h.month = month;
+ ray->h.day = day;
+ ray->h.hour = ray_info->hour;
+ ray->h.minute = ray_info->minute;
+ ray->h.sec = ray_info->second + ray_info->msec / 1000.;
+ ray->h.azimuth = ray_info->azimuth;
+ ray->h.ray_num = iray + 1;
+ ray->h.elev = ray_info->elevation;
+ ray->h.elev_num = ray_info->sweep_num;
+ ray->h.unam_rng = rd->unambiguous_range;
+ ray->h.nyq_vel = rd->unambiguous_velocity;
+ ray->h.azim_rate = ray_info->scan_rate;
+ /* TODO
+ ray->h.fix_angle = ;
+ */
+ ray->h.range_bin1 = range_bin1;
+ ray->h.gate_size = gate_size;
+ platform_info = dray->platform_info;
+ ray->h.pitch = platform_info->pitch;
+ ray->h.roll = platform_info->roll;
+ ray->h.heading = platform_info->heading;
+ ray->h.pitch_rate = platform_info->pitch_rate;
+ ray->h.heading_rate = platform_info->heading_rate;
+ ray->h.lat = platform_info->latitude;
+ ray->h.lon = platform_info->longitude;
+ ray->h.alt = platform_info->altitude;
+ ray->h.vel_east = platform_info->ew_speed;
+ ray->h.vel_north = platform_info->ns_speed;
+ ray->h.vel_up = platform_info->v_speed;
+ /* TODO: Does DORADE have Doppler velocity res?
+ A: No.
+ ray->h.vel_res = N/A
+ */
+ ray->h.beam_width = rd->horizontal_beam_width;
+ /* TODO
+ ray->h.frequency = @Radar Descriptor
+ Get parm_desc.xmit_freq
+ Find first set-bit, use frequency from corresponding index in rad_desc frequenies.
+
+ ray->h.pulse_count = Is it number samples used? @Param desc. Probably not.
+
+ * The following are DONE *
+ ray->h.pulse_width = @Parameter Descriptor--not same--it's in meters--
+ we use micro-sec. They're using pulse length.
+ pulse width = pulse length/c
+ ray->h.wavelength = Can compute using Nyquist velocity and PRF or PRT:
+ wl = 4*vmax/PRF or wl = 4*vmax*PRT
+ ray->h.prf = Can compute if nothing else: prf = c/(2*Rmax)
+ */
+ /* DORADE is actually using pulse length. Convert to pulse width. */
+ ray->h.pulse_width = (sd[0]->p_desc[iparam]->pulse_width/
+ RSL_SPEED_OF_LIGHT) * 1e6;
+ prf = RSL_SPEED_OF_LIGHT/(2.*rd->unambiguous_range*1000.);
+ ray->h.prf = prf;
+ ray->h.wavelength = (4.*rd->unambiguous_velocity)/prf;
+ ray->h.nbins = nbins;
+ ray->h.f = RSL_f_list[iv];
+ ray->h.invf = RSL_invf_list[iv];
+
+ /* Copy the ray data into the RSL ray. */
/* .... fill here .... */
+
+ /* Assign pointer to data.
+ * Get datum using word size and proper cast.
+ * Convert and store in rsl ray->range.
+ * Increment pointer to data based on word size.
+ */
+
+ ptr2bytes = (short *) pd->data;
+ ptr4bytes = (int *) pd->data;
+ scale = sd[0]->p_desc[iparam]->scale_factor;
+ offset = sd[0]->p_desc[iparam]->offset_factor;
+ missing_data_flag = sd[0]->p_desc[iparam]->missing_data_flag;
+
+ for (i=0; i<nbins; i++) {
+ if (word_size == 2) datum = *ptr2bytes++;
+ else datum = *ptr4bytes++;
+ /*
+ TODO: If there is a threshold parameter for this parameter then
+ apply threshold value. I think threshold works like this: If there's a
+ threshold parameter, as with VT (threshold parm = NCP), then use
+ threshold value from that parameter unless it is the missing data value.
+ */
+ if (datum != missing_data_flag) {
+ value = ((float)datum - offset) / scale;
}
+ else value = BADVAL;
+
+ ray->range[i] = ray->h.invf(value);
}
- nsweep++;
- if (radar_verbose_flag) fprintf(stderr, "______NEW SWEEP__<%d>____\n", nsweep);
- /* Save for loading into volume structure. */
- dorade_free_sweep(sr);
+
+ if (iray == 0) {
+ radar->v[iv]->h.nsweeps = nsweep + 1;
+ radar->v[iv]->h.f = RSL_f_list[iv];
+ radar->v[iv]->h.invf = RSL_invf_list[iv];
+ }
+ }
+ }
+ nsweep++;
+ if (radar_verbose_flag) fprintf(stderr, "______NEW SWEEP__<%d>____\n", nsweep);
+ /* Save for loading into volume structure. */
+ dorade_free_sweep(sr);
}
/* The following avoids a broken pipe message, since a VOLD at the end
return radar;
}
-
INCLUDES = -I$(prefix)/include
LOCAL_LIB = ../.libs/librsl.a
LDADD = @LIBS@ $(LOCAL_LIB)
-bin_PROGRAMS = any_to_gif any_to_uf
+bin_PROGRAMS = any_to_gif any_to_uf qlook
any_to_gif_LDFLAGS = -static
any_to_uf_LDFLAGS = -static
# Additional program to build but not install
noinst_PROGRAMS = any_to_ppm any_to_ufgz bscan \
cappi_image dorade_main killer_sweep \
kwaj_subtract_one_day lassen_to_gif print_hash_table \
- print_header_info qlook sector test_get_win \
+ print_header_info sector test_get_win \
wsr88d_to_gif wsr_hist_uf_test
#include <string.h>
#include <stdlib.h>
#include <unistd.h>
+#include <libgen.h>
+#define USE_RSL_VARS
#include "rsl.h"
#define ZDR_WIDTH 10
int verbose;
-void rebin_ray(Ray *r, int width);
-void rebin_sweep(Sweep *s, int width);
-void rebin_volume(Volume *v, int width);
+void RSL_rebin_ray(Ray *r, int width);
+void RSL_rebin_sweep(Sweep *s, int width);
+void RSL_rebin_volume(Volume *v, int width);
void adjust_gate_size(Radar *radar, float gate_size_adjustment);
void process_args(int argc, char **argv, char *in_file, int *verbose,
- char *site_id, char *tape_id,
- int *qc_reflectivity, int *total_reflectivity,
- int *differential_reflectivity,
- int *velocity, int *spectral_width,
- int *make_gif, int *make_pgm, int *make_bscan, int *make_uf,
- int *num_sweeps, float *dbz_offset,
- int *xdim, int *ydim, float *range,
- float *gate_size_adjustment, int *print_azim);
+ char *site_id, char *tape_id,
+ int *qc_reflectivity, int *total_reflectivity,
+ int *differential_reflectivity,
+ int *velocity, int *spectral_width,
+ int *make_gif, int *make_pgm, int *make_bscan, int *make_uf,
+ int *num_sweeps, float *dbz_offset,
+ int *xdim, int *ydim, float *range,
+ float *gate_size_adjustment, int *print_azim,
+ char *gifdir, char *pgmdir, char *ufdir);
+
+void make_pathname(char *filename, char *dir, char *pathname)
+{
+ /* Make pathname by combining directory name, if given, with filename. */
+
+ if (*dir == NULL) {
+ strcpy(pathname, filename);
+ }
+ else {
+ strcpy(pathname, dir);
+ strcat(pathname, "/");
+ strcat(pathname, filename);
+ }
+}
/***********************************************************************/
/* This program uses the NASA TRMM Office Radar Software Library (RSL) */
/***********************************************************************/
main(int argc, char **argv) {
- Radar *radar;
- Volume *dz_volume, *vr_volume;
- Volume *dr_volume, *zt_volume;
- Volume *sw_volume, *qc_volume, *volume;
-
- Sweep *sweep;
- Sweep *dz_sweep, *qc_sweep;
- Sweep *dr_sweep, *zt_sweep;
- Sweep *vr_sweep, *sw_sweep;
-
- Ray *ray;
- int reflectivity, qc_reflectivity, nvolumes;
- int differential_reflectivity, total_reflectivity;
- int velocity, spectral_width;
- int make_catalog, make_gif, make_pgm;
- int make_uf, make_bscan;
- int xdim, ydim, num_sweeps;
- int i,j,k,l,n, verbose, kk;
- int month, day, year, hour, min;
- int ncbins, width;
- int print_azim;
- float scale, dbz_offset;
- float latitude, longitude;
-
- float sec;
- float maxr;
- float nyquist, max_range, gate_size_adjustment;
-
- char tape_id[100];
- char in_file[100], site_id[100];
- char filename[100], outfile[100], nexfile[100];
- char command[100], file_prefix[100], file_suffix[3];
- char dir_string[100], red[120], grn[120], blu[120];
- char time_string[14], site_string[10],img_base[20];
- FILE *fp;
+ Radar *radar;
+ Volume *dz_volume, *vr_volume;
+ Volume *dr_volume, *zt_volume;
+ Volume *sw_volume, *qc_volume, *volume;
+
+ Sweep *sweep;
+ Sweep *dz_sweep, *qc_sweep;
+ Sweep *dr_sweep, *zt_sweep;
+ Sweep *vr_sweep, *sw_sweep;
+
+ Ray *ray;
+ int reflectivity, qc_reflectivity, nvolumes;
+ int differential_reflectivity, total_reflectivity;
+ int velocity, spectral_width;
+ int make_catalog, make_gif, make_pgm;
+ int make_uf, make_bscan;
+ int xdim, ydim, num_sweeps;
+ int i,j,k,l,n, verbose, kk;
+ int month, day, year, hour, min;
+ int ncbins, width;
+ int print_azim;
+ float scale, dbz_offset;
+ float latitude, longitude;
+
+ float sec;
+ float maxr;
+ float nyquist, max_range, gate_size_adjustment;
+
+ char tape_id[100];
+ char in_file[100], site_id[100];
+ char filename[100], outfile[100], nexfile[100];
+ char command[100], file_prefix[100], file_suffix[3];
+ char dir_string[100], red[120], grn[120], blu[120];
+ char time_string[14], site_string[10],img_base[20];
+ char pathname[256], gifdir[100], pgmdir[100], ufdir[100];
+ char *inpath;
+ FILE *fp;
/* Set default values */
- strcpy(site_id, "KMLB");
- strcpy(tape_id, "WSR88D");
-
- verbose = FALSE;
- reflectivity = TRUE;
- qc_reflectivity = FALSE;
- total_reflectivity = FALSE;
- differential_reflectivity = FALSE;
- velocity = FALSE;
- spectral_width = FALSE;
- make_gif = FALSE;
- make_pgm = FALSE;
- make_catalog = FALSE;
- make_uf = FALSE;
- print_azim = FALSE;
- num_sweeps = 1;
- xdim = ydim = 400;
- maxr = 200.;
- dbz_offset = 0.0;
- gate_size_adjustment = 1.0;
-
+ /*strcpy(site_id, "KMLB");*/
+ strcpy(tape_id, "WSR88D");
+
+ verbose = FALSE;
+ reflectivity = TRUE;
+ qc_reflectivity = FALSE;
+ total_reflectivity = FALSE;
+ differential_reflectivity = FALSE;
+ velocity = FALSE;
+ spectral_width = FALSE;
+ make_gif = FALSE;
+ make_pgm = FALSE;
+ make_catalog = FALSE;
+ make_uf = FALSE;
+ print_azim = FALSE;
+ num_sweeps = 1;
+ xdim = ydim = 400;
+ maxr = 200.;
+ dbz_offset = 0.0;
+ gate_size_adjustment = 1.0;
+ *gifdir = *pgmdir = *ufdir = '\0';
+ *site_id = '\0';
+
/* Process command_line arguments */
- process_args(argc, argv, in_file, &verbose,
- site_id, tape_id,
- &qc_reflectivity, &total_reflectivity,
- &differential_reflectivity,
- &velocity, &spectral_width,
- &make_gif, &make_pgm, &make_bscan, &make_uf,
- &num_sweeps, &dbz_offset,
- &xdim, &ydim, &maxr, &gate_size_adjustment,
- &print_azim);
-
+ process_args(argc, argv, in_file, &verbose,
+ site_id, tape_id,
+ &qc_reflectivity, &total_reflectivity,
+ &differential_reflectivity,
+ &velocity, &spectral_width,
+ &make_gif, &make_pgm, &make_bscan, &make_uf,
+ &num_sweeps, &dbz_offset,
+ &xdim, &ydim, &maxr, &gate_size_adjustment,
+ &print_azim, gifdir, pgmdir, ufdir);
+
+/* If site_id empty, use first 4 characters of file name. */
+
+ if (*site_id == 0) {
+ inpath = strdup(in_file);
+ strncpy(site_id, basename(inpath), 4);
+ site_id[4] = '\0';
+ if (strcmp(site_id, "KWA0") == 0) strcpy(site_id, "KWAJ");
+ }
+
/* Be a chatty Kathy? */
if(verbose)
/* Read data into radar */
- if(verbose) printf("Calling any_format_to_radar\n");
- radar = RSL_anyformat_to_radar(in_file, site_id);
- if(verbose) printf("Called any_format_to_radar\n");
- if(radar==NULL) {
- if (verbose)
- printf("No radar loaded, bye\n");
- exit(-1);
- }
+ if(verbose) printf("Calling any_format_to_radar\n");
+ radar = RSL_anyformat_to_radar(in_file, site_id);
+ if(verbose) printf("Called any_format_to_radar\n");
+ if(radar==NULL) {
+ if (verbose)
+ printf("No radar loaded, bye\n");
+ exit(-1);
+ }
/* Print command line parameters */
- if(verbose) {
- printf("Site ID = %s\n",site_id);
- printf("Tape ID = %s\n",tape_id);
- printf("Do reflectivity = %d\n",reflectivity);
- printf("Do qc_reflectivity = %d\n",qc_reflectivity);
- printf("Do differential_reflectivity = %d\n",
- differential_reflectivity);
- printf("Do total_reflectivity = %d\n",
- total_reflectivity);
- printf("Do qc_reflectivity = %d\n",qc_reflectivity);
- printf("Do velocity = %d\n",velocity);
- printf("Do spectral_width = %d\n",spectral_width);
- printf("Make gif = %d\n",make_gif);
- printf("Make pgm = %d\n",make_pgm);
- printf("Make UF file = %d\n",make_uf);
- printf("dBZ Offset = %.2f\n",dbz_offset);
- printf("Gate Size Adjust = %.2f\n",gate_size_adjustment);
- printf("Print Azimuths = %d\n",print_azim);
- }
+ if(verbose) {
+ printf("Site ID = %s\n",site_id);
+ printf("Tape ID = %s\n",tape_id);
+ printf("Do reflectivity = %d\n",reflectivity);
+ printf("Do qc_reflectivity = %d\n",qc_reflectivity);
+ printf("Do differential_reflectivity = %d\n",
+ differential_reflectivity);
+ printf("Do total_reflectivity = %d\n",
+ total_reflectivity);
+ printf("Do qc_reflectivity = %d\n",qc_reflectivity);
+ printf("Do velocity = %d\n",velocity);
+ printf("Do spectral_width = %d\n",spectral_width);
+ printf("Make gif = %d\n",make_gif);
+ printf("Make pgm = %d\n",make_pgm);
+ printf("Make UF file = %d\n",make_uf);
+ if (ufdir != NULL) printf("UF output dir = %s\n",ufdir);
+ if (gifdir != NULL) printf("GIF output dir = %s\n",gifdir);
+ if (pgmdir != NULL) printf("PGM output dir = %s\n",pgmdir);
+ printf("dBZ Offset = %.2f\n",dbz_offset);
+ printf("Gate Size Adjust = %.2f\n",gate_size_adjustment);
+ printf("Print Azimuths = %d\n",print_azim);
+ }
/*
If Gate Size Adjustment is not unity, then we must change the
following:
old_gs = radar->v[i]->sweep[sweepIndex]->ray[rayIndex=]->h.gate_size
- radar->v[i]->sweep[sweepIndex]->ray[rayIndex=]->h.gate_size =
- old_gs*gate_size_adjustment
+ radar->v[i]->sweep[sweepIndex]->ray[rayIndex=]->h.gate_size =
+ old_gs*gate_size_adjustment
Here are some comments from Sandra Yuter on the necessity of this fix.
> I dug into the revelant code and it looks like we can do a relatively
- > simple workaround for the SIGMET raw product file range bin size
- > errors for the RHB data pulses widths of 0.5 usec and 2.0 usec as follows.
- >
- > Since we are all converting from sigmet to UF I suggest we resize
- > the range bin size values in the ray headers in the qlook step
- > where the sigmet to UF conversion occurs.
- >
- > The resize requires only 1 additional line of code (I have included
- > a few others for context) in qlook.c
- >
- > rangeToFirstGate = 0.001 *
- > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.range_bin1;
- > gateSize = 0.001 *
- > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.gate_size;
- > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.gate_size=
- > gateSize*0.6*1000;
- >
- > I have used 0.6 adjustment factor since that is 75/125 which corresponds
- > to the error in my 0.5 usec data, for the SUR scans, this adjustment
- > factor is 133.33/125 or 1.067.
-
- The following is from Joe Holmes from SIGMET
-
- >
- > I think you have experienced a problem with the RVP7 range resolution
- > configuration. Both in IRIS and in the RVP7 you manually type in
- > the range resolution. The RVP7 allows a separate resolution for
- > each pulsewidth, while IRIS only allows 1 value. There is no feedback
- > if these values are not typed in the same. Based on setup information
- > we have here from the RH Brown from Oct 23, 1998, you had the following
- > configuration:
- >
- > RVP7:
- > 0 0.50 usec 75.0m
- > 1 0.80 usec 125.0m
- > 2 1.39 usec 125.0m
- > 3 2.00 usec 133.3m
- >
- > IRIS: 125.0 meters
- >
- > I think the error at PW#0 was corrected a while back, but
- > the error in PW#3 was never corrected. Next time someone is
- > at the ship, they should check this, fix the long pulse, and remake
- > the bandpass filter for the long pulse.
- >
- > In the short term, you can correct the error by taking all your
- > long pulse data and changing the header to correctly document the
- > range resolution. We have a program to do this, it is called "change_raw".
- > The source is on any IRIS system, which was installed with the
- > source, headers, and objects turned on. It is in the
- > ${IRIS_ROOT}utils/examples directory. We can supply you with
- > a compiled version of this program, if you want. Available platforms
- > are Linux, HP-UX, and IRIX.
+ > simple workaround for the SIGMET raw product file range bin size
+ > errors for the RHB data pulses widths of 0.5 usec and 2.0 usec as follows.
+ >
+ > Since we are all converting from sigmet to UF I suggest we resize
+ > the range bin size values in the ray headers in the qlook step
+ > where the sigmet to UF conversion occurs.
+ >
+ > The resize requires only 1 additional line of code (I have included
+ > a few others for context) in qlook.c
+ >
+ > rangeToFirstGate = 0.001 *
+ > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.range_bin1;
+ > gateSize = 0.001 *
+ > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.gate_size;
+ > radar->v[i]->sweep[sweepIndex]->ray[rayIndex]->h.gate_size=
+ > gateSize*0.6*1000;
+ >
+ > I have used 0.6 adjustment factor since that is 75/125 which corresponds
+ > to the error in my 0.5 usec data, for the SUR scans, this adjustment
+ > factor is 133.33/125 or 1.067.
+
+ The following is from Joe Holmes from SIGMET
+
+ >
+ > I think you have experienced a problem with the RVP7 range resolution
+ > configuration. Both in IRIS and in the RVP7 you manually type in
+ > the range resolution. The RVP7 allows a separate resolution for
+ > each pulsewidth, while IRIS only allows 1 value. There is no feedback
+ > if these values are not typed in the same. Based on setup information
+ > we have here from the RH Brown from Oct 23, 1998, you had the following
+ > configuration:
+ >
+ > RVP7:
+ > 0 0.50 usec 75.0m
+ > 1 0.80 usec 125.0m
+ > 2 1.39 usec 125.0m
+ > 3 2.00 usec 133.3m
+ >
+ > IRIS: 125.0 meters
+ >
+ > I think the error at PW#0 was corrected a while back, but
+ > the error in PW#3 was never corrected. Next time someone is
+ > at the ship, they should check this, fix the long pulse, and remake
+ > the bandpass filter for the long pulse.
+ >
+ > In the short term, you can correct the error by taking all your
+ > long pulse data and changing the header to correctly document the
+ > range resolution. We have a program to do this, it is called "change_raw".
+ > The source is on any IRIS system, which was installed with the
+ > source, headers, and objects turned on. It is in the
+ > ${IRIS_ROOT}utils/examples directory. We can supply you with
+ > a compiled version of this program, if you want. Available platforms
+ > are Linux, HP-UX, and IRIX.
*/
-
- if(gate_size_adjustment != 1.0) {
- printf("Adjusting Gate Size by Factor: %.3f\n",gate_size_adjustment);
- adjust_gate_size(radar, gate_size_adjustment);
- }
-
+
+ if(gate_size_adjustment != 1.0) {
+ printf("Adjusting Gate Size by Factor: %.3f\n",gate_size_adjustment);
+ adjust_gate_size(radar, gate_size_adjustment);
+ }
+
/*
Create the filename prefix. Consists of the Site ID,
and time string (YYMMDD_hhmm). The file suffix is
like MIME type (e.g, .gif, .pgm, .uf, etc.)
*/
- sprintf(time_string,"%4d/%2.2d%2.2d %2.2d:%2.2d UTC",
- radar->h.year, radar->h.month, radar->h.day,
- radar->h.hour, radar->h.minute);
+ sprintf(time_string,"%4d/%2.2d%2.2d %2.2d:%2.2d UTC",
+ radar->h.year, radar->h.month, radar->h.day,
+ radar->h.hour, radar->h.minute);
/*
Determine the location (lat/lon) of the radar.
*/
- latitude = radar->h.latd + radar->h.latm/60. + radar->h.lats/3600.;
- longitude = radar->h.lond + radar->h.lonm/60. + radar->h.lons/3600.;
+ latitude = radar->h.latd + radar->h.latm/60. + radar->h.lats/3600.;
+ longitude = radar->h.lond + radar->h.lonm/60. + radar->h.lons/3600.;
- printf("%s %s %s %.6f %.6f \n",
- in_file, radar->h.radar_name, time_string, longitude, latitude);
+ printf("%s %s %s %.6f %.6f \n",
+ in_file, radar->h.radar_name, time_string, longitude, latitude);
- sprintf(time_string,"%4d_%2.2d%2.2d_%2.2d%2.2d",
- radar->h.year, radar->h.month, radar->h.day,
- radar->h.hour, radar->h.minute);
+ sprintf(time_string,"%4d_%2.2d%2.2d_%2.2d%2.2d",
+ radar->h.year, radar->h.month, radar->h.day,
+ radar->h.hour, radar->h.minute);
/*
Print the radar/volume info.
*/
- if(verbose) {
- for(i=0; i< radar->h.nvolumes; i++) {
- if(radar->v[i] != NULL) {
- printf("Vol[%2.2d] has %d sweeps\n",i,radar->v[i]->h.nsweeps);
- } else {
- printf("Vol[%2.2d] == NULL\n",i);
- }
- }
- printf("--------------------------------------------\n");
- printf("Number of volumes in radar: %d\n",radar->h.nvolumes);
- }
-
- /* DZ_INDEX */
- if(radar->v[DZ_INDEX] == NULL) {
- printf("DZ_INDEX == NULL\n");
- reflectivity = FALSE;
- } else {
- dz_volume = radar->v[DZ_INDEX];
- if(verbose) printf("Number of sweeps in dz_volume = %d\n",
- dz_volume->h.nsweeps);
- }
-
- /* CZ_INDEX */
- if(radar->v[CZ_INDEX] == NULL) {
- if(verbose) printf("CZ_INDEX == NULL\n");
- qc_reflectivity = FALSE;
- } else {
- qc_volume = radar->v[CZ_INDEX];
- if(verbose) printf("Number of sweeps in qc_volume = %d\n",
- qc_volume->h.nsweeps);
- }
-
- /* ZT_INDEX */
- if(radar->v[ZT_INDEX] == NULL) {
- if(verbose) printf("ZT_INDEX == NULL\n");
- total_reflectivity = FALSE;
- } else {
- zt_volume = radar->v[ZT_INDEX];
- if(verbose) printf("Number of sweeps in zt_volume = %d\n",
- zt_volume->h.nsweeps);
- }
- /* ZD_INDEX */
- if(radar->v[ZD_INDEX] == NULL) {
- if(verbose) printf("ZD_INDEX == NULL\n");
- differential_reflectivity = FALSE;
- } else {
- dr_volume = radar->v[ZD_INDEX];
- if(verbose) printf("Number of sweeps in dr_volume = %d\n",
- dr_volume->h.nsweeps);
- }
- /* VR_INDEX */
- if(radar->v[VR_INDEX] == NULL) {
- if(verbose) printf("VR_INDEX == NULL\n");
- velocity = FALSE;
- } else {
- vr_volume = radar->v[VR_INDEX];
- if(verbose) printf("Number of sweeps in vr_volume = %d\n",
- vr_volume->h.nsweeps);
- }
-
- /* SW_INDEX */
- if(radar->v[SW_INDEX] == NULL) {
- if(verbose) printf("SW_INDEX == NULL\n");
- spectral_width = FALSE;
- } else {
- sw_volume = radar->v[SW_INDEX];
- if(verbose) printf("Number of sweeps in sw_volume = %d\n",
- sw_volume->h.nsweeps);
- }
- if(verbose) printf("--------------------------------------------\n");
+ if(verbose) {
+ for(i=0; i< radar->h.nvolumes; i++) {
+ if(radar->v[i] != NULL) {
+ printf("Vol[%2.2d] has %d sweeps\n",i,radar->v[i]->h.nsweeps);
+ } else {
+ printf("Vol[%2.2d] == NULL\n",i);
+ }
+ }
+ printf("--------------------------------------------\n");
+ printf("Number of volumes in radar: %d\n",radar->h.nvolumes);
+ }
+
+ /* DZ_INDEX */
+ if(radar->v[DZ_INDEX] == NULL) {
+ printf("DZ_INDEX == NULL\n");
+ reflectivity = FALSE;
+ } else {
+ dz_volume = radar->v[DZ_INDEX];
+ if(verbose) printf("Number of sweeps in dz_volume = %d\n",
+ dz_volume->h.nsweeps);
+ }
+
+ /* CZ_INDEX */
+ if(radar->v[CZ_INDEX] == NULL) {
+ if(verbose) printf("CZ_INDEX == NULL\n");
+ qc_reflectivity = FALSE;
+ } else {
+ qc_volume = radar->v[CZ_INDEX];
+ if(verbose) printf("Number of sweeps in qc_volume = %d\n",
+ qc_volume->h.nsweeps);
+ }
+
+ /* ZT_INDEX */
+ if(radar->v[ZT_INDEX] == NULL) {
+ if(verbose) printf("ZT_INDEX == NULL\n");
+ total_reflectivity = FALSE;
+ } else {
+ zt_volume = radar->v[ZT_INDEX];
+ if(verbose) printf("Number of sweeps in zt_volume = %d\n",
+ zt_volume->h.nsweeps);
+ }
+ /* DR_INDEX */
+ if(radar->v[DR_INDEX] == NULL) {
+ if(verbose) printf("DR_INDEX == NULL\n");
+ differential_reflectivity = FALSE;
+ } else {
+ dr_volume = radar->v[DR_INDEX];
+ if(verbose) printf("Number of sweeps in dr_volume = %d\n",
+ dr_volume->h.nsweeps);
+ }
+ /* VR_INDEX */
+ if(radar->v[VR_INDEX] == NULL) {
+ if(verbose) printf("VR_INDEX == NULL\n");
+ velocity = FALSE;
+ } else {
+ vr_volume = radar->v[VR_INDEX];
+ if(verbose) printf("Number of sweeps in vr_volume = %d\n",
+ vr_volume->h.nsweeps);
+ }
+
+ /* SW_INDEX */
+ if(radar->v[SW_INDEX] == NULL) {
+ if(verbose) printf("SW_INDEX == NULL\n");
+ spectral_width = FALSE;
+ } else {
+ sw_volume = radar->v[SW_INDEX];
+ if(verbose) printf("Number of sweeps in sw_volume = %d\n",
+ sw_volume->h.nsweeps);
+ }
+ if(verbose) printf("--------------------------------------------\n");
/*
Print out the elevation angles
*/
- if(verbose) {
- if(reflectivity) {
- printf("Reflectivity Tilts\n");
- printf("----------------------\n");
- if(dz_volume != NULL) {
- for(i=0; i<dz_volume->h.nsweeps; i++) {
- sweep = dz_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- printf("Tilt %2d, Elev=%4.1f\n",i,sweep->h.elev);
- }
- printf("----------------------\n");
- }
- }
- }
- /*
- Print out the values of the rays in each sweep requsted
- */
-
- if(print_azim) {
- printf("Ray angles\n");
- if(reflectivity) {
- for(i=0; i<dz_volume->h.nsweeps; i++) {
- if(dz_volume->sweep[i] != NULL) sweep = dz_volume->sweep[i];
- printf("Elevation angle: %f\n",sweep->h.elev);
- printf("Number of rays in sweep[%d] = %d\n",i,sweep->h.nrays);
-
- for(j=0; j<sweep->h.nrays-1; j++) {
- if(sweep->ray[j] != NULL) {
- ray = sweep->ray[j];
- printf("%d: %7.2f\n ",j,sweep->ray[j]->h.azimuth);
- }
- }
- }
- }
- }
+ if(verbose) {
+ if(reflectivity) {
+ printf("Reflectivity Tilts\n");
+ printf("----------------------\n");
+ if(dz_volume != NULL) {
+ for(i=0; i<dz_volume->h.nsweeps; i++) {
+ sweep = dz_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ printf("Tilt %2d, Elev=%4.1f\n",i,sweep->h.elev);
+ }
+ printf("----------------------\n");
+ }
+ }
+ }
+ /*
+ Print out the values of the rays in each sweep requsted
+ */
+
+ if(print_azim) {
+ printf("Ray angles\n");
+ if(reflectivity) {
+ for(i=0; i<dz_volume->h.nsweeps; i++) {
+ if(dz_volume->sweep[i] != NULL) sweep = dz_volume->sweep[i];
+ printf("Elevation angle: %f\n",sweep->h.elev);
+ printf("Number of rays in sweep[%d] = %d\n",i,sweep->h.nrays);
+
+ for(j=0; j<sweep->h.nrays-1; j++) {
+ if(sweep->ray[j] != NULL) {
+ ray = sweep->ray[j];
+ printf("%d: %7.2f\n ",j,sweep->ray[j]->h.azimuth);
+ }
+ }
+ }
+ }
+ }
/*
Write out some volume statistics
*/
- if(verbose) {
- printf("******************* Volume Statistics *******************\n");
- if(reflectivity) {
- sweep = RSL_get_first_sweep_of_volume(dz_volume);
- if(sweep != NULL) {
- printf("Number rays = %d\n", sweep->h.nrays);
- printf("Beam width = %.2f deg\n", sweep->h.beam_width);
- ray = RSL_get_first_ray_of_sweep(sweep);
- if(ray!= NULL) {
- max_range = ray->h.range_bin1/1000.0 +
- ray->h.nbins*ray->h.gate_size/1000.0;
- printf("Number of bins = %d\n",ray->h.nbins);
- printf("Max range = %.1f km\n", max_range);
- printf("Range to first bin = %d m\n",ray->h.range_bin1);
- printf("Gate size = %d m\n",ray->h.gate_size);
- printf("Pulse Rep. Freq. = %d Hz\n",ray->h.prf);
- printf("Pulse width = %.2f us\n",ray->h.pulse_width);
- printf("Wavelength = %.2f m\n",ray->h.wavelength);
- printf("Frequency = %.2f \n",ray->h.frequency);
- }
- }
- }
- }
+ if(verbose) {
+ printf("******************* Volume Statistics *******************\n");
+ if(reflectivity) {
+ sweep = RSL_get_first_sweep_of_volume(dz_volume);
+ if(sweep != NULL) {
+ printf("Number rays = %d\n", sweep->h.nrays);
+ printf("Beam width = %.2f deg\n", sweep->h.beam_width);
+ ray = RSL_get_first_ray_of_sweep(sweep);
+ if(ray!= NULL) {
+ max_range = ray->h.range_bin1/1000.0 +
+ ray->h.nbins*ray->h.gate_size/1000.0;
+ printf("Number of bins = %d\n",ray->h.nbins);
+ printf("Max range = %.1f km\n", max_range);
+ printf("Range to first bin = %d m\n",ray->h.range_bin1);
+ printf("Gate size = %d m\n",ray->h.gate_size);
+ printf("Pulse Rep. Freq. = %d Hz\n",ray->h.prf);
+ printf("Pulse width = %.2f us\n",ray->h.pulse_width);
+ printf("Wavelength = %.2f m\n",ray->h.wavelength);
+ printf("Frequency = %.2f \n",ray->h.frequency);
+ }
+ }
+ }
+ }
/*
- Add a dBZ offset if requested. The offset can be positive or negative.
- if(dbz_offset != 0.0) {
- printf("Adding dbz_offset to dz_volume: %.2f\n", dbz_offset);
- RSL_add_dbzoffset_to_volume(dz_volume, dbz_offset);
- printf("Added dbz_offset to dz_volume: %.2f\n", dbz_offset);
- exit(0);
- }
+ Add a dBZ offset if requested. The offset can be positive or negative.
+ if(dbz_offset != 0.0) {
+ printf("Adding dbz_offset to dz_volume: %.2f\n", dbz_offset);
+ RSL_add_dbzoffset_to_volume(dz_volume, dbz_offset);
+ printf("Added dbz_offset to dz_volume: %.2f\n", dbz_offset);
+ exit(0);
+ }
*/
/*
*/
- /* CZ_INDEX */
- if(qc_reflectivity) {
- if(verbose) printf("Loading refl colortable...\n");
- RSL_load_refl_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = qc_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"qc_%s_%2.2d.%s", time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"qc_%s_%2.2d.%s", time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
- /* DZ_INDEX */
- if(reflectivity) {
- if(verbose) printf("Loading refl colortable...\n");
- RSL_load_refl_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = dz_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"dz_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"dz_%s_%2.2d.%s", time_string,i,file_suffix);
-/*
- sprintf(filename,"dz_%s.%s.%s", time_string,in_file,file_suffix);
+ /* CZ_INDEX */
+ if(qc_reflectivity) {
+ if(verbose) printf("Loading refl colortable...\n");
+ RSL_load_refl_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = qc_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"qc_%s_%2.2d.%s", time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"qc_%s_%2.2d.%s", time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
+ /* DZ_INDEX */
+ if(reflectivity) {
+ if(verbose) printf("Loading refl colortable...\n");
+ RSL_load_refl_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = dz_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"dz_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"dz_%s_%2.2d.%s", time_string,i,file_suffix);
+/*
+ sprintf(filename,"dz_%s.%s.%s", time_string,in_file,file_suffix);
*/
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
- /* ZT_INDEX */
- if(total_reflectivity) {
- if(verbose) printf("Loading refl colortable...\n");
- RSL_load_refl_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = zt_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"zt_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"zt_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
- /* DR_INDEX */
- if(differential_reflectivity) {
- scale = 0.5;
- ncbins = 21;
- width = 10;
- printf("Calling RSL_rebin, %d %d %.2f\n", width);
- RSL_rebin_volume(dr_volume, width);
- if(verbose) printf("Loading zdr colortable...\n");
- RSL_load_zdr_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = dr_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"dr_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"dr_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
-
- /* VR_INDEX */
- if(velocity) {
- if(verbose) printf("Loading vel colortable...\n");
- RSL_load_vel_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = vr_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"vr_%s_%2.2d.%s", time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"vr_%s_%2.2d.%s", time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
- /* SW_INDEX */
- if(spectral_width) {
- if(verbose) printf("Loading sw colortable...\n");
- RSL_load_sw_color_table();
- for(i=0; i<num_sweeps; i++) {
- sweep = sw_volume->sweep[i];
- if(sweep == NULL) {
- printf("sweep[%d]==NULL\n",i);
- continue;
- }
- if(make_pgm) {
- sprintf(file_suffix,"pgm");
- sprintf(filename,"sw_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_pgm(sweep, filename, xdim, ydim, maxr);
- }
- if(make_gif) {
- sprintf(file_suffix,"gif");
- sprintf(filename,"sw_%s_%2.2d.%s",
- time_string,i,file_suffix);
- printf("Creating: %s\n", filename);
- RSL_sweep_to_gif(sweep,filename,xdim, ydim, maxr);
- }
- }
- }
-
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
+ /* ZT_INDEX */
+ if(total_reflectivity) {
+ if(verbose) printf("Loading refl colortable...\n");
+ RSL_load_refl_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = zt_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"zt_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"zt_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
+ /* DR_INDEX */
+ if(differential_reflectivity) {
+ scale = 0.5;
+ ncbins = 21;
+ width = 10;
+ printf("Calling RSL_rebin, %d %d %.2f\n", width);
+ RSL_rebin_volume(dr_volume, width);
+ if(verbose) printf("Loading zdr colortable...\n");
+ RSL_load_zdr_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = dr_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"dr_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"dr_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
+
+ /* VR_INDEX */
+ if(velocity) {
+ if(verbose) printf("Loading vel colortable...\n");
+ RSL_load_vel_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = vr_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"vr_%s_%2.2d.%s", time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"vr_%s_%2.2d.%s", time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
+ /* SW_INDEX */
+ if(spectral_width) {
+ if(verbose) printf("Loading sw colortable...\n");
+ RSL_load_sw_color_table();
+ for(i=0; i<num_sweeps; i++) {
+ sweep = sw_volume->sweep[i];
+ if(sweep == NULL) {
+ printf("sweep[%d]==NULL\n",i);
+ continue;
+ }
+ if(make_pgm) {
+ sprintf(file_suffix,"pgm");
+ sprintf(filename,"sw_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, pgmdir, pathname);
+ RSL_sweep_to_pgm(sweep, pathname, xdim, ydim, maxr);
+ }
+ if(make_gif) {
+ sprintf(file_suffix,"gif");
+ sprintf(filename,"sw_%s_%2.2d.%s",
+ time_string,i,file_suffix);
+ printf("Creating: %s\n", filename);
+ make_pathname(filename, gifdir, pathname);
+ RSL_sweep_to_gif(sweep,pathname,xdim, ydim, maxr);
+ }
+ }
+ }
+
/*
Write uf file if requested
*/
- if(make_uf) {
- sprintf(file_suffix,"uf.gz");
- sprintf(filename,"%s_%s.%s",site_id, time_string,file_suffix);
- printf("Creating UF file: %s\n", filename);
- RSL_radar_to_uf_gzip(radar, filename);
- }
-
- printf("-->> FIELDS: [ ");
- if(radar->v[0] != NULL) printf("DZ ");
- if(radar->v[1] != NULL) printf("VR ");
- if(radar->v[2] != NULL) printf("SW ");
- if(radar->v[3] != NULL) printf("CZ ");
- if(radar->v[4] != NULL) printf("ZT ");
- if(radar->v[5] != NULL) printf("DR ");
- if(radar->v[6] != NULL) printf("LR ");
- if(radar->v[7] != NULL) printf("ZD ");
- if(radar->v[8] != NULL) printf("DM ");
- if(radar->v[9] != NULL) printf("RH ");
- if(radar->v[10] != NULL) printf("PH ");
- if(radar->v[11] != NULL) printf("XZ ");
- if(radar->v[12] != NULL) printf("CR ");
- if(radar->v[13] != NULL) printf("MZ ");
- if(radar->v[14] != NULL) printf("MR ");
- if(radar->v[15] != NULL) printf("ZE ");
- if(radar->v[16] != NULL) printf("VE ");
- if(radar->v[17] != NULL) printf("KD ");
- if(radar->v[18] != NULL) printf("TI ");
- if(radar->v[19] != NULL) printf("DX ");
- if(radar->v[20] != NULL) printf("CH ");
- if(radar->v[21] != NULL) printf("AH ");
- if(radar->v[22] != NULL) printf("CV ");
- if(radar->v[23] != NULL) printf("AV ");
- if(radar->v[24] != NULL) printf("SQ ");
- printf("] \n\n");
+ if(make_uf) {
+ sprintf(file_suffix,"uf.gz");
+ sprintf(filename,"%s_%s.%s",site_id, time_string,file_suffix);
+ printf("Creating UF file: %s\n", filename);
+ make_pathname(filename, ufdir, pathname);
+ RSL_radar_to_uf_gzip(radar, pathname);
+ }
+
+ printf("-->> FIELDS: [ ");
+ /* Modified to use RSL_ftype from rsl.h (#define USE_RSL_VARS) and to
+ * loop through volumes. 10/16/2009, BLK
+ */
+ for (i=0; i < MAX_RADAR_VOLUMES; i++)
+ if (radar->v[i] != NULL) printf("%s ", RSL_ftype[i]);
+ /*
+ if(radar->v[0] != NULL) printf("DZ ");
+ if(radar->v[1] != NULL) printf("VR ");
+ if(radar->v[2] != NULL) printf("SW ");
+ if(radar->v[3] != NULL) printf("CZ ");
+ if(radar->v[4] != NULL) printf("ZT ");
+ if(radar->v[5] != NULL) printf("DR ");
+ if(radar->v[6] != NULL) printf("LR ");
+ if(radar->v[7] != NULL) printf("ZD ");
+ if(radar->v[8] != NULL) printf("DM ");
+ if(radar->v[9] != NULL) printf("RH ");
+ if(radar->v[10] != NULL) printf("PH ");
+ if(radar->v[11] != NULL) printf("XZ ");
+ if(radar->v[12] != NULL) printf("CR ");
+ if(radar->v[13] != NULL) printf("MZ ");
+ if(radar->v[14] != NULL) printf("MR ");
+ if(radar->v[15] != NULL) printf("ZE ");
+ if(radar->v[16] != NULL) printf("VE ");
+ if(radar->v[17] != NULL) printf("KD ");
+ if(radar->v[18] != NULL) printf("TI ");
+ if(radar->v[19] != NULL) printf("DX ");
+ if(radar->v[20] != NULL) printf("CH ");
+ if(radar->v[21] != NULL) printf("AH ");
+ if(radar->v[22] != NULL) printf("CV ");
+ if(radar->v[23] != NULL) printf("AV ");
+ if(radar->v[24] != NULL) printf("SQ ");
+ */
+ printf("] \n\n");
/*
Wrap it up!
*/
if(verbose)
- printf("Finished!\n");
+ printf("Finished!\n");
exit (0);
} /* End of main */
#include <unistd.h>
void process_args(int argc, char **argv, char *in_file, int *verbose,
- char *site_id, char *tape_id,
- int *qc_reflectivity, int *total_reflectivity,
- int *differential_reflectivity,
- int *velocity, int *spectral_width,
- int *make_gif, int *make_pgm, int *make_bscan, int *make_uf,
- int *num_sweeps, float *dbz_offset,
- int *xdim, int *ydim, float *range,
- float *gate_size_adjustment, int *print_azim)
+ char *site_id, char *tape_id,
+ int *qc_reflectivity, int *total_reflectivity,
+ int *differential_reflectivity,
+ int *velocity, int *spectral_width,
+ int *make_gif, int *make_pgm, int *make_bscan, int *make_uf,
+ int *num_sweeps, float *dbz_offset,
+ int *xdim, int *ydim, float *range,
+ float *gate_size_adjustment, int *print_azim,
+ char *gifdir, char *pgmdir, char *ufdir)
{
- extern char *optarg;
+ extern char *optarg;
extern int optind, optopt;
- char c;
+ char c;
- while ((c = getopt(argc, argv, "vgpus:t:n:x:y:r:o:a:ADCQTWV")) != -1) {
+ while ((c = getopt(argc, argv, "vgpus:t:n:x:y:r:o:a:ADCQTWVG:P:U:")) != -1) {
- switch(c) {
+ switch(c) {
/*
RSL Verbose flag
*/
- case 'v': *verbose = TRUE; break;
+ case 'v': *verbose = TRUE; break;
/*
s: First file or call sign
*/
- case 's': strcpy(site_id, optarg); break;
- case 't': strcpy(tape_id, optarg); break;
+ case 's': strcpy(site_id, optarg); break;
+ case 't': strcpy(tape_id, optarg); break;
/*
x: x dimension
r: max range
z: zoom factor (km/pixel)
*/
- case 'x': *xdim = atoi(optarg); break;
- case 'y': *ydim = atoi(optarg); break;
- case 'r': *range = atof(optarg); break;
- case 'a': *gate_size_adjustment = atof(optarg); break;
-
+ case 'x': *xdim = atoi(optarg); break;
+ case 'y': *ydim = atoi(optarg); break;
+ case 'r': *range = atof(optarg); break;
+ case 'a': *gate_size_adjustment = atof(optarg); break;
+
/* dBZ Offset
*/
- case 'o': *dbz_offset = atof(optarg); break;
+ case 'o': *dbz_offset = atof(optarg); break;
/*
T: Total reflectivity
Q: Do qc'd reflectivity
V: Do radial velocity
W: Do spectral width
*/
- case 'Q': *qc_reflectivity = TRUE; break;
- case 'T': *total_reflectivity = TRUE; break;
- case 'V': *velocity = TRUE; break;
- case 'W': *spectral_width = TRUE; break;
- case 'A': *print_azim = TRUE; break;
- case 'D': *differential_reflectivity = TRUE; break;
+ case 'Q': *qc_reflectivity = TRUE; break;
+ case 'T': *total_reflectivity = TRUE; break;
+ case 'V': *velocity = TRUE; break;
+ case 'W': *spectral_width = TRUE; break;
+ case 'A': *print_azim = TRUE; break;
+ case 'D': *differential_reflectivity = TRUE; break;
/*
g: Make gif images
p: Make pgm images
u: Make uf files
*/
- case 'g': *make_gif = TRUE; break;
- case 'p': *make_pgm = TRUE; break;
- case 'u': *make_uf = TRUE; break;
+ case 'g': *make_gif = TRUE; break;
+ case 'p': *make_pgm = TRUE; break;
+ case 'u': *make_uf = TRUE; break;
+
+/*
+ G: gif directory
+ P: pgm directory
+ U: uf directory
+*/
+ case 'G': strcpy(gifdir, optarg); break;
+ case 'P': strcpy(pgmdir, optarg); break;
+ case 'U': strcpy(ufdir, optarg); break;
/*
num_sweeps: Number of sweeps to make images of
*/
- case 'n': *num_sweeps = atoi(optarg); break;
+ case 'n': *num_sweeps = atoi(optarg); break;
/*
Deal with bad input
*/
- case '?': fprintf(stderr, "ERROR: option -%c is undefined\n", optopt);
- goto Usage;
- case ':': fprintf(stderr, "ERROR: option -%c requires an argument\n",optopt);
- goto Usage;
- default: break;
- }
- }
+ case '?': fprintf(stderr, "ERROR: option -%c is undefined\n", optopt);
+ goto Usage;
+ case ':': fprintf(stderr, "ERROR: option -%c requires an argument\n",optopt);
+ goto Usage;
+ default: break;
+ }
+ }
/*
Must have at the least a file listed on the command lines, everything
can be defaulted.
*/
- if (argc - optind != 1) {
+ if (argc - optind != 1) {
Usage:
- fprintf(stderr,"ERROR:::\n");
- fprintf(stderr,"%s [options] input_file:",argv[0]);
- fprintf(stderr,"\n[options]: ");
- fprintf(stderr,"\n\t[-v verbose_flag?] ");
- fprintf(stderr,"\n\t[-s First file or call sign?] ");
- fprintf(stderr,"\n\t[-t Tape ID] ");
- fprintf(stderr,"\n\t[-u Make UF file]");
- fprintf(stderr,"\n\t[-g Make GIF images?]");
- fprintf(stderr,"\n\t[-p Make PGM images?]");
- fprintf(stderr,"\n\t[-x X dimension]");
- fprintf(stderr,"\n\t[-y Y dimension]");
- fprintf(stderr,"\n\t[-r max range]");
- fprintf(stderr,"\n\t[-n Number of sweeps to make images]");
- fprintf(stderr,"\n\t[-Q Do qc reflectivity]");
- fprintf(stderr,"\n\t[-T Do total reflectivity]");
- fprintf(stderr,"\n\t[-V Do velocity]");
- fprintf(stderr,"\n\t[-W Do spectral_width]");
- fprintf(stderr,"\n\t[-D Do differential reflectivity");
- fprintf(stderr,"\n\t[-o Apply dBZ offset");
- fprintf(stderr,":::\n");
- exit(-1);
- }
-
- strcpy(in_file, argv[optind]);
+ fprintf(stderr,"ERROR:::\n");
+ fprintf(stderr,"%s [options] input_file:",argv[0]);
+ fprintf(stderr,"\n[options]: ");
+ fprintf(stderr,"\n\t[-v verbose_flag?] ");
+ fprintf(stderr,"\n\t[-s First file or call sign?] ");
+ fprintf(stderr,"\n\t[-t Tape ID] ");
+ fprintf(stderr,"\n\t[-u Make UF file]");
+ fprintf(stderr,"\n\t[-g Make GIF images?]");
+ fprintf(stderr,"\n\t[-p Make PGM images?]");
+ fprintf(stderr,"\n\t[-U Directory for UF output files]");
+ fprintf(stderr,"\n\t[-G Directory for GIF output files]");
+ fprintf(stderr,"\n\t[-P Directory for PGM output files]");
+ fprintf(stderr,"\n\t[-x X dimension]");
+ fprintf(stderr,"\n\t[-y Y dimension]");
+ fprintf(stderr,"\n\t[-r max range]");
+ fprintf(stderr,"\n\t[-n Number of sweeps to make images]");
+ fprintf(stderr,"\n\t[-Q Do qc reflectivity]");
+ fprintf(stderr,"\n\t[-T Do total reflectivity]");
+ fprintf(stderr,"\n\t[-V Do velocity]");
+ fprintf(stderr,"\n\t[-W Do spectral_width]");
+ fprintf(stderr,"\n\t[-D Do differential reflectivity");
+ fprintf(stderr,"\n\t[-o Apply dBZ offset");
+ fprintf(stderr,":::\n");
+ exit(-1);
+ }
+
+ strcpy(in_file, argv[optind]);
}
#endif
#include <stdlib.h>
#include <string.h>
+#include <unistd.h>
#include "rsl.h"
usage()
{
- fprintf(stderr,"Usage: wsr_hist_uf_test infile\n");
+ fprintf(stderr,"Usage: wsr_hist_uf_test infile [-s site_id]\n");
exit(-1);
}
-process_args(int argc, char **argv, char **in_file)
+process_args(int argc, char **argv, char **in_file, char **site)
{
- if (argc == 2) *in_file = strdup(argv[1]);
+ int c;
+
+ while ((c = getopt(argc, argv, "s:")) != -1)
+ switch (c) {
+ case 's': *site = strdup(optarg); break;
+ case '?': usage(argv); break;
+ default: break;
+ }
+ if (argc - optind == 1) *in_file = strdup(argv[optind]);
else usage();
}
main(int argc, char **argv)
{
char *infile;
+ char *site = NULL;
Radar *radar;
Histogram *histogram = NULL;
- process_args(argc, argv, &infile);
+ process_args(argc, argv, &infile, &site);
RSL_radar_verbose_on();
- if ((radar = RSL_anyformat_to_radar(infile, "KMLB")) == NULL) {
+ if ((radar = RSL_anyformat_to_radar(infile, site)) == NULL) {
/* RSL_wsr88d_to_radar writes an error message to stdout. */
exit(-1);
}
--- /dev/null
+#! /bin/sh
+# mkinstalldirs --- make directory hierarchy
+# Author: Noah Friedman <friedman@prep.ai.mit.edu>
+# Created: 1993-05-16
+# Public domain
+
+# $Id: mkinstalldirs,v 1.1 1999/12/06 05:21:06 merritt Exp $
+
+errstatus=0
+
+for file
+do
+ set fnord `echo ":$file" | sed -ne 's/^:\//#/;s/^://;s/\// /g;s/^#/\//;p'`
+ shift
+
+ pathcomp=
+ for d
+ do
+ pathcomp="$pathcomp$d"
+ case "$pathcomp" in
+ -* ) pathcomp=./$pathcomp ;;
+ esac
+
+ if test ! -d "$pathcomp"; then
+ echo "mkdir $pathcomp"
+
+ mkdir "$pathcomp" || lasterr=$?
+
+ if test ! -d "$pathcomp"; then
+ errstatus=$lasterr
+ fi
+ fi
+
+ pathcomp="$pathcomp/"
+ done
+done
+
+exit $errstatus
+
+# mkinstalldirs ends here
static float (*f)(Range x);
static Range (*invf)(float x);
-FILE *file;
+extern FILE *file;
void get_extended_header_info(NSIG_Sweep **nsig_sweep, int xh_size, int iray,
int nparams,
float second;
float pw;
float bin_space;
- float prf, wave, beam_width;
+ float prf, prf2, wave, beam_width;
+ int prf_mode;
+ float prf_modes[] = {1.0, 2.0/3.0, 3.0/4.0, 4.0/5.0};
float vert_half_bw, horz_half_bw;
float rng_last_bin;
float rng_first_bin, freq;
NSIG_Sweep **nsig_sweep;
NSIG_Ray *ray_p;
int itype, ifield;
- unsigned short nsig_u2byte; /* New for 2-byte data types, Aug 2009 */
+ unsigned short nsig_2byte; /* New for 2-byte data types, Aug 2009 */
+ twob nsig_twob;
Sweep *sweep;
int msec;
float azm, elev, pitch, roll, heading, azm_rate, elev_rate,
num_rays = 0;
pw = (NSIG_I4(prod_file->rec1.prod_end.pulse_wd))/100.0; /* pulse width */
prf = NSIG_I4(prod_file->rec1.prod_end.prf); /* pulse repetition frequency */
+ prf_mode = NSIG_I2(prod_file->rec2.task_config.dsp_info.prf_mode);
+ prf2 = prf * prf_modes[prf_mode];
wave = (NSIG_I4(prod_file->rec1.prod_end.wavelen))/100.0; /* wavelength (cm) */
rsl_kdp_wavelen = wave; /* EXTERNAL (volume.c) This sets KD_F and KD_INVF
* to operate with the proper wavelength.
num_samples = NSIG_I2(prod_file->rec1.prod_end.num_samp);
sweep_rate = 3.0; /** Approximate value -- info not stored **/
azim_rate = sweep_rate*360.0/60.0;
- max_vel = wave*prf/(100.0*4.0);
+ if (prf_mode != 0)
+ {
+ float max_vel1 = wave*prf/(100.0*4.0);
+ float max_vel2 = wave*prf2/(100.0*4.0);
+
+ max_vel = (max_vel1 * max_vel2)/(max_vel1-max_vel2);
+ }
+ else
+ {
+ max_vel = wave*prf/(100.0*4.0);
+ }
+
freq = (299793000.0/wave)*1.0e-4; /** freq in MHZ **/
sqi = NSIG_I2(prod_file->rec2.task_config.calib_info.sqi)/256.0;
}
if (radar_verbose_flag)
- fprintf(stderr, "vel: %f prf: %f\n", max_vel, prf);
+ fprintf(stderr, "vel: %f prf: %f prf2: %f\n", max_vel, prf, prf2);
/** Extracting Latitude and Longitude from nsig file **/
lat = nsig_from_fourb_ang(prod_file->rec2.ingest_head.lat_rad);
ray->h.unam_rng = 299793000.0 / (2.0 * prf * 1000.0); /* km */
else
ray->h.unam_rng = 0.0;
- ray->h.fix_angle = (float)sweep->h.elev;
+ ray->h.prf2 = (int) prf2;
+ ray->h.fix_angle = (float)sweep->h.elev;
ray->h.azim_rate = azim_rate;
- ray->h.pulse_count = (float)num_samples;
+ ray->h.pulse_count = num_samples;
ray->h.pulse_width = pw;
ray->h.beam_width = beam_width;
ray->h.frequency = freq / 1000.0; /* GHz */
case NSIG_DTB_VELC2:
case NSIG_DTB_ZDR2:
case NSIG_DTB_KDP2:
- memmove(&nsig_u2byte, &ray_p->range[2*k], 2);
- nsig_u2byte = NSIG_I2(&nsig_u2byte);
- if (nsig_u2byte == 0 || nsig_u2byte == 65535)
+ memmove(nsig_twob, &ray_p->range[2*k], 2);
+ nsig_2byte = NSIG_I2(nsig_twob);
+ if (nsig_2byte == 0 || nsig_2byte == 65535)
ray_data = NSIG_NO_ECHO2;
- else ray_data = (float)(nsig_u2byte-32768)/100.;
+ else ray_data = (float)(nsig_2byte-32768)/100.;
break;
case NSIG_DTB_WID2:
- memmove(&nsig_u2byte, &ray_p->range[2*k], 2);
- nsig_u2byte = NSIG_I2(&nsig_u2byte);
- if (nsig_u2byte == 0 || nsig_u2byte == 65535)
+ memmove(nsig_twob, &ray_p->range[2*k], 2);
+ nsig_2byte = NSIG_I2(nsig_twob);
+ if (nsig_2byte == 0 || nsig_2byte == 65535)
ray_data = NSIG_NO_ECHO2;
- else ray_data = (float)nsig_u2byte/100.;
+ else ray_data = (float)nsig_2byte/100.;
break;
case NSIG_DTB_PHIDP2:
- memmove(&nsig_u2byte, &ray_p->range[2*k], 2);
- nsig_u2byte = NSIG_I2(&nsig_u2byte);
- if (nsig_u2byte == 0 || nsig_u2byte == 65535)
+ memmove(nsig_twob, &ray_p->range[2*k], 2);
+ nsig_2byte = NSIG_I2(nsig_twob);
+ if (nsig_2byte == 0 || nsig_2byte == 65535)
ray_data = NSIG_NO_ECHO;
else
- ray_data = 360.*(nsig_u2byte-1)/65534.;
+ ray_data = 360.*(nsig_2byte-1)/65534.;
break;
case NSIG_DTB_SQI2:
case NSIG_DTB_RHOHV2:
- memmove(&nsig_u2byte, &ray_p->range[2*k], 2);
- nsig_u2byte = NSIG_I2(&nsig_u2byte);
- if (nsig_u2byte == 0 || nsig_u2byte == 65535)
+ memmove(nsig_twob, &ray_p->range[2*k], 2);
+ nsig_2byte = NSIG_I2(nsig_twob);
+ if (nsig_2byte == 0 || nsig_2byte == 65535)
ray_data = NSIG_NO_ECHO2;
- else ray_data = (float)(nsig_u2byte-1)/65533.;
+ else ray_data = (float)(nsig_2byte-1)/65533.;
break;
case NSIG_DTB_HCLASS2:
- memmove(&nsig_u2byte, &ray_p->range[2*k], 2);
- nsig_u2byte = NSIG_I2(&nsig_u2byte);
- if (nsig_u2byte == 0 || nsig_u2byte == 65535)
+ memmove(nsig_twob, &ray_p->range[2*k], 2);
+ nsig_2byte = NSIG_I2(nsig_twob);
+ if (nsig_2byte == 0 || nsig_2byte == 65535)
ray_data = NSIG_NO_ECHO2;
else
- ray_data = nsig_u2byte;
+ ray_data = nsig_2byte;
}
if (ray_data == NSIG_NO_ECHO || ray_data == NSIG_NO_ECHO2)
#include "rsl.h"
extern int radar_verbose_flag;
-/* Define global variable for pruning and the functions to set or unset it.
- * Added by Bart Kelley, SSAI, August 26, 2009
- */
-int prune_radar = 1;
-
-void RSL_prune_radar_on()
-{
- prune_radar = 1;
-}
-
-void RSL_prune_radar_off()
-{
- prune_radar = 0;
-}
-
Ray *RSL_prune_ray(Ray *ray)
{
if (ray == NULL) return NULL;
if (s == NULL) return NULL;
if (s->h.nrays == 0) {
- RSL_free_sweep(s);
- return NULL;
+ RSL_free_sweep(s);
+ return NULL;
}
/*
* Squash out all dataless rays. 'j' is the index for the squashed (pruned)
* rays.
*/
for (i=0,j=0; i<s->h.nrays; i++)
- if ((s->ray[i] = RSL_prune_ray(s->ray[i])))
- s->ray[j++] = s->ray[i]; /* Keep this ray. */
+ if ((s->ray[i] = RSL_prune_ray(s->ray[i])))
+ s->ray[j++] = s->ray[i]; /* Keep this ray. */
if (j==0) {
- RSL_free_sweep(s);
- return NULL; /* All rays were pruned. */
+ RSL_free_sweep(s);
+ return NULL; /* All rays were pruned. */
}
for (i=j; i<s->h.nrays; i++) s->ray[i] = NULL;
s->h.nrays = j;
if (v == NULL) return NULL;
if (v->h.nsweeps == 0) {
- RSL_free_volume(v);
- return NULL;
+ RSL_free_volume(v);
+ return NULL;
}
/*
* Squash out all dataless sweeps. 'j' is the index for sweep containing data.
*/
for (i=0,j=0; i<v->h.nsweeps; i++)
- if ((v->sweep[i] = RSL_prune_sweep(v->sweep[i])))
- v->sweep[j++] = v->sweep[i]; /* Keep this sweep. */
+ if ((v->sweep[i] = RSL_prune_sweep(v->sweep[i])))
+ v->sweep[j++] = v->sweep[i]; /* Keep this sweep. */
if (j==0) {
- RSL_free_volume(v);
- return NULL; /* All sweeps were pruned. */
+ RSL_free_volume(v);
+ return NULL; /* All sweeps were pruned. */
}
for (i=j; i<v->h.nsweeps; i++) v->sweep[i] = NULL;
v->h.nsweeps = j;
int i;
/* Volume indexes are fixed so we just prune the substructures. */
if (radar == NULL) return NULL;
- if (prune_radar)
- for (i=0; i<radar->h.nvolumes; i++)
- radar->v[i] = RSL_prune_volume(radar->v[i]);
- else if (radar_verbose_flag) fprintf(stderr,
- "RSL_prune_radar: No pruning done. prune_radar = %d\n", prune_radar);
+ for (i=0; i<radar->h.nvolumes; i++)
+ radar->v[i] = RSL_prune_volume(radar->v[i]);
return radar;
}
r = (Radar *) calloc(1, sizeof(Radar));
r->v = (Volume **) calloc(nvolumes, sizeof(Volume *));
r->h.nvolumes = nvolumes;
+ r->h.scan_mode = PPI; /* default PPI is enum constant defined in rsl.h */
return r;
}
#include <time.h>
#include <stdlib.h>
+#define USE_RSL_VARS
#include "rsl.h"
extern int radar_verbose_flag;
/* Missing data flag : -32768 when a signed short. */
#define UF_NO_DATA 0X8000
-/* Any convensions may be observed. */
+/* Field names. Any convensions may be observed. */
/* Typically:
* DZ = Reflectivity (dBZ).
* VR = Radial Velocity.
* KD = KDP wavelength*deg/km
* TI = TIME (units unknown).
* These fields may appear in any order in the UF file.
+ * There are more fields than appear here. See rsl.h.
*/
-char *UF_field_name[] = {"DZ", "VR", "SW", "CZ", "ZT", "DR", "LR",
- "ZD", "DM", "RH", "PH", "XZ", "CD", "MZ",
- "MD", "ZE", "VE", "KD", "TI", "DX", "CH",
- "AH", "CV", "AV", "SQ"
-};
-
typedef short UF_buffer[16384]; /* Bigger than documented 4096. */
int rec_len, save_rec_len;
int nfield;
float vr_az;
+ int max_field_names;
struct tm *tm;
time_t the_time;
int degree, minute;
float second;
+ int uf_sweep_mode = 1; /* default PPI */
/* Here are the arrays for each field type. Each dimension is the number
* of fields in the radar structure. I do this because the radar organization
float x;
if (r == NULL) {
- fprintf(stderr, "radar_to_uf_fp: radar pointer NULL\n");
- return;
+ fprintf(stderr, "radar_to_uf_fp: radar pointer NULL\n");
+ return;
}
-
/* Do all the headers first time around. Then, prune OP and LU. */
q_op = q_lu = q_dh = q_fh = 1;
sweep_num = ray_num = rec_num = 0;
true_nvolumes = nvolumes = maxsweeps = nrays = 0;
+ /* PPI and RHI are enum constants defined in rsl.h */
+ if (r->h.scan_mode == PPI) uf_sweep_mode = 1;
+ else if (r->h.scan_mode == RHI) uf_sweep_mode = 3;
+
/*
* The organization of the Radar structure is by volumes, then sweeps, then
* rays, then gates. This is different from the UF file organization.
* the main controlling loop variable.
*/
for (i=0; i<nvolumes; i++) {
- volume[i] = r->v[i];
- if(volume[i]) {
- nsweeps[i] = volume[i]->h.nsweeps;
- if (nsweeps[i] > maxsweeps) maxsweeps = nsweeps[i];
- true_nvolumes++;
- }
+ volume[i] = r->v[i];
+ if(volume[i]) {
+ nsweeps[i] = volume[i]->h.nsweeps;
+ if (nsweeps[i] > maxsweeps) maxsweeps = nsweeps[i];
+ true_nvolumes++;
+ }
}
if (radar_verbose_flag) {
- fprintf(stderr,"True number of volumes for UF is %d\n", true_nvolumes);
- fprintf(stderr,"Maximum # of volumes for UF is %d\n", nvolumes);
+ fprintf(stderr,"True number of volumes for UF is %d\n", true_nvolumes);
+ fprintf(stderr,"Maximum # of volumes for UF is %d\n", nvolumes);
}
+
+ max_field_names = sizeof(RSL_ftype) / 4;
+
/*--------
* LOOP for all sweeps (typically 11 or 16 for wsr88d data.
*
*/
for (i=0; i<maxsweeps; i++) {
/* Get the array of volume and sweep pointers; one for each field type. */
- nrays = 0;
- for (k=0; k<nvolumes; k++) {
- if (volume[k]) sweep[k] = volume[k]->sweep[i];
-
- /* Check if we really can access this sweep. Paul discovered that
- * if the actual number of sweeps is less than the maximum that we
- * could be chasing a bad pointer (a NON-NULL garbage pointer).
- */
- if (i >= nsweeps[k]) sweep[k] = NULL;
-
- if (sweep[k]) if (sweep[k]->h.nrays > nrays) nrays = sweep[k]->h.nrays;
- }
-
- sweep_num++; /* I guess it will be ok to count NULL sweeps. */
- ray_num = 0;
+ nrays = 0;
+ for (k=0; k<nvolumes; k++) {
+ if (volume[k]) sweep[k] = volume[k]->sweep[i];
+
+ /* Check if we really can access this sweep. Paul discovered that
+ * if the actual number of sweeps is less than the maximum that we
+ * could be chasing a bad pointer (a NON-NULL garbage pointer).
+ */
+ if (i >= nsweeps[k]) sweep[k] = NULL;
+
+ if (sweep[k]) if (sweep[k]->h.nrays > nrays) nrays = sweep[k]->h.nrays;
+ }
+
+ sweep_num++; /* I guess it will be ok to count NULL sweeps. */
+ ray_num = 0;
if (radar_verbose_flag)
- fprintf(stderr,"Processing sweep %d for %d rays.", i, nrays);
+ fprintf(stderr,"Processing sweep %d for %d rays.", i, nrays);
if (radar_verbose_flag)
- if (little_endian()) fprintf(stderr," ... On Little endian.\n");
- else fprintf(stderr,"\n");
+ if (little_endian()) fprintf(stderr," ... On Little endian.\n");
+ else fprintf(stderr,"\n");
/* Now LOOP for all rays within this particular sweep (i).
* Get all the field types together for the ray, see ray[k], and
* fill the UF data buffer appropriately.
*/
- for (j=0; j<nrays; j++) {
- memset(uf, 0, sizeof(uf));
- nfield = 0;
- ray_num++; /* And counting, possibly, NULL rays. */
- current_fh_index = 0;
-
-
- /* Find any ray for header information. It does not matter which
- * ray, since the information for the MANDITORY, OPTIONAL, and LOCAL
- * USE headers is common to any field type ray.
- */
- ray = NULL;
- for (k=0; k<nvolumes; k++) {
- if (sweep[k])
- if (j < sweep[k]->h.nrays)
- if (sweep[k]->ray)
- if ((ray = sweep[k]->ray[j])) break;
- }
-
- /* If there is no such ray, then continue on to the next ray. */
- if (ray) {
+ for (j=0; j<nrays; j++) {
+ memset(uf, 0, sizeof(uf));
+ nfield = 0;
+ ray_num++; /* And counting, possibly, NULL rays. */
+ current_fh_index = 0;
+
+
+ /* Find any ray for header information. It does not matter which
+ * ray, since the information for the MANDITORY, OPTIONAL, and LOCAL
+ * USE headers is common to any field type ray.
+ */
+ ray = NULL;
+ for (k=0; k<nvolumes; k++) {
+ if (sweep[k])
+ if (j < sweep[k]->h.nrays)
+ if (sweep[k]->ray)
+ if ((ray = sweep[k]->ray[j])) break;
+ }
+
+ /* If there is no such ray, then continue on to the next ray. */
+ if (ray) {
/*
- fprintf(stderr,"Ray: %.4d, Time: %2.2d:%2.2d:%f %.2d/%.2d/%.4d\n", ray_num, ray->h.hour, ray->h.minute, ray->h.sec, ray->h.month, ray->h.day, ray->h.year);
+ fprintf(stderr,"Ray: %.4d, Time: %2.2d:%2.2d:%f %.2d/%.2d/%.4d\n", ray_num, ray->h.hour, ray->h.minute, ray->h.sec, ray->h.month, ray->h.day, ray->h.year);
*/
- /*
- * ---- Begining of MANDITORY HEADER BLOCK.
- */
- uf_ma = uf;
- memcpy(&uf_ma[0], "UF", 2);
- if (little_endian()) memcpy(&uf_ma[0], "FU", 2);
- uf_ma[1] = 0; /* Not known yet. */
- uf_ma[2] = 0; /* Not known yet. Really, I do. */
- uf_ma[3] = 0; /* Not known yet. */
- uf_ma[4] = 0; /* Not known yet. */
-
- uf_ma[6] = 1;
- uf_ma[7] = ray_num;
- uf_ma[8 ] = 1;
- uf_ma[9 ] = sweep_num;
- memcpy(&uf_ma[10], r->h.radar_name, 8);
- if (little_endian()) swap2(&uf_ma[10], 8/2);
- memcpy(&uf_ma[14], r->h.name, 8);
- if (little_endian()) swap2(&uf_ma[14], 8/2);
- /* Convert decimal lat/lon to d:m:s */
-
- if (ray->h.lat != 0.0) {
- degree = (int)ray->h.lat;
- minute = (int)((ray->h.lat - degree) * 60);
- second = (ray->h.lat - degree - minute/60.0) * 3600.0;
- } else {
- degree = r->h.latd;
- minute = r->h.latm;
- second = r->h.lats;
- }
- uf_ma[18] = degree;
- uf_ma[19] = minute;
- if (second > 0.0) uf_ma[20] = second*64 + 0.5;
- else uf_ma[20] = second*64 - 0.5;
-
- if (ray->h.lon != 0.0) {
- degree = (int)ray->h.lon;
- minute = (int)((ray->h.lon - degree) * 60);
- second = (ray->h.lon - degree - minute/60.0) * 3600.0;
- } else {
- degree = r->h.lond;
- minute = r->h.lonm;
- second = r->h.lons;
- }
- uf_ma[21] = degree;
- uf_ma[22] = minute;
- if (second > 0.0) uf_ma[23] = second*64 + 0.5;
- else uf_ma[23] = second*64 - 0.5;
- if (ray->h.alt != 0)
- uf_ma[24] = ray->h.alt;
- else
- uf_ma[24] = r->h.height;
-
- uf_ma[25] = ray->h.year % 100; /* By definition: not year 2000 compliant. */
- uf_ma[26] = ray->h.month;
- uf_ma[27] = ray->h.day;
- uf_ma[28] = ray->h.hour;
- uf_ma[29] = ray->h.minute;
- uf_ma[30] = ray->h.sec;
- memcpy(&uf_ma[31], "UT", 2);
- if (little_endian()) memcpy(&uf_ma[31], "TU", 2);
- if (ray->h.azimuth > 0) uf_ma[32] = ray->h.azimuth*64 + 0.5;
- else uf_ma[32] = ray->h.azimuth*64 - 0.5;
- uf_ma[33] = ray->h.elev*64 + 0.5;
- uf_ma[34] = 1; /* Sweep mode: PPI = 1 */
- if (ray->h.fix_angle != 0.)
- uf_ma[35] = ray->h.fix_angle*64.0 + 0.5;
- else uf_ma[35] = sweep[k]->h.elev*64.0 + 0.5;
- uf_ma[36] = ray->h.sweep_rate*(360.0/60.0)*64.0 + 0.5;
-
- the_time = time(NULL);
- tm = gmtime(&the_time);
-
- uf_ma[37] = tm->tm_year % 100; /* Same format as data year */
- uf_ma[38] = tm->tm_mon+1;
- uf_ma[39] = tm->tm_mday;
- memcpy(&uf_ma[40], "RSL" RSL_VERSION_STR, 8);
- if (little_endian()) swap2(&uf_ma[40], 8/2);
- uf_ma[44] = (signed short)UF_NO_DATA;
- len_ma = 45;
- uf_ma[2] = len_ma+1;
- /*
- * ---- End of MANDITORY HEADER BLOCK.
- */
-
- /* ---- Begining of OPTIONAL HEADER BLOCK. */
- len_op = 0;
- if (q_op) {
- q_op = 0; /* Only once. */
- uf_op = uf+len_ma;
- memcpy(&uf_op[0], "TRMMGVUF", 8);
- if (little_endian()) swap2(&uf_op[0], 8/2);
- uf_op[4] = (signed short)UF_NO_DATA;
- uf_op[5] = (signed short)UF_NO_DATA;
- uf_op[6] = ray->h.hour;
- uf_op[7] = ray->h.minute;
- uf_op[8] = ray->h.sec;
- memcpy(&uf_op[9], "RADAR_UF", 8);
- if (little_endian()) swap2(&uf_op[9], 8/2);
- uf_op[13] = 2;
- len_op = 14;
- }
- /* ---- End of OPTIONAL HEADER BLOCK. */
-
- /* ---- Begining of LOCAL USE HEADER BLOCK. */
- /* If we have DZ and VR, check to see if their azimuths are
- * different. If they are, we store VR azimuth in Local Use
- * Header. These differences occur with WSR-88D radars, which
- * run separate sweeps for DZ and VR at low elevations.
- */
- q_lu = 0;
- if (sweep[DZ_INDEX] && sweep[VR_INDEX]) {
- if (sweep[DZ_INDEX]->ray[j] && sweep[VR_INDEX]->ray[j]) {
- vr_az = sweep[VR_INDEX]->ray[j]->h.azimuth;
- if (sweep[DZ_INDEX]->ray[j]->h.azimuth != vr_az)
- q_lu = 1; /* Set to use Local Use Header block. */
- }
- }
- len_lu = 0;
- if (q_lu) {
- /* Store azimuth for WSR-88D VR ray in Local Use Header. */
- uf_lu = uf+len_ma+len_op;
- memcpy(&uf_lu[0], "AZ", 2);
- if (little_endian()) memcpy(&uf_lu[0], "ZA", 2);
- if (vr_az > 0) uf_lu[1] = vr_az*64 + 0.5;
- else uf_lu[1] = vr_az*64 - 0.5;
- len_lu = 2;
- }
- /* ---- End of LOCAL USE HEADER BLOCK. */
-
-
- /* Here is where we loop on each field type. We need to keep
- * track of how many FIELD HEADER and FIELD DATA sections, one
- * for each field type, we fill. The variable that tracks this
- * index into 'uf' is 'current_fh_index'. It is bumped by
- * the length of the FIELD HEADER and FIELD DATA for each field
- * type encountered. Field types expected are: Reflectivity,
- * Velocity, and Spectrum width; this is a typicial list but it
- * is not restricted to it.
- */
-
- for (k=0; k<nvolumes; k++) {
- if (sweep[k])
- if (sweep[k]->ray)
- ray = sweep[k]->ray[j];
- else
- ray = NULL;
- else ray = NULL;
-
- if (k >= sizeof(UF_field_name)) {
- ray = NULL;
- fprintf(stderr,"RSL_uf_to_radar: Unknown field type encountered.\n");
- fprintf(stderr,"The field type index in Radar exceeds the number of known UF field types.\n");
- }
-
- if (ray) {
- /* ---- Begining of DATA HEADER. */
- nfield++;
- if (q_dh) {
- len_dh = 2*true_nvolumes + 3;
- uf_dh = uf+len_ma+len_op+len_lu;
- uf_dh[0] = nfield;
- uf_dh[1] = 1;
- uf_dh[2] = nfield;
- /* 'nfield' indexes the field number.
- * 'k' indexes the particular field from the volume.
- */
- memcpy(&uf_dh[3+2*(nfield-1)], UF_field_name[k], 2);
- if (little_endian()) swap2(&uf_dh[3+2*(nfield-1)], 2/2);
- if (current_fh_index == 0) current_fh_index = len_ma+len_op+len_lu+len_dh;
- uf_dh[4+2*(nfield-1)] = current_fh_index + 1;
- }
- /* ---- End of DATA HEADER. */
-
- /* ---- Begining of FIELD HEADER. */
- if (q_fh) {
- uf_fh = uf+current_fh_index;
- uf_fh[1] = scale_factor = 100;
- uf_fh[2] = ray->h.range_bin1/1000.0;
- uf_fh[3] = ray->h.range_bin1 - (1000*uf_fh[2]);
- uf_fh[4] = ray->h.gate_size;
- uf_fh[5] = ray->h.nbins;
- uf_fh[6] = ray->h.pulse_width*(RSL_SPEED_OF_LIGHT/1.0e6);
- uf_fh[7] = sweep[k]->h.beam_width*64.0 + 0.5;
- uf_fh[8] = sweep[k]->h.beam_width*64.0 + 0.5;
- uf_fh[9] = ray->h.frequency*64.0 + 0.5; /* Bandwidth (mHz). */
- uf_fh[10] = 0; /* Horizontal polarization. */
- uf_fh[11] = ray->h.wavelength*64.0*100.0; /* m to cm. */
- uf_fh[12] = ray->h.pulse_count;
- memcpy(&uf_fh[13], " ", 2);
- uf_fh[14] = (signed short)UF_NO_DATA;
- uf_fh[15] = (signed short)UF_NO_DATA;
- if (k == DZ_INDEX || k == ZT_INDEX) {
- uf_fh[16] = volume[k]->h.calibr_const*100.0 + 0.5;
- }
- else {
- memcpy(&uf_fh[16], " ", 2);
- }
- if (ray->h.prf != 0)
- uf_fh[17] = 1.0/ray->h.prf*1000000.0; /* Pulse repetition time(msec) = 1/prf */
- else
- uf_fh[17] = (signed short)UF_NO_DATA; /* Pulse repetition time = 1/prf */
- uf_fh[18] = 16;
- if (VR_INDEX == k || VE_INDEX == k) {
- uf_fh[19] = scale_factor*ray->h.nyq_vel;
- uf_fh[20] = 1;
- len_fh = 21;
- } else {
- len_fh = 19;
- }
-
- uf_fh[0] = current_fh_index + len_fh + 1;
- /* ---- End of FIELD HEADER. */
-
- /* ---- Begining of FIELD DATA. */
- uf_data = uf+len_fh+current_fh_index;
- len_data = ray->h.nbins;
- for (m=0; m<len_data; m++) {
- x = ray->h.f(ray->range[m]);
- if (x == BADVAL || x == RFVAL || x == APFLAG || x == NOECHO)
- uf_data[m] = (signed short)UF_NO_DATA;
- else
- uf_data[m] = scale_factor * x;
- }
-
- current_fh_index += (len_fh+len_data);
- }
- }
- /* ---- End of FIELD DATA. */
- }
- /* Fill in some infomation we didn't know. Like, buffer length,
- * record number, etc.
- */
- rec_num++;
- uf_ma[1] = current_fh_index;
- uf_ma[3] = len_ma + len_op + 1;
- uf_ma[4] = len_ma + len_op + len_lu + 1;
- uf_ma[5] = rec_num;
-
- /* WRITE the UF buffer. */
- rec_len =(int)uf_ma[1]*2;
- save_rec_len = rec_len; /* We destroy 'rec_len' when making it
- big endian on a little endian machine. */
- if (little_endian()) swap_4_bytes(&rec_len);
- (void)fwrite(&rec_len, sizeof(int), 1, fp);
- if (little_endian()) swap_uf_buffer(uf);
- (void)fwrite(uf, sizeof(char), save_rec_len, fp);
- (void)fwrite(&rec_len, sizeof(int), 1, fp);
- } /* if (ray) */
- }
+ /*
+ * ---- Begining of MANDITORY HEADER BLOCK.
+ */
+ uf_ma = uf;
+ memcpy(&uf_ma[0], "UF", 2);
+ if (little_endian()) memcpy(&uf_ma[0], "FU", 2);
+ uf_ma[1] = 0; /* Not known yet. */
+ uf_ma[2] = 0; /* Not known yet. Really, I do. */
+ uf_ma[3] = 0; /* Not known yet. */
+ uf_ma[4] = 0; /* Not known yet. */
+
+ uf_ma[6] = 1;
+ uf_ma[7] = ray_num;
+ uf_ma[8 ] = 1;
+ uf_ma[9 ] = sweep_num;
+ memcpy(&uf_ma[10], r->h.radar_name, 8);
+ if (little_endian()) swap2(&uf_ma[10], 8/2);
+ memcpy(&uf_ma[14], r->h.name, 8);
+ if (little_endian()) swap2(&uf_ma[14], 8/2);
+ /* Convert decimal lat/lon to d:m:s */
+
+ if (ray->h.lat != 0.0) {
+ degree = (int)ray->h.lat;
+ minute = (int)((ray->h.lat - degree) * 60);
+ second = (ray->h.lat - degree - minute/60.0) * 3600.0;
+ } else {
+ degree = r->h.latd;
+ minute = r->h.latm;
+ second = r->h.lats;
+ }
+ uf_ma[18] = degree;
+ uf_ma[19] = minute;
+ if (second > 0.0) uf_ma[20] = second*64 + 0.5;
+ else uf_ma[20] = second*64 - 0.5;
+
+ if (ray->h.lon != 0.0) {
+ degree = (int)ray->h.lon;
+ minute = (int)((ray->h.lon - degree) * 60);
+ second = (ray->h.lon - degree - minute/60.0) * 3600.0;
+ } else {
+ degree = r->h.lond;
+ minute = r->h.lonm;
+ second = r->h.lons;
+ }
+ uf_ma[21] = degree;
+ uf_ma[22] = minute;
+ if (second > 0.0) uf_ma[23] = second*64 + 0.5;
+ else uf_ma[23] = second*64 - 0.5;
+ if (ray->h.alt != 0)
+ uf_ma[24] = ray->h.alt;
+ else
+ uf_ma[24] = r->h.height;
+
+ uf_ma[25] = ray->h.year % 100; /* By definition: not year 2000 compliant. */
+ uf_ma[26] = ray->h.month;
+ uf_ma[27] = ray->h.day;
+ uf_ma[28] = ray->h.hour;
+ uf_ma[29] = ray->h.minute;
+ uf_ma[30] = ray->h.sec;
+ memcpy(&uf_ma[31], "UT", 2);
+ if (little_endian()) memcpy(&uf_ma[31], "TU", 2);
+ if (ray->h.azimuth > 0) uf_ma[32] = ray->h.azimuth*64 + 0.5;
+ else uf_ma[32] = ray->h.azimuth*64 - 0.5;
+ uf_ma[33] = ray->h.elev*64 + 0.5;
+ uf_ma[34] = uf_sweep_mode;
+ if (ray->h.fix_angle != 0.)
+ uf_ma[35] = ray->h.fix_angle*64.0 + 0.5;
+ else uf_ma[35] = sweep[k]->h.elev*64.0 + 0.5;
+ uf_ma[36] = ray->h.sweep_rate*(360.0/60.0)*64.0 + 0.5;
+
+ the_time = time(NULL);
+ tm = gmtime(&the_time);
+
+ uf_ma[37] = tm->tm_year % 100; /* Same format as data year */
+ uf_ma[38] = tm->tm_mon+1;
+ uf_ma[39] = tm->tm_mday;
+ memcpy(&uf_ma[40], "RSL" RSL_VERSION_STR, 8);
+ if (little_endian()) swap2(&uf_ma[40], 8/2);
+ uf_ma[44] = (signed short)UF_NO_DATA;
+ len_ma = 45;
+ uf_ma[2] = len_ma+1;
+ /*
+ * ---- End of MANDITORY HEADER BLOCK.
+ */
+
+ /* ---- Begining of OPTIONAL HEADER BLOCK. */
+ len_op = 0;
+ if (q_op) {
+ q_op = 0; /* Only once. */
+ uf_op = uf+len_ma;
+ memcpy(&uf_op[0], "TRMMGVUF", 8);
+ if (little_endian()) swap2(&uf_op[0], 8/2);
+ uf_op[4] = (signed short)UF_NO_DATA;
+ uf_op[5] = (signed short)UF_NO_DATA;
+ uf_op[6] = ray->h.hour;
+ uf_op[7] = ray->h.minute;
+ uf_op[8] = ray->h.sec;
+ memcpy(&uf_op[9], "RADAR_UF", 8);
+ if (little_endian()) swap2(&uf_op[9], 8/2);
+ uf_op[13] = 2;
+ len_op = 14;
+ }
+ /* ---- End of OPTIONAL HEADER BLOCK. */
+
+ /* ---- Begining of LOCAL USE HEADER BLOCK. */
+ q_lu = 0;
+
+ /* Note: Code within "#ifdef LUHDR_VR_AZ" below is retained for testing
+ * and is not normally compiled. It was used to deal with azimuth
+ * differences between DZ and VR in WSR-88D split-cuts, now handled in
+ * ingest routine.
+ */
+/* 5/18/2010 Temporarily define LUHDR_VR_AZ until merge_split_cuts is
+ completed. */
+#define LUHDR_VR_AZ
+#ifdef LUHDR_VR_AZ
+ /* If DZ and VR azimuths are different, store VR azimuth in Local Use
+ * Header. This is done for WSR-88D split cuts.
+ */
+ if (sweep[DZ_INDEX] && sweep[VR_INDEX]) {
+ if (sweep[DZ_INDEX]->ray[j] && sweep[VR_INDEX]->ray[j]) {
+ vr_az = sweep[VR_INDEX]->ray[j]->h.azimuth;
+ if (sweep[DZ_INDEX]->ray[j]->h.azimuth != vr_az)
+ q_lu = 1; /* Set to use Local Use Header block. */
+ }
+ }
+#endif
+ len_lu = 0;
+ if (q_lu) {
+ /* Store azimuth for WSR-88D VR ray in Local Use Header. */
+ uf_lu = uf+len_ma+len_op;
+ memcpy(&uf_lu[0], "AZ", 2);
+ if (little_endian()) memcpy(&uf_lu[0], "ZA", 2);
+ if (vr_az > 0) uf_lu[1] = vr_az*64 + 0.5;
+ else uf_lu[1] = vr_az*64 - 0.5;
+ len_lu = 2;
+ }
+ /* ---- End of LOCAL USE HEADER BLOCK. */
+
+
+ /* Here is where we loop on each field type. We need to keep
+ * track of how many FIELD HEADER and FIELD DATA sections, one
+ * for each field type, we fill. The variable that tracks this
+ * index into 'uf' is 'current_fh_index'. It is bumped by
+ * the length of the FIELD HEADER and FIELD DATA for each field
+ * type encountered. Field types expected are: Reflectivity,
+ * Velocity, and Spectrum width; this is a typicial list but it
+ * is not restricted to it.
+ */
+
+ for (k=0; k<nvolumes; k++) {
+ if (sweep[k])
+ if (j < sweep[k]->h.nrays && sweep[k]->ray[j])
+ ray = sweep[k]->ray[j];
+ else
+ ray = NULL;
+ else ray = NULL;
+
+ if (ray) {
+ /* ---- Begining of DATA HEADER. */
+ nfield++;
+ if (q_dh) {
+ len_dh = 2*true_nvolumes + 3;
+ uf_dh = uf+len_ma+len_op+len_lu;
+ uf_dh[0] = nfield;
+ uf_dh[1] = 1;
+ uf_dh[2] = nfield;
+ /* 'nfield' indexes the field number.
+ * 'k' indexes the particular field from the volume.
+ * RSL_ftype contains field names and is defined in rsl.h.
+ */
+ if (k > max_field_names-1) {
+ fprintf(stderr,
+ "RSL_uf_to_radar: No field name for volume index %d\n", k);
+ fprintf(stderr,"RSL_ftype must be updated in rsl.h for new field.\n");
+ fprintf(stderr,"Quitting now.\n");
+ return;
+ }
+ memcpy(&uf_dh[3+2*(nfield-1)], RSL_ftype[k], 2);
+ if (little_endian()) swap2(&uf_dh[3+2*(nfield-1)], 2/2);
+ if (current_fh_index == 0) current_fh_index = len_ma+len_op+len_lu+len_dh;
+ uf_dh[4+2*(nfield-1)] = current_fh_index + 1;
+ }
+ /* ---- End of DATA HEADER. */
+
+ /* ---- Begining of FIELD HEADER. */
+ if (q_fh) {
+ uf_fh = uf+current_fh_index;
+ uf_fh[1] = scale_factor = 100;
+ uf_fh[2] = ray->h.range_bin1/1000.0;
+ uf_fh[3] = ray->h.range_bin1 - (1000*uf_fh[2]);
+ uf_fh[4] = ray->h.gate_size;
+ uf_fh[5] = ray->h.nbins;
+ uf_fh[6] = ray->h.pulse_width*(RSL_SPEED_OF_LIGHT/1.0e6);
+ uf_fh[7] = sweep[k]->h.beam_width*64.0 + 0.5;
+ uf_fh[8] = sweep[k]->h.beam_width*64.0 + 0.5;
+ uf_fh[9] = ray->h.frequency*64.0 + 0.5; /* Bandwidth (mHz). */
+ uf_fh[10] = 0; /* Horizontal polarization. */
+ uf_fh[11] = ray->h.wavelength*64.0*100.0; /* m to cm. */
+ uf_fh[12] = ray->h.pulse_count;
+ memcpy(&uf_fh[13], " ", 2);
+ uf_fh[14] = (signed short)UF_NO_DATA;
+ uf_fh[15] = (signed short)UF_NO_DATA;
+ if (k == DZ_INDEX || k == ZT_INDEX) {
+ uf_fh[16] = volume[k]->h.calibr_const*100.0 + 0.5;
+ }
+ else {
+ memcpy(&uf_fh[16], " ", 2);
+ }
+ if (ray->h.prf != 0)
+ uf_fh[17] = 1.0/ray->h.prf*1000000.0; /* Pulse repetition time(msec) = 1/prf */
+ else
+ uf_fh[17] = (signed short)UF_NO_DATA; /* Pulse repetition time = 1/prf */
+ uf_fh[18] = 16;
+ if (VR_INDEX == k || VE_INDEX == k) {
+ uf_fh[19] = scale_factor*ray->h.nyq_vel;
+ uf_fh[20] = 1;
+ len_fh = 21;
+ } else {
+ len_fh = 19;
+ }
+
+ uf_fh[0] = current_fh_index + len_fh + 1;
+ /* ---- End of FIELD HEADER. */
+
+ /* ---- Begining of FIELD DATA. */
+ uf_data = uf+len_fh+current_fh_index;
+ len_data = ray->h.nbins;
+ for (m=0; m<len_data; m++) {
+ x = ray->h.f(ray->range[m]);
+ if (x == BADVAL || x == RFVAL || x == APFLAG || x == NOECHO)
+ uf_data[m] = (signed short)UF_NO_DATA;
+ else
+ uf_data[m] = scale_factor * x;
+ }
+
+ current_fh_index += (len_fh+len_data);
+ }
+ }
+ /* ---- End of FIELD DATA. */
+ }
+ /* Fill in some infomation we didn't know. Like, buffer length,
+ * record number, etc.
+ */
+ rec_num++;
+ uf_ma[1] = current_fh_index;
+ uf_ma[3] = len_ma + len_op + 1;
+ uf_ma[4] = len_ma + len_op + len_lu + 1;
+ uf_ma[5] = rec_num;
+
+ /* WRITE the UF buffer. */
+ rec_len =(int)uf_ma[1]*2;
+ save_rec_len = rec_len; /* We destroy 'rec_len' when making it
+ big endian on a little endian machine. */
+ if (little_endian()) swap_4_bytes(&rec_len);
+ (void)fwrite(&rec_len, sizeof(int), 1, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+ (void)fwrite(uf, sizeof(char), save_rec_len, fp);
+ (void)fwrite(&rec_len, sizeof(int), 1, fp);
+ } /* if (ray) */
+ }
}
}
{
FILE *fp;
if (r == NULL) {
- fprintf(stderr, "radar_to_uf: radar pointer NULL\n");
- return;
+ fprintf(stderr, "radar_to_uf: radar pointer NULL\n");
+ return;
}
if ((fp = fopen(outfile, "w")) == NULL) {
- perror(outfile);
- return;
+ perror(outfile);
+ return;
}
RSL_radar_to_uf_fp(r, fp);
{
FILE *fp;
if (r == NULL) {
- fprintf(stderr, "radar_to_uf_gzip: radar pointer NULL\n");
- return;
+ fprintf(stderr, "radar_to_uf_gzip: radar pointer NULL\n");
+ return;
}
if ((fp = fopen(outfile, "w")) == NULL) {
- perror(outfile);
- return;
+ perror(outfile);
+ return;
}
fp = compress_pipe(fp);
#include "config.h"
#endif
-#define RSL_VERSION_STR "v1.40"
+#define RSL_VERSION_STR "v1.40.3"
/**********************************************************************/
/* Configure: Define USE_TWO_BYTE_PRECISION to have RSL store internal*/
int minute;/* Date for this ray; minute (0-59).*/
float sec; /* Date for this ray; second + fraction of second. */
float unam_rng; /* Unambiguous range. (KM). */
- float azimuth; /* Azimuth angle. (degrees). Must be positive
- * 0=North, 90=east, -90/270=west.
+ float azimuth; /* Azimuth angle (degrees). Must be positive.
+ * 0=North, 90=east, -90/270=west.
* This angle is the mean azimuth for the whole ray.
- * Eg. for NSIG the beginning and end azimuths are
- * averaged.
- */
+ * E.g. for NSIG the beginning and end azimuths are averaged.
+ */
int ray_num; /* Ray no. within elevation scan. */
float elev; /* Elevation angle. (degrees). */
int elev_num; /* Elevation no. within volume scan. */
float sweep_rate; /* Sweep rate. Full sweeps/min. */
int prf; /* Pulse repetition frequency, in Hz. */
+ int prf2; /* Second PRF, for Sigmet dual PRF */
float azim_rate; /* Sweep rate in degrees/sec. */
float fix_angle; /* Elevation angle for the sweep. (degrees). */
float pitch; /* Pitch angle. */
float lon; /* Longitude (degrees) */
int alt; /* Altitude (m) */
float rvc; /* Radial velocity correction (m/sec) */
- float vel_east; /* Platform velocity to the east (m/sec) */
- float vel_north; /* Platform velocity to the north (m/sec) */
- float vel_up; /* Platform velocity toward up (m/sec) */
+ float vel_east; /* Platform velocity to the east (negative for west) (m/sec) */
+ float vel_north; /* Platform velocity to the north (negative south) (m/sec) */
+ float vel_up; /* Platform velocity toward up (negative down) (m/sec) */
int pulse_count; /* Pulses used in a single dwell time. */
float pulse_width; /* Pulse width (micro-sec). */
float beam_width; /* Beamwidth in degrees. */
typedef struct {
int sweep_num; /* Integer sweep number. This may be redundant, since
* this will be the same as the Volume.sweep array index.*/
- float elev; /* Elevation angle (mean) for the sweep. */
+ float elev; /* Elevation angle (mean) for the sweep. Value is -999 for
+ * RHI. */
+ float azimuth; /* Azimuth for the sweep (RHI). Value is -999 for PPI. */
float beam_width; /* This is in the ray header too. */
float vert_half_bw; /* Vertical beam width divided by 2 */
float horz_half_bw; /* Horizontal beam width divided by 2 */
int *data;
} Histogram;
+enum scan_mode {PPI, RHI};
+
typedef struct {
int month, day, year;
int hour, minute;
int height; /* height of site in meters above sea level*/
int spulse; /* length of short pulse (ns)*/
int lpulse; /* length of long pulse (ns) */
+ int scan_mode; /* 0 = PPI, 1 = RHI */
int vcp; /* Volume Coverage Pattern (for WSR-88D only) */
} Radar_header;
* 9 = RH_INDEX = RhoHV: Horz-Vert power corr coeff
*10 = PH_INDEX = PhiDP: Differential phase angle
*11 = XZ_INDEX = X-band reflectivity.
- *12 = CR_INDEX = Corrected DR.
+ *12 = CD_INDEX = Corrected DR.
*13 = MZ_INDEX = DZ mask for 1C-51 HDF.
*14 = MR_INDEX = DR mask for 1C-51 HDF.
*15 = ZE_INDEX = Edited reflectivity.
*22 = CV_INDEX
*23 = AV_INDEX
*24 = SQ_INDEX = Signal Quality Index (Sigmet)
- */
+ *25 = VS_INDEX = Radial Velocity Combined (DORADE)
+ *26 = VL_INDEX = Radial Velocity Combined (DORADE)
+ *27 = VG_INDEX = Radial Velocity Combined (DORADE)
+ *28 = VT_INDEX = Radial Velocity Combined (DORADE)
+ *29 = NP_INDEX = Normalized Coherent Power (DORADE)
+ *30 = HC_INDEX = HydroClass (Sigmet)
+ *31 = VC_INDEX = Radial Velocity Corrected (Sigmet)
+ *32 = V2_INDEX = Radial Velocity for VCP 121 second Doppler cut.
+ *33 = S2_INDEX = Spectrum Width for VCP 121 second Doppler cut.
+ *34 = V3_INDEX = Radial Velocity for VCP 121 third Doppler cut.
+ *35 = S3_INDEX = Spectrum Width for VCP 121 third Doppler cut.
+ */
} Radar;
* VE Edited Velocity. VE_INDEX
*
* KD KDP (deg/km) Differencial Phase KD_INDEX
- * [Sigmet, Lassen]
+ * (Sigmet, Lassen)
*
* TI TIME (unknown) for MCTEX data. TI_INDEX
- * SQ SQI: Signal Quality Index. [Sigmet] SQ_INDEX
+ *
+ * SQ SQI: Signal Quality Index. (Sigmet) SQ_INDEX
* Decimal fraction from 0 to 1, where 0
* is noise, 1 is noiseless.
+ *
+ * VS Radial Velocity, Short PRT (m/s) (DORADE) VS_INDEX
+ *
+ * VL Radial Velocity, Long PRT (m/s) (DORADE) VL_INDEX
+ *
+ * VG Radial Velocity, combined (m/s) (DORADE) VG_INDEX
+ *
+ * VT Radial Velocity, combined (m/s) (DORADE) VT_INDEX
+ *
+ * NP Normalized Coherent Power. (DORADE) NP_INDEX
+ *
+ * HC HydroClass: enumerated class. (Sigmet) HC_INDEX
+ *
+ * VC Radial Velocity corrected for (Sigmet) VC_INDEX
+ * Nyquist unfolding.
*/
/*
* The number of *_INDEX must never exceed MAX_RADAR_VOLUMES.
- * Increase MAX_RADAR_VOLUMES appropriately, for new ingest formats.
+ * Increase MAX_RADAR_VOLUMES appropriately, for new ingest formats.
+ *
+ * Also, when adding new *_INDEXes, you must update the following three arrays
+ * located near the end of this file: RSL_ftype, RSL_f_list, and RSL_invf_list.
+ * You also need to modify volume.c, updating the initialization of array
+ * rsl_qfield by adding a '1' for each new volume index.
*/
-#define MAX_RADAR_VOLUMES 25
+
+#define MAX_RADAR_VOLUMES 42
#define DZ_INDEX 0
#define VR_INDEX 1
#define CV_INDEX 22
#define AV_INDEX 23
#define SQ_INDEX 24
+#define VS_INDEX 25
+#define VL_INDEX 26
+#define VG_INDEX 27
+#define VT_INDEX 28
+#define NP_INDEX 29
+#define HC_INDEX 30
+#define VC_INDEX 31
+#define V2_INDEX 32
+#define S2_INDEX 33
+#define V3_INDEX 34
+#define S3_INDEX 35
+#define CR_INDEX 36
+#define CC_INDEX 37
+#define PR_INDEX 38
+#define SD_INDEX 39
+#define ZZ_INDEX 40
+#define RD_INDEX 41
/* Prototypes for functions. */
void RSL_print_histogram(Histogram *histogram, int min_range, int max_range,
char *filename);
void RSL_print_version();
+void RSL_prune_radar_on();
+void RSL_prune_radar_off();
void RSL_radar_to_uf(Radar *r, char *outfile);
void RSL_radar_to_uf_gzip(Radar *r, char *outfile);
void RSL_radar_verbose_off(void);
float CV_F(Range x);
float AV_F(Range x);
float SQ_F(Range x);
+float VS_F(Range x);
+float VL_F(Range x);
+float VG_F(Range x);
+float VT_F(Range x);
+float NP_F(Range x);
+float HC_F(Range x);
+float VC_F(Range x);
+float SD_F(Range x);
Range DZ_INVF(float x);
Range VR_INVF(float x);
Range CV_INVF(float x);
Range AV_INVF(float x);
Range SQ_INVF(float x);
+Range VS_INVF(float x);
+Range VL_INVF(float x);
+Range VG_INVF(float x);
+Range VT_INVF(float x);
+Range NP_INVF(float x);
+Range HC_INVF(float x);
+Range VC_INVF(float x);
+Range SD_INVF(float x);
/* If you like these variables, you can use them in your application
*/
#ifdef USE_RSL_VARS
static char *RSL_ftype[] = {"DZ", "VR", "SW", "CZ", "ZT", "DR",
- "LR", "ZD", "DM", "RH", "PH", "XZ",
- "CD", "MZ", "MD", "ZE", "VE", "KD",
- "TI", "DX", "CH", "AH", "CV", "AV",
- "SQ"};
+ "LR", "ZD", "DM", "RH", "PH", "XZ",
+ "CD", "MZ", "MD", "ZE", "VE", "KD",
+ "TI", "DX", "CH", "AH", "CV", "AV",
+ "SQ", "VS", "VL", "VG", "VT", "NP",
+ "HC", "VC", "V2", "S2", "V3", "S3",
+ "CR", "CC", "PR", "SD", "ZZ", "RD"};
static float (*RSL_f_list[])(Range x) = {DZ_F, VR_F, SW_F, CZ_F, ZT_F, DR_F,
- LR_F, ZD_F, DM_F, RH_F, PH_F, XZ_F,
- CD_F, MZ_F, MD_F, ZE_F, VE_F, KD_F,
- TI_F, DX_F, CH_F, AH_F, CV_F, AV_F,
- SQ_F};
+ LR_F, ZD_F, DM_F, RH_F, PH_F, XZ_F,
+ CD_F, MZ_F, MD_F, ZE_F, VE_F, KD_F,
+ TI_F, DX_F, CH_F, AH_F, CV_F, AV_F,
+ SQ_F, VS_F, VL_F, VG_F, VT_F, NP_F,
+ HC_F, VC_F, VR_F, SW_F, VR_F, SW_F,
+ DZ_F, CZ_F, PH_F, SD_F, DZ_F, DZ_F};
static Range (*RSL_invf_list[])(float x)
- = {DZ_INVF, VR_INVF, SW_INVF, CZ_INVF, ZT_INVF, DR_INVF,
- LR_INVF, ZD_INVF, DM_INVF, RH_INVF, PH_INVF, XZ_INVF,
- CD_INVF, MZ_INVF, MD_INVF, ZE_INVF, VE_INVF, KD_INVF,
- TI_INVF, DX_INVF, CH_INVF, AH_INVF, CV_INVF, AV_INVF,
- SQ_INVF};
+ = {DZ_INVF, VR_INVF, SW_INVF, CZ_INVF, ZT_INVF, DR_INVF,
+ LR_INVF, ZD_INVF, DM_INVF, RH_INVF, PH_INVF, XZ_INVF,
+ CD_INVF, MZ_INVF, MD_INVF, ZE_INVF, VE_INVF, KD_INVF,
+ TI_INVF, DX_INVF, CH_INVF, AH_INVF, CV_INVF, AV_INVF,
+ SQ_INVF, VS_INVF, VL_INVF, VG_INVF, VT_INVF, NP_INVF,
+ HC_INVF, VC_INVF, VR_INVF, SW_INVF, VR_INVF, SW_INVF,
+ DZ_INVF, CZ_INVF, PH_INVF, SD_INVF, DZ_INVF, DZ_INVF};
#endif
/* Secret routines that are quite useful and useful to developers. */
void radar_load_date_time(Radar *radar);
#include <stdlib.h>
#include <unistd.h>
+/* This allows us to use RSL_ftype, RSL_f_list, RSL_invf_list from rsl.h. */
+#define USE_RSL_VARS
#include "rsl.h"
extern int radar_verbose_flag;
typedef short UF_buffer[16384]; /* Some UF files are bigger than 4096
- * that the UF doc's specify.
- */
+ * that the UF doc's specify.
+ */
#define UF_MORE 0
#define UF_DONE 1
int i;
if (volume == NULL) return NULL;
for (i=volume->h.nsweeps; i>0; i--)
- if (volume->sweep[i-1] != NULL) {
- volume->h.nsweeps = i;
- break;
- }
+ if (volume->sweep[i-1] != NULL) {
+ volume->h.nsweeps = i;
+ break;
+ }
return volume;
}
Radar *reset_nsweeps_in_all_volumes(Radar *radar)
int i;
if (radar == NULL) return NULL;
for (i=0; i<radar->h.nvolumes; i++)
- radar->v[i] = reset_nsweeps_in_volume(radar->v[i]);
+ radar->v[i] = reset_nsweeps_in_volume(radar->v[i]);
return radar;
}
new_volume->h = old_volume->h;
new_volume->h.nsweeps = nsweeps;
for (i=0; i<old_volume->h.nsweeps; i++)
- new_volume->sweep[i] = old_volume->sweep[i]; /* Just copy pointers. */
+ new_volume->sweep[i] = old_volume->sweep[i]; /* Just copy pointers. */
/* Free the old sweep array, not the pointers to sweeps. */
free(old_volume->sweep);
return new_volume;
short *end_addr;
end_addr = buf + n;
while (buf < end_addr) {
- swap_2_bytes(buf);
- buf++;
+ swap_2_bytes(buf);
+ buf++;
}
}
int current_fh_index;
float scale_factor;
- int nfields, isweep, ifield, iray, i, m;
+ int nfields, isweep, ifield, iray, i, j, m;
static int pulled_time_from_first_ray = 0;
+ static int need_scan_mode = 1;
char *field_type; /* For printing the field type upon error. */
- short proj_name[4];
+ short proj_name[4];
Ray *ray;
Sweep *sweep;
Radar *radar;
isweep = uf_ma[9] - 1;
if (rsl_qsweep != NULL) {
- if (isweep > rsl_qsweep_max) return UF_DONE;
- if (rsl_qsweep[isweep] == 0) return UF_MORE;
+ if (isweep > rsl_qsweep_max) return UF_DONE;
+ if (rsl_qsweep[isweep] == 0) return UF_MORE;
}
/* Sticky solution here. */
#else
if (radar == NULL) {
- radar = RSL_new_radar(MAX_RADAR_VOLUMES);
- *the_radar = radar;
- pulled_time_from_first_ray = 0;
- for (i=0; i<MAX_RADAR_VOLUMES; i++)
- if (rsl_qfield[i]) /* See RSL_select_fields in volume.c */
- radar->v[i] = RSL_new_volume(20);
+ radar = RSL_new_radar(MAX_RADAR_VOLUMES);
+ *the_radar = radar;
+ pulled_time_from_first_ray = 0;
+ for (i=0; i<MAX_RADAR_VOLUMES; i++)
+ if (rsl_qfield[i]) /* See RSL_select_fields in volume.c */
+ radar->v[i] = RSL_new_volume(20);
}
#endif
+
+ if (need_scan_mode) {
+ /* PPI and RHI are enum constants defined in rsl.h */
+ if (uf_ma[34] == 1) radar->h.scan_mode = PPI;
+ else if (uf_ma[34] == 3) radar->h.scan_mode = RHI;
+ else {
+ fprintf(stderr,"Warning: UF sweep mode = %d\n", uf_ma[34]);
+ fprintf(stderr," Expected 1 or 3 (PPI or RHI)\n");
+ fprintf(stderr," Setting radar->h.scan_mode to PPI\n");
+ radar->h.scan_mode = PPI;
+ }
+ need_scan_mode = 0;
+ }
+
/* For LITTLE ENDIAN:
* WE "UNSWAP" character strings. Because there are so few of them,
* it is easier to catch them here. The entire UF buffer is swapped prior
*/
for (i=0; i<nfields; i++) {
- if (little_endian()) swap_2_bytes(&uf_dh[3+2*i]); /* Unswap. */
- if (strncmp((char *)&uf_dh[3+2*i], "DZ", 2) == 0) ifield = DZ_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "VR", 2) == 0) ifield = VR_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "SW", 2) == 0) ifield = SW_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "CZ", 2) == 0) ifield = CZ_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "UZ", 2) == 0) ifield = ZT_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "ZT", 2) == 0) ifield = ZT_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "DR", 2) == 0) ifield = DR_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "ZD", 2) == 0) ifield = ZD_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "DM", 2) == 0) ifield = DM_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "RH", 2) == 0) ifield = RH_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "PH", 2) == 0) ifield = PH_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "XZ", 2) == 0) ifield = XZ_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "CD", 2) == 0) ifield = CD_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "MZ", 2) == 0) ifield = MZ_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "MD", 2) == 0) ifield = MD_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "ZE", 2) == 0) ifield = ZE_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "VE", 2) == 0) ifield = VE_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "KD", 2) == 0) ifield = KD_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "TI", 2) == 0) ifield = TI_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "DX", 2) == 0) ifield = DX_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "CH", 2) == 0) ifield = CH_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "AH", 2) == 0) ifield = AH_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "CV", 2) == 0) ifield = CV_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "AV", 2) == 0) ifield = AV_INDEX;
- else if (strncmp((char *)&uf_dh[3+2*i], "SQ", 2) == 0) ifield = SQ_INDEX;
- else { /* etc. DON'T know how to handle this yet. */
- field_type = (char *)&uf_dh[3+2*i];
- fprintf(stderr, "Unknown field type %c%c\n", (char)field_type[0], (char)field_type[1]);
- continue;
- }
- switch (ifield) {
- case DZ_INDEX: f = DZ_F; invf = DZ_INVF; break;
- case VR_INDEX: f = VR_F; invf = VR_INVF; break;
- case SW_INDEX: f = SW_F; invf = SW_INVF; break;
- case CZ_INDEX: f = CZ_F; invf = CZ_INVF; break;
- case ZT_INDEX: f = ZT_F; invf = ZT_INVF; break;
- case DR_INDEX: f = DR_F; invf = DR_INVF; break;
- case ZD_INDEX: f = ZD_F; invf = ZD_INVF; break;
- case DM_INDEX: f = DM_F; invf = DM_INVF; break;
- case RH_INDEX: f = RH_F; invf = RH_INVF; break;
- case PH_INDEX: f = PH_F; invf = PH_INVF; break;
- case XZ_INDEX: f = XZ_F; invf = XZ_INVF; break;
- case CD_INDEX: f = CD_F; invf = CD_INVF; break;
- case MZ_INDEX: f = MZ_F; invf = MZ_INVF; break;
- case MD_INDEX: f = MD_F; invf = MD_INVF; break;
- case ZE_INDEX: f = ZE_F; invf = ZE_INVF; break;
- case VE_INDEX: f = VE_F; invf = VE_INVF; break;
- case KD_INDEX: f = KD_F; invf = KD_INVF; break;
- case TI_INDEX: f = TI_F; invf = TI_INVF; break;
- case DX_INDEX: f = DX_F; invf = DX_INVF; break;
- case CH_INDEX: f = CH_F; invf = CH_INVF; break;
- case AH_INDEX: f = AH_F; invf = AH_INVF; break;
- case CV_INDEX: f = CV_F; invf = CV_INVF; break;
- case AV_INDEX: f = AV_F; invf = AV_INVF; break;
- case SQ_INDEX: f = SQ_F; invf = SQ_INVF; break;
- default: f = DZ_F; invf = DZ_INVF; break;
- }
-
- /* Do we place the data into this volume? */
- if (radar->v[ifield] == NULL) continue; /* Nope. */
-
- if (isweep >= radar->v[ifield]->h.nsweeps) { /* Exceeded sweep limit.
- * Allocate more sweeps.
- * Copy all previous sweeps.
- */
- if (radar_verbose_flag)
- fprintf(stderr,"Exceeded sweep allocation of %d. Adding 20 more.\n", isweep);
- new_volume = RSL_new_volume(radar->v[ifield]->h.nsweeps+20);
- new_volume = copy_sweeps_into_volume(new_volume, radar->v[ifield]);
- radar->v[ifield] = new_volume;
- }
-
- if (radar->v[ifield]->sweep[isweep] == NULL) {
- if (radar_verbose_flag)
- fprintf(stderr,"Allocating new sweep for field %d, isweep %d\n", ifield, isweep);
- radar->v[ifield]->sweep[isweep] = RSL_new_sweep(1000);
- radar->v[ifield]->sweep[isweep]->h.nrays = 0; /* Increment this for each
- * ray encountered.
- */
- radar->v[ifield]->h.f = f;
- radar->v[ifield]->h.invf = invf;
- radar->v[ifield]->sweep[isweep]->h.f = f;
- radar->v[ifield]->sweep[isweep]->h.invf = invf;
- radar->v[ifield]->sweep[isweep]->h.sweep_num = uf_ma[9];
- radar->v[ifield]->sweep[isweep]->h.elev = uf_ma[35] / 64.0;
- }
-
-
-
- current_fh_index = uf_dh[4+2*i];
- uf_fh = uf + current_fh_index - 1;
- sweep = radar->v[ifield]->sweep[isweep];
- iray = sweep->h.nrays;
- nbins = uf_fh[5];
- radar->v[ifield]->sweep[isweep]->ray[iray] = RSL_new_ray(nbins);
- ray = radar->v[ifield]->sweep[isweep]->ray[iray];
- sweep->h.nrays += 1;
-
-
- if (ray) {
- /*
- * ---- Beginning of MANDATORY HEADER BLOCK.
- */
- ray->h.ray_num = uf_ma[7];
- if (little_endian()) swap2(&uf_ma[10], 8);
- memcpy(radar->h.radar_name, &uf_ma[10], 8);
- if (little_endian()) swap2(&uf_ma[10], 8/2);
- memcpy(radar->h.name, &uf_ma[14], 8);
- if (little_endian()) swap2(&uf_ma[14], 8/2);
-
+ if (little_endian()) swap_2_bytes(&uf_dh[3+2*i]); /* Unswap. */
+ ifield = -1;
+ field_type = (char *)&uf_dh[3+2*i];
+ for (j=0; j<MAX_RADAR_VOLUMES; j++) {
+ if (strncmp(field_type, RSL_ftype[j], 2) == 0) {
+ ifield = j;
+ break;
+ }
+ }
+ if (ifield < 0) { /* DON'T know how to handle this yet. */
+ fprintf(stderr, "Unknown field type %c%c\n", (char)field_type[0],
+ (char)field_type[1]);
+ continue;
+ }
+
+ f = RSL_f_list[ifield];
+ invf = RSL_invf_list[ifield];
+
+ /* Do we place the data into this volume? */
+ if (radar->v[ifield] == NULL) continue; /* Nope. */
+
+ if (isweep >= radar->v[ifield]->h.nsweeps) { /* Exceeded sweep limit.
+ * Allocate more sweeps.
+ * Copy all previous sweeps.
+ */
+ if (radar_verbose_flag)
+ fprintf(stderr,"Exceeded sweep allocation of %d. Adding 20 more.\n", isweep);
+ new_volume = RSL_new_volume(radar->v[ifield]->h.nsweeps+20);
+ new_volume = copy_sweeps_into_volume(new_volume, radar->v[ifield]);
+ radar->v[ifield] = new_volume;
+ }
+
+ if (radar->v[ifield]->sweep[isweep] == NULL) {
+ if (radar_verbose_flag)
+ fprintf(stderr,"Allocating new sweep for field %d, isweep %d\n", ifield, isweep);
+ radar->v[ifield]->sweep[isweep] = RSL_new_sweep(1000);
+ radar->v[ifield]->sweep[isweep]->h.nrays = 0; /* Increment this for each
+ * ray encountered.
+ */
+ radar->v[ifield]->h.f = f;
+ radar->v[ifield]->h.invf = invf;
+ radar->v[ifield]->sweep[isweep]->h.f = f;
+ radar->v[ifield]->sweep[isweep]->h.invf = invf;
+ radar->v[ifield]->sweep[isweep]->h.sweep_num = uf_ma[9];
+ radar->v[ifield]->sweep[isweep]->h.elev = uf_ma[35] / 64.0;
+ }
+
+
+
+ current_fh_index = uf_dh[4+2*i];
+ uf_fh = uf + current_fh_index - 1;
+ sweep = radar->v[ifield]->sweep[isweep];
+ iray = sweep->h.nrays;
+ nbins = uf_fh[5];
+ radar->v[ifield]->sweep[isweep]->ray[iray] = RSL_new_ray(nbins);
+ ray = radar->v[ifield]->sweep[isweep]->ray[iray];
+ sweep->h.nrays += 1;
+
+
+ if (ray) {
+ /*
+ * ---- Beginning of MANDATORY HEADER BLOCK.
+ */
+ ray->h.ray_num = uf_ma[7];
+ if (little_endian()) swap2(&uf_ma[10], 8);
+ memcpy(radar->h.radar_name, &uf_ma[10], 8);
+ if (little_endian()) swap2(&uf_ma[10], 8/2);
+ memcpy(radar->h.name, &uf_ma[14], 8);
+ if (little_endian()) swap2(&uf_ma[14], 8/2);
+
/* All components of lat/lon are the same sign. If not, then
* what ever wrote the UF was in error. A simple RSL program
* can repair the damage, however, not here.
*/
- ray->h.lat = uf_ma[18] + uf_ma[19]/60.0 + uf_ma[20]/64.0/3600;
- ray->h.lon = uf_ma[21] + uf_ma[22]/60.0 + uf_ma[23]/64.0/3600;
- ray->h.alt = uf_ma[24];
- ray->h.year = uf_ma[25];
- if (ray->h.year < 1900) {
- ray->h.year += 1900;
- if (ray->h.year < 1980) ray->h.year += 100; /* Year >= 2000. */
- }
- ray->h.month = uf_ma[26];
- ray->h.day = uf_ma[27];
- ray->h.hour = uf_ma[28];
- ray->h.minute = uf_ma[29];
- ray->h.sec = uf_ma[30];
- ray->h.azimuth = uf_ma[32] / 64.0;
-
- /* If Local Use Header is present and contains azimuth, use that
- * azimuth for VR and SW. This is for WSR-88D, which runs separate
- * scans for DZ and VR/SW at the lower elevations, which means DZ
- * VR/SW and have different azimuths in the "same" ray.
- */
- len_lu = uf_ma[4] - uf_ma[3];
- if (len_lu == 2 && (ifield == VR_INDEX || ifield == SW_INDEX)) {
- if (strncmp((char *)uf_lu,"ZA",2) == 0 ||
- strncmp((char *)uf_lu,"AZ",2) == 0)
- ray->h.azimuth = uf_lu[1] / 64.0;
- }
- if (ray->h.azimuth < 0.) ray->h.azimuth += 360.; /* make it 0 to 360. */
- ray->h.elev = uf_ma[33] / 64.0;
- ray->h.elev_num = sweep->h.sweep_num;
- ray->h.fix_angle = sweep->h.elev = uf_ma[35] / 64.0;
- ray->h.azim_rate = uf_ma[36] / 64.0;
- ray->h.sweep_rate = ray->h.azim_rate * (60.0/360.0);
- missing_data = uf_ma[44];
-
- if (pulled_time_from_first_ray == 0) {
- radar->h.height = uf_ma[24];
- radar->h.latd = uf_ma[18];
- radar->h.latm = uf_ma[19];
- radar->h.lats = uf_ma[20] / 64.0;
- radar->h.lond = uf_ma[21];
- radar->h.lonm = uf_ma[22];
- radar->h.lons = uf_ma[23] / 64.0;
- radar->h.year = ray->h.year;
- radar->h.month = ray->h.month;
- radar->h.day = ray->h.day;
- radar->h.hour = ray->h.hour;
- radar->h.minute = ray->h.minute;
- radar->h.sec = ray->h.sec;
- strcpy(radar->h.radar_type, "uf");
- pulled_time_from_first_ray = 1;
- }
- /*
- * ---- End of MANDATORY HEADER BLOCK.
- */
-
- /* ---- Optional header used for MCTEX files. */
- /* If this is a MCTEX file, the first 4 words following the
- mandatory header contain the string 'MCTEX'. */
- memcpy(proj_name, (short *)(uf + uf_ma[2] - 1), 8);
- if (little_endian()) swap2(proj_name, 4);
-
-
- /* ---- Local Use Header (if present) was checked during Mandatory
- * Header processing above.
- */
-
- /* ---- Begining of FIELD HEADER. */
- uf_fh = uf+current_fh_index - 1;
- scale_factor = uf_fh[1];
- ray->h.range_bin1 = uf_fh[2] * 1000.0 + uf_fh[3];
- ray->h.gate_size = uf_fh[4];
-
- ray->h.nbins = uf_fh[5];
- ray->h.pulse_width = uf_fh[6]/(RSL_SPEED_OF_LIGHT/1.0e6);
-
- if (strncmp((char *)proj_name, "MCTEX", 5) == 0) /* MCTEX? */
- {
- /* The beamwidth values are not correct in Mctex UF files. */
- ray->h.beam_width = 1.0;
- sweep->h.beam_width = ray->h.beam_width;
- sweep->h.horz_half_bw = ray->h.beam_width/2.0;
- sweep->h.vert_half_bw = ray->h.beam_width/2.0;
- }
- else /* Not MCTEX */
- {
- ray->h.beam_width = uf_fh[7] / 64.0;
- sweep->h.beam_width = uf_fh[7] / 64.0;
- sweep->h.horz_half_bw = uf_fh[7] / 128.0; /* DFF 4/4/95 */
- sweep->h.vert_half_bw = uf_fh[8] / 128.0; /* DFF 4/4/95 */
- }
- /* fprintf (stderr, "uf_fh[7] = %d, [8] = %d\n", (int)uf_fh[7], (int)uf_fh[8]); */
- if((int)uf_fh[7] == -32768) {
- ray->h.beam_width = 1;
- sweep->h.beam_width = 1;
- sweep->h.horz_half_bw = .5;
- sweep->h.vert_half_bw = .5;
- }
-
- ray->h.frequency = uf_fh[9] / 64.0;
- ray->h.wavelength = uf_fh[11] / 64.0 / 100.0; /* cm to m. */
- ray->h.pulse_count = uf_fh[12];
- if (ifield == DZ_INDEX || ifield == ZT_INDEX) {
- radar->v[ifield]->h.calibr_const = uf_fh[16] / 100.0;
- /* uf value scaled by 100 */
- }
- else {
- radar->v[ifield]->h.calibr_const = 0.0;
- }
- if (uf_fh[17] == (short)UF_NO_DATA) x = 0;
- else x = uf_fh[17] / 1000000.0; /* PRT in seconds. */
- if (x != 0) {
- ray->h.prf = 1/x;
- ray->h.unam_rng = RSL_SPEED_OF_LIGHT / (2.0 * ray->h.prf * 1000.0);
- }
- else {
- ray->h.prf = 0.0;
- ray->h.unam_rng = 0.0;
- }
-
- if (VR_INDEX == ifield || VE_INDEX == ifield) {
- ray->h.nyq_vel = uf_fh[19] / scale_factor;
- }
-
- /* ---- End of FIELD HEADER. */
-
- ray->h.f = f;
- ray->h.invf = invf;
-
- /* ---- Begining of FIELD DATA. */
- uf_data = uf+uf_fh[0] - 1;
-
- len_data = ray->h.nbins; /* Known because of RSL_new_ray. */
- for (m=0; m<len_data; m++) {
- if (uf_data[m] == (short)UF_NO_DATA)
- ray->range[m] = invf(BADVAL); /* BADVAL */
- else {
- if(uf_data[m] == missing_data)
- ray->range[m] = invf(NOECHO); /* NOECHO */
- else
- ray->range[m] = invf((float)uf_data[m]/scale_factor);
- }
- }
- }
+ ray->h.lat = uf_ma[18] + uf_ma[19]/60.0 + uf_ma[20]/64.0/3600;
+ ray->h.lon = uf_ma[21] + uf_ma[22]/60.0 + uf_ma[23]/64.0/3600;
+ ray->h.alt = uf_ma[24];
+ ray->h.year = uf_ma[25];
+ if (ray->h.year < 1900) {
+ ray->h.year += 1900;
+ if (ray->h.year < 1980) ray->h.year += 100; /* Year >= 2000. */
+ }
+ ray->h.month = uf_ma[26];
+ ray->h.day = uf_ma[27];
+ ray->h.hour = uf_ma[28];
+ ray->h.minute = uf_ma[29];
+ ray->h.sec = uf_ma[30];
+ ray->h.azimuth = uf_ma[32] / 64.0;
+
+ /* If Local Use Header is present and contains azimuth, use that
+ * azimuth for VR and SW. This is for WSR-88D, which runs separate
+ * scans for DZ and VR/SW at the lower elevations, which means DZ
+ * VR/SW and have different azimuths in the "same" ray.
+ */
+ len_lu = uf_ma[4] - uf_ma[3];
+ if (len_lu == 2 && (ifield == VR_INDEX || ifield == SW_INDEX)) {
+ if (strncmp((char *)uf_lu,"ZA",2) == 0 ||
+ strncmp((char *)uf_lu,"AZ",2) == 0)
+ ray->h.azimuth = uf_lu[1] / 64.0;
+ }
+ if (ray->h.azimuth < 0.) ray->h.azimuth += 360.; /* make it 0 to 360. */
+ ray->h.elev = uf_ma[33] / 64.0;
+ ray->h.elev_num = sweep->h.sweep_num;
+ ray->h.fix_angle = sweep->h.elev = uf_ma[35] / 64.0;
+ ray->h.azim_rate = uf_ma[36] / 64.0;
+ ray->h.sweep_rate = ray->h.azim_rate * (60.0/360.0);
+ missing_data = uf_ma[44];
+
+ if (pulled_time_from_first_ray == 0) {
+ radar->h.height = uf_ma[24];
+ radar->h.latd = uf_ma[18];
+ radar->h.latm = uf_ma[19];
+ radar->h.lats = uf_ma[20] / 64.0;
+ radar->h.lond = uf_ma[21];
+ radar->h.lonm = uf_ma[22];
+ radar->h.lons = uf_ma[23] / 64.0;
+ radar->h.year = ray->h.year;
+ radar->h.month = ray->h.month;
+ radar->h.day = ray->h.day;
+ radar->h.hour = ray->h.hour;
+ radar->h.minute = ray->h.minute;
+ radar->h.sec = ray->h.sec;
+ strcpy(radar->h.radar_type, "uf");
+ pulled_time_from_first_ray = 1;
+ }
+ /*
+ * ---- End of MANDATORY HEADER BLOCK.
+ */
+
+ /* ---- Optional header used for MCTEX files. */
+ /* If this is a MCTEX file, the first 4 words following the
+ mandatory header contain the string 'MCTEX'. */
+ memcpy(proj_name, (short *)(uf + uf_ma[2] - 1), 8);
+ if (little_endian()) swap2(proj_name, 4);
+
+
+ /* ---- Local Use Header (if present) was checked during Mandatory
+ * Header processing above.
+ */
+
+ /* ---- Begining of FIELD HEADER. */
+ uf_fh = uf+current_fh_index - 1;
+ scale_factor = uf_fh[1];
+ ray->h.range_bin1 = uf_fh[2] * 1000.0 + uf_fh[3];
+ ray->h.gate_size = uf_fh[4];
+
+ ray->h.nbins = uf_fh[5];
+ ray->h.pulse_width = uf_fh[6]/(RSL_SPEED_OF_LIGHT/1.0e6);
+
+ if (strncmp((char *)proj_name, "MCTEX", 5) == 0) /* MCTEX? */
+ {
+ /* The beamwidth values are not correct in Mctex UF files. */
+ ray->h.beam_width = 1.0;
+ sweep->h.beam_width = ray->h.beam_width;
+ sweep->h.horz_half_bw = ray->h.beam_width/2.0;
+ sweep->h.vert_half_bw = ray->h.beam_width/2.0;
+ }
+ else /* Not MCTEX */
+ {
+ ray->h.beam_width = uf_fh[7] / 64.0;
+ sweep->h.beam_width = uf_fh[7] / 64.0;
+ sweep->h.horz_half_bw = uf_fh[7] / 128.0; /* DFF 4/4/95 */
+ sweep->h.vert_half_bw = uf_fh[8] / 128.0; /* DFF 4/4/95 */
+ }
+ /* fprintf (stderr, "uf_fh[7] = %d, [8] = %d\n", (int)uf_fh[7], (int)uf_fh[8]); */
+ if((int)uf_fh[7] == -32768) {
+ ray->h.beam_width = 1;
+ sweep->h.beam_width = 1;
+ sweep->h.horz_half_bw = .5;
+ sweep->h.vert_half_bw = .5;
+ }
+
+ ray->h.frequency = uf_fh[9] / 64.0;
+ ray->h.wavelength = uf_fh[11] / 64.0 / 100.0; /* cm to m. */
+ ray->h.pulse_count = uf_fh[12];
+ if (ifield == DZ_INDEX || ifield == ZT_INDEX) {
+ radar->v[ifield]->h.calibr_const = uf_fh[16] / 100.0;
+ /* uf value scaled by 100 */
+ }
+ else {
+ radar->v[ifield]->h.calibr_const = 0.0;
+ }
+ if (uf_fh[17] == (short)UF_NO_DATA) x = 0;
+ else x = uf_fh[17] / 1000000.0; /* PRT in seconds. */
+ if (x != 0) {
+ ray->h.prf = 1/x;
+ ray->h.unam_rng = RSL_SPEED_OF_LIGHT / (2.0 * ray->h.prf * 1000.0);
+ }
+ else {
+ ray->h.prf = 0.0;
+ ray->h.unam_rng = 0.0;
+ }
+
+ if (VR_INDEX == ifield || VE_INDEX == ifield) {
+ ray->h.nyq_vel = uf_fh[19] / scale_factor;
+ }
+
+ /* ---- End of FIELD HEADER. */
+
+ ray->h.f = f;
+ ray->h.invf = invf;
+
+ /* ---- Begining of FIELD DATA. */
+ uf_data = uf+uf_fh[0] - 1;
+
+ len_data = ray->h.nbins; /* Known because of RSL_new_ray. */
+ for (m=0; m<len_data; m++) {
+ if (uf_data[m] == (short)UF_NO_DATA)
+ ray->range[m] = invf(BADVAL); /* BADVAL */
+ else {
+ if(uf_data[m] == missing_data)
+ ray->range[m] = invf(NOECHO); /* NOECHO */
+ else
+ ray->range[m] = invf((float)uf_data[m]/scale_factor);
+ }
+ }
+ }
}
return UF_MORE;
}
addr_end = uf + sizeof(UF_buffer)/sizeof(short);
while (uf < addr_end)
- swap_2_bytes(uf++);
+ swap_2_bytes(uf++);
}
enum UF_type {NOT_UF, TRUE_UF, TWO_BYTE_UF, FOUR_BYTE_UF};
Radar *RSL_uf_to_radar_fp(FILE *fp)
{
union {
- char buf[6];
- short sword;
- int word;
+ char buf[6];
+ short sword;
+ int word;
} magic;
Radar *radar;
int nbytes;
switch (uf_type) {
case FOUR_BYTE_UF:
- if (radar_verbose_flag) fprintf(stderr,"UF file with 4 byte FORTRAN record delimeters.\n");
- /* Handle first record specially, since we needed magic information. */
- nbytes = magic.word;
- if (little_endian()) swap_4_bytes(&nbytes);
- memcpy(uf, &magic.buf[4], 2);
- (void)fread(&uf[1], sizeof(char), nbytes-2, fp);
- if (little_endian()) swap_uf_buffer(uf);
- (void)fread(&nbytes, sizeof(int), 1, fp);
- if (uf_into_radar(uf, &radar) == UF_DONE) break;
- /* Now the rest of the file. */
- while(fread(&nbytes, sizeof(int), 1, fp) > 0) {
- if (little_endian()) swap_4_bytes(&nbytes);
-
- (void)fread(uf, sizeof(char), nbytes, fp);
- if (little_endian()) swap_uf_buffer(uf);
-
- (void)fread(&nbytes, sizeof(int), 1, fp);
-
- if (uf_into_radar(uf, &radar) == UF_DONE) break;
- }
- break;
+ if (radar_verbose_flag) fprintf(stderr,"UF file with 4 byte FORTRAN record delimeters.\n");
+ /* Handle first record specially, since we needed magic information. */
+ nbytes = magic.word;
+ if (little_endian()) swap_4_bytes(&nbytes);
+ memcpy(uf, &magic.buf[4], 2);
+ (void)fread(&uf[1], sizeof(char), nbytes-2, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+ (void)fread(&nbytes, sizeof(int), 1, fp);
+ if (uf_into_radar(uf, &radar) == UF_DONE) break;
+ /* Now the rest of the file. */
+ while(fread(&nbytes, sizeof(int), 1, fp) > 0) {
+ if (little_endian()) swap_4_bytes(&nbytes);
+
+ (void)fread(uf, sizeof(char), nbytes, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+
+ (void)fread(&nbytes, sizeof(int), 1, fp);
+
+ if (uf_into_radar(uf, &radar) == UF_DONE) break;
+ }
+ break;
case TWO_BYTE_UF:
- if (radar_verbose_flag) fprintf(stderr,"UF file with 2 byte FORTRAN record delimeters.\n");
- /* Handle first record specially, since we needed magic information. */
- sbytes = magic.sword;
- if (little_endian()) swap_2_bytes(&sbytes);
- memcpy(uf, &magic.buf[2], 4);
- (void)fread(&uf[2], sizeof(char), sbytes-4, fp);
- if (little_endian()) swap_uf_buffer(uf);
- (void)fread(&sbytes, sizeof(short), 1, fp);
- uf_into_radar(uf, &radar);
- /* Now the rest of the file. */
- while(fread(&sbytes, sizeof(short), 1, fp) > 0) {
- if (little_endian()) swap_2_bytes(&sbytes);
-
- (void)fread(uf, sizeof(char), sbytes, fp);
- if (little_endian()) swap_uf_buffer(uf);
-
- (void)fread(&sbytes, sizeof(short), 1, fp);
-
- if (uf_into_radar(uf, &radar) == UF_DONE) break;
- }
- break;
+ if (radar_verbose_flag) fprintf(stderr,"UF file with 2 byte FORTRAN record delimeters.\n");
+ /* Handle first record specially, since we needed magic information. */
+ sbytes = magic.sword;
+ if (little_endian()) swap_2_bytes(&sbytes);
+ memcpy(uf, &magic.buf[2], 4);
+ (void)fread(&uf[2], sizeof(char), sbytes-4, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+ (void)fread(&sbytes, sizeof(short), 1, fp);
+ uf_into_radar(uf, &radar);
+ /* Now the rest of the file. */
+ while(fread(&sbytes, sizeof(short), 1, fp) > 0) {
+ if (little_endian()) swap_2_bytes(&sbytes);
+
+ (void)fread(uf, sizeof(char), sbytes, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+
+ (void)fread(&sbytes, sizeof(short), 1, fp);
+
+ if (uf_into_radar(uf, &radar) == UF_DONE) break;
+ }
+ break;
case TRUE_UF:
- if (radar_verbose_flag) fprintf(stderr,"UF file with no FORTRAN record delimeters. Good.\n");
- /* Handle first record specially, since we needed magic information. */
- memcpy(&sbytes, &magic.buf[2], 2); /* Record length is in word #2. */
- if (little_endian()) swap_2_bytes(&sbytes); /* # of 2 byte words. */
-
- memcpy(uf, &magic.buf[0], 6);
- (void)fread(&uf[3], sizeof(short), sbytes-3, fp);
- if (little_endian()) swap_uf_buffer(uf);
- uf_into_radar(uf, &radar);
- /* Now the rest of the file. */
- while(fread(uf, sizeof(short), 2, fp) > 0) {
- memcpy(&sbytes, &uf[1], 2); /* Record length is in word #2. */
- if (little_endian()) swap_2_bytes(&sbytes);
-
- (void)fread(&uf[2], sizeof(short), sbytes-2, fp); /* Have words 1,2. */
- if (little_endian()) swap_uf_buffer(uf);
-
- if (uf_into_radar(uf, &radar) == UF_DONE) break;
- }
- break;
-
+ if (radar_verbose_flag) fprintf(stderr,"UF file with no FORTRAN record delimeters. Good.\n");
+ /* Handle first record specially, since we needed magic information. */
+ memcpy(&sbytes, &magic.buf[2], 2); /* Record length is in word #2. */
+ if (little_endian()) swap_2_bytes(&sbytes); /* # of 2 byte words. */
+
+ memcpy(uf, &magic.buf[0], 6);
+ (void)fread(&uf[3], sizeof(short), sbytes-3, fp);
+ if (little_endian()) swap_uf_buffer(uf);
+ uf_into_radar(uf, &radar);
+ /* Now the rest of the file. */
+ while(fread(uf, sizeof(short), 2, fp) > 0) {
+ memcpy(&sbytes, &uf[1], 2); /* Record length is in word #2. */
+ if (little_endian()) swap_2_bytes(&sbytes);
+
+ (void)fread(&uf[2], sizeof(short), sbytes-2, fp); /* Have words 1,2. */
+ if (little_endian()) swap_uf_buffer(uf);
+
+ if (uf_into_radar(uf, &radar) == UF_DONE) break;
+ }
+ break;
+
case NOT_UF: return NULL; break;
}
radar = reset_nsweeps_in_all_volumes(radar);
radar = NULL;
if (infile == NULL) {
- int save_fd;
- save_fd = dup(0);
- fp = fdopen(save_fd, "r");
+ int save_fd;
+ save_fd = dup(0);
+ fp = fdopen(save_fd, "r");
} else if ((fp = fopen(infile, "r")) == NULL) {
- perror(infile);
- return radar;
+ perror(infile);
+ return radar;
}
fp = uncompress_pipe(fp); /* Transparently gunzip. */
radar = RSL_uf_to_radar_fp(fp);
rsl_pclose(fp);
-
+
return radar;
}
float DZ_F(Range x) {
if (x >= F_OFFSET) /* This test works when Range is unsigned. */
- return (((float)x-F_OFFSET)/F_FACTOR - F_DZ_RANGE_OFFSET); /* Default wsr88d. */
+ return (((float)x-F_OFFSET)/F_FACTOR - F_DZ_RANGE_OFFSET); /* Default wsr88d. */
if (x == 0) return BADVAL;
if (x == 1) return RFVAL;
if (x == 2) return APFLAG;
float VR_F(Range x) {
float val;
if (x >= F_OFFSET) { /* This test works when Range is unsigned. */
- val = (((float)x-F_OFFSET)/F_FACTOR - F_VR_OFFSET); /* Default wsr88d coding. */
- /* fprintf(stderr, "x=%d, val=%f\n", x, val); */
- return val;
+ val = (((float)x-F_OFFSET)/F_FACTOR - F_VR_OFFSET); /* Default wsr88d coding. */
+ /* fprintf(stderr, "x=%d, val=%f\n", x, val); */
+ return val;
}
if (x == 0) return BADVAL;
if (x == 1) return RFVAL;
float DR_F(Range x) { /* Differential reflectivity */
float val;
if (x >= F_OFFSET) { /* This test works when Range is unsigned. */
- val = (((float)x-F_OFFSET)/F_DR_FACTOR - F_DR_OFFSET);
- return val;
+ val = (((float)x-F_OFFSET)/F_DR_FACTOR - F_DR_OFFSET);
+ return val;
}
if (x == 0) return BADVAL;
if (x == 1) return RFVAL;
float LR_F(Range x) {/* From MCTEX */
if (x >= F_OFFSET) /* This test works when Range is unsigned. */
- return (float) (x - 250.)/6.;
+ return (float) (x - 250.)/6.;
if (x == 0) return BADVAL;
if (x == 1) return RFVAL;
if (x == 2) return APFLAG;
return (x-32768.)/100.;
}
-float NP_F(Range x) { /* Normalized Coherent Power (DORADE) */
+/* Normalized Coherent Power (DORADE) */
+float NP_F(Range x)
+{
+ if (x == 0) return BADVAL;
+ return (float)(x - 1) / 100.;
+}
+
+/* Standard Deviation (for Dual-pole QC testing.) */
+float SD_F(Range x)
+{
if (x == 0) return BADVAL;
return (float)x / 100.;
}
/* Signal Quality Index */
-float SQ_F(Range x) {
+float SQ_F(Range x)
+{
if (x == 0) return BADVAL;
return (float)(x-1) / 65533.;
}
/* KD_INVF for 1 or 2 byte data. */
Range KD_INVF(float x) {
if (x == BADVAL) return (Range)0;
-/****** Commented-out code for 1-byte Sigmet native data format.
+ return (Range)(x * 100. + 32768. + 0.5);
+/****** Old code for 1-byte Sigmet native data format commented-out:
if (x == RFVAL) return (Range)1;
if (x == APFLAG) return (Range)2;
if (x == NOECHO) return (Range)3;
if (rsl_kdp_wavelen == 0.0) return (Range)0;
if (x < 0) {
- x = 127 -
- 126 * (log((double)-x) - log((double)(0.25/rsl_kdp_wavelen))) /
- log((double)600.0) +
- 0.5;
+ x = 127 -
+ 126 * (log((double)-x) - log((double)(0.25/rsl_kdp_wavelen))) /
+ log((double)600.0) +
+ 0.5;
} else if (x > 0) {
- x = 129 +
- 126 * (log((double)x) - log((double)0.25/rsl_kdp_wavelen)) /
- log((double)600.0) +
- 0.5;
+ x = 129 +
+ 126 * (log((double)x) - log((double)0.25/rsl_kdp_wavelen)) /
+ log((double)600.0) +
+ 0.5;
} else {
- x = 128;
+ x = 128;
}
x += F_OFFSET;
******/
- return (Range)(x * 100. + 32768. + 0.5);
-
+}
+
+/* Standard Deviation (for Dual-pole QC testing.) */
+Range SD_INVF(float x)
+{
+ if (x == BADVAL) return (Range)0;
+ return (Range)(x * 100.);
}
/* Signal Quality Index */
return (Range)(x * 65533. + 1. +.5);
}
-
-Range NP_INVF(float x) /* Normalized Coherent Power (DORADE) */
+/* Normalized Coherent Power (DORADE) */
+Range NP_INVF(float x)
{
if (x == BADVAL) return (0);
- return (Range)(x * 100.);
+ return (Range)(x * 100. + 1.);
}
Range TI_INVF(float x) /* MCTEX */
#define STATIC
STATIC int RSL_max_sweeps = 0; /* Initial allocation for sweep_list.
- * RSL_new_sweep will allocate the space first
- * time around.
- */
+ * RSL_new_sweep will allocate the space first
+ * time around.
+ */
STATIC int RSL_nsweep_addr = 0; /* A count of sweeps in the table. */
STATIC Sweep_list *RSL_sweep_list = NULL;
STATIC int RSL_nextents = 0;
int i;
if (table == NULL) return;
for (i=0; i<table->nindexes; i++)
- FREE_HASH_NODE(table->indexes[i]); /* A possible linked list of Rays. */
+ FREE_HASH_NODE(table->indexes[i]); /* A possible linked list of Rays. */
free(table->indexes);
free(table);
}
int j;
/* Find where it goes, split the list and slide the tail down one. */
for (i=0; i<RSL_nsweep_addr; i++)
- if (s == RSL_sweep_list[i].s_addr) break;
+ if (s == RSL_sweep_list[i].s_addr) break;
if (i == RSL_nsweep_addr) return; /* Not found. */
/* This sweep is at 'i'. */
RSL_nsweep_addr--;
for (j=i; j<RSL_nsweep_addr; j++)
- RSL_sweep_list[j] = RSL_sweep_list[j+1];
+ RSL_sweep_list[j] = RSL_sweep_list[j+1];
RSL_sweep_list[RSL_nsweep_addr].s_addr = NULL;
RSL_sweep_list[RSL_nsweep_addr].hash = NULL;
int i,j;
if (RSL_nsweep_addr >= RSL_max_sweeps) { /* Current list is too small. */
- RSL_nextents++;
- new_list = (Sweep_list *) calloc(100*RSL_nextents, sizeof(Sweep_list));
- if (new_list == NULL) {
- perror("INSERT_SWEEP");
- exit(2);
- }
+ RSL_nextents++;
+ new_list = (Sweep_list *) calloc(100*RSL_nextents, sizeof(Sweep_list));
+ if (new_list == NULL) {
+ perror("INSERT_SWEEP");
+ exit(2);
+ }
/* Copy the old list to the new one. */
- for (i=0; i<RSL_max_sweeps; i++) new_list[i] = RSL_sweep_list[i];
- RSL_max_sweeps = 100*RSL_nextents;
- free(RSL_sweep_list);
- RSL_sweep_list = new_list;
+ for (i=0; i<RSL_max_sweeps; i++) new_list[i] = RSL_sweep_list[i];
+ RSL_max_sweeps = 100*RSL_nextents;
+ free(RSL_sweep_list);
+ RSL_sweep_list = new_list;
}
/* Find where it goes, split the list and slide the tail down one. */
for (i=0; i<RSL_nsweep_addr; i++)
- if (s < RSL_sweep_list[i].s_addr) break;
-
+ if (s < RSL_sweep_list[i].s_addr) break;
+
/* This sweep goes at 'i'. But first we must split the list. */
for (j=RSL_nsweep_addr; j>i; j--)
- RSL_sweep_list[j] = RSL_sweep_list[j-1];
+ RSL_sweep_list[j] = RSL_sweep_list[j-1];
RSL_sweep_list[i].s_addr = s;
RSL_sweep_list[i].hash = NULL;
/* Simple linear search; but this will be a binary search. */
int i;
for (i=0; i<RSL_nsweep_addr; i++)
- if (s == RSL_sweep_list[i].s_addr) return i;
+ if (s == RSL_sweep_list[i].s_addr) return i;
return -1;
}
s->ray = (Ray **) calloc(max_rays, sizeof(Ray*));
if (s->ray == NULL) perror("RSL_new_sweep, Ray*");
s->h.nrays = max_rays; /* A default setting. */
+ s->h.elev = -999.;
+ s->h.azimuth = -999.;
return s;
}
int i;
if (s == NULL) return s;
for (i=0; i<s->h.nrays; i++) {
- RSL_clear_ray(s->ray[i]);
+ RSL_clear_ray(s->ray[i]);
}
return s;
}
int i;
if (v == NULL) return v;
for (i=0; i<v->h.nsweeps; i++) {
- RSL_clear_sweep(v->sweep[i]);
+ RSL_clear_sweep(v->sweep[i]);
}
return v;
}
int i;
if (s == NULL) return;
for (i=0; i<s->h.nrays; i++) {
- RSL_free_ray(s->ray[i]);
+ RSL_free_ray(s->ray[i]);
}
if (s->ray) free(s->ray);
REMOVE_SWEEP(s); /* Remove from internal Sweep list. */
if (v == NULL) return;
for (i=0; i<v->h.nsweeps; i++)
- {
- RSL_free_sweep(v->sweep[i]);
- }
+ {
+ RSL_free_sweep(v->sweep[i]);
+ }
if (v->sweep) free(v->sweep);
free(v);
}
n_sweep->h = s->h;
for (i=0; i<s->h.nrays; i++) {
- n_sweep->ray[i] = RSL_copy_ray(s->ray[i]);
+ n_sweep->ray[i] = RSL_copy_ray(s->ray[i]);
}
return n_sweep;
}
new_vol->h = v->h;
for (i=0; i<v->h.nsweeps; i++) {
- new_vol->sweep[i] = RSL_copy_sweep(v->sweep[i]);
+ new_vol->sweep[i] = RSL_copy_sweep(v->sweep[i]);
}
return new_vol;
}
double cwise_angle_diff(float x,float y)
{
/* Returns the clockwise angle difference of x to y.
- * If x = 345 and y = 355 return 10.
- * If x = 345 and y = 335 return 350
- */
+ * If x = 345 and y = 355 return 10.
+ * If x = 345 and y = 335 return 350
+ */
double d;
d = (double)(y - x);
double ccwise_angle_diff(float x,float y)
{
/* Returns the counterclockwise angle differnce of x to y.
- * If x = 345 and y = 355 return 350.
- * If x = 345 and y = 335 return 10
- */
+ * If x = 345 and y = 355 return 350.
+ * If x = 345 and y = 335 return 10
+ */
double d;
d = (double)(x - y);
if (hash == NULL) return NULL;
/* Set low pointer to hash index with ray angle just below
- * requested angle and high pointer to just above requested
- * angle.
- */
+ * requested angle and high pointer to just above requested
+ * angle.
+ */
/* set low and high pointers to initial search locations*/
low = hash;
high = hash->ray_high;
/* Search until clockwise angle to high is less then clockwise
- * angle to low.
- */
+ * angle to low.
+ */
clow = cwise_angle_diff(ray_angle,low->ray->h.azimuth);
chigh = cwise_angle_diff(ray_angle,high->ray->h.azimuth);
cclow = ccwise_angle_diff(ray_angle,low->ray->h.azimuth);
while((chigh > clow) && (clow != 0))
- {
- if (clow < cclow)
- {
- low = low->ray_low;
- high = low->ray_high; /* Not the same low as line before ! */
- }
- else
- {
- low = low->ray_high;
- high = low->ray_high; /* Not the same low as line before ! */
- }
-
- clow = cwise_angle_diff(ray_angle,low->ray->h.azimuth);
- chigh = cwise_angle_diff(ray_angle,high->ray->h.azimuth);
- cclow = ccwise_angle_diff(ray_angle,low->ray->h.azimuth);
- }
+ {
+ if (clow < cclow)
+ {
+ low = low->ray_low;
+ high = low->ray_high; /* Not the same low as line before ! */
+ }
+ else
+ {
+ low = low->ray_high;
+ high = low->ray_high; /* Not the same low as line before ! */
+ }
+
+ clow = cwise_angle_diff(ray_angle,low->ray->h.azimuth);
+ chigh = cwise_angle_diff(ray_angle,high->ray->h.azimuth);
+ cclow = ccwise_angle_diff(ray_angle,low->ray->h.azimuth);
+ }
if(chigh <= cclow)
- {
- return high;
- }
+ {
+ return high;
+ }
else
- {
- return low;
- }
+ {
+ return low;
+ }
}
* why bother?
*/
while(table->indexes[hash] == NULL) {
- hash++;
- if(hash >= table->nindexes) hash = 0;
+ hash++;
+ if(hash >= table->nindexes) hash = 0;
}
return hash;
i = SWEEP_INDEX(s);
if (i==-1) { /* Obviously, an unregistered sweep. Most likely the
- * result of pointer assignments.
- */
- i = INSERT_SWEEP(s);
+ * result of pointer assignments.
+ */
+ i = INSERT_SWEEP(s);
}
if (RSL_sweep_list[i].hash == NULL) { /* First time. Construct the table. */
- RSL_sweep_list[i].hash = construct_sweep_hash_table(s);
+ RSL_sweep_list[i].hash = construct_sweep_hash_table(s);
}
return RSL_sweep_list[i].hash;
/*
* Return closest Ray in Sweep within limit (angle) specified
* in parameter list. Assume PPI mode.
- */
+ */
int hindex;
Hash_table *hash_table;
Azimuth_hash *closest;
closest = the_closest_hash(hash_table->indexes[hindex],ray_angle);
/* Is closest ray within limit parameter ? If
- * so return ray, else return NULL.
- */
+ * so return ray, else return NULL.
+ */
close_diff = angle_diff(ray_angle,closest->ray->h.azimuth);
if (ray == NULL) return BADVAL;
if(ray->h.gate_size == 0)
- {
- if(radar_verbose_flag)
- {
- fprintf(stderr,"RSL_get_value_from_ray: ray->h.gate_size == 0\n");
- }
- return BADVAL;
- }
+ {
+ if(radar_verbose_flag)
+ {
+ fprintf(stderr,"RSL_get_value_from_ray: ray->h.gate_size == 0\n");
+ }
+ return BADVAL;
+ }
/* range_bin1 is range to center of first bin */
bin_index = (int)(((rm - ray->h.range_bin1)/ray->h.gate_size) + 0.5);
i++;
while( i < v->h.nsweeps)
- {
- if(v->sweep[i] != NULL) break;
- i++;
- }
+ {
+ if(v->sweep[i] != NULL) break;
+ i++;
+ }
if(i >= v->h.nsweeps) return NULL;
i--;
while( i >= 0)
- {
- if(v->sweep[i] != NULL) break;
- i--;
- }
+ {
+ if(v->sweep[i] != NULL) break;
+ i--;
+ }
if(i < 0) return NULL;
smallest_ray_num = 9999999;
if (s == NULL) return r;
for (i=0; i<s->h.nrays; i++)
- if (s->ray[i]) {
- if (s->ray[i]->h.ray_num <= 1) return s->ray[i];
- if (s->ray[i]->h.ray_num < smallest_ray_num) {
- r = s->ray[i];
- smallest_ray_num = r->h.ray_num;
- }
- }
+ if (s->ray[i]) {
+ if (s->ray[i]->h.ray_num <= 1) return s->ray[i];
+ if (s->ray[i]->h.ray_num < smallest_ray_num) {
+ r = s->ray[i];
+ smallest_ray_num = r->h.ray_num;
+ }
+ }
return r;
}
int i;
if (v == NULL) return NULL;
for (i=0; i<v->h.nsweeps; i++)
- if (RSL_get_first_ray_of_sweep(v->sweep[i])) return v->sweep[i];
+ if (RSL_get_first_ray_of_sweep(v->sweep[i])) return v->sweep[i];
return NULL;
}
1, 1, 1, 1, 1,
1, 1, 1, 1, 1,
1, 1, 1, 1, 1,
+ 1, 1, 1, 1, 1,
+ 1, 1, 1, 1, 1,
1, 1
};
if (radar_verbose_flag) fprintf(stderr,"Selected fields for ingest:");
while (c_field) {
- /* CHECK EACH FIELD. This is a fancier case statement than C provides. */
- if (radar_verbose_flag) fprintf(stderr," %s", c_field);
- if (strcasecmp(c_field, "all") == 0) {
- for (i=0; i<MAX_RADAR_VOLUMES; i++) rsl_qfield[i] = 1;
- } else if (strcasecmp(c_field, "none") == 0) {
- for (i=0; i<MAX_RADAR_VOLUMES; i++) rsl_qfield[i] = 0;
- } else {
-
- for (i=0; i<MAX_RADAR_VOLUMES; i++)
- if (strcasecmp(c_field, RSL_ftype[i]) == 0) {
- rsl_qfield[i] = 1;
- break; /* Break the for loop. */
- }
-
- if (i == MAX_RADAR_VOLUMES) {
- if (radar_verbose_flag)
- fprintf(stderr, "\nRSL_select_fields: Invalid field name <<%s>> specified.\n", c_field);
- }
- }
- c_field = va_arg(ap, char *);
+ /* CHECK EACH FIELD. This is a fancier case statement than C provides. */
+ if (radar_verbose_flag) fprintf(stderr," %s", c_field);
+ if (strcasecmp(c_field, "all") == 0) {
+ for (i=0; i<MAX_RADAR_VOLUMES; i++) rsl_qfield[i] = 1;
+ } else if (strcasecmp(c_field, "none") == 0) {
+ for (i=0; i<MAX_RADAR_VOLUMES; i++) rsl_qfield[i] = 0;
+ } else {
+
+ for (i=0; i<MAX_RADAR_VOLUMES; i++)
+ if (strcasecmp(c_field, RSL_ftype[i]) == 0) {
+ rsl_qfield[i] = 1;
+ break; /* Break the for loop. */
+ }
+
+ if (i == MAX_RADAR_VOLUMES) {
+ if (radar_verbose_flag)
+ fprintf(stderr, "\nRSL_select_fields: Invalid field name <<%s>> specified.\n", c_field);
+ }
+ }
+ c_field = va_arg(ap, char *);
}
if (radar_verbose_flag) fprintf(stderr,"\n");
RSL_invf_list[0] = RSL_invf_list[0];
for (i=0; i<MAX_RADAR_VOLUMES; i++)
- if (strcasecmp(c_field, RSL_ftype[i]) == 0) {
- rsl_qfield[i] = 1;
- break; /* Break the for loop. */
- }
+ if (strcasecmp(c_field, RSL_ftype[i]) == 0) {
+ rsl_qfield[i] = 1;
+ break; /* Break the for loop. */
+ }
if (i == MAX_RADAR_VOLUMES) { /* We should never see this message for
- * properly written ingest code.
- */
- fprintf(stderr, "rsl_query_field: Invalid field name <<%s>> specified.\n", c_field);
+ * properly written ingest code.
+ */
+ fprintf(stderr, "rsl_query_field: Invalid field name <<%s>> specified.\n", c_field);
}
/* 'i' is the index. Is it set? */
/* Could be static and force use of 'rsl_query_sweep' */
int *rsl_qsweep = NULL; /* If NULL, then read all sweeps. Otherwise,
- * read what is on the list.
- */
+ * read what is on the list.
+ */
#define RSL_MAX_QSWEEP 500 /* It'll be rediculious to have more. :-) */
int rsl_qsweep_max = RSL_MAX_QSWEEP;
if (radar_verbose_flag) fprintf(stderr,"Selected sweeps for ingest:");
- for (;c_sweep; c_sweep = va_arg(ap, char *))
- {
- /* CHECK EACH FIELD. This is a fancier case statement than C provides. */
- if (radar_verbose_flag) fprintf(stderr," %s", c_sweep);
- if (strcasecmp(c_sweep, "all") == 0) {
- for (i=0; i<RSL_MAX_QSWEEP; i++) rsl_qsweep[i] = 1;
- rsl_qsweep_max = RSL_MAX_QSWEEP;
- } else if (strcasecmp(c_sweep, "none") == 0) {
- /* Commented this out to save runtime -GJW
- * rsl_qsweep[] already initialized to 0 above.
- *
- * for (i=0; i<RSL_MAX_QSWEEP; i++) rsl_qsweep[i] = 0;
- * rsl_qsweep_max = -1;
- */
- } else {
- i = sscanf(c_sweep,"%d", &isweep);
- if (i == 0) { /* No match, bad argument. */
- if (radar_verbose_flag) fprintf(stderr,"\nRSL_read_these_sweeps: bad parameter %s. Ignoring.\n", c_sweep);
- continue;
- }
-
- if (isweep < 0 || isweep > RSL_MAX_QSWEEP) {
- if (radar_verbose_flag) fprintf(stderr,"\nRSL_read_these_sweeps: parameter %s not in [0,%d). Ignoring.\n", c_sweep, RSL_MAX_QSWEEP);
- continue;
- }
-
- if (isweep > rsl_qsweep_max) rsl_qsweep_max = isweep;
- rsl_qsweep[isweep] = 1;
- }
+ for (;c_sweep; c_sweep = va_arg(ap, char *))
+ {
+ /* CHECK EACH FIELD. This is a fancier case statement than C provides. */
+ if (radar_verbose_flag) fprintf(stderr," %s", c_sweep);
+ if (strcasecmp(c_sweep, "all") == 0) {
+ for (i=0; i<RSL_MAX_QSWEEP; i++) rsl_qsweep[i] = 1;
+ rsl_qsweep_max = RSL_MAX_QSWEEP;
+ } else if (strcasecmp(c_sweep, "none") == 0) {
+ /* Commented this out to save runtime -GJW
+ * rsl_qsweep[] already initialized to 0 above.
+ *
+ * for (i=0; i<RSL_MAX_QSWEEP; i++) rsl_qsweep[i] = 0;
+ * rsl_qsweep_max = -1;
+ */
+ } else {
+ i = sscanf(c_sweep,"%d", &isweep);
+ if (i == 0) { /* No match, bad argument. */
+ if (radar_verbose_flag) fprintf(stderr,"\nRSL_read_these_sweeps: bad parameter %s. Ignoring.\n", c_sweep);
+ continue;
+ }
+
+ if (isweep < 0 || isweep > RSL_MAX_QSWEEP) {
+ if (radar_verbose_flag) fprintf(stderr,"\nRSL_read_these_sweeps: parameter %s not in [0,%d). Ignoring.\n", c_sweep, RSL_MAX_QSWEEP);
+ continue;
+ }
+
+ if (isweep > rsl_qsweep_max) rsl_qsweep_max = isweep;
+ rsl_qsweep[isweep] = 1;
+ }
}
if (radar_verbose_flag) fprintf(stderr,"\n");
if (r == NULL) return;
for (ibin=0; ibin<r->h.nbins; ibin++)
{
- val = r->h.f(r->range[ibin]);
- if ( val >= (float)NOECHO ) continue; /* Invalid value */
- r->range[ibin] = r->h.invf(val + dbz_offset);
+ val = r->h.f(r->range[ibin]);
+ if ( val >= (float)NOECHO ) continue; /* Invalid value */
+ r->range[ibin] = r->h.invf(val + dbz_offset);
}
}
int iray;
if (s == NULL) return;
for (iray=0; iray<s->h.nrays; iray++)
- RSL_add_dbz_offset_to_ray(s->ray[iray], dbz_offset);
+ RSL_add_dbz_offset_to_ray(s->ray[iray], dbz_offset);
}
/*********************************************************************/
int isweep;
if (v == NULL) return;
for (isweep=0; isweep<v->h.nsweeps; isweep++)
- RSL_add_dbz_offset_to_sweep(v->sweep[isweep], dbz_offset);
+ RSL_add_dbz_offset_to_sweep(v->sweep[isweep], dbz_offset);
}
if (ray->vol_cpat == 31) return 31;
if (ray->vol_cpat == 32) return 32;
if (ray->vol_cpat == 121) return 121;
+ if (ray->vol_cpat == 211) return 211;
+ if (ray->vol_cpat == 212) return 212;
+ if (ray->vol_cpat == 221) return 221;
return 0;
}
53833 KOHX NASHVILLE TN 36 14 50 -86 33 45 176
94703 KOKX NEW_YORK_CITY NY 40 51 56 -72 51 50 26
4106 KOTX SPOKANE WA 47 40 49 -117 37 36 727
+3948 KOUN NORMAN OK 35 14 10 -97 27 44 390
3816 KPAH PADUCAH KY 37 4 6 -88 46 19 119
4832 KPBZ PITTSBURGH PA 40 31 54 -80 13 6 361
24155 KPDT PENDLETON OR 45 41 26 -118 51 10 462
#include "rsl.h"
#include "wsr88d.h"
+#include <string.h>
/* Data descriptions in the following data structures are from the "Interface
* Control Document for the RDA/RPG", Build 10.0 Draft, WSR-88D Radar
*/
typedef struct {
- short rpg[6]; /* Unused. Not really part of message header, but is
- * inserted by RPG Communications Manager for each message.
- */
- unsigned short msg_size; /* for this segment, in halfwords */
+ short rpg[6]; /* 12 bytes inserted by RPG Communications Mgr. Ignored. */
+ unsigned short msg_size; /* Message size for this segment, in halfwords */
unsigned char channel; /* RDA Redundant Channel */
- unsigned char msg_type; /* Message Type */
- unsigned short id_seq; /* I.d. Seq = 0 to 7FFF, then roll over to 0 */
- unsigned short msg_date; /* modified Julian date from 1/1/70 */
- unsigned int msg_time; /* packet generation time in ms past midnite */
- unsigned short num_segs;
- unsigned short seg_num;
+ unsigned char msg_type; /* Message type. For example, 31 */
+ unsigned short id_seq; /* Msg seq num = 0 to 7FFF, then roll over to 0 */
+ unsigned short msg_date; /* Modified Julian date from 1/1/70 */
+ unsigned int msg_time; /* Packet generation time in ms past midnight */
+ unsigned short num_segs; /* Number of segments for this message */
+ unsigned short seg_num; /* Number of this segment */
} Wsr88d_msg_hdr;
typedef struct {
float offset;
} Data_moment_hdr;
-typedef struct {
- Data_moment_hdr data_hdr;
- /* TODO: data type will need to be changed to "void" to handle Dual Pole:
- * some data are in 2-byte words.
- */
- unsigned char *data;
-} Data_moment;
+#define MAX_RADIAL_LENGTH 14288
typedef struct {
Ray_header_m31 ray_hdr;
- short unamb_rng;
- short nyq_vel;
- /* Pointers to data */
- /* TODO: Replace data moment pointers with single pointer to radial.
- * Can simply make unsigned char array as in MSG1 version.
- * Maximum radial length is 14288 bytes.
- */
- Data_moment *ref;
- Data_moment *vel;
- Data_moment *sw;
+ float unamb_rng;
+ float nyq_vel;
+ unsigned char data[MAX_RADIAL_LENGTH];
} Wsr88d_ray_m31;
-#define MAXRAYS_M31 800
-#define MAXSWEEPS 20
-enum {START_OF_ELEV, INTERMED_RADIAL, END_OF_ELEV, BEGIN_VOS, END_VOS};
+enum radial_status {START_OF_ELEV, INTERMED_RADIAL, END_OF_ELEV, BEGIN_VOS,
+ END_VOS};
void wsr88d_swap_m31_hdr(Wsr88d_msg_hdr *msghdr)
void wsr88d_swap_m31_ray(Ray_header_m31 *wsr88d_ray)
{
- int *fullword;
+ int *data_ptr;
swap_4_bytes(&wsr88d_ray->ray_time);
swap_2_bytes(&wsr88d_ray->ray_date);
swap_2_bytes(&wsr88d_ray->radial_len);
swap_4_bytes(&wsr88d_ray->elev);
swap_2_bytes(&wsr88d_ray->data_block_count);
- fullword = (int *) &wsr88d_ray->dbptr_vol_const;
- for (; fullword <= (int *) &wsr88d_ray->dbptr_rho; fullword++)
- swap_4_bytes(fullword);
+ data_ptr = (int *) &wsr88d_ray->dbptr_vol_const;
+ for (; data_ptr <= (int *) &wsr88d_ray->dbptr_rho; data_ptr++)
+ swap_4_bytes(data_ptr);
}
-void wsr88d_swap_data_moment(Data_moment_hdr *this_field)
+void wsr88d_swap_data_hdr(Data_moment_hdr *this_field)
{
short *halfword;
halfword = (short *) &this_field->ngates;
}
-void testprt(Ray_header_m31 ray_hdr)
-{
- /* For testing: Print some values from data block header. */
- printf("\nray_time: %d\n",ray_hdr.ray_time);
- printf("ray_date: %d\n",ray_hdr.ray_date);
- printf("azm_num: %d\n",ray_hdr.azm_num);
- printf("azm: %f\n",ray_hdr.azm);
- printf("radial_len: %d\n",ray_hdr.radial_len);
- printf("elev: %f\n",ray_hdr.elev);
- printf("data block count: %d\n",ray_hdr.data_block_count);
- printf("dbptr_vol_const: %d\n",ray_hdr.dbptr_vol_const);
- printf("dbptr_elev_const: %d\n",ray_hdr.dbptr_elev_const);
- printf("dbptr_radial_const: %d\n",ray_hdr.dbptr_radial_const);
- printf("dbptr_ref: %d\n",ray_hdr.dbptr_ref);
- printf("dbptr_vel: %d\n",ray_hdr.dbptr_vel);
- printf("dbptr_sw: %d\n",ray_hdr.dbptr_sw);
-}
-
-
float wsr88d_get_angle(short bitfield)
{
short mask = 1;
}
-Data_moment *read_data_moment(unsigned char *radial, int begin_data_block)
+void get_wsr88d_unamb_and_nyq_vel(Wsr88d_ray_m31 *wsr88d_ray, float *unamb_rng,
+ float *nyq_vel)
{
- Data_moment *this_field;
- unsigned char* data;
- unsigned char* dptr;
- int hdr_size, n;
-
- /*
- printf("read_data_moment: ");
- printf("sizeof(Data_moment) = %d\n", sizeof(Data_moment));
- printf("sizeof(Data_moment_hdr) = %d\n", sizeof(Data_moment_hdr));
- */
- this_field = malloc(sizeof(Data_moment));
- hdr_size = sizeof(Data_moment_hdr);
- /* printf("hdr_size = %d\n", hdr_size); */
-
- dptr = radial + begin_data_block;
- memcpy(&this_field->data_hdr, dptr, hdr_size);
- dptr += hdr_size;
-
- if (little_endian()) wsr88d_swap_data_moment(&this_field->data_hdr);
- data = (unsigned char *) malloc(this_field->data_hdr.ngates);
- memcpy(data, dptr, this_field->data_hdr.ngates);
- this_field->data = data;
- return this_field;
+ int dbptr, found, ray_hdr_len;
+ short nyq_vel_sh, unamb_rng_sh;
+
+ found = 0;
+ ray_hdr_len = sizeof(wsr88d_ray->ray_hdr);
+ dbptr = wsr88d_ray->ray_hdr.dbptr_radial_const - ray_hdr_len;
+ if (strncmp((char *) &wsr88d_ray->data[dbptr], "RRAD", 4) == 0) found = 1;
+ else {
+ dbptr = wsr88d_ray->ray_hdr.dbptr_elev_const - ray_hdr_len;
+ if (strncmp((char *) &wsr88d_ray->data[dbptr], "RRAD", 4) == 0)
+ found = 1;
+ else {
+ dbptr = wsr88d_ray->ray_hdr.dbptr_vol_const - ray_hdr_len;
+ if (strncmp((char *) &wsr88d_ray->data[dbptr], "RRAD", 4) == 0)
+ found = 1;
+ }
+ }
+ if (found) {
+ memcpy(&unamb_rng_sh, &wsr88d_ray->data[dbptr+6], 2);
+ memcpy(&nyq_vel_sh, &wsr88d_ray->data[dbptr+16], 2);
+ if (little_endian()) {
+ swap_2_bytes(&unamb_rng_sh);
+ swap_2_bytes(&nyq_vel_sh);
+ }
+ *unamb_rng = unamb_rng_sh / 10.;
+ *nyq_vel = nyq_vel_sh / 100.;
+ } else {
+ *unamb_rng = 0.;
+ *nyq_vel = 0.;
+ }
}
Wsr88d_ray_m31 *wsr88d_ray)
{
Ray_header_m31 ray_hdr;
- int n, begin_data_block, bytes_to_read, ray_hdr_len;
- unsigned char *radial;
+ int n, bytes_to_read, ray_hdr_len;
+ float nyq_vel, unamb_rng;
/* Read ray data header block */
n = fread(&wsr88d_ray->ray_hdr, sizeof(Ray_header_m31), 1, wf->fptr);
fprintf(stderr,"read_wsr88d_ray_m31: Read failed.\n");
return 0;
}
+
+ /* Byte swap header if needed. */
if (little_endian()) wsr88d_swap_m31_ray(&wsr88d_ray->ray_hdr);
ray_hdr = wsr88d_ray->ray_hdr;
-
- /*
- if (ray_hdr.azm_num == 1) testprt(wsr88d_ray->ray_hdr);
- */
-
ray_hdr_len = sizeof(ray_hdr);
- /* Read in radial
- * Use header_size offset with data pointers
- * Copy data into structures
- */
+ /* Read data portion of radial. */
bytes_to_read = msg_size - ray_hdr_len;
- radial = (unsigned char *) malloc(bytes_to_read);
- n = fread(radial, bytes_to_read, 1, wf->fptr);
+ n = fread(wsr88d_ray->data, bytes_to_read, 1, wf->fptr);
if (n < 1) {
fprintf(stderr,"read_wsr88d_ray_m31: Read failed.\n");
return 0;
}
-
- if (ray_hdr.dbptr_radial_const != 0) {
- begin_data_block = ray_hdr.dbptr_radial_const - ray_hdr_len;
- memcpy(&wsr88d_ray->unamb_rng, &radial[begin_data_block+6], 2);
- memcpy(&wsr88d_ray->nyq_vel, &radial[begin_data_block+16], 2);
- if (little_endian()) {
- swap_2_bytes(&wsr88d_ray->unamb_rng);
- swap_2_bytes(&wsr88d_ray->nyq_vel);
- }
- }
- else {
- wsr88d_ray->unamb_rng = 0;
- wsr88d_ray->nyq_vel = 0;
- }
-
- if (ray_hdr.dbptr_ref != 0) {
- begin_data_block = ray_hdr.dbptr_ref - ray_hdr_len;
- wsr88d_ray->ref = read_data_moment(radial, begin_data_block);
- }
- if (ray_hdr.dbptr_vel != 0) {
- begin_data_block = ray_hdr.dbptr_vel - ray_hdr_len;
- wsr88d_ray->vel = read_data_moment(radial, begin_data_block);
- }
- if (ray_hdr.dbptr_sw != 0) {
- begin_data_block = ray_hdr.dbptr_sw - ray_hdr_len;
- wsr88d_ray->sw = read_data_moment(radial, begin_data_block);
- }
-
-
- /* For testing: print reflectivity data. */
- int prtdata = 0;
- { int i;
- if (prtdata) {
- for (i=0; i < wsr88d_ray->ref->data_hdr.ngates; i++) {
- if (i % 10 == 0) printf("\n");
- printf(" %d", wsr88d_ray->ref->data[i]);
- }
- printf("\n");
- }
- }
-
- free(radial);
- return 1;
-}
-
-
-void wsr88d_load_ray_data(Data_moment *data_block, Ray *ray)
-{
- Data_moment_hdr data_hdr;
- int ngates;
- int i;
- float value, scale, offset;
- unsigned char *data;
- Range (*invf)(float x);
- float (*f)(Range x);
-
- data_hdr = data_block->data_hdr;
- data = data_block->data;
-
- /*
- printf("wsr88d_load_ray_data: ");
- printf(" DataName: %s\n", data_hdr.dataname);
- */
-
- ngates = data_hdr.ngates;
- /*
- printf(" ngates = %d\n", ngates);
- printf(" scale = %f\n", data_hdr.scale);
- printf(" offset = %f\n", data_hdr.offset);
- */
- offset = data_hdr.offset;
- scale = data_hdr.scale;
-
- /* Note: data range is 2-255. 0 means signal is below threshold, and 1
- * means range folded.
+ /* We retrieve unambiguous range and Nyquist velocity here so that we don't
+ * have to do it repeatedly for each data moment later.
*/
+ get_wsr88d_unamb_and_nyq_vel(wsr88d_ray, &unamb_rng, &nyq_vel);
+ wsr88d_ray->unamb_rng = unamb_rng;
+ wsr88d_ray->nyq_vel = nyq_vel;
- /* Test print
- for (i=0; i < 10; i++) {
- value = (data[i] - offset) / scale;
- printf(" bytevalue = %d, value = %f\n", data[i], value);
- }
- */
-
- /* Reflectivity */
- if (strncmp(data_hdr.dataname, "DREF", 4) == 0) {
- /* convert data to float
- * convert float to range and put in ray.range
- */
- f = DZ_F;
- invf = DZ_INVF;
- for (i = 0; i < ngates; i++) {
- if (data[i] > 1)
- value = (data[i] - offset) / scale;
- else value = (data[i] == 0) ? BADVAL : RFVAL;
- ray->range[i] = invf(value);
- ray->h.f = f;
- ray->h.invf = invf;
- }
- }
-
- /* Velocity */
- if (strncmp(data_hdr.dataname, "DVEL", 4) == 0) {
- /* convert data to float
- * convert float to range and put in ray.range
- */
- f = VR_F;
- invf = VR_INVF;
- for (i = 0; i < ngates; i++) {
- if (data[i] > 1)
- value = (data[i] - offset) / scale;
- else value = (data[i] == 0) ? BADVAL : RFVAL;
- ray->range[i] = invf(value);
- ray->h.f = f;
- ray->h.invf = invf;
- }
- }
-
- /* Spectrum Width */
- if (strncmp(data_hdr.dataname, "DSW", 3) == 0) {
- /* convert data to float
- * convert float to range and put in ray.range
- */
- f = SW_F;
- invf = SW_INVF;
- for (i = 0; i < ngates; i++) {
- if (data[i] > 1)
- value = (data[i] - offset) / scale;
- else value = (data[i] == 0) ? BADVAL : RFVAL;
- ray->range[i] = invf(value);
- ray->h.f = f;
- ray->h.invf = invf;
- }
- }
- ray->h.range_bin1 = data_hdr.range_first_gate;
- ray->h.gate_size = data_hdr.range_samp_interval;
- ray->h.nbins = ngates;
+ return 1;
}
-
void wsr88d_load_ray_hdr(Wsr88d_ray_m31 wsr88d_ray, Ray *ray)
{
int month, day, year, hour, minute, sec;
ray->h.ray_num = ray_hdr.azm_num;
ray->h.elev = ray_hdr.elev;
ray->h.elev_num = ray_hdr.elev_num;
- ray->h.unam_rng = wsr88d_ray.unamb_rng / 10.;
- ray->h.nyq_vel = wsr88d_ray.nyq_vel / 100.;
+ ray->h.unam_rng = wsr88d_ray.unamb_rng;
+ ray->h.nyq_vel = wsr88d_ray.nyq_vel;
int elev_index;
elev_index = ray_hdr.elev_num - 1;
ray->h.azim_rate = vcp_data.azim_rate[elev_index];
ray->h.fix_angle = vcp_data.fixed_angle[elev_index];
ray->h.vel_res = vcp_data.vel_res;
+ if (ray_hdr.azm_res != 1)
+ ray->h.beam_width = 1.0;
+ else ray->h.beam_width = 0.5;
- /* Get some values using message type 1 routines.
- * First load VCP and elevation numbers into msg 1 ray.
+ /* For convenience, use message type 1 routines to get some values.
+ * First load VCP and elevation numbers into a msg 1 ray.
*/
m1_ray.vol_cpat = vcp_data.vcp;
m1_ray.elev_num = ray_hdr.elev_num;
- m1_ray.unam_rng = wsr88d_ray.unamb_rng;
- if (ray_hdr.azm_res != 1)
- ray->h.beam_width = 1.0;
- else ray->h.beam_width = 0.5;
+ m1_ray.unam_rng = (short) (wsr88d_ray.unamb_rng * 10.);
+ /* Get values from message type 1 routines. */
ray->h.frequency = wsr88d_get_frequency(&m1_ray);
ray->h.pulse_width = wsr88d_get_pulse_width(&m1_ray);
ray->h.pulse_count = wsr88d_get_pulse_count(&m1_ray);
- ray->h.prf = wsr88d_get_prf(&m1_ray);
+ ray->h.prf = (int) wsr88d_get_prf(&m1_ray);
ray->h.wavelength = 0.1071;
}
int wsr88d_get_vol_index(char* dataname)
{
- int vol_index;
+ int vol_index = -1;
if (strncmp(dataname, "DREF", 4) == 0) vol_index = DZ_INDEX;
if (strncmp(dataname, "DVEL", 4) == 0) vol_index = VR_INDEX;
- if (strncmp(dataname, "DSW", 3) == 0) vol_index = SW_INDEX;
- /* TODO: Add the other data moments. */
+ if (strncmp(dataname, "DSW", 3) == 0) vol_index = SW_INDEX;
+ if (strncmp(dataname, "DZDR", 3) == 0) vol_index = DR_INDEX;
+ if (strncmp(dataname, "DPHI", 3) == 0) vol_index = PH_INDEX;
+ if (strncmp(dataname, "DRHO", 3) == 0) vol_index = RH_INDEX;
return vol_index;
}
-void wsr88d_load_ray_into_radar(Wsr88d_ray_m31 wsr88d_ray, int isweep, int iray,
- Radar *radar)
+#define MAXRAYS_M31 800
+#define MAXSWEEPS 20
+
+void wsr88d_load_ray(Wsr88d_ray_m31 wsr88d_ray, int data_ptr,
+ int isweep, int iray, Radar *radar)
{
- Ray *ray;
+ /* Load data into ray structure for this field or data moment. */
+
+ Data_moment_hdr data_hdr;
+ int ngates;
+ int i, hdr_size;
+ float value, scale, offset;
+ unsigned char *data;
Range (*invf)(float x);
float (*f)(Range x);
+ Ray *ray;
int vol_index, waveform;
- enum waveforms {surveillance=1, doppler_ambres, doppler_no_ambres, batch};
+ enum waveforms {surveillance=1, doppler_w_amb_res, doppler_no_amb_res,
+ batch};
-/* for each data moment do
- get the name of the data moment
- get new volume if needed
- get new sweep if needed
- get new ray
- put this data moment in ray
- endfor each data moment
-*/
-
- /* Note: The data moment type can only be determined by the name within the
- * individual data moment block. The name of the pointer is not reliable.
- * For example, in legacy resolution mode, dbptr_ref points to the velocity
- * data moment in velocity split cuts. Apparently the location of the first
- * data moment is always stored in the reflectivity pointer (dbptr_ref), and
- * sometimes this is velocity. When this occurs, the velocity data pointer
- * (dbptr_vel) then points to spectrum width.
- * With super resolution, there actually is reflectivity in the velocity
- * split cut. It's there for use with the recombination algorithm.
- */
+ /* Get data moment header. */
+ hdr_size = sizeof(data_hdr);
+ memcpy(&data_hdr, &wsr88d_ray.data[data_ptr], hdr_size);
+ data_ptr += hdr_size;
+ if (little_endian()) wsr88d_swap_data_hdr(&data_hdr);
- if (wsr88d_ray.ray_hdr.dbptr_ref > 0) {
- vol_index = wsr88d_get_vol_index(wsr88d_ray.ref->data_hdr.dataname);
- switch (vol_index) {
- case DZ_INDEX: f = DZ_F; invf = DZ_INVF; break;
- case VR_INDEX: f = VR_F; invf = VR_INVF; break;
- case SW_INDEX: f = SW_F; invf = SW_INVF; break;
- default: f = DZ_F; invf = DZ_INVF; break;
- }
- /* If this is reflectivity, check the waveform type to make sure
- * it isn't from a Doppler split cut.
- * We only keep reflectivity if the waveform type is surveillance or
- * batch, or the elevation is above the split cut elevations.
- */
- waveform = vcp_data.waveform[isweep];
- if (vol_index != DZ_INDEX ||
- (waveform == surveillance || waveform == batch ||
- vcp_data.fixed_angle[isweep] >= 6.0))
- {
- if (radar->v[vol_index] == NULL) {
- radar->v[vol_index] = RSL_new_volume(MAXSWEEPS);
- radar->v[vol_index]->h.f = f;
- radar->v[vol_index]->h.invf = invf;
- }
- if (radar->v[vol_index]->sweep[isweep] == NULL) {
- radar->v[vol_index]->sweep[isweep] = RSL_new_sweep(MAXRAYS_M31);
- radar->v[vol_index]->sweep[isweep]->h.f = f;
- radar->v[vol_index]->sweep[isweep]->h.invf = invf;
- }
- ray = RSL_new_ray(wsr88d_ray.ref->data_hdr.ngates);
- wsr88d_load_ray_data(wsr88d_ray.ref, ray);
- wsr88d_load_ray_hdr(wsr88d_ray, ray);
- radar->v[vol_index]->sweep[isweep]->ray[iray] = ray;
- radar->v[vol_index]->sweep[isweep]->h.nrays = iray+1;
- }
+ vol_index = wsr88d_get_vol_index(data_hdr.dataname);
+ if (vol_index < 0) {
+ fprintf(stderr,"wsr88d_load_ray: Unknown dataname %s. isweep = %d, "
+ "iray = %d.\n", data_hdr.dataname, isweep, iray);
+ return;
+ }
+ switch (vol_index) {
+ case DZ_INDEX: f = DZ_F; invf = DZ_INVF; break;
+ case VR_INDEX: f = VR_F; invf = VR_INVF; break;
+ case SW_INDEX: f = SW_F; invf = SW_INVF; break;
+ case DR_INDEX: f = DR_F; invf = DR_INVF; break;
+ case PH_INDEX: f = PH_F; invf = PH_INVF; break;
+ case RH_INDEX: f = RH_F; invf = RH_INVF; break;
}
- if (wsr88d_ray.ray_hdr.dbptr_vel > 0) {
- vol_index = wsr88d_get_vol_index(wsr88d_ray.vel->data_hdr.dataname);
- switch (vol_index) {
- case DZ_INDEX: f = DZ_F; invf = DZ_INVF; break;
- case VR_INDEX: f = VR_F; invf = VR_INVF; break;
- case SW_INDEX: f = SW_F; invf = SW_INVF; break;
- default: f = DZ_F; invf = DZ_INVF; break;
- }
+ waveform = vcp_data.waveform[isweep];
+ /* The following somewhat complicated if-statement says we'll do the
+ * body of the statement if:
+ * a) the data moment is not reflectivity, or
+ * b) the data moment is reflectivity and one of the following is true:
+ * - waveform is surveillance
+ * - waveform is batch
+ * - elevation is greater than highest split cut (i.e., 6 deg.)
+ */
+ if (vol_index != DZ_INDEX || (waveform == surveillance ||
+ waveform == batch || vcp_data.fixed_angle[isweep] >= 6.0)) {
if (radar->v[vol_index] == NULL) {
radar->v[vol_index] = RSL_new_volume(MAXSWEEPS);
radar->v[vol_index]->h.f = f;
radar->v[vol_index]->sweep[isweep]->h.f = f;
radar->v[vol_index]->sweep[isweep]->h.invf = invf;
}
- ray = RSL_new_ray(wsr88d_ray.vel->data_hdr.ngates);
- wsr88d_load_ray_data(wsr88d_ray.vel, ray);
- wsr88d_load_ray_hdr(wsr88d_ray, ray);
- radar->v[vol_index]->sweep[isweep]->ray[iray] = ray;
- radar->v[vol_index]->sweep[isweep]->h.nrays = iray+1;
- }
+ ngates = data_hdr.ngates;
+ ray = RSL_new_ray(ngates);
- if (wsr88d_ray.ray_hdr.dbptr_sw > 0) {
- vol_index = wsr88d_get_vol_index(wsr88d_ray.sw->data_hdr.dataname);
- switch (vol_index) {
- case DZ_INDEX: f = DZ_F; invf = DZ_INVF; break;
- case VR_INDEX: f = VR_F; invf = VR_INVF; break;
- case SW_INDEX: f = SW_F; invf = SW_INVF; break;
- default: f = DZ_F; invf = DZ_INVF; break;
- }
- if (radar->v[vol_index] == NULL) {
- radar->v[vol_index] = RSL_new_volume(MAXSWEEPS);
- radar->v[vol_index]->h.f = f;
- radar->v[vol_index]->h.invf = invf;
- }
- if (radar->v[vol_index]->sweep[isweep] == NULL) {
- radar->v[vol_index]->sweep[isweep] = RSL_new_sweep(MAXRAYS_M31);
- radar->v[vol_index]->sweep[isweep]->h.f = f;
- radar->v[vol_index]->sweep[isweep]->h.invf = invf;
+ /* Convert data to float, then use range function to store in ray.
+ * Note: data range is 2-255. 0 means signal is below threshold, and 1
+ * means range folded.
+ */
+
+ offset = data_hdr.offset;
+ scale = data_hdr.scale;
+ if (data_hdr.scale == 0) scale = 1.0;
+ data = &wsr88d_ray.data[data_ptr];
+ for (i = 0; i < ngates; i++) {
+ if (data[i] > 1)
+ value = (data[i] - offset) / scale;
+ else value = (data[i] == 0) ? BADVAL : RFVAL;
+ ray->range[i] = invf(value);
+ ray->h.f = f;
+ ray->h.invf = invf;
}
- ray = RSL_new_ray(wsr88d_ray.sw->data_hdr.ngates);
- wsr88d_load_ray_data(wsr88d_ray.sw, ray);
wsr88d_load_ray_hdr(wsr88d_ray, ray);
+ ray->h.range_bin1 = data_hdr.range_first_gate;
+ ray->h.gate_size = data_hdr.range_samp_interval;
+ ray->h.nbins = ngates;
radar->v[vol_index]->sweep[isweep]->ray[iray] = ray;
radar->v[vol_index]->sweep[isweep]->h.nrays = iray+1;
}
-
}
+void wsr88d_load_ray_into_radar(Wsr88d_ray_m31 wsr88d_ray, int isweep, int iray,
+ Radar *radar)
+{
+ int *data_ptr, hdr_size;
+ int i, ndatablocks, nconstblocks = 3;
+
+ hdr_size = sizeof(wsr88d_ray.ray_hdr);
+
+ ndatablocks = wsr88d_ray.ray_hdr.data_block_count;
+ data_ptr = (int *) &wsr88d_ray.ray_hdr.dbptr_ref;
+ for (i=0; i < ndatablocks-nconstblocks; i++) {
+ wsr88d_load_ray(wsr88d_ray, *data_ptr-hdr_size, isweep, iray, radar);
+ data_ptr++;
+ }
+}
+
-void wsr88d_load_sweep_header(Radar *radar, int isweep,
- Wsr88d_ray_m31 wsr88d_ray)
+void wsr88d_load_sweep_header(Radar *radar, int isweep)
{
- int ivolume;
- Ray_header_m31 ray_hdr;
+ int ivolume, nrays;
Sweep *sweep;
- int vcp;
-
- ray_hdr = wsr88d_ray.ray_hdr;
+ Ray *last_ray;
for (ivolume=0; ivolume < MAX_RADAR_VOLUMES; ivolume++) {
- if (radar->v[ivolume] != NULL && radar->v[ivolume]->sweep[isweep] != NULL) {
+ if (radar->v[ivolume] != NULL &&
+ radar->v[ivolume]->sweep[isweep] != NULL) {
sweep = radar->v[ivolume]->sweep[isweep];
- radar->v[ivolume]->sweep[isweep]->h.sweep_num = ray_hdr.elev_num;
- radar->v[ivolume]->sweep[isweep]->h.elev =
- vcp_data.fixed_angle[isweep];
- if (ray_hdr.azm_res != 1)
- sweep->h.beam_width = 1.0;
- else sweep->h.beam_width = 0.5;
+ nrays = sweep->h.nrays;
+ if (nrays == 0) continue;
+ last_ray = sweep->ray[nrays-1];
+ sweep->h.sweep_num = last_ray->h.elev_num;
+ sweep->h.elev = vcp_data.fixed_angle[isweep];
+ sweep->h.beam_width = last_ray->h.beam_width;
sweep->h.vert_half_bw = sweep->h.beam_width / 2.;
sweep->h.horz_half_bw = sweep->h.beam_width / 2.;
}
int msg_hdr_size, msg_size, n;
Radar *radar = NULL;
-/* Message type 31 is a variable length message, whereas all other message
- * types are made up of 2432 byte segments. To handle this, we read the
- * message header and check the message type. If it is not 31, then we read
- * the remainder of the constant size segment. If message type is 31, we read
- * the remainder of the message by the size given in message header.
- * When reading the message header, we must include 12 bytes inserted
- * by RPG, which we ignore, followed by the 8 halfwords (short integers) which
- * make up the actual message header.
- * For more information, see the "Interface Control Document for the RDA/RPG"
- * at the WSR-88D Radar Operations Center web site.
- */
+ /* Message type 31 is a variable length message. All other message types
+ * are made up of 2432-byte segments. To handle all types, we read the
+ * message header and check the message type. If it is not 31, we simply
+ * read the remainder of the 2432-byte segment. If message type is 31, we
+ * use the size given in message header to determine how many bytes to
+ * read. For more information, see the "Interface Control Document for the
+ * RDA/RPG" at the WSR-88D Radar Operations Center web site.
+ */
n = fread(&msghdr, sizeof(Wsr88d_msg_hdr), 1, wf->fptr);
/* printf("msgtype = %d\n", msghdr.msg_type); */
msg_hdr_size = sizeof(Wsr88d_msg_hdr) - sizeof(msghdr.rpg);
-
radar = RSL_new_radar(MAX_RADAR_VOLUMES);
-
+
while (! end_of_vos) {
if (msghdr.msg_type == 31) {
if (little_endian()) wsr88d_swap_m31_hdr(&msghdr);
- /* Get size in bytes of message following message header.
- * The size given in message header is in halfwords, so double it.
+ /* Get size of the remainder of message. The given size is in
+ * halfwords, but we want it in bytes, so double it.
*/
msg_size = (int) msghdr.msg_size * 2 - msg_hdr_size;
n = read_wsr88d_ray_m31(wf, msg_size, &wsr88d_ray);
+ /* Assume error message was issued from read routine */
if (n <= 0) return NULL;
- /* Load this ray into radar structure ray. */
+ /* We need to check radial status for start of new elevation.
+ * Sometimes this occurs without an end-of-elevation flag for the
+ * previous sweep, which is the trigger for loading the sweep
+ * header. "iray" is set to zero after end-of-elev is received,
+ * so that's why we check it. We issue a warning because when this
+ * occurs, the number of rays in the previous sweep is usually less
+ * than expected.
+ */
+ if (wsr88d_ray.ray_hdr.radial_status == START_OF_ELEV &&
+ iray != 0) {
+ fprintf(stderr,"Warning: Radial status is Start-of-Elevation, "
+ "but End-of-Elevation was not\n"
+ "issued for elevation number %d. Number of rays = %d"
+ "\n", isweep+1, iray);
+ wsr88d_load_sweep_header(radar, isweep);
+ isweep++;
+ iray = 0;
+ }
+
+ /* Load ray into radar structure. */
wsr88d_load_ray_into_radar(wsr88d_ray, isweep, iray, radar);
iray++;
+ if (iray >= MAXRAYS_M31) {
+ fprintf(stderr,"Error: iray = %d, equals or exceeds MAXRAYS_M31"
+ " (%d)\n", iray, MAXRAYS_M31);
+ fprintf(stderr,"isweep = %d\n", isweep);
+ RSL_free_radar(radar);
+ return NULL;
+ }
}
else { /* msg_type not 31 */
n = fread(&non31_seg_remainder, sizeof(non31_seg_remainder), 1,
fprintf(stderr,"Warning: load_wsr88d_m31_into_radar: ");
if (feof(wf->fptr) != 0)
fprintf(stderr, "Unexpected end of file.\n");
- else fprintf(stderr,"Read failed.\n");
- fprintf(stderr,"Current sweep number: %d\n"
- "Last ray read: %d\n", isweep+1, iray);
- wsr88d_load_sweep_header(radar, isweep, wsr88d_ray);
+ else
+ fprintf(stderr,"Read failed.\n");
+ fprintf(stderr,"Current sweep index: %d\n"
+ "Last ray read: %d\n", isweep, iray);
+ wsr88d_load_sweep_header(radar, isweep);
return radar;
}
if (msghdr.msg_type == 5) {
/* Check for end of sweep */
if (wsr88d_ray.ray_hdr.radial_status == END_OF_ELEV) {
- wsr88d_load_sweep_header(radar, isweep, wsr88d_ray);
+ wsr88d_load_sweep_header(radar, isweep);
isweep++;
iray = 0;
}
+ /* If not at end of volume scan, read next message header. */
if (wsr88d_ray.ray_hdr.radial_status != END_VOS) {
n = fread(&msghdr, sizeof(Wsr88d_msg_hdr), 1, wf->fptr);
if (n < 1) {
if (feof(wf->fptr) != 0) fprintf(stderr,
"Unexpected end of file.\n");
else fprintf(stderr,"Failed reading msghdr.\n");
- fprintf(stderr,"Current sweep number: %d\n"
- "Last ray read: %d\n", isweep+1, iray);
- wsr88d_load_sweep_header(radar, isweep, wsr88d_ray);
+ fprintf(stderr,"Current sweep index: %d\n"
+ "Last ray read: %d\n", isweep, iray);
+ wsr88d_load_sweep_header(radar, isweep);
return radar;
}
}
else {
end_of_vos = 1;
- wsr88d_load_sweep_header(radar, isweep, wsr88d_ray);
+ wsr88d_load_sweep_header(radar, isweep);
}
if (feof(wf->fptr) != 0) end_of_vos = 1;
- }
+ } /* while not end of vos */
return radar;
}
Wsr88d_tape_header wsr88d_tape_header;
int n;
int nsweep;
+ int i;
int iv;
int nvolumes;
int volume_mask[] = {WSR88D_DZ, WSR88D_VR, WSR88D_SW};
*/
if (n > 0) {
strncpy(version, wsr88d_file_header.title.filename, 8);
- if (strncmp(version,"AR2V0004",8) == 0 ||
- strncmp(version,"AR2V0003",8) ==0 ||
- strncmp(version,"AR2V0002",8) == 0) {
- expected_msgtype = 31;
+ if (strncmp(version,"AR2V0006",8) == 0 ||
+ strncmp(version,"AR2V0004",8) == 0 ||
+ strncmp(version,"AR2V0003",8) == 0 ||
+ strncmp(version,"AR2V0002",8) == 0) {
+ expected_msgtype = 31;
}
else if (strncmp(version,"ARCHIVE2",8) == 0 ||
- strncmp(version,"AR2V0001",8) == 0) {
- expected_msgtype = 1;
+ strncmp(version,"AR2V0001",8) == 0) {
+ expected_msgtype = 1;
}
}
}
}
}
+ if (nsweep == 0) {
+ /* Get Volume Coverage Pattern number for radar header. */
+ i=0;
+ while (i < MAX_RAYS_IN_SWEEP && wsr88d_sweep.ray[i] == NULL) i++;
+ if (i < MAX_RAYS_IN_SWEEP) radar->h.vcp = wsr88d_get_volume_coverage(
+ wsr88d_sweep.ray[i]);
+ }
+
free_and_clear_sweep(&wsr88d_sweep, 0, MAX_RAYS_IN_SWEEP);
}