]> git.mxchange.org Git - flightgear.git/blob - src/FDM/JSBSim/math/FGTable.cpp
Merge branch 'jmt/units-fix' into maint
[flightgear.git] / src / FDM / JSBSim / math / FGTable.cpp
1 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2
3  Module:       FGTable.cpp
4  Author:       Jon S. Berndt
5  Date started: 1/9/2001
6  Purpose:      Models a lookup table
7
8  ------------- Copyright (C) 2001  Jon S. Berndt (jon@jsbsim.org) -------------
9
10  This program is free software; you can redistribute it and/or modify it under
11  the terms of the GNU Lesser General Public License as published by the Free Software
12  Foundation; either version 2 of the License, or (at your option) any later
13  version.
14
15  This program is distributed in the hope that it will be useful, but WITHOUT
16  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17  FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
18  details.
19
20  You should have received a copy of the GNU Lesser General Public License along with
21  this program; if not, write to the Free Software Foundation, Inc., 59 Temple
22  Place - Suite 330, Boston, MA  02111-1307, USA.
23
24  Further information about the GNU Lesser General Public License can also be found on
25  the world wide web at http://www.gnu.org.
26
27 FUNCTIONAL DESCRIPTION
28 --------------------------------------------------------------------------------
29 Models a lookup table
30
31 HISTORY
32 --------------------------------------------------------------------------------
33 JSB  1/9/00          Created
34
35 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
36 INCLUDES
37 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
38
39 #include "FGTable.h"
40 #include "input_output/FGXMLElement.h"
41 #include "input_output/FGPropertyManager.h"
42 #include <iostream>
43 #include <sstream>
44 #include <cstdlib>
45
46 using namespace std;
47
48 namespace JSBSim {
49
50 static const char *IdSrc = "$Id$";
51 static const char *IdHdr = ID_TABLE;
52
53 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
54 CLASS IMPLEMENTATION
55 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
56
57 FGTable::FGTable(int NRows) : nRows(NRows), nCols(1), PropertyManager(0)
58 {
59   Type = tt1D;
60   colCounter = 0;
61   rowCounter = 1;
62   nTables = 0;
63
64   Data = Allocate();
65   Debug(0);
66   lastRowIndex=lastColumnIndex=2;
67 }
68
69 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
70
71 FGTable::FGTable(const FGTable& t) : PropertyManager(t.PropertyManager)
72 {
73   Type = t.Type;
74   colCounter = t.colCounter;
75   rowCounter = t.rowCounter;
76   tableCounter = t.tableCounter;
77   nRows = t.nRows;
78   nCols = t.nCols;
79   nTables = t.nTables;
80   dimension = t.dimension;
81   internal = t.internal;
82   Name = t.Name;
83   lookupProperty[0] = t.lookupProperty[0];
84   lookupProperty[1] = t.lookupProperty[1];
85   lookupProperty[2] = t.lookupProperty[2];
86
87   Tables = t.Tables;
88   Data = Allocate();
89   for (unsigned int r=0; r<=nRows; r++) {
90     for (unsigned int c=0; c<=nCols; c++) {
91       Data[r][c] = t.Data[r][c];
92     }
93   }
94   lastRowIndex = t.lastRowIndex;
95   lastColumnIndex = t.lastColumnIndex;
96   lastTableIndex = t.lastTableIndex;
97 }
98
99 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
100
101 FGTable::FGTable(FGPropertyManager* propMan, Element* el) : PropertyManager(propMan)
102 {
103   unsigned int i;
104
105   stringstream buf;
106   string property_string;
107   string lookup_axis;
108   string call_type;
109   string parent_type;
110   string brkpt_string;
111   FGPropertyManager* node;
112   Element *tableData=0;
113   Element *parent_element=0;
114   Element *axisElement=0;
115   string operation_types = "function, product, sum, difference, quotient,"
116                            "pow, abs, sin, cos, asin, acos, tan, atan, table";
117
118   nTables = 0;
119
120   // Is this an internal lookup table?
121
122   internal = false;
123   Name = el->GetAttributeValue("name"); // Allow this table to be named with a property
124   call_type = el->GetAttributeValue("type");
125   if (call_type == string("internal")) {
126     parent_element = el->GetParent();
127     parent_type = parent_element->GetName();
128     if (operation_types.find(parent_type) == string::npos) {
129       internal = true;
130     } else {
131       // internal table is a child element of a restricted type
132       cerr << endl << fgred << "  An internal table cannot be nested within another type," << endl;
133       cerr << "  such as a function. The 'internal' keyword is ignored." << fgdef << endl << endl;
134     }
135   } else if (!call_type.empty()) {
136     cerr << endl << fgred << "  An unknown table type attribute is listed: " << call_type
137                  << ". Execution cannot continue." << fgdef << endl << endl;
138     abort();
139   }
140
141   // Determine and store the lookup properties for this table unless this table
142   // is part of a 3D table, in which case its independentVar property indexes will
143   // be set by a call from the owning table during creation
144
145   dimension = 0;
146
147   axisElement = el->FindElement("independentVar");
148   if (axisElement) {
149
150     // The 'internal' attribute of the table element cannot be specified
151     // at the same time that independentVars are specified.
152     if (internal) {
153       cerr << endl << fgred << "  This table specifies both 'internal' call type" << endl;
154       cerr << "  and specific lookup properties via the 'independentVar' element." << endl;
155       cerr << "  These are mutually exclusive specifications. The 'internal'" << endl;
156       cerr << "  attribute will be ignored." << fgdef << endl << endl;
157       internal = false;
158     }
159
160     for (i=0; i<3; i++) lookupProperty[i] = 0;
161
162     while (axisElement) {
163       property_string = axisElement->GetDataLine();
164       // The property string passed into GetNode() must have no spaces or tabs.
165       node = PropertyManager->GetNode(property_string);
166
167       if (node == 0) {
168         cerr << "IndependenVar property, " << property_string
169              << " in Table definition is not defined." << endl;
170         abort();
171       }
172
173       lookup_axis = axisElement->GetAttributeValue("lookup");
174       if (lookup_axis == string("row")) {
175         lookupProperty[eRow] = node;
176       } else if (lookup_axis == string("column")) {
177         lookupProperty[eColumn] = node;
178       } else if (lookup_axis == string("table")) {
179         lookupProperty[eTable] = node;
180       } else { // assumed single dimension table; row lookup
181         lookupProperty[eRow] = node;
182       }
183       dimension++;
184       axisElement = el->FindNextElement("independentVar");
185     }
186
187   } else if (internal) { // This table is an internal table
188
189   // determine how many rows, columns, and tables in this table (dimension).
190
191     if (el->GetNumElements("tableData") > 1) {
192       dimension = 3; // this is a 3D table
193     } else {
194       tableData = el->FindElement("tableData");
195       string test_line = tableData->GetDataLine(1);  // examine second line in table for dimension
196       if (FindNumColumns(test_line) == 2) dimension = 1;    // 1D table
197       else if (FindNumColumns(test_line) > 2) dimension = 2; // 2D table
198       else {
199         cerr << "Invalid number of columns in table" << endl;
200       }
201     }
202
203   } else {
204     brkpt_string = el->GetAttributeValue("breakPoint");
205     if (brkpt_string.empty()) {
206      // no independentVars found, and table is not marked as internal, nor is it a 3D table
207       cerr << endl << fgred << "No independent variable found for table."  << fgdef << endl << endl;
208       abort();
209     }
210   }
211   // end lookup property code
212
213   if (brkpt_string.empty()) {                  // Not a 3D table "table element"
214     tableData = el->FindElement("tableData");
215   } else {                                     // This is a table in a 3D table
216     tableData = el;
217     dimension = 2;                             // Currently, infers 2D table
218   }
219
220   for (i=0; i<tableData->GetNumDataLines(); i++) {
221     buf << tableData->GetDataLine(i) << string(" ");
222   }
223   switch (dimension) {
224   case 1:
225     nRows = tableData->GetNumDataLines();
226     nCols = 1;
227     Type = tt1D;
228     colCounter = 0;
229     rowCounter = 1;
230     Data = Allocate();
231     Debug(0);
232     lastRowIndex = lastColumnIndex = 2;
233     *this << buf;
234     break;
235   case 2:
236     nRows = tableData->GetNumDataLines()-1;
237
238     if (nRows >= 2) nCols = FindNumColumns(tableData->GetDataLine(0));
239     else {
240       cerr << endl << fgred << "Not enough rows in this table." << fgdef << endl;
241       abort();
242     }
243
244     Type = tt2D;
245     colCounter = 1;
246     rowCounter = 0;
247
248     Data = Allocate();
249     lastRowIndex = lastColumnIndex = 2;
250     *this << buf;
251     break;
252   case 3:
253     nTables = el->GetNumElements("tableData");
254     nRows = nTables;
255     nCols = 1;
256     Type = tt3D;
257     colCounter = 1;
258     rowCounter = 1;
259     lastRowIndex = lastColumnIndex = 2;
260
261     Data = Allocate(); // this data array will contain the keys for the associated tables
262     Tables.reserve(nTables); // necessary?
263     tableData = el->FindElement("tableData");
264     for (i=0; i<nTables; i++) {
265       Tables.push_back(new FGTable(PropertyManager, tableData));
266       Data[i+1][1] = tableData->GetAttributeValueAsNumber("breakPoint");
267       Tables[i]->SetRowIndexProperty(lookupProperty[eRow]);
268       Tables[i]->SetColumnIndexProperty(lookupProperty[eColumn]);
269       tableData = el->FindNextElement("tableData");
270     }
271
272     Debug(0);
273     break;
274   default:
275     cout << "No dimension given" << endl;
276     break;
277   }
278
279   bind();
280
281   if (debug_lvl & 1) Print();
282 }
283
284 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
285
286 double** FGTable::Allocate(void)
287 {
288   Data = new double*[nRows+1];
289   for (unsigned int r=0; r<=nRows; r++) {
290     Data[r] = new double[nCols+1];
291     for (unsigned int c=0; c<=nCols; c++) {
292       Data[r][c] = 0.0;
293     }
294   }
295   return Data;
296 }
297
298 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
299
300 FGTable::~FGTable()
301 {
302   if (nTables > 0) {
303     for (unsigned int i=0; i<nTables; i++) delete Tables[i];
304     Tables.clear();
305   }
306   for (unsigned int r=0; r<=nRows; r++) delete[] Data[r];
307   delete[] Data;
308
309   Debug(1);
310 }
311
312 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
313
314 unsigned int FGTable::FindNumColumns(const string& test_line)
315 {
316   // determine number of data columns in table (first column is row lookup - don't count)
317   size_t position=0;
318   unsigned int nCols=0;
319   while ((position = test_line.find_first_not_of(" \t", position)) != string::npos) {
320     nCols++;
321     position = test_line.find_first_of(" \t", position);
322   }
323   return nCols;
324 }
325
326 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
327
328 double FGTable::GetValue(void) const
329 {
330   double temp = 0;
331   double temp2 = 0;
332
333   switch (Type) {
334   case tt1D:
335     temp = lookupProperty[eRow]->getDoubleValue();
336     temp2 = GetValue(temp);
337     return temp2;
338   case tt2D:
339     return GetValue(lookupProperty[eRow]->getDoubleValue(),
340                     lookupProperty[eColumn]->getDoubleValue());
341   case tt3D:
342     return GetValue(lookupProperty[eRow]->getDoubleValue(),
343                     lookupProperty[eColumn]->getDoubleValue(),
344                     lookupProperty[eTable]->getDoubleValue());
345   default:
346     cerr << "Attempted to GetValue() for invalid/unknown table type" << endl;
347     throw(string("Attempted to GetValue() for invalid/unknown table type"));
348   }
349 }
350
351 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
352
353 double FGTable::GetValue(double key) const
354 {
355   double Factor, Value, Span;
356   unsigned int r = lastRowIndex;
357
358   //if the key is off the end of the table, just return the
359   //end-of-table value, do not extrapolate
360   if( key <= Data[1][0] ) {
361     lastRowIndex=2;
362     //cout << "Key underneath table: " << key << endl;
363     return Data[1][1];
364   } else if ( key >= Data[nRows][0] ) {
365     lastRowIndex=nRows;
366     //cout << "Key over table: " << key << endl;
367     return Data[nRows][1];
368   }
369
370   // the key is somewhere in the middle, search for the right breakpoint
371   // The search is particularly efficient if 
372   // the correct breakpoint has not changed since last frame or
373   // has only changed very little
374
375   while (r > 2     && Data[r-1][0] > key) { r--; }
376   while (r < nRows && Data[r][0]   < key) { r++; }
377
378   lastRowIndex=r;
379   // make sure denominator below does not go to zero.
380
381   Span = Data[r][0] - Data[r-1][0];
382   if (Span != 0.0) {
383     Factor = (key - Data[r-1][0]) / Span;
384     if (Factor > 1.0) Factor = 1.0;
385   } else {
386     Factor = 1.0;
387   }
388
389   Value = Factor*(Data[r][1] - Data[r-1][1]) + Data[r-1][1];
390
391   return Value;
392 }
393
394 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
395
396 double FGTable::GetValue(double rowKey, double colKey) const
397 {
398   double rFactor, cFactor, col1temp, col2temp, Value;
399   unsigned int r = lastRowIndex;
400   unsigned int c = lastColumnIndex;
401
402   while(r > 2     && Data[r-1][0] > rowKey) { r--; }
403   while(r < nRows && Data[r]  [0] < rowKey) { r++; }
404
405   while(c > 2     && Data[0][c-1] > colKey) { c--; }
406   while(c < nCols && Data[0][c]   < colKey) { c++; }
407
408   lastRowIndex=r;
409   lastColumnIndex=c;
410
411   rFactor = (rowKey - Data[r-1][0]) / (Data[r][0] - Data[r-1][0]);
412   cFactor = (colKey - Data[0][c-1]) / (Data[0][c] - Data[0][c-1]);
413
414   if (rFactor > 1.0) rFactor = 1.0;
415   else if (rFactor < 0.0) rFactor = 0.0;
416
417   if (cFactor > 1.0) cFactor = 1.0;
418   else if (cFactor < 0.0) cFactor = 0.0;
419
420   col1temp = rFactor*(Data[r][c-1] - Data[r-1][c-1]) + Data[r-1][c-1];
421   col2temp = rFactor*(Data[r][c] - Data[r-1][c]) + Data[r-1][c];
422
423   Value = col1temp + cFactor*(col2temp - col1temp);
424
425   return Value;
426 }
427
428 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
429
430 double FGTable::GetValue(double rowKey, double colKey, double tableKey) const
431 {
432   double Factor, Value, Span;
433   unsigned int r = lastRowIndex;
434
435   //if the key is off the end  (or before the beginning) of the table,
436   // just return the boundary-table value, do not extrapolate
437
438   if( tableKey <= Data[1][1] ) {
439     lastRowIndex=2;
440     return Tables[0]->GetValue(rowKey, colKey);
441   } else if ( tableKey >= Data[nRows][1] ) {
442     lastRowIndex=nRows;
443     return Tables[nRows-1]->GetValue(rowKey, colKey);
444   }
445
446   // the key is somewhere in the middle, search for the right breakpoint
447   // The search is particularly efficient if 
448   // the correct breakpoint has not changed since last frame or
449   // has only changed very little
450
451   while(r > 2     && Data[r-1][1] > tableKey) { r--; }
452   while(r < nRows && Data[r]  [1] < tableKey) { r++; }
453
454   lastRowIndex=r;
455   // make sure denominator below does not go to zero.
456
457   Span = Data[r][1] - Data[r-1][1];
458   if (Span != 0.0) {
459     Factor = (tableKey - Data[r-1][1]) / Span;
460     if (Factor > 1.0) Factor = 1.0;
461   } else {
462     Factor = 1.0;
463   }
464
465   Value = Factor*(Tables[r-1]->GetValue(rowKey, colKey) - Tables[r-2]->GetValue(rowKey, colKey))
466                               + Tables[r-2]->GetValue(rowKey, colKey);
467
468   return Value;
469 }
470
471 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
472
473 void FGTable::operator<<(istream& in_stream)
474 {
475   int startRow=0;
476   int startCol=0;
477
478 // In 1D table, no pseudo-row of column-headers (i.e. keys):
479   if (Type == tt1D) startRow = 1;
480
481   for (unsigned int r=startRow; r<=nRows; r++) {
482     for (unsigned int c=startCol; c<=nCols; c++) {
483       if (r != 0 || c != 0) {
484         in_stream >> Data[r][c];
485       }
486     }
487   }
488 }
489
490 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
491
492 FGTable& FGTable::operator<<(const double n)
493 {
494   Data[rowCounter][colCounter] = n;
495   if (colCounter == (int)nCols) {
496     colCounter = 0;
497     rowCounter++;
498   } else {
499     colCounter++;
500   }
501   return *this;
502 }
503
504 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
505
506 FGTable& FGTable::operator<<(const int n)
507 {
508   *this << (double)n;
509   return *this;
510 }
511
512 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
513
514 void FGTable::Print(void)
515 {
516   int startRow=0;
517   int startCol=0;
518
519   if (Type == tt1D || Type == tt3D) startRow = 1;
520   if (Type == tt3D) startCol = 1;
521
522 #if defined (sgi) && !defined(__GNUC__) && (_COMPILER_VERSION < 740)
523   unsigned long flags = cout.setf(ios::fixed);
524 #else
525   ios::fmtflags flags = cout.setf(ios::fixed); // set up output stream
526 #endif
527
528   switch(Type) {
529     case tt1D:
530       cout << "    1 dimensional table with " << nRows << " rows." << endl;
531       break;
532     case tt2D:
533       cout << "    2 dimensional table with " << nRows << " rows, " << nCols << " columns." << endl;
534       break;
535     case tt3D:
536       cout << "    3 dimensional table with " << nRows << " rows, "
537                                           << nCols << " columns "
538                                           << nTables << " tables." << endl;
539       break;
540   }
541   cout.precision(4);
542   for (unsigned int r=startRow; r<=nRows; r++) {
543     cout << "   ";
544     for (unsigned int c=startCol; c<=nCols; c++) {
545       if (r == 0 && c == 0) {
546         cout << "       ";
547       } else {
548         cout << Data[r][c] << " ";
549         if (Type == tt3D) {
550           cout << endl;
551           Tables[r-1]->Print();
552         }
553       }
554     }
555     cout << endl;
556   }
557   cout.setf(flags); // reset
558 }
559
560 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
561
562 void FGTable::bind(void)
563 {
564   typedef double (FGTable::*PMF)(void) const;
565   if ( !Name.empty() && !internal) {
566     string tmp = PropertyManager->mkPropertyName(Name, false); // Allow upper
567     PropertyManager->Tie( tmp, this, (PMF)&FGTable::GetValue);
568   }
569 }
570 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
571 //    The bitmasked value choices are as follows:
572 //    unset: In this case (the default) JSBSim would only print
573 //       out the normally expected messages, essentially echoing
574 //       the config files as they are read. If the environment
575 //       variable is not set, debug_lvl is set to 1 internally
576 //    0: This requests JSBSim not to output any messages
577 //       whatsoever.
578 //    1: This value explicity requests the normal JSBSim
579 //       startup messages
580 //    2: This value asks for a message to be printed out when
581 //       a class is instantiated
582 //    4: When this value is set, a message is displayed when a
583 //       FGModel object executes its Run() method
584 //    8: When this value is set, various runtime state variables
585 //       are printed out periodically
586 //    16: When set various parameters are sanity checked and
587 //       a message is printed out when they go out of bounds
588
589 void FGTable::Debug(int from)
590 {
591   if (debug_lvl <= 0) return;
592
593   if (debug_lvl & 1) { // Standard console startup message output
594     if (from == 0) { // Constructor
595
596     }
597   }
598   if (debug_lvl & 2 ) { // Instantiation/Destruction notification
599     if (from == 0) cout << "Instantiated: FGTable" << endl;
600     if (from == 1) cout << "Destroyed:    FGTable" << endl;
601   }
602   if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
603   }
604   if (debug_lvl & 8 ) { // Runtime state variables
605   }
606   if (debug_lvl & 16) { // Sanity checking
607   }
608   if (debug_lvl & 64) {
609     if (from == 0) { // Constructor
610       cout << IdSrc << endl;
611       cout << IdHdr << endl;
612     }
613   }
614 }
615 }