]> git.mxchange.org Git - flightgear.git/blob - src/FDM/JSBSim/JSBSim.cxx
5a1be0a29b25d1d77b8570e90fb4bb05e099277e
[flightgear.git] / src / FDM / JSBSim / JSBSim.cxx
1 // JSBsim.cxx -- interface to the JSBsim flight model
2 //
3 // Written by Curtis Olson, started February 1999.
4 //
5 // Copyright (C) 1999  Curtis L. Olson  - curt@flightgear.org
6 //
7 // This program is free software; you can redistribute it and/or
8 // modify it under the terms of the GNU General Public License as
9 // published by the Free Software Foundation; either version 2 of the
10 // License, or (at your option) any later version.
11 //
12 // This program is distributed in the hope that it will be useful, but
13 // WITHOUT ANY WARRANTY; without even the implied warranty of
14 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
15 // 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., 675 Mass Ave, Cambridge, MA 02139, USA.
20 //
21 // $Id$
22
23
24 #ifdef HAVE_CONFIG_H
25 #  include <config.h>
26 #endif
27
28 #include <simgear/compiler.h>
29
30 #include <stdio.h>    //    size_t
31 #include <string>
32
33 #include <simgear/constants.h>
34 #include <simgear/debug/logstream.hxx>
35 #include <simgear/math/sg_geodesy.hxx>
36 #include <simgear/misc/sg_path.hxx>
37 #include <simgear/structure/commands.hxx>
38
39 #include <FDM/flight.hxx>
40
41 #include <Aircraft/aircraft.hxx>
42 #include <Aircraft/controls.hxx>
43 #include <Main/globals.hxx>
44 #include <Main/fg_props.hxx>
45
46 #include "JSBSim.hxx"
47 #include <FDM/JSBSim/FGFDMExec.h>
48 #include <FDM/JSBSim/FGJSBBase.h>
49 #include <FDM/JSBSim/FGState.h>
50 #include <FDM/JSBSim/initialization/FGInitialCondition.h>
51 #include <FDM/JSBSim/initialization/FGTrim.h>
52 #include <FDM/JSBSim/models/FGModel.h>
53 #include <FDM/JSBSim/models/FGAircraft.h>
54 #include <FDM/JSBSim/models/FGFCS.h>
55 #include <FDM/JSBSim/models/FGPropagate.h>
56 #include <FDM/JSBSim/models/FGAuxiliary.h>
57 #include <FDM/JSBSim/models/FGInertial.h>
58 #include <FDM/JSBSim/models/FGAtmosphere.h>
59 #include <FDM/JSBSim/models/FGMassBalance.h>
60 #include <FDM/JSBSim/models/FGAerodynamics.h>
61 #include <FDM/JSBSim/models/FGLGear.h>
62 #include <FDM/JSBSim/models/propulsion/FGEngine.h>
63 #include <FDM/JSBSim/models/propulsion/FGPiston.h>
64 #include <FDM/JSBSim/models/propulsion/FGTurbine.h>
65 #include <FDM/JSBSim/models/propulsion/FGTurboProp.h>
66 #include <FDM/JSBSim/models/propulsion/FGRocket.h>
67 #include <FDM/JSBSim/models/propulsion/FGElectric.h>
68 #include <FDM/JSBSim/models/propulsion/FGNozzle.h>
69 #include <FDM/JSBSim/models/propulsion/FGPropeller.h>
70 #include <FDM/JSBSim/models/propulsion/FGRotor.h>
71 #include <FDM/JSBSim/models/propulsion/FGTank.h>
72 #include <FDM/JSBSim/input_output/FGPropertyManager.h>
73 #include <FDM/JSBSim/input_output/FGGroundCallback.h>
74
75 using namespace JSBSim;
76
77 static inline double
78 FMAX (double a, double b)
79 {
80   return a > b ? a : b;
81 }
82
83 class FGFSGroundCallback : public FGGroundCallback {
84 public:
85   FGFSGroundCallback(FGJSBsim* ifc) : mInterface(ifc) {}
86   virtual ~FGFSGroundCallback() {}
87
88   /** Get the altitude above sea level depenent on the location. */
89   virtual double GetAltitude(const FGLocation& l) const {
90     double pt[3] = { SG_FEET_TO_METER*l(eX),
91                      SG_FEET_TO_METER*l(eY),
92                      SG_FEET_TO_METER*l(eZ) };
93     double lat, lon, alt;
94     sgCartToGeod( pt, &lat, &lon, &alt);
95     return alt * SG_METER_TO_FEET;
96   }
97
98   /** Compute the altitude above ground. */
99   virtual double GetAGLevel(double t, const FGLocation& l,
100                             FGLocation& cont,
101                             FGColumnVector3& n, FGColumnVector3& v) const {
102     double loc_cart[3] = { l(eX), l(eY), l(eZ) };
103     double contact[3], normal[3], vel[3], agl = 0;
104     mInterface->get_agl_ft(t, loc_cart, SG_METER_TO_FEET*2, contact, normal,
105                            vel, &agl);
106     n = l.GetTec2l()*FGColumnVector3( normal[0], normal[1], normal[2] );
107     v = l.GetTec2l()*FGColumnVector3( vel[0], vel[1], vel[2] );
108     cont = FGColumnVector3( contact[0], contact[1], contact[2] );
109     return agl;
110   }
111 private:
112   FGJSBsim* mInterface;
113 };
114
115 /******************************************************************************/
116
117 FGJSBsim::FGJSBsim( double dt )
118   : FGInterface(dt), got_wire(false)
119 {
120     bool result;
121                                 // Set up the debugging level
122                                 // FIXME: this will not respond to
123                                 // runtime changes
124
125                                 // if flight is excluded, don't bother
126     if ((logbuf::get_log_classes() & SG_FLIGHT) != 0) {
127
128                                 // do a rough-and-ready mapping to
129                                 // the levels documented in FGFDMExec.h
130         switch (logbuf::get_log_priority()) {
131         case SG_BULK:
132             FGJSBBase::debug_lvl = 0x1f;
133             break;
134         case SG_DEBUG:
135             FGJSBBase::debug_lvl = 0x0f;
136         case SG_INFO:
137             FGJSBBase::debug_lvl = 0x01;
138             break;
139         case SG_WARN:
140         case SG_ALERT:
141             FGJSBBase::debug_lvl = 0x00;
142             break;
143         }
144     }
145
146     fdmex = new FGFDMExec( (FGPropertyManager*)globals->get_props() );
147
148     // begin ugly hack
149     // Untie the write-state-file property to avoid creating an initfile.xml
150     // file on each FlightGear reset.
151     fgGetNode("/fdm/jsbsim/simulation/write-state-file")->untie();
152     fgGetNode("/fdm/jsbsim/simulation")->removeChild("write-state-file", false);
153     // end ugly hack
154
155     // Register ground callback.
156     fdmex->SetGroundCallback( new FGFSGroundCallback(this) );
157
158     State           = fdmex->GetState();
159     Atmosphere      = fdmex->GetAtmosphere();
160     FCS             = fdmex->GetFCS();
161     MassBalance     = fdmex->GetMassBalance();
162     Propulsion      = fdmex->GetPropulsion();
163     Aircraft        = fdmex->GetAircraft();
164     Propagate        = fdmex->GetPropagate();
165     Auxiliary       = fdmex->GetAuxiliary();
166     Inertial        = fdmex->GetInertial();
167     Aerodynamics    = fdmex->GetAerodynamics();
168     GroundReactions = fdmex->GetGroundReactions();
169
170     fgic=fdmex->GetIC();
171     needTrim=true;
172
173     SGPath aircraft_path( fgGetString("/sim/aircraft-dir") );
174
175     SGPath engine_path( fgGetString("/sim/aircraft-dir") );
176     engine_path.append( "Engine" );
177
178     SGPath systems_path( fgGetString("/sim/aircraft-dir") );
179     systems_path.append( "Systems" );
180
181     State->Setdt( dt );
182
183     result = fdmex->LoadModel( aircraft_path.str(),
184                                engine_path.str(),
185                                systems_path.str(),
186                                fgGetString("/sim/aero"), false );
187
188     if (result) {
189       SG_LOG( SG_FLIGHT, SG_INFO, "  loaded aero.");
190     } else {
191       SG_LOG( SG_FLIGHT, SG_INFO,
192               "  aero does not exist (you may have mis-typed the name).");
193       throw(-1);
194     }
195
196     SG_LOG( SG_FLIGHT, SG_INFO, "" );
197     SG_LOG( SG_FLIGHT, SG_INFO, "" );
198     SG_LOG( SG_FLIGHT, SG_INFO, "After loading aero definition file ..." );
199
200     int Neng = Propulsion->GetNumEngines();
201     SG_LOG( SG_FLIGHT, SG_INFO, "num engines = " << Neng );
202
203     if ( GroundReactions->GetNumGearUnits() <= 0 ) {
204         SG_LOG( SG_FLIGHT, SG_ALERT, "num gear units = "
205                 << GroundReactions->GetNumGearUnits() );
206         SG_LOG( SG_FLIGHT, SG_ALERT, "This is a very bad thing because with 0 gear units, the ground trimming");
207         SG_LOG( SG_FLIGHT, SG_ALERT, "routine (coming up later in the code) will core dump.");
208         SG_LOG( SG_FLIGHT, SG_ALERT, "Halting the sim now, and hoping a solution will present itself soon!");
209         exit(-1);
210     }
211
212     init_gear();
213
214     // Set initial fuel levels if provided.
215     for (unsigned int i = 0; i < Propulsion->GetNumTanks(); i++) {
216       SGPropertyNode * node = fgGetNode("/consumables/fuel/tank", i, true);
217       if (node->getChild("level-gal_us", 0, false) != 0) {
218         Propulsion->GetTank(i)->SetContents(node->getDoubleValue("level-gal_us") * 6.6);
219       } else {
220         node->setDoubleValue("level-lbs", Propulsion->GetTank(i)->GetContents());
221         node->setDoubleValue("level-gal_us", Propulsion->GetTank(i)->GetContents() / 6.6);
222       }
223       node->setDoubleValue("capacity-gal_us",
224                            Propulsion->GetTank(i)->GetCapacity() / 6.6);
225     }
226     Propulsion->SetFuelFreeze((fgGetNode("/sim/freeze/fuel",true))->getBoolValue());
227
228     fgSetDouble("/fdm/trim/pitch-trim", FCS->GetPitchTrimCmd());
229     fgSetDouble("/fdm/trim/throttle",   FCS->GetThrottleCmd(0));
230     fgSetDouble("/fdm/trim/aileron",    FCS->GetDaCmd());
231     fgSetDouble("/fdm/trim/rudder",     FCS->GetDrCmd());
232
233     startup_trim = fgGetNode("/sim/presets/trim", true);
234
235     trimmed = fgGetNode("/fdm/trim/trimmed", true);
236     trimmed->setBoolValue(false);
237
238     pitch_trim = fgGetNode("/fdm/trim/pitch-trim", true );
239     throttle_trim = fgGetNode("/fdm/trim/throttle", true );
240     aileron_trim = fgGetNode("/fdm/trim/aileron", true );
241     rudder_trim = fgGetNode("/fdm/trim/rudder", true );
242
243     stall_warning = fgGetNode("/sim/alarms/stall-warning",true);
244     stall_warning->setDoubleValue(0);
245
246
247     flap_pos_pct=fgGetNode("/surface-positions/flap-pos-norm",true);
248     elevator_pos_pct=fgGetNode("/surface-positions/elevator-pos-norm",true);
249     left_aileron_pos_pct
250         =fgGetNode("/surface-positions/left-aileron-pos-norm",true);
251     right_aileron_pos_pct
252         =fgGetNode("/surface-positions/right-aileron-pos-norm",true);
253     rudder_pos_pct=fgGetNode("/surface-positions/rudder-pos-norm",true);
254     speedbrake_pos_pct
255         =fgGetNode("/surface-positions/speedbrake-pos-norm",true);
256     spoilers_pos_pct=fgGetNode("/surface-positions/spoilers-pos-norm",true);
257     tailhook_pos_pct=fgGetNode("/gear/tailhook/position-norm",true);
258     wing_fold_pos_pct=fgGetNode("surface-positions/wing-fold-pos-norm",true);
259
260     elevator_pos_pct->setDoubleValue(0);
261     left_aileron_pos_pct->setDoubleValue(0);
262     right_aileron_pos_pct->setDoubleValue(0);
263     rudder_pos_pct->setDoubleValue(0);
264     flap_pos_pct->setDoubleValue(0);
265     speedbrake_pos_pct->setDoubleValue(0);
266     spoilers_pos_pct->setDoubleValue(0);
267
268     temperature = fgGetNode("/environment/temperature-degc",true);
269     pressure = fgGetNode("/environment/pressure-inhg",true);
270     density = fgGetNode("/environment/density-slugft3",true);
271     turbulence_gain = fgGetNode("/environment/turbulence/magnitude-norm",true);
272     turbulence_rate = fgGetNode("/environment/turbulence/rate-hz",true);
273
274     wind_from_north= fgGetNode("/environment/wind-from-north-fps",true);
275     wind_from_east = fgGetNode("/environment/wind-from-east-fps" ,true);
276     wind_from_down = fgGetNode("/environment/wind-from-down-fps" ,true);
277
278     slaved = fgGetNode("/sim/slaved/enabled", true);
279
280     for (unsigned int i = 0; i < Propulsion->GetNumEngines(); i++) {
281       SGPropertyNode * node = fgGetNode("engines/engine", i, true);
282       Propulsion->GetEngine(i)->GetThruster()->SetRPM(node->getDoubleValue("rpm") /
283                      Propulsion->GetEngine(i)->GetThruster()->GetGearRatio());
284     }
285
286     hook_root_struct = FGColumnVector3(
287         fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-offset-x-in", 196),
288         fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-offset-y-in", 0),
289         fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-offset-z-in", -16));
290
291     crashed = false;
292 }
293
294 /******************************************************************************/
295 FGJSBsim::~FGJSBsim(void)
296 {
297   delete fdmex;
298 }
299
300 /******************************************************************************/
301
302 // Initialize the JSBsim flight model, dt is the time increment for
303 // each subsequent iteration through the EOM
304
305 void FGJSBsim::init()
306 {
307     double tmp;
308
309     SG_LOG( SG_FLIGHT, SG_INFO, "Starting and initializing JSBsim" );
310
311     // Explicitly call the superclass's
312     // init method first.
313
314     if (fgGetBool("/environment/params/control-fdm-atmosphere")) {
315       Atmosphere->UseExternal();
316       Atmosphere->SetExTemperature(
317                   9.0/5.0*(temperature->getDoubleValue()+273.15) );
318       Atmosphere->SetExPressure(pressure->getDoubleValue()*70.726566);
319       Atmosphere->SetExDensity(density->getDoubleValue());
320       Atmosphere->SetTurbType(FGAtmosphere::ttCulp);
321       Atmosphere->SetTurbGain(turbulence_gain->getDoubleValue());
322       Atmosphere->SetTurbRate(turbulence_rate->getDoubleValue());
323
324     } else {
325       Atmosphere->UseInternal();
326     }
327
328     fgic->SetVNorthFpsIC( -wind_from_north->getDoubleValue() );
329     fgic->SetVEastFpsIC( -wind_from_east->getDoubleValue() );
330     fgic->SetVDownFpsIC( -wind_from_down->getDoubleValue() );
331
332     //Atmosphere->SetExTemperature(get_Static_temperature());
333     //Atmosphere->SetExPressure(get_Static_pressure());
334     //Atmosphere->SetExDensity(get_Density());
335     SG_LOG(SG_FLIGHT,SG_INFO,"T,p,rho: " << fdmex->GetAtmosphere()->GetTemperature()
336      << ", " << fdmex->GetAtmosphere()->GetPressure()
337      << ", " << fdmex->GetAtmosphere()->GetDensity() );
338
339     if (fgGetBool("/sim/presets/running")) {
340           for (unsigned int i=0; i < Propulsion->GetNumEngines(); i++) {
341             SGPropertyNode * node = fgGetNode("engines/engine", i, true);
342             node->setBoolValue("running", true);
343             Propulsion->GetEngine(i)->SetRunning(true);
344           }
345     }
346
347     FCS->SetDfPos( ofNorm, globals->get_controls()->get_flaps() );
348
349     common_init();
350
351     copy_to_JSBsim();
352     fdmex->RunIC();     //loop JSBSim once w/o integrating
353     copy_from_JSBsim(); //update the bus
354
355     SG_LOG( SG_FLIGHT, SG_INFO, "  Initialized JSBSim with:" );
356
357     switch(fgic->GetSpeedSet()) {
358     case setned:
359         SG_LOG(SG_FLIGHT,SG_INFO, "  Vn,Ve,Vd= "
360                << Propagate->GetVel(FGJSBBase::eNorth) << ", "
361                << Propagate->GetVel(FGJSBBase::eEast) << ", "
362                << Propagate->GetVel(FGJSBBase::eDown) << " ft/s");
363     break;
364     case setuvw:
365         SG_LOG(SG_FLIGHT,SG_INFO, "  U,V,W= "
366                << Propagate->GetUVW(1) << ", "
367                << Propagate->GetUVW(2) << ", "
368                << Propagate->GetUVW(3) << " ft/s");
369     break;
370     case setmach:
371         SG_LOG(SG_FLIGHT,SG_INFO, "  Mach: "
372                << Auxiliary->GetMach() );
373     break;
374     case setvc:
375     default:
376         SG_LOG(SG_FLIGHT,SG_INFO, "  Indicated Airspeed: "
377                << Auxiliary->GetVcalibratedKTS() << " knots" );
378     break;
379     }
380
381     stall_warning->setDoubleValue(0);
382
383     SG_LOG( SG_FLIGHT, SG_INFO, "  Bank Angle: "
384             << Propagate->GetEuler(FGJSBBase::ePhi)*RADTODEG << " deg" );
385     SG_LOG( SG_FLIGHT, SG_INFO, "  Pitch Angle: "
386             << Propagate->GetEuler(FGJSBBase::eTht)*RADTODEG << " deg" );
387     SG_LOG( SG_FLIGHT, SG_INFO, "  True Heading: "
388             << Propagate->GetEuler(FGJSBBase::ePsi)*RADTODEG << " deg" );
389     SG_LOG( SG_FLIGHT, SG_INFO, "  Latitude: "
390             << Propagate->GetLocation().GetLatitudeDeg() << " deg" );
391     SG_LOG( SG_FLIGHT, SG_INFO, "  Longitude: "
392             << Propagate->GetLocation().GetLongitudeDeg() << " deg" );
393     SG_LOG( SG_FLIGHT, SG_INFO, "  Altitude: "
394             << Propagate->Geth() << " feet" );
395     SG_LOG( SG_FLIGHT, SG_INFO, "  loaded initial conditions" );
396
397     SG_LOG( SG_FLIGHT, SG_INFO, "  set dt" );
398
399     SG_LOG( SG_FLIGHT, SG_INFO, "Finished initializing JSBSim" );
400
401     SG_LOG( SG_FLIGHT, SG_INFO, "FGControls::get_gear_down()= " <<
402                                   globals->get_controls()->get_gear_down() );
403 }
404
405 /******************************************************************************/
406
407 // Run an iteration of the EOM (equations of motion)
408
409 void FGJSBsim::update( double dt )
410 {
411     if(crashed) {
412       if(!fgGetBool("/sim/crashed"))
413         fgSetBool("/sim/crashed", true);
414       return;
415     }
416
417     if (is_suspended())
418       return;
419
420     int multiloop = _calc_multiloop(dt);
421
422     int i;
423
424     // Compute the radius of the aircraft. That is the radius of a ball
425     // where all gear units are in. At the moment it is at least 10ft ...
426     double acrad = 10.0;
427     int n_gears = GroundReactions->GetNumGearUnits();
428     for (i=0; i<n_gears; ++i) {
429       FGColumnVector3 bl = GroundReactions->GetGearUnit(i)->GetBodyLocation();
430       double r = bl.Magnitude();
431       if (acrad < r)
432         acrad = r;
433     }
434
435     // Compute the potential movement of this aircraft and query for the
436     // ground in this area.
437     double groundCacheRadius = acrad + 2*dt*Propagate->GetUVW().Magnitude();
438     double alt, slr, lat, lon;
439     FGColumnVector3 cart = Auxiliary->GetLocationVRP();
440     if ( needTrim && startup_trim->getBoolValue() ) {
441       alt = fgic->GetAltitudeFtIC();
442       slr = fgic->GetSeaLevelRadiusFtIC();
443       lat = fgic->GetLatitudeDegIC() * SGD_DEGREES_TO_RADIANS;
444       lon = fgic->GetLongitudeDegIC() * SGD_DEGREES_TO_RADIANS;
445       cart = FGLocation(lon, lat, alt+slr);
446     }
447     double cart_pos[3] = { cart(1), cart(2), cart(3) };
448     double t0 = State->Getsim_time();
449     bool cache_ok = prepare_ground_cache_ft( t0, t0 + dt, cart_pos,
450                                              groundCacheRadius );
451     if (!cache_ok) {
452       SG_LOG(SG_FLIGHT, SG_WARN,
453              "FGInterface is being called without scenery below the aircraft!");
454
455       alt = fgic->GetAltitudeFtIC();
456       SG_LOG(SG_FLIGHT, SG_WARN, "altitude         = " << alt);
457
458       slr = fgic->GetSeaLevelRadiusFtIC();
459       SG_LOG(SG_FLIGHT, SG_WARN, "sea level radius = " << slr);
460
461       lat = fgic->GetLatitudeDegIC() * SGD_DEGREES_TO_RADIANS;
462       SG_LOG(SG_FLIGHT, SG_WARN, "latitude         = " << lat);
463
464       lon = fgic->GetLongitudeDegIC() * SGD_DEGREES_TO_RADIANS;
465       SG_LOG(SG_FLIGHT, SG_WARN, "longitude        = " << lon);
466       //return;
467     }
468
469     copy_to_JSBsim();
470
471     trimmed->setBoolValue(false);
472
473     if ( needTrim ) {
474       if ( startup_trim->getBoolValue() ) {
475         double contact[3], d[3], agl;
476         get_agl_ft(State->Getsim_time(), cart_pos, SG_METER_TO_FEET*2, contact,
477                    d, d, &agl);
478         double terrain_alt = sqrt(contact[0]*contact[0] + contact[1]*contact[1]
479              + contact[2]*contact[2]) - fgic->GetSeaLevelRadiusFtIC();
480
481         SG_LOG(SG_FLIGHT, SG_INFO,
482           "Ready to trim, terrain altitude is: "
483             << terrain_alt * SG_METER_TO_FEET );
484
485         fgic->SetTerrainAltitudeFtIC( terrain_alt );
486         do_trim();
487       } else {
488         fdmex->RunIC();  //apply any changes made through the set_ functions
489       }
490       needTrim = false;
491     }
492
493     for ( i=0; i < multiloop; i++ ) {
494       fdmex->Run();
495       update_external_forces(State->Getsim_time() + i * State->Getdt());      
496     }
497
498     FGJSBBase::Message* msg;
499     while (fdmex->SomeMessages()) {
500       msg = fdmex->ProcessMessage();
501       switch (msg->type) {
502       case FGJSBBase::Message::eText:
503         if (msg->text == "Crash Detected: Simulation FREEZE.")
504           crashed = true;
505         SG_LOG( SG_FLIGHT, SG_INFO, msg->messageId << ": " << msg->text );
506         break;
507       case FGJSBBase::Message::eBool:
508         SG_LOG( SG_FLIGHT, SG_INFO, msg->messageId << ": " << msg->text << " " << msg->bVal );
509         break;
510       case FGJSBBase::Message::eInteger:
511         SG_LOG( SG_FLIGHT, SG_INFO, msg->messageId << ": " << msg->text << " " << msg->iVal );
512         break;
513       case FGJSBBase::Message::eDouble:
514         SG_LOG( SG_FLIGHT, SG_INFO, msg->messageId << ": " << msg->text << " " << msg->dVal );
515         break;
516       default:
517         SG_LOG( SG_FLIGHT, SG_INFO, "Unrecognized message type." );
518         break;
519       }
520     }
521
522     // translate JSBsim back to FG structure so that the
523     // autopilot (and the rest of the sim can use the updated values
524     copy_from_JSBsim();
525 }
526
527 /******************************************************************************/
528
529 // Convert from the FGInterface struct to the JSBsim generic_ struct
530
531 bool FGJSBsim::copy_to_JSBsim()
532 {
533     double tmp;
534     unsigned int i;
535
536     // copy control positions into the JSBsim structure
537
538     FCS->SetDaCmd( globals->get_controls()->get_aileron());
539     FCS->SetRollTrimCmd( globals->get_controls()->get_aileron_trim() );
540     FCS->SetDeCmd( globals->get_controls()->get_elevator());
541     FCS->SetPitchTrimCmd( globals->get_controls()->get_elevator_trim() );
542     FCS->SetDrCmd( -globals->get_controls()->get_rudder() );
543     FCS->SetYawTrimCmd( -globals->get_controls()->get_rudder_trim() );
544     FCS->SetDsCmd( globals->get_controls()->get_rudder() );
545     FCS->SetDfCmd( globals->get_controls()->get_flaps() );
546     FCS->SetDsbCmd( globals->get_controls()->get_speedbrake() );
547     FCS->SetDspCmd( globals->get_controls()->get_spoilers() );
548
549         // Parking brake sets minimum braking
550         // level for mains.
551     double parking_brake = globals->get_controls()->get_brake_parking();
552     FCS->SetLBrake(FMAX(globals->get_controls()->get_brake_left(), parking_brake));
553     FCS->SetRBrake(FMAX(globals->get_controls()->get_brake_right(), parking_brake));
554     FCS->SetCBrake( 0.0 );
555     // FCS->SetCBrake( globals->get_controls()->get_brake(2) );
556
557     FCS->SetGearCmd( globals->get_controls()->get_gear_down());
558     for (i = 0; i < Propulsion->GetNumEngines(); i++) {
559       SGPropertyNode * node = fgGetNode("engines/engine", i, true);
560
561       FCS->SetThrottleCmd(i, globals->get_controls()->get_throttle(i));
562       FCS->SetMixtureCmd(i, globals->get_controls()->get_mixture(i));
563       FCS->SetPropAdvanceCmd(i, globals->get_controls()->get_prop_advance(i));
564       FCS->SetFeatherCmd(i, globals->get_controls()->get_feather(i));
565
566       switch (Propulsion->GetEngine(i)->GetType()) {
567       case FGEngine::etPiston:
568         { // FGPiston code block
569         FGPiston* eng = (FGPiston*)Propulsion->GetEngine(i);
570         eng->SetMagnetos( globals->get_controls()->get_magnetos(i) );
571         break;
572         } // end FGPiston code block
573       case FGEngine::etTurbine:
574         { // FGTurbine code block
575         FGTurbine* eng = (FGTurbine*)Propulsion->GetEngine(i);
576         eng->SetAugmentation( globals->get_controls()->get_augmentation(i) );
577         eng->SetReverse( globals->get_controls()->get_reverser(i) );
578         //eng->SetInjection( globals->get_controls()->get_water_injection(i) );
579         eng->SetCutoff( globals->get_controls()->get_cutoff(i) );
580         eng->SetIgnition( globals->get_controls()->get_ignition(i) );
581         break;
582         } // end FGTurbine code block
583       case FGEngine::etRocket:
584         { // FGRocket code block
585         FGRocket* eng = (FGRocket*)Propulsion->GetEngine(i);
586         break;
587         } // end FGRocket code block
588       case FGEngine::etTurboprop:
589         { // FGTurboProp code block
590         FGTurboProp* eng = (FGTurboProp*)Propulsion->GetEngine(i);
591         eng->SetReverse( globals->get_controls()->get_reverser(i) );
592         eng->SetCutoff( globals->get_controls()->get_cutoff(i) );
593         eng->SetIgnition( globals->get_controls()->get_ignition(i) );
594
595         eng->SetGeneratorPower( globals->get_controls()->get_generator_breaker(i) );
596         eng->SetCondition( globals->get_controls()->get_condition(i) );
597         break;
598         } // end FGTurboProp code block
599       default:
600         break;
601       }
602
603       { // FGEngine code block
604       FGEngine* eng = Propulsion->GetEngine(i);
605
606       eng->SetStarter( globals->get_controls()->get_starter(i) );
607       eng->SetRunning( node->getBoolValue("running") );
608       } // end FGEngine code block
609     }
610
611
612     Propagate->SetSeaLevelRadius( get_Sea_level_radius() );
613
614     Atmosphere->SetExTemperature(
615                   9.0/5.0*(temperature->getDoubleValue()+273.15) );
616     Atmosphere->SetExPressure(pressure->getDoubleValue()*70.726566);
617     Atmosphere->SetExDensity(density->getDoubleValue());
618
619     tmp = turbulence_gain->getDoubleValue();
620     //Atmosphere->SetTurbGain(tmp * tmp * 100.0);
621
622     tmp = turbulence_rate->getDoubleValue();
623     //Atmosphere->SetTurbRate(tmp);
624
625     Atmosphere->SetWindNED( -wind_from_north->getDoubleValue(),
626                             -wind_from_east->getDoubleValue(),
627                             -wind_from_down->getDoubleValue() );
628 //    SG_LOG(SG_FLIGHT,SG_INFO, "Wind NED: "
629 //                  << get_V_north_airmass() << ", "
630 //                  << get_V_east_airmass()  << ", "
631 //                  << get_V_down_airmass() );
632
633     for (i = 0; i < Propulsion->GetNumTanks(); i++) {
634       SGPropertyNode * node = fgGetNode("/consumables/fuel/tank", i, true);
635       FGTank * tank = Propulsion->GetTank(i);
636       tank->SetContents(node->getDoubleValue("level-gal_us") * 6.6);
637 //       tank->SetContents(node->getDoubleValue("level-lbs"));
638     }
639
640     Propulsion->SetFuelFreeze((fgGetNode("/sim/freeze/fuel",true))->getBoolValue());
641     fdmex->SetSlave(slaved->getBoolValue());
642
643     return true;
644 }
645
646 /******************************************************************************/
647
648 // Convert from the JSBsim generic_ struct to the FGInterface struct
649
650 bool FGJSBsim::copy_from_JSBsim()
651 {
652     unsigned int i, j;
653 /*
654     _set_Inertias( MassBalance->GetMass(),
655                    MassBalance->GetIxx(),
656                    MassBalance->GetIyy(),
657                    MassBalance->GetIzz(),
658                    MassBalance->GetIxz() );
659 */
660     _set_CG_Position( MassBalance->GetXYZcg(1),
661                       MassBalance->GetXYZcg(2),
662                       MassBalance->GetXYZcg(3) );
663
664     _set_Accels_Body( Aircraft->GetBodyAccel(1),
665                       Aircraft->GetBodyAccel(2),
666                       Aircraft->GetBodyAccel(3) );
667
668     _set_Accels_CG_Body_N ( Aircraft->GetNcg(1),
669                             Aircraft->GetNcg(2),
670                             Aircraft->GetNcg(3) );
671
672     _set_Accels_Pilot_Body( Auxiliary->GetPilotAccel(1),
673                             Auxiliary->GetPilotAccel(2),
674                             Auxiliary->GetPilotAccel(3) );
675
676     _set_Nlf( Aircraft->GetNlf() );
677
678     // Velocities
679
680     _set_Velocities_Local( Propagate->GetVel(FGJSBBase::eNorth),
681                            Propagate->GetVel(FGJSBBase::eEast),
682                            Propagate->GetVel(FGJSBBase::eDown) );
683
684     _set_Velocities_Wind_Body( Propagate->GetUVW(1),
685                                Propagate->GetUVW(2),
686                                Propagate->GetUVW(3) );
687
688     // Make the HUD work ...
689     _set_Velocities_Ground( Propagate->GetVel(FGJSBBase::eNorth),
690                             Propagate->GetVel(FGJSBBase::eEast),
691                             -Propagate->GetVel(FGJSBBase::eDown) );
692
693     _set_V_rel_wind( Auxiliary->GetVt() );
694
695     _set_V_equiv_kts( Auxiliary->GetVequivalentKTS() );
696
697     _set_V_calibrated_kts( Auxiliary->GetVcalibratedKTS() );
698
699     _set_V_ground_speed( Auxiliary->GetVground() );
700
701     _set_Omega_Body( Propagate->GetPQR(FGJSBBase::eP),
702                      Propagate->GetPQR(FGJSBBase::eQ),
703                      Propagate->GetPQR(FGJSBBase::eR) );
704
705     _set_Euler_Rates( Auxiliary->GetEulerRates(FGJSBBase::ePhi),
706                       Auxiliary->GetEulerRates(FGJSBBase::eTht),
707                       Auxiliary->GetEulerRates(FGJSBBase::ePsi) );
708
709     _set_Mach_number( Auxiliary->GetMach() );
710
711     // Positions of Visual Reference Point
712     FGLocation l = Auxiliary->GetLocationVRP();
713     _updateGeocentricPosition( l.GetLatitude(), l.GetLongitude(),
714                                l.GetRadius() - get_Sea_level_radius() );
715
716     _set_Altitude_AGL( Propagate->GetDistanceAGL() );
717     {
718       double loc_cart[3] = { l(FGJSBBase::eX), l(FGJSBBase::eY), l(FGJSBBase::eZ) };
719       double contact[3], d[3], sd, t;
720       is_valid_m(&t, d, &sd);
721       get_agl_ft(t, loc_cart, SG_METER_TO_FEET*2, contact, d, d, &sd);
722       double rwrad
723         = FGColumnVector3( contact[0], contact[1], contact[2] ).Magnitude();
724       _set_Runway_altitude( rwrad - get_Sea_level_radius() );
725     }
726
727     _set_Euler_Angles( Propagate->GetEuler(FGJSBBase::ePhi),
728                        Propagate->GetEuler(FGJSBBase::eTht),
729                        Propagate->GetEuler(FGJSBBase::ePsi) );
730
731     _set_Alpha( Auxiliary->Getalpha() );
732     _set_Beta( Auxiliary->Getbeta() );
733
734
735     _set_Gamma_vert_rad( Auxiliary->GetGamma() );
736
737     _set_Earth_position_angle( Inertial->GetEarthPositionAngle() );
738
739     _set_Climb_Rate( Propagate->Gethdot() );
740
741     const FGMatrix33& Tl2b = Propagate->GetTl2b();
742     for ( i = 1; i <= 3; i++ ) {
743         for ( j = 1; j <= 3; j++ ) {
744             _set_T_Local_to_Body( i, j, Tl2b(i,j) );
745         }
746     }
747
748     // Copy the engine values from JSBSim.
749     for ( i=0; i < Propulsion->GetNumEngines(); i++ ) {
750       SGPropertyNode * node = fgGetNode("engines/engine", i, true);
751       char buf[30];
752       sprintf(buf, "engines/engine[%d]/thruster", i);
753       SGPropertyNode * tnode = fgGetNode(buf, true);
754       FGThruster * thruster = Propulsion->GetEngine(i)->GetThruster();
755
756       switch (Propulsion->GetEngine(i)->GetType()) {
757       case FGEngine::etPiston:
758         { // FGPiston code block
759         FGPiston* eng = (FGPiston*)Propulsion->GetEngine(i);
760         node->setDoubleValue("egt-degf", eng->getExhaustGasTemp_degF());
761         node->setDoubleValue("oil-temperature-degf", eng->getOilTemp_degF());
762         node->setDoubleValue("oil-pressure-psi", eng->getOilPressure_psi());
763         node->setDoubleValue("mp-osi", eng->getManifoldPressure_inHg());
764         // NOTE: mp-osi is not in ounces per square inch.
765         // This error is left for reasons of backwards compatibility with
766         // existing FlightGear sound and instrument configurations.
767         node->setDoubleValue("mp-inhg", eng->getManifoldPressure_inHg());
768         node->setDoubleValue("cht-degf", eng->getCylinderHeadTemp_degF());
769         node->setDoubleValue("rpm", eng->getRPM());
770         } // end FGPiston code block
771         break;
772       case FGEngine::etRocket:
773         { // FGRocket code block
774         FGRocket* eng = (FGRocket*)Propulsion->GetEngine(i);
775         } // end FGRocket code block
776         break;
777       case FGEngine::etTurbine:
778         { // FGTurbine code block
779         FGTurbine* eng = (FGTurbine*)Propulsion->GetEngine(i);
780         node->setDoubleValue("n1", eng->GetN1());
781         node->setDoubleValue("n2", eng->GetN2());
782         node->setDoubleValue("egt_degf", 32 + eng->GetEGT()*9/5);
783         node->setBoolValue("augmentation", eng->GetAugmentation());
784         node->setBoolValue("water-injection", eng->GetInjection());
785         node->setBoolValue("ignition", eng->GetIgnition());
786         node->setDoubleValue("nozzle-pos-norm", eng->GetNozzle());
787         node->setDoubleValue("inlet-pos-norm", eng->GetInlet());
788         node->setDoubleValue("oil-pressure-psi", eng->getOilPressure_psi());
789         node->setBoolValue("reversed", eng->GetReversed());
790         node->setBoolValue("cutoff", eng->GetCutoff());
791         node->setDoubleValue("epr", eng->GetEPR());
792         globals->get_controls()->set_reverser(i, eng->GetReversed() );
793         globals->get_controls()->set_cutoff(i, eng->GetCutoff() );
794         globals->get_controls()->set_water_injection(i, eng->GetInjection() );
795         globals->get_controls()->set_augmentation(i, eng->GetAugmentation() );
796         } // end FGTurbine code block
797         break;
798       case FGEngine::etTurboprop:
799         { // FGTurboProp code block
800         FGTurboProp* eng = (FGTurboProp*)Propulsion->GetEngine(i);
801         node->setDoubleValue("n1", eng->GetN1());
802         //node->setDoubleValue("n2", eng->GetN2());
803         node->setDoubleValue("itt_degf", 32 + eng->GetITT()*9/5);
804         node->setBoolValue("ignition", eng->GetIgnition());
805         node->setDoubleValue("nozzle-pos-norm", eng->GetNozzle());
806         node->setDoubleValue("inlet-pos-norm", eng->GetInlet());
807         node->setDoubleValue("oil-pressure-psi", eng->getOilPressure_psi());
808         node->setBoolValue("reversed", eng->GetReversed());
809         node->setBoolValue("cutoff", eng->GetCutoff());
810         node->setBoolValue("starting", eng->GetEngStarting());
811         node->setBoolValue("generator-power", eng->GetGeneratorPower());
812         node->setBoolValue("damaged", eng->GetCondition());
813         node->setBoolValue("ielu-intervent", eng->GetIeluIntervent());
814         node->setDoubleValue("oil-temperature-degf", eng->getOilTemp_degF());
815 //        node->setBoolValue("onfire", eng->GetFire());
816         globals->get_controls()->set_reverser(i, eng->GetReversed() );
817         globals->get_controls()->set_cutoff(i, eng->GetCutoff() );
818         } // end FGTurboProp code block
819         break;
820       case FGEngine::etElectric:
821         { // FGElectric code block
822         FGElectric* eng = (FGElectric*)Propulsion->GetEngine(i);
823         node->setDoubleValue("rpm", eng->getRPM());
824         } // end FGElectric code block
825         break;
826       }
827
828       { // FGEngine code block
829       FGEngine* eng = Propulsion->GetEngine(i);
830       node->setDoubleValue("fuel-flow-gph", eng->getFuelFlow_gph());
831       node->setDoubleValue("thrust_lb", thruster->GetThrust());
832       node->setDoubleValue("fuel-flow_pph", eng->getFuelFlow_pph());
833       node->setBoolValue("running", eng->GetRunning());
834       node->setBoolValue("starter", eng->GetStarter());
835       node->setBoolValue("cranking", eng->GetCranking());
836       globals->get_controls()->set_starter(i, eng->GetStarter() );
837       } // end FGEngine code block
838
839       switch (thruster->GetType()) {
840       case FGThruster::ttNozzle:
841         { // FGNozzle code block
842         FGNozzle* noz = (FGNozzle*)thruster;
843         } // end FGNozzle code block
844         break;
845       case FGThruster::ttPropeller:
846         { // FGPropeller code block
847         FGPropeller* prop = (FGPropeller*)thruster;
848         tnode->setDoubleValue("rpm", thruster->GetRPM());
849         tnode->setDoubleValue("pitch", prop->GetPitch());
850         tnode->setDoubleValue("torque", prop->GetTorque());
851         tnode->setBoolValue("feathered", prop->GetFeather());
852         } // end FGPropeller code block
853         break;
854       case FGThruster::ttRotor:
855         { // FGRotor code block
856         FGRotor* rotor = (FGRotor*)thruster;
857         } // end FGRotor code block
858         break;
859       case FGThruster::ttDirect:
860         { // Direct code block
861         } // end Direct code block
862         break;
863       }
864
865     }
866
867     // Copy the fuel levels from JSBSim if fuel
868     // freeze not enabled.
869     if ( ! Propulsion->GetFuelFreeze() ) {
870       for (i = 0; i < Propulsion->GetNumTanks(); i++) {
871         SGPropertyNode * node = fgGetNode("/consumables/fuel/tank", i, true);
872         FGTank* tank = Propulsion->GetTank(i);
873         double contents = tank->GetContents();
874         double temp = tank->GetTemperature_degC();
875         node->setDoubleValue("level-gal_us", contents/6.6);
876         node->setDoubleValue("level-lbs", contents);
877         if (temp != -9999.0) node->setDoubleValue("temperature_degC", temp);
878       }
879     }
880
881     update_gear();
882
883     stall_warning->setDoubleValue( Aerodynamics->GetStallWarn() );
884
885     elevator_pos_pct->setDoubleValue( FCS->GetDePos(ofNorm) );
886     left_aileron_pos_pct->setDoubleValue( FCS->GetDaLPos(ofNorm) );
887     right_aileron_pos_pct->setDoubleValue( FCS->GetDaRPos(ofNorm) );
888     rudder_pos_pct->setDoubleValue( -1*FCS->GetDrPos(ofNorm) );
889     flap_pos_pct->setDoubleValue( FCS->GetDfPos(ofNorm) );
890     speedbrake_pos_pct->setDoubleValue( FCS->GetDsbPos(ofNorm) );
891     spoilers_pos_pct->setDoubleValue( FCS->GetDspPos(ofNorm) );
892     tailhook_pos_pct->setDoubleValue( FCS->GetTailhookPos() );
893     wing_fold_pos_pct->setDoubleValue( FCS->GetWingFoldPos() );
894
895     // force a sim crashed if crashed (altitude AGL < 0)
896     if (get_Altitude_AGL() < -100.0) {
897          State->SuspendIntegration();
898          crashed = true;
899     }
900
901     return true;
902 }
903
904
905 bool FGJSBsim::ToggleDataLogging(void)
906 {
907   // ToDo: handle this properly
908   fdmex->DisableOutput();
909   return false;
910 }
911
912
913 bool FGJSBsim::ToggleDataLogging(bool state)
914 {
915     if (state) {
916       fdmex->EnableOutput();
917       return true;
918     } else {
919       fdmex->DisableOutput();
920       return false;
921     }
922 }
923
924
925 //Positions
926 void FGJSBsim::set_Latitude(double lat)
927 {
928     static SGConstPropertyNode_ptr altitude = fgGetNode("/position/altitude-ft");
929     double alt;
930     double sea_level_radius_meters, lat_geoc;
931
932     // In case we're not trimming
933     FGInterface::set_Latitude(lat);
934
935     if ( altitude->getDoubleValue() > -9990 ) {
936       alt = altitude->getDoubleValue();
937     } else {
938       alt = 0.0;
939     }
940
941     update_ic();
942     SG_LOG(SG_FLIGHT,SG_INFO,"FGJSBsim::set_Latitude: " << lat );
943     SG_LOG(SG_FLIGHT,SG_INFO," cur alt (ft) =  " << alt );
944
945     sgGeodToGeoc( lat, alt * SG_FEET_TO_METER,
946                       &sea_level_radius_meters, &lat_geoc );
947     _set_Sea_level_radius( sea_level_radius_meters * SG_METER_TO_FEET  );
948     fgic->SetSeaLevelRadiusFtIC( sea_level_radius_meters * SG_METER_TO_FEET  );
949     fgic->SetLatitudeRadIC( lat_geoc );
950     needTrim=true;
951 }
952
953
954 void FGJSBsim::set_Longitude(double lon)
955 {
956     SG_LOG(SG_FLIGHT,SG_INFO,"FGJSBsim::set_Longitude: " << lon );
957
958     // In case we're not trimming
959     FGInterface::set_Longitude(lon);
960
961     update_ic();
962     fgic->SetLongitudeRadIC( lon );
963     needTrim=true;
964 }
965
966 void FGJSBsim::set_Altitude(double alt)
967 {
968     static SGConstPropertyNode_ptr latitude = fgGetNode("/position/latitude-deg");
969
970     double sea_level_radius_meters,lat_geoc;
971
972     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Altitude: " << alt );
973     SG_LOG(SG_FLIGHT,SG_INFO, "  lat (deg) = " << latitude->getDoubleValue() );
974
975     // In case we're not trimming
976     FGInterface::set_Altitude(alt);
977
978     update_ic();
979     sgGeodToGeoc( latitude->getDoubleValue() * SGD_DEGREES_TO_RADIANS, alt,
980                   &sea_level_radius_meters, &lat_geoc);
981     _set_Sea_level_radius( sea_level_radius_meters * SG_METER_TO_FEET  );
982     fgic->SetSeaLevelRadiusFtIC( sea_level_radius_meters * SG_METER_TO_FEET );
983     SG_LOG(SG_FLIGHT, SG_INFO,
984           "Terrain altitude: " << cur_fdm_state->get_Runway_altitude() * SG_METER_TO_FEET );
985     fgic->SetLatitudeRadIC( lat_geoc );
986     fgic->SetAltitudeFtIC(alt);
987     needTrim=true;
988 }
989
990 void FGJSBsim::set_V_calibrated_kts(double vc)
991 {
992     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_V_calibrated_kts: " <<  vc );
993
994     // In case we're not trimming
995     FGInterface::set_V_calibrated_kts(vc);
996
997     update_ic();
998     fgic->SetVcalibratedKtsIC(vc);
999     needTrim=true;
1000 }
1001
1002 void FGJSBsim::set_Mach_number(double mach)
1003 {
1004     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Mach_number: " <<  mach );
1005
1006     // In case we're not trimming
1007     FGInterface::set_Mach_number(mach);
1008
1009     update_ic();
1010     fgic->SetMachIC(mach);
1011     needTrim=true;
1012 }
1013
1014 void FGJSBsim::set_Velocities_Local( double north, double east, double down )
1015 {
1016     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Velocities_Local: "
1017        << north << ", " <<  east << ", " << down );
1018
1019     // In case we're not trimming
1020     FGInterface::set_Velocities_Local(north, east, down);
1021
1022     update_ic();
1023     fgic->SetVNorthFpsIC(north);
1024     fgic->SetVEastFpsIC(east);
1025     fgic->SetVDownFpsIC(down);
1026     needTrim=true;
1027 }
1028
1029 void FGJSBsim::set_Velocities_Wind_Body( double u, double v, double w)
1030 {
1031     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Velocities_Wind_Body: "
1032        << u << ", " <<  v << ", " <<  w );
1033
1034     // In case we're not trimming
1035     FGInterface::set_Velocities_Wind_Body(u, v, w);
1036
1037     update_ic();
1038     fgic->SetUBodyFpsIC(u);
1039     fgic->SetVBodyFpsIC(v);
1040     fgic->SetWBodyFpsIC(w);
1041     needTrim=true;
1042 }
1043
1044 //Euler angles
1045 void FGJSBsim::set_Euler_Angles( double phi, double theta, double psi )
1046 {
1047     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Euler_Angles: "
1048        << phi << ", " << theta << ", " << psi );
1049
1050     // In case we're not trimming
1051     FGInterface::set_Euler_Angles(phi, theta, psi);
1052
1053     update_ic();
1054     fgic->SetThetaRadIC(theta);
1055     fgic->SetPhiRadIC(phi);
1056     fgic->SetPsiRadIC(psi);
1057     needTrim=true;
1058 }
1059
1060 //Flight Path
1061 void FGJSBsim::set_Climb_Rate( double roc)
1062 {
1063     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Climb_Rate: " << roc );
1064
1065     // In case we're not trimming
1066     FGInterface::set_Climb_Rate(roc);
1067
1068     update_ic();
1069     //since both climb rate and flight path angle are set in the FG
1070     //startup sequence, something is needed to keep one from cancelling
1071     //out the other.
1072     if( !(fabs(roc) > 1 && fabs(fgic->GetFlightPathAngleRadIC()) < 0.01) ) {
1073       fgic->SetClimbRateFpsIC(roc);
1074     }
1075     needTrim=true;
1076 }
1077
1078 void FGJSBsim::set_Gamma_vert_rad( double gamma)
1079 {
1080     SG_LOG(SG_FLIGHT,SG_INFO, "FGJSBsim::set_Gamma_vert_rad: " << gamma );
1081
1082     update_ic();
1083     if( !(fabs(gamma) < 0.01 && fabs(fgic->GetClimbRateFpsIC()) > 1) ) {
1084       fgic->SetFlightPathAngleRadIC(gamma);
1085     }
1086     needTrim=true;
1087 }
1088
1089 void FGJSBsim::init_gear(void )
1090 {
1091     FGGroundReactions* gr=fdmex->GetGroundReactions();
1092     int Ngear=GroundReactions->GetNumGearUnits();
1093     for (int i=0;i<Ngear;i++) {
1094       FGLGear *gear = gr->GetGearUnit(i);
1095       SGPropertyNode * node = fgGetNode("gear/gear", i, true);
1096       node->setDoubleValue("xoffset-in", gear->GetBodyLocation()(1));
1097       node->setDoubleValue("yoffset-in", gear->GetBodyLocation()(2));
1098       node->setDoubleValue("zoffset-in", gear->GetBodyLocation()(3));
1099       node->setBoolValue("wow", gear->GetWOW());
1100       node->setDoubleValue("rollspeed-ms", gear->GetWheelRollVel()*0.3043);
1101       node->setBoolValue("has-brake", gear->GetBrakeGroup() > 0);
1102       node->setDoubleValue("position-norm", gear->GetGearUnitPos());
1103       node->setDoubleValue("tire-pressure-norm", gear->GetTirePressure());
1104       node->setDoubleValue("compression-norm", gear->GetCompLen());
1105       node->setDoubleValue("compression-ft", gear->GetCompLen());
1106       if ( gear->GetSteerable() )
1107         node->setDoubleValue("steering-norm", gear->GetSteerNorm());
1108     }
1109 }
1110
1111 void FGJSBsim::update_gear(void)
1112 {
1113     FGGroundReactions* gr=fdmex->GetGroundReactions();
1114     int Ngear=GroundReactions->GetNumGearUnits();
1115     for (int i=0;i<Ngear;i++) {
1116       FGLGear *gear = gr->GetGearUnit(i);
1117       SGPropertyNode * node = fgGetNode("gear/gear", i, true);
1118       node->getChild("wow", 0, true)->setBoolValue( gear->GetWOW());
1119       node->getChild("rollspeed-ms", 0, true)->setDoubleValue(gear->GetWheelRollVel()*0.3043);
1120       node->getChild("position-norm", 0, true)->setDoubleValue(gear->GetGearUnitPos());
1121       gear->SetTirePressure(node->getDoubleValue("tire-pressure-norm"));
1122       node->setDoubleValue("compression-norm", gear->GetCompLen());
1123       node->setDoubleValue("compression-ft", gear->GetCompLen());
1124       if ( gear->GetSteerable() )
1125         node->setDoubleValue("steering-norm", gear->GetSteerNorm());
1126     }
1127 }
1128
1129 void FGJSBsim::do_trim(void)
1130 {
1131   FGTrim *fgtrim;
1132
1133   if ( fgGetBool("/sim/presets/onground") )
1134   {
1135     fgtrim = new FGTrim(fdmex,tGround);
1136   } else {
1137     fgtrim = new FGTrim(fdmex,tLongitudinal);
1138   }
1139
1140   if ( !fgtrim->DoTrim() ) {
1141     fgtrim->Report();
1142     fgtrim->TrimStats();
1143   } else {
1144     trimmed->setBoolValue(true);
1145   }
1146 //  if (FGJSBBase::debug_lvl > 0)
1147 //      State->ReportState();
1148
1149   delete fgtrim;
1150
1151   pitch_trim->setDoubleValue( FCS->GetPitchTrimCmd() );
1152   throttle_trim->setDoubleValue( FCS->GetThrottleCmd(0) );
1153   aileron_trim->setDoubleValue( FCS->GetDaCmd() );
1154   rudder_trim->setDoubleValue( FCS->GetDrCmd() );
1155
1156   globals->get_controls()->set_elevator_trim(FCS->GetPitchTrimCmd());
1157   globals->get_controls()->set_elevator(FCS->GetDeCmd());
1158   globals->get_controls()->set_throttle(FGControls::ALL_ENGINES,
1159   FCS->GetThrottleCmd(0));
1160
1161   globals->get_controls()->set_aileron(FCS->GetDaCmd());
1162   globals->get_controls()->set_rudder( FCS->GetDrCmd());
1163
1164   SG_LOG( SG_FLIGHT, SG_INFO, "  Trim complete" );
1165 }
1166
1167 void FGJSBsim::update_ic(void)
1168 {
1169    if ( !needTrim ) {
1170      fgic->SetLatitudeRadIC(get_Lat_geocentric() );
1171      fgic->SetLongitudeRadIC( get_Longitude() );
1172      fgic->SetAltitudeFtIC( get_Altitude() );
1173      fgic->SetVcalibratedKtsIC( get_V_calibrated_kts() );
1174      fgic->SetThetaRadIC( get_Theta() );
1175      fgic->SetPhiRadIC( get_Phi() );
1176      fgic->SetPsiRadIC( get_Psi() );
1177      fgic->SetClimbRateFpsIC( get_Climb_Rate() );
1178    }
1179 }
1180
1181 bool
1182 FGJSBsim::get_agl_ft(double t, const double pt[3], double alt_off,
1183                      double contact[3], double normal[3], double vel[3],
1184                      double *agl)
1185 {
1186    double angularVel[3];
1187    const SGMaterial* material;
1188    simgear::BVHNode::Id id;
1189    if (!FGInterface::get_agl_ft(t, pt, alt_off, contact, normal, vel,
1190                                 angularVel, material, id))
1191        return false;
1192    SGGeod geodPt = SGGeod::fromCart(SG_FEET_TO_METER*SGVec3d(pt));
1193    SGQuatd hlToEc = SGQuatd::fromLonLat(geodPt);
1194    *agl = dot(hlToEc.rotate(SGVec3d(0, 0, 1)), SGVec3d(contact) - SGVec3d(pt));
1195    return true;
1196 }
1197
1198 inline static double dot3(const FGColumnVector3& a, const FGColumnVector3& b)
1199 {
1200     return a(1) * b(1) + a(2) * b(2) + a(3) * b(3);
1201 }
1202
1203 inline static double sqr(double x)
1204 {
1205     return x * x;
1206 }
1207
1208 static double angle_diff(double a, double b)
1209 {
1210     double diff = fabs(a - b);
1211     if (diff > 180) diff = 360 - diff;
1212     
1213     return diff;
1214 }
1215
1216 static void check_hook_solution(const FGColumnVector3& ground_normal_body, double E, double hook_length, double sin_fi_guess, double cos_fi_guess, double* sin_fis, double* cos_fis, double* fis, int* points)
1217 {
1218     FGColumnVector3 tip(-hook_length * cos_fi_guess, 0, hook_length * sin_fi_guess);
1219     double dist = dot3(tip, ground_normal_body);
1220     if (fabs(dist + E) < 0.0001) {
1221         sin_fis[*points] = sin_fi_guess;
1222         cos_fis[*points] = cos_fi_guess;
1223         fis[*points] = atan2(sin_fi_guess, cos_fi_guess) * SG_RADIANS_TO_DEGREES;
1224         (*points)++;
1225     } 
1226 }
1227
1228
1229 static void check_hook_solution(const FGColumnVector3& ground_normal_body, double E, double hook_length, double sin_fi_guess, double* sin_fis, double* cos_fis, double* fis, int* points)
1230 {
1231     if (sin_fi_guess >= -1 && sin_fi_guess <= 1) {
1232         double cos_fi_guess = sqrt(1 - sqr(sin_fi_guess));
1233         check_hook_solution(ground_normal_body, E, hook_length, sin_fi_guess, cos_fi_guess, sin_fis, cos_fis, fis, points);
1234         if (fabs(cos_fi_guess) > SG_EPSILON) {
1235             check_hook_solution(ground_normal_body, E, hook_length, sin_fi_guess, -cos_fi_guess, sin_fis, cos_fis, fis, points);
1236         }
1237     }
1238 }
1239
1240 void FGJSBsim::update_external_forces(double t_off)
1241 {
1242     const FGMatrix33& Tb2l = Propagate->GetTb2l();
1243     const FGMatrix33& Tl2b = Propagate->GetTl2b();
1244     const FGLocation& Location = Propagate->GetLocation();
1245     const FGMatrix33& Tec2l = Location.GetTec2l();
1246         
1247     double hook_area[4][3];
1248     
1249     FGColumnVector3 hook_root_body = MassBalance->StructuralToBody(hook_root_struct);
1250     FGColumnVector3 hook_root = Location.LocalToLocation(Tb2l *   hook_root_body);
1251     hook_area[1][0] = hook_root(1);
1252     hook_area[1][1] = hook_root(2);
1253     hook_area[1][2] = hook_root(3);
1254     
1255     hook_length = fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-length-ft", 6.75);
1256     double fi_min = fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-pos-min-deg", -18);
1257     double fi_max = fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-pos-max-deg", 30);
1258     double fi = fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-pos-norm") * (fi_max - fi_min) + fi_min;
1259     double cos_fi = cos(fi * SG_DEGREES_TO_RADIANS);
1260     double sin_fi = sin(fi * SG_DEGREES_TO_RADIANS);
1261
1262     FGColumnVector3 hook_tip_body = hook_root_body;
1263     hook_tip_body(1) -= hook_length * cos_fi;
1264     hook_tip_body(3) += hook_length * sin_fi;    
1265     bool hook_tip_valid = true;
1266     
1267     double contact[3];
1268     double ground_normal[3];
1269     double ground_vel[3];
1270     double root_agl_ft;
1271
1272     if (!got_wire) {
1273         bool got = get_agl_ft(t_off, hook_area[1], 0, contact, ground_normal, ground_vel, &root_agl_ft);
1274         if (got && root_agl_ft > 0 && root_agl_ft < hook_length) {
1275             FGColumnVector3 ground_normal_body = Tl2b * (Tec2l * FGColumnVector3(ground_normal[0], ground_normal[1], ground_normal[2]));
1276             FGColumnVector3 contact_body = Tl2b * Location.LocationToLocal(FGColumnVector3(contact[0], contact[1], contact[2]));
1277             double D = -dot3(contact_body, ground_normal_body);
1278
1279             // check hook tip agl against same ground plane
1280             double hook_tip_agl_ft = dot3(hook_tip_body, ground_normal_body) + D;
1281             if (hook_tip_agl_ft < 0) {
1282
1283                 // hook tip: hx - l cos, hy, hz + l sin
1284                 // on ground:  - n0 l cos + n2 l sin + E = 0
1285
1286                 double E = D + dot3(hook_root_body, ground_normal_body);
1287
1288                 // substitue x = sin fi, cos fi = sqrt(1 - x * x)
1289                 // and rearrange to get a quadratic with coeffs:
1290                 double a = sqr(hook_length) * (sqr(ground_normal_body(1)) + sqr(ground_normal_body(3)));
1291                 double b = 2 * E * ground_normal_body(3) * hook_length;
1292                 double c = sqr(E) - sqr(ground_normal_body(1) * hook_length);   
1293
1294                 double disc = sqr(b) - 4 * a * c;
1295                 if (disc >= 0) {
1296                     double delta = sqrt(disc) / (2 * a);
1297                 
1298                     // allow 4 solutions for safety, should never happen
1299                     double sin_fis[4];
1300                     double cos_fis[4];
1301                     double fis[4];
1302                     int points = 0;
1303                 
1304                     double sin_fi_guess = -b / (2 * a) - delta;
1305                     check_hook_solution(ground_normal_body, E, hook_length, sin_fi_guess, sin_fis, cos_fis, fis, &points);
1306                     check_hook_solution(ground_normal_body, E, hook_length, sin_fi_guess + 2 * delta, sin_fis, cos_fis, fis, &points);
1307                 
1308                     if (points == 2) {
1309                         double diff1 = angle_diff(fi, fis[0]);
1310                         double diff2 = angle_diff(fi, fis[1]);
1311                         int point = diff1 < diff2 ? 0 : 1;
1312                         fi = fis[point];
1313                         sin_fi = sin_fis[point];
1314                         cos_fi = cos_fis[point];
1315                         hook_tip_body(1) = hook_root_body(1) - hook_length * cos_fi;
1316                         hook_tip_body(3) = hook_root_body(3) + hook_length * sin_fi;
1317                     }
1318                 }
1319             }
1320         }
1321     } else {
1322         FGColumnVector3 hook_root_vel = Propagate->GetVel() + (Tb2l * (Propagate->GetPQR() *  hook_root_body));
1323         double wire_ends_ec[2][3];
1324         double wire_vel_ec[2][3];
1325         get_wire_ends_ft(t_off, wire_ends_ec, wire_vel_ec);
1326         FGColumnVector3 wire_vel_1 = Tec2l * FGColumnVector3(wire_vel_ec[0][0], wire_vel_ec[0][1], wire_vel_ec[0][2]);
1327         FGColumnVector3 wire_vel_2 = Tec2l * FGColumnVector3(wire_vel_ec[1][0], wire_vel_ec[1][1], wire_vel_ec[1][2]);
1328         FGColumnVector3 rel_vel = hook_root_vel - (wire_vel_1 + wire_vel_2) / 2;
1329         if (rel_vel.Magnitude() < 3) {
1330             got_wire = false;
1331             release_wire();
1332             fgSetDouble("/fdm/jsbsim/external_reactions/hook/magnitude", 0.0);
1333         } else {
1334             FGColumnVector3 wire_end1_body = Tl2b * Location.LocationToLocal(FGColumnVector3(wire_ends_ec[0][0], wire_ends_ec[0][1], wire_ends_ec[0][2])) - hook_root_body;
1335             FGColumnVector3 wire_end2_body = Tl2b * Location.LocationToLocal(FGColumnVector3(wire_ends_ec[1][0], wire_ends_ec[1][1], wire_ends_ec[1][2])) - hook_root_body;
1336             FGColumnVector3 force_plane_normal = wire_end1_body * wire_end2_body;
1337             force_plane_normal.Normalize();
1338             cos_fi = dot3(force_plane_normal, FGColumnVector3(0, 0, 1));
1339             if (cos_fi < 0) cos_fi = -cos_fi;
1340             sin_fi = sqrt(1 - sqr(cos_fi));
1341             fi = atan2(sin_fi, cos_fi) * SG_RADIANS_TO_DEGREES;
1342         
1343             fgSetDouble("/fdm/jsbsim/external_reactions/hook/x", -cos_fi);
1344             fgSetDouble("/fdm/jsbsim/external_reactions/hook/y", 0);
1345             fgSetDouble("/fdm/jsbsim/external_reactions/hook/z", sin_fi);
1346             fgSetDouble("/fdm/jsbsim/external_reactions/hook/magnitude", fgGetDouble("/fdm/jsbsim/systems/hook/force"));
1347         }
1348     }
1349
1350     FGColumnVector3 hook_tip = Location.LocalToLocation(Tb2l * hook_tip_body);
1351
1352     hook_area[0][0] = hook_tip(1);
1353     hook_area[0][1] = hook_tip(2);
1354     hook_area[0][2] = hook_tip(3);
1355
1356     if (!got_wire) {
1357         // The previous positions.
1358         hook_area[2][0] = last_hook_root[0];
1359         hook_area[2][1] = last_hook_root[1];
1360         hook_area[2][2] = last_hook_root[2];
1361         hook_area[3][0] = last_hook_tip[0];
1362         hook_area[3][1] = last_hook_tip[1];
1363         hook_area[3][2] = last_hook_tip[2];
1364
1365         // Check if we caught a wire.
1366         // Returns true if we caught one.
1367         if (caught_wire_ft(t_off, hook_area)) {
1368                 got_wire = true;
1369         }
1370     }
1371     
1372     // save actual position as old position ...
1373     last_hook_tip[0] = hook_area[0][0];
1374     last_hook_tip[1] = hook_area[0][1];
1375     last_hook_tip[2] = hook_area[0][2];
1376     last_hook_root[0] = hook_area[1][0];
1377     last_hook_root[1] = hook_area[1][1];
1378     last_hook_root[2] = hook_area[1][2];
1379     
1380     fgSetDouble("/fdm/jsbsim/systems/hook/tailhook-pos-deg", fi);
1381 }
1382