simgear/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  - http://www.flightgear.org/~curt
   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 General Public License
  20 // along with this program; if not, write to the Free Software
  21 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
  22 //
  23 // $Id$
  24
  25
  26 #include <simgear/compiler.h>
  27
  28 #include STL_STRING
  29
  30 #include <simgear/debug/logstream.hxx>
  31 #include <simgear/misc/sgstream.hxx>
  32
  33 #include "interpolater.hxx"
  34
  35 SG_USING_STD(string);
  36
  37 // Constructor -- starts with an empty table.
  38 SGInterpTable::SGInterpTable()
  39 {
  40 }
  41
  42
  43 // Constructor -- loads the interpolation table from the specified
  44 // file
  45 SGInterpTable::SGInterpTable( const string& file )
  46 {
  47     SG_LOG( SG_MATH, SG_INFO, "Initializing Interpolator for " << file );
  48
  49     sg_gzifstream in( file );
  50     if ( !in.is_open() ) {
  51         SG_LOG( SG_GENERAL, SG_ALERT, "Cannot open file: " << file );
  52         return;
  53     }
  54
  55     in >> skipcomment;
  56     while ( in ) {
  57       double ind, dep;
  58       in >> ind >> dep;
  59       in >> skipws;
  60       _table[ind] = dep;
  61     }
  62 }
  63
  64
  65 // Add an entry to the table.
  66 void SGInterpTable::addEntry (double ind, double dep)
  67 {
  68   _table[ind] = dep;
  69 }
  70
  71 // Given an x value, linearly interpolate the y value from the table
  72 double SGInterpTable::interpolate(double x) const
  73 {
  74   // Empty table??
  75   if (_table.empty())
  76     return 0;
  77
  78   // Find the table bounds for the requested input.
  79   Table::const_iterator upBoundIt = _table.upper_bound(x);
  80   // points to a value outside the map. That is we are out of range.
  81   // use the last entry
  82   if (upBoundIt == _table.end())
  83     return _table.rbegin()->second;
  84
  85   // points to the first key must be lower
  86   // use the first entry
  87   if (upBoundIt == _table.begin())
  88     return upBoundIt->second;
  89
  90   // we know that we do not stand at the beginning, so it is safe to do so
  91   Table::const_iterator loBoundIt = upBoundIt;
  92   --loBoundIt;
  93
  94   // Just do linear interpolation.
  95   double loBound = loBoundIt->first;
  96   double upBound = upBoundIt->first;
  97   double loVal = loBoundIt->second;
  98   double upVal = upBoundIt->second;
  99
 100   // division by zero should not happen since the std::map
 101   // has sorted out duplicate entries before. Also since we have a
 102   // map, we know that we get different first values for different iterators
 103   return loVal + (upVal - loVal)*(x - loBound)/(upBound - loBound);
 104 }
 105
 106
 107 // Destructor
 108 SGInterpTable::~SGInterpTable() {
 109 }
 110
 111