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

   1 /*
   2     NASA/TRMM, Code 910.1.
   3     This is the TRMM Office Radar Software Library.
   4     Copyright (C) 1996, 1997
   5             David B. Wolff
   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 "rsl.h"
  24 /*
  25  * Author: David B. Wolff
  26  * Date:   8/4/94
  27  */
  28 /**************************************************************************
  29         Ideally this would be done via a CAPPI, but could also be done via
  30         the PPI. Here we will send a volume to produce a CAPPI. From the CAPPI(Z)
  31         we compute the CAPPI(R) given a Z-R relationship. Then convert that to
  32         a cartesian array for imaging or later fiddling.
  33
  34     NOTE:
  35         In keep with the object oriented approach to the library, the following
  36     approach is used to deal with the Z-R conversion:
  37
  38         Four functions are produced: RSL_volume_Z_to_R, RSL_sweep_Z_to_R,
  39     RSL_ray_Z_to_r, z_to_r. RSL_volume_Z_to_R call RSL_sweep_Z_to_R for each
  40     sweep, RSL_sweep_Z_to_R calls RSL__ray_Z_to_R for each ray, RSL__ray_Z_to_R
  41     calls z_to_r for each range bin.
  42
  43 ***************************************************************************/
  44
  45 Volume  *RSL_volume_z_to_r(Volume *z_volume, float k, float a);
  46 Sweep   *RSL_sweep_z_to_r(Sweep *z_sweep, float k, float a);
  47 Ray     *RSL_ray_z_to_r(Ray *z_ray, float k, float a);
  48 float   RSL_z_to_r(float z, float k, float a);
  49
  50
  51 Volume *RSL_volume_z_to_r(Volume *z_volume, float k, float a)
  52 {
  53         Volume  *r_volume;
  54         int             i;
  55         if(z_volume == NULL) return NULL;
  56         r_volume = RSL_new_volume(z_volume->h.nsweeps);
  57         r_volume->h = z_volume->h;
  58         for(i=0; i<z_volume->h.nsweeps; i++) {
  59                 r_volume->sweep[i] = RSL_sweep_z_to_r(z_volume->sweep[i], k, a);
  60         }
  61         return r_volume;
  62 }
  63
  64 Sweep *RSL_sweep_z_to_r(Sweep *z_sweep, float k, float a)
  65 {
  66         Sweep   *r_sweep;
  67         int             i;
  68         if(z_sweep == NULL) return NULL;
  69         r_sweep = RSL_new_sweep(z_sweep->h.nrays);
  70         r_sweep->h = z_sweep->h;
  71         for(i=0; i<z_sweep->h.nrays; i++) {
  72                 r_sweep->ray[i] = RSL_ray_z_to_r(z_sweep->ray[i], k, a);
  73         }
  74         return r_sweep;
  75 }
  76
  77 Ray *RSL_ray_z_to_r(Ray *z_ray, float k, float a)
  78 {
  79         Ray     *r_ray;
  80         int             i;
  81         Range (*invf)(float x);
  82         float (*f)(Range x);
  83
  84         if (z_ray == NULL) return NULL;
  85         r_ray = RSL_new_ray(z_ray->h.nbins);
  86         r_ray->h = z_ray->h;
  87         f = r_ray->h.f;
  88         invf = r_ray->h.invf;
  89         for(i=0; i<z_ray->h.nbins; i++) {
  90                 r_ray->range[i] = invf(RSL_z_to_r( f(z_ray->range[i]), k, a));
  91         }
  92         return r_ray;
  93
  94 }
  95
  96 #include <math.h>
  97
  98 float RSL_z_to_r(float dbz, float k, float a) {
  99
 100         float   dbr, dbk, r;
 101
 102         if (dbz >= NOECHO) return dbz; /* Be careful, NOECHO is the smallest,
 103                                                                           so far, of the reserved numbers. */
 104
 105         dbk = (float)10*log10((double)k);
 106         dbr = 1./a*(dbz- dbk);
 107         r   = (float)pow((double)10., (double)(dbr/10.) );
 108 /*      fprintf(stderr,"dbz=%.1f; \tdbr=%.1f \tRate= %.3f\n",dbz,dbr,r);  */
 109         return r;
 110 }
 111