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

   1 /*
   2     NASA/TRMM, Code 910.1.
   3     This is the TRMM Office Radar Software Library.
   4     Copyright (C) 1996, 1997
   5             Mike Kolander
   6             Space Applications Corporation
   7             Vienna, Virginia
   8
   9     This library is free software; you can redistribute it and/or
  10     modify it under the terms of the GNU Library General Public
  11     License as published by the Free Software Foundation; either
  12     version 2 of the License, or (at your option) any later version.
  13
  14     This library is distributed in the hope that it will be useful,
  15     but WITHOUT ANY WARRANTY; without even the implied warranty of
  16     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17     Library General Public License for more details.
  18
  19     You should have received a copy of the GNU Library General Public
  20     License along with this library; if not, write to the Free
  21     Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  22 */
  23 #ifdef HAVE_CONFIG_H
  24 #include "config.h"
  25 #endif
  26 #ifdef HAVE_LIBTSDISTK
  27 /******************************************************************
  28
  29          Subroutines to write one VOS from a RSL radar structure into one
  30          1B-51/1C-51 HDF file.
  31
  32          A 1B-51/1C-51 HDF file contains multiple VOS's recorded by a radar
  33          site during a 1-hour time period, and the HDF file is named using
  34          the date/hour of the constituent VOS's.
  35
  36          All functions defined herein build the components of the TSDIS
  37          toolkit 'L1B_1C_GV' structure using the data from a RSL radar
  38          structure. These functions are executed via a call from the top-level
  39          RSL function 'RSL_radar_to_hdf()', defined in RSL file
  40          'radar_to_hdf_1.c'.
  41
  42   -----------------------------------------------------------------
  43          Libraries required for execution of this code :
  44       -ltsdistk                    : TSDIS toolkit
  45       -lmfhdf -ldf -ljpeg -lz      : HDF
  46       -lrsl                        : rsl
  47       -lm                          : C math
  48
  49   -----------------------------------------------------------------
  50 *******************************************************************/
  51
  52 #include <math.h>
  53 #include <time.h>
  54 #include <string.h>
  55 #include <stdio.h>
  56 #include <stdlib.h>
  57 #include <unistd.h>
  58 #include <sys/types.h>
  59 #include <sys/stat.h>
  60
  61 /* TSDIS toolkit function and structure definitions. */
  62 #include "IO.h"
  63 #include "IO_GV.h"
  64 /* RSL function and structure definitions. */
  65 #include "rsl.h"
  66 /* Parameter definitions for 1B-51 and 1C-51 HDF
  67          file handling applications using the TSDIS toolkit. */
  68 #include "toolkit_1BC-51_appl.h"
  69
  70
  71 /*************************************************************/
  72 /*                                                           */
  73 /*                    Function Prototypes                    */
  74 /*                                                           */
  75 /*************************************************************/
  76 static int julian(int year, int mo, int day);
  77 int8 ***parmData1byteBuild(Volume *v, int pindex, VosSize *vs);
  78 int16 ***parmData2byteBuild(Radar *radar, PARAMETER_DESCRIPTOR *parmDesc,
  79                                                                                                                 int vindex, int pindex, VosSize *vs);
  80 void cellRangeVectorFill(CELL_RANGE_VECTOR *cellRangeVector,
  81                                                                                                  int vindex, int pindex, VosSize *vs);
  82 void parmDescFill(PARAMETER_DESCRIPTOR *parmDesc, Radar *radar, int vindex);
  83 PARAMETER *parmBuild(Radar *radar, VosSize *vs, int vindex,
  84                                                                                  int pindex);
  85 void rayInfoFill(int32 rayInfoInteger[MAX_SWEEP][MAX_RAY][7],
  86                                                                  float32 rayInfoFloat[MAX_SWEEP][MAX_RAY][4],
  87                                                                  Radar *radar, VosSize *vs);
  88 void sweepInfoFill(SWEEP_INFO sweepInfo[MAX_SWEEP], Radar *radar, VosSize *vs);
  89 void radarDescFill(RADAR_DESCRIPTOR *radarDesc, Radar *radar, int vindex,
  90                                                                          VosSize *vs);
  91 void sensorFill(SENSORS *sensor, Radar *radar, VosSize *vs, int productID);
  92 struct tm *timeUTC(void);
  93 void volDesFill(VOLUME_DESCRIPTORS *volDes, Radar_header *h, VosSize *vs);
  94 void commentsFill(char *comments, VosSize *vs, Radar *radar,
  95                                                                         float *qcParm, int productID);
  96 L1B_1C_GV *gvl1Build(Radar *radar, float *qcParm, VosSize *vs,
  97                                                                                  int productID);
  98
  99 extern int nextVolume(Radar *radar, int last_volume);
 100 extern int radar_verbose_flag;
 101 extern Ray *first_ray_in_volume[MAX_RADAR_VOLUMES];
 102
 103 /*************************************************************/
 104 /*                                                           */
 105 /*                          julian                           */
 106 /*                                                           */
 107 /*************************************************************/
 108 static int daytab[2][13] = {
 109   {0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334, 365},
 110   {0, 31, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335, 366}
 111 };
 112
 113 static int julian(int year, int mo, int day)
 114 {
 115 /* Converts a calendar date (month, day, year) to a Julian date.
 116          Returns:
 117            Julian day.
 118 */
 119   int leap;
 120
 121   leap = (year%4 == 0 && year%100 != 0) || year%400 == 0;
 122   return(day + daytab[leap][mo-1]);
 123 }
 124
 125 /*************************************************************/
 126 /*                                                           */
 127 /*                   parmData1byteBuild                      */
 128 /*                                                           */
 129 /*************************************************************/
 130 int8 ***parmData1byteBuild(Volume *v, int pindex, VosSize *vs)
 131 /* Move all data from one mask volume of the RSL structure into
 132          the array parmData1byte[][][] .
 133 */
 134 {
 135         int sindex, rindex, bindex;  /* Indices for rsl arrays. */
 136         int tk_sindex, tk_rindex;    /* Indices for toolkit arrays. */
 137         int ncell;
 138         int8 ***data1byte;
 139         Ray *ray;
 140
 141         ncell = vs->rsl.ncell[pindex][0];
 142         data1byte = (int8 ***)TKnewParmData1byte(vs->tk.nsweep, vs->rsl.maxNray,
 143                                                                                                                                                                          ncell);
 144         /* Move data values from all non-NULL rsl sweeps into the
 145                  'parmData1byte' array. */
 146         sindex = -1;
 147         for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 148         {
 149                 if (vs->tk.ncell[tk_sindex][pindex] == 0)
 150                 {
 151                         /* No data from this parm in this physical sweep.
 152                                  Fill bins of this 'parmData1byte' sweep with 0 */
 153                         for (rindex=0; rindex<vs->rsl.maxNray; rindex++)
 154                     for (bindex=0; bindex<ncell; bindex++)
 155                             data1byte[tk_sindex][rindex][bindex] = (int8) 0;
 156                         continue;
 157                 }
 158                 sindex++;
 159                 tk_rindex = -1;
 160                 /* Move data values from all non-NULL rsl rays into the
 161                          'parmData1byte' array. */
 162                 for (rindex=0; rindex<v->sweep[sindex]->h.nrays; rindex++)
 163                 {
 164                         if (v->sweep[sindex]->ray[rindex] == NULL) continue;
 165                         tk_rindex++;
 166                         ray = v->sweep[sindex]->ray[rindex];
 167                         /* Move the rsl bin values which exist into the 'parmData1byte' ray.*/
 168                         for (bindex=0; bindex<vs->tk.ncell[tk_sindex][pindex]; bindex++)
 169                         {
 170                           if (bindex >= ray->h.nbins)  /* Short rsl ray? */
 171                                 data1byte[tk_sindex][tk_rindex][bindex] = (int8) 0;
 172                           else  /* Valid bin */
 173                             data1byte[tk_sindex][tk_rindex][bindex] = (int8) ray->range[bindex];
 174                         } /* end for (bindex=0... */
 175                         /* Fill all remaining bins of 'parmData1byte' ray with 0. */
 176                         for (bindex=vs->tk.ncell[tk_sindex][pindex]; bindex<ncell; bindex++)
 177                           data1byte[tk_sindex][tk_rindex][bindex] = (int8) 0;
 178                 } /* end (rindex=0;... */
 179
 180                 /* Fill bins of all remaining rays of this sweep with 0. */
 181                 for (rindex=tk_rindex+1; rindex<vs->rsl.maxNray; rindex++)
 182                   for (bindex=0; bindex<ncell; bindex++)
 183                           data1byte[tk_sindex][rindex][bindex] = (int8) 0;
 184         } /* end for (tk_sindex=0;...*/
 185
 186         return(data1byte);
 187 }
 188
 189 /*************************************************************/
 190 /*                                                           */
 191 /*                   parmData2byteBuild                      */
 192 /*                                                           */
 193 /*************************************************************/
 194 int16 ***parmData2byteBuild(Radar *radar, PARAMETER_DESCRIPTOR *parmDesc,
 195                                                                                                                  int vindex, int pindex, VosSize *vs)
 196 /* Move all ray data from one volume of the RSL structure into
 197          the array parmData2byte[][][] . Data values placed into
 198          parmData2byte[][][] are scaled.
 199 */
 200 {
 201   int sindex, rindex, bindex;  /* Indices for rsl arrays. */
 202   int tk_sindex, tk_rindex;    /* Indices for toolkit arrays. */
 203   int ncell;
 204   float value;
 205   Volume *v;
 206   Ray *ray;
 207   int16 ***data2byte;
 208
 209   ncell = vs->rsl.ncell[pindex][0];
 210   data2byte = (int16 ***)TKnewParmData2byte(vs->tk.nsweep, vs->rsl.maxNray, ncell);
 211   v = radar->v[vindex];
 212   /* Move data values from all non-NULL rsl sweeps into the
 213          'parmData2byte' array. */
 214   sindex = -1;
 215   for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 216   {
 217         if (vs->tk.ncell[tk_sindex][pindex] == 0)
 218         {
 219           /* No data from this parm in this physical sweep.
 220                  Fill bins of this 'parmData2byte' sweeps with NO_VALUE. */
 221           for (rindex=0; rindex<vs->rsl.maxNray; rindex++)
 222                 for (bindex=0; bindex<ncell; bindex++)
 223                   data2byte[tk_sindex][rindex][bindex] = (int16) NO_VALUE;
 224           continue;
 225         }
 226         sindex++;
 227         tk_rindex = -1;
 228         /* Move data values from all non-NULL rsl rays into the
 229            'parmData2byte' array. */
 230         for (rindex=0; rindex<v->sweep[sindex]->h.nrays; rindex++)
 231         {
 232           if (v->sweep[sindex]->ray[rindex] == NULL) continue;
 233           tk_rindex++;
 234           ray = v->sweep[sindex]->ray[rindex];
 235           /* Move the rsl bin values which exist into the 'parmData2byte' ray.*/
 236           for (bindex=0; bindex<vs->tk.ncell[tk_sindex][pindex]; bindex++)
 237           {
 238                 /* if short rsl ray, fill cell with NO_VALUE */
 239                 if (bindex >= ray->h.nbins)
 240                 {
 241                   data2byte[tk_sindex][tk_rindex][bindex] = (int16) NO_VALUE;
 242                   continue;
 243                 }
 244                 value = v->h.f(ray->range[bindex]);
 245                 if (value >= NOECHO)  /* Handle anomalous condition flags */
 246                 {
 247                   if (value == BADVAL)
 248                         data2byte[tk_sindex][tk_rindex][bindex] = (int16) NO_VALUE;
 249                   else if (value == RFVAL)
 250                         data2byte[tk_sindex][tk_rindex][bindex] = (int16) RNG_AMBIG_VALUE;
 251                   else if (value == APFLAG)
 252                         data2byte[tk_sindex][tk_rindex][bindex] = (int16) AP_VALUE;
 253                   else
 254                         data2byte[tk_sindex][tk_rindex][bindex] = (int16) NOECHO_VALUE;
 255                 }
 256                 else  /* Valid rsl data */
 257                 {
 258                   if (vindex == CZ_INDEX)  /* Corrected Z data */
 259                         {
 260                           /* CZ = DZ + dzCal - mask_val * X  ... From Ferrier memo. */
 261                           value = (value + v->h.calibr_const -
 262                                            radar->v[MZ_INDEX]->sweep[sindex]->ray[rindex]->range[bindex] * X);
 263                         }
 264                   else if (vindex == CD_INDEX)  /* Corrected differential Z data */
 265                         {
 266                           /* CD = ZD - mask_val * X */
 267                           value = value -
 268                                 radar->v[MD_INDEX]->sweep[sindex]->ray[rindex]->range[bindex] * X;
 269                         }
 270                   /* Apply scale and offset factors, then store value in
 271                          parmData2byte structure. */
 272                   data2byte[tk_sindex][tk_rindex][bindex] = (int16)
 273                         (value * parmDesc->scaleFactor + parmDesc->offsetFactor);
 274                 } /* end else Valid rsl data */
 275           } /* end for (bindex=0;... */
 276
 277           /* Fill all remaining bins of 'parmData2byte' ray with NO_VALUE. */
 278           for (bindex=vs->tk.ncell[tk_sindex][pindex]; bindex<ncell; bindex++)
 279                 data2byte[tk_sindex][tk_rindex][bindex] = (int16) NO_VALUE;
 280         } /* end (rindex=0;... */
 281
 282         /* Fill bins of all remaining rays of this sweep with NO_VALUE. */
 283         for (rindex=tk_rindex+1; rindex<vs->rsl.maxNray; rindex++)
 284           for (bindex=0; bindex<ncell; bindex++)
 285                 data2byte[tk_sindex][rindex][bindex] = (int16) NO_VALUE;
 286   } /* end for (sindex=0;... */
 287
 288   return(data2byte);
 289 }
 290
 291 /*************************************************************/
 292 /*                                                           */
 293 /*                   cellRangeVectorFill                     */
 294 /*                                                           */
 295 /*************************************************************/
 296 void cellRangeVectorFill(CELL_RANGE_VECTOR *cellRangeVector,
 297                                                                                                  int vindex, int pindex, VosSize *vs)
 298 {
 299 /* Fill in the cellRangeVector values.
 300          cellRangeVector[j] = distance(m) from radar to center of cell[j] .
 301 */
 302         int j;
 303         Ray_header *ray_head;
 304
 305         cellRangeVector->numOfCells = (int32) vs->rsl.ncell[pindex][0];
 306         ray_head = &first_ray_in_volume[vindex]->h;
 307
 308         /* Do first cell (m) */
 309         cellRangeVector->distanceToCell[0] = (float32) (ray_head->range_bin1 +
 310                                                                                                                                                                                                         0.5 * ray_head->gate_size);
 311
 312         /* Do remaining cells. Just add gate_size to previous cellRange value. */
 313         for (j=1; j<cellRangeVector->numOfCells; j++)
 314         {
 315                 cellRangeVector->distanceToCell[j] = (float32)
 316                      (cellRangeVector->distanceToCell[j-1] + ray_head->gate_size);
 317         }
 318 }
 319
 320 /*************************************************************/
 321 /*                                                           */
 322 /*                      parmDescFill                         */
 323 /*                                                           */
 324 /*************************************************************/
 325 void parmDescFill(PARAMETER_DESCRIPTOR *parmDesc, Radar *radar, int vindex)
 326 {
 327         static char *parm_list[20][3] =
 328         {
 329                 { "Z",      "dBz",    "Reflectivity"                     },
 330                 { "V",      "m/s",    "Radial Velocity"                  },
 331                 { "SW",     "m2/s2",  "Spectral Width"                   },
 332                 { "QCZ",    "dBz",    "QC'ed Reflectivity"               },
 333                 { "ZT",     "dBz",    "Total Reflectivity"               },
 334                 { "DR",     "?",      "Differential reflectivity"        },
 335                 { "LR",     "?",      "Differential reflectivity"        },
 336                 { "ZDR",    "dB",     "Differential Reflectivity"        },
 337                 { "DM",     "dBm",    "Received power"                   },
 338                 { "RH",     "-",      "Correlation Coefficient"          },
 339                 { "PH",     "?",      "Phi"                              },
 340                 { "XZ",     "dBz",    "X-band Reflectivity"              },
 341                 { "QCZDR",  "dB",     "QC'ed Differential Reflectivity"  },
 342                 { "QCMZ",   "-",      "Z Mask"                           },
 343                 { "QCMZDR", "-",      "ZDR Mask"                         },
 344                 { "ZE",     "-",      "Edited Reflectivity"              },
 345                 { "VE",     "-",      "Edited Velocity"                  },
 346                 { "--",     "-",      "*******"                          },
 347                 { "--",     "-",      "*******"                          },
 348                 { "--",     "-",      "*******"                          }
 349         };
 350
 351   strncpy(parmDesc->parmName, parm_list[vindex][0], 7);
 352         strncpy(parmDesc->parmDesc, parm_list[vindex][2], 39);
 353         strncpy(parmDesc->parmUnits, parm_list[vindex][1], 7);
 354         parmDesc->interPulsePeriod = (int16) 0;
 355         parmDesc->transFreq = (int16) 1;
 356         /* Receiver Bandwidth (MHz) */
 357         parmDesc->receiverBandwidth = (float32) 0.0;
 358         /* Pulse width (m) */
 359   parmDesc->pulseWidth = (int16) (300.0 *
 360                                                                                                   first_ray_in_volume[vindex]->h.pulse_width);
 361         parmDesc->polarTransWave = (int16) 0;
 362         parmDesc->numOfsamples = (int16) 0;
 363         /* No thresholding done in 1B-51, 1C-51 HDF files. */
 364         strncpy(parmDesc->thresholdField, "NONE", 8);
 365         parmDesc->thresholdValue = (float32) 0.0;
 366         parmDesc->offsetFactor = (float32) 0.0;
 367         if ((vindex == MZ_INDEX) || (vindex == MD_INDEX))
 368         {
 369           parmDesc->parmDataType = (int16) 1;   /* 1_byte mask value. */
 370                 parmDesc->scaleFactor = (float32) 1.0;
 371                 parmDesc->deletedOrMissDataFlag = (int32) 0;
 372         }
 373         else  /* 2_byte data value. */
 374         {
 375           parmDesc->parmDataType = (int16) 2;
 376                 parmDesc->scaleFactor = (float32) SCALE_FACTOR;
 377                 parmDesc->deletedOrMissDataFlag = (int32) NO_VALUE;
 378         }
 379 }
 380
 381 /*************************************************************/
 382 /*                                                           */
 383 /*                       parmBuild                           */
 384 /*                                                           */
 385 /*************************************************************/
 386 PARAMETER *parmBuild(Radar *radar, VosSize *vs, int vindex,
 387                                                                                  int pindex)
 388 {
 389         PARAMETER *parm;
 390
 391         /* Allocate memory for a new parameter structure. */
 392         parm = (PARAMETER *)TKnewGVL1parm();
 393
 394         parmDescFill(&parm->parmDesc, radar, vindex);
 395         cellRangeVectorFill(&parm->cellRangeVector, vindex, pindex, vs);
 396         if ((vindex == MZ_INDEX) || (vindex == MD_INDEX)) /* Mask? */
 397           parm->parmData1byte = (int8 ***)parmData1byteBuild(radar->v[vindex],
 398                                                                                                                                                                                                                         pindex, vs);
 399         else
 400           parm->parmData2byte = (int16 ***)parmData2byteBuild(radar, &parm->parmDesc,
 401                                                                                                                                                                                                                          vindex, pindex, vs);
 402         return(parm);
 403 }
 404
 405 /*************************************************************/
 406 /*                                                           */
 407 /*                         rayInfoFill                       */
 408 /*                                                           */
 409 /*************************************************************/
 410 void rayInfoFill(int32 rayInfoInteger[MAX_SWEEP][MAX_RAY][7],
 411                                                                  float32 rayInfoFloat[MAX_SWEEP][MAX_RAY][4],
 412                                                                  Radar *radar, VosSize *vs)
 413 /* For each ray in the rsl structure, move ray header info into the
 414          arrays rayInfoInteger[][][] and rayInfoFloat[][][].
 415 */
 416 {
 417         int tk_sindex, rindex, tk_rindex;
 418         double second;
 419         Ray_header *ray_head;
 420         Sweep *sweep;
 421         static int32 julday;
 422         static int day=-1;
 423
 424         /* For each physical sweep...*/
 425         for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 426         {
 427                 sweep = vs->rsl.sweep[tk_sindex];
 428                 tk_rindex = -1;
 429                 for (rindex=0; rindex<sweep->h.nrays; rindex++)
 430                 {
 431                         if (sweep->ray[rindex] == NULL) continue;
 432                         tk_rindex++;
 433                                 /***********  Fill in ray info fields.  ***************/
 434                         ray_head = &sweep->ray[rindex]->h;
 435                         /* No. of sweep which contains this ray. */
 436                         rayInfoInteger[tk_sindex][tk_rindex][0] = (int32)(tk_sindex + 1);
 437                         /* Compute Julian Day. Usually do only once per vos. */
 438                         if (day != ray_head->day)
 439                         {
 440                                 day = ray_head->day;  /* Note day & julday are static. */
 441                                 julday = (int32)
 442                                                  julian(ray_head->year, ray_head->month, ray_head->day);
 443                         }
 444                         rayInfoInteger[tk_sindex][tk_rindex][1] = (int32) julday;
 445                         rayInfoInteger[tk_sindex][tk_rindex][2] = (int32) ray_head->hour;
 446                         rayInfoInteger[tk_sindex][tk_rindex][3] = (int32) ray_head->minute;
 447                         rayInfoInteger[tk_sindex][tk_rindex][5] = (int32) (1000.0 *
 448                                                                 modf(ray_head->sec, &second));
 449                         rayInfoInteger[tk_sindex][tk_rindex][4] = (int32) second;
 450                         /* Ray status. 0:Normal , 1:Tansition , 2:Bad , 3:Questionable */
 451 /*                      rayInfoInteger[tk_sindex][tk_rindex][6] = (int32) 0; */
 452
 453                         rayInfoFloat[tk_sindex][tk_rindex][0] = (float32) ray_head->azimuth;
 454                         rayInfoFloat[tk_sindex][tk_rindex][1] = (float32) ray_head->elev;
 455                         /* Store num_of_samples here instead of power. */
 456                         rayInfoFloat[tk_sindex][tk_rindex][2] = (float32) ray_head->pulse_count;
 457                         /* Store prf here instead of Sweep Rate */
 458                         rayInfoFloat[tk_sindex][tk_rindex][3] = (float32) ray_head->prf;
 459                 } /* end (rindex=0;... */
 460         } /* end for (tk_sindex=0;...*/
 461 }
 462
 463 /*************************************************************/
 464 /*                                                           */
 465 /*                      sweepInfoFill                        */
 466 /*                                                           */
 467 /*************************************************************/
 468 void sweepInfoFill(SWEEP_INFO sweepInfo[MAX_SWEEP], Radar *radar, VosSize *vs)
 469 {
 470         int rindex, tk_sindex;
 471         Sweep *sweep;
 472
 473         /* For each physical sweep...*/
 474         for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 475         {
 476                 sweep = vs->rsl.sweep[tk_sindex];
 477                 /* Fill in sweep info from the 1st non-NULL rsl ray we find. */
 478                 for (rindex=0; rindex<sweep->h.nrays; rindex++)
 479                 {
 480                         if (sweep->ray[rindex] == NULL) continue;
 481                         strncpy(sweepInfo[tk_sindex].radarName, radar->h.radar_name, 8);
 482                         /* 1st sweep number (for tsdis structures) is 1 */
 483                         sweepInfo[tk_sindex].sweepNum = (int32) (tk_sindex + 1);
 484                         /* No. of rays in sweep */
 485                         sweepInfo[tk_sindex].numRays = (int32) vs->tk.nray[tk_sindex];
 486                         sweepInfo[tk_sindex].trueStartAngle = (float32)
 487                                                                sweep->ray[rindex]->h.azimuth;
 488                         if (sweep->ray[sweep->h.nrays-1] != NULL)
 489                           sweepInfo[tk_sindex].trueStopAngle = (float32)
 490                                              sweep->ray[sweep->h.nrays-1]->h.azimuth;
 491                         else
 492                           sweepInfo[tk_sindex].trueStopAngle = (float32)
 493                                                               sweepInfo[tk_sindex].trueStartAngle;
 494                         /* degrees. Only for PPI scans. */
 495                         sweepInfo[tk_sindex].fixedAngle = (float32) sweep->h.elev;
 496                         /* Filter Flag. 0: No filtering, 1: filtered (descr in comment block) */
 497 /*                      sweepInfo[tk_sindex].filterFlag = (int32) 0; */
 498                         break;
 499                 }
 500         } /* end for (tk_sindex=0;... */
 501 }
 502
 503 /*************************************************************/
 504 /*                                                           */
 505 /*                      radarDescFill                        */
 506 /*                                                           */
 507 /*************************************************************/
 508 void radarDescFill(RADAR_DESCRIPTOR *radarDesc, Radar *radar, int vindex,
 509                                                                          VosSize *vs)
 510 {
 511         strncpy(radarDesc->radarName, radar->h.name, 8);
 512 /*radarDesc->radarConstant = (float32) 0.0;
 513         radarDesc->nomPeakPower = (float32) 0.0;
 514         radarDesc->nomNoisePower = (float32) 0.0;
 515         radarDesc->receiverGain = (float32) 0.0;
 516         radarDesc->antennaGain = (float32) 0.0;
 517         radarDesc->radarSystemGain = (float32) 0.0;
 518 */
 519         radarDesc->horBeamWidth = (float32) first_ray_in_volume[vindex]->h.beam_width;
 520         radarDesc->verBeamWidth = (float32) first_ray_in_volume[vindex]->h.beam_width;
 521         radarDesc->radarType = (int16) 0;            /* 0: Ground-based radar */
 522         radarDesc->scanMode = (int16) 1;             /* 1: PPI  ,  3: RHI  */
 523         /* radar sweep rate (deg/sec) */
 524         radarDesc->nomScanRate = (float32)
 525                    (first_ray_in_volume[vindex]->h.sweep_rate * 6.0);
 526         /* Following holds only for PPI scans. ??? */
 527         radarDesc->nomStartAngle = (float32)
 528                    first_ray_in_volume[vindex]->h.azimuth;
 529         radarDesc->nomStopAngle = (float32) radarDesc->nomStartAngle;
 530         radarDesc->numParmDesc = (int16) vs->tk.nparm;
 531         radarDesc->numDesc = (int16) vs->tk.nparm;
 532         /* Data compression. Always 0 . Data compression done by HDF.*/
 533         radarDesc->dataComp = (int16) 0;
 534         /* Data reduction algorithm. */
 535         radarDesc->dataReductAlg = (int16) 0;          /* 0: No reduction */
 536         radarDesc->dataReductParm1 = (float32) 4.0;    /* TBD */
 537         radarDesc->dataReductParm2 = (float32) 4.0;    /* TBD */
 538         radarDesc->radarLon = (float32) (radar->h.lond + radar->h.lonm/60.0 +
 539                                                                                                                                          radar->h.lons/3600.0);
 540         radarDesc->radarLat = (float32) (radar->h.latd + radar->h.latm/60.0 +
 541                                                                                                                                          radar->h.lats/3600.0);
 542         /* altitude above Mean Sea Level (km) */
 543         radarDesc->radarAlt = (float32)((float)(radar->h.height) / 1000.0);
 544         /* Effective unambiguous velocity (m/s) Leave blank. See range comments
 545            below. */
 546 /*
 547         if (radar->v[VR_INDEX] != NULL)
 548           radarDesc->velocity = (float32)
 549                                  first_ray_in_volume[VR_INDEX]->h.nyq_vel;
 550         else
 551           radarDesc->velocity = (float32)0.0;
 552 */
 553         /* Effective unambiguous range (km). Leave blank. For wsr88d,
 554                  unambig_range varies with sweep. See range_info_float block. */
 555 /*
 556         radarDesc->range = (float32)
 557                             first_ray_in_volume[vindex]->h.unam_rng;
 558 */
 559         /* No. of transmitted freqs */
 560         radarDesc->numTransfreqency = (int16) 1;
 561 /*      radarDesc->numInterPulsePeriods = (int16) 0; */
 562         /* Freq. GHz. */
 563         radarDesc->frequency1 = (float32) first_ray_in_volume[vindex]->h.frequency;
 564 }
 565
 566 /*************************************************************/
 567 /*                                                           */
 568 /*                      sensorFill                           */
 569 /*                                                           */
 570 /*************************************************************/
 571 void sensorFill(SENSORS *sensor, Radar *radar, VosSize *vs, int productID)
 572 /* Fill the substructures of the sensor data structure using
 573          data from the radar structure. */
 574 {
 575         int pindex, tk_sindex, vindex;
 576
 577         vindex = nextVolume(radar, -1);  /* Find 1st non-NULL rsl volume. */
 578
 579         radarDescFill(&sensor->radarDesc, radar, vindex, vs);
 580         sweepInfoFill(sensor->sweepInfo, radar, vs);
 581         rayInfoFill(sensor->rayInfoInteger, sensor->rayInfoFloat, radar, vs);
 582
 583         /* Move data from each of the radar structure volumes into a
 584                  corresponding L1GV parameter structure. */
 585         if (productID == TK_L1B_GV)
 586         {
 587           for (pindex=0; pindex<vs->tk.nparm; pindex++)
 588                 {
 589                         sensor->parm[pindex] = (PARAMETER       *)parmBuild(radar, vs, vindex, pindex);
 590                         /* Find the next non-NULL volume in radar structure. */
 591                         vindex = nextVolume(radar, vindex);
 592                 }  /* end for (pindex=0; ... */
 593         }
 594         else  /* 1C-51 */
 595         {
 596                 /* This is a hatchet job to conform with newest toolkit. The toolkit
 597                          arbitrarily assumes that the mask volume precedes the corresponding
 598                          data volume.
 599                 */
 600           for (pindex=0; pindex<vs->tk.nparm/2; pindex++)
 601                 {
 602                         if (vindex == CZ_INDEX)
 603                           sensor->parm[pindex] = (PARAMETER *)parmBuild(radar, vs, MZ_INDEX,
 604                                                                                                                                                                                                                  pindex);
 605                         else if (vindex == CD_INDEX)
 606                           sensor->parm[pindex*2] = (PARAMETER *)parmBuild(radar, vs, MD_INDEX,
 607                                                                                                                                                                                                                          pindex);
 608                         else
 609                           continue;
 610                         sensor->parm[pindex*2+1] = (PARAMETER *)parmBuild(radar, vs, vindex,
 611                                                                                                                                                                                                                          pindex);
 612                         /* Find the next non-NULL volume in radar structure. */
 613                         vindex = nextVolume(radar, vindex);
 614                 }  /* end for (pindex=0... */
 615         }  /* end else 1C-51 */
 616
 617         if (radar_verbose_flag)
 618         {
 619                 for (pindex=0; pindex<vs->tk.nparm; pindex++)
 620                 {
 621                         fprintf(stderr, "Toolkit parameter type : %s '%s'\n",
 622                                                         sensor->parm[pindex]->parmDesc.parmDesc,
 623                                                         sensor->parm[pindex]->parmDesc.parmName);
 624                   for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 625                           if (vs->tk.ncell[tk_sindex][pindex] == 0)
 626                                   fprintf(stderr, "  tk_sweep[%.2d]  elev=%4.1f  nrays=%3d  cells/ray=%d\n",
 627                                                                         tk_sindex, sensor->sweepInfo[tk_sindex].fixedAngle,
 628                                                                         (int)0,
 629                                                                         vs->tk.ncell[tk_sindex][pindex]);
 630                                 else
 631                                   fprintf(stderr, "  tk_sweep[%.2d]  elev=%4.1f  nrays=%3d  cells/ray=%d\n",
 632                                                                         tk_sindex, sensor->sweepInfo[tk_sindex].fixedAngle,
 633                                                                         (int)sensor->sweepInfo[tk_sindex].numRays,
 634                                                                         vs->tk.ncell[tk_sindex][pindex]);
 635                 } /* end for (pindex=0;... */
 636         } /* end if (radar_verbose_flag) */
 637 }
 638
 639 /*************************************************************/
 640 /*                                                           */
 641 /*                  timeUTC                                  */
 642 /*                                                           */
 643 /*************************************************************/
 644 struct tm *timeUTC(void)
 645 {
 646 /* Find the current time (UTC).
 647          If success: return pointer to filled time_t structure.
 648          If failure: return NULL.
 649 */
 650         time_t time_current;
 651
 652         /* Get the current system clock time */
 653         time_current = (time_t) time(NULL);
 654         if (time_current != -1)      /* valid time? */
 655                 return(gmtime(&time_current));  /* Convert to UTC and return. */
 656         else
 657           return(NULL);         /* Couldn't get the current time. */
 658 }
 659
 660
 661 /*************************************************************/
 662 /*                                                           */
 663 /*                      volDesFill                           */
 664 /*                                                           */
 665 /*************************************************************/
 666 void volDesFill(VOLUME_DESCRIPTORS *volDes, Radar_header *h, VosSize *vs)
 667 {
 668         struct tm *time_utc;
 669         int32 max_ncell;
 670         int pindex;
 671
 672         /* Version no. of DORADE specifications used. Currently 1 */
 673         volDes->verNum = (int16) 1;
 674         /* No. of this volume scan in granule */
 675         volDes->volNum = (int16) (vs->vos_num + 1);
 676         /* Max size of DORADE data record in this VOS. 2bytes x max(ncell) */
 677         max_ncell = 0;
 678         for (pindex=0; pindex<vs->tk.nparm; pindex++)
 679           if (vs->tk.ncell[0][pindex] > max_ncell)
 680             max_ncell = vs->tk.ncell[0][pindex];
 681         volDes->sizeDataRec = (int32) (2 * max_ncell);
 682         strncpy(volDes->projectName, "TRMM GV", 20);
 683         volDes->year = (int16) h->year;             /* Year of volume scan */
 684         volDes->month = (int16) h->month;    /* Month of volume scan */
 685         volDes->day = (int16) h->day;        /* Day ... */
 686         volDes->hour = (int16) h->hour;      /* Hour ... */
 687         volDes->minute = (int16) h->minute;  /* Minute ... */
 688         volDes->second = (int16) floor((double)h->sec);  /* Second ... */
 689         /* Flight no. for airborne radar, or IOP no. for ground radar */
 690         strncpy(volDes->flightNum, "***", 8);
 691         /* Data product generation facility name */
 692         strncpy(volDes->facName, "TSDIS", 8);
 693         /* Get the current time; ie, the time of creation of this hdf file. */
 694         time_utc = timeUTC();
 695         if (time_utc != NULL)   /* Valid time? */
 696         {
 697                 volDes->genYear = (int16) (1900 + time_utc->tm_year);
 698                 volDes->genMonth = (int16) (time_utc->tm_mon + 1);
 699                 volDes->genDay = (int16) time_utc->tm_mday;
 700         }
 701         else  /* Couldn't get valid time. */
 702         {
 703                 volDes->genYear = (int16) 0;
 704                 volDes->genMonth = (int16) 0;
 705                 volDes->genDay = (int16) 0;
 706         }
 707         /* No. of sensor descriptors in this volume scan */
 708         volDes->numSensorDesc = (int16) 1;     /* 1 for ground-based radar */
 709 }
 710
 711 /*************************************************************/
 712 /*                                                           */
 713 /*                      commentsFill                         */
 714 /*                                                           */
 715 /*************************************************************/
 716 void commentsFill(char *comments, VosSize *vs, Radar *radar,
 717                                                                         float *qcParm, int productID)
 718 {
 719 /* Write the following into the comments field of the gvl1
 720          structure:
 721            1. VOS comment_field header line.
 722            2. cell/ray/sweep count.
 723            3. 1C-51 QC parameters, if we're writing a 1C-51 HDF file.
 724 */
 725         char buf[256];
 726         int pindex, tk_sindex;
 727
 728         /* Write out the dimensions of the toolkit structure which contains
 729                  this VOS. */
 730         sprintf(buf, "nSweep=%d\n", vs->tk.nsweep);
 731         strcat(comments, buf);
 732         for (tk_sindex=0; tk_sindex<vs->tk.nsweep; tk_sindex++)
 733         {
 734                 sprintf(buf, "sweep[%.2d]--\n nRay=%d\n", tk_sindex,
 735                                                 vs->tk.nray[tk_sindex]);
 736                 strcat(comments, buf);
 737           for (pindex=0; pindex<vs->tk.nparm; pindex++)
 738                 {
 739                         sprintf(buf, "  nCell_parm[%d]=%d\n", pindex,
 740                                                         vs->tk.ncell[tk_sindex][pindex]);
 741                         strcat(comments, buf);
 742                 }
 743                 strcat(comments, "\n");
 744         } /* end for (pindex=0 ... */
 745
 746         strcat(comments, "********\n");
 747         /* If 1C-51 file, write out the QC parameters. */
 748         if (productID == TK_L1C_GV)
 749         {
 750                 sprintf(buf, "-hThresh1 %.2f  -hThresh2 %.2f  -hThresh3 %.2f  -zThresh0 %.2f  -zThresh1 %.2f  -zThresh2 %.2f  -zThresh3 %.2f  -hFreeze %.2f  -dbzNoise %.2f  -zCal %.2f\n\n",
 751                                                 qcParm[HTHRESH1], qcParm[HTHRESH2], qcParm[HTHRESH3],
 752                                                 qcParm[ZTHRESH0], qcParm[ZTHRESH1], qcParm[ZTHRESH2], qcParm[ZTHRESH3],
 753                                                 qcParm[HFREEZE], qcParm[DBZNOISE], qcParm[ZCAL]);
 754                 strcat(comments, buf);
 755         } /* end if (productID == TK_L1C_GV) */
 756 }
 757
 758 /*************************************************************/
 759 /*                                                           */
 760 /*                        gvl1Build                          */
 761 /*                                                           */
 762 /*************************************************************/
 763 L1B_1C_GV *gvl1Build(Radar *radar, float *qcParm, VosSize *vs,
 764                                                                                  int productID)
 765 {
 766 /* Build the components of the Toolkit 'L1B_1C_GV' structure using
 767          the data from the RSL radar structure.
 768 */
 769         L1B_1C_GV *gvl1;
 770
 771         /* Allocate memory for a TSDIS level_1 structure. */
 772         gvl1 = (L1B_1C_GV *)TKnewGVL1();
 773
 774         /* Fill the structure, using data from the radar structure. */
 775         commentsFill(gvl1->comments, vs, radar, qcParm, productID);
 776         volDesFill(&gvl1->volDes, &radar->h, vs);
 777         sensorFill(&gvl1->sensor, radar, vs, productID);
 778         return(gvl1);
 779 }
 780
 781 #endif