1 /*******************************************************************************
7 --------- Copyright (C) 1999 Anthony K. Peden (apeden@earthlink.net) ---------
9 This program is free software; you can redistribute it and/or modify it under
10 the terms of the GNU General Public License as published by the Free Software
11 Foundation; either version 2 of the License, or (at your option) any later
14 This program is distributed in the hope that it will be useful, but WITHOUT
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
19 You should have received a copy of the GNU General Public License along with
20 this program; if not, write to the Free Software Foundation, Inc., 59 Temple
21 Place - Suite 330, Boston, MA 02111-1307, USA.
23 Further information about the GNU General Public License can also be found on
24 the world wide web at http://www.gnu.org.
28 --------------------------------------------------------------------------------
31 ********************************************************************************
33 *******************************************************************************/
38 #include "FGFDMExec.h"
39 #include "FGAtmosphere.h"
40 #include "FGInitialCondition.h"
41 #include "FGTrimAxis.h"
42 #include "FGAircraft.h"
44 static const char *IdSrc = "$Header$";
45 static const char *IdHdr = ID_TRIMAXIS;
47 /*****************************************************************************/
49 FGTrimAxis::FGTrimAxis(FGFDMExec* fdex, FGInitialCondition* ic, Accel acc,
50 Control ctrl, float ff) {
60 its_to_stable_value=0;
62 total_stability_iterations=0;
73 control_min=-30*DEGTORAD;
74 control_max=30*DEGTORAD;
75 control_convert=RADTODEG;
78 control_min=fdmex->GetAircraft()->GetAlphaCLMin();
79 control_max=fdmex->GetAircraft()->GetAlphaCLMax();
80 if(control_max <= control_min) {
81 control_max=20*DEGTORAD;
82 control_min=-5*DEGTORAD;
84 control_value= (control_min+control_max)/2;
85 control_convert=RADTODEG;
86 solver_eps=tolerance/100;
96 accel_convert=RADTODEG;
97 solver_eps=tolerance/100;
102 control_value=fdmex->GetPosition()->GetDistanceAGL();
103 solver_eps=tolerance/100;
106 control_min=-10*DEGTORAD;
107 control_max=10*DEGTORAD;
108 accel_convert=RADTODEG;
111 control_min=-30*DEGTORAD;
112 control_max=30*DEGTORAD;
113 accel_convert=RADTODEG;
116 solver_eps=tolerance/100;
117 control_min=-80*DEGTORAD;
118 control_max=80*DEGTORAD;
119 control_convert=RADTODEG;
125 /*****************************************************************************/
127 FGTrimAxis::~FGTrimAxis() {}
129 /*****************************************************************************/
131 void FGTrimAxis::getAccel(void) {
133 case tUdot: accel_value=fdmex -> GetTranslation()->GetUVWdot()(1); break;
134 case tVdot: accel_value=fdmex -> GetTranslation()->GetUVWdot()(2); break;
135 case tWdot: accel_value=fdmex -> GetTranslation()->GetUVWdot()(3); break;
136 case tQdot: accel_value=fdmex -> GetRotation()->GetPQRdot()(2);break;
137 case tPdot: accel_value=fdmex -> GetRotation()->GetPQRdot()(1); break;
138 case tRdot: accel_value=fdmex -> GetRotation()->GetPQRdot()(3); break;
142 /*****************************************************************************/
144 //Accels are not settable
146 void FGTrimAxis::getControl(void) {
148 case tThrottle: control_value=fdmex->GetFCS()->GetThrottleCmd(0); break;
149 case tBeta: control_value=fdmex->GetTranslation()->Getalpha(); break;
150 case tAlpha: control_value=fdmex->GetTranslation()->Getbeta(); break;
151 case tPitchTrim: control_value=fdmex->GetFCS() -> GetPitchTrimCmd(); break;
152 case tElevator: control_value=fdmex->GetFCS() -> GetDeCmd(); break;
154 case tAileron: control_value=fdmex->GetFCS() -> GetDaCmd(); break;
156 case tRudder: control_value=fdmex->GetFCS() -> GetDrCmd(); break;
157 case tAltAGL: control_value=fdmex->GetPosition()->GetDistanceAGL();break;
158 case tTheta: control_value=fdmex->GetRotation()->Gettht(); break;
159 case tPhi: control_value=fdmex->GetRotation()->Getphi(); break;
160 case tGamma: control_value=fdmex->GetPosition()->GetGamma();break;
164 /*****************************************************************************/
167 void FGTrimAxis::setControl(void) {
169 case tThrottle: setThrottlesPct(); break;
170 case tBeta: fgic->SetBetaRadIC(control_value); break;
171 case tAlpha: fgic->SetAlphaRadIC(control_value); break;
172 case tPitchTrim: fdmex->GetFCS() -> SetPitchTrimCmd(control_value); break;
173 case tElevator: fdmex-> GetFCS() -> SetDeCmd(control_value); break;
175 case tAileron: fdmex-> GetFCS() -> SetDaCmd(control_value); break;
177 case tRudder: fdmex-> GetFCS() -> SetDrCmd(control_value); break;
178 case tAltAGL: fgic->SetAltitudeAGLFtIC(control_value); break;
179 case tTheta: fgic->SetPitchAngleRadIC(control_value); break;
180 case tPhi: fgic->SetRollAngleRadIC(control_value); break;
181 case tGamma: fgic->SetFlightPathAngleRadIC(control_value); break;
185 /*****************************************************************************/
187 // the aircraft center of rotation is no longer the cg once the gear
188 // contact the ground so the altitude needs to be changed when pitch
189 // and roll angle are adjusted. Instead of attempting to calculate the
190 // new center of rotation, pick a gear unit as a reference and use its
191 // location vector to calculate the new height change. i.e. new altitude =
192 // earth z component of that vector (which is in body axes )
193 void FGTrimAxis::SetThetaOnGround(float ff) {
196 // favor an off-center unit so that the same one can be used for both
197 // pitch and roll. An on-center unit is used (for pitch)if that's all
198 // that's in contact with the ground.
199 i=0; ref=-1; center=-1;
200 while( (ref < 0) && (i < fdmex->GetAircraft()->GetNumGearUnits()) ) {
201 if(fdmex->GetAircraft()->GetGearUnit(i)->GetWOW()) {
202 if(fabs(fdmex->GetAircraft()->GetGearUnit(i)->GetBodyLocation()(2)) > 0.01)
209 if((ref < 0) && (center >= 0)) {
212 cout << "SetThetaOnGround ref gear: " << ref << endl;
214 float sp=fdmex->GetRotation()->GetSinphi();
215 float cp=fdmex->GetRotation()->GetCosphi();
216 float lx=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(1);
217 float ly=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(2);
218 float lz=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(3);
219 float hagl = -1*lx*sin(ff) +
223 fgic->SetAltitudeAGLFtIC(hagl);
224 cout << "SetThetaOnGround new alt: " << hagl << endl;
226 fgic->SetPitchAngleRadIC(ff);
227 cout << "SetThetaOnGround new theta: " << ff << endl;
230 /*****************************************************************************/
232 void FGTrimAxis::SetPhiOnGround(float ff) {
236 //must have an off-center unit here
237 while( (ref < 0) && (i < fdmex->GetAircraft()->GetNumGearUnits()) ) {
238 if( (fdmex->GetAircraft()->GetGearUnit(i)->GetWOW()) &&
239 (fabs(fdmex->GetAircraft()->GetGearUnit(i)->GetBodyLocation()(2)) > 0.01))
244 float st=fdmex->GetRotation()->GetSintht();
245 float ct=fdmex->GetRotation()->GetCostht();
246 float lx=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(1);
247 float ly=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(2);
248 float lz=fdmex->GetAircraft()->GetGearUnit(ref)->GetBodyLocation()(3);
249 float hagl = -1*lx*st +
253 fgic->SetAltitudeAGLFtIC(hagl);
255 fgic->SetRollAngleRadIC(ff);
259 /*****************************************************************************/
261 void FGTrimAxis::Run(void) {
263 float last_accel_value;
266 //cout << "FGTrimAxis::Run: " << control_value << endl;
271 last_accel_value=accel_value;
275 if((fabs(last_accel_value - accel_value) < tolerance) || (i >= 100) )
280 its_to_stable_value=i;
281 total_stability_iterations+=its_to_stable_value;
285 /*****************************************************************************/
287 void FGTrimAxis::setThrottlesPct(void) {
289 for(unsigned i=0;i<fdmex->GetAircraft()->GetNumEngines();i++) {
290 tMin=fdmex->GetAircraft()->GetEngine(i)->GetThrottleMin();
291 tMax=fdmex->GetAircraft()->GetEngine(i)->GetThrottleMax();
292 //cout << "setThrottlespct: " << i << ", " << control_min << ", " << control_max << ", " << control_value;
293 fdmex -> GetFCS() -> SetThrottleCmd(i,tMin+control_value*(tMax-tMin));
298 /*****************************************************************************/
301 void FGTrimAxis::AxisReport(void) {
304 sprintf(out," %20s: %6.2f %5s: %9.2e Tolerance: %3.0e\n",
305 GetControlName().c_str(), GetControl()*control_convert,
306 GetAccelName().c_str(), GetAccel(), GetTolerance());
312 /*****************************************************************************/
314 float FGTrimAxis::GetAvgStability( void ) {
315 if(total_iterations > 0) {
316 return float(total_stability_iterations)/float(total_iterations);