/*
* ##################################################################
* ##################################################################
* ###### ######
* ###### 88D2ARPS ######
* ###### ######
* ###### Developed by ######
* ###### Center for Analysis and Prediction of Storms ######
* ###### University of Oklahoma ######
* ###### ######
* ##################################################################
* ##################################################################
*
* 88d2arps.c
*
*#######################################################################
*
* PURPOSE:
*
* Reads WSR 88D radar data from circular buffer using NSSL
* Realtime/Archive II radar reading routines and passes the data
* to ARPS remapping routines.
*
* Usage:
* 88d2arps Kxxx reads data from live data
* 88d2arps Kxxx -f /dev/rmt0 reads data from tape drive rmt0
* 88d2arps_fake Kxxx test code using fakio.c analytic data
*
* Advertised to work for live WSR-88D radar data as well
* as archive II taped data.
*
* Requires NSSL libraries: -la2rt for real-time
* -la2tp and -ltpio for tape playback
*
*#######################################################################
*
* AUTHOR: Keith Brewster, CAPS June 1995
* June 1995, Based on CAPS remap.c written for FSL LAPS
*
* MODIFICATIONS:
* 10/15/1996 (K. Brewster)
* Updated for changes to remaptilt
*
* 12/11/1997 (K. Brewster)
* Ranamed from remaparps and set-up for arps4.3.1
*
* 05/25/2000 (K. Brewster)
* Modified for changes to the setlookup argument list.
* Added file read capability to tape and realtime.
*
* 08/31/2000 (Gene Bassett)
* Converted to use array allocation.
*
* 03/19/2001 (K. Brewster)
* Merged updated versions from KB and GB.
*
* 05/10/2001 (K. Brewster)
* Added new interpolation schemes.
*
* 05/16/2003 (K. Brewster)
* Added writing of min,max range and elevation to datafile header.
*
*
*#######################################################################
*
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
/* When VERBOSE is defined, the program generates verbose output
* that is useful for testing but annoying for operations */
#define VERBOSE
#undef VERBOSE
/* When UNFOLD_DBG is defined, the program generates files containing
* raw az-ran data for checking the unfolding algorithms */
#define UNFOLD_DBG
#undef UNFOLD_DBG
/* Files may be compressed using the program identified by REMAP_COMPRESS.
* Currently supported REMAP_COMPRESS options are: gzip, compress.
*
* When the proper environment variables are defined
* the program will try to send the remapped files
* to directory REMAP_DEST on machine REMAP_REMOTE owned by REMAP_USER,
* using transmission mode REMAP_XMIT.
*
* currently supported REMAP_XMIT methods are
* rcp, ftp, mv
*
* The REMAP_ variables are environment variables of the parent shell.
*
* Sample:
* setenv REMAP_COMPRESS gzip
* setenv REMAP_XMIT rcp
* setenv REMAP_USER kbrews
* setenv REMAP_REMOTE cirrus.gcn.uoknor.edu
* setenv REMAP_DEST /scratch/cirrus/kbrews/arps40
*/
/* Note:
* The radar name is determined from the environment variable
* RADARNAME
*
* Sample:
* setenv RADARNAME KTLX
*
* A file containing the lat,lon and elevation of the radar sites
* should be in the directory of remap execution.
* This file is named radarinfo.dat.
*/
/* The following defines the circular buffer to use
* "RAW" is RAW
* "E" is unfolded
*/
#define RT_BUFFER "RAW"
#define N_AZIM_TILT 400
#define N_REF_GATES 460
#define N_VEL_GATES 920
#define MAX_ELEV 20
#define MAX_SORT 3000
/* The following are for comparison to return
values from the radar data access libraries */
#define GOOD_STATUS 1
#define BAD_STATUS 0
#define MAX_BAD_STAT 1000
#define BLANK ""
#define TRUE 1
#define FALSE 0
#define IORDER_REF 2
#define IORDER_VEL 2
#define REF_CHEK -20.0
#define VEL_CHEK -90.0
#define REF_MISS -999.0
#define VEL_MISS -999.0
#define REF_DAZLIM 360.0
#define VEL_DAZLIM 30.0
#define REF_MEDLIM 20.0
#define VEL_MEDLIM 15.0
#define DELTA_AZIM 1.0
#define RANGE_MIN 10000.0
#define RANGE_MAX 230000.0
#define PLAYBACK 0
#define REALTIME 1
#define MODE_SET 0
#define NZ_SND 200
#define NTS 1
#define APWINSZ_RAD 7000.0
#define APWINSZ_AZIM 10.0
#define UNFBKGOPT 1
#define UNFSHROPT 1
#define RFRACTOPT 1
#define RANGE_VAD 30000.
#define RANGE_WIDTH 1000.
#define MIN_COUNT 20
#ifdef UNDERSCORE
#define ctim2abss ctim2abss_
#define initpara initpara_
#define mkradfnm mkradfnm_
#define mkvadfnm mkvadfnm_
#define radcoord radcoord_
#define rmpinit rmpinit_
#define despekl despekl_
#define unfnqc unfnqc_
#define volbuild volbuild_
#define vadvol vadvol_
#define remapvol remapvol_
#define remap2d remap2d_
#define wrtrad88 wrtrad88_
#define wtradcol wtradcol_
#define wtvadprf wtvadprf_
#define wrtvel wrtvel_
#define set_lbcopt set_lbcopt_
#define get_dims_from_data get_dims_from_data_
#define initgrdvar initgrdvar_
#define extenvprf extenvprf_
#define uv2vr uv2vr_
#define quadfit quadfit_
#define apdetect apdetect_
#define wrttilts wrttilts_
#endif
#define allocate(s) (s= (char *) malloc(80))
/* compr_file_sub controls
compression of the remapped data file */
void compr_file_sub(file_ptr,compr_ptr,i_status)
char *file_ptr,*compr_ptr;
int i_status;
{
char *string1;
FILE *comp_file;
int comp_stat;
i_status=0;
allocate(string1);
if(strcmp(compr_ptr,"gzip") == 0)
{
strcpy(string1,"gzip -f ");
strcat(string1,file_ptr);
comp_file= popen(string1,"r");
if(comp_file == (FILE *) NULL) {
printf("\n Error occurred starting compression..\n");
i_status=-1;
}
else {
comp_stat=pclose(comp_file);
if( comp_stat == 0) {
strcat(file_ptr,".gz");
printf("\n File Successfully Compressed\n");
}
else {
printf("\n Error occured while compressing %s\n",file_ptr);
printf("\n Returned status = %i\n",comp_stat);
}
i_status=comp_stat;
}
}
else if(strcmp(compr_ptr,"compress") == 0)
{
strcpy(string1,"compress -f ");
strcat(string1,file_ptr);
comp_file= popen(string1,"r");
strcat(file_ptr,".Z");
fclose(comp_file);
}
else
{
printf(" %s Compression Unsupported\n",compr_ptr);
printf(" Set REMAP_COMPRESS to gzip or compress\n");
i_status=-2;
}
}
/* send_file_sub controls
transmission of the remapped data file */
void send_file_sub(file_ptr,xmit_ptr,
user_ptr,host_ptr,dest_ptr,i_status)
char *file_ptr,*xmit_ptr,*user_ptr,*host_ptr,*dest_ptr;
int i_status;
{
char *string1;
FILE *p1;
int xmit_stat;
int slash;
char buf[200];
i_status=0;
allocate(string1);
/* Strip any directory info from the file name so all
* directory info in the destination file name is from
* dest_ptr, environment variable REMAP_DEST.
* The source file name is unchanged */
slash=47;
string1=(char *)strrchr(file_ptr,slash);
if (strncmp(xmit_ptr,"mv",2) == 0) {
/* Create mv command */
if (string1 == NULL)
sprintf(buf,"mv %s %s/%s",file_ptr,
dest_ptr,file_ptr);
else
sprintf(buf,"mv %s %s%s",file_ptr,
dest_ptr,string1);
printf("\n Renaming file.....\n");
}
else if (strncmp(xmit_ptr,"rcp",3) == 0) {
/* Create the rcp command */
if (string1 == NULL)
sprintf(buf,"rcp %s %s@%s:%s/%s",file_ptr,user_ptr,host_ptr,
dest_ptr,file_ptr);
else
sprintf(buf,"rcp %s %s@%s:%s%s",file_ptr,user_ptr,host_ptr,
dest_ptr,string1);
printf("\n Copying file to remote host.....\n");
}
else if (xmit_ptr != NULL) {
/* Create the custom user command
For this option, a command line is created from the
REMAP environment variables and the long and short filenames
The REMAP_XMIT variable (represented here by xmit_ptr) is the
name of a user-programmed script or executable that accepts
the arguments as passed here. These arguments should give
the script or program enough info to work with.
The command line created consists of:
REMAP_XMIT long_filename short_filename REMAP_USER REMAP_REMOTE REMAP_DEST
*/
if (string1 == NULL)
sprintf(buf,"%s %s %s %s %s %s",xmit_ptr,file_ptr,file_ptr,
user_ptr,host_ptr,dest_ptr);
else
sprintf(buf,"%s %s %s %s %s %s",xmit_ptr,file_ptr,string1,
user_ptr,host_ptr,dest_ptr);
printf("\n Issuing command to user-defined data management program...\n");
}
fprintf(stdout,buf,"%s","\n");
p1= popen(buf,"r");
if (p1 == (FILE *) NULL)
printf("\n Error occurred starting file transmission..\n");
else {
xmit_stat = pclose(p1);
if( xmit_stat == 0 )
printf("\n File Successfully Transmitted.\n");
else {
printf("\n Error occured during transmission %s\n",file_ptr);
printf("\n Returned status = %i\n",xmit_stat);
}
}
}
/******************************************************
* *
* *
* MAIN REMAPPING DRIVER *
* *
* *
*******************************************************/
main(argc,argv)
int argc;
char *argv[];
{
/* Switches for controlling remapping */
int i_last_scan,i_first_scan;
/* Variables used for data access and in remapping */
float *ref_ptr, *ref0_ptr;
float *vel_ptr, *vel0_ptr;
float *spw_ptr, *spw0_ptr;
float *unfv_ptr;
float *bkgv_ptr;
float *tmp1_ptr;
float *azi0_ptr;
float *elv0_ptr;
int *tim0_ptr;
int *bgate_ptr;
int *egate_ptr;
int krelv, *kntrelv;
int *kntrazm;
int *kntrgat;
int kvelv, *kntvelv;
int *kntvazm;
int *kntvgat;
int ref_index,vel_index,spw_index,io_stat,i_status;
int ref_ok,vel_ok;
int gsp_ref,gsp_vel; /* gate spacing (m), ref, vel */
int rfrst_ref,rfrst_vel; /* range to 1st gate (m), ref, vel */
int maxgate, maxazim, maxelv, maxsort;
int n_azim, i_scan, i_tilt, n_ref_gates, n_vel_gates;
float dazim,rngmin,rngmax;
float refmiss,velmiss,refchek,velchek;
float refdazl,veldazl,refmedl,velmedl;
/* File transfer variables */
char *compr_name;
char *xmit_method, *remote_user, *remote_host, *dest_dir;
static char *blnk = BLANK;
/* Radar Location variables */
int i_lat,i_lon,i_alt;
float radar_lat, *rlat_ptr;
float radar_lon, *rlon_ptr;
float radar_alt, *ralt_ptr;
float radarx,radary,envavgr;
float refelvmin,refelvmax,refrngmax,refrngmin;
/* Diagnostic write variables */
/* Misc Local variables */
char *sw;
char source[256], *src_ptr;
char refid[6] = "reflct", *refid_ptr;
char refname[20] = "Reflectivity", *refnam_ptr;
char ref2did[6] = "refl2d", *ref2did_ptr;
char ref2dname[20] = "Reflect-2d", *ref2dnam_ptr;
char refunits[20] = "dBZ", *refunit_ptr;
char velid[6] = "radvel", *velid_ptr;
char velname[20] = "Radial Velocity",*velnam_ptr;
char vel2did[6] = "radv2d", *vel2did_ptr;
char vel2dname[20] = "Radial Vel 2D",*vel2dnam_ptr;
char velunits[20] = "m/s", *velunit_ptr;
char dir_name[80], *dir_ptr;
char radar_name[5], *rname_ptr;
char vad_fname[100], *vadfnm_ptr;
char full_fname[100], *fnm_ptr;
int iyr,iyear,iday,imon,ihour,imin,isec,itime; /* time variables */
int iyr1,iday1,imon1,ihour1,imin1,isec1,itime1; /* time variables */
int *timevol;
int initial_ray; /* flag for first ray in tilt*/
int initvol_ray; /* flag for first ray in volume*/
int alls_well, alloc_vol, read_next, knt_bad_stat;
int i_angle, past_angle, ng_ref;
int past_scan, past_tilt;
int len_dir, len_fname, len_vfname;
int lbc_one, settim, setoff, seton;
int compr_on, xmit_on, qc_on, write_and_exit;
int wrttilts_on;
int velproc, refdump, veldump, ref2d, vel2d, vad;
int dmpfmt, hdf4cmpr;
int i, i_mode, i_arg, istatus;
int i_vcp,i_src,irad_fmt;
int nx,ny,nz,nzsoil,nstyps;
int iordref,iordvel;
float v_nyq, rmax;
float *nyqvvol;
float *elvrvol,*azmrvol,*rngrvol,*refvol;
float *elvvvol,*azmvvol,*rngvvol,*velvol;
float *rxvol,*ryvol,*rzvol;
float *varsort, *elvmean;
int nt,exbcbufsz;
float *x,*y,*z,*zp,*xs,*ys,*zps,*hterain,*mapfct;
float *tem3d1,*tem3d2,*tem3d3,*tem2d,*tem3dsoil;
float *tem2dyz,*tem2dxz,*tem2dns,*tem4dsoilns;
float *u,*v,*qv,*uprt,*vprt,*ptprt,*pprt;
float *ubar,*vbar,*ptbar,*pbar,*qvbar;
float *rhostr,*j1,*j2,*j3;
float *trigs1,*trigs2;
float *wsave1,*wsave2,*vwork1,*vwork2;
float *stypfrct;
int *soiltyp,*vegtyp,*ifax;
float *qcumsrc;
float *prcrate,*exbcbuf;
float *zsnd,*usnd,*vsnd,*rfrsnd,*ktsnd;
float *gridvel,*gridref,*gridnyq,*gridtim;
float *tem1d1,*tem1d2,*tem1d3,*tem1d4,*tem1d5,*tem1d6;
float *vadhgt,*vaddir,*vadspd;
float winszrad,winszazim;
float dzsnd;
int nzsnd;
int bkgopt,shropt,rfropt;
float rngvad,rngwid;
int minknt;
#ifdef VERBOSE
int ng_vel; /* number of gates vel returned */
#endif
#ifdef UNFOLD_DBG
int p_angle,p_tilt;
float elev;
char varid[6];
#endif
/* Beginning of Executable Code */
/* Some initializations */
n_vel_gates = N_VEL_GATES;
n_ref_gates = N_REF_GATES;
if ( n_ref_gates > n_vel_gates )
maxgate = n_ref_gates;
else
maxgate = n_vel_gates;
maxazim = N_AZIM_TILT;
iordref = IORDER_REF;
iordvel = IORDER_VEL;
dazim = DELTA_AZIM;
rngmin = RANGE_MIN;
rngmax = RANGE_MAX;
envavgr = rngmax/1.4142;
refchek = REF_CHEK;
velchek = VEL_CHEK;
refmiss = REF_MISS;
velmiss = VEL_MISS;
refdazl = REF_DAZLIM;
veldazl = VEL_DAZLIM;
refmedl = REF_MEDLIM;
velmedl = VEL_MEDLIM;
nzsnd = NZ_SND;
rngvad = RANGE_VAD;
rngwid = RANGE_WIDTH;
minknt = MIN_COUNT;
bkgopt = UNFBKGOPT;
shropt = UNFSHROPT;
rfropt = RFRACTOPT;
nt = NTS;
winszrad = APWINSZ_RAD;
winszazim = APWINSZ_AZIM;
maxsort = MAX_SORT;
fnm_ptr = &full_fname[0];
vadfnm_ptr = &vad_fname[0];
kntrelv = &krelv;
kntvelv = &kvelv;
rlat_ptr=&radar_lat;
rlon_ptr=&radar_lon;
ralt_ptr=&radar_alt;
refid_ptr=refid;
refnam_ptr=refname;
ref2did_ptr=ref2did;
ref2dnam_ptr=ref2dname;
refunit_ptr=refunits;
velid_ptr=velid;
velnam_ptr=velname;
vel2did_ptr=vel2did;
vel2dnam_ptr=vel2dname;
velunit_ptr=velunits;
refelvmin=999.;
refelvmax=-999.;
refrngmax=0.;
refrngmin=999999.;
lbc_one = 1;
setoff = 0;
seton = 1;
irad_fmt = 1;
nt = 1;
/* Get Radar name variable */
dir_ptr=dir_name;
dir_ptr=getenv("REMAP_DIR");
if ( dir_ptr == NULL || strcmp(dir_ptr,blnk) == 0 )
dir_ptr = "./";
len_dir = (int) strlen(dir_ptr);
i_arg=0;
rname_ptr=radar_name;
if (argc > 1 ) {
sw = argv[1];
printf (" First argument: %s\n",sw);
if (index(sw,'-') != NULL ) {
rname_ptr=getenv("RADARNAME");
if ( rname_ptr == NULL || strcmp(rname_ptr,blnk) == 0 ) {
printf (" Radar name not first in command line and\n");
printf (" RADARNAME environmental variable not defined\n");
printf (" Usage %s Kxxx -f /dev/rmt/0n [-reffile] [-velfile]\n",
argv[0]);
exit (1);
} else {
printf (" Radar name obtained from env variable %s\n",rname_ptr);
}
}
else {
i_arg=1;
rname_ptr = sw;
printf (" Radar name from command line %s\n",rname_ptr);
}
}
if (rname_ptr == NULL )
{
printf (" Couldn't evaluate RADARNAME environment variable.\n");
printf (" Set the 4-character radar name using:\n");
printf (" setenv RADARNAME Kxxx\n");
printf (" before running 88d2arps.\n");
printf (" Or enter radar name as first argument on command line.\n");
exit(1);
}
set_radar_name(rname_ptr);
/* Get REMAP environment variables */
veldump = FALSE;
refdump = FALSE;
vel2d = FALSE;
ref2d = FALSE;
vad = FALSE;
velproc = TRUE;
compr_on = TRUE;
xmit_on = TRUE;
qc_on = TRUE;
wrttilts_on = FALSE;
allocate(compr_name);
compr_name=getenv("REMAP_COMPRESS");
allocate(xmit_method);
xmit_method=getenv("REMAP_XMIT");
allocate(remote_user);
remote_user=getenv("REMAP_USER");
allocate(remote_host);
remote_host=getenv("REMAP_REMOTE");
allocate(dest_dir);
dest_dir=getenv("REMAP_DEST");
if ( remote_user == NULL || strcmp(remote_user,blnk) == 0 ||
remote_host == NULL || strcmp(remote_host,blnk) == 0 ||
dest_dir == NULL || strcmp(dest_dir,blnk) == 0 ) {
xmit_on = FALSE;
printf ("\n\n\n REMAP environment variables\n");
printf (" which specify file transfer options not set.\n");
printf (" Will assume no file transferring desired.\n\n");
printf ("\n To enable data transfer set the following before 88d2arps\n");
printf (" setenv REMAP_USER user\n");
printf (" setenv REMAP_HOST host.domain.ext\n");
printf (" setenv REMAP_DEST /destination/directory\n\n");
}
if (xmit_method == NULL || strcmp(xmit_method,blnk) == 0 )
xmit_method = "rcp";
if (compr_name == NULL || strcmp(compr_name,blnk) == 0 )
compr_on = FALSE;
if ( xmit_on == TRUE ) {
printf (" REMAP_COMPRESS: %s\n",compr_name);
printf (" REMAP_XMIT: %s\n",xmit_method);
printf (" REMAP_USER: %s\n",remote_user);
printf (" REMAP_REMOTE: %s\n",remote_host);
printf (" REMAP_DEST: %s\n",dest_dir);
}
/* Read and process command line switches */
dmpfmt=1;
hdf4cmpr=0;
i_mode=MODE_SET;
if ( strstr(argv[0],"88d2arps_fake") != NULL ) i_mode=REALTIME;
src_ptr=source;
if ((argc-i_arg) < 3) {
src_ptr = RT_BUFFER;
if ( i_mode == REALTIME )
{
printf (" Reading real-time datastream from radar %s\n",rname_ptr);
printf (" Using circular buffer %s\n",source);
}
else
{
printf (" Remapper compiled with tape/file playback libraries\n");
printf (" Provide file or tape drive name,\n");
printf (" or use real-time version.\n");
printf (" Usage: 88d2arps_rt KTLX OR\n");
printf (" Usage: 88d2arps_a2 KTLX -f tape_device OR\n");
printf (" Usage: 88d2arps_a2 KTLX -diskf filename \n");
exit(1);
}
} else {
if ( i_mode == PLAYBACK ) {
while (i_arg < (argc-1)) {
i_arg++;
sw = argv[i_arg];
if (strcmp(sw,"-diskf") == 0) {
if (i_arg < argc ) {
i_arg++;
src_ptr = argv[i_arg];
printf(" Command line switch: %s, diskfile name %s\n",
sw,src_ptr);
printf (" Archive data is from radar %s\n",rname_ptr);
printf (" Reading file %s\n",src_ptr);
}
else {
printf (" Switch %s requires an argument, try again\n",
sw);
exit(2);
}
}
else if (strcmp(sw,"-f") == 0) {
if (i_arg < argc ) {
i_arg++;
src_ptr = argv[i_arg];
printf(" Command line switch: %s, drive name %s\n",
sw,src_ptr);
printf (" Archive data is from radar %s\n",rname_ptr);
printf (" Reading from tape drive %s\n",src_ptr);
}
else {
printf (" Switch %s requires an argument, try again\n",
sw);
exit(3);
}
}
else if (strcmp(sw,"-dir") == 0) {
if (i_arg < argc ) {
i_arg++;
dir_ptr = argv[i_arg];
len_dir = (int) strlen(dir_ptr);
printf(" Command line switch: %s, directory name %s\n",
sw,dir_ptr);
printf (" Writing remap data to directory: %s\n",dir_ptr);
}
else {
printf (" Switch %s requires an argument, try again\n",
sw);
exit(4);
}
}
else if (strcmp(sw,"-hdf") == 0 ) {
printf (" Using hdf format for output file.\n");
dmpfmt=2;
hdf4cmpr=0;
if( i_arg < argc ) {
i_arg++;
hdf4cmpr=atoi(argv[i_arg]);
if( hdf4cmpr < 0 ) hdf4cmpr=0;
if( hdf4cmpr > 7 ) hdf4cmpr=7;
}
printf (" Using hdf compression level: %i.\n",hdf4cmpr);
}
else if (strcmp(sw,"-binary") == 0 ) {
printf (" Using binary format for output file.\n");
dmpfmt=1;
hdf4cmpr=0;
}
else if (strcmp(sw,"-reffile") == 0 ) {
printf (" 3-D reflectivity file generated for arpsplt\n");
refdump = TRUE;
}
else if (strcmp(sw,"-velfile") == 0 ) {
printf (" 3-D velocity file generated for arpsplt\n");
veldump = TRUE;
}
else if (strcmp(sw,"-vad") == 0 ) {
printf (" VAD wind profile file generated for ADAS\n");
vad = TRUE;
}
else if (strcmp(sw,"-ref2d") == 0 ) {
printf (" 2-D low-level reflectivity file generated for arpsplt\n");
ref2d = TRUE;
}
else if (strcmp(sw,"-vel2d") == 0 ) {
printf (" 2-D low-level velocity file generated for arpsplt\n");
vel2d = TRUE;
}
else if (strcmp(sw,"-novel") == 0 ) {
printf (" Witholding processing of velocity\n");
velproc = FALSE;
}
else if (strcmp(sw,"-noqc") == 0 ) {
printf (" Witholding quality control processing\n");
qc_on = FALSE;
}
else if (strcmp(sw,"-wrttilts_on") == 0 ) {
printf (" Write out radar scan.\n");
wrttilts_on = TRUE;
}
else {
printf (" Unknown command line switch: %s, try again.\n",sw);
exit(1);
}
}
}
else {
printf (" Remapper compiled with real-time libraries\n");
printf (" Remove tape drive name, or use real-time version.\n");
printf (" Usage: 88d2arps_rt OR\n");
printf (" Usage: 8dd2arps_a2 -f tape_device OR\n");
printf (" Usage: 8dd2arps_a2 -diskf filename \n");
}
}
if ( veldump == TRUE ) velproc = TRUE;
if ( velproc == TRUE ) {
settim = setoff;
} else {
settim = seton;
}
printf (" Evaluated dir name as %s Length: %i\n",dir_ptr,len_dir);
/* call Archive II initialization routine */
radar_init(src_ptr);
i_alt=get_altitude();
i_lat=get_latitude();
i_lon=get_longitude();
if ( i_lat == 0 && i_lon == 0 ) {
printf (" Radar not found in radarinfo.dat file.\n");
printf (" Or radarinfo.dat missing in local directory.\n");
printf (" radarinfo.dat may be copied from ./data/adas\n");
exit(1);
}
radar_alt= (float) i_alt;
radar_lat= 0.00001 * (float) i_lat;
radar_lon= -0.00001 * (float) i_lon;
printf ("\n Radar altitude (m): %f\n",radar_alt);
printf (" Radar latitude (degrees): %f\n",radar_lat);
printf (" Radar longitude (degrees): %f\n\n",radar_lon);
/* call FORTRAN routine to set up ARPS variables */
initpara(&nx,&ny,&nz,&nzsoil,&nstyps);
printf (" Back from initpara\n");
printf (" nx: %d , ny: %d , nz: %d , nz soil: %d , nstyps: %d \n",
nx,ny,nz,nzsoil,nstyps);
/* Set the lbcopt to override input to save memory in initgrdvar */
set_lbcopt(&lbc_one);
exbcbufsz=1;
/* allocate memory for ARPS grids needed for environmental wind profile */
printf (" Allocating ARPS grid memory\n");
x = (float *)calloc((size_t)(nx),sizeof(float));
y = (float *)calloc((size_t)(ny),sizeof(float));
z = (float *)calloc((size_t)(nz),sizeof(float));
zp = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
hterain = (float *)calloc((size_t)(nx*ny),sizeof(float));
mapfct = (float *)calloc((size_t)(8*nx*ny),sizeof(float));
j1 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
j2 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
j3 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
tem2d = (float *)calloc((size_t)(nx*ny),sizeof(float));
tem2dyz = (float *)calloc((size_t)(ny*nz),sizeof(float));
tem2dxz = (float *)calloc((size_t)(nx*nz),sizeof(float));
tem2dns = (float *)calloc((size_t)(nx*ny*(nstyps+1)),sizeof(float));
tem3d1 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
tem3d2 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
tem3d3 = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
tem3dsoil = (float *)calloc((size_t)(nx*ny*nzsoil),sizeof(float));
tem4dsoilns=(float *)calloc((size_t)(nx*ny*nzsoil*(nstyps+1)),sizeof(float));
u = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
v = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
qv = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
ubar = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
vbar = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
uprt = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
vprt = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
ptprt = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
pprt = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
ptbar = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
pbar = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
qvbar = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
rhostr = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
trigs1 = (float *)calloc((size_t)(3*(nx-1)/2+1),sizeof(float));
trigs2 = (float *)calloc((size_t)(3*(ny-1)/2+1),sizeof(float));
ifax = (int *)calloc((size_t)(13),sizeof(int));
wsave1 = (float *)calloc((size_t)(3*(ny-1)+15),sizeof(float));
wsave2 = (float *)calloc((size_t)(3*(nx-1)+15),sizeof(float));
vwork1 = (float *)calloc((size_t)((nx+1)*(ny+1)),sizeof(float));
vwork2 = (float *)calloc((size_t)(ny*(nx+1)),sizeof(float));
soiltyp = (int *)calloc((size_t)(nx*ny*nstyps),sizeof(int));
stypfrct = (float *)calloc((size_t)(nx*ny*nstyps),sizeof(float));
vegtyp = (int *)calloc((size_t)(nx*ny),sizeof(int));
qcumsrc = (float *)calloc((size_t)(nx*ny*nz*5),sizeof(float));
prcrate = (float *)calloc((size_t)(nx*ny*4),sizeof(float));
exbcbuf = (float *)calloc((size_t)(exbcbufsz),sizeof(float));
initgrdvar(&nx,&ny,&nz,&nzsoil,&nt,&nstyps,&exbcbufsz,
x,y,z,zp,tem3dsoil,hterain,mapfct,
j1,j2,j3,tem3dsoil,tem3d1,tem3d1,tem3d1,tem3d1,tem3dsoil,
u,v,tem3d1,tem3d1,ptprt,pprt,
qv,tem3d1,tem3d1,tem3d1,tem3d1,tem3d1,tem3d1,
tem2dyz,tem2dyz,tem2dxz,tem2dxz,
tem2dyz,tem2dyz,tem2dxz,tem2dxz,
trigs1,trigs2,ifax,ifax,
wsave1,wsave2,vwork1,vwork2,
ubar,vbar,ptbar,pbar,tem3d1,tem3d1,
rhostr,tem3d1,qvbar,tem3d1,tem3d1,
soiltyp,stypfrct,vegtyp,tem2d,tem2d,tem2d,
tem4dsoilns,tem4dsoilns,tem2dns,tem2d,tem2dns,
tem3d1,qcumsrc,tem3d1,vegtyp,tem2d,
tem2d,tem2d,tem2d,
tem2d,tem2d,prcrate,exbcbuf,
tem3d1,tem2d,tem2d,tem2d,tem2d,
tem2d,tem2d,tem2d,tem2d,
tem3dsoil,tem3dsoil,tem3dsoil,tem3dsoil,tem3dsoil,
uprt,vprt,tem3d1,tem3d1,tem3d1,
tem3d1,tem3d1,tem3d1,tem3d1);
printf( " Back from initgrdvar \n");
/* Deallocate unneeded ARPS grid arrays */
free((void *) tem2dyz);
free((void *) tem2dxz);
free((void *) tem2dns);
free((void *) tem4dsoilns);
free((void *) tem3dsoil);
free((void *) ubar);
free((void *) vbar);
free((void *) uprt);
free((void *) vprt);
free((void *) qvbar);
free((void *) rhostr);
free((void *) j1);
free((void *) j2);
free((void *) j3);
free((void *) hterain);
free((void *) mapfct);
free((void *) trigs1);
free((void *) trigs2);
free((void *) ifax);
free((void *) wsave1);
free((void *) wsave2);
free((void *) vwork1);
free((void *) vwork2);
free((void *) soiltyp);
free((void *) stypfrct);
free((void *) vegtyp);
free((void *) qcumsrc);
free((void *) prcrate);
free((void *) exbcbuf);
/* Allocate arrays for computing scalar coordinate system */
xs = (float *)calloc((size_t)(nx),sizeof(float));
ys = (float *)calloc((size_t)(ny),sizeof(float));
zps = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
radcoord(&nx,&ny,&nz,
x,y,z,zp,xs,ys,zps,
&radar_lat,&radar_lon,&radarx,&radary);
/* Allocate arrays for computing environmental wind profile */
zsnd = (float *)calloc((size_t)(nzsnd),sizeof(float));
usnd = (float *)calloc((size_t)(nzsnd),sizeof(float));
vsnd = (float *)calloc((size_t)(nzsnd),sizeof(float));
rfrsnd = (float *)calloc((size_t)(nzsnd),sizeof(float));
ktsnd = (float *)calloc((size_t)(nzsnd),sizeof(float));
dzsnd = z[nz-1]/(nzsnd-1);
for(i=0;i<nzsnd;i++) zsnd[i]=i*dzsnd;
printf(" Environmental wind averaging radius: %f m\n",envavgr);
extenvprf(&nx,&ny,&nz,&nzsnd,
x,y,zp,xs,ys,zps,
u,v,ptprt,pprt,qv,ptbar,pbar,tem3d1,tem3d2,tem3d3,
&radarx,&radary,&radar_alt,&envavgr,
zsnd,ktsnd,usnd,vsnd,rfrsnd);
printf( " Back from extenvprf \n");
free((void *) tem3d1);
free((void *) tem3d2);
free((void *) tem3d3);
free((void *) u);
free((void *) v);
free((void *) ptprt);
free((void *) pprt);
free((void *) qv);
free((void *) ptbar);
free((void *) pbar);
free((void *) ktsnd);
/* printf(" lvl zsnd usnd vsnd rfrsnd"); */
/* for(i=0;i<nzsnd;i++)
printf(" %d %f %f %f %f\n",
i,zsnd[i],usnd[i],vsnd[i],rfrsnd[i]); */
/* Allocate arrays for base radar data in polar coords */
ref0_ptr = (float *)calloc((size_t)(maxazim*n_ref_gates),sizeof(float));
vel0_ptr = (float *)calloc((size_t)(maxazim*n_vel_gates),sizeof(float));
spw0_ptr = (float *)calloc((size_t)(maxazim*n_vel_gates),sizeof(float));
unfv_ptr = (float *)calloc((size_t)(maxazim*n_vel_gates),sizeof(float));
bkgv_ptr = (float *)calloc((size_t)(maxazim*n_vel_gates),sizeof(float));
tmp1_ptr = (float *)calloc((size_t)(maxazim*n_vel_gates),sizeof(float));
azi0_ptr = (float *)calloc((size_t)(maxazim),sizeof(float));
elv0_ptr = (float *)calloc((size_t)(maxazim),sizeof(float));
tim0_ptr = (int *)calloc((size_t)(maxazim),sizeof(int));
bgate_ptr = (int *)calloc((size_t)(maxazim),sizeof(int));
egate_ptr = (int *)calloc((size_t)(maxazim),sizeof(int));
/* Allocate arrays for gridded radar data */
gridvel = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
gridref = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
gridnyq = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
gridtim = (float *)calloc((size_t)(nx*ny*nz),sizeof(float));
tem1d1 = (float *)calloc((size_t)(nz),sizeof(float));
tem1d2 = (float *)calloc((size_t)(nz),sizeof(float));
tem1d3 = (float *)calloc((size_t)(nz),sizeof(float));
tem1d4 = (float *)calloc((size_t)(nz),sizeof(float));
tem1d5 = (float *)calloc((size_t)(nz),sizeof(float));
tem1d6 = (float *)calloc((size_t)(nz),sizeof(float));
varsort = (float *)calloc((size_t)(maxsort),sizeof(float));
/* get data indices needed for other a2io library calls */
ref_index=get_field_num("DBZ");
printf (" Retrieved reflectivity index as %i\n",ref_index);
vel_index=get_field_num("VEL");
printf (" Retrieved velocity index as %i\n",vel_index);
spw_index=get_field_num("SPW");
printf (" Retrieved spectrum width index as %i\n",spw_index);
/* Misc initializations */
n_azim=0;
v_nyq=-999.;
initvol_ray = TRUE;
initial_ray = TRUE;
ref_ok = FALSE;
vel_ok = FALSE;
i_first_scan = TRUE;
i_last_scan = FALSE;
write_and_exit = FALSE;
read_next = TRUE;
alls_well = TRUE;
alloc_vol = TRUE;
/* Begin infinite loop to continuously read radar data */
while (alls_well) {
/* Begin loop to fill buffer arrays with data from the circular buffer.
Call remap routines and reset pointers at the end of a volume scan */
if (read_next)
{
# ifdef VERBOSE
printf( " Calling read_radial \n" );
# endif
io_stat=read_radial();
if(io_stat == 1)
{
printf( " Read_radial returned double eof \n");
write_and_exit = TRUE;
}
# ifdef VERBOSE
printf( " Back from read_radial \n");
# endif
}
else
read_next = TRUE;
/* Test for existence of velocity or reflectivity data */
if ( get_status(ref_index) == GOOD_STATUS ||
get_status(vel_index) == GOOD_STATUS ) {
knt_bad_stat = 0;
i_angle = get_fixed_angle();
i_scan = get_scan();
i_tilt = get_tilt();
# ifdef VERBOSE
printf( " Good status received\n");
printf( " i_angle = %i i_tilt = %i i_scan = %i\n", i_angle, i_tilt,i_scan);
# endif
if ( initial_ray == TRUE ) {
ref_ptr = ref0_ptr;
vel_ptr = vel0_ptr;
spw_ptr = spw0_ptr;
past_scan = i_scan;
past_tilt = i_tilt;
past_angle = i_angle;
iyr = get_year();
imon = get_month();
iday = get_day();
ihour = get_hour();
imin = get_min();
isec = get_sec();
iyear=iyr+1900;
if ( iyear < 1960 )
iyear = iyear + 100;
rfrst_ref = get_first_gate(ref_index);
gsp_ref = get_gate_spacing(ref_index);
ng_ref = get_number_of_gates(ref_index);
if ( rfrst_ref < refrngmin ) refrngmin=rfrst_ref;
rmax=rfrst_ref+ng_ref*gsp_ref;
if ( rmax > refrngmax ) refrngmax=rmax;
i_vcp=get_vcp();
i_src=1;
v_nyq=-999.;
if( initvol_ray == TRUE ) {
iyr1 = iyr,
imon1 = imon;
iday1 = iday;
ihour1 = ihour;
imin1 = imin;
isec1 = isec;
ctim2abss(&iyear,&imon,&iday,&ihour,&imin,&isec,&itime1);
initvol_ray = FALSE;
}
if( alloc_vol == TRUE) {
switch (i_vcp)
{
case 11 :
maxelv = 16;
break;
case 21 :
maxelv = 11;
break;
case 31 :
maxelv = 7;
break;
case 32 :
maxelv = 7;
break;
default :
maxelv = MAX_ELEV;
break;
}
/* allocate radar volume arrays */
printf (" Allocating radar memory\n");
kntrgat = (int *)calloc((size_t)(maxelv*maxazim),sizeof(int));
kntrazm = (int *)calloc((size_t)(maxelv),sizeof(int));
kntvgat = (int *)calloc((size_t)(maxelv*maxazim),sizeof(int));
kntvazm = (int *)calloc((size_t)(maxelv),sizeof(int));
nyqvvol = (float *)calloc((size_t)(maxelv),sizeof(float));
if ( nyqvvol == NULL ) printf( " Error allocating nyqvvol \n");
timevol = (int *)calloc((size_t)(maxazim*maxelv),sizeof(float));
if ( timevol == NULL ) printf( " Error allocating timevol \n");
rngrvol = (float *)calloc((size_t)(n_ref_gates*maxelv),sizeof(float));
if ( rngrvol == NULL ) printf( " Error allocating rngrvol \n");
azmrvol = (float *)calloc((size_t)(maxazim*maxelv),sizeof(float));
if ( azmrvol == NULL ) printf( " Error allocating azmrvol \n");
elvrvol = (float *)calloc((size_t)(maxazim*maxelv),sizeof(float));
if ( elvrvol == NULL ) printf( " Error allocating elvrvol \n");
rngvvol = (float *)calloc((size_t)(n_vel_gates*maxelv),sizeof(float));
if ( rngvvol == NULL ) printf( " Error allocating rngvvol \n");
azmvvol = (float *)calloc((size_t)(maxazim*maxelv),sizeof(float));
if ( azmvvol == NULL ) printf( " Error allocating azmvol \n");
elvvvol = (float *)calloc((size_t)(maxazim*maxelv),sizeof(float));
if ( elvvvol == NULL ) printf( " Error allocating elvvol \n");
refvol = (float *)calloc((size_t)(n_ref_gates*maxazim*maxelv),sizeof(float));
if ( refvol == NULL ) printf( " Error allocating refvol \n");
velvol = (float *)calloc((size_t)(n_vel_gates*maxazim*maxelv),sizeof(float));
if ( velvol == NULL ) printf( " Error allocating velvol \n");
rxvol = (float *)calloc((size_t)(maxgate*maxazim*maxelv),sizeof(float));
if ( rxvol == NULL ) printf( " Error allocating rxvol \n");
ryvol = (float *)calloc((size_t)(maxgate*maxazim*maxelv),sizeof(float));
if ( ryvol == NULL ) printf( " Error allocating ryvol \n");
rzvol = (float *)calloc((size_t)(maxgate*maxazim*maxelv),sizeof(float));
if ( rzvol == NULL ) printf( " Error allocating rzvol \n");
elvmean = (float *)calloc((size_t)(maxelv),sizeof(float));
if ( elvmean == NULL ) printf( " Error allocating elvmean \n");
/* Allocate arrays for VAD wind profile */
vadhgt = (float *)calloc((size_t)(maxelv),sizeof(float));
if ( vadhgt == NULL ) printf( " Error allocating vadhgt \n");
vaddir = (float *)calloc((size_t)(maxelv),sizeof(float));
if ( vaddir == NULL ) printf( " Error allocating vaddir \n");
vadspd = (float *)calloc((size_t)(maxelv),sizeof(float));
if ( vadspd == NULL ) printf( " Error allocating vadspd \n");
alloc_vol = FALSE;
}
# ifdef VERBOSE
printf (" iyr = %i imon = %i iday = %i\n",iyr,imon,iday);
printf (" VCP number: %i Max tilts: %i\n",i_vcp,maxelv);
printf ("ihour = %i imin = %i isec = %i\n",ihour,imin,isec);
# endif
initial_ray = FALSE;
} /* initial_ray = TRUE */
if( i_tilt == past_tilt && i_scan == past_scan &&
n_azim < N_AZIM_TILT )
{
n_azim ++;
iyr = get_year();
imon = get_month();
iday = get_day();
ihour = get_hour();
imin = get_min();
isec = get_sec();
iyear=iyr+1900;
if ( iyear < 1960 )
iyear = iyear + 100;
ctim2abss(&iyear,&imon,&iday,&ihour,&imin,&isec,&itime);
*(azi0_ptr+n_azim-1) = 0.01 * (float) get_azi();
*(elv0_ptr+n_azim-1) = 0.01 * (float) get_fixed_angle();
*(tim0_ptr+n_azim-1) = itime;
v_nyq = 0.01 * (float) get_nyquist();
# ifdef VERBOSE
printf( " INFO FOR n_azim = %i \n",n_azim);
printf( " ref_ptr = %i vel_ptr = %i\n",
ref_ptr,vel_ptr);
# endif
if ( (n_azim-1) % 60 == 0)
printf( " eleva = %f azim = %f Nyqst = %f\n",
*(elv0_ptr+n_azim-1),*(azi0_ptr+n_azim-1),v_nyq);
gsp_ref = get_gate_spacing(ref_index);
rfrst_ref = get_first_gate(ref_index);
gsp_vel = get_gate_spacing(vel_index);
rfrst_vel = get_first_gate(vel_index);
# ifdef VERBOSE
ng_ref = get_number_of_gates(ref_index);
ng_vel = get_number_of_gates(vel_index);
printf( " ref: Number of gates = %i, first gate = %i\n",
ng_ref,rfrst_ref);
printf( " vel: Number of gates = %i, first gate = %i\n",
ng_vel,rfrst_vel);
# endif
io_stat = get_data_field(ref_index, ref_ptr, n_ref_gates);
io_stat = get_data_field(vel_index, vel_ptr, n_vel_gates);
io_stat = get_data_field(spw_index, spw_ptr, n_vel_gates);
# ifdef VERBOSE
printf( " sample reflectivities %f %f %f %f\n",
*(ref_ptr+20),*(ref_ptr+40),
*(ref_ptr+60),*(ref_ptr+80));
printf( " sample velocities %f %f %f %f\n",
*(vel_ptr+20),*(vel_ptr+40),
*(vel_ptr+60),*(vel_ptr+80));
printf( " sample spectrum wdths %f %f %f %f\n",
*(spw_ptr+20),*(spw_ptr+40),
*(spw_ptr+60),*(spw_ptr+80));
# endif
if( get_status(ref_index) == GOOD_STATUS ) ref_ok = TRUE;
if( get_status(vel_index) == GOOD_STATUS ) vel_ok = TRUE;
ref_ptr += N_REF_GATES;
vel_ptr += N_VEL_GATES;
spw_ptr += N_VEL_GATES;
}
else
{
if( (i_angle - past_angle ) < -50 || i_scan != past_scan )
i_last_scan = TRUE;
ctim2abss(&iyear,&imon,&iday,&ihour,&imin,&isec,&itime);
if(gsp_vel < 0 ) gsp_vel = 250;
if(rfrst_vel < 0 ) rfrst_vel = 250;
if( i_first_scan )
{
printf( " Calling rmpinit i_first_scan = %d\n",
i_first_scan);
rmpinit(&nx,&ny,&nz,
&n_ref_gates,&n_vel_gates,&maxazim,&maxelv,
kntrgat,kntrazm,kntrelv,
kntvgat,kntvazm,kntvelv,
nyqvvol,timevol,
rngrvol,azmrvol,elvrvol,refvol,
rngvvol,azmvvol,elvvvol,velvol,
gridvel,gridref,gridnyq,gridtim);
printf( " Back from rmpinit \n");
}
if ( ref_ok == TRUE )
{
if ( qc_on == TRUE ) {
printf( " Calling despekl 1 for reflectivity\n");
despekl(&n_ref_gates,&maxazim,&n_ref_gates,&n_azim,&refchek,
ref0_ptr,tmp1_ptr);
}
for ( i = 0; i < n_azim; ++i) {
if (*(elv0_ptr+i) > refelvmax) refelvmax = *(elv0_ptr+i);
if (*(elv0_ptr+i) < refelvmin) refelvmin = *(elv0_ptr+i);
}
printf( " Calling volbuild 1 for reflectivity\n");
volbuild(&n_ref_gates,&maxazim,&maxelv,&n_ref_gates,&n_azim,
&setoff,&settim,
kntrgat,kntrazm,kntrelv,
&gsp_ref,&rfrst_ref,&refchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,ref0_ptr,
nyqvvol,timevol,
rngrvol,azmrvol,elvrvol,refvol);
printf( " Back from volbuild \n");
# ifdef VERBOSE
printf( "Write ref at tilt = %d\n",i_tilt);
# endif
# ifdef UNFOLD_DBG
strcpy(varid,"rref__");
p_tilt = i_tilt-1;
elev = *(elv0_ptr+5);
p_angle = (int)(100*elev);
wrtvel(&p_angle,&p_tilt,varid,
&iyr,&imon,&iday,&ihour,&imin,&isec,
&gsp_ref,&rfrst_ref,
&n_ref_gates,&maxazim,&n_ref_gates,&n_azim,
azi0_ptr,elv0_ptr,ref0_ptr);
# endif
}
if ( vel_ok == TRUE ) {
if ( qc_on == TRUE ) {
printf( " Calling desepkl 1 for velocity\n");
despekl(&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,&velchek,
vel0_ptr,tmp1_ptr);
# ifdef VERBOSE
printf( "Write vel at tilt = %d\n",i_tilt);
# endif
# ifdef UNFOLD_DBG
p_tilt = i_tilt-1;
elev = *(elv0_ptr+5);
p_angle = (int)(100*elev);
strcpy(varid,"rvel__");
wrtvel(&p_angle,&p_tilt,varid,
&iyr,&imon,&iday,&ihour,&imin,&isec,
&gsp_vel,&rfrst_vel,
&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
azi0_ptr,elv0_ptr,vel0_ptr);
# endif
printf( " Calling unfnqc 1\n");
unfnqc(&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
&nzsnd, zsnd,usnd,vsnd,rfrsnd,
&bkgopt,&shropt,&rfropt,&gsp_vel,&rfrst_vel,
rlat_ptr,rlon_ptr,ralt_ptr,
vel0_ptr,spw0_ptr,elv0_ptr,azi0_ptr,&v_nyq,
unfv_ptr,bkgv_ptr,bgate_ptr,egate_ptr,tmp1_ptr);
# ifdef UNFOLD_DBG
strcpy(varid,"unfvel");
wrtvel(&p_angle,&p_tilt,varid,
&iyr,&imon,&iday,&ihour,&imin,&isec,
&gsp_vel,&rfrst_vel,
&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
azi0_ptr,elv0_ptr,unfv_ptr);
# endif
printf( " Calling volbuild 1-1 for velocity\n");
volbuild(&n_vel_gates,&maxazim,&maxelv,&n_vel_gates,&n_azim,
&seton,&seton,
kntvgat,kntvazm,kntvelv,
&gsp_vel,&rfrst_vel,&velchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,unfv_ptr,
nyqvvol,timevol,
rngvvol,azmvvol,elvvvol,velvol);
printf( " Back from volbuild \n");
} else {
printf( " Calling volbuild 1-2 for velocity\n");
volbuild(&n_vel_gates,&maxazim,&maxelv,&n_vel_gates,&n_azim,
&seton,&seton,
kntvgat,kntvazm,kntvelv,
&gsp_vel,&rfrst_vel,&velchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,vel0_ptr,
nyqvvol,timevol,
rngvvol,azmvvol,elvvvol,velvol);
printf( " Back from volbuild \n");
}
}
ref_ok = FALSE;
vel_ok = FALSE;
if ( i_last_scan ) {
if ( qc_on == TRUE ) {
printf( " Calling APDETECT 1 \n");
apdetect(&n_ref_gates,&n_vel_gates,&maxazim,&maxelv,
kntrgat,kntrazm,kntrelv,
kntvgat,kntvazm,kntvelv,
&refchek,&velchek,
&gsp_ref,&gsp_vel,
&winszrad,&winszazim,&i_vcp,
rngrvol,azmrvol,elvrvol,
rngvvol,azmvvol,elvvvol,
refvol,velvol,tmp1_ptr,
&istatus);
printf( " Back from APDETECT 1\n");
}
printf( " Calling mkradfnm \n");
mkradfnm(&dmpfmt,dir_ptr,&len_dir,rname_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,
fnm_ptr,&len_fname);
*(fnm_ptr+len_fname)='\0';
printf ("\n Filename for this volume 1: %s\n\n",fnm_ptr);
printf( " Calling remapvol 1 for reflectivity\n");
printf( " Reflectivity tilt count: %d\n", krelv);
remapvol(&n_ref_gates,&maxazim,&maxelv,&nx,&ny,&nz,&nzsnd,&maxsort,
&refdump,&setoff,&settim,&rfropt,
refid_ptr,refnam_ptr,refunit_ptr,
&refchek,&refmiss,&refmedl,&refdazl,&iordref,
kntrgat,kntrazm,kntrelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,&itime1,
rngrvol,azmrvol,elvrvol,
refvol,timevol,nyqvvol,rxvol,ryvol,rzvol,
xs,ys,zps,zsnd,rfrsnd,varsort,elvmean,
gridref,gridtim,gridnyq,&istatus);
printf( " Back from remapvol\n");
if ( ref2d == TRUE ) {
printf( " Calling remap2d 1 for reflectivity\n");
remap2d(&n_ref_gates,&maxazim,&maxelv,&nx,&ny,&nzsnd,&maxsort,
&ref2d,&rfropt,
ref2did_ptr,ref2dnam_ptr,refunit_ptr,
&refchek,&refmiss,&refmedl,&refdazl,&iordref,
kntrgat,kntrazm,kntrelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,
rngrvol,azmrvol,elvrvol,
refvol,rxvol,ryvol,
xs,ys,zsnd,rfrsnd,varsort,tem2d,&istatus);
printf( " Back from remap2d\n");
}
if ( velproc == TRUE || vad == TRUE ) {
if( qc_on == TRUE ) {
printf( " Calling quadfit 1 for velocity\n");
quadfit(&n_vel_gates,&maxazim,&maxelv,&nzsnd,&maxsort,
&velchek,&velmedl,&veldazl,&iordvel,&rfropt,
kntvgat,kntvazm,kntvelv,
&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,zsnd,rfrsnd,
rngvvol,azmvvol,elvvvol,rxvol,ryvol,rzvol,
varsort,elvmean,velvol,
&istatus);
printf( " Back from quadfit 1\n");
}
}
if ( vad == TRUE ) {
mkvadfnm(dir_ptr,&len_dir,rname_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,
vadfnm_ptr,&len_vfname);
*(vadfnm_ptr+len_vfname)='\0';
printf ("\n Filename for VAD profile: %s\n\n",vadfnm_ptr);
vadvol(&n_vel_gates,&maxazim,&maxelv,
&velchek,&rngvad,&rngwid,&minknt,
kntvgat,kntvazm,kntvelv,
rngvvol,azmvvol,elvvvol,velvol,
vadhgt,vaddir,vadspd);
wtvadprf(&maxelv,vadfnm_ptr,rname_ptr,
rlat_ptr,rlon_ptr,ralt_ptr,
vadhgt,vaddir,vadspd);
}
if ( velproc == TRUE ) {
printf( " Calling remapvol 1 for velocity\n");
printf( " Velocity tilt count: %d\n", kvelv);
remapvol(&n_vel_gates,&maxazim,&maxelv,
&nx,&ny,&nz,&nzsnd,&maxsort,
&veldump,&seton,&seton,&rfropt,
velid_ptr,velnam_ptr,velunit_ptr,
&velchek,&velmiss,&velmedl,&veldazl,&iordvel,
kntvgat,kntvazm,kntvelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,&itime1,
rngvvol,azmvvol,elvvvol,
velvol,timevol,nyqvvol,rxvol,ryvol,rzvol,
xs,ys,zps,zsnd,rfrsnd,varsort,elvmean,
gridvel,gridtim,gridnyq,&istatus);
printf( " Back from remapvol\n");
if ( vel2d == TRUE ) {
printf( " Calling remap2d 1 for velocity\n");
remap2d(&n_vel_gates,&maxazim,&maxelv,&nx,&ny,&nzsnd,&maxsort,
&vel2d,&rfropt,
vel2did_ptr,vel2dnam_ptr,velunit_ptr,
&velchek,&velmiss,&velmedl,&veldazl,&iordvel,
kntvgat,kntvazm,kntvelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,
rngvvol,azmvvol,elvvvol,
velvol,rxvol,ryvol,
xs,ys,zsnd,rfrsnd,varsort,tem2d,&istatus);
printf( " Back from remap2d\n");
}
} /* velproc */
/* wrtrad88(1,fnm_ptr,rname_ptr,rlat_ptr,rlon_ptr,ralt_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,&i_vcp,
&nx,&ny,&nz,gridvel,gridref,gridnyq,gridtim,xs,ys,zps);
*/
if (refrngmin < rngmin) refrngmin = rngmin;
if (refrngmax > rngmax) refrngmax = rngmax;
printf( " Calling wtradcol\n");
wtradcol(&nx,&ny,&nz,&dmpfmt,&irad_fmt,&hdf4cmpr,
fnm_ptr,rname_ptr,rlat_ptr,rlon_ptr,ralt_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,&i_vcp,&i_src,
&refelvmin,&refelvmax,&refrngmin,&refrngmax,
xs,ys,zps,gridvel,gridref,gridnyq,gridtim,
tem1d1,tem1d2,tem1d3,tem1d4,tem1d5,tem1d6);
printf( " Back from wtradcol\n");
if(compr_on == TRUE)
compr_file_sub(fnm_ptr,compr_name,i_status);
if(xmit_on == TRUE)
send_file_sub(fnm_ptr,xmit_method,
remote_user,remote_host,dest_dir,i_status);
initvol_ray = TRUE;
/* Deallocate radar volume arrays */
printf( " Freeing radar volume array memory\n");
free ((void *)kntrgat);
free ((void *)kntrazm);
free ((void *)kntvgat);
free ((void *)kntvazm);
free ((void *)rngrvol);
free ((void *)azmrvol);
free ((void *)elvrvol);
free ((void *)rngvvol);
free ((void *)azmvvol);
free ((void *)elvvvol);
free ((void *)refvol);
free ((void *)velvol);
free ((void *)rxvol);
free ((void *)ryvol);
free ((void *)elvmean);
free ((void *)vadhgt);
free ((void *)vaddir);
free ((void *)vadspd);
alloc_vol = TRUE;
} /* if i_last_scan */
i_last_scan = FALSE;
i_first_scan = FALSE;
if( (i_angle - past_angle) < -50 || i_scan != past_scan )
{
i_first_scan = TRUE;
past_angle= i_angle;
}
n_azim = 0;
initial_ray = TRUE;
read_next = FALSE;
} /* i_tilt=past_tilt end of tilt-down or n_ray overflow block */
} /* GOOD STATUS */
/* For bad status, increment bad status counter and try again. */
else if ( knt_bad_stat < MAX_BAD_STAT &&
write_and_exit == FALSE )
{
# ifdef VERBOSE
printf( " Bad status received for data\n");
# endif
knt_bad_stat ++;
/* Once 1000 consecutive bad stati have been received, assume end
of data and dump what might be in the buffer. */
}
else
{
printf ( " %i bad read status reports received \n",
knt_bad_stat);
knt_bad_stat = 0;
if ( n_azim > 0)
{
printf ( " Transferring %i available radials\n", n_azim);
if( (i_angle - past_angle) < -50 || write_and_exit == TRUE )
i_last_scan = TRUE;
/* call the FORTRAN remapper module */
if( i_first_scan ) {
printf( " Calling rmpinit 2 i_first_scan = %d\n",
i_first_scan);
rmpinit(&nx,&ny,&nz,
&n_ref_gates,&n_vel_gates,&maxazim,&maxelv,
kntrgat,kntrazm,kntrelv,
kntvgat,kntvazm,kntvelv,
nyqvvol,timevol,
rngrvol,azmrvol,elvrvol,refvol,
rngvvol,azmvvol,elvvvol,velvol,
gridvel,gridref,gridnyq,gridtim);
printf( " Back from rmpinit \n");
}
if ( ref_ok == TRUE ) {
if ( qc_on == TRUE ) {
printf( " Calling despekl 2 for reflectivity\n");
despekl(&n_ref_gates,&maxazim,&n_ref_gates,&n_azim,&refchek,
ref0_ptr,tmp1_ptr);
}
printf( " Calling volbuild 2-1 for reflectivity\n");
volbuild(&n_ref_gates,&maxazim,&maxelv,&n_ref_gates,&n_azim,
&setoff,&settim,
kntrgat,kntrazm,kntrelv,
&gsp_ref,&rfrst_ref,&refchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,ref0_ptr,
nyqvvol,timevol,
rngrvol,azmrvol,elvrvol,refvol);
printf( " Back from volbuild \n");
}
if ( vel_ok == TRUE ) {
if ( qc_on == TRUE ) {
printf( " Calling despekl 1 for velocity\n");
despekl(&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,&velchek,
vel0_ptr,tmp1_ptr);
# ifdef UNFOLD_DBG
p_tilt = i_tilt-1;
elev = *(elv0_ptr+5);
p_angle = (int)(100*elev);
strcpy(varid,"rvel__");
wrtvel(&p_angle,&p_tilt,varid,
&iyr,&imon,&iday,&ihour,&imin,&isec,
&gsp_vel,&rfrst_vel,
&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
azi0_ptr,elv0_ptr,vel0_ptr);
# endif
printf( " Calling unfnqc 1\n");
unfnqc(&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
&nzsnd,zsnd,usnd,vsnd,rfrsnd,
&bkgopt,&shropt,&rfropt,&gsp_vel,&rfrst_vel,
rlat_ptr,rlon_ptr,ralt_ptr,
vel0_ptr,spw0_ptr,elv0_ptr,azi0_ptr,&v_nyq,
unfv_ptr,bkgv_ptr,bgate_ptr,egate_ptr,tmp1_ptr);
# ifdef UNFOLD_DBG
strcpy(varid,"unfvel");
wrtvel(&p_angle,&p_tilt,varid,
&iyr,&imon,&iday,&ihour,&imin,&isec,
&gsp_vel,&rfrst_vel,
&n_vel_gates,&maxazim,&n_vel_gates,&n_azim,
azi0_ptr,elv0_ptr,unfv_ptr);
# endif
printf( " Calling volbuild 2-1 for velocity\n");
volbuild(&n_vel_gates,&maxazim,&maxelv,&n_vel_gates,&n_azim,
&seton,&seton,
kntvgat,kntvazm,kntvelv,
&gsp_vel,&rfrst_vel,&velchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,unfv_ptr,
nyqvvol,timevol,
rngvvol,azmvvol,elvvvol,velvol);
printf( " Back from volbuild \n");
} else {
printf( " Calling volbuild 2-2 for velocity\n");
volbuild(&n_vel_gates,&maxazim,&maxelv,&n_vel_gates,&n_azim,
&seton,&seton,
kntvgat,kntvazm,kntvelv,
&gsp_vel,&rfrst_vel,&velchek,
&v_nyq,tim0_ptr,
azi0_ptr,elv0_ptr,vel0_ptr,
nyqvvol,timevol,
rngvvol,azmvvol,elvvvol,velvol);
printf( " Back from volbuild \n");
}
}
ref_ok = FALSE;
vel_ok = FALSE;
if ( i_last_scan ) {
mkradfnm(&dmpfmt,dir_ptr,&len_dir,rname_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,
fnm_ptr,&len_fname);
*(fnm_ptr+len_fname)='\0';
printf ("\n Filename for this volume: %s\n\n",fnm_ptr);
if ( qc_on == TRUE ) {
printf( " Calling APDETECT 2 \n");
apdetect(&n_ref_gates,&n_vel_gates,&maxazim,&maxelv,
kntrgat,kntrazm,kntrelv,
kntvgat,kntvazm,kntvelv,
&refchek,&velchek,
&gsp_ref,&gsp_vel,
&winszrad,&winszazim,&i_vcp,
rngrvol,azmrvol,elvrvol,
rngvvol,azmvvol,elvvvol,
refvol,velvol,tmp1_ptr,
&istatus);
printf( " Back from APDETECT 2\n");
}
printf( " Calling remapvol 2 for reflectivity\n");
printf( " Reflectivity tilt count: %d\n", krelv);
remapvol(&n_ref_gates,&maxazim,&maxelv,
&nx,&ny,&nz,&nzsnd,&maxsort,
&refdump,&setoff,&settim,&rfropt,
refid_ptr,refnam_ptr,refunit_ptr,
&refchek,&refmiss,&refmedl,&refdazl,&iordref,
kntrgat,kntrazm,kntrelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,&itime1,
rngrvol,azmrvol,elvrvol,
refvol,timevol,nyqvvol,rxvol,ryvol,rzvol,
xs,ys,zps,zsnd,rfrsnd,varsort,elvmean,
gridref,gridtim,gridnyq,&istatus);
printf( " Back from remapvol\n");
if ( ref2d == TRUE ) {
printf( " Calling remap2d 2 for reflectivity\n");
remap2d(&n_ref_gates,&maxazim,&maxelv,&nx,&ny,&nzsnd,&maxsort,
&ref2d,&rfropt,
ref2did_ptr,ref2dnam_ptr,refunit_ptr,
&refchek,&refmiss,&refmedl,&refdazl,&iordref,
kntrgat,kntrazm,kntrelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,
rngrvol,azmrvol,elvrvol,
refvol,rxvol,ryvol,
xs,ys,zsnd,rfrsnd,varsort,tem2d,&istatus);
printf( " Back from remap2d\n");
}
if ( velproc == TRUE || vad == TRUE ) {
if ( qc_on == TRUE ) {
printf( " Calling quadfit 2 for velocity\n");
quadfit(&n_vel_gates,&maxazim,&maxelv,&nzsnd,&maxsort,
&velchek,&velmedl,&veldazl,&iordvel,&rfropt,
kntvgat,kntvazm,kntvelv,
&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,zsnd,rfrsnd,
rngvvol,azmvvol,elvvvol,rxvol,ryvol,rzvol,
varsort,elvmean,velvol,
&istatus);
printf( " Back from quadfit 2\n");
}
}
if ( vad == TRUE ) {
mkvadfnm(dir_ptr,&len_dir,rname_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,
vadfnm_ptr,&len_vfname);
*(vadfnm_ptr+len_vfname)='\0';
printf ("\n Filename for VAD profile: %s\n\n",vadfnm_ptr);
vadvol(&n_vel_gates,&maxazim,&maxelv,
&velchek,&rngvad,&rngwid,&minknt,
kntvgat,kntvazm,kntvelv,
rngvvol,azmvvol,elvvvol,velvol,
vadhgt,vaddir,vadspd);
wtvadprf(&maxelv,vadfnm_ptr,rname_ptr,
rlat_ptr,rlon_ptr,ralt_ptr,
vadhgt,vaddir,vadspd);
}
if ( velproc == TRUE ) {
printf( " Calling remapvol 2 for velocity\n");
printf( " Velocity tilt count: %d\n", kvelv);
remapvol(&n_vel_gates,&maxazim,&maxelv,
&nx,&ny,&nz,&nzsnd,&maxsort,
&veldump,&seton,&seton,&rfropt,
velid_ptr,velnam_ptr,velunit_ptr,
&velchek,&velmiss,&refmedl,&veldazl,&iordvel,
kntvgat,kntvazm,kntvelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,&itime1,
rngvvol,azmvvol,elvvvol,
velvol,timevol,nyqvvol,rxvol,ryvol,rzvol,
xs,ys,zps,zsnd,rfrsnd,varsort,elvmean,
gridvel,gridtim,gridnyq,&istatus);
printf( " Back from remapvol\n");
if ( vel2d == TRUE ) {
printf( " Calling remap2d 2 for velocity\n");
remap2d(&n_vel_gates,&maxazim,&maxelv,&nx,&ny,&nzsnd,&maxsort,
&vel2d,&rfropt,
vel2did_ptr,vel2dnam_ptr,velunit_ptr,
&velchek,&velmiss,&velmedl,&veldazl,&iordvel,
kntvgat,kntvazm,kntvelv,
rlat_ptr,rlon_ptr,&radarx,&radary,ralt_ptr,&dazim,
&rngmin,&rngmax,
rngvvol,azmvvol,elvvvol,
velvol,rxvol,ryvol,
xs,ys,zsnd,rfrsnd,varsort,tem2d,&istatus);
printf( " Back from remap2d\n");
}
} /* velproc */
if (refrngmin < rngmin) refrngmin = rngmin;
if (refrngmax > rngmax) refrngmax = rngmax;
printf(" Calling wtradcol\n");
wtradcol(&nx,&ny,&nz,&dmpfmt,&irad_fmt,&hdf4cmpr,
fnm_ptr,rname_ptr,rlat_ptr,rlon_ptr,ralt_ptr,
&iyr1,&imon1,&iday1,&ihour1,&imin1,&isec1,&i_vcp,&i_src,
&refelvmin,&refelvmax,&refrngmin,&refrngmax,
xs,ys,zps,gridvel,gridref,gridnyq,gridtim,
tem1d1,tem1d2,tem1d3,tem1d4,tem1d5,tem1d6);
printf(" Back from wtradcol\n");
/*
* Added by Ming Hu
*/
if (wrttilts_on == TRUE ) {
wrttilts(fnm_ptr,&n_vel_gates,&n_ref_gates,&maxazim,&maxelv,
rngvvol,azmvvol,elvvvol,velvol,
rngrvol,azmrvol,elvrvol,refvol,
&radar_alt,&radar_lat,&radar_lon,
&iyear,&imon,&iday,&ihour,&imin,&isec );
}
if(compr_on == TRUE)
compr_file_sub(fnm_ptr,compr_name,i_status);
if(xmit_on == TRUE)
send_file_sub(fnm_ptr,xmit_method,
remote_user,remote_host,dest_dir,i_status);
initvol_ray = TRUE;
/* Deallocate radar volume arrays */
printf( " Freeing data array memory\n");
free ((void *)rngrvol);
free ((void *)azmrvol);
free ((void *)elvrvol);
free ((void *)rngvvol);
free ((void *)azmvvol);
free ((void *)elvvvol);
free ((void *)refvol);
free ((void *)velvol);
free ((void *)rxvol);
free ((void *)ryvol);
free ((void *)rzvol);
free ((void *)elvmean);
free ((void *)vadhgt);
free ((void *)vaddir);
free ((void *)vadspd);
alloc_vol = TRUE;
}
if(write_and_exit == TRUE) exit(0);
i_last_scan = FALSE;
i_first_scan = FALSE;
printf( "i_angle %d past_angle %d\n",i_angle,past_angle);
if( (i_angle - past_angle) < -50 ) {
i_first_scan = TRUE;
past_angle = i_angle;
}
n_azim = 0;
initial_ray = TRUE;
printf(" Deallocating radar volume arrays\n");
/* Deallocate radar volume arrays */
free ((void *)rngrvol);
free ((void *)azmrvol);
free ((void *)elvrvol);
free ((void *)rngvvol);
free ((void *)azmvvol);
free ((void *)elvvvol);
free ((void *)refvol);
free ((void *)velvol);
free ((void *)rxvol);
free ((void *)ryvol);
free ((void *)rzvol);
} /* close n_azim > 0 block */
if(write_and_exit == TRUE) exit(0);
} /* close velocity status block */
if(write_and_exit == TRUE) exit(0);
} /* close infinite while loop */
} /* end of main */