/* Changes for RSL
*
*---------------------------------------------------------------------
+ * v1.41 In progress
+ *
+ * 1. wsr88d_m31.c: Simplified the WSR-88D ray structure and supporting code.
+ * Removed function read_data_moment.
+ * 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()
+ *---------------------------------------------------------------------
* v1.40 Released 10/10/2008
*
* 1. wsr88d_m31.c (load_wsr88d_m31_into_radar): Corrected a potential problem
/* Byte swapped ? */
if (strncmp(" P A B", magic, 6) == 0) return MCGILL_FILE;
if (strncmp("Volume", magic, 6) == 0) return EDGE_FILE;
+ if (strncmp("SSWB", magic, 4) == 0) return DORADE_FILE;
if (strncmp("VOLD", magic, 4) == 0) return DORADE_FILE;
return UNKNOWN;
dnl Process this file with autoconf to produce a configure script.
-AC_INIT(rsl, v1.40)
+AC_INIT(rsl, v1.40.1)
AC_CONFIG_SRCDIR(volume.c)
AM_INIT_AUTOMAKE
dorade_verbose = 0;
}
+static int do_swap = 0;
+
+/**********************************************************************/
+/* */
+/* read_extra_bytes */
+/* */
+/**********************************************************************/
+int read_extra_bytes(int nbytes, FILE *in)
+{
+ char *extra;
+ int nread;
+ /*
+ * Read and discard nbytes bytes. The return value is the byte count
+ * returned by fread, unless there is an error, then it is 0.
+ */
+ extra = (char *) malloc(nbytes);
+ if (!extra) {
+ perror("\nError: read_extra_bytes");
+ fprintf(stderr,"Tried to allocate %d bytes\n", nbytes);
+ return 0;
+ }
+ nread = fread(extra, sizeof(char), nbytes, in);
+ free(extra);
+ return nread;
+}
+
/**********************************************************************/
/* */
/* dorade_read_comment_block */
return NULL;
}
fread(cb->code, sizeof(cb->code), 1, in);
- fread(&cb->len, sizeof(cb->len), 1, in); cb->len = ntohl(cb->len);
- cb->comment = (char *) calloc(cb->len, sizeof(char));
+ fread(&cb->len, sizeof(cb->len), 1, in);
+
+ /* Check for big endian data on little endian platform. The smallest value
+ * cb->len could have is 8 (length of cb->code + cb->len), so we put that in
+ * first byte of test value, and also test for negative, since sign bit could
+ * be set in a larger value.
+ */
+ if (cb->len > 0x08000000 || cb->len < 0) do_swap = 1;
+ if (do_swap) cb->len = ntohl(cb->len);
+
+ /* Length of cb->comment is cb->len-8 since cb->code and cb->len have
+ * already been read.
+ */
+ cb->comment = (char *) calloc(cb->len-8, sizeof(char));
if (cb->comment == NULL) {
perror("dorade_read_comment_block: cb->comment");
return cb;
}
- fread(cb->comment, sizeof(char), cb->len, in);
+ fread(cb->comment, sizeof(char), cb->len-8, in);
return cb;
}
fread(vd, sizeof(Volume_desc), 1, in);
/* Now, convert from Big Endian. */
- vd->len = ntohl(vd->len);
- vd->version = ntohs(vd->version);
- vd->volume_number = ntohs(vd->volume_number);
- vd->max_bytes = ntohl(vd->max_bytes);
- vd->year = ntohs(vd->year);
- vd->month = ntohs(vd->month);
- vd->day = ntohs(vd->day);
- vd->hour = ntohs(vd->hour);
- vd->minute = ntohs(vd->minute);
- vd->second = ntohs(vd->second);
- vd->gen_year = ntohs(vd->gen_year);
- vd->gen_month = ntohs(vd->gen_month);
- vd->gen_day = ntohs(vd->gen_day);
- vd->nsensors = ntohs(vd->nsensors);
+ if (do_swap) {
+ vd->len = ntohl(vd->len);
+ vd->version = ntohs(vd->version);
+ vd->volume_number = ntohs(vd->volume_number);
+ vd->max_bytes = ntohl(vd->max_bytes);
+ vd->year = ntohs(vd->year);
+ vd->month = ntohs(vd->month);
+ vd->day = ntohs(vd->day);
+ vd->hour = ntohs(vd->hour);
+ vd->minute = ntohs(vd->minute);
+ vd->second = ntohs(vd->second);
+ vd->gen_year = ntohs(vd->gen_year);
+ vd->gen_month = ntohs(vd->gen_month);
+ vd->gen_day = ntohs(vd->gen_day);
+ vd->nsensors = ntohs(vd->nsensors);
+ }
return vd;
}
fread(rd, sizeof(Radar_desc), 1, in);
/* Now, convert from Big Endian. */
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&rd->len);
swap_4_bytes(&rd->radar_constant); /* Yes, even the ieee floating values. */
swap_4_bytes(&rd->peak_power);
swap_4_bytes(&rd->period[i]);
}
}
+ /* If RADD block is longer than structure, read through extra bytes.
+ * This sometimes happens.
+ */
+ if (rd->len > sizeof(Radar_desc)) {
+ if (read_extra_bytes(rd->len - sizeof(Radar_desc), in) <= 0)
+ fprintf(stderr,"Called from %s, line: %d\n", __FILE__, __LINE__ - 1);
+ }
return rd;
}
/**********************************************************************/
fread(pd, sizeof(Parameter_desc), 1, in);
/* Now, convert from Big Endian. */
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&pd->len);
swap_2_bytes(&pd->ipp);
swap_2_bytes(&pd->xmit_freq);
swap_4_bytes(&pd->offset_factor);
swap_4_bytes(&pd->missing_data_flag);
}
+ /* If the descriptor block is longer than the structure, read past the extra bytes.
+ */
+ if (pd->len > sizeof(Parameter_desc)) {
+ if (read_extra_bytes(pd->len - sizeof(Parameter_desc), in) <= 0)
+ fprintf(stderr,"Called from %s, line: %d\n", __FILE__, __LINE__ - 1);
+ }
return pd;
}
fread(&cv->code, sizeof(cv->code), 1, in);
fread(&cv->len, sizeof(cv->len), 1, in);
fread(&cv->ncells, sizeof(cv->ncells), 1, in);
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&cv->len);
swap_4_bytes(&cv->ncells);
}
}
fread(cv->range_cell, sizeof(float), cv->ncells, in);
- if (little_endian()) {
+ if (do_swap) {
for (i=0; i<cv->ncells; i++)
swap_4_bytes(&cv->range_cell[i]);
}
Correction_factor_desc *dorade_read_correction_factor_desc(FILE *in)
{
Correction_factor_desc *cf;
+ char *remaining;
+ int is_cfac = 0;
cf = (Correction_factor_desc *) calloc(1, sizeof(Correction_factor_desc));
if(!cf) {
return NULL;
}
- fread(cf, sizeof(Correction_factor_desc), 1, in);
+ /* Make sure we have Correction Factor Descriptor. */
+ while (!is_cfac) {
+ fread(cf->code, sizeof(cf->code), 1, in);
+ if (strncmp(cf->code, "CFAC", 4) == 0)
+ is_cfac = 1;
+ else {
+ fread(&cf->len, sizeof(cf->len), 1, in);
+ if (do_swap) swap_4_bytes(&cf->len);
+ remaining = (char *) malloc(cf->len-8);
+ if (!remaining) {
+ perror("\ndorade_read_correction_factor_desc");
+ fprintf(stderr,"cf->len = %d\n\n", cf->len);
+ return NULL;
+ }
+ fread(remaining, sizeof(char), cf->len-8, in);
+ free(remaining);
+ }
+ }
+ fread(&cf->len, sizeof(Correction_factor_desc)-4, 1, in);
/* Now, convert from Big Endian. */
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&cf->len);
swap_4_bytes(&cf->azimuth);
swap_4_bytes(&cf->elevation);
}
/* Now, convert from Big Endian. */
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&si->len);
swap_4_bytes(&si->sweep_num);
swap_4_bytes(&si->nrays);
fread(ri, sizeof(Ray_info), 1, in);
/* Now, convert from Big Endian. */
- if (little_endian()) {
+ if (do_swap) {
swap_4_bytes(&ri->len);
swap_4_bytes(&ri->sweep_num);
swap_4_bytes(&ri->jday);
Platform_info *dorade_read_platform_info (FILE *in)
{
Platform_info *pi;
+ int len_first_two;
pi = (Platform_info *) calloc(1, sizeof(Platform_info));
if(!pi) {
return NULL;
}
- fread(pi, sizeof(Platform_info), 1, in);
- /* Now, convert from Big Endian. */
- if (little_endian()) {
- swap_4_bytes(&pi->len);
- swap_4_bytes(&pi->longitude);
- swap_4_bytes(&pi->latitude);
- swap_4_bytes(&pi->altitude);
- swap_4_bytes(&pi->height);
- swap_4_bytes(&pi->ew_speed);
- swap_4_bytes(&pi->ns_speed);
- swap_4_bytes(&pi->v_speed);
- swap_4_bytes(&pi->heading);
- swap_4_bytes(&pi->roll);
- swap_4_bytes(&pi->pitch);
- swap_4_bytes(&pi->drift);
- swap_4_bytes(&pi->rotation);
- swap_4_bytes(&pi->tilt);
- swap_4_bytes(&pi->ew_wind_speed);
- swap_4_bytes(&pi->ns_wind_speed);
- swap_4_bytes(&pi->v_wind_speed);
- swap_4_bytes(&pi->heading_rate);
- swap_4_bytes(&pi->pitch_rate);
+ /* Read the id code to make sure we have "ASIB" for platform info. If
+ * id is ASIB, then read data into the Platform_info structure. If it is
+ * XSTF, read and discard remainder of block, which will have a different
+ * size than Platform_info. XSTF is undocumented, but apparently it takes
+ * the place of ASIB when radar is grounded.
+ */
+
+ fread(pi->code, sizeof(pi->code), 1, in);
+ fread(&pi->len, sizeof(pi->len), 1, in);
+ if (do_swap) swap_4_bytes(&pi->len);
+ len_first_two = sizeof(pi->code) + sizeof(pi->len);
+
+ if (strncmp(pi->code, "ASIB", 4) == 0) {
+ fread(&pi->longitude, sizeof(Platform_info)-len_first_two, 1, in);
+ /* Read past any extra bytes. */
+ if (pi->len > sizeof(Platform_info)) {
+ if (read_extra_bytes(pi->len - sizeof(Platform_info), in) <= 0)
+ fprintf(stderr,"Called from %s, line: %d\n",__FILE__,__LINE__-1);
+ }
+ /* Now, convert from Big Endian. */
+ if (do_swap) {
+ swap_4_bytes(&pi->longitude);
+ swap_4_bytes(&pi->latitude);
+ swap_4_bytes(&pi->altitude);
+ swap_4_bytes(&pi->height);
+ swap_4_bytes(&pi->ew_speed);
+ swap_4_bytes(&pi->ns_speed);
+ swap_4_bytes(&pi->v_speed);
+ swap_4_bytes(&pi->heading);
+ swap_4_bytes(&pi->roll);
+ swap_4_bytes(&pi->pitch);
+ swap_4_bytes(&pi->drift);
+ swap_4_bytes(&pi->rotation);
+ swap_4_bytes(&pi->tilt);
+ swap_4_bytes(&pi->ew_wind_speed);
+ swap_4_bytes(&pi->ns_wind_speed);
+ swap_4_bytes(&pi->v_wind_speed);
+ swap_4_bytes(&pi->heading_rate);
+ swap_4_bytes(&pi->pitch_rate);
+ }
+ } else if (strncmp(pi->code, "XSTF", 4) == 0) {
+ /* Read to end of XSTF block. */
+ if (read_extra_bytes(pi->len - len_first_two, in) <= 0)
+ fprintf(stderr,"Called from %s, line: %d\n", __FILE__, __LINE__ - 1);
+ } else {
+ fprintf(stderr,"Unexpected block id: \"%s\"."
+ " Expected \"ASIB\" or \"XSTF\"\n", pi->code);
}
return pi;
}
+
/**********************************************************************/
/* */
-/* dorade_read_parameter_info */
+/* dorade_read_parameter_data */
/* */
/**********************************************************************/
fread(&pd->code, sizeof(pd->code), 1, in);
fread(&pd->len, sizeof(pd->len), 1, in);
fread(&pd->name, sizeof(pd->name), 1, in);
- if (little_endian()) swap_4_bytes(&pd->len);
+ if (do_swap) swap_4_bytes(&pd->len);
/* Length is in parameter data block? or calculate if from pd->len. */
len = pd->len /* Use pd->len for now. */
parameter_desc[j]->parameter_type == 4)
sr->data_ray[i]->word_size[j] = 4; /* 4 bytes per word */
- if (little_endian()) { /* Numbers were read big-endian. */
+ if (do_swap) { /* Numbers were read big-endian. */
if (sr->data_ray[i]->word_size[j] == 2)
for (k=0; k<len; k+=2)
swap_2_bytes(&pd->data[k]);
return sr;
}
-/**********************************************************************/
-/* */
-/* dorade_read_ray */
-/* */
-/**********************************************************************/
-Data_ray *dorade_read_ray (FILE *in);
-
/* MEMORY MANAGEMENT ROUTINES */
short minute; /* Minute 0-59 (2 byte int) */
short second; /* Second 0-59 (2 byte int) */
char flight_num[8]; /* Flight number (8 ASCII) for airborne raar or IOP number
- for ground based radars. */
+ for ground based radars. */
char facility_name[8]; /* Generation facility. (8 ASCII) */
short gen_year; /* Generation year YYYY. (2 byte int) */
short gen_month; /* Generation month 1-12. (2 byte int) */
float horizontal_beam_width; /* [deg] */
float vertical_beam_width; /* [deg] */
short radar_type; /* 0--Ground
- 1--Airborne fore
- 2--Airborne aft
- 3--Airborne tail
- 4--Airborne lower fuselage
- 5--Shipborne
- */
+ 1--Airborne fore
+ 2--Airborne aft
+ 3--Airborne tail
+ 4--Airborne lower fuselage
+ 5--Shipborne
+ */
short scan_mode; /* 0--Calibration
- 1--PPI (Constant elevation)
- 2--Coplane
- 3--RHI (Constant azimuth)
- 4--Vertical pointing
- 5--Target (Stationary, not vertical pointing)
- 6--Manual
- 7--Idle (out of control)
- 8--Surveillance
- 9--Vertical sweep (rotation axis parallels the fuselage)
- */
+ 1--PPI (Constant elevation)
+ 2--Coplane
+ 3--RHI (Constant azimuth)
+ 4--Vertical pointing
+ 5--Target (Stationary, not vertical pointing)
+ 6--Manual
+ 7--Idle (out of control)
+ 8--Surveillance
+ 9--Vertical sweep (rotation axis parallels the fuselage)
+ */
float scan_rate; /* Nominal scan rate. [deg/sec] */
float start_angle; /* Nominal start angle. [deg] */
float stop_angle; /* Nominal stop angle. [deg] */
short nparam_desc; /* Total number of parameter descriptors for this radar. (2 byte int) */
short ndesc; /* Total number of descriptors for this radar. (2 byte int) */
short compress_code; /* Data compression format code. (2 byte int):
- 0--no compression
- 1--data compression (compression algorithm described
- in the ASCII file athe the begining of the file.
- */
+ 0--no compression
+ 1--data compression (compression algorithm described
+ in the ASCII file athe the begining of the file.
+ */
short compress_algo; /* Data reduction algorithm:
- 0--No data reduction.
- 1--Data recorded between two rotation angles.
- 2--Data recorded between two concentic angles.
- 3--Data recorded between two altitudes.
- 4-N--Other types of data reduction.
- */
+ 0--No data reduction.
+ 1--Data recorded between two rotation angles.
+ 2--Data recorded between two concentic angles.
+ 3--Data recorded between two altitudes.
+ 4-N--Other types of data reduction.
+ */
float data_reduction_param1; /* Data reduction specific parameter #1
- 1--Smallest positive angle [deg].
- 2--Inner circle diameter [km].
- 3--Minimum altitude [km].
- 4-N--Will be defined if other types created.
- */
+ 1--Smallest positive angle [deg].
+ 2--Inner circle diameter [km].
+ 3--Minimum altitude [km].
+ 4-N--Will be defined if other types created.
+ */
float data_reduction_param2; /* Data reduction specific parameter #2
- 1--Largest positive angle [deg].
- 2--Outer circle diameter [km].
- 3--Maximum altitude [km].
- 4-N--Will be defined if other types created.
- */
+ 1--Largest positive angle [deg].
+ 2--Outer circle diameter [km].
+ 3--Maximum altitude [km].
+ 4-N--Will be defined if other types created.
+ */
float longitude; /* Radar longitude [deg]. If airborne, airport longitude. */
float latitude; /* Radar latitude [deg]. If airborne, airport latitude. */
float altitude; /* Radar altitude of mean sea level (msl) [km].
- If airborne, airport altitude. */
+ If airborne, airport altitude. */
float unambiguous_velocity; /* Effective unambiguous velocity [m/s]. */
float unambiguous_range; /* Effective unambiguous range [km]. */
short nfreq; /* Number of freqencies transmitted (2 byte int). */
short npulse_periods; /* Number of different inter-pulse periods (IPP's) transmitted.
- (2 byte int). */
+ (2 byte int). */
float freq[5]; /* Frequency 1..5 [GHz] (float) */
float period[5]; /* Interpulse Period (IPP) 1..5 [ms] (float) */
} Radar_desc;
char description[40]; /* Description of the parameter. */
char units[8]; /* Units (8 ASCII) */
short ipp; /* Inter-pulse periods.
- Bit 0 set to 1 indicates IPP#1 is used in this parameter.
- Similiarly for bits 1,2,3 and 4 and IPP's 2,3,4 and 5.
- */
+ Bit 0 set to 1 indicates IPP#1 is used in this parameter.
+ Similiarly for bits 1,2,3 and 4 and IPP's 2,3,4 and 5.
+ */
short xmit_freq; /* Transmittd frequencies.
- Bit 0 set to 1 indicates Frequency#1 is used in this parameter.
- Similiarly for bits 1,2,3 and 4 and Frequencies 2,3,4 and 5.
- */
+ Bit 0 set to 1 indicates Frequency#1 is used in this parameter.
+ Similiarly for bits 1,2,3 and 4 and Frequencies 2,3,4 and 5.
+ */
float rcvr_bandwidth;/* [MHz] */
short pulse_width; /* [m] */
short polarization; /* 0--Horizontal
- 1--Vertical
- 2--Circular, Right handed
- 3--Elliptical
- 4--Circular, Left handed
- 5--Dual polarization
- */
+ 1--Vertical
+ 2--Circular, Right handed
+ 3--Elliptical
+ 4--Circular, Left handed
+ 5--Dual polarization
+ */
short nsamp_in_dwell_time; /* Number of samples in dwell time. */
short parameter_type; /* 1--8 bit integer
- 2--16 bit integer
- 3--32 bit integer
- 4--floating point (32 bit IEEE)
- */
+ 2--16 bit integer
+ 3--32 bit integer
+ 4--floating point (32 bit IEEE)
+ */
char threshold_field[8];
float threshold_value; /* Units depend on the threshold field. */
float scale_factor; /* Scale factor. */
- float offset_factor; /* meteorological val = (recorded val - offset factor / scale factor */
+ float offset_factor; /* meteorological val = (recorded val - offset factor) / scale factor */
int missing_data_flag; /* Deleted or missing data flag. 256for bytes, -999 for all others. */
} Parameter_desc;
float stop_angle; /* True stop angle [deg]. */
float fixed_angle; /* Fixed angle [deg]. */
int filter_flag; /* Filter flag:
- 0--No filtering in use.
- 1--ON (Algorithm described in ASCII file at beginning of the file.)
- */
+ 0--No filtering in use.
+ 1--ON (Algorithm described in ASCII file at beginning of the file.)
+ */
} Sweep_info;
typedef struct {
short msec; /* Millisecond. */
float azimuth; /* [deg] */
float elevation; /* [deg] */
- float peak_power;/* [kw] */
- float scan_rate; /* Peak transmitted power [deg/sec]. */
+ float peak_power;/* Peak transmitted power [kw] */
+ float scan_rate; /* [deg/sec]. */
int status; /* Ray status:
- 0--Normal
- 1--Transition (antenna repositioning)
- 2--Bad
- 3--Questionable
- */
+ 0--Normal
+ 1--Transition (antenna repositioning)
+ 2--Bad
+ 3--Questionable
+ */
} Ray_info;
typedef struct { /* Especially for moving radars. */
float latitude; /* Radar latitude [deg]. */
float altitude; /* Radar pressure altitude (msl) [km]. */
float height; /* Radar above ground altitude (agl) [km]. */
- float ew_speed; /* Platform ground speed (E->W) [m/s]. */
- float ns_speed; /* Platform ground speed (N->S) [m/s]. */
+ float ew_speed; /* Platform ground speed (E (positive) or W) [m/s]. */
+ float ns_speed; /* Platform ground speed (N (positive) or S) [m/s]. */
float v_speed; /* Platform vertical velocity [m/s]. */
float heading; /* Platform heading [deg]. */
float roll; /* Platform roll [deg]. */
float drift; /* Platform drift [deg]. */
float rotation; /* Platform rotation angle [deg]. */
float tilt; /* Platform tilt [deg]. */
- float ew_wind_speed; /* Horizontal wind speed at radar (East->West) [m/s]. */
- float ns_wind_speed; /* Horizontal wind speed at radar (North->South) [m/s]. */
- float v_wind_speed; /* Vertical wind speed at radar (East->West) [m/s]. */
+ float ew_wind_speed; /* Horizontal wind speed at radar (toward East positive) [m/s]. */
+ float ns_wind_speed; /* Horizontal wind speed at radar (toward North positive) [m/s]. */
+ float v_wind_speed; /* Vertical wind speed at radar (up is positive) [m/s]. */
float heading_rate; /* Heading change rate [deg/sec]. */
float pitch_rate; /* Pitch change rate [deg/sec]. */
} Platform_info;
Platform_info *platform_info;
int nparam;
int *data_len; /* 0..nparam-1
- * Length of *parameter_data[i] in words.
- * Words can be 2 or 4 bytes.
- */
+ * Length of *parameter_data[i] in bytes.
+ * This is length of data portion.
+ */
int *word_size; /* 0..nparam-1
- * Size of each word in *parameter_data[i].
- */
+ * Size of each word in *parameter_data[i].
+ */
Parameter_data **parameter_data; /* 0..nparam-1 */
} Data_ray;
printf(">\n");
printf("len = %d\n", d->len);
printf("d->ncells = %d\n", d->ncells);
+ /*
for (i=0; i<d->ncells; i++)
printf("d->range_cell[%d] = %f\n", i, d->range_cell[i]);
+ */
}
void dorade_print_parameter_desc(Parameter_desc *d)
dorade_print_parameter_desc(s->p_desc[i]);
printf("=====================================================\n");
}
- /* dorade_print_cell_range_vector(s->cell_range_vector); */
+ dorade_print_cell_range_vector(s->cell_range_vector);
dorade_print_correction_factor_desc(s->correction_factor_desc);
}
/* Open input file */
if (file_name == NULL) { /* Use stdin */
- save_fd = dup(0);
- fp = fdopen(save_fd, "r");
+ save_fd = dup(0);
+ fp = fdopen(save_fd, "r");
} else if((fp = fopen(file_name,"r")) == NULL) {
- perror(file_name);
- return fp;
+ perror(file_name);
+ return fp;
}
fp = uncompress_pipe(fp); /* Transparently gunzip. */
if (feof(fp)) return -1;
nbytes = 0;
- while((n = fread(&buf[nbytes], sizeof(char),NSIG_BLOCK-nbytes, fp)) > 0) {
- nbytes += n;
+ while((n = fread(&buf[nbytes], sizeof(char), NSIG_BLOCK-nbytes, fp)) > 0) {
+ nbytes += n;
}
return nbytes;
}
*********************************************************/
void nsig_close(FILE *fp)
{
- rsl_pclose(fp);
+ rsl_pclose(fp);
}
static int do_swap;
printf("id = %d %d\n", (int)rec1->struct_head.id[0], (int)rec1->struct_head.id[1]);
*/
if (rec1->struct_head.id[0] == 0) { /* Possible little-endian */
- if (rec1->struct_head.id[1] >= 20)
- /* This is a big-endian file. Typically, version 2. */
- do_swap = little_endian();
- else
- do_swap = big_endian();
+ if (rec1->struct_head.id[1] >= 20)
+ /* This is a big-endian file. Typically, version 2. */
+ do_swap = little_endian();
+ else
+ do_swap = big_endian();
} else if ((rec1->struct_head.id[1] == 0)) { /* Possible big-endian */
- if (rec1->struct_head.id[0] <= 7)
- /* This is a little-endian file. Version 1. */
- do_swap = big_endian();
+ if (rec1->struct_head.id[0] <= 7)
+ /* This is a little-endian file. Version 1. */
+ do_swap = big_endian();
}
/*
printf("DO SWAP = %d\n", do_swap);
int i=0,itype;
if (s == NULL) return;
for (itype=0; itype<s[0]->nparams; itype++) {
- if (s[itype] == NULL) continue;
-
- if (s[itype]->idh.data_type == NSIG_DTB_EXH)
- free(s[itype]->ray[i]);
- else {
- for (i=0; i<NSIG_I2(s[itype]->idh.num_rays_act); i++)
- nsig_free_ray(s[itype]->ray[i]);
- }
- free(s[itype]->ray);
- free(s[itype]);
+ if (s[itype] == NULL) continue;
+
+ if (s[itype]->idh.data_type == NSIG_DTB_EXH)
+ free(s[itype]->ray[i]);
+ else {
+ for (i=0; i<NSIG_I2(s[itype]->idh.num_rays_act); i++)
+ nsig_free_ray(s[itype]->ray[i]);
+ }
+ free(s[itype]->ray);
+ free(s[itype]);
}
free(s);
}
#define Vprint
#undef Vprint
while(the_code != 1) {
- if (feof(fp)) return -1;
- if (ipos == sizeof(data)) { /* the_code is in the next chunk */
+ if (feof(fp)) return -1;
+ if (ipos == sizeof(data)) { /* the_code is in the next chunk */
#ifdef Vprint
- printf("Exceeded block size looking for the_code. Get it from next buffer.\n");
+ printf("Exceeded block size looking for the_code. Get it from next buffer.\n");
#endif
- n = nsig_read_record(fp, (char *)data);
- if (n <= 0) return n; /* Problem. */
+ n = nsig_read_record(fp, (char *)data);
+ if (n <= 0) return n; /* Problem. */
#ifdef Vprint
- printf("Read %d bytes.\n", n);
- printf("Resetting ipos\n");
+ printf("Read %d bytes.\n", n);
+ printf("Resetting ipos\n");
#endif
- ipos = sizeof(NSIG_Raw_prod_bhdr);
- }
-
+ ipos = sizeof(NSIG_Raw_prod_bhdr);
+ }
+
#ifdef Vprint
- printf("ipos = %d -- ", ipos);
+ printf("ipos = %d -- ", ipos);
#endif
- the_code = NSIG_I2(&data[ipos]);
+ the_code = NSIG_I2(&data[ipos]);
#ifdef Vprint
- printf("the_code = %d (%d) -- ", the_code, (unsigned short)the_code);
+ printf("the_code = %d (%d) -- ", the_code, (unsigned short)the_code);
#endif
- ipos += sizeof(twob);
-
- if (the_code < 0) { /* THIS IS DATA */
- nwords = the_code & 0x7fff;
+ ipos += sizeof(twob);
+
+ if (the_code < 0) { /* THIS IS DATA */
+ nwords = the_code & 0x7fff;
#ifdef Vprint
- printf("#data words (2-bytes) is %d\n", nwords);
+ printf("#data words (2-bytes) is %d\n", nwords);
#endif
- if (ipos + sizeof(twob)*nwords > sizeof(data)) {
+ if (ipos + sizeof(twob)*nwords > sizeof(data)) {
#ifdef Vprint
- printf("Exceeded block size... transferring and reading new (i=%d).\n", i);
+ printf("Exceeded block size... transferring and reading new (i=%d).\n", i);
#endif
- /* Need another phyical block. */
-
- /* But, first transfer the remainder to the output chunk. */
- /* And, transfer begining of next block */
- /* Transfer end of current buffer. */
- end_nwords = (NSIG_BLOCK - ipos)/sizeof(twob);
- memmove(&chunk[i], &data[ipos], sizeof(twob)*end_nwords);
- i += end_nwords * sizeof(twob);
-
- n = nsig_read_record(fp, (char *)data);
- if (n <= 0) return n; /* Problem. */
- /* New ipos */
- nwords -= end_nwords;
- ipos = sizeof(NSIG_Raw_prod_bhdr);
-
- /* Transfer beginning of new buffer */
- if (i+nwords * sizeof(twob) > NSIG_BLOCK) return -1;
- memmove(&chunk[i], &data[ipos], sizeof(twob) * nwords);
- i += nwords * sizeof(twob);
+ /* Need another phyical block. */
+
+ /* But, first transfer the remainder to the output chunk. */
+ /* And, transfer begining of next block */
+ /* Transfer end of current buffer. */
+ end_nwords = (NSIG_BLOCK - ipos)/sizeof(twob);
+ memmove(&chunk[i], &data[ipos], sizeof(twob)*end_nwords);
+ i += end_nwords * sizeof(twob);
+
+ n = nsig_read_record(fp, (char *)data);
+ if (n <= 0) return n; /* Problem. */
+ /* New ipos */
+ nwords -= end_nwords;
+ ipos = sizeof(NSIG_Raw_prod_bhdr);
+
+ /* Transfer beginning of new buffer */
+ if (i+nwords * sizeof(twob) > NSIG_BLOCK) return -1;
+ memmove(&chunk[i], &data[ipos], sizeof(twob) * nwords);
+ i += nwords * sizeof(twob);
ipos += nwords * sizeof(twob);
-
+
#ifdef Vprint
- printf("Words to transfer (at end of block) is %d\n", end_nwords);
- printf("Transfer %d words from beginning of next buffer.\n", nwords);
- printf("ipos in new buffer is %d\n", ipos);
+ printf("Words to transfer (at end of block) is %d\n", end_nwords);
+ printf("Transfer %d words from beginning of next buffer.\n", nwords);
+ printf("ipos in new buffer is %d\n", ipos);
#endif
- } else { /* Normal situation. Position to end of data.
- * But, first transfer it to the chunk.
- */
- if (i+nwords * sizeof(twob) > NSIG_BLOCK) return -1;
- memmove(&chunk[i], &data[ipos], sizeof(twob) * nwords);
- i += nwords * sizeof(twob);
- ipos += sizeof(twob) * nwords;
- }
-
- } else if (the_code == 1) { /* END OF THE RAY. */
+ } else { /* Normal situation. Position to end of data.
+ * But, first transfer it to the chunk.
+ */
+ if (i+nwords * sizeof(twob) > NSIG_BLOCK) return -1;
+ memmove(&chunk[i], &data[ipos], sizeof(twob) * nwords);
+ i += nwords * sizeof(twob);
+ ipos += sizeof(twob) * nwords;
+ }
+
+ } else if (the_code == 1) { /* END OF THE RAY. */
#ifdef Vprint
- printf("------------------------------> Reached end of ray.\n");
+ printf("------------------------------> Reached end of ray.\n");
#endif
- break; /* or continue; */
+ break; /* or continue; */
- } else if (the_code == 0) { /* UNKNOWN */
- break;
- } else { /* NUMBER OF ZERO's */
+ } else if (the_code == 0) { /* UNKNOWN */
+ break;
+ } else { /* NUMBER OF ZERO's */
#ifdef Vprint
- printf("#000000000000 to skip is %d (i=%d)\n", the_code, i);
+ printf("#000000000000 to skip is %d (i=%d)\n", the_code, i);
#endif
- if (i+the_code * sizeof(twob) > NSIG_BLOCK) return -1;
- memset(&chunk[i], 0, the_code*sizeof(twob));
- i += the_code * sizeof(twob);
- }
-
- if (ipos >= sizeof(data)) {
+ if (i+the_code * sizeof(twob) > NSIG_BLOCK) return -1;
+ memset(&chunk[i], 0, the_code*sizeof(twob));
+ i += the_code * sizeof(twob);
+ }
+
+ if (ipos >= sizeof(data)) {
#ifdef Vprint
- printf("Exceeded block size ... ipos = %d\n", ipos);
- printf("This should be right at the end of the block.\n");
+ printf("Exceeded block size ... ipos = %d\n", ipos);
+ printf("This should be right at the end of the block.\n");
#endif
- n = nsig_read_record(fp, (char *)data);
- if (n <= 0) return n; /* Problem. */
- ipos = sizeof(NSIG_Raw_prod_bhdr);
- }
-
- } /* End while. */
+ n = nsig_read_record(fp, (char *)data);
+ if (n <= 0) return n; /* Problem. */
+ ipos = sizeof(NSIG_Raw_prod_bhdr);
+ }
+
+ } /* End while. */
return i;
}
if (n == 0) return NULL; /* Silent error. */
if (n < 0) {
- fprintf(stderr, "nsig_read_ray: chunk return code = %d.\n", n);
- return NULL;
+ fprintf(stderr, "nsig_read_ray: chunk return code = %d.\n", n);
+ return NULL;
}
if (n > NSIG_BLOCK) { /* Whoa! */
- fprintf(stderr, "nsig_read_ray: chunk bigger than buffer. n = %d,\
+ fprintf(stderr, "nsig_read_ray: chunk bigger than buffer. n = %d,\
maximum block size allowed is %d\n", n, NSIG_BLOCK);
- return NULL;
+ return NULL;
}
ray = (NSIG_Ray *) calloc(1, sizeof(NSIG_Ray));
#ifdef Vprint
printf(" rayh.num_bins = %d (nbins %d, n %d)\n", NSIG_I2(rayh.num_bins), nbins, n);
#endif
+ ray->range = (unsigned char *)calloc(n, sizeof(unsigned char));
+ /* Changed calloc nbins to calloc n for 2-byte data.
ray->range = (unsigned char *)calloc(nbins, sizeof(unsigned char));
memmove(ray->range, &chunk[sizeof(NSIG_Ray_header)], nbins);
+ Can remove this commented-out code once we know changes work.
+ */
+ memmove(ray->range, &chunk[sizeof(NSIG_Ray_header)], n);
return ray;
}
NSIG_Sweep **nsig_read_sweep(FILE *fp, NSIG_Product_file *prod_file)
{
NSIG_Sweep **s;
- int i, n;
+ int i, j, n;
static NSIG_Ingest_data_header **idh = NULL;
static NSIG_Raw_prod_bhdr *bhdr = NULL;
NSIG_Ray *nsig_ray;
int data_mask, iray, nrays[12], max_rays;
+ int masks[5];
int nparams;
int is_new_ray;
int idtype[12];
(void)nsig_endianess(&prod_file->rec1);
/* Setup the array of ingest data headers [0..nparams-1] */
+#ifdef NSIG_VER2
+ memmove(&masks[0], prod_file->rec2.task_config.dsp_info.data_mask_cur.mask_word_0,
+ sizeof(fourb));
+ memmove(&masks[1], &prod_file->rec2.task_config.dsp_info.data_mask_cur.mask_word_1,
+ 4*sizeof(fourb));
+ nparams = 0;
+ for (j=0; j < 5; j++) {
+ data_mask = masks[j];
+ for (i=0; i<32; i++)
+ nparams += (data_mask >> i) & 0x1;
+ }
+#else
memmove(&data_mask, prod_file->rec2.task_config.dsp_info.data_mask, sizeof(fourb));
+ for (nparams=i=0; i<32; i++)
+ nparams += (data_mask >> i) & 0x1;
#ifdef Vprint
printf("data_mask %x\n", data_mask);
#endif
- for (nparams=i=0; i<32; i++)
- nparams += (data_mask >> i) & 0x1;
+#endif
/* Number of sweeps */
#ifdef Vprint
{int nsweeps;
if (idh == NULL) {
-
- idh = (NSIG_Ingest_data_header **)calloc(nparams, sizeof(NSIG_Ingest_data_header *));
- ipos = 0;
- for (i=0; i<nparams; i++) {
- idh[i] = (NSIG_Ingest_data_header *)&data[sizeof(NSIG_Raw_prod_bhdr) + i*sizeof(NSIG_Ingest_data_header)];
- }
- bhdr = (NSIG_Raw_prod_bhdr *)prod_file->data;
+
+ idh = (NSIG_Ingest_data_header **)calloc(nparams, sizeof(NSIG_Ingest_data_header *));
+ ipos = 0;
+ for (i=0; i<nparams; i++) {
+ idh[i] = (NSIG_Ingest_data_header *)&data[sizeof(NSIG_Raw_prod_bhdr) + i*sizeof(NSIG_Ingest_data_header)];
+ }
+ bhdr = (NSIG_Raw_prod_bhdr *)prod_file->data;
}
xh_size = NSIG_I2(prod_file->rec2.ingest_head.size_ext_ray_headers);
/* Now pointers to all possible rays. */
for (i=0; i<nparams; i++) {
/* Load sweep headers */
- s[i] = (NSIG_Sweep *) calloc (nparams, sizeof(NSIG_Sweep));
- s[i]->nparams = nparams;
- memmove(&s[i]->bhdr, &bhdr, sizeof(NSIG_Raw_prod_bhdr));
- memmove(&s[i]->idh, idh[i], sizeof(NSIG_Ingest_data_header));
- s[i]->ray = (NSIG_Ray **) calloc (max_rays, sizeof(NSIG_Ray *));
- }
+ s[i] = (NSIG_Sweep *) calloc (nparams, sizeof(NSIG_Sweep));
+ s[i]->nparams = nparams;
+ memmove(&s[i]->bhdr, &bhdr, sizeof(NSIG_Raw_prod_bhdr));
+ memmove(&s[i]->idh, idh[i], sizeof(NSIG_Ingest_data_header));
+ s[i]->ray = (NSIG_Ray **) calloc (max_rays, sizeof(NSIG_Ray *));
+ }
/* Process this sweep. Keep track of the end of the ray. */
ipos = sizeof(NSIG_Raw_prod_bhdr); /* Position in the 'data' array */
max_rays = 0;
for (i=0; i<nparams; i++) {
- idtype[i] = NSIG_I2(idh[i]->data_type);
- nrays[i] = (int)NSIG_I2(idh[i]->num_rays_act);
- if (nrays[i] > max_rays) max_rays = nrays[i];
+ idtype[i] = NSIG_I2(idh[i]->data_type);
+ nrays[i] = (int)NSIG_I2(idh[i]->num_rays_act);
+ if (nrays[i] > max_rays) max_rays = nrays[i];
#ifdef Vprint
- printf("New ray: parameter %d has idtype=%d\n", i, idtype[i]);
- printf("Number of expected rays in sweep %d is %d\n", isweep, (int)NSIG_I2(idh[i]->num_rays_exp));
- printf("Number of actual rays in sweep %d is %d\n", isweep, (int)NSIG_I2(idh[i]->num_rays_act));
+ printf("New ray: parameter %d has idtype=%d\n", i, idtype[i]);
+ printf("Number of expected rays in sweep %d is %d\n", isweep, (int)NSIG_I2(idh[i]->num_rays_exp));
+ printf("Number of actual rays in sweep %d is %d\n", isweep, (int)NSIG_I2(idh[i]->num_rays_act));
#endif
}
if (is_new_sweep)
- ipos += nparams * sizeof(NSIG_Ingest_data_header);
+ ipos += nparams * sizeof(NSIG_Ingest_data_header);
- /* ipos = sizeof(NSIG_Raw_prod_bhdr) + nparams*sizeof(NSIG_Ingest_data_header); */
+ /* ipos = sizeof(NSIG_Raw_prod_bhdr) + nparams*sizeof(NSIG_Ingest_data_header); */
/* 'iray' is the true ray index into 's', whereas, 'nsig_iray' is what
* the NSIG file says it is. I'll trust 'iray'
*
do {
#ifdef Vprint
- printf("---------------------- New Ray <%d> --------------------\n", iray);
+ printf("---------------------- New Ray <%d> --------------------\n", iray);
#endif
- if (feof(fp)) { /* Premature eof */
- return NULL; /* This will have to do. */
- }
- /* For all parameters present. */
- is_new_ray = 0;
- for (i=0; i<nparams; i++) {
-
- /* Keep track of the cursor within the buffer, and, possibly
+ if (feof(fp)) { /* Premature eof */
+ return NULL; /* This will have to do. */
+ }
+ /* For all parameters present. */
+ is_new_ray = 0;
+ for (i=0; i<nparams; i++) {
+
+ /* Keep track of the cursor within the buffer, and, possibly
* read another buffer when a ray is split across two NSIG blocks
- */
- nsig_ray = NULL;
- if (idtype[i] != 0) { /* Not an extended header. */
- nsig_ray = nsig_read_ray(fp);
-
- } else { /* Check extended header version. */
- if (xh_size <= 20) {
- exh0 = nsig_read_ext_header_ver0(fp);
- if (exh0) {
- nsig_ray = (NSIG_Ray *)calloc(1, sizeof(NSIG_Ray));
- nsig_ray->range = (unsigned char *)exh0;
- }
- } else {
- exh1 = nsig_read_ext_header_ver1(fp);
- if (exh1) {
- nsig_ray = (NSIG_Ray *)calloc(1, sizeof(NSIG_Ray));
- nsig_ray->range = (unsigned char *)exh1;
- }
- }
- }
- if (nsig_ray) is_new_ray = 1;
- if (iray > nrays[i]) break;
- s[i]->ray[iray] = nsig_ray;
-
- } /* End for */
- if (is_new_ray) iray++;
-
+ */
+ nsig_ray = NULL;
+ if (idtype[i] != 0) { /* Not an extended header. */
+ nsig_ray = nsig_read_ray(fp);
+
+ } else { /* Check extended header version. */
+ if (xh_size <= 20) {
+ exh0 = nsig_read_ext_header_ver0(fp);
+ if (exh0) {
+ nsig_ray = (NSIG_Ray *)calloc(1, sizeof(NSIG_Ray));
+ nsig_ray->range = (unsigned char *)exh0;
+ }
+ } else {
+ exh1 = nsig_read_ext_header_ver1(fp);
+ if (exh1) {
+ nsig_ray = (NSIG_Ray *)calloc(1, sizeof(NSIG_Ray));
+ nsig_ray->range = (unsigned char *)exh1;
+ }
+ }
+ }
+ if (nsig_ray) is_new_ray = 1;
+ if (iray > nrays[i]) break;
+ s[i]->ray[iray] = nsig_ray;
+
+ } /* End for */
+ if (is_new_ray) iray++;
+
} while (iray < max_rays);
#ifdef Vprint
printf("iray = %d\n", iray);
#define NSIG_DTB_VEL 3
#define NSIG_DTB_WID 4
#define NSIG_DTB_ZDR 5
+#define NSIG_DTB_UCR2 8
+#define NSIG_DTB_CR2 9
+#define NSIG_DTB_VEL2 10
+#define NSIG_DTB_WID2 11
+#define NSIG_DTB_ZDR2 12
#define NSIG_DTB_KDP 14
+#define NSIG_DTB_KDP2 15
#define NSIG_DTB_PHIDP 16
+#define NSIG_DTB_VELC 17
#define NSIG_DTB_SQI 18
#define NSIG_DTB_RHOHV 19
+#define NSIG_DTB_RHOHV2 20
+#define NSIG_DTB_VELC2 22
+#define NSIG_DTB_SQI2 23
+#define NSIG_DTB_PHIDP2 24
+#define NSIG_DTB_HCLASS 55
+#define NSIG_DTB_HCLASS2 56
/* Product type code ,value for byte 12 in product configuration
* struct, III-35
* Used for warning, shearline, and track only */
twob iresults_count;
-#define PROD_END_PAD (20)
+#define PROD_END_PAD 20
char ipad_end[PROD_END_PAD];
#else
char color_unscan; /* Color used for unscanned area. */
/*============================================================*/
/*============================================================*/
+#ifdef NSIG_VER2
+typedef struct {
+ twob low_prf; /* Hertz */
+ twob low_prf_frac; /* Fraction part, scaled by 2**-16 */
+ twob low_prf_sample_size;
+ twob low_prf_range_averaging; /* In bins */
+ twob thresh_refl_unfolding; /* Threshold for reflectivity unfolding in 1/100 dB */
+ twob thresh_vel_unfolding; /* Threshold for velocity unfolding in 1/100 dB */
+ twob thresh_sw_unfolding; /* Threshold for width unfolding in 1/100 dB */
+ char spare[18];
+} NSIG_Task_dsp_mode_batch;
+#endif
+/*============================================================*/
+/*============================================================*/
+
+#ifdef NSIG_VER2
+typedef struct {
+ fourb mask_word_0;
+ fourb ext_hdr_type;
+ fourb mask_word_1;
+ fourb mask_word_2;
+ fourb mask_word_3;
+ fourb mask_word_4;
+} NSIG_Dsp_data_mask;
+#endif
+/*============================================================*/
+/*============================================================*/
+
/* Task dsp info 3.3.2.2, page III-22 */
typedef struct {
#ifdef NSIG_VER2
twob dsp_num;
fourb dsp_type;
+ NSIG_Dsp_data_mask data_mask_cur;
+ NSIG_Dsp_data_mask data_mask_orig;
+ NSIG_Task_dsp_mode_batch task_dsp_mode;
+ char spare[52];
#else
twob dsp_num;
twob dsp_type;
-#endif
fourb data_mask;
fourb aux_data_def[32];
+#endif
fourb prf;
fourb pwid;
twob prf_mode;
#endif
twob atten_gain;
#ifdef NSIG_VER2
- twob igas_atten ; /* 100000 * db/km */
-#define TASK_DSP_INFO_PAD (148)
+ twob igas_atten; /* 100000 * db/km */
+ twob clutter_map_flag;
+ twob xmt_phase_seq;
+ fourb ray_hdr_mask;
+ twob time_series_playback;
+ twob spare2;
+ char custom_ray_hdr_name[16];
+#define TASK_DSP_INFO_PAD 120
char ipad_end[TASK_DSP_INFO_PAD];
#else
char spare[150];
#define MAX_NSIG_SWEEPS 30
#define MAX_NSIG_RAYS 400
#define NSIG_NO_ECHO -32.0
+#define NSIG_NO_ECHO2 -999.0
static float (*f)(Range x);
static Range (*invf)(float x);
FILE *file;
void get_extended_header_info(NSIG_Sweep **nsig_sweep, int xh_size, int iray,
- int nparams,
- int *msec, float *azm, float *elev,
- float *pitch, float *roll, float *heading,
- float *azm_rate, float *elev_rate,
- float *pitch_rate, float *roll_rate, float *heading_rate,
- float *lat, float *lon, int *alt, float *rvc,
- float *vel_east, float *vel_north, float *vel_up)
+ int nparams,
+ int *msec, float *azm, float *elev,
+ float *pitch, float *roll, float *heading,
+ float *azm_rate, float *elev_rate,
+ float *pitch_rate, float *roll_rate,
+ float *heading_rate,
+ float *lat, float *lon, int *alt, float *rvc,
+ float *vel_east, float *vel_north, float *vel_up)
{
static NSIG_Ext_header_ver1 xh;
int data_type, itype;
/* Determine where 'itype' for extended header is. */
for (itype = 0; itype<nparams; itype++) {
data_type = NSIG_I2(nsig_sweep[itype]->idh.data_type);
- if (data_type == NSIG_DTB_EXH) break;
+ if (data_type == NSIG_DTB_EXH) break;
}
/* printf("...extended header itype=%d, nparams=%d\n", itype, nparams); */
if (itype == nparams) return; /* No extended header. */
*msec = NSIG_I4(xh.msec);
/* printf("...extended header msec= %d\n", *msec); */
if (xh_size <= 20) /* Stop, only have version 0. */
- return;
+ return;
/* Version 1 processing. */
*azm = nsig_from_bang(xh.azm);
#endif
(char *filename)
{
- FILE *fp;
- /* RSL structures */
- Radar *radar;
- Ray *ray;
-
- int i, j, k, n;
- int year, month, day;
- int hour, minute, sec;
- int numbins, numsweep;
- int num_rays, sea_lvl_hgt;
- int radar_number, num_samples;
- int latd, latm, lats, lond, lonm, lons;
- int data_type;
- int bin_num;
- int sweep_year, sweep_day, sweep_month;
- int sweep_hour, sweep_minute, sweep_second;
- int sweep_sec;
- int z_flag_unc, z_flag_cor, v_flag, w_flag, speckle;
- int ant_scan_mode;
- float second;
- float pw;
- float bin_space;
- float prf, wave, beam_width;
- float vert_half_bw, horz_half_bw;
- float rng_last_bin;
- float rng_first_bin, freq;
- float max_vel, sweep_rate, azim_rate;
- float ray_data;
- float az1, az2;
- double tmp;
- float sqi, log, csr, sig, cal_dbz;
- char radar_type[50], state[2], city[15];
- char site_name[16];
+ FILE *fp;
+ /* RSL structures */
+ Radar *radar;
+ Ray *ray;
+
+ int i, j, k, n;
+ int year, month, day;
+ int hour, minute, sec;
+ int numbins, numsweep;
+ int num_rays, sea_lvl_hgt;
+ int radar_number, num_samples;
+ int latd, latm, lats, lond, lonm, lons;
+ int data_type;
+ int bin_num;
+ int sweep_year, sweep_day, sweep_month;
+ int sweep_hour, sweep_minute, sweep_second;
+ int sweep_sec;
+ int z_flag_unc, z_flag_cor, v_flag, w_flag, speckle;
+ int ant_scan_mode;
+ float second;
+ float pw;
+ float bin_space;
+ float prf, wave, beam_width;
+ float vert_half_bw, horz_half_bw;
+ float rng_last_bin;
+ float rng_first_bin, freq;
+ float max_vel, sweep_rate, azim_rate;
+ float ray_data;
+ float az1, az2;
+ double tmp;
+ float sqi, log, csr, sig, cal_dbz;
+ char radar_type[50], state[2], city[15];
+ char site_name[16];
NSIG_Product_file *prod_file;
short id;
int data_mask, nrays;
+ int masks[5];
int nparams, nsweeps;
NSIG_Sweep **nsig_sweep;
NSIG_Ray *ray_p;
int itype, ifield;
+ unsigned short nsig_u2byte; /* New for 2-byte data types, Aug 2009 */
Sweep *sweep;
int msec;
float azm, elev, pitch, roll, heading, azm_rate, elev_rate,
- pitch_rate, roll_rate, heading_rate,
- lat, lon;
+ pitch_rate, roll_rate, heading_rate,
+ lat, lon;
int alt; /* Altitude */
float rvc; /* Radial correction velocity m/s */
float vel_east, vel_north, vel_up; /* Platform velocity vectors m/sec */
int xh_size;
+ float incr;
extern int *rsl_qsweep; /* See RSL_read_these_sweeps in volume.c */
extern int rsl_qsweep_max;
extern float rsl_kdp_wavelen;
radar = NULL;
if (radar_verbose_flag)
- fprintf(stderr, "open file: %s\n", filename);
+ fprintf(stderr, "open file: %s\n", filename);
/** Opening nsig file **/
if((fp = nsig_open(filename)) == NULL) return NULL;
n = nsig_read_record(fp, (char *)&prod_file->rec1);
nsig_endianess(&prod_file->rec1);
if (radar_verbose_flag)
- fprintf(stderr, "Read %d bytes for rec1.\n", n);
+ fprintf(stderr, "Read %d bytes for rec1.\n", n);
id = NSIG_I2(prod_file->rec1.struct_head.id);
if (radar_verbose_flag)
- fprintf(stderr, "ID = %d\n", (int)id);
+ fprintf(stderr, "ID = %d\n", (int)id);
if (id != 7 && id != 27) { /* testing: Use 27 for Version 2 data */
- fprintf(stderr, "File is not a SIGMET version 1 nor version 2 raw product file.\n");
- return NULL;
+ fprintf(stderr, "File is not a SIGMET version 1 nor version 2 raw product file.\n");
+ return NULL;
}
n = nsig_read_record(fp, (char *)&prod_file->rec2);
if (radar_verbose_flag)
- fprintf(stderr, "Read %d bytes for rec2.\n", n);
+ fprintf(stderr, "Read %d bytes for rec2.\n", n);
/** Test for scan mode -- If scan is a RHI will return NULL **/
/** because RSL can't handle RHI's. In the future, replace **/
/** NULL will a routine to convert RHI's to RSL Format **/
ant_scan_mode =NSIG_I2(prod_file->rec2.task_config.scan_info.ant_scan_mode);
if(ant_scan_mode == 2)
- {
- if (radar_verbose_flag)
- fprintf(stderr, "RHI scan detected. Unable to process, returning NULL.\n");
- /* return NULL; */
- }
+ {
+ if (radar_verbose_flag)
+ fprintf(stderr, "RHI scan detected. Unable to process, returning NULL.\n");
+ /* return NULL; */
+ }
/* Count the bits set in 'data_mask' to determine the number
* of parameters present.
xh_size = NSIG_I2(prod_file->rec2.ingest_head.size_ext_ray_headers);
nrays = NSIG_I2(prod_file->rec2.ingest_head.num_rays);
if (radar_verbose_flag)
- fprintf(stderr, "Expecting %d rays in each sweep.\n", nrays);
+ fprintf(stderr, "Expecting %d rays in each sweep.\n", nrays);
+#ifdef NSIG_VER2
+ memmove(&masks[0], prod_file->rec2.task_config.dsp_info.data_mask_cur.mask_word_0,
+ sizeof(fourb));
+ memmove(&masks[1], &prod_file->rec2.task_config.dsp_info.data_mask_cur.mask_word_1,
+ 4*sizeof(fourb));
+ nparams = 0;
+ for (j=0; j < 5; j++) {
+ data_mask = masks[j];
+ for (i=0; i<32; i++)
+ nparams += (data_mask >> i) & 0x1;
+ }
+#else
memmove(&data_mask, prod_file->rec2.task_config.dsp_info.data_mask, sizeof(fourb));
for (nparams=i=0; i<32; i++)
- nparams += (data_mask >> i) & 0x1;
+ nparams += (data_mask >> i) & 0x1;
+#endif
/* Number of sweeps */
nsweeps = NSIG_I2(prod_file->rec2.task_config.scan_info.num_swp);
memmove(site_name, prod_file->rec1.prod_end.site_name, sizeof(prod_file->rec1.prod_end.site_name));
site_name[sizeof(site_name)-1] = '\0';
if (radar_verbose_flag) {
- fprintf(stderr, "nparams = %d, nsweeps = %d\n", nparams, nsweeps);
- fprintf(stderr, "Site name = <%s>\n", site_name);
+ fprintf(stderr, "nparams = %d, nsweeps = %d\n", nparams, nsweeps);
+ fprintf(stderr, "Site name = <%s>\n", site_name);
}
/* nsig_sweep = nsig_read_sweep(fp, prod_file)
sea_lvl_hgt = NSIG_I2(prod_file->rec1.prod_end.grnd_sea_ht);
if (radar_verbose_flag)
- fprintf(stderr, "sea: %d\n", sea_lvl_hgt);
+ fprintf(stderr, "sea: %d\n", sea_lvl_hgt);
if (radar_verbose_flag)
- fprintf(stderr, "site_name: %s", site_name);
+ fprintf(stderr, "site_name: %s", site_name);
/** Determine beamwidth from input variables (not saved in nsig file) **/
if(strncmp(site_name,"mit",3) == 0 || strncmp(site_name,"MIT",3) == 0)
- beam_width = MIT_BEAMWIDTH;
+ beam_width = MIT_BEAMWIDTH;
else if(strncmp(site_name,"tog",3) == 0 || strncmp(site_name,"TOG",3) == 0)
- beam_width = TOG_BEAMWIDTH;
+ beam_width = TOG_BEAMWIDTH;
else if(strncmp(site_name,"kwa",3) == 0 || strncmp(site_name,"KWA",3) == 0)
- beam_width = KWA_BEAMWIDTH;
+ beam_width = KWA_BEAMWIDTH;
else
- beam_width = DEFAULT_BEAMWIDTH;
+ beam_width = DEFAULT_BEAMWIDTH;
if (radar_verbose_flag)
- fprintf(stderr, "beamwidth: %f\n", beam_width);
+ fprintf(stderr, "beamwidth: %f\n", beam_width);
vert_half_bw = beam_width/2.0;
horz_half_bw = beam_width/2.0;
prf = NSIG_I4(prod_file->rec1.prod_end.prf); /* pulse repetition frequency */
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.
- */
+ * to operate with the proper wavelength.
+ */
numbins = NSIG_I4(prod_file->rec1.prod_end.num_bin); /* # bins in ray */
rng_first_bin = (float)NSIG_I4(prod_file->rec1.prod_end.rng_f_bin)/100.0;
rng_last_bin = (float)NSIG_I4(prod_file->rec1.prod_end.rng_l_bin)/100.0;
/** Verbose calibration information **/
if (radar_verbose_flag)
- {
- fprintf(stderr, "LOG = %5.2f\n", log);
- fprintf(stderr, "SQI = %5.2f\n", sqi);
- fprintf(stderr, "CSR = %5.2f\n", csr);
- fprintf(stderr, "SIG = %5.2f\n", sig);
- fprintf(stderr, "Calibration reflectivity: %5.2f dBZ\n", cal_dbz);
- fprintf(stderr, "ZT flags: %d\n", z_flag_unc); /** can find these **/
- fprintf(stderr, "DZ flags: %d\n", z_flag_cor); /** defn in the **/
- fprintf(stderr, "VR flags: %d\n", v_flag); /** SIGMET Doc **/
- fprintf(stderr, "SW flags: %d\n", w_flag);
- fprintf(stderr, "Flags: -3856 = SQI thresholding\n");
- fprintf(stderr, " -21846 = LOG thresholding\n");
- fprintf(stderr, " -24416 = LOG & SQI thresholding\n");
- fprintf(stderr, " -24516 = LOG & SQI & SIG thresholding\n");
- fprintf(stderr, "speckle remover: %d\n", speckle);
- }
+ {
+ fprintf(stderr, "LOG = %5.2f\n", log);
+ fprintf(stderr, "SQI = %5.2f\n", sqi);
+ fprintf(stderr, "CSR = %5.2f\n", csr);
+ fprintf(stderr, "SIG = %5.2f\n", sig);
+ fprintf(stderr, "Calibration reflectivity: %5.2f dBZ\n", cal_dbz);
+ fprintf(stderr, "ZT flags: %d\n", z_flag_unc); /** can find these **/
+ fprintf(stderr, "DZ flags: %d\n", z_flag_cor); /** defn in the **/
+ fprintf(stderr, "VR flags: %d\n", v_flag); /** SIGMET Doc **/
+ fprintf(stderr, "SW flags: %d\n", w_flag);
+ fprintf(stderr, "Flags: -3856 = SQI thresholding\n");
+ fprintf(stderr, " -21846 = LOG thresholding\n");
+ fprintf(stderr, " -24416 = LOG & SQI thresholding\n");
+ fprintf(stderr, " -24516 = LOG & SQI & SIG thresholding\n");
+ fprintf(stderr, "speckle remover: %d\n", speckle);
+ }
if (radar_verbose_flag)
- fprintf(stderr, "vel: %f prf: %f\n", max_vel, prf);
+ fprintf(stderr, "vel: %f prf: %f\n", max_vel, prf);
/** Extracting Latitude and Longitude from nsig file **/
lat = nsig_from_fourb_ang(prod_file->rec2.ingest_head.lat_rad);
if(lat > 180.0) lat -= 360.0;
if(lon > 180.0) lon -= 360.0;
if (radar_verbose_flag)
- fprintf(stderr, "nsig_to_radar: lat %f, lon %f\n", lat, lon);
+ fprintf(stderr, "nsig_to_radar: lat %f, lon %f\n", lat, lon);
/** Latitude deg, min, sec **/
latd = (int)lat;
tmp = (lat - latd) * 60.0;
/** Allocating memory for radar structure **/
radar = RSL_new_radar(MAX_RADAR_VOLUMES);
if (radar == NULL)
- {
- fprintf(stderr, "nsig_to_radar: radar is NULL\n");
- return NULL;
- }
+ {
+ fprintf(stderr, "nsig_to_radar: radar is NULL\n");
+ return NULL;
+ }
/** Filling Radar Header **/
radar->h.month = month;
if (radar_verbose_flag) {
#ifdef NSIG_VER2
- fprintf(stderr, "\nSIGMET version 2 raw product file.\n");
+ fprintf(stderr, "\nSIGMET version 2 raw product file.\n");
#else
- fprintf(stderr, "\nSIGMET version 1 raw product file.\n");
+ fprintf(stderr, "\nSIGMET version 1 raw product file.\n");
#endif
- fprintf(stderr, "Date: %2.2d/%2.2d/%4.4d %2.2d:%2.2d:%f\n",
- radar->h.month, radar->h.day, radar->h.year,
- radar->h.hour, radar->h.minute, radar->h.sec);
- fprintf(stderr, "Name: ");
- for (i=0; i<sizeof(radar->h.name); i++)
- fprintf(stderr, "%c", radar->h.name[i]);
- fprintf(stderr, "\n");
- fprintf(stderr, "Lat/lon (%d %d' %d'', %d %d' %d'')\n",
- radar->h.latd, radar->h.latm, radar->h.lats,
- radar->h.lond, radar->h.lonm, radar->h.lons);
+ fprintf(stderr, "Date: %2.2d/%2.2d/%4.4d %2.2d:%2.2d:%f\n",
+ radar->h.month, radar->h.day, radar->h.year,
+ radar->h.hour, radar->h.minute, radar->h.sec);
+ fprintf(stderr, "Name: ");
+ for (i=0; i<sizeof(radar->h.name); i++)
+ fprintf(stderr, "%c", radar->h.name[i]);
+ fprintf(stderr, "\n");
+ fprintf(stderr, "Lat/lon (%d %d' %d'', %d %d' %d'')\n",
+ radar->h.latd, radar->h.latm, radar->h.lats,
+ radar->h.lond, radar->h.lonm, radar->h.lons);
}
/** Converting data **/
if (radar_verbose_flag) fprintf(stderr, "Expecting %d sweeps.\n", numsweep);
for(i = 0; i < numsweep; i++)
- {
- nsig_sweep = nsig_read_sweep(fp, prod_file);
- if (nsig_sweep == NULL) { /* EOF possibility */
- if (feof(fp)) break;
- else continue;
- }
- if (rsl_qsweep != NULL) {
- if (i > rsl_qsweep_max) break;
- if (rsl_qsweep[i] == 0) continue;
- }
- if (radar_verbose_flag)
- fprintf(stderr, "Read sweep # %d\n", i);
- /* The whole sweep is 'nsig_sweep' ... pretty slick.
+ {
+ nsig_sweep = nsig_read_sweep(fp, prod_file);
+ if (nsig_sweep == NULL) { /* EOF possibility */
+ if (feof(fp)) break;
+ else continue;
+ }
+ if (rsl_qsweep != NULL) {
+ if (i > rsl_qsweep_max) break;
+ if (rsl_qsweep[i] == 0) continue;
+ }
+ if (radar_verbose_flag)
+ fprintf(stderr, "Read sweep # %d\n", i);
+ /* The whole sweep is 'nsig_sweep' ... pretty slick.
*
* nsig_sweep[itype] -- [0..nparams], if non-null.
*/
- for (itype=0; itype<nparams; itype++) {
- if (nsig_sweep[itype] == NULL) continue;
-
- /** Reading date and time **/
- sweep_month = NSIG_I2(nsig_sweep[itype]->idh.time.month);
- sweep_year = NSIG_I2(nsig_sweep[itype]->idh.time.year);
- sweep_day = NSIG_I2(nsig_sweep[itype]->idh.time.day);
- sweep_sec = NSIG_I4(nsig_sweep[itype]->idh.time.sec);
+ for (itype=0; itype<nparams; itype++) {
+ if (nsig_sweep[itype] == NULL) continue;
+
+ /** Reading date and time **/
+ sweep_month = NSIG_I2(nsig_sweep[itype]->idh.time.month);
+ sweep_year = NSIG_I2(nsig_sweep[itype]->idh.time.year);
+ sweep_day = NSIG_I2(nsig_sweep[itype]->idh.time.day);
+ sweep_sec = NSIG_I4(nsig_sweep[itype]->idh.time.sec);
#ifdef NSIG_VER2
- msec = NSIG_I2(nsig_sweep[itype]->idh.time.msec);
- /* printf("....... msec == %d\n", msec); */
+ msec = NSIG_I2(nsig_sweep[itype]->idh.time.msec);
+ /* printf("....... msec == %d\n", msec); */
#endif
- /* converting seconds since mid to time of day */
- tmp = sweep_sec/3600.0;
- sweep_hour = (int)tmp;
- tmp = (tmp - sweep_hour) * 60.0;
- sweep_minute = (int)tmp;
- sweep_second = sweep_sec - (sweep_hour*3600 + sweep_minute*60);
-
- num_rays = NSIG_I2(nsig_sweep[itype]->idh.num_rays_exp);
-
- data_type = NSIG_I2(nsig_sweep[itype]->idh.data_type);
-
- ifield = 0;
- switch (data_type) {
- case NSIG_DTB_EXH:
- ifield = -1;
+ /* converting seconds since mid to time of day */
+ tmp = sweep_sec/3600.0;
+ sweep_hour = (int)tmp;
+ tmp = (tmp - sweep_hour) * 60.0;
+ sweep_minute = (int)tmp;
+ sweep_second = sweep_sec - (sweep_hour*3600 + sweep_minute*60);
+
+ num_rays = NSIG_I2(nsig_sweep[itype]->idh.num_rays_exp);
+
+ data_type = NSIG_I2(nsig_sweep[itype]->idh.data_type);
+
+ ifield = 0;
+ switch (data_type) {
+ case NSIG_DTB_EXH:
+ ifield = -1;
+ break;
+ case NSIG_DTB_UCR:
+ case NSIG_DTB_UCR2:
+ ifield = ZT_INDEX;
+ f = ZT_F;
+ invf = ZT_INVF;
+ break;
+ case NSIG_DTB_CR:
+ case NSIG_DTB_CR2:
+ ifield = DZ_INDEX;
+ f = DZ_F;
+ invf = DZ_INVF;
+ break;
+ case NSIG_DTB_VEL:
+ case NSIG_DTB_VEL2:
+ ifield = VR_INDEX;
+ f = VR_F;
+ invf = VR_INVF;
+ break;
+ case NSIG_DTB_WID:
+ case NSIG_DTB_WID2:
+ ifield = SW_INDEX;
+ f = SW_F;
+ invf = SW_INVF;
+ break;
+ case NSIG_DTB_ZDR:
+ case NSIG_DTB_ZDR2:
+ ifield = DR_INDEX;
+ f = DR_F;
+ invf = DR_INVF;
+ break;
+ case NSIG_DTB_KDP:
+ ifield = KD_INDEX;
+ f = KD_F;
+ invf = KD_INVF;
+ break;
+ case NSIG_DTB_PHIDP: /* SRB 990127 */
+ ifield = PH_INDEX;
+ f = PH_F;
+ invf = PH_INVF;
+ break;
+ case NSIG_DTB_RHOHV: /* SRB 000414 */
+ ifield = RH_INDEX;
+ f = RH_F;
+ invf = RH_INVF;
+ break;
+ case NSIG_DTB_VELC:
+ case NSIG_DTB_VELC2:
+ ifield = VC_INDEX;
+ f = VC_F;
+ invf = VC_INVF;
+ break;
+ case NSIG_DTB_KDP2:
+ ifield = KD_INDEX;
+ f = KD_F;
+ invf = KD_INVF;
+ break;
+ case NSIG_DTB_PHIDP2:
+ ifield = PH_INDEX;
+ f = PH_F;
+ invf = PH_INVF;
+ break;
+ case NSIG_DTB_RHOHV2:
+ ifield = RH_INDEX;
+ f = RH_F;
+ invf = RH_INVF;
+ break;
+ case NSIG_DTB_SQI:
+ case NSIG_DTB_SQI2:
+ ifield = SQ_INDEX;
+ f = SQ_F;
+ invf = SQ_INVF;
+ break;
+ case NSIG_DTB_HCLASS:
+ case NSIG_DTB_HCLASS2:
+ ifield = HC_INDEX;
+ f = HC_F;
+ invf = HC_INVF;
+ break;
+ default:
+ fprintf(stderr,"Unknown field type: %d Skipping it.\n", data_type);
+ continue;
+ }
+
+ if (radar_verbose_flag)
+ fprintf(stderr, " nsig_sweep[%d], data_type = %d, rays(expected) = %d, nrays(actual) = %d\n", itype, data_type, num_rays, NSIG_I2(nsig_sweep[itype]->idh.num_rays_act));
+
+ if (data_type != NSIG_DTB_EXH) {
+ if ((radar->v[ifield] == NULL)) {
+ if (rsl_qfield[ifield]) {
+ radar->v[ifield] = RSL_new_volume(numsweep);
+ radar->v[ifield]->h.f = f;
+ radar->v[ifield]->h.invf = invf;
+ } else {
+ /* Skip this field, because, the user does not want it. */
+ continue;
+ }
+ }
+ if (radar->v[ifield]->sweep[i] == NULL)
+ radar->v[ifield]->sweep[i] = RSL_new_sweep(num_rays);
+ }
+ else
+ continue; /* Skip the actual extended header processing.
+ * This is different than getting it, so that
+ * the information is available for the other
+ * fields when filling the RSL ray headers.
+ */
+
+ /** DATA conversion time **/
+ sweep = radar->v[ifield]->sweep[i];
+ sweep->h.f = f;
+ sweep->h.invf = invf;
+ sweep->h.sweep_num = i;
+ sweep->h.beam_width = beam_width;
+ sweep->h.vert_half_bw = vert_half_bw;
+ sweep->h.horz_half_bw = horz_half_bw;
+ elev = nsig_from_bang(nsig_sweep[itype]->idh.fix_ang);
+ sweep->h.elev = elev;
+
+ for(j = 0; j < num_rays; j++)
+ {
+ ray_p = nsig_sweep[itype]->ray[j];
+ if (ray_p == NULL) continue;
+ bin_num = NSIG_I2(ray_p->h.num_bins);
+
+ /* Load extended header information, if available.
+ * We need to pass the entire nsig_sweep and search for
+ * the extended header field (it may not be data_type==0).
+ */
+ get_extended_header_info(nsig_sweep, xh_size, j, nparams,
+ &msec, &azm, &elev,
+ &pitch, &roll, &heading,
+ &azm_rate, &elev_rate,
+ &pitch_rate, &roll_rate, &heading_rate,
+ &lat, &lon, &alt, &rvc,
+ &vel_east, &vel_north, &vel_up);
+
+
+ if (radar->v[ifield]->sweep[i]->ray[j] == NULL)
+ radar->v[ifield]->sweep[i]->ray[j] = RSL_new_ray(bin_num);
+ ray = radar->v[ifield]->sweep[i]->ray[j];
+ ray->h.f = f;
+ ray->h.invf = invf;
+ /** Ray is at nsig_sweep[itype].ray->... **/
+ /** Loading nsig data into data structure **/
+
+ ray->h.month = sweep_month;
+ ray->h.day = sweep_day;
+ ray->h.year = sweep_year; /* Year 2000 compliant. */
+ ray->h.hour = sweep_hour;
+ ray->h.minute = sweep_minute;
+ if (msec == 0) { /* No extended header */
+ ray->h.sec = NSIG_I2(ray_p->h.sec) + sweep_second;
+ elev = sweep->h.elev;
+ } else
+ ray->h.sec = sweep_second + msec/1000.0;
+
+ /* add time ... handles end of min,hour,month,year and century. */
+ if (ray->h.sec >= 60) /* Should I fix the time no matter what? */
+ RSL_fix_time(ray); /* Repair second overflow. */
+
+ ray->h.ray_num = j;
+ ray->h.elev_num = i;
+ ray->h.range_bin1 = (int)rng_first_bin;
+ ray->h.gate_size = (int)(bin_space+.5); /* Nearest int */
+ ray->h.vel_res = bin_space;
+ ray->h.sweep_rate = sweep_rate;
+ ray->h.prf = (int)prf;
+ if (prf != 0)
+ 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.azim_rate = azim_rate;
+ ray->h.pulse_count = (float)num_samples;
+ ray->h.pulse_width = pw;
+ ray->h.beam_width = beam_width;
+ ray->h.frequency = freq / 1000.0; /* GHz */
+ ray->h.wavelength = wave/100.0; /* meters */
+ ray->h.nyq_vel = max_vel; /* m/s */
+ if (elev == 0.) elev = sweep->h.elev;
+ ray->h.elev = elev;
+ /* Compute mean azimuth angle for ray. */
+ az1 = nsig_from_bang(ray_p->h.beg_azm);
+ az2 = nsig_from_bang(ray_p->h.end_azm);
+ /* printf("az1, %f, az2 %f\n", az1, az2); */
+ if(az1 > az2)
+ if((az1 - az2) > 180.0) az2 += 360.0;
+ else
+ ;
+ else
+ if((az2 - az1) > 180.0) az1 += 360.0;
+
+ az1 = (az1 + az2) / 2.0;
+ if (az1 > 360) az1 -= 360;
+ ray->h.azimuth = az1;
+
+ /* From the extended header information, we learn the following. */
+ ray->h.pitch = pitch;
+ ray->h.roll = roll;
+ ray->h.heading = heading;
+ ray->h.pitch_rate = pitch_rate;
+ ray->h.roll_rate = roll_rate;
+ ray->h.heading_rate = heading_rate;
+ ray->h.lat = lat;
+ ray->h.lon = lon;
+ ray->h.alt = alt;
+ ray->h.rvc = rvc;
+ ray->h.vel_east = vel_east;
+ ray->h.vel_north = vel_north;
+ ray->h.vel_up = vel_up;
+
+ /* printf("Processing sweep[%d]->ray[%d]: %d %f %f %f %f %f %f %f %f %d nbins=%d, bin1=%d gate=%d\n",
+ i, j, msec, ray->h.sec, ray->h.azimuth, ray->h.elev, ray->h.pitch, ray->h.roll, ray->h.heading, ray->h.lat, ray->h.lon, ray->h.alt, ray->h.nbins, ray->h.range_bin1, ray->h.gate_size);
+ */
+ /* TODO: ingest data header contains a value for bits-per-bin.
+ * This might be of use to allocate an array for ray->range with
+ * either 1-byte or 2-byte elements. Then there's no need for
+ * memmove() whenever we need 2 bytes.
+ */
+
+ if (data_type == NSIG_DTB_EXH) continue;
+ ray_data = 0;
+ for(k = 0; k < bin_num; k++) {
+ switch(data_type) {
+ case NSIG_DTB_UCR:
+ case NSIG_DTB_CR:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else ray_data = (float)((ray_p->range[k]-64.0)/2.0);
+ break;
+ /* Simplified the velocity conversion for NSIG_DTB_VEL, using
+ * formula from IRIS Programmer's Manual. BLK, Oct 9 2009.
+ */
+ case NSIG_DTB_VEL:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else ray_data = (float)((ray_p->range[k]-128.0)/127.0)*max_vel;
+ break;
+
+ case NSIG_DTB_WID:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else ray_data =(float)((ray_p->range[k])/256.0)*max_vel;
+ break;
+
+ case NSIG_DTB_ZDR:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else ray_data = (float)((ray_p->range[k]-128.0)/16.0);
+ break;
+
+ case NSIG_DTB_KDP:
+ if (ray_p->range[k] == 0 || ray_p->range[k] == 255 ||
+ rsl_kdp_wavelen == 0.0) {
+ ray_data = NSIG_NO_ECHO;
break;
- case NSIG_DTB_UCR:
- ifield = ZT_INDEX;
- f = ZT_F;
- invf = ZT_INVF;
- break;
- case NSIG_DTB_CR:
- ifield = DZ_INDEX;
- f = DZ_F;
- invf = DZ_INVF;
- break;
- case NSIG_DTB_VEL:
- ifield = VR_INDEX;
- f = VR_F;
- invf = VR_INVF;
- break;
- case NSIG_DTB_WID:
- ifield = SW_INDEX;
- f = SW_F;
- invf = SW_INVF;
- break;
- case NSIG_DTB_ZDR:
- ifield = DR_INDEX;
- f = DR_F;
- invf = DR_INVF;
- break;
- case NSIG_DTB_KDP:
- ifield = KD_INDEX;
- f = KD_F;
- invf = KD_INVF;
- break;
- case NSIG_DTB_PHIDP: /* SRB 990127 */
- ifield = PH_INDEX;
- f = PH_F;
- invf = PH_INVF;
- break;
- case NSIG_DTB_RHOHV: /* SRB 000414 */
- ifield = RH_INDEX;
- f = RH_F;
- invf = RH_INVF;
- break;
- case NSIG_DTB_SQI:
- ifield = SQ_INDEX;
- f = SQ_F;
- invf = SQ_INVF;
- break;
- default:
- fprintf(stderr,"Unknown field type: %d Skipping it.\n", data_type);
- continue;
- }
-
- if (radar_verbose_flag)
- fprintf(stderr, " nsig_sweep[%d], data_type = %d, rays(expected) = %d, nrays(actual) = %d\n", itype, data_type, num_rays, NSIG_I2(nsig_sweep[itype]->idh.num_rays_act));
-
- if (data_type != NSIG_DTB_EXH) {
- if ((radar->v[ifield] == NULL)) {
- if (rsl_qfield[ifield]) {
- radar->v[ifield] = RSL_new_volume(numsweep);
- radar->v[ifield]->h.f = f;
- radar->v[ifield]->h.invf = invf;
- } else {
- /* Skip this field, because, the user does not want it. */
- continue;
- }
- }
- if (radar->v[ifield]->sweep[i] == NULL)
- radar->v[ifield]->sweep[i] = RSL_new_sweep(num_rays);
- }
- else
- continue; /* Skip the actual extended header processing.
- * This is different than getting it, so that
- * the information is available for the other
- * fields when filling the RSL ray headers.
- */
-
- /** DATA conversion time **/
- sweep = radar->v[ifield]->sweep[i];
- sweep->h.f = f;
- sweep->h.invf = invf;
- sweep->h.sweep_num = i;
- sweep->h.beam_width = beam_width;
- sweep->h.vert_half_bw = vert_half_bw;
- sweep->h.horz_half_bw = horz_half_bw;
- elev = nsig_from_bang(nsig_sweep[itype]->idh.fix_ang);
- sweep->h.elev = elev;
-
- for(j = 0; j < num_rays; j++)
- {
- ray_p = nsig_sweep[itype]->ray[j];
- if (ray_p == NULL) continue;
- bin_num = NSIG_I2(ray_p->h.num_bins);
-
- /* Load extended header information, if available.
- * We need to pass the entire nsig_sweep and search for
- * the extended header field (it may not be data_type==0).
- */
- get_extended_header_info(nsig_sweep, xh_size, j, nparams,
- &msec, &azm, &elev,
- &pitch, &roll, &heading,
- &azm_rate, &elev_rate,
- &pitch_rate, &roll_rate, &heading_rate,
- &lat, &lon, &alt, &rvc,
- &vel_east, &vel_north, &vel_up);
-
-
- if (radar->v[ifield]->sweep[i]->ray[j] == NULL)
- radar->v[ifield]->sweep[i]->ray[j] = RSL_new_ray(bin_num);
- ray = radar->v[ifield]->sweep[i]->ray[j];
- ray->h.f = f;
- ray->h.invf = invf;
- /** Ray is at nsig_sweep[itype].ray->... **/
- /** Loading nsig data into data structure **/
-
- ray->h.month = sweep_month;
- ray->h.day = sweep_day;
- ray->h.year = sweep_year; /* Year 2000 compliant. */
- ray->h.hour = sweep_hour;
- ray->h.minute = sweep_minute;
- if (msec == 0) { /* No extended header */
- ray->h.sec = NSIG_I2(ray_p->h.sec) + sweep_second;
- elev = sweep->h.elev;
- } else
- ray->h.sec = sweep_second + msec/1000.0;
-
- /* add time ... handles end of min,hour,month,year and century. */
- if (ray->h.sec >= 60) /* Should I fix the time no matter what? */
- RSL_fix_time(ray); /* Repair second overflow. */
-
- ray->h.ray_num = j;
- ray->h.elev_num = i;
- ray->h.range_bin1 = (int)rng_first_bin;
- ray->h.gate_size = (int)(bin_space+.5); /* Nearest int */
- ray->h.vel_res = bin_space;
- ray->h.sweep_rate = sweep_rate;
- ray->h.prf = (int)prf;
- if (prf != 0)
- 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.azim_rate = azim_rate;
- ray->h.pulse_count = (float)num_samples;
- ray->h.pulse_width = pw;
- ray->h.beam_width = beam_width;
- ray->h.frequency = freq / 1000.0; /* GHz */
- ray->h.wavelength = wave/100.0; /* meters */
- ray->h.nyq_vel = max_vel; /* m/s */
- if (elev == 0.) elev = sweep->h.elev;
- ray->h.elev = elev;
- /* Compute mean azimuth angle for ray. */
- az1 = nsig_from_bang(ray_p->h.beg_azm);
- az2 = nsig_from_bang(ray_p->h.end_azm);
- /* printf("az1, %f, az2 %f\n", az1, az2); */
- if(az1 > az2)
- if((az1 - az2) > 180.0) az2 += 360.0;
- else
- ;
- else
- if((az2 - az1) > 180.0) az1 += 360.0;
-
- az1 = (az1 + az2) / 2.0;
- if (az1 > 360) az1 -= 360;
- ray->h.azimuth = az1;
-
- /* From the extended header information, we learn the following. */
- ray->h.pitch = pitch;
- ray->h.roll = roll;
- ray->h.heading = heading;
- ray->h.pitch_rate = pitch_rate;
- ray->h.roll_rate = roll_rate;
- ray->h.heading_rate = heading_rate;
- ray->h.lat = lat;
- ray->h.lon = lon;
- ray->h.alt = alt;
- ray->h.rvc = rvc;
- ray->h.vel_east = vel_east;
- ray->h.vel_north = vel_north;
- ray->h.vel_up = vel_up;
-
- /* printf("Processing sweep[%d]->ray[%d]: %d %f %f %f %f %f %f %f %f %d nbins=%d, bin1=%d gate=%d\n",
- i, j, msec, ray->h.sec, ray->h.azimuth, ray->h.elev, ray->h.pitch, ray->h.roll, ray->h.heading, ray->h.lat, ray->h.lon, ray->h.alt, ray->h.nbins, ray->h.range_bin1, ray->h.gate_size);
- */
- if (data_type == NSIG_DTB_EXH) continue;
- ray_data = 0;
- for(k = 0; k < bin_num; k++) {
- switch(data_type) {
- case NSIG_DTB_UCR:
- case NSIG_DTB_CR:
- if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
- else ray_data = (float)((ray_p->range[k]-64.0)/2.0);
- break;
-
- case NSIG_DTB_VEL:
- if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
- else ray_data = (float)((ray_p->range[k]*max_vel/127.0)+
- max_vel*(1.0-255.0/127.0));
- break;
-
- case NSIG_DTB_WID:
- if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
- else ray_data =(float)((ray_p->range[k])/256.0)*max_vel;
- break;
-
- case NSIG_DTB_ZDR:
- if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
- else ray_data = (float)((ray_p->range[k]-128.0)/16.0);
- break;
-
- /*
- * Special optimization note:
- * For KDP, PHIDP, RHOHV we skip the float conversion,
- * and carry the native sigmet data values into RSL storage.
- */
- case NSIG_DTB_KDP:
- ray_data = ray_p->range[k];
- /* F_OFFSET *must* match the definition in volume.c */
-#define F_OFFSET 4
- if (ray_data == 0 || ray_data == 255) ray_data = NSIG_NO_ECHO;
- else ray_data += F_OFFSET;
- break;
-
- case NSIG_DTB_RHOHV:
- case NSIG_DTB_PHIDP:
- if (ray_p->range[k] <= 0 || ray_p->range[k] >= 255)
- ray_data = NSIG_NO_ECHO;
- else
- ray_data = ray_p->range[k];
- break;
- case NSIG_DTB_SQI:
- if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
- else ray_data =
- (float)sqrt((ray_p->range[k]-1.0)/253.0);
- }
-
- if (ray_data == NSIG_NO_ECHO)
- ray->range[k] = ray->h.invf(BADVAL);
- else
- if ( (data_type == NSIG_DTB_KDP) ||
- (data_type == NSIG_DTB_RHOHV) ||
- (data_type == NSIG_DTB_PHIDP) )
- ray->range[k] = (Range)ray_data;
- else
- ray->range[k] = ray->h.invf(ray_data);
-
- /*
- if (data_type == NSIG_DTB_KDP)
- printf("v[%d]->sweep[%d]->ray[%d]->range[%d] = %f, %d, %f\n",
- ifield, i, j, k, ray->h.f(ray->range[k]),
- (int)ray_p->range[k], ray_data);
- */
- }
}
- }
+ if (ray_p->range[k] < 128)
+ ray_data = (-0.25 *
+ pow((double)600.0,(double)((127-ray_p->range[k])/126.0))) /
+ rsl_kdp_wavelen;
+ else if (ray_p->range[k] > 128)
+ ray_data = (0.25 *
+ pow((double)600.0,(double)((ray_p->range[k]-129)/126.0))) /
+ rsl_kdp_wavelen;
+ else
+ ray_data = 0.0;
+ break;
+
+ case NSIG_DTB_PHIDP:
+ if (ray_p->range[k] == 0 || ray_p->range[k] == 255)
+ ray_data = NSIG_NO_ECHO;
+ else
+ ray_data = 180.0*((ray_p->range[k]-1.0)/254.0);
+ break;
+
+ case NSIG_DTB_RHOHV:
+ if (ray_p->range[k] == 0 || ray_p->range[k] == 255)
+ ray_data = NSIG_NO_ECHO;
+ else
+ ray_data = sqrt((double)((ray_p->range[k]-1.0)/253.0));
+ break;
+
+ case NSIG_DTB_HCLASS:
+ if (ray_p->range[k] == 0 || ray_p->range[k] == 255)
+ ray_data = NSIG_NO_ECHO;
+ else
+ ray_data = ray_p->range[k];
+ break;
+
+ case NSIG_DTB_SQI:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else ray_data = (float)sqrt((ray_p->range[k]-1.0)/253.0);
+ break;
+
+ case NSIG_DTB_VELC:
+ if (ray_p->range[k] == 0) ray_data = NSIG_NO_ECHO;
+ else {
+ incr=75./127.; /* (+|- 75m/s) / 254 values */
+ ray_data = (float)(ray_p->range[k]-128)*incr;
+ }
+ break;
+
+ case NSIG_DTB_UCR2:
+ case NSIG_DTB_CR2:
+ case NSIG_DTB_VEL2:
+ 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)
+ ray_data = NSIG_NO_ECHO2;
+ else ray_data = (float)(nsig_u2byte-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)
+ ray_data = NSIG_NO_ECHO2;
+ else ray_data = (float)nsig_u2byte/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)
+ ray_data = NSIG_NO_ECHO;
+ else
+ ray_data = 360.*(nsig_u2byte-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)
+ ray_data = NSIG_NO_ECHO2;
+ else ray_data = (float)(nsig_u2byte-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)
+ ray_data = NSIG_NO_ECHO2;
+ else
+ ray_data = nsig_u2byte;
+ }
+
+ if (ray_data == NSIG_NO_ECHO || ray_data == NSIG_NO_ECHO2)
+ ray->range[k] = ray->h.invf(BADVAL);
+ else
+ ray->range[k] = ray->h.invf(ray_data);
+
+ /*
+ if (data_type == NSIG_DTB_KDP)
+ printf("v[%d]->sweep[%d]->ray[%d]->range[%d] = %f, %d, %f\n",
+ ifield, i, j, k, ray->h.f(ray->range[k]),
+ (int)ray_p->range[k], ray_data);
+ */
+ }
+ }
+ }
nsig_free_sweep(nsig_sweep);
- }
+ }
/* Do not reset radar->h.nvolumes. It is already set properly. */
if (radar_verbose_flag)
- fprintf(stderr, "Max index of radar->v[0..%d]\n", radar->h.nvolumes);
+ fprintf(stderr, "Max index of radar->v[0..%d]\n", radar->h.nvolumes);
/** close nsig file **/
#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;
- for (i=0; i<radar->h.nvolumes; i++)
- radar->v[i] = RSL_prune_volume(radar->v[i]);
+ 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);
return radar;
}
#define F_DZ_RANGE_OFFSET 32
#endif
+/* #define F_VR_OFFSET 63.5 */
+#define F_VR_OFFSET 127.0
+#define F_DR_OFFSET 12.0
+
/* IMPORTANT: This is the offset from reserved values. This
* number must be exactly (or >=) the number of
* reserved values in XX_F and XX_INVF.
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 - 63.5); /* Default wsr88d coding. */
+ val = (((float)x-F_OFFSET)/F_FACTOR - F_VR_OFFSET); /* Default wsr88d coding. */
/* fprintf(stderr, "x=%d, val=%f\n", x, val); */
return val;
}
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 - 12.0);
+ val = (((float)x-F_OFFSET)/F_DR_FACTOR - F_DR_OFFSET);
return val;
}
if (x == 0) return BADVAL;
return BADVAL;
}
+float HC_F(Range x) { /* HydroClass (Sigmet) */
+ if (x == 0) return BADVAL;
+ return (float)x;
+}
+
/****************************
Sigmet RhoHV : one_byte values
> RohHV = sqrt((N-1)/253)
*******************************/
float RH_F(Range x) {
if (x == 0) return BADVAL;
- return (float)(sqrt((double)((x-1.0)/253.0)));
- }
+ /* return (float)(sqrt((double)((x-1.0)/253.0))); */
+ return (float)(x-1) / 65533.;
+}
/*****************************
Sigmet PhiDP : one_byte values
******************************/
float PH_F(Range x) {
if (x == 0) return BADVAL;
- return (float)(180.0*((x-1.0)/254.0));
+ /*return (float)(180.0*((x-1.0)/254.0));*/
+ return (360.*(x-1.))/65534.;
}
+/* TODO: Should this be 5. cm. instead of 0.5? Or maybe 10. cm.? */
float rsl_kdp_wavelen = 0.5; /* Default radar wavelen = .5 cm. See
* nsig_to_radar.c for the code that sets this.
- */
+ */
+/* KD_F for 1 or 2 byte. */
float KD_F(Range x)
{
+/****** Commented-out code for 1-byte Sigmet native data format.
if (x >= F_OFFSET) {
- x -= F_OFFSET;
- if (rsl_kdp_wavelen == 0.0) return BADVAL;
- if (x < 128)
- return (float)(
- -0.25 * pow((double)600.0,(double)((127-x)/126.0))
- )/rsl_kdp_wavelen;
- else if (x > 128)
- return (float)(
- 0.25 * pow((double)600.0,(double)((x-129)/126.0))
- )/rsl_kdp_wavelen;
- else
- return 0.0;
+ x -= F_OFFSET;
+ if (rsl_kdp_wavelen == 0.0) return BADVAL;
+ if (x < 128)
+ return (float)(
+ -0.25 * pow((double)600.0,(double)((127-x)/126.0))
+ )/rsl_kdp_wavelen;
+ else if (x > 128)
+ return (float)(
+ 0.25 * pow((double)600.0,(double)((x-129)/126.0))
+ )/rsl_kdp_wavelen;
+ else
+ return 0.0;
}
- if (x == 0) return BADVAL;
if (x == 1) return RFVAL;
if (x == 2) return APFLAG;
if (x == 3) return NOECHO;
return BADVAL;
+******/
+
+ if (x == 0) return BADVAL;
+ return (x-32768.)/100.;
+}
+
+float NP_F(Range x) { /* Normalized Coherent Power (DORADE) */
+ if (x == 0) return BADVAL;
+ return (float)x / 100.;
}
+/* Signal Quality Index */
float SQ_F(Range x) {
- if (x >= F_OFFSET) return (float)((x - F_OFFSET)/10000.);
if (x == 0) return BADVAL;
- if (x == 1) return RFVAL;
- if (x == 2) return APFLAG;
- if (x == 3) return NOECHO;
- return BADVAL; /* Can't get here, but quiets the compiler. */
+ return (float)(x-1) / 65533.;
}
float TI_F(Range x)
float AH_F(Range x) { return DZ_F(x); }
float CV_F(Range x) { return DZ_F(x); }
float AV_F(Range x) { return DZ_F(x); }
+float VS_F(Range x) { return VR_F(x); }
+float VL_F(Range x) { return VR_F(x); }
+float VG_F(Range x) { return VR_F(x); }
+float VT_F(Range x) { return VR_F(x); }
+float VC_F(Range x) { return VR_F(x); }
if (x == RFVAL) return (Range)1;
if (x == APFLAG) return (Range)2;
if (x == NOECHO) return (Range)3;
+ if (x < -F_DZ_RANGE_OFFSET) return (Range)0;
return (Range)(F_FACTOR*(x+F_DZ_RANGE_OFFSET)+.5 + F_OFFSET); /* Default wsr88d. */
}
if (x == RFVAL) return (Range)1;
if (x == APFLAG) return (Range)2;
if (x == NOECHO) return (Range)3;
- return (Range)(F_FACTOR*(x+63.5)+.5 + F_OFFSET); /* Default wsr88d coding. */
+ if (x < -F_VR_OFFSET) return (Range)0;
+ return (Range)(F_FACTOR*(x+F_VR_OFFSET)+.5 + F_OFFSET); /* Default wsr88d coding. */
}
Range DR_INVF(float x) /* Differential reflectivity */
if (x == RFVAL) return (Range)1;
if (x == APFLAG) return (Range)2;
if (x == NOECHO) return (Range)3;
- return (Range)(F_DR_FACTOR*(x + 12.0) + F_OFFSET + 0.5);
+ if (x < -F_DR_OFFSET) return (Range)0;
+ return (Range)(F_DR_FACTOR*(x + F_DR_OFFSET) + F_OFFSET + 0.5);
+}
+
+Range HC_INVF(float x) /* HydroClass (Sigmet) */
+{
+ if (x == BADVAL) return (Range)0;
+ return (Range)(x + 0.5); /* Round */
}
Range LR_INVF(float x) /* MCTEX */
> 253 0.9980
> 254 1.0000
****************************/
+/* RH_INVF for 1 or 2 byte data. */
Range RH_INVF(float x) {
if (x == BADVAL) return (Range)0;
- return (Range)(x * x * 253.0 + 1.0 + 0.5);
+ /* return (Range)(x * x * 253.0 + 1.0 + 0.5); */
+ return (Range)(x * 65533. + 1. +.5);
}
/******************************
*******************************/
Range PH_INVF(float x) {
if (x == BADVAL) return (Range)0;
- return (Range)((x / 180.0) * 254.0 + 1.0 + 0.5);
+ /* return (Range)((x / 180.0) * 254.0 + 1.0 + 0.5); */
+ return (Range)(x*65534./360. + 1.0 + 0.5);
}
+
+/* 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.
if (x == RFVAL) return (Range)1;
if (x == APFLAG) return (Range)2;
if (x == NOECHO) return (Range)3;
x = 128;
}
x += F_OFFSET;
- return (Range) x;
+******/
+ return (Range)(x * 100. + 32768. + 0.5);
}
-Range SQ_INVF(float x) /* Signal Quality Index */
+/* Signal Quality Index */
+Range SQ_INVF(float x)
{
if (x == BADVAL) return (Range)0;
- if (x == RFVAL) return (Range)1;
- if (x == APFLAG) return (Range)2;
- if (x == NOECHO) return (Range)3;
- return (Range)(x * 10000. + F_OFFSET);
+ return (Range)(x * 65533. + 1. +.5);
+}
+
+
+Range NP_INVF(float x) /* Normalized Coherent Power (DORADE) */
+{
+ if (x == BADVAL) return (0);
+ return (Range)(x * 100.);
}
Range TI_INVF(float x) /* MCTEX */
Range AH_INVF(float x) { return DZ_INVF(x); }
Range CV_INVF(float x) { return DZ_INVF(x); }
Range AV_INVF(float x) { return DZ_INVF(x); }
+Range VS_INVF(float x) { return VR_INVF(x); }
+Range VL_INVF(float x) { return VR_INVF(x); }
+Range VG_INVF(float x) { return VR_INVF(x); }
+Range VT_INVF(float x) { return VR_INVF(x); }
+Range VC_INVF(float x) { return VR_INVF(x); }
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, 1, 1, 1,
+ 1, 1
};