- return vWindNED.Magnitude();
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-void FGAtmosphere::SetWindPsi(double dir)
-{
- double mag = GetWindspeed();
- psiw = dir;
- SetWindspeed(mag);
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-void FGAtmosphere::Turbulence(void)
-{
- switch (turbType) {
- case ttStandard: {
- TurbGain = TurbGain * TurbGain * 100.0;
-
- vDirectiondAccelDt(eX) = 1 - 2.0*(double(rand())/double(RAND_MAX));
- vDirectiondAccelDt(eY) = 1 - 2.0*(double(rand())/double(RAND_MAX));
- vDirectiondAccelDt(eZ) = 1 - 2.0*(double(rand())/double(RAND_MAX));
-
- MagnitudedAccelDt = 1 - 2.0*(double(rand())/double(RAND_MAX)) - Magnitude;
- // Scale the magnitude so that it moves
- // away from the peaks
- MagnitudedAccelDt = ((MagnitudedAccelDt - Magnitude) /
- (1 + fabs(Magnitude)));
- MagnitudeAccel += MagnitudedAccelDt*rate*TurbRate*State->Getdt();
- Magnitude += MagnitudeAccel*rate*State->Getdt();
- Magnitude = fabs(Magnitude);
-
- vDirectiondAccelDt.Normalize();
-
- // deemphasise non-vertical forces
- vDirectiondAccelDt(eX) = square_signed(vDirectiondAccelDt(eX));
- vDirectiondAccelDt(eY) = square_signed(vDirectiondAccelDt(eY));
-
- vDirectionAccel += vDirectiondAccelDt*rate*TurbRate*State->Getdt();
- vDirectionAccel.Normalize();
- vDirection += vDirectionAccel*rate*State->Getdt();
-
- vDirection.Normalize();
-
- // Diminish turbulence within three wingspans
- // of the ground
- vTurbulenceNED = TurbGain * Magnitude * vDirection;
- double HOverBMAC = Auxiliary->GetHOverBMAC();
- if (HOverBMAC < 3.0)
- vTurbulenceNED *= (HOverBMAC / 3.0) * (HOverBMAC / 3.0);
-
- // I don't believe these next two statements calculate the proper gradient over
- // the aircraft body. One reason is because this has no relationship with the
- // orientation or velocity of the aircraft, which it must have. What is vTurbulenceGrad
- // supposed to represent? And the direction and magnitude of the turbulence can change,
- // so both accelerations need to be accounted for, no?
-
- // Need to determine the turbulence change in body axes between two time points.
-
- vTurbulenceGrad = TurbGain*MagnitudeAccel * vDirection;
- vBodyTurbGrad = Propagate->GetTl2b()*vTurbulenceGrad;
-
- if (Aircraft->GetWingSpan() > 0) {
- vTurbPQR(eP) = vBodyTurbGrad(eY)/Aircraft->GetWingSpan();
- } else {
- vTurbPQR(eP) = vBodyTurbGrad(eY)/30.0;
- }
-// if (Aircraft->GetHTailArm() != 0.0)
-// vTurbPQR(eQ) = vBodyTurbGrad(eZ)/Aircraft->GetHTailArm();
-// else
-// vTurbPQR(eQ) = vBodyTurbGrad(eZ)/10.0;
-
- if (Aircraft->GetVTailArm() > 0)
- vTurbPQR(eR) = vBodyTurbGrad(eX)/Aircraft->GetVTailArm();
- else
- vTurbPQR(eR) = vBodyTurbGrad(eX)/10.0;
-
- // Clear the horizontal forces
- // actually felt by the plane, now
- // that we've used them to calculate
- // moments.
- // Why? (JSB)
-// vTurbulenceNED(eX) = 0.0;
-// vTurbulenceNED(eY) = 0.0;