jDttSP/chap.c

   1 /* chap.c
   2
   3 This file is part of a program that implements a Software-Defined Radio.
   4
   5 Copyright (C) 2004 by Frank Brickle, AB2KT and Bob McGwier, N4HY
   6
   7 This program is free software; you can redistribute it and/or modify
   8 it under the terms of the GNU General Public License as published by
   9 the Free Software Foundation; either version 2 of the License, or
  10 (at your option) any later version.
  11
  12 This program is distributed in the hope that it will be useful,
  13 but WITHOUT ANY WARRANTY; without even the implied warranty of
  14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15 GNU General Public License for more details.
  16
  17 You should have received a copy of the GNU General Public License
  18 along with this program; if not, write to the Free Software
  19 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  20
  21 The authors can be reached by email at
  22
  23 ab2kt@arrl.net
  24 or
  25 rwmcgwier@comcast.net
  26
  27 or by paper mail at
  28
  29 The DTTS Microwave Society
  30 6 Kathleen Place
  31 Bridgewater, NJ 08807
  32 */
  33
  34 #include <chap.h>
  35
  36 REAL
  37 ChAp_eval(ChAp ca, REAL x, BOOLEAN *err) {
  38   int i;
  39   REAL d, e, u, v;
  40
  41   if (x < ChApLob(ca)) {
  42     x = ChApLob(ca);
  43     if (err) *err = TRUE;
  44   } else if (x > ChApHib(ca)) {
  45     x = ChApHib(ca);
  46     if (err) *err = TRUE;
  47   } else if (err) *err = FALSE;
  48
  49   d = e = 0.0;
  50   u = (2.0 * x - ChApDif(ca)) / ChApDif(ca);
  51   v = 2.0 * u;
  52   for (i = ChApLen(ca) - 1; i > 0; --i) {
  53     REAL sv = d;
  54     d = v * d - e + ChApCoef(ca, i);
  55     e = sv;
  56   }
  57
  58   return u * d - e + 0.5 * ChApCoef(ca, 0);
  59 }
  60
  61 ChAp
  62 ChAp_fit(ChAp ca) {
  63   int i, j, n = ChApLen(ca);
  64   REAL (*func)(REAL) = ChApFunc(ca),
  65        bma = 0.5 * (ChApHib(ca) - ChApLob(ca)),
  66        bpa = 0.5 * (ChApHib(ca) + ChApLob(ca)),
  67        fac = 2.0 / n,
  68       *tbl = newvec_REAL(n, ChApTag(ca));
  69   for (j = 0; j < n; j++) {
  70     REAL y = cos(M_PI * (j + 0.5) / n);
  71     tbl[j] = (*func)(y * bma + bpa);
  72   }
  73   for (i = 0; i < n; i++) {
  74     REAL sum = 0.0;
  75     for (j = 0; j < n; j++)
  76       sum += tbl[j] * cos((M_PI * i) * ((j + 0.5) / n));
  77     ChApCoef(ca, i) = fac * sum;
  78   }
  79   delvec_REAL(tbl);
  80   return ca;
  81 }
  82
  83 ChAp
  84 newChAp(REAL (*func)(REAL arg),
  85         REAL lo,
  86         REAL hi,
  87         int len,
  88         char *tag) {
  89   ChAp ca = (ChAp) safealloc(1, sizeof(ChApDesc), tag);
  90   ChApFunc(ca) = func;
  91   ChApLen(ca) = len;
  92   ChApLob(ca) = lo;
  93   ChApHib(ca) = hi;
  94   ChApDif(ca) = hi - lo;
  95   ChApCoefBase(ca) = newvec_REAL(len, tag);
  96   ChApTag(ca) = strdup(tag);
  97   return ChAp_fit(ca);
  98 }
  99
 100 void
 101 delChAp(ChAp ca) {
 102   if (ca) {
 103     safefree((char *) ChApCoefBase(ca));
 104     safefree(ChApTag(ca));
 105     safefree((char *) ca);
 106   }
 107 }
 108
 109