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

   1 /*
   2     NASA/TRMM, Code 910.1.
   3     This is the TRMM Office Radar Software Library.
   4     Copyright (C) 1997
   5             John H. Merritt
   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 #include <stdio.h>
  24 #include <stdlib.h>
  25 #include <unistd.h>
  26 #include "rsl.h"
  27 #include "africa.h"
  28
  29 static void ymd(int jday, int yy, int *mm, int *dd);
  30
  31 static int daytab[2][13] = {
  32   {0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334, 365},
  33   {0, 31, 60, 91, 121, 152, 182, 213, 244, 274, 305, 335, 366}
  34 };
  35
  36 static void ymd(int jday, int year, int *mm, int *dd)
  37 {
  38   /*  Input: jday, yyyy */
  39   /* Output: mm, dd */
  40   int leap;
  41   int i;
  42
  43   leap = (year%4 == 0 && year%100 != 0) || year%400 == 0;
  44   for (i=0; daytab[leap][i]<jday; i++) continue;
  45   *mm = i;
  46   i--;
  47   *dd = jday - daytab[leap][i];
  48 }
  49
  50 extern int radar_verbose_flag;
  51
  52 Radar *RSL_africa_to_radar(char *infile)
  53 {
  54
  55   FILE *fp;
  56   Africa_ray *ray;
  57   Africa_sweep *sweep;
  58   int n;
  59   int hour, min=0, sec=0, month, day, year, jday;
  60   float elev, azim, s_azim, e_azim;
  61   int iray, ielev, isite;
  62
  63   Radar *radar;
  64   Volume *v;
  65   Sweep  *s;
  66   Ray    *r;
  67   int     i,ibin;
  68   int save_fd;
  69
  70   if (infile == NULL) {
  71         save_fd = dup(0);
  72         fp = fdopen(save_fd, "r");
  73   } else {
  74         fp = fopen(infile, "r");
  75   }
  76   fp = uncompress_pipe(fp);
  77   n = 0;
  78
  79   radar = RSL_new_radar(MAX_RADAR_VOLUMES);
  80   radar->v[DZ_INDEX] = RSL_new_volume(20);
  81   v = radar->v[DZ_INDEX];
  82
  83   for(i=0; (sweep = africa_read_sweep(fp)); i++) {
  84
  85         /* Load the sweep into the radar volume */
  86         v->sweep[i] = RSL_new_sweep((int)sweep->nrays);
  87         s = v->sweep[i];
  88         if (radar_verbose_flag) printf("NUMBER OF RAYS: %d\n", sweep->nrays);
  89         for (n=0; n < sweep->nrays; n++) {
  90           ray = sweep->ray[n];
  91           if (ray == NULL) continue;
  92           year = ray->yearjday/512;
  93           jday = ray->yearjday & 0x1ff;
  94           year += 1900;
  95           if (year < 1970) year += 100; /* >=2000 */
  96           ymd(jday, year, &month, &day);
  97           hour = ray->hour;
  98           min  = ray->minute_sec/60;
  99           sec  = ray->minute_sec - min*60;
 100
 101           elev = africa_bcd_convert (ray->bcd_elevation);
 102           azim = africa_bcd_convert (ray->bcd_azimuth);
 103           s_azim = africa_bcd_convert (ray->bcd_start_azim);
 104           e_azim = africa_bcd_convert (ray->bcd_end_azim);
 105
 106           /* Values that I CANNOT trust.
 107            *    e_azim
 108            *
 109            * DATA ORGANIZATION:
 110            *   Usually, 10 ray groups for each field type.  Field type is in
 111            *            ray->xmit_power_site (Doc wrong?)
 112            *   iray == 0 when at the end of the sweep.  (Don't need lookahead)
 113            */
 114
 115           ielev = ray->raycount/512;
 116           iray  = ray->raycount - ielev*512;
 117           isite = ray->xmit_power_site & 0x1f;
 118           printf("Record %d, time = %.2d:%.2d:%.2d %.2d/%.2d/%.2d  --> elev,azim = %f, %f.  Start: %f, %f, iray/ielev %d %d, site=%d\n", n,
 119                          hour, min, sec, month, day, year-1900, elev, azim, s_azim, e_azim, iray, ielev, isite);
 120
 121
 122           if (isite != 22) continue;
 123           /* ONLY LOAD from isite==22, for now. */
 124           r = RSL_new_ray(224);
 125           s->ray[n]  = r;
 126           r->h.month = month;
 127           r->h.day   = day;
 128           r->h.year  = year;
 129           r->h.hour  = hour;
 130           r->h.minute  = min;
 131           r->h.sec   = sec;
 132           r->h.azimuth = azim;
 133           r->h.ray_num = iray;
 134           r->h.elev    = elev;
 135           r->h.elev_num = ielev;
 136           r->h.beam_width = 1.0;
 137           r->h.gate_size = 1000;
 138
 139           /* What are the others? */
 140
 141           r->h.invf = DZ_INVF;
 142           r->h.f    = DZ_F;
 143           s->h.horz_half_bw = 0.5;
 144           s->h.vert_half_bw = 0.5;
 145           s->h.beam_width = r->h.beam_width;
 146
 147           for (ibin=0; ibin<r->h.nbins; ibin++) {
 148                 r->range[ibin] = r->h.invf(ray->bin[ibin]/8.0/100.0 + 0.5);
 149           }
 150         }
 151   }
 152   fclose(fp);
 153
 154   sprintf(radar->h.radar_type, "south_africa");
 155   sprintf(radar->h.name, "SAFRICA");
 156   sprintf(radar->h.radar_name, "SAFRICA");
 157   sprintf(radar->h.city, "I don't know");
 158   sprintf(radar->h.state, "??");
 159   sprintf(radar->h.country, "South Africa");
 160
 161   radar_load_date_time(radar);
 162   return radar;
 163 }