Math/interpolater.cxx

   1 //
   2 // interpolater.cxx -- routines to handle linear interpolation from a table of
   3 //                     x,y   The table must be sorted by "x" in ascending order
   4 //
   5 // Written by Curtis Olson, started April 1998.
   6 //
   7 // Copyright (C) 1998  Curtis L. Olson  - curt@me.umn.edu
   8 //
   9 // This program is free software; you can redistribute it and/or
  10 // modify it under the terms of the GNU General Public License as
  11 // published by the Free Software Foundation; either version 2 of the
  12 // License, or (at your option) any later version.
  13 //
  14 // This program is distributed in the hope that it will be useful, but
  15 // WITHOUT ANY WARRANTY; without even the implied warranty of
  16 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17 // General Public License for more details.
  18 //
  19 // You should have received a copy of the GNU General Public License
  20 // along with this program; if not, write to the Free Software
  21 // Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  22 //
  23 // $Id$
  24 // (Log is kept at end of this file)
  25
  26
  27 #include <string.h>
  28
  29 #include <Debug/fg_debug.h>
  30 #include <zlib/zlib.h>
  31
  32 #include "interpolater.hxx"
  33
  34
  35 // Constructor -- loads the interpolation table from the specified
  36 // file
  37 fgINTERPTABLE::fgINTERPTABLE( char *file ) {
  38     char gzfile[256], line[256];
  39     gzFile fd;
  40
  41     fgPrintf( FG_MATH, FG_INFO, "Initializing Interpolator for %s\n", file);
  42
  43     // First try "file.gz"
  44     strcpy(gzfile, file);
  45     strcat(gzfile, ".gz");
  46     if ( (fd = gzopen(gzfile, "r")) == NULL ) {
  47         // Next try "path"
  48         if ( (fd = gzopen(file, "r")) == NULL ) {
  49             fgPrintf(FG_MATH, FG_ALERT, "Cannot open file: %s\n", file);
  50         }
  51     }
  52
  53     size = 0;
  54     while ( gzgets(fd, line, 250) != NULL ) {
  55         if ( size < MAX_TABLE_SIZE ) {
  56             sscanf(line, "%lf %lf\n", &(table[size][0]), &(table[size][1]));
  57             size++;
  58         } else {
  59             fgPrintf( FG_MATH, FG_ALERT,
  60                       "fgInterpolateInit(): Exceed max table size = %d\n",
  61                       MAX_TABLE_SIZE );
  62         }
  63     }
  64
  65     gzclose(fd);
  66 }
  67
  68
  69 // Given an x value, linearly interpolate the y value from the table
  70 double fgINTERPTABLE::interpolate(double x) {
  71     int i;
  72     double y;
  73
  74     i = 0;
  75
  76     while ( (x > table[i][0]) && (i < size) ) {
  77         i++;
  78     }
  79
  80     printf ("i = %d ", i);
  81
  82     if ( (i == 0) && (x < table[0][0]) ) {
  83         fgPrintf( FG_MATH, FG_ALERT,
  84                   "fgInterpolateInit(): lookup error, x to small = %.2f\n", x);
  85         return(0.0);
  86     }
  87
  88     if ( x > table[i][0] ) {
  89         fgPrintf( FG_MATH, FG_ALERT,
  90                   "fgInterpolateInit(): lookup error, x to big = %.2f\n",  x);
  91         return(0.0);
  92     }
  93
  94     // y = y1 + (y0 - y1)(x - x1) / (x0 - x1)
  95     y = table[i][1] +
  96         ( (table[i-1][1] - table[i][1]) *
  97           (x - table[i][0]) ) /
  98         (table[i-1][0] - table[i][0]);
  99
 100     return(y);
 101 }
 102
 103
 104 // Destructor
 105 fgINTERPTABLE::~fgINTERPTABLE( void ) {
 106 }
 107
 108
 109 // $Log$
 110 // Revision 1.1  1998/04/21 19:14:23  curt
 111 // Modified Files:
 112 //     Makefile.am Makefile.in
 113 // Added Files:
 114 //     interpolater.cxx interpolater.hxx
 115 //
 116 //