src/FDM/JSBSim/FGAircraft.cpp

   1 /*******************************************************************************
   2
   3  Module:       FGAircraft.cpp
   4  Author:       Jon S. Berndt
   5  Date started: 12/12/98
   6  Purpose:      Encapsulates an aircraft
   7  Called by:    FGFDMExec
   8
   9  ------------- Copyright (C) 1999  Jon S. Berndt (jsb@hal-pc.org) -------------
  10
  11  This program is free software; you can redistribute it and/or modify it under
  12  the terms of the GNU General Public License as published by the Free Software
  13  Foundation; either version 2 of the License, or (at your option) any later
  14  version.
  15
  16  This program is distributed in the hope that it will be useful, but WITHOUT
  17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
  18  FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
  19  details.
  20
  21  You should have received a copy of the GNU General Public License along with
  22  this program; if not, write to the Free Software Foundation, Inc., 59 Temple
  23  Place - Suite 330, Boston, MA  02111-1307, USA.
  24
  25  Further information about the GNU General Public License can also be found on
  26  the world wide web at http://www.gnu.org.
  27
  28 FUNCTIONAL DESCRIPTION
  29 --------------------------------------------------------------------------------
  30 Models the aircraft reactions and forces. This class is instantiated by the
  31 FGFDMExec class and scheduled as an FDM entry. LoadAircraft() is supplied with a
  32 name of a valid, registered aircraft, and the data file is parsed.
  33
  34 HISTORY
  35 --------------------------------------------------------------------------------
  36 12/12/98   JSB   Created
  37 04/03/99   JSB   Changed Aero() method to correct body axis force calculation
  38                  from wind vector. Fix provided by Tony Peden.
  39 05/03/99   JSB   Changed (for the better?) the way configurations are read in.
  40 9/17/99     TP   Combined force and moment functions. Added aero reference
  41                  point to config file. Added calculations for moments due to
  42                                  difference in cg and aero reference point
  43
  44 ********************************************************************************
  45 COMMENTS, REFERENCES,  and NOTES
  46 ********************************************************************************
  47 [1] Cooke, Zyda, Pratt, and McGhee, "NPSNET: Flight Simulation Dynamic Modeling
  48       Using Quaternions", Presence, Vol. 1, No. 4, pp. 404-420  Naval Postgraduate
  49       School, January 1994
  50 [2] D. M. Henderson, "Euler Angles, Quaternions, and Transformation Matrices",
  51       JSC 12960, July 1977
  52 [3] Richard E. McFarland, "A Standard Kinematic Model for Flight Simulation at
  53       NASA-Ames", NASA CR-2497, January 1975
  54 [4] Barnes W. McCormick, "Aerodynamics, Aeronautics, and Flight Mechanics",
  55       Wiley & Sons, 1979 ISBN 0-471-03032-5
  56 [5] Bernard Etkin, "Dynamics of Flight, Stability and Control", Wiley & Sons,
  57       1982 ISBN 0-471-08936-2
  58
  59 The aerodynamic coefficients used in this model are:
  60
  61 Longitudinal
  62   CL0 - Reference lift at zero alpha
  63   CD0 - Reference drag at zero alpha
  64   CDM - Drag due to Mach
  65   CLa - Lift curve slope (w.r.t. alpha)
  66   CDa - Drag curve slope (w.r.t. alpha)
  67   CLq - Lift due to pitch rate
  68   CLM - Lift due to Mach
  69   CLadt - Lift due to alpha rate
  70
  71   Cmadt - Pitching Moment due to alpha rate
  72   Cm0 - Reference Pitching moment at zero alpha
  73   Cma - Pitching moment slope (w.r.t. alpha)
  74   Cmq - Pitch damping (pitch moment due to pitch rate)
  75   CmM - Pitch Moment due to Mach
  76
  77 Lateral
  78   Cyb - Side force due to sideslip
  79   Cyr - Side force due to yaw rate
  80
  81   Clb - Dihedral effect (roll moment due to sideslip)
  82   Clp - Roll damping (roll moment due to roll rate)
  83   Clr - Roll moment due to yaw rate
  84   Cnb - Weathercocking stability (yaw moment due to sideslip)
  85   Cnp - Rudder adverse yaw (yaw moment due to roll rate)
  86   Cnr - Yaw damping (yaw moment due to yaw rate)
  87
  88 Control
  89   CLDe - Lift due to elevator
  90   CDDe - Drag due to elevator
  91   CyDr - Side force due to rudder
  92   CyDa - Side force due to aileron
  93
  94   CmDe - Pitch moment due to elevator
  95   ClDa - Roll moment due to aileron
  96   ClDr - Roll moment due to rudder
  97   CnDr - Yaw moment due to rudder
  98   CnDa - Yaw moment due to aileron
  99
 100 ********************************************************************************
 101 INCLUDES
 102 *******************************************************************************/
 103
 104 #include <sys/stat.h>
 105 #include <sys/types.h>
 106
 107 #ifdef FGFS
 108 # ifndef __BORLANDC__
 109 #  include <Include/compiler.h>
 110 # endif
 111 #  ifdef FG_HAVE_STD_INCLUDES
 112 #    include <cmath>
 113 #  else
 114 #    include <math.h>
 115 #  endif
 116 #else
 117 #  include <cmath>
 118 #endif
 119
 120 #include "FGAircraft.h"
 121 #include "FGTranslation.h"
 122 #include "FGRotation.h"
 123 #include "FGAtmosphere.h"
 124 #include "FGState.h"
 125 #include "FGFDMExec.h"
 126 #include "FGFCS.h"
 127 #include "FGPosition.h"
 128 #include "FGAuxiliary.h"
 129 #include "FGOutput.h"
 130
 131 /*******************************************************************************
 132 ************************************ CODE **************************************
 133 *******************************************************************************/
 134
 135 FGAircraft::FGAircraft(FGFDMExec* fdmex) : FGModel(fdmex)
 136 {
 137   int i;
 138
 139   Name = "FGAircraft";
 140
 141   for (i=0;i<6;i++) coeff_ctr[i] = 0;
 142 }
 143
 144
 145 FGAircraft::~FGAircraft(void)
 146 {
 147 }
 148
 149 bool FGAircraft::LoadAircraft(string aircraft_path, string engine_path, string fname)
 150 {
 151   string path;
 152   string fullpath;
 153   string filename;
 154   string aircraftDef;
 155   string tag;
 156   string holding_string;
 157   char scratch[128];
 158   ifstream coeffInFile;
 159   streampos gpos;
 160   int axis;
 161   string axis_descript;
 162
 163   axis = -1;
 164   aircraftDef = aircraft_path + "/" + fname + "/" + fname + ".cfg";
 165   ifstream aircraftfile(aircraftDef.c_str());
 166   cout << "Reading Aircraft Configuration File: " << aircraftDef << endl;
 167   Output->SocketStatusOutput("Reading Aircraft Configuration File: " + aircraftDef);
 168
 169   numTanks = numEngines = 0;
 170   numSelectedOxiTanks = numSelectedFuelTanks = 0;
 171
 172   while (!aircraftfile.fail()) {
 173     holding_string.erase();
 174     aircraftfile >> holding_string;
 175 #if defined(__BORLANDC__) || defined(FG_HAVE_NATIVE_SGI_COMPILERS) || defined(_MSC_VER)
 176     if (holding_string.compare(0, 2, "//") != 0) {
 177 #else
 178     if (holding_string.compare("//",0,2) != 0) {
 179 #endif
 180       if (holding_string == "CFG_VERSION") {
 181         aircraftfile >> CFGVersion;
 182         cout << "Config file version: " << CFGVersion << endl;
 183         if (CFGVersion < NEEDED_CFG_VERSION) {
 184           cout << endl << "YOU HAVE AN OLD CFG FILE FOR THIS AIRCRAFT."
 185                           " RESULTS WILL BE UNPREDICTABLE !!" << endl << endl;
 186         }
 187       } else if (holding_string == "AIRCRAFT") {
 188         cout << "Reading in Aircraft parameters ..." << endl;
 189       } else if (holding_string == "AERODYNAMICS") {
 190         cout << "Reading in Aerodynamic parameters ..." << endl;
 191       } else if (holding_string == "AC_NAME") {
 192         aircraftfile >> AircraftName;   // String with no embedded spaces
 193         cout << "Aircraft Name: " << AircraftName << endl;
 194       } else if (holding_string == "AC_WINGAREA") {
 195         aircraftfile >> WingArea;
 196         cout << "Aircraft Wing Area: " << WingArea << endl;
 197       } else if (holding_string == "AC_WINGSPAN") {
 198         aircraftfile >> WingSpan;
 199         cout << "Aircraft WingSpan: " << WingSpan << endl;
 200       } else if (holding_string == "AC_CHORD") {
 201         aircraftfile >> cbar;
 202         cout << "Aircraft Chord: " << cbar << endl;
 203       } else if (holding_string == "AC_IXX") {
 204         aircraftfile >> baseIxx;
 205         cout << "Aircraft Base Ixx: " << baseIxx << endl;
 206       } else if (holding_string == "AC_IYY") {
 207         aircraftfile >> baseIyy;
 208         cout << "Aircraft Base Iyy: " << baseIyy << endl;
 209       } else if (holding_string == "AC_IZZ") {
 210         aircraftfile >> baseIzz;
 211         cout << "Aircraft Base Izz: " << baseIzz << endl;
 212       } else if (holding_string == "AC_IXZ") {
 213         aircraftfile >> baseIxz;
 214         cout << "Aircraft Base Ixz: " << baseIxz << endl;
 215       } else if (holding_string == "AC_EMPTYWT") {
 216         aircraftfile >> EmptyWeight;
 217         EmptyMass = EmptyWeight / GRAVITY;
 218         cout << "Aircraft Empty Weight: " << EmptyWeight << endl;
 219       } else if (holding_string == "AC_AERORP") {
 220         aircraftfile >> Xrp >> Yrp >> Zrp;
 221         cout << "Aerodynamic Reference Point: " << Xrp << " " << Yrp << " " << Zrp << endl;
 222       } else if (holding_string == "AC_CGLOC") {
 223         aircraftfile >> baseXcg >> baseYcg >> baseZcg;
 224         cout << "Aircraft Base C.G.: " << baseXcg << " " << baseYcg << " " << baseZcg << endl;
 225       } else if (holding_string == "AC_EYEPTLOC") {
 226         aircraftfile >> Xep >> Yep >> Zep;
 227         cout << "Pilot Eyepoint: " << Xep << " " << Yep << " " << Zep << endl;
 228       } else if (holding_string == "AC_TANK") {
 229         Tank[numTanks] = new FGTank(aircraftfile);
 230         switch(Tank[numTanks]->GetType()) {
 231         case FGTank::ttFUEL:
 232           numSelectedFuelTanks++;
 233           cout << "Reading in Fuel Tank #" << numSelectedFuelTanks << " parameters ..." << endl;
 234           break;
 235         case FGTank::ttOXIDIZER:
 236           numSelectedOxiTanks++;
 237           cout << "Reading in Oxidizer Tank #" << numSelectedOxiTanks << " parameters ..." << endl;
 238           break;
 239         }
 240         numTanks++;
 241
 242       } else if (holding_string == "AC_GEAR") {
 243
 244         lGear.push_back(new FGLGear(aircraftfile));
 245
 246       } else if (holding_string == "AC_ENGINE") {
 247
 248         aircraftfile >> tag;
 249         cout << "Reading in " << tag << " Engine parameters ..." << endl;
 250         Engine[numEngines] = new FGEngine(FDMExec, engine_path, tag, numEngines);
 251         numEngines++;
 252
 253       } else if (holding_string == "}") {
 254
 255       } else if (holding_string == "{") {
 256
 257       } else if (holding_string == "LIFT") {
 258
 259         axis_descript = "   Lift Coefficients ...";
 260         axis = LiftCoeff;
 261
 262       } else if (holding_string == "DRAG") {
 263
 264         axis_descript = "   Drag Coefficients ...";
 265         axis = DragCoeff;
 266
 267       } else if (holding_string == "SIDE") {
 268
 269         axis_descript = "   Side Coefficients ...";
 270         axis = SideCoeff;
 271
 272       } else if (holding_string == "ROLL") {
 273
 274         axis_descript = "   Roll Coefficients ...";
 275         axis = RollCoeff;
 276
 277       } else if (holding_string == "PITCH") {
 278
 279         axis_descript = "   Pitch Coefficients ...";
 280         axis = PitchCoeff;
 281
 282       } else if (holding_string == "YAW") {
 283
 284         axis_descript = "   Yaw Coefficients ...";
 285         axis = YawCoeff;
 286
 287       }
 288
 289       if (axis >= 0) {
 290         cout << axis_descript << endl;
 291         aircraftfile >> tag;
 292         gpos = aircraftfile.tellg();
 293         aircraftfile >> tag;
 294         if ( !(tag == "}") ) {
 295           while ( !(tag == "}") ) {
 296             aircraftfile.seekg(gpos);
 297             Coeff[axis][coeff_ctr[axis]] = new FGCoefficient(FDMExec, aircraftfile);
 298             coeff_ctr[axis]++;
 299             aircraftfile >> tag;
 300             gpos = aircraftfile.tellg();
 301             aircraftfile >> tag;
 302           }
 303         } else {
 304           cout << "      None found ..." << endl;
 305         }
 306       }
 307       axis = -1;
 308
 309     } else {
 310       aircraftfile.getline(scratch, 127);
 311     }
 312   }
 313   cout << "End of Configuration File Parsing." << endl;
 314
 315   return true;
 316 }
 317
 318
 319 bool FGAircraft::Run(void)
 320 {
 321   if (!FGModel::Run()) {                 // if false then execute this Run()
 322     GetState();
 323
 324     for (int i = 0; i < 3; i++)  Forces[i] = Moments[i] = 0.0;
 325
 326     MassChange();
 327
 328     FMProp(); FMAero(); FMGear(); FMMass();
 329
 330     PutState();
 331   } else {                               // skip Run() execution this time
 332   }
 333   return false;
 334 }
 335
 336
 337 void FGAircraft::MassChange()
 338 {
 339   float Xt, Xw, Yt, Yw, Zt, Zw, Tw;
 340   float IXXt, IYYt, IZZt, IXZt;
 341   int t;
 342
 343   // UPDATE TANK CONTENTS
 344   //
 345   // For each engine, cycle through the tanks and draw an equal amount of
 346   // fuel (or oxidizer) from each active tank. The needed amount of fuel is
 347   // determined by the engine in the FGEngine class. If more fuel is needed
 348   // than is available in the tank, then that amount is considered a shortage,
 349   // and will be drawn from the next tank. If the engine cannot be fed what it
 350   // needs, it will be considered to be starved, and will shut down.
 351
 352   float Oshortage, Fshortage;
 353
 354   for (int e=0; e<numEngines; e++) {
 355     Fshortage = Oshortage = 0.0;
 356     for (t=0; t<numTanks; t++) {
 357       switch(Engine[e]->GetType()) {
 358       case FGEngine::etRocket:
 359
 360         switch(Tank[t]->GetType()) {
 361         case FGTank::ttFUEL:
 362           if (Tank[t]->GetSelected()) {
 363             Fshortage = Tank[t]->Reduce((Engine[e]->CalcFuelNeed()/
 364                                          numSelectedFuelTanks)*(dt*rate) + Fshortage);
 365           }
 366           break;
 367         case FGTank::ttOXIDIZER:
 368           if (Tank[t]->GetSelected()) {
 369             Oshortage = Tank[t]->Reduce((Engine[e]->CalcOxidizerNeed()/
 370                                          numSelectedOxiTanks)*(dt*rate) + Oshortage);
 371           }
 372           break;
 373         }
 374         break;
 375
 376       case FGEngine::etPiston:
 377       case FGEngine::etTurboJet:
 378       case FGEngine::etTurboProp:
 379
 380         if (Tank[t]->GetSelected()) {
 381           Fshortage = Tank[t]->Reduce((Engine[e]->CalcFuelNeed()/
 382                                        numSelectedFuelTanks)*(dt*rate) + Fshortage);
 383         }
 384         break;
 385       }
 386     }
 387     if ((Fshortage <= 0.0) || (Oshortage <= 0.0)) Engine[e]->SetStarved();
 388     else Engine[e]->SetStarved(false);
 389   }
 390
 391   Weight = EmptyWeight;
 392   for (t=0; t<numTanks; t++)
 393     Weight += Tank[t]->GetContents();
 394
 395   Mass = Weight / GRAVITY;
 396
 397   // Calculate new CG here.
 398
 399   Xt = Yt = Zt = Tw = 0;
 400   Xw = Yw = Zw = 0;
 401   for (t=0; t<numTanks; t++) {
 402     Xt += Tank[t]->GetX()*Tank[t]->GetContents();
 403     Yt += Tank[t]->GetY()*Tank[t]->GetContents();
 404     Zt += Tank[t]->GetZ()*Tank[t]->GetContents();
 405
 406     Tw += Tank[t]->GetContents();
 407   }
 408
 409   Xcg = (Xt + EmptyWeight*baseXcg) / (Tw + EmptyWeight);
 410   Ycg = (Yt + EmptyWeight*baseYcg) / (Tw + EmptyWeight);
 411   Zcg = (Zt + EmptyWeight*baseZcg) / (Tw + EmptyWeight);
 412
 413   // Calculate new moments of inertia here
 414
 415   IXXt = IYYt = IZZt = IXZt = 0.0;
 416   for (t=0; t<numTanks; t++) {
 417     IXXt += ((Tank[t]->GetX()-Xcg)/12.0)*((Tank[t]->GetX() - Xcg)/12.0)*Tank[t]->GetContents()/GRAVITY;
 418     IYYt += ((Tank[t]->GetY()-Ycg)/12.0)*((Tank[t]->GetY() - Ycg)/12.0)*Tank[t]->GetContents()/GRAVITY;
 419     IZZt += ((Tank[t]->GetZ()-Zcg)/12.0)*((Tank[t]->GetZ() - Zcg)/12.0)*Tank[t]->GetContents()/GRAVITY;
 420     IXZt += ((Tank[t]->GetX()-Xcg)/12.0)*((Tank[t]->GetZ() - Zcg)/12.0)*Tank[t]->GetContents()/GRAVITY;
 421   }
 422
 423   Ixx = baseIxx + IXXt;
 424   Iyy = baseIyy + IYYt;
 425   Izz = baseIzz + IZZt;
 426   Ixz = baseIxz + IXZt;
 427
 428 }
 429
 430
 431 void FGAircraft::FMAero(void)
 432 {
 433   float F[3];
 434   float Fxaero,Fyaero,Fzaero;
 435   float dxcg,dycg,dzcg;
 436   int axis_ctr,ctr;
 437   F[0] = F[1] = F[2] = 0.0;
 438
 439   for (axis_ctr = 0; axis_ctr < 3; axis_ctr++)
 440     for (ctr=0; ctr < coeff_ctr[axis_ctr]; ctr++)
 441       F[axis_ctr] += Coeff[axis_ctr][ctr]->TotalValue();
 442
 443   Fxaero = - F[DragCoeff]*cos(alpha)*cos(beta)
 444            - F[SideCoeff]*cos(alpha)*sin(beta)
 445            + F[LiftCoeff]*sin(alpha);
 446   Fyaero =   F[DragCoeff]*sin(beta)
 447              + F[SideCoeff]*cos(beta);
 448   Fzaero = - F[DragCoeff]*sin(alpha)*cos(beta)
 449            - F[SideCoeff]*sin(alpha)*sin(beta)
 450            - F[LiftCoeff]*cos(alpha);
 451
 452   Forces[0] += - F[DragCoeff]*cos(alpha)*cos(beta)
 453                - F[SideCoeff]*cos(alpha)*sin(beta)
 454                + F[LiftCoeff]*sin(alpha);
 455   Forces[1] +=   F[DragCoeff]*sin(beta)
 456                  + F[SideCoeff]*cos(beta);
 457   Forces[2] += - F[DragCoeff]*sin(alpha)*cos(beta)
 458                - F[SideCoeff]*sin(alpha)*sin(beta)
 459                - F[LiftCoeff]*cos(alpha);
 460
 461   // The d*cg distances below, given in inches, are the distances FROM the c.g.
 462   // TO the reference point. Since the c.g. and ref point are given in inches in
 463   // the structural system (X positive rearwards) and the body coordinate system
 464   // is given with X positive out the nose, the dxcg and dzcg values are
 465   // *rotated* 180 degrees about the Y axis.
 466
 467   dxcg = -(Xrp - Xcg)/12; //cg and rp values are in inches
 468   dycg =  (Yrp - Ycg)/12;
 469   dzcg = -(Zrp - Zcg)/12;
 470
 471   Moments[0] +=  Fzaero*dycg - Fyaero*dzcg; //rolling moment
 472   Moments[1] +=  Fxaero*dzcg - Fzaero*dxcg; //pitching moment
 473   Moments[2] += -Fxaero*dycg + Fyaero*dxcg; //yawing moment
 474
 475   for (axis_ctr = 0; axis_ctr < 3; axis_ctr++) {
 476     for (ctr = 0; ctr < coeff_ctr[axis_ctr+3]; ctr++) {
 477       Moments[axis_ctr] += Coeff[axis_ctr+3][ctr]->TotalValue();
 478     }
 479   }
 480 }
 481
 482
 483 void FGAircraft::FMGear(void)
 484 {
 485   if (GearUp) {
 486     // crash routine
 487   } else {
 488     for (int i=0;i<lGear.size();i++) {
 489       //      lGear[i].
 490     }
 491   }
 492 }
 493
 494
 495 void FGAircraft::FMMass(void)
 496 {
 497   Forces[0] += -GRAVITY*sin(tht) * Mass;
 498   Forces[1] +=  GRAVITY*sin(phi)*cos(tht) * Mass;
 499   Forces[2] +=  GRAVITY*cos(phi)*cos(tht) * Mass;
 500 }
 501
 502
 503 void FGAircraft::FMProp(void)
 504 {
 505   for (int i=0;i<numEngines;i++) {
 506     Forces[0] += Engine[i]->CalcThrust();
 507   }
 508 }
 509
 510 void FGAircraft::GetState(void)
 511 {
 512   dt = State->Getdt();
 513
 514   alpha = Translation->Getalpha();
 515   beta = Translation->Getbeta();
 516   phi = Rotation->Getphi();
 517   tht = Rotation->Gettht();
 518   psi = Rotation->Getpsi();
 519 }
 520
 521
 522 void FGAircraft::PutState(void)
 523 {
 524 }
 525