]> git.mxchange.org Git - flightgear.git/blob - src/FDM/JSBSim/models/FGMassBalance.cpp
Merge branch 'next' of gitorious.org:fg/flightgear into next
[flightgear.git] / src / FDM / JSBSim / models / FGMassBalance.cpp
1 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2
3  Module:       FGMassBalance.cpp
4  Author:       Jon S. Berndt
5  Date started: 09/12/2000
6  Purpose:      This module models weight and balance
7
8  ------------- Copyright (C) 2000  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
30 This class models the change in weight and balance of the aircraft due to fuel
31 burnoff, etc.
32
33 HISTORY
34 --------------------------------------------------------------------------------
35 09/12/2000  JSB  Created
36
37 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
38 INCLUDES
39 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
40
41 #include "FGMassBalance.h"
42 #include "FGPropulsion.h"
43 #include "propulsion/FGTank.h"
44 #include "FGBuoyantForces.h"
45 #include "input_output/FGPropertyManager.h"
46 #include <iostream>
47 #include <iomanip>
48 #include <cstdlib>
49
50 using namespace std;
51
52 namespace JSBSim {
53
54 static const char *IdSrc = "$Id: FGMassBalance.cpp,v 1.31 2010/02/19 00:30:00 jberndt Exp $";
55 static const char *IdHdr = ID_MASSBALANCE;
56
57 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
58 CLASS IMPLEMENTATION
59 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
60
61
62 FGMassBalance::FGMassBalance(FGFDMExec* fdmex) : FGModel(fdmex)
63 {
64   Name = "FGMassBalance";
65   Weight = EmptyWeight = Mass = 0.0;
66
67   vbaseXYZcg.InitMatrix(0.0);
68   vXYZcg.InitMatrix(0.0);
69   vLastXYZcg.InitMatrix(0.0);
70   vDeltaXYZcg.InitMatrix(0.0);
71   baseJ.InitMatrix();
72   mJ.InitMatrix();
73   mJinv.InitMatrix();
74   pmJ.InitMatrix();
75
76   bind();
77
78   Debug(0);
79 }
80
81 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
82
83 FGMassBalance::~FGMassBalance()
84 {
85   for (unsigned int i=0; i<PointMasses.size(); i++) delete PointMasses[i];
86   PointMasses.clear();
87
88   Debug(1);
89 }
90
91 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
92
93 bool FGMassBalance::InitModel(void)
94 {
95   if (!FGModel::InitModel()) return false;
96
97   vLastXYZcg.InitMatrix(0.0);
98   vDeltaXYZcg.InitMatrix(0.0);
99
100   return true;
101 }
102
103 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
104
105 bool FGMassBalance::Load(Element* el)
106 {
107   Element *element;
108   string element_name = "";
109   double bixx, biyy, bizz, bixy, bixz, biyz;
110
111   FGModel::Load(el); // Perform base class Load.
112
113   bixx = biyy = bizz = bixy = bixz = biyz = 0.0;
114   if (el->FindElement("ixx"))
115     bixx = el->FindElementValueAsNumberConvertTo("ixx", "SLUG*FT2");
116   if (el->FindElement("iyy"))
117     biyy = el->FindElementValueAsNumberConvertTo("iyy", "SLUG*FT2");
118   if (el->FindElement("izz"))
119     bizz = el->FindElementValueAsNumberConvertTo("izz", "SLUG*FT2");
120   if (el->FindElement("ixy"))
121     bixy = el->FindElementValueAsNumberConvertTo("ixy", "SLUG*FT2");
122   if (el->FindElement("ixz"))
123     bixz = el->FindElementValueAsNumberConvertTo("ixz", "SLUG*FT2");
124   if (el->FindElement("iyz"))
125     biyz = el->FindElementValueAsNumberConvertTo("iyz", "SLUG*FT2");
126   SetAircraftBaseInertias(FGMatrix33(  bixx,  -bixy,  bixz,
127                                       -bixy,  biyy,  -biyz,
128                                        bixz,  -biyz,  bizz ));
129   if (el->FindElement("emptywt")) {
130     EmptyWeight = el->FindElementValueAsNumberConvertTo("emptywt", "LBS");
131   }
132
133   element = el->FindElement("location");
134   while (element) {
135     element_name = element->GetAttributeValue("name");
136     if (element_name == "CG") vbaseXYZcg = element->FindElementTripletConvertTo("IN");
137     element = el->FindNextElement("location");
138   }
139
140 // Find all POINTMASS elements that descend from this METRICS branch of the
141 // config file.
142
143   element = el->FindElement("pointmass");
144   while (element) {
145     AddPointMass(element);
146     element = el->FindNextElement("pointmass");
147   }
148
149   double ChildFDMWeight = 0.0;
150   for (int fdm=0; fdm<FDMExec->GetFDMCount(); fdm++) {
151     if (FDMExec->GetChildFDM(fdm)->mated) ChildFDMWeight += FDMExec->GetChildFDM(fdm)->exec->GetMassBalance()->GetWeight();
152   }
153
154   Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetTotalPointMassWeight()
155     + BuoyantForces->GetGasMass()*slugtolb + ChildFDMWeight;
156
157   Mass = lbtoslug*Weight;
158
159   FGModel::PostLoad(el);
160
161   Debug(2);
162   return true;
163 }
164
165 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
166
167 bool FGMassBalance::Run(void)
168 {
169   double denom, k1, k2, k3, k4, k5, k6;
170   double Ixx, Iyy, Izz, Ixy, Ixz, Iyz;
171
172   if (FGModel::Run()) return true;
173   if (FDMExec->Holding()) return false;
174
175   RunPreFunctions();
176
177   double ChildFDMWeight = 0.0;
178   for (int fdm=0; fdm<FDMExec->GetFDMCount(); fdm++) {
179     if (FDMExec->GetChildFDM(fdm)->mated) ChildFDMWeight += FDMExec->GetChildFDM(fdm)->exec->GetMassBalance()->GetWeight();
180   }
181
182   Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetTotalPointMassWeight()
183     + BuoyantForces->GetGasMass()*slugtolb + ChildFDMWeight;
184
185   Mass = lbtoslug*Weight;
186
187 // Calculate new CG
188
189   vXYZcg = (Propulsion->GetTanksMoment() + EmptyWeight*vbaseXYZcg
190             + GetPointMassMoment()
191             + BuoyantForces->GetGasMassMoment()) / Weight;
192
193   // Track frame-by-frame delta CG, and move the EOM-tracked location
194   // by this amount.
195   if (vLastXYZcg.Magnitude() == 0.0) vLastXYZcg = vXYZcg;
196   vDeltaXYZcg = vXYZcg - vLastXYZcg;
197   vDeltaXYZcgBody = StructuralToBody(vLastXYZcg) - StructuralToBody(vXYZcg);
198   vLastXYZcg = vXYZcg;
199   Propagate->NudgeBodyLocation(vDeltaXYZcgBody);
200
201 // Calculate new total moments of inertia
202
203   // At first it is the base configuration inertia matrix ...
204   mJ = baseJ;
205   // ... with the additional term originating from the parallel axis theorem.
206   mJ += GetPointmassInertia( lbtoslug * EmptyWeight, vbaseXYZcg );
207   // Then add the contributions from the additional pointmasses.
208   mJ += CalculatePMInertias();
209   mJ += Propulsion->CalculateTankInertias();
210   mJ += BuoyantForces->GetGasMassInertia();
211
212   Ixx = mJ(1,1);
213   Iyy = mJ(2,2);
214   Izz = mJ(3,3);
215   Ixy = -mJ(1,2);
216   Ixz = -mJ(1,3);
217   Iyz = -mJ(2,3);
218
219 // Calculate inertia matrix inverse (ref. Stevens and Lewis, "Flight Control & Simulation")
220
221   k1 = (Iyy*Izz - Iyz*Iyz);
222   k2 = (Iyz*Ixz + Ixy*Izz);
223   k3 = (Ixy*Iyz + Iyy*Ixz);
224
225   denom = 1.0/(Ixx*k1 - Ixy*k2 - Ixz*k3 );
226   k1 = k1*denom;
227   k2 = k2*denom;
228   k3 = k3*denom;
229   k4 = (Izz*Ixx - Ixz*Ixz)*denom;
230   k5 = (Ixy*Ixz + Iyz*Ixx)*denom;
231   k6 = (Ixx*Iyy - Ixy*Ixy)*denom;
232
233   mJinv.InitMatrix( k1, k2, k3,
234                     k2, k4, k5,
235                     k3, k5, k6 );
236
237   RunPostFunctions();
238
239   Debug(0);
240
241   return false;
242 }
243
244 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
245
246 void FGMassBalance::AddPointMass(Element* el)
247 {
248   double radius=0, length=0;
249   Element* loc_element = el->FindElement("location");
250   string pointmass_name = el->GetAttributeValue("name");
251   if (!loc_element) {
252     cerr << "Pointmass " << pointmass_name << " has no location." << endl;
253     exit(-1);
254   }
255
256   double w = el->FindElementValueAsNumberConvertTo("weight", "LBS");
257   FGColumnVector3 vXYZ = loc_element->FindElementTripletConvertTo("IN");
258
259   PointMass *pm = new PointMass(w, vXYZ);
260   pm->SetName(pointmass_name);
261
262   Element* form_element = el->FindElement("form");
263   if (form_element) {
264     string shape = form_element->GetAttributeValue("shape");
265     Element* radius_element = form_element->FindElement("radius");
266     Element* length_element = form_element->FindElement("length");
267     if (radius_element) radius = form_element->FindElementValueAsNumberConvertTo("radius", "FT");
268     if (length_element) length = form_element->FindElementValueAsNumberConvertTo("length", "FT");
269     if (shape == "tube") {
270       pm->SetPointMassShapeType(PointMass::esTube);
271       pm->SetRadius(radius);
272       pm->SetLength(length);
273       pm->CalculateShapeInertia();
274     } else if (shape == "cylinder") {
275       pm->SetPointMassShapeType(PointMass::esCylinder);
276       pm->SetRadius(radius);
277       pm->SetLength(length);
278       pm->CalculateShapeInertia();
279     } else if (shape == "sphere") {
280       pm->SetPointMassShapeType(PointMass::esSphere);
281       pm->SetRadius(radius);
282       pm->CalculateShapeInertia();
283     } else {
284     }
285   }
286
287   pm->bind(PropertyManager, PointMasses.size());
288   PointMasses.push_back(pm);
289 }
290
291 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
292
293 double FGMassBalance::GetTotalPointMassWeight(void)
294 {
295   double PM_total_weight = 0.0;
296
297   for (unsigned int i=0; i<PointMasses.size(); i++) {
298     PM_total_weight += PointMasses[i]->Weight;
299   }
300   return PM_total_weight;
301 }
302
303 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
304
305 FGColumnVector3& FGMassBalance::GetPointMassMoment(void)
306 {
307   PointMassCG.InitMatrix();
308
309   for (unsigned int i=0; i<PointMasses.size(); i++) {
310     PointMassCG += PointMasses[i]->Weight*PointMasses[i]->Location;
311   }
312   return PointMassCG;
313 }
314
315 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
316
317 FGMatrix33& FGMassBalance::CalculatePMInertias(void)
318 {
319   unsigned int size;
320
321   size = PointMasses.size();
322   if (size == 0) return pmJ;
323
324   pmJ = FGMatrix33();
325
326   for (unsigned int i=0; i<size; i++) {
327     pmJ += GetPointmassInertia( lbtoslug * PointMasses[i]->Weight, PointMasses[i]->Location );
328     pmJ += PointMasses[i]->GetPointMassInertia();
329   }
330
331   return pmJ;
332 }
333
334 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
335
336 FGColumnVector3 FGMassBalance::StructuralToBody(const FGColumnVector3& r) const
337 {
338   // Under the assumption that in the structural frame the:
339   //
340   // - X-axis is directed afterwards,
341   // - Y-axis is directed towards the right,
342   // - Z-axis is directed upwards,
343   //
344   // (as documented in http://jsbsim.sourceforge.net/JSBSimCoordinates.pdf)
345   // we have to subtract first the center of gravity of the plane which
346   // is also defined in the structural frame:
347   //
348   //   FGColumnVector3 cgOff = r - vXYZcg;
349   //
350   // Next, we do a change of units:
351   //
352   //   cgOff *= inchtoft;
353   //
354   // And then a 180 degree rotation is done about the Y axis so that the:
355   //
356   // - X-axis is directed forward,
357   // - Y-axis is directed towards the right,
358   // - Z-axis is directed downward.
359   //
360   // This is needed because the structural and body frames are 180 degrees apart.
361
362   return FGColumnVector3(inchtoft*(vXYZcg(1)-r(1)),
363                          inchtoft*(r(2)-vXYZcg(2)),
364                          inchtoft*(vXYZcg(3)-r(3)));
365 }
366
367 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
368
369 void FGMassBalance::bind(void)
370 {
371   typedef double (FGMassBalance::*PMF)(int) const;
372   PropertyManager->Tie("inertia/mass-slugs", this,
373                        &FGMassBalance::GetMass);
374   PropertyManager->Tie("inertia/weight-lbs", this,
375                        &FGMassBalance::GetWeight);
376   PropertyManager->Tie("inertia/empty-weight-lbs", this,
377                        &FGMassBalance::GetEmptyWeight);
378   PropertyManager->Tie("inertia/cg-x-in", this,1,
379                        (PMF)&FGMassBalance::GetXYZcg);
380   PropertyManager->Tie("inertia/cg-y-in", this,2,
381                        (PMF)&FGMassBalance::GetXYZcg);
382   PropertyManager->Tie("inertia/cg-z-in", this,3,
383                        (PMF)&FGMassBalance::GetXYZcg);
384 }
385
386 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
387 //
388 // This function binds properties for each pointmass object created.
389 //
390 void FGMassBalance::PointMass::bind(FGPropertyManager* PropertyManager, int num) {
391   string tmp = CreateIndexedPropertyName("inertia/pointmass-weight-lbs", num);
392   PropertyManager->Tie( tmp.c_str(), this, &PointMass::GetPointMassWeight,
393                                        &PointMass::SetPointMassWeight);
394
395   tmp = CreateIndexedPropertyName("inertia/pointmass-location-X-inches", num);
396   PropertyManager->Tie( tmp.c_str(), this, eX, &PointMass::GetPointMassLocation,
397                                            &PointMass::SetPointMassLocation);
398   tmp = CreateIndexedPropertyName("inertia/pointmass-location-Y-inches", num);
399   PropertyManager->Tie( tmp.c_str(), this, eY, &PointMass::GetPointMassLocation,
400                                            &PointMass::SetPointMassLocation);
401   tmp = CreateIndexedPropertyName("inertia/pointmass-location-Z-inches", num);
402   PropertyManager->Tie( tmp.c_str(), this, eZ, &PointMass::GetPointMassLocation,
403                                            &PointMass::SetPointMassLocation);
404 }
405
406 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
407
408 void FGMassBalance::GetMassPropertiesReport(void) const
409 {
410   cout << endl << fgblue << highint 
411        << "  Mass Properties Report (English units: lbf, in, slug-ft^2)"
412        << reset << endl;
413   cout << "                                  " << underon << "    Weight    CG-X    CG-Y"
414        << "    CG-Z         Ixx         Iyy         Izz" << underoff << endl;
415   cout.precision(1);
416   cout << highint << setw(34) << left << "    Base Vehicle " << normint
417        << right << setw(10) << EmptyWeight << setw(8) << vbaseXYZcg(eX) << setw(8)
418        << vbaseXYZcg(eY) << setw(8) << vbaseXYZcg(eZ) << setw(12) << baseJ(1,1)
419        << setw(12) << baseJ(2,2) << setw(12) << baseJ(3,3) << endl;
420
421   for (unsigned int i=0;i<PointMasses.size();i++) {
422     PointMass* pm = PointMasses[i];
423     double pmweight = pm->GetPointMassWeight();
424     cout << highint << left << setw(4) << i << setw(30) << pm->GetName() << normint
425          << right << setw(10) << pmweight << setw(8) << pm->GetLocation()(eX)
426          << setw(8) << pm->GetLocation()(eY) << setw(8) << pm->GetLocation()(eZ)
427          << setw(12) << pm->GetPointMassMoI(1,1) << setw(12) << pm->GetPointMassMoI(2,2)
428          << setw(12) << pm->GetPointMassMoI(3,3) << endl;
429   }
430
431   for (unsigned int i=0;i<Propulsion->GetNumTanks() ;i++) {
432     FGTank* tank = Propulsion->GetTank(i);
433     string tankname="";
434     if (tank->GetType() == FGTank::ttFUEL && tank->GetGrainType() != FGTank::gtUNKNOWN) {
435       tankname = "Solid Fuel";
436     } else if (tank->GetType() == FGTank::ttFUEL) {
437       tankname = "Fuel";
438     } else if (tank->GetType() == FGTank::ttOXIDIZER) {
439       tankname = "Oxidizer";
440     } else {
441       tankname = "(Unknown tank type)";
442     }
443     cout << highint << left << setw(4) << i << setw(30) << tankname << normint
444       << right << setw(10) << tank->GetContents() << setw(8) << tank->GetXYZ(eX)
445          << setw(8) << tank->GetXYZ(eY) << setw(8) << tank->GetXYZ(eZ)
446          << setw(12) << "*" << setw(12) << "*"
447          << setw(12) << "*" << endl;
448   }
449
450   cout << underon << setw(104) << " " << underoff << endl;
451   cout << highint << left << setw(30) << "    Total: " << right << setw(14) << Weight 
452        << setw(8) << vXYZcg(eX)
453        << setw(8) << vXYZcg(eY)
454        << setw(8) << vXYZcg(eZ)
455        << setw(12) << mJ(1,1)
456        << setw(12) << mJ(2,2)
457        << setw(12) << mJ(3,3)
458        << normint << endl;
459
460   cout.setf(ios_base::fixed);
461 }
462
463 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
464 //    The bitmasked value choices are as follows:
465 //    unset: In this case (the default) JSBSim would only print
466 //       out the normally expected messages, essentially echoing
467 //       the config files as they are read. If the environment
468 //       variable is not set, debug_lvl is set to 1 internally
469 //    0: This requests JSBSim not to output any messages
470 //       whatsoever.
471 //    1: This value explicity requests the normal JSBSim
472 //       startup messages
473 //    2: This value asks for a message to be printed out when
474 //       a class is instantiated
475 //    4: When this value is set, a message is displayed when a
476 //       FGModel object executes its Run() method
477 //    8: When this value is set, various runtime state variables
478 //       are printed out periodically
479 //    16: When set various parameters are sanity checked and
480 //       a message is printed out when they go out of bounds
481
482 void FGMassBalance::Debug(int from)
483 {
484   if (debug_lvl <= 0) return;
485
486   if (debug_lvl & 1) { // Standard console startup message output
487     if (from == 2) { // Loading
488       cout << endl << "  Mass and Balance:" << endl;
489       cout << "    baseIxx: " << baseJ(1,1) << " slug-ft2" << endl;
490       cout << "    baseIyy: " << baseJ(2,2) << " slug-ft2" << endl;
491       cout << "    baseIzz: " << baseJ(3,3) << " slug-ft2" << endl;
492       cout << "    baseIxy: " << baseJ(1,2) << " slug-ft2" << endl;
493       cout << "    baseIxz: " << baseJ(1,3) << " slug-ft2" << endl;
494       cout << "    baseIyz: " << baseJ(2,3) << " slug-ft2" << endl;
495       cout << "    Empty Weight: " << EmptyWeight << " lbm" << endl;
496       cout << "    CG (x, y, z): " << vbaseXYZcg << endl;
497       // ToDo: Need to add point mass outputs here
498       for (unsigned int i=0; i<PointMasses.size(); i++) {
499         cout << "    Point Mass Object: " << PointMasses[i]->Weight << " lbs. at "
500                    << "X, Y, Z (in.): " << PointMasses[i]->Location(eX) << "  "
501                    << PointMasses[i]->Location(eY) << "  "
502                    << PointMasses[i]->Location(eZ) << endl;
503       }
504     }
505   }
506   if (debug_lvl & 2 ) { // Instantiation/Destruction notification
507     if (from == 0) cout << "Instantiated: FGMassBalance" << endl;
508     if (from == 1) cout << "Destroyed:    FGMassBalance" << endl;
509   }
510   if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
511   }
512   if (debug_lvl & 8 ) { // Runtime state variables
513   }
514   if (debug_lvl & 16) { // Sanity checking
515     if (from == 2) {
516       if (EmptyWeight <= 0.0 || EmptyWeight > 1e9)
517         cout << "MassBalance::EmptyWeight out of bounds: " << EmptyWeight << endl;
518       if (Weight <= 0.0 || Weight > 1e9)
519         cout << "MassBalance::Weight out of bounds: " << Weight << endl;
520       if (Mass <= 0.0 || Mass > 1e9)
521         cout << "MassBalance::Mass out of bounds: " << Mass << endl;
522     }
523   }
524   if (debug_lvl & 64) {
525     if (from == 0) { // Constructor
526       cout << IdSrc << endl;
527       cout << IdHdr << endl;
528     }
529   }
530 }
531 }