Accurate Trajectory Prediction for Typical Artillery Projectile

0 downloads 0 Views 266KB Size Report
Jun 2, 2018 - dynamics give its angular velocity as a function of the corresponding ... 0 sin sec cos sec t t t p q r. (3). Where, pt, qt, and rt are the angular acceleration taking in ... Where, sin and cos are abbreviated with s and c, respectively.
Proceedings of the 33rd Chinese Control Conference July 28-30, 2014, Nanjing, China

Accurate Trajectory Prediction for Typical Artillery Projectile A. Elsaadany1, YI Wen-jun2 1. Nanjing University of Science and Technology, Nanjing 210094, China E-mail: [email protected] 2. Nanjing University of Science and Technology, Nanjing 210094, China E-mail: [email protected] Abstract: A full six degree-of-freedom nonlinear model is proposed for the accurate prediction at short and long range trajectories of high and low spin stabilized projectiles via atmospheric flight. The model includes Earth’s rotation and ellipsoidal shape, Magnus effect and the atmospheric wind. Furthermore, a modified standard atmospheric model to simulate air density and the speed of sound is used. The computation is performed in time marching scheme. The fourth-order Runge Kutta method is employed for integration. Aerodynamic forces and moments are calculated using aerodynamic coefficients of the projectile that predicted using PRODAS program. The computational flight analysis takes into consideration all the aerodynamics variations by means of the variable aerodynamic coefficients. The computational results of the proposed synthesized analysis give satisfactory agreement. The Static stability, also called gyroscopic stability, for the developed projectile model, is examined to ensure the projectile would remain gyroscopically stable in flight. In addition, the correction capability of range using a drag ring brake module is presented. The simulation results show that the impact accuracy of a conventional projectile using this drag brake module can be improved. The drag ring brake is found to be highly capable for range correction. The deploying of drag brake in early stage of trajectory results in large range correction. The correction occasion time can be predefined depending on required correction of range. Key Words: 6-DOF Model, Trajectory prediction, Aerodynamic coefficients, Gyroscopic stability, Range correction.



1

generally also require the most input data and as such have the most parameter and initial state error.

Introduction

The field of guided artillery projectiles has attracted increased interest over the last few years due mostly to increased requirements concerning impact accuracy in order to decrease dispersion as well as. For the purpose of regulating a projectile trajectory through the synthesis of appropriate control laws, the availability of a sufficiently accurate nonlinear projectile model is of crucial importance. Therefore, the accurate model can play a significant role in the efficient tracking control system design for current and future control actions applied to fin and spin-stabilized projectiles. A different trajectory models are available with varying level of fidelity and data requirements. These models split into three categories: point-mass (PM), modified point-mass (MPM) and rigid body which have been developed for computing trajectory elements and each of which has varying degree of complexity. As the name suggests, PM model assumes that the projectile is a point mass that poetesses three position degrees of freedom. The MPM model assumes that the projectile is a point mass, but the roll dynamic equation is appended to the model along with side force to better replicate cross range of spin stabilized projectiles [1]. The rigid body model assumes that the projectile is a rigid body that possesses three position and three orientation degrees of freedom [2-6]. While rigid body models are considered the most accurate with the least model error, they

6368

In this study, a predicted trajectory calculation using an accurate Six Degree-of-Freedom (6-DOF) nonlinear model is employed. The proposed nonlinear model is developed, based on Newton’s equations of motion, taking into consideration the influence of the Earth’s rotation and ellipsoidal shape, Magnus effect and the atmospheric wind. Furthermore, a modified standard atmospheric model to simulate air density and the speed of sound is used. The computation is performed in time marching scheme. The fourth-order Runge Kutta method was employed for integration. Aerodynamic forces and moments are calculated using aerodynamic coefficients that predicted using PRODAS program, which a well known aerodynamics prediction code. The computational flight analysis takes into consideration the Mach number and total angle of attack effects by means of the variable aerodynamic coefficients. The model is applicable to predict the velocity magnitude, position, time of flight and angle of inclination at any point at the trajectory of the projectile. The maximum height, maximum range and cross range of the projectile trajectory can also be calculated with this simulation. The computational results of the proposed synthesized analysis give satisfactory agreement. The developed simulated model features a very quick and very cheap prediction opportunity for the trajectory parameter calculations of the unguided projectile. Also it is applicable for different projectiles calibers and different flight regimes. The Static stability, also called gyroscopic stability, is examined for the developed projectile model. In addition,

the correction capability of range using a drag ring brake module is presented. The simulation results show that the impact accuracy of a conventional projectile using this drag brake module can be improved. The drag ring brake is found to be highly capable for range correction. The deploying of drag brake in early stage of trajectory results in large range correction. The correction occasion time can be predefined depending on required correction of range.

2

ªc T c\ « c T s\ « «¬  s T

J BE

O VN P VE

The motion of a projectile in atmospheric flight is assumed to be adequately described with six rigid body degrees of freedom comprised of three inertial position body coordinates as well as three Euler angle body attitudes. This section presents the equation of motion that model the projectile’s atmospheric trajectory according to a set of initial launch conditions. As well as the methodology used to compute forces and moments acting on the projectile body [2, 7]. The projectile is assumed to be both rigid (non-flexible), and rotationally symmetric about its spin axis. The translation dynamics give the projectile linear velocity as a result of the externally applied forces, whereas the rotation dynamics give its angular velocity as a function of the corresponding moments. The translational dynamic equations are obtained using a force balance equation on the projectile written in the body-fixed frame, given by

­axb ½ ° b° ®a y ¾ °a b ° ¯ z¿

­ FX ½ 1° ° ® FY ¾ m° ° ¯ FZ ¿

(1)

ª ­ L ½ ª 0 r 1 « ° > I @ « ® M ¾°  «« r 0 «¬ °¯ N °¿ «¬  q p

q º ­ p ½º ° °»  p »» > I @ ® q ¾» 0 ¼» ¯° r ¿°»¼

(2)

Note that, the terms L, M, and N are the sum of steady aerodynamic, unsteady aerodynamic and Magnus moments resolved in the body-fixed reference frame. The rotational kinematic equations relate derivatives of Euler angles to angular velocity states are given by

­M ½ ° ° ®T ¾ °\ ° ¯ ¿

ª1 sin M tan T «0 cos M « ¬«0 sin M sec T

cos M tan T º ­ pt ½ ° °  sin M »» ® qt ¾ cos M sec T ¼» ¯° rt ¿°

(3)

Where, pt, qt, and rt are the angular acceleration taking in account the error resulting from the Earth’s rotation that expressed in terms of the vehicle-carried North East Down (NED) velocity and given by

­M ½ ° ° ®T ¾ °\ ° ¯ ¿

ª Z e  P cos O º ­ p½ « » ° ° O » ® q ¾  J BE « « » °r ° e ¯ ¿ «¬  Z  P sin O ¼»

(4)

Rmer  h Rnorm  h cos O

(6) (7)

P  2Z sin O V  P  2Z cos O V V O V  P  2Z cos O V  a

VE

e

e

N

D

 aE

e

D

N

E

D

(10) (11)

The translational kinematic equations, relating vehicle-carried NED acceleration states to body-fixed acceleration states, are given by

­ aN ½ ° ° ® aE ¾ °a ° ¯ D¿

­axb ½ 1 ° b ° J BE ®a y ¾ °a b ° ¯ z¿

(12)

1 J BE is the inverse of the transformation matrix

which rotates from the body frame to the vehicle-carried NED frame. The parameters of equations from 6 to 11 are defined and derived based on the WGS 84 (world geodetic system 84, which was originally proposed in 1984 and lastly updated in 2004) [8]. 2.1

Force and moment model for projectile body

During flight there are two kinds of forces acting on Projectile motion. They are weight and resultant body aerodynamic forces. When the aerodynamic forces do not pass through the center of gravity, moments are created. The aerodynamic forces and moments are calculated within the flight simulation using the aerodynamic coefficients that predicted using PRODAS program. The total forces acting on the projectile can be expressed as

­ FX ½ ° ° ® FY ¾ °F ° ¯ Z¿

­ FWX ½ ­ FAX ½ ­ FMX ½ ° ° ° ° ° ° ® FWY ¾  ® FAY ¾  ® FMY ¾ °F ° °F ° °F ° ¯ WZ ¿ ¯ AZ ¿ ¯ MZ ¿

(13)

­ FWX ½ ª  sin T º ° ° « sin M cos T » F mg ® WY ¾ « » °F ° ¯ WZ ¿ ¬«cos M cos T ¼»

(14)

­ FAX ½ ª CD º ° ° « » ® FAY ¾ q S «CND E » °F ° «¬CND D »¼ ¯ AZ ¿

(15)

­ FMX ½ ° ° ® FMY ¾ °F ° ¯ MZ ¿

6369

(5)

(8) h VD The derivatives of the vehicle-carried NED velocity components taking in account Earth’s rotation are respectively given by VN  P  2Z e sin O VE  O VD  a N (9)

Where,

Note that, the terms FX, FY, and FZ are the sum of weight, Magnus, and aerodynamic forces resolved in the body-fixed reference frame. The rotational dynamic equations are obtained using a moment equation about the projectile mass center and written in the body-fixed frame, given by

­ p ½ ° ° ® q ¾ ° r ° ¯ ¿

s M s T s\  c M c\ s M cT

c M s T c\  s M s\ º c M s T s\  s M c\ »» »¼ cM cT

Where, sin and cos are abbreviated with s and c, respectively. The derivative of the geodetic position can be expressed in terms of the vehicle-carried NED velocity and given by

Projectile Dynamic Model

­ u ½ ° ° ® v ¾ ° w ° ¯ ¿

s M s T c\  c M s\

ª 0 º § PD · « » qS¨ ¸ « CYpD D » © Vt ¹ « C E » ¬ YpD ¼

(16)

Where, Vt is the body velocity magnitude and given by

T

(17)

The lateral and longitudinal aerodynamic angles of attack, respectively, are computed as

E sin 1 v Vt , D

tan 1 w u

(18)

These forces are expressed in the body-fixed frame and split into contributions due to weight (W), body aerodynamic forces, respectively. The body aerodynamic force is split into a standard aerodynamic (A) and Magnus (M) forces, respectively. The total moments acting on the projectile can be expressed by

­L½ ° ° ®M ¾ °N ° ¯ ¿ ­ LSA ½ ° ° ® M SA ¾ °N ° ¯ SA ¿

­ LSA ½ ­ LUA ½ ­ LM ½ ° ° ° ° ° ° ® M SA ¾  ® M UA ¾  ® M M ¾ °N ° °N ° °N ° ¯ SA ¿ ¯ UA ¿ ¯ M ¿ ª 0 «R « CPZ ¬«  RCPY

 RCPZ 0 RCPX

(19)

RCPY º ­ FAX ½ ° °  RCPX »» ® FAY ¾ 0 ¼» ¯° FAZ ¿°

­ LUA ½ ° ° ® M UA ¾ °N ° ¯ UA ¿

ª P D Vt Clp º « » q S D « q D Vt Cmq  CmD D » «¬ r D Vt Cnr  CmD E »¼ ª 0 º ­ LM ½ § PD · « ° ° » ®M M ¾ q S D ¨ ¸ «CnpD E » V © t ¹ «C D » °N ° ¯ M¿ ¬ npD ¼

§ g

a

· 1¸ ¹

e

g Z

trop  Z



RT



J RT

3 2.5 2 3

Density [Kg/m ] Speed of sound/100 [m/sec] temperature/100 [K] pressure/100 [KPa]

1.5 1

0

(21)

0

0.5

1 Altitude [m]

1.5

2 x 10

4

Fig. 1: Atmospheric parameters versus altitude

4 (22)

Variations in meteorological conditions have an effect on the projectile traveling through the atmosphere and hence affect its trajectory. The artillery projectile typically has peak altitudes of about 20 kilometers which is within the troposphere and is thus subjected to air density and drag. With increasing altitude, air properties such as density, temperature, pressure and air viscosity change. Therefore, this changing is taken into the consideration during the trajectory calculation to get an accurate prediction. For this purpose, a standard atmosphere model is developed based on the International Standard Atmosphere (ISA) [9]. Expressions for air density, ȡ, and speed of sound, a, as a function of altitude, Z, can be derived and are given by

U 0 T T0 ¨© R L

3.5

0.5

Atmospheric model

U

(25)

Where, ȡ0 and T0 are the air density and air temperature at the sea level, respectively. g, R, L, Ztrop, and J are gravity acceleration, real gas constant for air, temperature lapse rate, tropopause altitude, and adiabatic gas constant, respectively. Fig. 1 depicts the variation of the air density, sound speed, air temperature, and air pressure as a function of altitude.

(20)

These moments contain steady aerodynamic (SA), unsteady aerodynamic (UA) and Magnus (M), respectively. The steady aerodynamic moments is computed with a cross product between the distance vector from the center of gravity to the body center of pressure and the body aerodynamic force vector in equation 15. The aerodynamic coefficients and aerodynamic center distance are all a function of local Mach number at the mass center of the projectile. Computationally, these Mach number dependent parameters are obtained by a table lookup scheme using linear interpolation.

3

T0  L Z

(23) (24)

Where, T is the air temperature and is given by

6370

Wind model

The Earth's atmosphere generally fosters air velocity perturbations that can significantly modify the trajectory of a projectile. It is common to separate the atmospheric air velocity perturbations into steady and fluctuating components, typically called mean wind and atmospheric turbulence, respectively. The atmospheric wind model adopted in the current work includes only a steady wind component. Fig. 2 shows a recorded data for wind speed (in m/s) and direction (in degrees) as a function of altitude that used in this work. Wind speed and direction 30

280

20

260

10

240

Wind direction [deg]

2

Atmospheric parameters

u v w 2

Wind speed [m/sec]

Vt

2

Wind speed Wind direction 0

0

2

4

6

8

220 10

Altitude [Km]

Fig. 2: Wind speed and direction versus altitude

5

Projectile Configuration

The developed nonlinear 6-DOF model in this work uses the nominal physical properties of a typical 155mm spin- stabilized projectile. These physical properties are listed in Table 1. This is followed with a range correction capability. The range correction is obtained using a drag ring brake module that fitted

onto the projectile within a fuze. This fuze screws into the forward portion of the projectile (see Fig. 4). 6

Simulation conditions

Trajectory simulations have been performed in the following flight conditions: Muzzle velocity, V0 = 910 m/sec and the initial elevation angle, ș0 = 450. 6.1

Muzzle spin rate estimation

According to McCoy definition [2], the muzzle spin rate can be estimated by

p0

2 S V0 K D

(26)

Where, Ș is the rifling twist rate at the gun muzzle (caliber per turn).

7

Aerodynamic coefficients

The body non-dimensional aerodynamic coefficients and distance that predicted using PRODAS program, based on the dimensions of the nominal 155mm projectile, are shown in Fig. 3. The total drag coefficient, in case of range correction (body + opened drag ring), is also shown in Fig. 3.

0.2

0

1

2

1.5

3

0

1

2

2

CDt [nd] 0

1

2

3

0

1

2

1

2

3

0

1

2

3

0

1

2

3

-8 -10 -12

3

0.35

0.5

0

0

-6

1

-0.2

-1

-1.5

3

4

3.5

3

0

CnpD [nd]

1

Cmq [nd]

-0.02

-0.4

0

4.5

CmD [nd]

Clp [nd]

0

-0.04

2

RCPx [m]

0

-0.5

CNpD [nd]

2.5

CND [nd]

CD [nd]

0.4

mind that the forces are attempting to raise the projectile’s axis of rotation. Gyroscopic stability ensures that a projectile will not immediately tumble upon leaving the muzzle of the rifle [2]. In Fig. 5, two cases of static stability are demonstrated: in the top figure, CP lies behind the CG so that a clockwise (restoring) moment is produced. This case tends to reduce the yaw angle and return the body to its trajectory, therefore statically stable. Conversely, the lower figure, with CP ahead of CG, produces an anti-clockwise (overturning) moment, which increases a further and is therefore statically unstable. It also possible to have a neutral case in which CP and CG are coincident whereby no moment is produced. There is clearly an important correspondence in the distance between the center of pressure and the center of gravity and the center of the round. This distance is called the static margin. By definition, it is positive for positive static stability, zero for neutral stability and negative for negative stability. The projectile is initially flying at zero angle of incidence along its flight trajectory and is then struck by the wind gust (see Fig. 6). So that the nose is deflected upwards, producing an incidence angle (Į). The response of the projectile to this disturbed incidence angle determines its stability characteristics. In particular, the initial response, determines whether it is statically stable or unstable.

0

2

4

0.3

0.25

Fig. 5: Static stability/instability conditions

Mach [nd]

Fig. 6: Gust producing incidence angle

Fig. 3: Body aerodynamic parameters versus Mach number

If the initial response is to move the nose back towards zero incidence (i.e. reducing the incidence angle) then it is statically stable. If the incidence angle initially increases as a response then it is statically unstable while if the disturbed incidence angle is retained the projectile is neutral. It is therefore clear that it is the direction or sign of the resulting pitching moment generated that defines the static stability. This depends upon the aerodynamic force produced on the body due to the incidence angle and, in particular, upon the normal force (i.e. the aerodynamic force component perpendicular to the body axis) and the position at which it acts along the body’s axis. McCoy [2] defines the gyroscopic stability factor, S g , in the following generalized form:

Fig. 4: Basic geometrical data of the proposed projectile Table 1: Projectile physical properties Parameters Caliber Length Total Mass Center of Gravity from the Nose Axial Moment of Inertia Lateral Moment of Inertia

8

Value D = 155 mm l = 902 mm m = 45 kg XCG = 558 mm IXX = 0.162 Kg. m2 IYY = IZZ =1.763 Kg. m2

Sg

2 I XX P2 2 U IYY S DVt 2 CmD

(27)

Where, the condition for stable projectile flight in the initial section of the trajectory is as follow

Gyroscopic stability

Any spinning object will have gyroscopic properties. In spin stabilized projectile, the center of pressure is located in front of the center of gravity. Hence, as the projectile leaves the muzzle it experiences an overturning movement caused by air force acting about the center of mass. It must be kept in

6371

Sg ! 1

9

(28)

Simulation results

For the purpose of verifying the developed nonlinear 6-DOF model of the given projectile, Monte Carlo

12

10

60 40 20

Pitch angle [deg]

simulations are obtained. Fig. 7 shows the projectile altitude versus range. It is observed that, the maximum range is slightly over 30 km, and the maximum altitude, at the summit time, is nearly 10 km, where the summit time is approximately 41 sec. Fig. 8 shows the axial acceleration acting on the projectile as a function of flight time. At the firing point, the axial deceleration is nearly 4.73 g. And due to, the aerodynamic axial force acting on the projectile and the gravitational acceleration component in the opposite direction of its flight, this deceleration will be decreased. Fig. 9 shows the velocity magnitude of the projectile versus flight time from firing to impact point, where the velocity at the firing point is 910 m/sec and it will be decreased as well as the projectile goes up, until the summit point, then the projectile will go down to increase the projectile velocity due to a gravitational acceleration component in the direction of the velocity vector, and reached to 363 m/sec at the impact point. Fig. 10 shows the variation of the elevation angle of the projectile versus flight time, where starting from 450 and decreases until reached to 00 at the summit point, then decreases again until reached to -640 at the impact point.

0 -20 -40 -60 -80 0

10

20

30

40

50

60

70

80

90

100

Flight time [s]

Fig. 10: Pitch angle vs. flight time

Fig. 11 shows the spin rate of the given projectile as a function of flight time, where the spin rate at the firing point is 1845 rad/sec and then decreases all flight due to the friction that acting on the projectile body until reached to 1377 rad/sec at the impact point. After the damping of the initial transient motion, the stability factor for 450 has increased from 2.18 at muzzle to 40 at the summit point and then decreased to value of 7.35 at final impact point as shown in Fig. 12. The computational result presented in Fig. 12 shows the projectile has been made gyroscopically stable by appropriate selection of the spin rate ( S g ! 1 ).

Altitude [km]

8

1850 6

1750

Spin rate [rad/s]

4

2

0 0

5

10

15

20

25

30

35

Range [km]

Fig. 7: Altitude vs. range

1550

1450

1

0

Axial acceleration [gees]

1650

1350 0

10

20

30

40

50

60

70

80

90

100

Flight time [s] -1

Fig. 11: Spin rate vs. flight time 45

-2

40 35

-3

30

-5

Sg [nd]

-4

0

10

20

30

40

50

60

70

80

90

25 20

100

Flight time [s]

15

Fig. 8: Axial acceleration vs. flight time

10

1000 5

Flight velocity [m/s]

900

0

0

0.5

1

1.5

2

2.5

Range [m]

800

3

3.5 4

x 10

Fig. 12: Gyroscopic stability vs. range 700

Figs.13 and 14 depict comparative predicted trajectories of the 155mm projectile at pitch angle of 450, for both the nominal (no-wind) and wind-based cases. It is observed that the drift in the nominal trajectory case is approximately 806 m to the right at the impact point, while the drift in the wind-based trajectory case is nearly 6 m to the right (i.e., the drift is reduced by nearly 800 m).

600 500 400 300 0

10

20

30

40

50

60

70

80

90

100

Flight time [s]

Fig. 9: Flight velocity vs. flight time

6372

closer to the intended target. The results show that the drag ring brake is found to be highly capable for range correction. The deploying of the drag brake in early stage of trajectory results in maximum range correction (see Fig. 16). Hence, maximum range correction is observed at the earliest occasion time, where at occasion time of 20 seconds gives correction of 1609 m (i.e., 5.34% of the total travelled range). The differences in velocities due to drag brake are shown in Fig. 17. It is observed that deployment of drag brake throughout the flight reduces the remaining velocity, by nearly 44 m/sec., which gives the required correction in range. From the results listed in Table 2, it can conclude that the occasion time can be predefined depending on required correction of range.

Altitude [km]

15

without wind with wind

10

5

0 -500 40

0

30 20

500 10 1000

Drift [km]

0

Range [km]

Fig. 13: Flight trajectory path 900

without wind with wind

800

Table 2: Correction vs. occasion time

700

Ctrl start time [sec]

Drift [m]

600 500

Correction [m]

t = 20 t = 30 t = 40 t = 50

400 300 200

1609 1084 780 480

100 12000

0 -100 0

10

20

30

40

50

60

70

80

90

No correction At time=50 s At time=40 s At time=30 s At time=20 s

10000

100

Time [sec] 8000

Altitude [m]

Fig. 14: Drift vs. flight time

The Unguided projectiles exhibit large miss distances on impact, even at relatively short ranges. As shown in Fig. 15, the projectile is fired at intended target B, but due to wind and other meteorological conditions, muzzle velocity error, aiming error, etc., the projectile actually impact point A. For this reason, a several concepts and modifications, for trajectory correction, have been proposed. In this work, a range correction is obtained using a drag ring brake module, shown in Fig. 4. The concept of drag ring brake has discussed by Hollis [10], which replace the standard fuze of conventional artillery projectile with trajectory corrector module, without any changes within the ogive shape of the artillery projectile. Drag brake is designed to fit onto a spin-stabilized projectile within the fuze which screws into the forward portion of the projectile (see Fig. 4).

6000

4000

2000

0 1.6

1.8

2

2.2

2.4

2.6

2.8

3

3.2 4

Range [m]

x 10

Fig. 16: Altitude vs. range 480

No correction At time=50 s At time=40 s At time=30 s At time=20 s

460 440

Flight velocity [m/sec]

10 Range correction

420 400 380 360 340 320 300 280 20

30

40

50

60

70

80

90

100

Time [sec]

Fig. 17: Flight velocity vs. flight time

11 Conclusions

Fig. 15: The intended and actual paths of the projectile

During the course correction phase, semi-circular plates will deploy from the module. The plates create a blunt cross-sectional area in front of the projectile, thus creating more drag and effectively slowing the projectile. The analysis is carried out for deploying the drag brake at various stages of trajectory and once the drag ring module is deployed, it will remain open throughout the trajectory, and then the projectile will slow down, ultimately bringing it

6373

A complicated 6-DOF dynamic nonlinear model is developed for the accurate prediction at short and long-range trajectories of high and low spin stabilized projectiles. The computational results of the proposed synthesized analysis give satisfactory agreement. Criteria and analysis of gyroscopic stability, for the developed projectile model, is examined to ensure the projectile would remain gyroscopically stable in flight. In addition, the correction in range for different occasion times, using a drag ring brake module, is predicted. The results show that the drag ring brake is found to be highly capable for range correction. The

deploying of the drag brake in early stage of trajectory results in maximum range correction. Hence, maximum range correction is observed at the earliest occasion time, where at occasion time of 20 seconds gives correction of 1609 m (i.e., 5.34% of the total travelled range). The occasion time can be predefined depending on required correction of range.

References [1]

Robert F. Lieske, Use of the Magnus Force in the Modified Point Mass Trajectory Model, Ballistic Research Laboratory, Aberdeen Proving Ground, Maryland, October 1990. [2] McCoy, R., Modern Exterior Ballistics, Schiffer, Attlen, PA, 1999. [3] G. R. Cooper and Mark Costello, Trajectory Prediction of Spin-Stabilized Projectiles with a Liquid Payload, Journal of Spacecraft and Rockets, Vol. 48, No. 4, July–August 2011. [4] M. Khalil, H. Abdalla, and O. Kamal, Dispersion Analysis for Spinning Artillery Projectile, 13th International Conference on Aerospace Science and Aviation Technology, MTC, May 26-28 2009. [5] Joseph K., Costello M., and Jubaraj S., Generating an Aerodynamic Model for Projectile Flight Simulation Using Unsteady Time Accurate Computational Fluid Dynamic Results, Army Research Laboratory, ARL-CR-577, 2006. [6] Amoruso, Michael J., Euler Angles and Quaternions in Six Degree of Freedom Simulations of Projectiles, US Army Armament, Munitions and Chemicals Command, Picatinny Arsenal, New Jersey, March, 1996. [7] Etkin, B., Dynamics of Atmospheric Flight, John Wiley and Sons, New York, 1972. [8] NIMA TR8350.2: Department of Defense World Geodetic System 1984, Its Definition and Relationship with Local Geodetic Systems. [9] Gyatt, Graham (2006-01-14): The Standard Atmosphere. A mathematical model of the 1976 U.S. Standard Atmosphere. [10] Michael S. L. Hollis, Fred J. Brandon, Range Correction Module for a Spin Stabilized Projectile, U.S. Patent Documents, US-5816531, Oct. 6, 1998.

6374