BASIC programs for battleship fire control and antiaircraft gun fire control (also a space battleship fire control program) By Ian W. Douglas, July 2011. Appendix (see end of article) added January 2012 and revised March 2012. Further revised December 2011: these further revisions involve the programs TICTAC1.BAS, BATCON1.BAS, ATICTAC.BAS, SPACBAT.BAS, and AIRBAT.BAS . I rewrote sections of these programs to enable the programs to better handle situations where the target ship or aircraft are following perfect 0 degree, 90 degree, 180 degree and 270 degree course azimuths, as measured by the (own) ship or antiaircraft gun tracking instrument and computer. Further revised March 2012: standardized the flow and logic between IOELLIP.BAS and IOELOID.BAS to make the flow and logic of these two programs as similar as possible, wherever possible. Fixed some bugs in these two programs (division by zero was sometimes possible in the earlier versions of these programs, in cases where a target or point was on an axis or had zero as one of the coordinates). Added programs GUNEOD1.BAS and EIDGUN2.BAS, together with explanatory notes for these programs; these programs determine whether a stationary airborne target is within range of an antiaircraft gun, assuming that the only forces acting on the shell are the gun's initial impulse and gravity. Further revised July 2012: clarified function (via remarks in program listing) of some of the parts of certain programs; parts of some of the programs convert a vector and starting point to an end point. I did not explain this clearly in earlier versions of programs. Further revised December 2012: added more discussion of the problem of whether the target is within range of an antiaircraft gun; added program UPARTR.BAS at end of Appendix; made some comments in program listings more consistent internally; and alerted the reader to weaknesses in UPART.BAS, AIRBAT.BAS, and UPARTR.BAS. I have written some programs in BASIC which deal with the problem of tracking an enemy ship from a moving ship, and for battleship gun fire control. The listings for the programs are included below on this page and file. I also wrote some programs (included on this page as well) which solve the same problem, except with aircraft tracking and antiaircraft gun fire control. Some of the REM (remark or comment) lines in the program listings refer to flowcharts I made in order to write the programs, and I would be willing to send text descriptions of the flowcharts to those who are interested (my flowcharts are hand drawn and hand written and take up several pages, so they would be difficult to send graphically). Readers who are familiar with flowcharting and programming can reconstruct my flowcharts from my program listings, and those who wish to can send me their guesses to me as to what my flowcharts look like. I should be able to send back my comments as to what is right or wrong about your guesses about my flowcharts. I have sprinkled many comments throughout the program listings so as to make understanding of the programs' functions and logic relatively easy to understand, assuming that the reader has some knowledge of BASIC programming. I have test run the programs successfully using the Microsoft GW-BASIC interpreter. The programs use a cross between the Dreyer plot (relative plot) and the Pollen plot (true plot). I used a modified true plot which assumes that our own ship's heading is always 0 degrees/360 degrees; however, given this, our own ship moves against a fixed background in a plot of the battle. If our own ship does not turn, then this is a good "Occam's razor" approach to the problem in which we assume a 2 dimensional Cartesian coordinate system for the two ships. The origin of the coordinate system I used puts our own ship at (50 nautical miles, 50 nautical miles) at the start of the problem so as to avoid using negative numbers to represent points in space occupied by the ships and projectile. The projectile does move in three dimensions; however, since gun and target are in 2 dimensions, the problem of the projectile's motion is--for intents and purposes inherent in ship-to-ship firing--more or less two dimensional. My program uses, and takes into account, the arc of trajectory. There are two programs I wrote which address the battleship fire control problem. The first I called "TICTAC1.BAS" because of the acronym TICTAC I used, which stands for "time interval compensating tachymetry." This program uses data from a range and bearing finder (assuming a flat sea, a ship which does not pitch, yaw, or roll, etc.) and our own ship's speed (negative speeds accepted by the program for our own ship going backwards) to determine the enemy ship's starting coordinates, course and speed, and end coordinates. The program allows for our own ship's speed when calculating the enemy ship's speed, direction, and absolute coordinate data. The second program I called "BATCON1.BAS" because it is the battleship fire control program. It uses data from TICTAC.BAS, together with own ship's speed (the program accepts negative speeds for a ship going in reverse, with forward heading still assumed to be 360 degrees) and muzzle velocity of the gun, to obtain a fire and hit solution. The program assumes a flat sea, no air resistance, no wind speed and direction, no Coriolis effect, acceleration due to gravity of 32 feet per second per second, and no rotating earth. If the target is forward of the beam of our own ship, and if our own ship is moving forwards, then the program uses the following method. The program uses trial and error to find a point uprange from the target ship to shoot at--uprange by a distance equal to the distance travelled by our own ship during the shell's time of flight, and along a line parallel to our own ship's course. At the outset of the problem we have a "chicken and egg" situation. How far uprange (on the sea's surface) from the target to aim for (I called this point "point C" in the program) depends on the shell's time of flight, but the shell's time of flight depends on where point C is. The program reiterates itself until there is convergence between time of shell's flight and a time to traverse from point C to the enemy at our own ship's speed. If the target is on, or abaft the beam of, our own ship, and our own ship is moving forwards, then point C is downrange of the target. Generally speaking, though, the above method still applies. To put it another way: let's say that one wishes to do a drive by shooting from a boat at a buoy. The boat is travelling due north at 10 feet per second, and when the buoy is abeam of the boat's gun, it will be due west of the boat and gun. At that point, the buoy will be 6840 feet off the port beam. The gun's muzzle velocity is 800 feet per second, and to make the projectile go those 6840 feet over the water's surface, the gun must be elevated to 20 degrees, and the projectile will take 17.1 seconds to reach the target once fired. (I did all these calculations using formulae used in the subsequent programs on this web page.) This means that one must fire the gun directly off the port beam (i.e. due west of the boat) when the boat is 17.1 seconds away from reaching a point where the buoy is off the port beam. Point "C" is the point on the water's surface that you fire at, and it is 10 feet per second * 17.1 seconds = 17.1 feet due south of the buoy. (Distance = rate * time). Point "C" is also 6840 feet due west of the boat at the instant the gun is fired. If necessary, I could scan a diagram I made of the problem and sent it to those who are interested as a GIF file. Also I would encourage readers to try the programs and see if they give answers which make sense. As I mentioned above, I also wrote programs to do aircraft tracking (ATICTAC.BAS) and antiaircraft gun fire control (AIRBAT.BAS). I also wrote what I call a "space battleship fire control program" (SPACBAT.BAS, listing at the end of this file). This is similar to the antiaircraft gun fire control program I wrote (I should note here that it even takes as input data the results of my aircraft tracking program) and generates an antiaircraft or anti spacecraft gun firing solution, assuming you are in deep space with no gravity. There is a strange phenomenon I have observed when running BATCON1.BAS. Let's say that our own ship is heading due north, and that the target ship is behind us, abaft the port beam, and heading due west. The target ship has a red flag on its bow and a blue flag on its stern. We have crossed the target ship's path and it is re- ceding into the distance abaft the port beam. We look at the centreline of our own wake, and this is a line extending due south into the distance. We sweep our gaze to the right (westward) and see the target ship. Our own ship is travelling at 20 knots and the target ship is travelling at 5 knots. Due to the changing relative positions of the target ship and own ship, and the geometry of the overall layout of the situation, we see the target ship appear to get smaller; but it also appears to go backwards, apparently moving towards the centreline of our own wake (you can verify this the next time you are out boating). Running BATCON1.BAS, though, I noticed that in order to get a firing solution, one must make a further step. You have to aim at a point apparently behind the target ship in order to hit it. That is, when you sweep your gaze to the left beginning at the target ship (from own ship's taffrail), you see (1) the target ship's red flag, (2) the target ship's blue flag, and (3) the point to aim at, or point "C," a considerable apparent distance behind the target ship on the open water. Acknowledgement: I used information from the "Hyper Physics Trajectories" web page by R. Nave, web site address is http://hyperphysics.phy-astr.gsu.edu/hbase/traj.html although the formulae on that page are in the public domain. My email address is idouglas1304@yahoo.ca . Program ARCSIN1.BAS. Microsoft GW-BASIC has no arc sine command or statement included, so I had to write one of my own in order to do gun angle of elevation calculations. I obtained the algorithm for this from "Van Nostrand's Scientific Encyclopedia," 6th edition, Considine et al editors, published 1983, page 1646a. 10 PRINT "Calculate ELANG = arc sin (RHEXPR)" 20 PRINT "ELANG means gun's angle of elevation; RHEXPR is input argument for the series approximation polynomial used on right hand side of ELANG = arc sin (RHEXPR)" 30 PRINT:INPUT"Input argument for which arc sine is to be calculated ";RHEXPR 40 IF RHEXP>1 THEN PRINT "Input argument is out of range. ":GOTO 200 50 ELANG=RHEXPR+((RHEXPR^3)/6)+(.075*(RHEXPR^5)) 60 ELANG=ELANG+(4.464285714285714D-02*(RHEXPR^7)) 70 ELANG=ELANG+(3.038194444444444D-02*(RHEXPR^9)) 80 ELANG=ELANG+(2.237215909090909D-02*(RHEXPR^11)) 90 ELANG=ELANG+(1.735276442307692D-02*(RHEXPR^13)) 100 ELANG=ELANG+(.01396484375#*(RHEXPR^15)) 110 ELANG=ELANG+(1.155180089613971D-02*(RHEXPR^17)) 120 ELANG=ELANG+(9.761609529194079D-03*(RHEXPR^19)) 130 ELANG=ELANG+(8.390335809616816D-03*(RHEXPR^21)) 140 ELANG=ELANG+(7.312525877035907D-03*(RHEXPR^23)) 150 ELANG=ELANG+(6.447210311379962D-03*(RHEXPR^25)) 160 ELANG=ELANG+(5.74003766956053D-03*(RHEXPR^27)) 170 ELANG=ELANG+(.00515330968258#*(RHEXPR^29)) 180 PRINT:PRINT "Arc sin of input argument (for this purpose, gun's angle of elevation):" 190 PRINT ELANG 200 PRINT:END 210 REM Did a loopback test in September 2009 and this program gave a reasonable 220 REM final answer given that there were only 7 decimal places in each answer 230 REM given by this program ARCSIN1.BAS. 240 REM I.e. did arcsin (0.523598775 radians) with this program; 250 REM got 0.5510696; did sin of this in GWBASIC using a direct statement; 260 REM got 0.5235989 radians via ARCSIN1.bas as the final answer. Program TICTAC1.BAS 10 REM Begin flowchart box 1 20 CLEAR:CLS:PRINT "TICTAC program, that is, Time Interval" 30 PRINT "Compensating Tachymetry program, or" 40 PRINT "program to take data from two observations" 50 PRINT "of a moving target ship taken from" 60 PRINT "a moving ship and use" 70 PRINT "the data to determine the target ship's" 80 PRINT "starting and end X Y coordinates, azimuthal course" 90 PRINT "bearing, and speed. This program" 100 PRINT "is for IBM PC or similar computers" 110 PRINT "which do trigonometric calculations in radians." 120 PRINT "X y coordinates of own ship at start" 130 PRINT "are always assumed to be 50 nautical miles" 140 PRINT "by 50 nautical miles. Own ship's" 150 PRINT "true direction always assumed to be a" 160 PRINT "bearing of 0 degrees or 360 degrees." 170 PRINT "by Ian Douglas, year 2009." 180 PRINT "This program works for any azimuthal" 190 PRINT "bearing of target ship both at beginning and" 200 PRINT "end of time interval between the two observations." 210 REM End flowchart box 1 220 REM Begin flowchart box 2 230 PRINT "Enter the azimuthal bearing off own" 240 PRINT "ship of the target ship for first" 250 INPUT "observation (decimal degrees) :";ABT1:ABT1=ABT1/57.29577951# 260 PRINT "Enter the range to the target ship" 270 INPUT "for the first observation (feet): ";RT1 280 PRINT "Enter the azimuthal bearing off own" 290 PRINT "ship of the target ship for second" 300 INPUT "observation (decimal degrees) :";ABT2:ABT2=ABT2/57.29577951# 310 PRINT "Enter the range to the target ship" 320 INPUT "for the second observation (feet): ";RT2 330 PRINT "Enter the time elapsed in seconds between" 340 INPUT "The two observations ";TEBO 350 PRINT "Enter own ship's speed in" 360 INPUT "feet per second ";OSS 370 J1X=303805!:J1Y=303805!:REM Own 380 REM ship's starting coordinates in feet, i.e. 390 REM 50 nautical miles by 50 nautical miles 400 REM End flowchart box 2 410 REM Begin flowchart box 3 420 DOST=OSS*TEBO:REM DOST is distance 430 REM own ship travelled between the two observations 440 REM End flowchart box 3 450 REM Begin flowchart box 4 460 REM Calculate the coordinates of point J2, i.e. 470 REM end coordinates of own ship. J2's x 480 REM coordinate is simply J1's x coordinate 490 REM (303805 feet). J2's y coordinate is 500 REM J1's y coordinate (303805 feet) plus DOST. 510 J2X=J1X:J2Y=J1Y+DOST 520 REM End flowchart box 4 530 REM Begin flowchart box 5 540 REM Do a coordinated inverse (convert a vector and starting point 550 REM to an end point) to get the coordinates 560 REM of point K1, using point J1's coordinates 570 REM as the starting coordinates 580 REM Calculate end coordinates of vector 590 REM from angle, distance, and starting coordinates. 600 REM Distance is range to target 610 REM from first observation; azimuthal 620 REM bearing is that from first observation. 630 K1X=J1X+(RT1*SIN(ABT1)):K1Y=J1Y+(RT1*COS(ABT1)) 640 REM End flowchart box 5 650 REM Begin flowchart box 6 660 REM Do another coordinate inverse (as before, 670 REM convert a vector and starting point to an end point) 680 REM to get the coordinates of point K2, 690 REM using the coordinates of point J2 700 REM as the starting coordinates and 710 REM observational data (range, bearing) 720 REM from point J2. 730 REM Calculate end coordinates of vector 740 REM from angle, distance, and starting coordinates. 750 REM Distance is range to target 760 REM from second observation; azimuthal 770 REM bearing is that from second observation. 780 K2X=J2X+(RT2*SIN(ABT2)):K2Y=J2Y+(RT2*COS(ABT2)) 790 REM End flowchart box 6 800 REM Begin flowchart box 7 810 REM Do a coordinate inverse ("Convert 2 820 REM Cartesian points to vector") to get the 830 REM bearing and distance from K1 to K2. 840 DXAA=K2X-K1X:DYAA=K2Y-K1Y 850 DTAR=SQR((DXAA^2)+(DYAA^2)):REM Distance travelled by target 860 REM between the two observations. 870 IF DXAA=0 THEN GOTO 930 880 REM Above line prevents division by zero in next line 890 SLOPE=DYAA/DXAA:REM Slope 900 ANG1=ATN(SLOPE):REM Angle referred to standard position, calculated from slope 910 IF DXAA<0 THEN ANG1=ANG1+3.141592654# 920 IF DXAA>0 AND DYAA<0 THEN ANG1=ANG1+6.283185307# 930 IF DXAA=0 AND DYAA=0 THEN ANG1=0 940 IF DXAA=0 AND DYAA>0 THEN ANG1=1.570796327# 950 IF DXAA<0 AND DYAA=0 THEN ANG1=3.141592654# 960 IF DXAA=0 AND DYAA<0 THEN ANG1=4.712388981# 970 IF DXAA>0 AND DYAA=0 THEN ANG1=6.283185307# 980 AN1Z=6.283185308#-ANG1+1.570796327# 990 IF AN1Z>6.283185308# THEN AN1Z=AN1Z-6.283185308# 1000 REM Previous calculation gives azimuth from 1010 REM K1 to K2; 360 degrees minus angle plus 90 degrees 1020 REM End flowchart box 7 1030 REM Begin flowchart box 8 1040 REM Calculate target ship's speed 1050 TSS=DTAR/TEBO 1060 REM Recall that DTAR is distance travelled 1070 REM by target ship over TEBO, 1080 REM which is time elapsed between observations. 1090 REM End flowchart box 8 1100 REM Begin flowchart box 9 1110 PRINT "Target ship's speed (feet/sec): ";TSS 1120 PRINT "Target ship's azimuthal course" 1130 PRINT "bearing: ";AN1Z*57.29577951#;" decimal degrees of arc" 1140 PRINT "Target ship's starting coordinates (feet): ";K1X,K1Y 1150 PRINT "Target ship's end coordinates (feet): ";K2X,K2Y 1160 REM End flowchart box 9 1170 REM Begin flowchart box 10 1180 END 1190 REM End flowchart box 10 Program BATCON1.BAS 10 REM Begin flowchart box 1 20 CLEAR:CLS:PRINT "Battleship fire control program" 30 PRINT "for IBM PC or similar computers" 40 PRINT "which do trigonometric calculations in radians" 50 PRINT "by Ian Douglas, year 2007" 60 PRINT "This program works for any" 70 PRINT "azimuthal bearing of target ship" 80 PRINT "both at beginning and end of" 90 PRINT "shell's time of flight." 100 PRINT "Also works for any" 110 PRINT "bearing of point to aim at," 120 PRINT "i.e. point 'C.'" 130 PRINT "X y coordinates of own ship at start" 140 PRINT "of problem are 50 nautical miles by 50 nautical miles." 150 PRINT "own ship's true direction always assumed to" 160 PRINT "be a bearing of 0 degrees or 360 degrees." 170 REM End flowchart box 1 180 REM Begin flowchart box 2 190 INPUT "Muzzle velocity in feet/sec ";MV 200 INPUT "Target ship's present range in feet ";TSR 210 INPUT "Target ship's true speed in feet/sec ";TSS 220 PRINT "Target ship's modified true plot direction of travel as an azimuthal bearing in" 230 INPUT "decimal degrees of arc "; TSD:TSD=TSD/57.29577951# 240 INPUT "Own ship's true speed in feet/sec";OSS 250 PRINT "Target ship's bearing off own ship" 260 INPUT "at time of gun's firing in decimal degrees of arc ";TSB:TSB=TSB/57.29577951# 270 G=32:REM acceleration due to gravity in feet per second per second 280 MAXRAN=((MV*MV)*SIN(1.570796327#))/G:REM Maximum theoretical range of gun 290 MAXRAN=MAXRAN*.85:REM Must use 85% (or thereabouts) of maximum theoretical range of gun 300 REM as maximum program-workable range of gun in order to prevent 310 REM possible looping, overflow, etc. errors associated with attempting 320 REM to fire gun at excessive ranges (albeit in some cases still under 330 REM maximum theoretical range of gun). 340 REM End flowchart box 2 350 REM Begin flowchart box 3 360 REM Calculate target ship's starting (i.e. initially sighted) coordinates in (feet, feet) 370 TSPX=TSR*SIN(TSB)+303805! 380 TSPY=TSR*COS(TSB)+303805! 390 REM Above addition of 303805 is number of feet in 50 nautical miles; 400 REM recall that own ship's position is (50 naut. miles, 50 naut. miles) 410 REM Calculate first approximation of gun's angle of elevation 420 THETA1=(TSR*G)/(MV*MV) 430 REM The following arc sine routine 'refines' THETA1 to an unnecessarily high 440 REM level of precision since all we need is an initial approximation 450 REM of the angle of the gun's elevation, and since we use a different 460 REM formula, later on, to calculate this angle; also the different formula 470 REM formula we use later on does not require calculation of y = arc sin (x). 480 THETA1=THETA1+((THETA1^3)/6)+(.075*(THETA1^5)) 490 THETA1=THETA1+(4.464285714285714D-02*(THETA1^7)) 500 THETA1=THETA1+(3.038194444444444D-02*(THETA1^9)) 510 THETA1=THETA1+(2.237215909090909D-02*(THETA1^11)) 520 THETA1=THETA1+(1.735276442307692D-02*(THETA1^13)) 530 THETA1=THETA1+(.01396484375#*(THETA1^15)) 540 THETA1=THETA1+(1.155180089613971D-02*(THETA1^17)) 550 THETA1=THETA1+(9.761609529194079D-03*(THETA1^19)) 560 THETA1=THETA1+(8.390335809616816D-03*(THETA1^21)) 570 THETA1=THETA1+(7.312525877035907D-03*(THETA1^23)) 580 THETA1=THETA1+(6.447210311379962D-03*(THETA1^25)) 590 THETA1=THETA1+(5.74003766956053D-03*(THETA1^27)) 600 THETA1=THETA1+(.00515330968258#*(THETA1^29)) 610 THETA1=THETA1/2 620 REM Angle of elevation is theta1; used formula on p. 12 of notes. 630 PRINT "Target ship's starting (i.e. initially sighted) position in (feet, feet): " 640 PRINT TSPX,TSPY 650 PRINT "First approximation of gun's angle of elevation" 660 PRINT "(decimal degrees): ";THETA1*57.29577951# 670 REM Previous 3 lines are for debugging and error isolation only; delete or rem them out later 680 REM End flowchart box 3 690 REM Begin flowchart box 4 700 TOF1=2*((MV*SIN(THETA1))/G) 710 REM TOF1 is first approximation of shell's time of flight 720 PRINT "First approximation of shell's time of flight in seconds:" 730 PRINT TOF1 740 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 750 REM End flowchart box 4 760 REM Begin flowchart box 5 770 BEEP:REM Gets computer to beep to show that the program is looping back properly 780 REM on first (presumably wrong) guess of time of flight; need to use 790 REM BEEP since I cannot get any of my computers to 'listen' to the 800 REM scroll lock key on the keyboard(s). BEEP statement also delays operation 810 REM of program so you can see on screen whether proper reiterative looping 820 REM to get progressively better answers is taking place. 830 REM First step in this flowchart box is to get distance target ship 840 REM travelled during TOF (time of flight) calculated above. 850 DTST=TSS*TOF1:REM DTST is distance target ship travels during TOF1 860 REM (time of flight) calculated above 870 REM Now calculate how far own ship travelled during time of shell's flight 880 DOST=OSS*TOF1:REM DOST is how far own ship travelled during shell's 890 REM time of flight 900 REM Next step is to use section B of coordinate inverse program pp. 75-76 910 REM to calculate the end coordinates of the target ship at end of time 920 REM of flight, using azimuth of target ship, distance travelled, and 930 REM target's starting coordinates. 940 REM Input is target ship's initial sighted coordinates (TSPX, TSPY), 950 REM target ship's modified true plot direction of travel azimuthal bearing 960 REM TSD and distance target ship travels during shell's time of flight DTST. 970 REM Output is target ship's endpoint coordinates (TSEPX, TSEPY). 980 REM Calculate end coordinates 990 REM of vector, from angle, distance 1000 REM and starting coordinates 1010 TSEPX=TSPX+(DTST*SIN(TSD)):TSEPY=TSPY+(DTST*COS(TSD)) 1020 REM End calculation of target ship's end coordinates 1030 REM Begin calculation of coordinates of point "C" (PCX, PCY) 1040 REM X coordinate of point "C" is simply x coordinate of target ship's end point 1050 REM at end of shell's time of flight, i.e. PCX=TSEPX 1060 PCX=TSEPX 1070 REM Y coordinate of point "C" is the y coordinate of target ship's endpoint at end of time of flight of shell TSEPY 1080 REM minus distance own ship travels during shell's time of flight 1090 PCY=TSEPY-DOST 1100 REM End calculation of coordinates of point "C" 1110 PRINT "Distance target ship travels during shell's time of flight (feet):" 1120 PRINT DTST 1130 PRINT "Distance own ship travels during shell's time of flight (feet):" 1140 PRINT DOST 1150 PRINT "Target ship's end coordinates (feet, feet):" 1160 PRINT TSEPX,TSEPY 1170 PRINT "Coordinates of point 'C' (feet, feet):" 1180 PRINT PCX,PCY 1190 REM Previous 8 lines are for debugging and error isolation only; delete or rem them out later 1200 REM End flowchart box 5 1210 REM Begin flowchart box 6 1220 REM Calculate straight line distance between own ship at time of gun's firing 1230 REM and point "C," DOSPC 1240 REM Recall that own ship's coordinates at time of gun's firing are 1250 REM (303805 feet, 303805 feet) (i.e. 50 nautical miles by 50 nautical miles) 1260 DOSPC=SQR(((PCX-303805!)^2)+((PCY-303805!)^2)) 1270 IF DOSPC>MAXRAN THEN PRINT "Target ship is out of range. ":GOTO 1760 1280 PRINT "Distance between own ship and point 'C' when gun is fired (feet):" 1290 PRINT DOSPC 1300 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 1310 REM End flowchart box 6 1320 REM Begin flowchart box 7 1330 THETA2=(DOSPC*G)/(MV*MV) 1340 TOF2=2*(MV*SIN(THETA2))/G 1350 PRINT "New elevation and time of flight (decimal degrees of arc and seconds of time):" 1360 PRINT THETA2*57.29577951#,TOF2 1370 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 1380 REM End flowchart box 7 1390 REM Begin flowchart box 8 1400 IF (ABS(TOF2-TOF1))<.001 THEN GOTO 1490 1410 REM Note that the above line of code gives 1/1000 of a second accuracy!! I have found through 1420 REM experimentation that seeking higher accuracy will in some cases cause 1430 REM infinite looping of the program in a search for unnecessarily high accuracy. 1440 REM It is debatable whether in this program you need even 1/100 of a second accuracy. 1450 REM End flowchart box 8 1460 REM Begin flowchart box 9 1470 TOF1=TOF2:GOTO 760 1480 REM End flowchart box 9 1490 REM Begin flowchart box 10 1500 REM End flowchart box 10 1510 REM Begin flowchart box 11 1520 REM Calculate gun's bearing using coordinates of point "C" 1530 REM From flowchart box 5, and using coordinates of own ship at time of gun's firing; 1540 REM Convert points to vector 1550 REM 1st point is coordinates of own ship at gun's time of firing (303805 feet, 303805 feet) 1560 REM 2nd point is coordinates of point "C," (PCX, PCY) 1570 DXAA=PCX-303805!:DYAA=PCY-303805! 1580 IF DXAA=0 THEN GOTO 1640 1590 REM Above line prevents division by zero in next line 1600 SLOPE=DYAA/DXAA:REM Slope 1610 ANGLE#=ATN(SLOPE):REM Angle referred to standard position, calculated from slope 1620 IF DXAA<0 THEN ANGLE#=ANGLE#+3.141592654# 1630 IF DXAA>0 AND DYAA<0 THEN ANGLE#=ANGLE#+6.283185307# 1640 IF DXAA=0 AND DYAA=0 THEN ANGLE#=0 1650 IF DXAA=0 AND DYAA>0 THEN ANGLE#=1.570796327# 1660 IF DXAA<0 AND DYAA=0 THEN ANGLE#=3.141592654# 1670 IF DXAA=0 AND DYAA<0 THEN ANGLE#=4.712388981# 1680 IF DXAA>0 AND DYAA=0 THEN ANGLE#=6.283185307# 1690 REM End flowchart box 11 1700 REM Begin flowchart box 12 1710 REM Gun's angle of elevation is previously calculated THETA2 1720 PRINT "Gun's azimuth to shoot along: ";GUNAZ#*57.29577951#;" decimal degrees of arc" 1730 PRINT "Gun's angle of elevation: ";THETA2*57.29577951#;" decimal degrees of arc" 1740 REM If latest value for time of flight is to be printed here, that value is TOF2. 1750 REM End flowchart box 12 1760 REM Begin flowchart box 13 1770 END 1780 REM End flowchart box 13 Program UPART.BAS: this program uses a trial and error method to calculate the angle of elevation of a gun needed to hit a target of higher elevation than the gun. The formula to do this is a transcendental function (similar to Kepler's Equation) and cannot be solved directly. This program assumes that the target is within range and will simply crash or go into an infinite loop if the target whose coordinates are entered into the program is out of range. See the Appendix for a discussion of how to determine whether a target of higher elevation than the gun is within the gun's range, and for a more advanced version of this program which traps such out-of-range-of-the-gun situations. Also note that this program goes into an infinite loop if the fire and hit solution involves the gun having an angle of elevation of greater than about 78 degrees (the true value is higher but I have not determined this exactly so to be safe I have given a lower and rounded off figure). For example, if your gun has a muzzle velocity of 611.8823 feet per second and the target is separated from the gun vertically and horizontally by 4900 feet and 1460 feet, then the program will go into an infinite loop; but if the horizontal separation increases to 1470 feet then the program will tell you to elevate the gun to 78.4921 decimal degrees above horizontal. In other words, the program will not work for targets directly overhead and above the gun or nearly so. I hope at some time to address this great weakness but for the time being the program's end user must beware of this flaw. On the other hand, note that this weakness creeps in when the gun is firing only 12 degrees of arc--or less--from vertical. I am still debating how to handle the last mentioned types of situation; such instances mean that the program has a blind spot where a fire and hit solution is physically possible but not practically computable. There are two ways to handle such instances. One is to get the program to end if it tries too many times to find an answer to the question of what the gun's angle of elevation should be to hit the target. The other method is to get the program to end if it starts to use, as trial angles of elevation, angles somewhere between 75 and 78 degrees. The first method is virtually guaranteed to work but might take some time to be processed by a program; the second method is potentially faster but could be dangerous, as the program might accidentally cross the threshold of "too high" an angle of elevation by too great an increment. This having occurred, the program would run unchecked into the domain where it cannot use practical means to find an answer--plunging itself into an infinite loop. At some point in the future I hope to make a decision and incorporate the decision into a revised program UPART.BAS and other programs on this web site related to it, e.g. AIRBAT.BAS. 10 REM Begin flowchart box E1 (begin/end marker) 20 REM End flowchart box E1; the flowchart numbers in this program correspond to flowchart on p. 251 of notes. 30 REM Begin flowchart box E2 40 PRINT "Program to use a variant of successive approximation to get angle of elevation" 50 PRINT "of a gun firing at a target separated from the gun significantly" 60 PRINT "separated from it both vertically and horizontally (e.g. an aircraft) using equation" 70 PRINT "211a on page 211 of notes; by Ian Douglas, year 2010. This program also works" 80 PRINT "for targets below gun. This program is for upward-only trajectories if target is above gun" 90 PRINT "and not for upward-to-peak-to-downward trajectories." 100 REM Program is intended to make the problem manageable by taking the maximum possible shell height the gun 110 REM can attain and obtaining a constant to multiply 90 degrees of arc (1.570796327 radians) (gun's maximum angle of 120 REM elevation) by to get any given possible height the shell can attain; this enables a resultant vertical 130 REM separation of the target from gun to be compared to an actual vertical separation attained by a trial angle 140 REM of elevation of gun, for comparison purposes; also this enables the successive approximation process to occur. 150 REM For this reason, I do not recommend using this program to compute data concerning guns with greater 160 REM than 10,000 feet per second muzzle velocities; but even the most powerful guns have about 1/4 that. 170 REM Also note that the program will go into an infinite loop if the target is directly above gun, even though 180 REM a fire and hit solution is possible in such a case. 190 PRINT "Vertical separation (Y) of target from gun (feet) " 200 PRINT "(negative for a target below gun)" 210 INPUT VESEP 220 INPUT "Horizontal separation (in XY plane) of target from gun (feet) ";HOSEP 230 INPUT "Gun's muzzle velocity (feet per second) ";MV 240 G=32:G1=G/2:REM Acceleration due to gravity in feet per second per second; also half this for use in 250 REM instances of equation 211a. 260 REM End flowchart box E2 270 REM Begin flowchart box E3 280 REM Calculate the shell's maximum attainable height. 290 MAXHT=(MV*MV)/(2*G) 300 REM Calculate the constant mentioned in remarks for flowchart box E2. 310 CONST=MAXHT/1.570796327# 320 REM Calculate the minimum angle of elevation of gun needed to reach HOSEP. 330 MIELAN=(HOSEP*G)/(MV*MV) 340 REM The following is an arc sine routine with the final angle divided by 2 for the purpose of the 350 REM above formula to get the minimum necessary angle of elevation given HOSEP, muzzle velocity, and G. 360 MIELAN=MIELAN+((MIELAN^3)/6)+(.075*(MIELAN^5)) 370 MIELAN=MIELAN+(4.464285714285714D-02*(MIELAN^7)) 380 MIELAN=MIELAN+(3.038194444444444D-02*(MIELAN^9)) 390 MIELAN=MIELAN+(2.237215909090909D-02*(MIELAN^11)) 400 MIELAN=MIELAN+(1.735276442307692D-02*(MIELAN^13)) 410 MIELAN=MIELAN+(.01396484375#*(MIELAN^15)) 420 MIELAN=MIELAN+(1.155180089613971D-02*(MIELAN^17)) 430 MIELAN=MIELAN+(9.761609529194079D-03*(MIELAN^19)) 440 MIELAN=MIELAN+(8.390335809616816D-03*(MIELAN^21)) 450 MIELAN=MIELAN+(7.312525877035907D-03*(MIELAN^23)) 460 MIELAN=MIELAN+(6.447210311379962D-03*(MIELAN^25)) 470 MIELAN=MIELAN+(5.74003766956053D-03*(MIELAN^27)) 480 MIELAN=MIELAN+(.00515330968258#*(MIELAN^29)) 490 MIELAN=MIELAN/2 500 REM Angle of elevation is MIELAN; used formula on p. 12 of notes. 510 REM Set shell's minimum possible height as zero feet; set gun's maximum angle of elevation as 90 degrees of arc. 520 REM End flowchart box E3 530 REM Begin flowchart box E4 540 REM Set first trial angle of elevation of gun as (45 degrees of arc minus minimum angle of elevation of gun 550 REM to reach HOSEP)/2; the reason I use 45 and not 90 degrees here is that as soon as you use an angle greater 560 REM than 45 degrees you meet up with the possibility of the shell falling short of HOSEP. 570 ELANG=(.785398163#-MIELAN)/2 580 REM End flowchart box E4 590 REM Begin flowchart box E5 600 REM Calculate trial vertical (Y) separation resulting from ELANG and inputted values from flowchart box E2 610 REM (except VESEP) using equation 211a; calculate difference between trial vertical separation and VESEP. 620 A1=MV*SIN(ELANG):A2=MV*COS(ELANG) 630 TRIZ=(HOSEP*A1)/A2 640 TRIZ=TRIZ-G1*((HOSEP*HOSEP)/(A2*A2)) 650 REM PRINT "Trial Z and trial gun's angle of elevation are ";TRIZ,ELANG*57.29577951#:REM For debugging only 660 REM End flowchart box E5 670 REM Begin flowchart box E6 680 IF ABS(TRIZ-VESEP)<.01 THEN GOTO 770:REM That is, go to flowchart box E8. Note that this represents a 1/100 690 REM of a foot accuracy of the shell's striking point in relation to the program's perception of the target, 700 REM where the target is an aircraft! Compare this with the size of any conventional aircraft. 710 REM I have tried to get the program to achieve a higher level of precision, but this only caused 720 REM endless looping due to the problem of attainable accuracy of TRIZ and ELANG, loss of digits, etc. 730 REM End flowchart box E6 740 REM Begin flowchart box E7 750 ELANG=ELANG-(((TRIZ-VESEP)/2)/CONST):GOTO 590:REM That is, go to flowchart box E5 760 REM End flowchart box E7 770 REM Begin flowchart box E8 780 PRINT:PRINT "Gun's angle of elevation from horizontal (depression below horizontal negative) = "; 790 PRINT ELANG*57.29577951#;" decimal degrees of arc":PRINT 800 REM End flowchart box E8 810 REM Begin flowchart box E9 (begin/end marker) 820 END 830 REM End flowchart box E9 Program ATICTAC.BAS 10 REM Begin flowchart box 1 20 CLEAR:CLS:PRINT "ATICTAC program, that is, Aircraft Time Interval" 30 PRINT "Compensating Tachymetry program, or" 40 PRINT "program to take data from two observations" 50 PRINT "of a moving target aircraft taken from" 60 PRINT "a moving ship and use" 70 PRINT "the data to determine the target aircraft's" 80 PRINT "starting and end XYZ coordinates, azimuthal course" 90 PRINT "bearing, speed, and rate of climb (descent negative). This program" 100 PRINT "is for IBM PC or similar computers" 110 PRINT "which do trigonometric calculations in radians." 120 REM Note that X and Y coordinates are 130 REM ground coordinates and Z 140 REM is altitiude above own ship. 150 PRINT:PRINT "X Y coordinates of own ship at start" 160 PRINT "are always assumed to be 50 nautical miles" 170 PRINT "by 50 nautical miles. Own ship's" 180 PRINT "true direction always assumed to be a" 190 PRINT "bearing of 0 degrees or 360 degrees." 200 PRINT:PRINT "by Ian Douglas, year 2010." 210 PRINT:PRINT "This program works for any azimuthal" 220 PRINT "bearing of target aircraft both at beginning and" 230 PRINT "end of time interval between the two observations." 240 REM This program assumes altitude of own ship always to be 0 feet and works only for positive altitudes of aircraft. 250 REM End flowchart box 1 260 REM Begin flowchart box 2 270 PRINT:PRINT "Enter the azimuthal bearing off own" 280 PRINT "ship of the target aircraft for first" 290 INPUT "observation (decimal degrees) :";ABT1:ABT1=ABT1/57.29577951# 300 PRINT "Enter the straight line of sight range to the target aircraft" 310 INPUT "for the first observation (feet): ";RT1 320 PRINT "Enter the observed angle of elevation of the target aircraft in decimal degrees" 330 INPUT "above horizontal for the first observation";ELANG1:ELANG1=ELANG1/57.29577951# 340 PRINT:PRINT "Enter the azimuthal bearing off own" 350 PRINT "ship of the target aircraft for second" 360 INPUT "observation (decimal degrees) :";ABT2:ABT2=ABT2/57.29577951# 370 PRINT "Enter the straight line of sight range to the target aircraft" 380 INPUT "for the second observation (feet): ";RT2 390 PRINT "Enter the observed angle of elevation of the target aircraft in decimal degrees" 400 INPUT "above horizontal for the second observation";ELANG2:ELANG2=ELANG2/57.29577951# 410 PRINT:PRINT "Enter the time elapsed in seconds between" 420 INPUT "The two observations ";TEBO 430 PRINT:PRINT "Enter own ship's speed in" 440 INPUT "feet per second ";OSS 450 J1X=303805!:J1Y=303805!:REM Own 460 REM Own ship's starting coordinates in feet, i.e. 470 REM 50 nautical miles by 50 nautical miles 480 REM End flowchart box 2 490 REM Begin flowchart box 3 500 ALT1=RT1*(SIN(ELANG1)):ALT2=RT2*(SIN(ELANG2)):REM Computes altitudes of aircraft for both observations 510 RCLMB=(ALT2-ALT1)/TEBO:REM rate of climb of aircraft (descent negative) 520 HRAN1=RT1*(COS(ELANG1)):HRAN2=RT2*(COS(ELANG2)):REM Computes horizontal ground track ranges to target aircraft for both observations 530 DOST=OSS*TEBO:REM DOST is distance 540 REM own ship travelled between the two observations 550 REM End flowchart box 3 560 REM Begin flowchart box 4 570 REM Calculate the coordinates of point J2, i.e. 580 REM end coordinates of own ship. J2's x 590 REM coordinate is simply J1's x coordinate 600 REM (303805 feet). J2's y coordinate is 610 REM J1's y coordinate (303805 feet) plus DOST. 620 J2X=J1X:J2Y=J1Y+DOST 630 REM End flowchart box 4 640 REM Begin flowchart box 5 650 REM Do a coordinate inverse (convert a vector and starting point 660 REM to an end point) to get the coordinates 670 REM of point K1, using point J1's coordinates 680 REM as the starting coordinates 690 REM Calculate end coordinates of vector 700 REM from angle, distance, and starting coordinates. 710 REM Distance is horizontal ground track range to target aircraft 720 REM from first observation; azimuthal 730 REM bearing is that from first observation. 740 K1X=J1X+(HRAN1*SIN(ABT1)):K1Y=J1Y+(HRAN1*COS(ABT1)):K1Z=ALT1 750 REM End flowchart box 5 760 REM Begin flowchart box 6 770 REM Do another coordinate inverse (as before, 780 REM convert a vector and a starting point to an end point) 790 REM to get the coordinates of point K2, 800 REM using the coordinates of point J2 810 REM as the starting coordinates and 820 REM observational data (range, bearing) 830 REM from point J2. 840 REM Calculate end coordinates of vector 850 REM from angle, distance, and starting coordinates. 860 REM Distance is horizontal ground track range to target 870 REM from second observation; azimuthal 880 REM bearing is that from second observation. 890 K2X=J2X+(HRAN2*SIN(ABT2)):K2Y=J2Y+(HRAN2*COS(ABT2)):K2Z=ALT2 900 REM End flowchart box 6 910 REM Begin flowchart box 7 920 REM Do a coordinate inverse ("Convert 2 930 REM Cartesian points to vector") to get the 940 REM bearing and distance from K1 to K2. 950 DXAA=K2X-K1X:DYAA=K2Y-K1Y:DZAA=K2Z-K1Z 960 DTAR=SQR((DXAA^2)+(DYAA^2)):REM Ground track distance travelled by target aircraft 970 REM between the two observations. 980 IF DXAA=0 THEN GOTO 1040 990 REM Above line prevents division by zero in next line 1000 SLOPE=DYAA/DXAA:REM Slope 1010 ANG1=ATN(SLOPE):REM Angle referred to standard position, calculated from slope 1020 IF DXAA<0 THEN ANG1=ANG1+3.141592654# 1030 IF DXAA>0 AND DYAA<0 THEN ANG1=ANG1+6.283185307# 1040 IF DXAA=0 AND DYAA=0 THEN ANG1=0 1050 IF DXAA=0 AND DYAA>0 THEN ANG1=1.570796327# 1060 IF DXAA<0 AND DYAA=0 THEN ANG1=3.141592654# 1070 IF DXAA=0 AND DYAA<0 THEN ANG1=4.712388981# 1080 IF DXAA>0 AND DYAA=0 THEN ANG1=6.283185307# 1090 AN1Z=6.283185308#-ANG1+1.570796327# 1100 IF AN1Z>6.283185308# THEN AN1Z=AN1Z-6.283185308# 1110 REM Previous calculation gives azimuth from 1120 REM K1 to K2; 360 degrees minus angle plus 90 degrees 1130 REM End flowchart box 7 1140 REM Begin flowchart box 8 1150 REM Calculate target aircraft's ground track speed 1160 TAS=DTAR/TEBO 1170 REM Recall that DTAR is ground track distance travelled 1180 REM by target aircraft over TEBO, 1190 REM which is time elapsed between observations. 1200 REM End flowchart box 8 1210 REM Begin flowchart box 9 1220 PRINT:PRINT "Target aircraft's ground track speed (feet/sec): ";TAS 1230 PRINT "Target aircraft's azimuthal course" 1240 PRINT "bearing: ";AN1Z*57.29577951#;" decimal degrees of arc" 1250 PRINT "Target aircraft's rate of" 1260 PRINT "climb (descent negative) (feet per second) ";RCLMB 1270 PRINT "Target aircraft's starting XYZ coordinates (feet): ";K1X,K1Y,K1Z 1280 PRINT "Target aircraft's end XYZ coordinates (feet): ";K2X,K2Y,K2Z 1290 REM End flowchart box 9 1300 REM Begin flowchart box 10 1310 END 1320 REM End flowchart box 10 Program AIRBAT.BAS : This program assumes that the target is within range both at the outset of the problem (initial acqisition of the target) and at the end of the problem (hitting the target). The program will simply crash or go into an infinite loop if the target whose coordinates are entered into the program is out of range. See the Appendix for a discussion of how to determine whether a target of higher elevation than the gun is within the gun's range. Also note that this program has a weakness also present in UPART.BAS, which is that it goes into an infinite loop if the fire and hit solution involves the gun having an angle of elevation of greater than about 78 degrees: see the introductory notes for the program UPART.BAS above for a fuller explanation. 10 REM Begin flowchart box F1 20 CLEAR:CLS:PRINT "Antiaircraft gun fire control program for firing at a moving aircraft from a moving ship." 30 PRINT "This program is for IBM PC or similar computers" 40 PRINT "which do trigonometric calculations in radians" 50 PRINT "by Ian Douglas, year 2010" 60 PRINT "This program works for any azimuthal bearing of target aircraft" 70 PRINT "both at beginning and end of shell's time of flight." 80 PRINT "Also works for any bearing of point to aim at, i.e. point 'C.'" 90 PRINT "Note that this program will not point the gun directly at point 'C' but to a point directly above it," 100 PRINT "not named or calculated in the program, and let the shell fall to point 'C' by gravity." 110 PRINT "X Y Z coordinates of own ship at start" 120 PRINT "of problem are 50 nautical miles by 50 nautical miles by 0." 130 PRINT "Own ship's true direction always assumed to" 140 PRINT "be a bearing of 0 degrees or 360 degrees." 150 REM This program assumes that it is working at or near the earth's surface, where the shell flies at the target 160 REM under gravitational influences. 170 PRINT "This program works for targets whose starting and/or ending at elevations are above or below own ship;" 180 PRINT "likewise for the shell's point of impact. If gun must be depressed to hit target, then gun's" 190 PRINT "angle of elevation is negative." 200 REM End flowchart box F1 210 REM Begin flowchart box F2 220 INPUT "Muzzle velocity in feet/sec ";MV 230 INPUT "Target aircraft's present straight line-of-sight range in feet ";TAR 240 INPUT "Target aircraft's true 'ground track' speed in feet/sec ";TAS 250 PRINT "Target aircraft's modified true plot 'ground track' direction of travel as an" 260 INPUT "azimuthal bearing in decimal degrees of arc "; TAD:TAD=TAD/57.29577951# 270 PRINT "Target aircraft's first observed angle of elevation above horizontal" 280 INPUT "in decimal degrees of arc ";TAAE:TAAE=TAAE/57.29577951# 290 PRINT "Target aircraft's rate of climb in feet/sec" 300 INPUT "(descent negative) ";TARC 310 INPUT "Own ship's true speed in feet/sec";OSS 320 PRINT "Target aircraft's bearing off own ship" 330 INPUT "at time of gun's firing in decimal degrees of arc ";TAB:TAB=TAB/57.29577951# 340 G=32:G1=G/2:REM Acceleration due to gravity in feet per second per second; also half this for use in 350 REM instances of equation 211a. 360 REM End flowchart box F2 370 REM Begin flowchart box F3 380 REM Calculate the shell's maximum attainable height for the purposes of calculating the value of the 390 REM constant mentioned in remarks for flowchart box E2 for subroutine 'A'; also do that second calculation. 400 REM I have moved the function of these calculations to here since subroutine 'A' 410 REM works inside a larger program, and since these calculations need only be done once. 420 MAXHT=(MV*MV)/(2*G):CONST=MAXHT/1.570796327# 430 REM Calculate target aircraft's flat range (straight line of sight range reduced to XY plane) and then 440 REM calculate beginning (i.e. initially sighted) coordinates in (feet, feet) 450 FLRANG=TAR*COS(TAAE) 460 TABX=(FLRANG*SIN(TAB))+303805! 470 TABY=(FLRANG*COS(TAB))+303805! 480 TABZ=TAR*SIN(TAAE) 490 REM Above addition of 303805 is number of feet in 50 nautical mile; 500 REM recall that own ship's position is (50 naut. mile, 50 naut. mile, 0) 510 PRINT "Target aircraft's starting (i.e. initially sighted) XYZ position in (feet, feet, feet): " 520 PRINT TABX,TABY,TABZ 530 REM Previous 3 lines are for debugging and error isolation only; delete or rem them out later 540 REM Calculate initial horizontal and vertical separation between own ship at time of gun's firing 550 REM and target aircraft (i.e. initial values of HOSEP and VESEP) 560 REM Recall that own ship's coordinates at time of gun's firing are 570 REM (303805 feet, 303805 feet, 0 feet) (i.e. 50 nautical miles by 50 nautical miles by 0) 580 HOSEP=SQR(((TABX-303805!)^2)+(TABY-303805!)^2) 590 VESEP=TABZ 600 PRINT "Horizontal and vertical distances between own ship and target aircraft:" 610 PRINT HOSEP,VESEP 620 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 630 GOSUB 1770:TOF1=TOF2:REM That is, go to subroutine A to get latest approximation of shell's time of flight TOF1, and 640 REM first approximation of gun's angle of elevation (this last not needed in first iterations of program, 650 REM but I wanted to have a single 'package' subroutine to do the necessary multiple tasks in later 660 REM iterations of the problem). Note that since the subroutine calls the latest time of flight TOF2, and 670 REM since we are in a position in the program where we need only the first approximation of time of flight 680 REM (and use this in the initial calculations immediately following this line of code) we name this result TOF1. 690 REM Later we will stick with the name the subroutine calls this time of flight, i.e. TOF2, unchanged. 700 PRINT "First approximation of shell's time of flight in seconds:" 710 PRINT TOF1 720 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 730 REM End flowchart box F3 740 REM Begin flowchart box F4 (off page connector) 750 REM End flowchart box F4 (off page connector) 760 REM Begin flowchart box F5 770 BEEP:REM Gets computer to beep to show that the program is looping back properly 780 REM on first (presumably wrong) guess of time of flight; need to use 790 REM BEEP since I cannot get any of my computers to 'listen' to the 800 REM scroll lock key on the keyboard(s). BEEP statement also delays operation 810 REM of program so you can see on screen whether proper reiterative looping 820 REM to get progressively better answers is taking place. 830 REM First step in this flowchart box is to get distance target aircraft 840 REM travelled during TOF (time of flight) calculated above. 850 DTAT=TAS*TOF1:REM DTAT is 'ground track' distance target aircraft travels during TOF1 860 REM (time of flight) calculated above 870 REM Now calculate how far own ship travelled during time of shell's flight 880 DOST=OSS*TOF1:REM DOST is 'ground track'distance own ship travelled during shell's 890 REM time of flight 900 REM Next step is to use section B of coordinate inverse program pp. 75-76 910 REM to calculate the end coordinates of the target aircraft at end of time 920 REM of flight, using azimuth of target aircraft, distance travelled, and 930 REM target aircraft's starting coordinates. 940 REM Input is target aircraft's initial sighted coordinates (TABX, TABY, TABZ), 950 REM target aircraft's modified true plot direction of travel azimuthal bearing 960 REM TAD and distance target aircraft travels during shell's time of flight DTAT. 970 REM Output is target aircraft's endpoint coordinates (TAEPX, TAEPY, TAEPZ). 980 REM Calculate end coordinates 990 REM of vector, from angle, distance 1000 REM and starting coordinates 1010 TAEPX=TABX+(DTAT*SIN(TAD)):TAEPY=TABY+(DTAT*COS(TAD)):TAEPZ=TABZ+(TARC*TOF1) 1020 REM End calculation of target aircraft's end coordinates 1030 REM Begin calculation of coordinates of point "C" (PCX, PCY, PCZ) 1040 REM X coordinate of point "C" is simply x coordinate of target aircraft's end point 1050 REM at end of shell's time of flight, i.e. PCX=TAEPX 1060 PCX=TAEPX 1070 REM Y coordinate of point "C" is the y coordinate of target aircraft's endpoint at end of time of flight of shell TAEPY 1080 REM minus distance own ship travels during shell's time of flight 1090 PCY=TAEPY-DOST 1100 REM Z coordinate of point "C" is simply target aircraft's altitude at end of shell's time of flight, i.e. TAEPZ 1110 PCZ=TAEPZ 1120 REM End calculation of coordinates of point "C" 1130 PRINT "'Ground track' distance target aircraft travels during shell's time of flight (feet):" 1140 PRINT DTAT 1150 PRINT "'Ground track' distance own ship travels during shell's time of flight (feet):" 1160 PRINT DOST 1170 PRINT "Target aircraft's XYZ end coordinates (feet, feet, feet):" 1180 PRINT TAEPX,TAEPY,TAEPZ 1190 PRINT "XYZ coordinates of point 'C' (feet, feet, feet):" 1200 PRINT PCX,PCY,PCZ 1210 REM Previous 8 lines are for debugging and error isolation only; delete or rem them out later 1220 REM End flowchart box F5 1230 REM Begin flowchart box F6 1240 REM Calculate horizontal and vertical separation between own ship at time of gun's firing 1250 REM and point "C," (i.e. HOSEP and VESEP) 1260 REM Recall that own ship's coordinates at time of gun's firing are 1270 REM (303805 feet, 303805 feet, 0 feet) (i.e. 50 nautical miles by 50 nautical miles by 0) 1280 HOSEP=SQR(((PCX-303805!)^2)+(PCY-303805!)^2) 1290 VESEP=PCZ 1300 PRINT "Horizontal and vertical distances between own ship and point 'C' when gun is fired (feet):" 1310 PRINT HOSEP,VESEP 1320 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 1330 GOSUB 1770:REM That is, go to subroutine A to get latest approximation of shell's time of flight TOF2, and 1340 REM latest approximation of gun's angle of elevation. 1350 REM End flowchart box F6 1360 REM Begin flowchart box F7 1370 IF (ABS(TOF2-TOF1))<.001 THEN GOTO 1470:REM That is, go to flowchart box F9 1380 REM Note that the above line of code gives 1/1000 of a second accuracy!! I have found through 1390 REM experimentation that seeking higher accuracy will in some cases cause 1400 REM infinite looping of the program in a search for unnecessarily high accuracy. 1410 REM It is debatable whether in this program you need even 1/100 of a second accuracy. To make the program get an 1420 REM answer in less time, reduce accuracy by ten or a hundred fold in the line of code 5 lines before this one. 1430 REM End flowchart box F7 1440 REM Begin flowchart box F8 1450 TOF1=TOF2:GOTO 760:REM That is, go to flowchart box F5 1460 REM End flowchart box F8 1470 REM Begin flowchart box F9 1480 REM Calculate gun's bearing (not angle of elevation yet) using only X Y coordinates of point "C" 1490 REM From flowchart box 4, and using coordinates of own ship at time of gun's firing; 1500 REM Convert points to vector 1510 REM 1st point is coordinates of own ship at gun's time of firing (303805 feet, 303805 feet) 1520 REM 2nd point is coordinates of point "C," (PCX, PCY) 1530 DXAA=PCX-303805!:DYAA=PCY-303805! 1540 IF DXAA=0 THEN GOTO 1600 1550 REM Above line prevents division by zero in next line 1560 SLOPE=DYAA/DXAA:REM Slope 1570 ANGLE#=ATN(SLOPE):REM Angle referred to standard position, calculated from slope 1580 IF DXAA<0 THEN ANGLE#=ANGLE#+3.141592654# 1590 IF DXAA>0 AND DYAA<0 THEN ANGLE#=ANGLE#+6.283185307# 1600 IF DXAA=0 AND DYAA=0 THEN ANGLE#=0 1610 IF DXAA=0 AND DYAA>0 THEN ANGLE#=1.570796327# 1620 IF DXAA<0 AND DYAA=0 THEN ANGLE#=3.141592654# 1630 IF DXAA=0 AND DYAA<0 THEN ANGLE#=4.712388981# 1640 IF DXAA>0 AND DYAA=0 THEN ANGLE#=6.283185307# 1650 GUNAZ#=6.283185307#-ANGLE#+1.570796327#:REM 360 degrees minus angle plus 90 degrees 1660 IF GUNAZ#>6.283185307# THEN GUNAZ#=GUNAZ#-6.283185307# 1670 REM End flowchart box F9 1680 REM Begin flowchart box F10 1690 REM Gun's angle of elevation is previously calculated ELANG 1700 PRINT "Gun's azimuth to shoot along: ";GUNAZ#*57.29577951#;" decimal degrees of arc" 1710 PRINT "Gun's angle of elevation: "; ELANG*57.29577951#;" decimal degrees of arc" 1720 REM If latest value for time of flight is to be printed here, that value is TOF2. 1730 REM End flowchart box F10 1740 REM Begin flowchart box F11 1750 END 1760 REM End flowchart box F11 1770 REM Begin flowchart box E1 (begin/end marker) 1780 REM End flowchart box E1; the flowchart numbers in this program correspond to flowchart on p. 251 of notes 1790 REM since this subroutine is simply (for the most part, except for the time of flight calculations 1800 REM at the end) a repeat of the program UPART.BAS. 1810 REM Begin flowchart box E2 1820 REM *** Subroutine 'A' ***, to use a variant of successive approximation to get angle of elevation 1830 REM of a gun firing at a target separated from the gun significantly 1840 REM separated from it both vertically and horizontally (e.g. an aircraft) using equation 1850 REM 211a on page 211 of notes; by Ian Douglas, year 2010. This program also works 1860 REM for targets below gun. This program is for upward-only trajectories if target is above gun 1870 REM and not for upward-to-peak-to-downward trajectories. This subroutine also calculates the 1880 REM shell's time of flight. 1890 REM Program is intended to make the problem manageable by taking the maximum possible shell height the gun 1900 REM can attain and obtaining a constant to multiply 90 degrees of arc (1.570796327 radians) (gun's maximum angle of 1910 REM elevation) by to get any given possible height the shell can attain; this enables a resultant vertical 1920 REM separation of the target from gun to be compared to an actual vertical separation attained by a trial angle 1930 REM of elevation of gun, for comparison purposes; also this enables the successive approximation process to occur. 1940 REM For this reason, I do not recommend using this program to compute data concerning guns with greater 1950 REM than 10,000 feet per second muzzle velocities; but even the most powerful guns have about 1/4 that. 1960 REM End flowchart box E2 1970 REM Begin flowchart box E3 1980 REM If this subroutine were working on its own, I would have put the calculations for the shell's maximum 1990 REM attainable height and for the constant needed in this subroutine here; but since this subroutine 2000 REM functions inside a larger program, and since these calculations need only be done once, I have 2010 REM left them to flowchart box F3. 2020 REM Calculate the minimum angle of elevation of gun needed to reach HOSEP. 2030 MIELAN=(HOSEP*G)/(MV*MV) 2040 REM The following is an arc sine routine with the final angle divided by 2 for the purpose of the 2050 REM above formula to get the minimum necessary angle of elevation given HOSEP, muzzle velocity, and G. 2060 MIELAN=MIELAN+((MIELAN^3)/6)+(.075*(MIELAN^5)) 2070 MIELAN=MIELAN+(4.464285714285714D-02*(MIELAN^7)) 2080 MIELAN=MIELAN+(3.038194444444444D-02*(MIELAN^9)) 2090 MIELAN=MIELAN+(2.237215909090909D-02*(MIELAN^11)) 2100 MIELAN=MIELAN+(1.735276442307692D-02*(MIELAN^13)) 2110 MIELAN=MIELAN+(.01396484375#*(MIELAN^15)) 2120 MIELAN=MIELAN+(1.155180089613971D-02*(MIELAN^17)) 2130 MIELAN=MIELAN+(9.761609529194079D-03*(MIELAN^19)) 2140 MIELAN=MIELAN+(8.390335809616816D-03*(MIELAN^21)) 2150 MIELAN=MIELAN+(7.312525877035907D-03*(MIELAN^23)) 2160 MIELAN=MIELAN+(6.447210311379962D-03*(MIELAN^25)) 2170 MIELAN=MIELAN+(5.74003766956053D-03*(MIELAN^27)) 2180 MIELAN=MIELAN+(.00515330968258#*(MIELAN^29)) 2190 MIELAN=MIELAN/2 2200 REM Angle of elevation is MIELAN; used formula on p. 12 of notes. 2210 REM Set shell's minimum possible height as zero feet; set gun's maximum angle of elevation as 90 degrees of arc. 2220 REM End flowchart box E3 2230 REM Begin flowchart box E4 2240 REM Set first trial angle of elevation of gun as (45 degrees of arc minus minimum angle of elevation of gun 2250 REM to reach HOSEP)/2; the reason I use 45 and not 90 degrees here is that as soon as you use an angle greater 2260 REM than 45 degrees you meet up with the possibility of the shell falling short of HOSEP. 2270 ELANG=(.785398163#-MIELAN)/2 2280 REM End flowchart box E4 2290 REM Begin flowchart box E5 2300 REM Calculate trial vertical (Z) separation resulting from ELANG and inputted values from flowchart box E2 2310 REM (except VESEP) using equation 211a; calculate difference between trial vertical separation and VESEP. 2320 A1=MV*SIN(ELANG):A2=MV*COS(ELANG) 2330 TRIZ=(HOSEP*A1)/A2 2340 TRIZ=TRIZ-G1*((HOSEP*HOSEP)/(A2*A2)) 2350 REM PRINT "Trial Z and trial gun's angle of elevation are ";TRIZ,ELANG*57.29577951#:REM For debugging only 2360 REM End flowchart box E5 2370 REM Begin flowchart box E6 2380 IF ABS(TRIZ-VESEP)<.01 THEN GOTO 2470:REM That is, go to flowchart box E8. Note that this represents a 1/100 2390 REM of a foot accuracy of the shell's striking point in relation to the program's perception of the target, 2400 REM where the target is an aircraft! Compare this with the size of any conventional aircraft. 2410 REM I have tried to get the program to achieve a higher level of precision, but this only caused 2420 REM endless looping due to the problem of attainable accuracy of TRIZ and ELANG, loss of digits, etc. 2430 REM End flowchart box E6 2440 REM Begin flowchart box E7 2450 ELANG=ELANG-(((TRIZ-VESEP)/2)/CONST):GOTO 2290:REM That is, go to flowchart box E5 2460 REM End flowchart box E7 2470 REM Begin flowchart box E8 2480 REM PRINT:PRINT "Gun's angle of elevation from horizontal (depression below horizontal negative) = "; 2490 REM PRINT ELANG*57.29577951#;" decimal degrees of arc":PRINT 2500 REM Previous two lines are for debugging only. 2510 REM This next section is not in the original flowchart box E8 or any other flowchart box; however, let 2520 REM us consider this section, for all intents and purposes, to be in flowchart box E8 since its result 2530 REM (i.e. time of flight) is dependent on angle of elevation and horizontal or vertical separation (these 2540 REM last depending on the gun's angle of elevation) via a direct algebraic formula. 2550 REM Calculate the shell's time of flight. 2560 IF ELANG>1.483529864# THEN GOTO 2620:REM That is, we must use a 'special technique' to calculate 2570 REM the shell's time of flight if the gun's angle of elevation is extremely high, i.e. between 2580 REM 85 and 90 degrees of arc. 2590 TOF2=HOSEP/A2:GOTO 2680:REM Time of flight for relatively low gun angles of elevation, that is, angles of 2600 REM elevation lower than 85 degrees of arc; the GOTO etc. command in the previous line of code gets 2610 REM the program to skip over the high-gun-angle-time-of-flight algorithm. 2620 REM Calculate time of flight for gun angles of elevation between 85 and 90 degrees of arc; this uses 2630 REM equation 254a from page 254 of notes. 2640 B1=A1/G:B2=(A1*A1)/(G*G):B3=(2*VESEP)/G:TOF2=B1-(SQR(B2-B3)):REM TOF2 is time of flight in seconds 2650 REM of time for high gun angles of elevation, that is, greater than 85 degrees of arc. 2660 REM End flowchart box E8 2670 REM Begin flowchart box E9 (begin/end marker) 2680 RETURN 2690 REM End flowchart box E9 Program SPACBAT.BAS 10 REM Begin flowchart box 1 20 CLEAR:CLS:PRINT "Space battleship fire control program" 30 PRINT "for IBM PC or similar computers" 40 PRINT "which do trigonometric calculations in radians" 50 PRINT "by Ian Douglas, year 2010" 60 PRINT "This program works for any" 70 PRINT "azimuthal bearing of target ship" 80 PRINT "both at beginning and end of" 90 PRINT "shell's time of flight." 100 PRINT "Also works for any" 110 PRINT "bearing of point to aim at," 120 PRINT "i.e. point 'C.'" 130 PRINT "X Y Z coordinates of own ship at start" 140 PRINT "of problem are 50 nautical miles by 50 nautical miles by 0." 150 PRINT "Own ship's true direction always assumed to" 160 PRINT "be a bearing of 0 degrees or 360 degrees." 170 REM This program assumes that it is working in deep space, where the shell flies straight at the target 180 REM under no gravitational influences. 190 PRINT "This program works for targets whose starting and/or ending at elevations are above or below own ship;" 200 PRINT "likewise for the shell's point of impact. If gun must be depressed to hit target, then gun's" 210 PRINT "angle of elevation is negative." 220 REM End flowchart box 1 230 REM Begin flowchart box 2 240 INPUT "Muzzle velocity in feet/sec ";MV 250 INPUT "Target spaceship's present straight line-of-sight range in feet ";TSR 260 INPUT "Target spaceship's true 'ground track' speed in feet/sec ";TSS 270 PRINT "Target spaceship's modified true plot 'ground track' direction of travel as an" 280 INPUT "azimuthal bearing in decimal degrees of arc "; TSD:TSD=TSD/57.29577951# 290 PRINT "Target spaceship's first observed angle of elevation above horizontal" 300 INPUT "in decimal degrees of arc ";TSAE:TSAE=TSAE/57.29577951# 310 PRINT "Target spaceship's rate of climb in feet/sec" 320 INPUT "(descent negative) ";TSRC 330 INPUT "Own ship's true speed in feet/sec";OSS 340 PRINT "Target spaceship's bearing off own ship" 350 INPUT "at time of gun's firing in decimal degrees of arc ";TSB:TSB=TSB/57.29577951# 360 REM End flowchart box 2 370 REM Begin flowchart box 3 380 REM Calculate target aircraft's flat range (straight line of sight range reduced to XY plane) and then 390 REM calculate beginning (i.e. initially sighted) coordinates in (feet, feet) 400 FLRANG=TSR*COS(TSAE) 410 TSPX=(FLRANG*SIN(TSB))+303805! 420 TSPY=(FLRANG*COS(TSB))+303805! 430 TSPZ=TSR*SIN(TSAE) 440 REM Above addition of 303805 is number of feet in 50 nautical mile; 450 REM recall that own ship's position is (50 naut. mile, 50 naut. mile, 0) 460 PRINT "Target spaceship's starting (i.e. initially sighted) XYZ position in (feet, feet, feet): " 470 PRINT TSPX,TSPY,TSPZ 480 REM Previous 3 lines are for debugging and error isolation only; delete or rem them out later 490 TOF1=TSR/MV 500 REM TOF1 is first approximation of shell's time of flight 510 PRINT "First approximation of shell's time of flight in seconds:" 520 PRINT TOF1 530 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 540 REM End flowchart box 3 550 REM Begin flowchart box 4 560 BEEP:REM Gets computer to beep to show that the program is looping back properly 570 REM on first (presumably wrong) guess of time of flight; need to use 580 REM BEEP since I cannot get any of my computers to 'listen' to the 590 REM scroll lock key on the keyboard(s). BEEP statement also delays operation 600 REM of program so you can see on screen whether proper reiterative looping 610 REM to get progressively better answers is taking place. 620 REM First step in this flowchart box is to get distance target spaceship 630 REM travelled during TOF (time of flight) calculated above. 640 DTST=TSS*TOF1:REM DTST is 'ground track' distance target spaceship travels during TOF1 650 REM (time of flight) calculated above 660 REM Now calculate how far own ship travelled during time of shell's flight 670 DOST=OSS*TOF1:REM DOST is 'ground track'distance own ship travelled during shell's 680 REM time of flight 690 REM Next step is to use section B of coordinate inverse program pp. 75-76 700 REM to calculate the end coordinates of the target spaceship at end of time 710 REM of flight, using azimuth of target spaceship, distance travelled, and 720 REM target spaceship's starting coordinates. 730 REM Input is target spaceship's initial sighted coordinates (TSPX, TSPY, TSPZ), 740 REM target spaceship's modified true plot direction of travel azimuthal bearing 750 REM TSD and distance target spaceship travels during shell's time of flight DTST. 760 REM Output is target spaceship's endpoint coordinates (TSEPX, TSEPY, TSEPZ). 770 REM Calculate end coordinates 780 REM of vector, from angle, distance 790 REM and starting coordinates 800 TSEPX=TSPX+(DTST*SIN(TSD)):TSEPY=TSPY+(DTST*COS(TSD)):TSEPZ=TSPZ+(TSRC*TOF1) 810 REM End calculation of target spaceship's end coordinates 820 REM Begin calculation of coordinates of point "C" (PCX, PCY, PCZ) 830 REM X coordinate of point "C" is simply x coordinate of target spaceship's end point 840 REM at end of shell's time of flight, i.e. PCX=TSEPX 850 PCX=TSEPX 860 REM Y coordinate of point "C" is the y coordinate of target spaceship's endpoint at end of time of flight of shell TSEPY 870 REM minus distance own ship travels during shell's time of flight 880 PCY=TSEPY-DOST 890 REM Z coordinate of point "C" is simply target spaceship's altitude at end of shell's time of flight, i.e. TSEPZ 900 PCZ=TSEPZ 910 REM End calculation of coordinates of point "C" 920 PRINT "'Ground track' distance target spaceship travels during shell's time of flight (feet):" 930 PRINT DTST 940 PRINT "'Ground track' distance own ship travels during shell's time of flight (feet):" 950 PRINT DOST 960 PRINT "Target spaceship's XYZ end coordinates (feet, feet, feet):" 970 PRINT TSEPX,TSEPY, TSEPZ 980 PRINT "XYZ coordinates of point 'C' (feet, feet, feet):" 990 PRINT PCX,PCY,PCZ 1000 REM Previous 8 lines are for debugging and error isolation only; delete or rem them out later 1010 REM End flowchart box 4 1020 REM Begin flowchart box 5 1030 REM Calculate straight line distance between own ship at time of gun's firing 1040 REM and point "C," DOSPC 1050 REM Recall that own ship's coordinates at time of gun's firing are 1060 REM (303805 feet, 303805 feet, 0 feet) (i.e. 50 nautical miles by 50 nautical miles by 0) 1070 DOSPC=SQR(((PCX-303805!)^2)+((PCY-303805!)^2)+(PCZ^2)) 1080 PRINT "Distance between own ship and point 'C' when gun is fired (feet):" 1090 PRINT DOSPC 1100 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 1110 REM End flowchart box 5 1120 REM Begin flowchart box 6 1130 TOF2=DOSPC/MV 1140 PRINT "New time of flight (seconds of time):" 1150 PRINT TOF2 1160 REM Previous 2 lines are for debugging and error isolation only; delete or rem them out later 1170 REM End flowchart box 6 1180 REM Begin flowchart box 7 1190 IF (ABS(TOF2-TOF1))<.001 THEN GOTO 1280 1200 REM Note that the above line of code gives 1/1000 of a second accuracy!! I have found through 1210 REM experimentation that seeking higher accuracy will in some cases cause 1220 REM infinite looping of the program in a search for unnecessarily high accuracy. 1230 REM It is debatable whether in this program you need even 1/100 of a second accuracy. 1240 REM End flowchart box 7 1250 REM Begin flowchart box 8 1260 TOF1=TOF2:GOTO 550 1270 REM End flowchart box 8 1280 REM Begin flowchart box 9 1290 REM End flowchart box 9 1300 REM Begin flowchart box 10 1310 REM Calculate gun's bearing (not angle of elevation yet) using only X Y coordinates of point "C" 1320 REM From flowchart box 4, and using coordinates of own ship at time of gun's firing; 1330 REM Convert points to vector 1340 REM 1st point is coordinates of own ship at gun's time of firing (303805 feet, 303805 feet) 1350 REM 2nd point is coordinates of point "C," (PCX, PCY) 1360 DXAA=PCX-303805!:DYAA=PCY-303805! 1370 IF DXAA=0 THEN GOTO 1430 1380 REM Above line prevents division by zero in next line 1390 SLOPE=DYAA/DXAA:REM Slope 1400 ANGLE#=ATN(SLOPE):REM Angle referred to standard position, calculated from slope 1410 IF DXAA<0 THEN ANGLE#=ANGLE#+3.141592654# 1420 IF DXAA>0 AND DYAA<0 THEN ANGLE#=ANGLE#+6.283185307# 1430 IF DXAA=0 AND DYAA=0 THEN ANGLE#=0 1440 IF DXAA=0 AND DYAA>0 THEN ANGLE#=1.570796327# 1450 IF DXAA<0 AND DYAA=0 THEN ANGLE#=3.141592654# 1460 IF DXAA=0 AND DYAA<0 THEN ANGLE#=4.712388981# 1470 IF DXAA>0 AND DYAA=0 THEN ANGLE#=6.283185307# 1480 GUNAZ#=6.283185307#-ANGLE#+1.570796327#:REM 360 degrees minus angle plus 90 degrees 1490 IF GUNAZ#>6.283185307# THEN GUNAZ#=GUNAZ#-6.283185307# 1500 REM Calculate distance between own ship and 'ground track' position of point 'C,' that is, 1510 REM 'top view' position of point 'C' 1520 DOCG=SQR((DXAA^2)+(DYAA^2)):REM Above distance is calculated via Pythagoras 1530 REM Calculate gun's angle of elevation using arctan of (PCZ/DOCG) 1540 THETA1=ATN(PCZ/DOCG) 1550 REM THETA1 is the gun's angle of elevation such that the gun points directly at point "C" 1560 REM when it is fired 1570 REM End flowchart box 10 1580 REM Begin flowchart box 11 1590 REM Gun's angle of elevation is previously calculated THETA1 1600 PRINT "Gun's azimuth to shoot along: ";GUNAZ#*57.29577951#;" decimal degrees of arc" 1610 PRINT "Gun's angle of elevation: "; THETA1*57.29577951#;" decimal degrees of arc" 1620 REM If latest value for time of flight is to be printed here, that value is TOF2. 1630 IF TOF2<301 THEN GOTO 1680 1640 PRINT:PRINT "Time of shell's flight is greater than 5 minutes; target spaceship" 1650 PRINT "is effectively out of range.":PRINT 1660 REM End flowchart box 11 1670 REM Begin flowchart box 12 1680 END 1690 REM End flowchart box 12 Appendix, January 2012: Is an aircraft or airborne target within range of an antiaircraft gun? A surface-to-surface artillery piece is surrounded by a 2 dimensional "circle of death" whose radius is equivalent to the gun's surface (of the ground or sea) range; this circle's plane is also on the ground or sea surface, level with the gun. Analysis using basic trajectory equations will show that the space within range of an antiaircraft gun's upward-only and upward-to-peak trajectories is an "ellipsoid of death," that is, a (3 dimensional) bubble-like ellipsoid of revolution. The semimajor axis of this ellipsoid is horizontal, and equivalent to half the gun's ground/sea range. The minor axis is vertical, and is equivalent to the shell's maximum vertically-fired attainable height. The gun's location just touches the bottom of the ellipsoid, at the bottom of the minor axis; the major axis's closest point to the gun is above the gun by a distance equal to half the shell's maximum attainable vertically-fired height. The major axis is exactly twice as long as the minor axis, and the "bubble" space of the ellipsoid of death has a flattened appearance. All the above statements hold true given assumptions I have made elsewhere in this article, namely, that there is no air resistance, that the earth is flat, gravity is the only force other than the gun's initial impulse in effect, the earth is not rotating, that there is no wind (or even air for that matter), that there are no Magnus or Coriolis effects, etc. I have added some program listings to this article to determine whether a target is within range of an antiaircraft gun. The geometry of ellipses (related to these types of problems) is discussed (besides many other sources) in Krafft Ehricke's book "Space Flight Vol. I, Environment and Celestial Mechanics" (Princeton, N.J: D. Van Nostrand Company, Inc., 1960, pp. 293 et seq.). At a future date, I hope to integrate such a program with UPART.BAS and AIRBAT.BAS. As the reader may have noticed by now, my programs UPART.BAS and AIRBAT.BAS do not indicate or trap situations where the target aircraft is out of range; in my programs at present, these situations may cause infinite looping or an alarm-like exit to the operating system. Below (program IOELLIP.BAS) is a program to determine whether a 2 dimensional point is inside or outside of a 2D ellipse whose major axis is on the x axis. Part of the algorithm used here is from Ehricke's book mentioned above. This program would be useful for problems where the target aircraft is flying on a straight ground track, this ground track passing directly over the site of the gun. Regarding coordinate systems, there are two in use here: the gun's coordinate system (implied in this situation but not used in the program below) which has the location of the gun at the origin; and the ellipse's coordinate system, which has the centre of the ellipse at its origin. In a combat situation, the gun and ellipse centre are separated vertically; the ellipse centre and origin are above the gun location and origin by a distance equal to the ellipse's semi minor axis length. That semi minor axis length is in turn equal to half the shell's maximum attainable vertically-fired height. So ultimately we would have to convert between coordinate systems thus: YE = ellipse Y coordinate YG = gun Y coordinate B = ellipse semi minor axis length (i.e. half the shell's maximum attainable vertically-fired height) A gun Y coordinate is always positive in these types of problem. To convert from YE to YG: if YE is positive, then YG = ABS(B) + YE if YE is negative, then YG = ABS(B) - ABS(YE) (works only when gun Y coordinate is positive, which applies to our problem) [ ABS (expression/variable/number) means "the absolute value of (expression/variable/number) ] To convert from YG to YE: If YG < ABS(B) then YE = - {[ABS(B)] - YG} If YG > ABS(B) then YE = YG - ABS(B) Program IOELLIP.BAS 10 REM Begin flowchart box 1 20 CLS:CLEAR:PRINT "Determination of whether a given xy point in the 2D" 30 PRINT "Cartesian plane is inside or outside an ellipse" 40 PRINT "of given major and minor axis lengths" 50 PRINT "with the ellipse centre being at (0,0)." 60 PRINT 70 REM End flowchart box 1 80 REM Begin flowchart box 2 90 PRINT "Enter the data for the ellipse." 100 INPUT "Semi major axis length ";MAJA 110 INPUT "Semi minor axis length ";MINA 120 PRINT:PRINT "Enter the coordinates for the point whose general location (i.e. inside or outside" 130 PRINT "the ellipse) we wish to determine." 140 PRINT:INPUT "Point's x coordinate ";XB 150 INPUT "Point's y coordinate ";YB 160 REM End flowchart box 2 170 REM Begin flowchart box 3 180 XA=0:YA=0:REM Note that the starting coordinates XA and YA are at the centre of the ellipse 190 REM which is also the origin (0,0). 200 DX=XB-XA:DY=YB-YA 210 D=SQR((DX^2)+(DY^2)):REM Distance 220 IF DX=0 THEN GOTO 280 230 REM Above line prevents division by zero in next line 240 S=DY/DX:REM Slope 250 AN=ATN(S) 260 IF DX<0 THEN AN=AN+3.141592654# 270 IF DX>0 AND DY<0 THEN AN=AN+6.283185307# 280 IF DX=0 AND DY=0 THEN AN=0 290 IF DX=0 AND DY>0 THEN AN=1.570796327# 300 IF DX<0 AND DY=0 THEN AN=3.141592654# 310 IF DX=0 AND DY<0 THEN AN=4.712388981# 320 IF DX>0 AND DY=0 THEN AN=6.283185307# 330 REM End flowchart box 3 340 REM Begin flowchart box 4 350 REM Use Ehricke Vol. 1 p. 293 eq. 4-92. 360 REM Plug in central polar anomaly 'psi' 370 REM (which is the variable AN obtained in flowchart 380 REM box 3) to get the central polar radius vector 390 REM of the ellipse 'rho.' Above mentioned 400 REM central polar anomaly is the same for both 410 REM the ellipse and the given point. 420 A2=MAJA^2:B2=MINA^2:COS2P=(COS(AN))^2:TERM1=A2*B2:TERM2=(A2-B2)*COS2P 430 TERM3=A2-TERM2:REM Must keep this line and the next two non-remmed-out lines separate from one another 440 REM or program will not work! Logically they should work but the interpreter will not recognize this! 450 IF TERM3=0 THEN TERM3=.0000001 460 RHO=TERM1/TERM3:RHO=SQR(RHO) 470 REM Line before previous line prevents division by zero later in that same line to obtain 'rho'. 480 PRINT "Rho (i.e. central polar radius vector of" 490 PRINT "the ellipse for the same central" 500 PRINT "polar anomaly as the given point) ";RHO 510 PRINT "Distance of given point from centre of" 520 PRINT "ellipse ";D 530 PRINT "Psi (i.e. central polar anomaly of the point and of the ellipse intersection point of that line) " 540 PRINT AN;" radians":PRINT 550 REM End flowchart box 4 560 REM Begin flowchart box 5 (off page connector) 570 REM 580 REM End flowchart box 5 (off page connector) 590 REM Begin flowchart box 6 (off page connector) 600 REM 610 REM End flowchart box 6 (off page connector) 620 REM Begin flowchart box 7 630 IF RHO>D THEN PRINT "Given point is inside ellipse ":GOTO 660:REM Latter part of this line is flowchart box 8 640 PRINT "Given point is on or outside ellipse":REM This is flowchart box 9 650 REM End flowchart box 7 660 REM Begin flowchart box 10 670 END 680 REM End flowchart box 10 Here is a program to determine whether a 3 dimensional XYZ point is inside or outside an ellipsoid of revolution, sort of a "geodesic ellipsoid." The program below works with the left hand orthogonal coordinate system, but probably works equally well with the right hand orthogonal system. The coordinate system being used must have Z being vertical, with the semi minor axis being on Z. The ellipsoid is like two frisbees glued together by their bottom surfaces, clamshell-like or flying-saucer-like; its major axes are on X and Y, and are both the same length. The ellipsoid used is not a "rugby football" or "torpedo" ellipsoid, with differing lengths for each ellipsoid axis. Program IOELOID.BAS 10 CLS:CLEAR:PRINT "Program to determine whether a given XYZ point" 20 PRINT "is inside or outside an ellipsoid with semi minor axis being on the Z axis" 30 PRINT "and semi major axes, both of same length, being on X and Y axes." 40 PRINT "Note that this is a 'clamshell flying saucer' ellipsoid and not a 'torpedo' or 'rugby football' ellipsoid." 50 PRINT "Program calculates central polar latitude and radius vector" 60 PRINT "given a point's XYZ coordinates and lengths of" 70 PRINT "major and minor axes. Coordinate system" 80 PRINT "used here is right handed, sort of a 'geodesic ellipsoid' with" 90 PRINT "'north pole' being on positive Z axis" 100 PRINT "and X and Y being 'equatorial radii.'" 110 REM For in-or-out-of-flat-ellipse calculations, see Ehricke, "Space Flight Vol. 1", page 293, eq. 4-92 120 PRINT:PRINT "Enter the given point's XYZ coordinates" 130 PRINT:INPUT "X = ";GPX 140 INPUT "Y = ";GPY 150 INPUT "Z = ";GPZ 160 PRINT:PRINT "Enter the ellipsoid's semi major and semi minor axis lengths." 170 PRINT:INPUT "Semi major axis length = ";SMAJA 180 INPUT "Semi minor axis length = ";SMINA 190 REM Calculate the non-reduced non-projected horizontal distance from the origin to the 200 REM given point (i.e. to the location of the given point as seen projected onto the XY plane as seen from 210 REM well above the ellipsoid's "north pole" (the "north pole" being the Z axis). 220 HDIST=SQR((GPX^2)+(GPY^2)) 230 REM Calculate the straight line distance from 240 REM the origin to the given point 250 SLD=SQR((GPX^2)+(GPY^2)+(GPZ^2)) 260 REM Calculate the central polar non-reduced non-projected anomaly "psi" from the XY plane 270 REM to the given point; this is not a full circle arctan situation, as the denominator argument 280 REM is always positive. 290 DENOM1=SQR((GPX^2)+(GPY^2)):IF DENOM1=0 THEN DENOM1=.0000001:REM Last statement prevents division by zero in next line. 300 PSI=ATN(GPZ/DENOM1) 310 REM Calculate central polar radius vector of ellipsoid 320 TERM1=SMAJA^2:TERM2=SMINA^2:TERM3=TERM1*TERM2:TERM4=TERM1-TERM2 330 TERM5=(COS(PSI))^2:TERM6=TERM4*TERM5 340 TERM7=TERM1-TERM6:REM Must keep this line and the next two non-remmed-out lines separate from one another 350 REM or program will not work! Logically they should work but the interpreter will not recognize this! 360 IF TERM7=0 THEN TERM7=.0000001 370 TERM8=TERM3/TERM7:RHO=SQR(TERM8) 380 REM Line before previous line prevents division by zero later in that same line. 390 PRINT:PRINT "Straight line distance from origin to given point: " 400 PRINT SLD 410 PRINT:PRINT "Central polar radius vector from origin to ellipsoid:" 420 PRINT RHO:PRINT 430 IF RHO>SLD THEN PRINT "Given point is inside ellipsoid":GOTO 450 440 PRINT "Given point is on or outside ellipsoid" 450 PRINT:END The following two programs (GUNEOD1.BAS and EIDGUN2.BAS) deal with the problem of taking a stationary airborne target's antiaircraft-gun-centred coordinates and converting them to ellipsoid-of-death-centred coordinates; then the programs determine whether the airborne target is within range of the antiaircraft gun. Upward-only trajectories are considered, not upward-to-peak-to-downward trajectories. As before, the only forces assumed to be acting on the shell are the gun's initial impulse and gravity. Gun-centred coordinates have the gun's location at the origin; ellipse of death centred coordinates, or ellipsoid of death centred coordinates have the centre of the ellipse or ellipsoid at the origin. GUNEOD1.BAS solves the problem in 2 dimensions, where x is horizontal and y is vertical; EIDGUN2.BAS solves the problem in 3 dimensions, where x and y are horizontal and z is vertical. In the 2D version of the problem, x is the same for both gun-centred and ellipse of death centred coordinates; the only conversion that needs to be made is in y (vertical). In the 3D version of the problem, x and y are the same for both gun-centred and ellipse of death centred coordinates; the only conversion that needs to be made is in z (vertical). Program GUNEOD1.BAS 10 REM Begin flowchart box 1 20 CLS:CLEAR:PRINT "Determination of whether a target's xy location in the 2D" 30 PRINT "Cartesian plane is inside or outside an AA gun's ellipse of death (EOD);" 40 PRINT "the 2D plane has its y axis as vertical. Gun-centred coordinates used below have the location" 50 PRINT "of the gun at the origin; EOD-centred coordinates used below have the centre of the ellipse at" 60 PRINT "the origin." 70 PRINT 80 REM End flowchart box 1 90 REM Begin flowchart box 2 100 PRINT "Enter the data for the ellipse." 110 INPUT "Semi major axis length ";MAJA 120 INPUT "Semi minor axis length ";MINA 130 PRINT:PRINT "Enter the coordinates for the target's point whose general location (i.e. inside or outside" 140 PRINT "the ellipse) we wish to determine." 150 PRINT:INPUT "Target's gun-centred x coordinate ";XBG:XB=XBG:REM Last statement is simply intended to keep program's 160 REM and logic easier to follow than would be otherwise. 170 INPUT "Target's gun-centred y coordinate ";YBG 180 REM Convert y (vertical) coordinate from gun centred coordinate to ellipse-of-death-centred coordinate 190 IF YBG=ABS(MINA) THEN YB=YBG-ABS(MINA) 210 IF YBG=0 THEN YB=-(ABS(MINA)) 220 REM XB, YB obtained above are EOD-centred x, y coordinates. 230 PRINT "Target's EOD-centred x coordinate ";XB 240 PRINT "Target's EOD-centred y coordinate ";YB 250 REM End flowchart box 2 260 REM Begin flowchart box 3 270 XA=0:YA=0:REM Note that the starting coordinates XA and YA are at the centre of the ellipse 280 REM which is also the origin (0,0). 290 DX=XB-XA:DY=YB-YA 300 D=SQR((DX^2)+(DY^2)):REM Distance 310 IF DX=0 THEN GOTO 370 320 REM Above line prevents division by zero in next line 330 S=DY/DX:REM Slope 340 AN=ATN(S) 350 IF DX<0 THEN AN=AN+3.141592654# 360 IF DX>0 AND DY<0 THEN AN=AN+6.283185307# 370 IF DX=0 AND DY=0 THEN AN=0 380 IF DX=0 AND DY>0 THEN AN=1.570796327# 390 IF DX<0 AND DY=0 THEN AN=3.141592654# 400 IF DX=0 AND DY<0 THEN AN=4.712388981# 410 IF DX>0 AND DY=0 THEN AN=6.283185307# 420 REM End flowchart box 3 430 REM Begin flowchart box 4 440 REM Use Ehricke Vol. 1 p. 293 eq. 4-92. 450 REM Plug in central polar anomaly 'psi' 460 REM (which is the variable AN obtained in flowchart 470 REM box 3) to get the central polar radius vector 480 REM of the ellipse 'rho.' Above mentioned 490 REM central polar anomaly is the same for both 500 REM the ellipse and the target. 510 A2=MAJA^2:B2=MINA^2:COS2P=(COS(AN))^2:TERM1=A2*B2:TERM2=(A2-B2)*COS2P 520 TERM3=A2-TERM2:REM Must keep this line and the next two non-remmed-out lines separate from one another 530 REM or program will not work! Logically they should work but the interpreter will not recognize this! 540 IF TERM3=0 THEN TERM3=.0000001 550 RHO=TERM1/TERM3:RHO=SQR(RHO) 560 REM Line before previous line prevents division by zero later in that same line to obtain 'rho'. 570 PRINT "Rho (i.e. central polar radius vector of" 580 PRINT "the ellipse for the same central" 590 PRINT "polar anomaly as the target) ";RHO 600 PRINT "Distance of target from centre of" 610 PRINT "ellipse ";D 620 PRINT "Psi (i.e. central polar anomaly of the point and of the ellipse intersection point of that line) " 630 PRINT AN;" radians":PRINT 640 REM End flowchart box 4 650 REM Begin flowchart box 5 (off page connector) 660 REM 670 REM End flowchart box 5 (off page connector) 680 REM Begin flowchart box 6 (off page connector) 690 REM 700 REM End flowchart box 6 (off page connector) 710 REM Begin flowchart box 7 720 IF RHO>D THEN PRINT "Target is inside ellipse of death":GOTO 750:REM Latter part of this line is flowchart box 8 730 PRINT "Target is on or outside ellipse of death":REM This is flowchart box 9 740 REM End flowchart box 7 750 REM Begin flowchart box 10 760 END 770 REM End flowchart box 10 Program EIDGUN2.BAS 10 CLS:CLEAR:PRINT "Program to determine whether a target's XYZ is inside or outside an" 20 PRINT "AA gun's ellipsoid of death (EOD) with semi minor axis being on the Z axis" 30 PRINT "and semi major axes, both of same length, being on X and Y axes." 40 PRINT "Note that this is a 'clamshell flying saucer' ellipsoid and not a 'torpedo' or 'rugby football' ellipsoid." 50 PRINT "Program calculates central polar latitude and radius vector" 60 PRINT "given a target's XYZ coordinates and lengths of major and minor axes. Coordinate system" 70 PRINT "used here is right handed, sort of a 'geodesic ellipsoid' with 'north pole' being on positive Z axis" 80 PRINT "and X and Y being 'equatorial radii.' Gun-centred coordinates used below have the location" 90 PRINT "of the gun at the origin; EOD-centred coordinates used below have the centre of the EOD" 100 PRINT "at the origin." 110 REM For in-or-out-of-flat-ellipse calculations, see Ehricke, "Space Flight Vol. 1", page 293, eq. 4-92 120 PRINT:PRINT "Enter the target's gun-centred XYZ coordinates" 130 PRINT:INPUT "X = ";GPGCX 140 INPUT "Y = ";GPGCY 150 INPUT "Z = ";GPGCZ 160 GPX=GPGCX:GPY=GPGCY:REM Last statement is simply intended to keep program's logic easier to follow 170 REM than would be otherwise. 180 PRINT:PRINT "Enter the ellipsoid's semi major and semi minor axis lengths." 190 PRINT:INPUT "Semi major axis length = ";SMAJA 200 INPUT "Semi minor axis length = ";SMINA 210 REM Convert z (vertical) coordinate from gun centred coordinate to EOD-centred coordinate 220 IF GPGCZ=ABS(SMINA) THEN GPZ=GPGCZ-ABS(SMINA) 240 IF GPGCZ=0 THEN GPZ=-(ABS(SMINA)) 250 REM GPX, GPY, GPZ obtained above are EOD-centred x, y, z coordinates. 260 PRINT:PRINT "Target's EOD-centred x, y, z coordinates:" 270 PRINT:PRINT "X = ";GPX 280 PRINT "Y = ";GPY 290 PRINT "Z = ";GPZ 300 REM Calculate the non-reduced non-projected horizontal distance from the origin to the 310 REM target (i.e. to the location of the target as seen projected onto the XY plane as seen from 320 REM well above the ellipsoid's "north pole") (the "north pole" being the Z axis). 330 HDIST=SQR((GPX^2)+(GPY^2)) 340 REM Calculate the straight line distance from 350 REM the origin to the target 360 SLD=SQR((GPX^2)+(GPY^2)+(GPZ^2)) 370 REM Calculate the central polar non-reduced non-projected anomaly "psi" from the XY plane 380 REM to the target; this is not a full circle arctan situation, as the denominator argument 390 REM is always zero or positive. 400 DENOM1=SQR((GPX^2)+(GPY^2)) 410 IF DENOM1=0 THEN DENOM1=.0000001:REM Last statement prevents division by zero in next line. 420 PSI=ATN(GPZ/DENOM1) 430 REM Calculate central polar radius vector of ellipsoid 440 TERM1=SMAJA^2:TERM2=SMINA^2:TERM3=TERM1*TERM2:TERM4=TERM1-TERM2 450 TERM5=(COS(PSI))^2:TERM6=TERM4*TERM5 460 TERM7=TERM1-TERM6:REM Must keep this line and the next two non-remmed-out lines separate from one another 470 REM or program will not work! Logically they should work but the interpreter will not recognize this! 480 IF TERM7=0 THEN TERM7=.0000001 490 TERM8=TERM3/TERM7:RHO=SQR(TERM8) 500 REM Line before previous line prevents division by zero later in that same line. 510 PRINT:PRINT "Straight line distance from origin to target: " 520 PRINT SLD 530 PRINT:PRINT "Central polar radius vector from origin to ellipsoid:" 540 PRINT RHO:PRINT 550 IF RHO>SLD THEN PRINT "Target is inside ellipsoid of death":GOTO 570 560 PRINT "Target is on or outside ellipsoid of death" 570 PRINT:END Program UPARTR.BAS :This program expands on program UPART.BAS so as to trap instances where the target is outside the EOD and therefore out of range. Also note that this program has a weakness also present in UPART.BAS, which is that it goes into an infinite loop if the fire and hit solution involves the gun having an angle of elevation of greater than about 78 degrees: see the introductory notes for the program UPART.BAS above for a fuller explanation. This program appears to solve the problem in UPART.BAS of the program simply crashing, or looping forever, when the target is out of range. 10 REM Begin flowchart box E1 (begin/end marker) 20 REM End flowchart box E1; the flowchart numbers in this program correspond to flowchart on p. 251 of notes. 30 REM Begin flowchart box E2 40 PRINT "Program to use a variant of successive approximation to get angle of elevation" 50 PRINT "of a gun firing at a target separated from the gun significantly" 60 PRINT "separated from it both vertically and horizontally (e.g. an aircraft) using equation" 70 PRINT "211a on page 211 of notes; by Ian Douglas, year 2010. This program also works" 80 PRINT "for targets below gun. This program is for upward-only trajectories if target is above gun" 90 PRINT "and not for upward-to-peak-to-downward trajectories." 100 REM Program is intended to make the problem manageable by taking the maximum possible shell height the gun 110 REM can attain and obtaining a constant to multiply 90 degrees of arc (1.570796327 radians) (gun's maximum angle of 120 REM elevation) by to get any given possible height the shell can attain; this enables a resultant vertical 130 REM separation of the target from gun to be compared to an actual vertical separation attained by a trial angle 140 REM of elevation of gun, for comparison purposes; also this enables the successive approximation process to occur. 150 REM For this reason, I do not recommend using this program to compute data concerning guns with greater 160 REM than 10,000 feet per second muzzle velocities; but even the most powerful guns have about 1/4 that. 170 REM Both the program and its subroutine GUNEOD2 assume that the problem is in a 2 dimensional plane, with 180 REM the plane's horizontal portion being on the line on the (perfectly flat and level) ground separating 190 REM the target from the gun horizontally, and the vertical portion of the plane being the vertical separation 200 REM of the gun from the target. 210 REM Also note that the program will go into an infinite loop if the target is directly above gun, even though 220 REM a fire and hit solution is possible in such a case. 230 PRINT "Vertical separation (Y) of target from gun (feet) " 240 PRINT "(negative for a target below gun)" 250 INPUT VESEP 260 INPUT "Horizontal separation (in XY plane) of target from gun (feet) ";HOSEP 270 INPUT "Gun's muzzle velocity (feet per second) ";MV 280 G=32:G1=G/2:REM Acceleration due to gravity in feet per second per second; also half this for use in 290 REM instances of equation 211a. 300 REM End flowchart box E2 310 REM Begin flowchart box E3 320 REM Calculate the shell's maximum attainable height. 330 MAXHT=(MV*MV)/(2*G) 340 REM Calculate the semi major axis length of the Ellipse of Death (EOD); this is half the gun's surface-to- 350 REM surface range. 360 MAJA=((MV*MV)/G)/2 370 REM Calculate the semi minor axis length of the Ellipse of Death (EOD); this is half the shell's maximum 380 REM attainable height, i.e. the maximum height attained by the shell when fired straight up. 390 MINA=MAXHT/2:GOSUB 970 400 REM Calculate the constant mentioned in remarks for flowchart box E2. 410 CONST=MAXHT/1.570796327# 420 REM Calculate the minimum angle of elevation of gun needed to reach HOSEP. 430 MIELAN=(HOSEP*G)/(MV*MV) 440 REM The following is an arc sine routine with the final angle divided by 2 for the purpose of the 450 REM above formula to get the minimum necessary angle of elevation given HOSEP, muzzle velocity, and G. 460 MIELAN=MIELAN+((MIELAN^3)/6)+(.075*(MIELAN^5)) 470 MIELAN=MIELAN+(4.464285714285714D-02*(MIELAN^7)) 480 MIELAN=MIELAN+(3.038194444444444D-02*(MIELAN^9)) 490 MIELAN=MIELAN+(2.237215909090909D-02*(MIELAN^11)) 500 MIELAN=MIELAN+(1.735276442307692D-02*(MIELAN^13)) 510 MIELAN=MIELAN+(.01396484375#*(MIELAN^15)) 520 MIELAN=MIELAN+(1.155180089613971D-02*(MIELAN^17)) 530 MIELAN=MIELAN+(9.761609529194079D-03*(MIELAN^19)) 540 MIELAN=MIELAN+(8.390335809616816D-03*(MIELAN^21)) 550 MIELAN=MIELAN+(7.312525877035907D-03*(MIELAN^23)) 560 MIELAN=MIELAN+(6.447210311379962D-03*(MIELAN^25)) 570 MIELAN=MIELAN+(5.74003766956053D-03*(MIELAN^27)) 580 MIELAN=MIELAN+(.00515330968258#*(MIELAN^29)) 590 MIELAN=MIELAN/2 600 REM Angle of elevation is MIELAN; used formula on p. 12 of notes. 610 REM Set shell's minimum possible height as zero feet; set gun's maximum angle of elevation as 90 degrees of arc. 620 REM End flowchart box E3 630 REM Begin flowchart box E4 640 REM Set first trial angle of elevation of gun as (45 degrees of arc minus minimum angle of elevation of gun 650 REM to reach HOSEP)/2; the reason I use 45 and not 90 degrees here is that as soon as you use an angle greater 660 REM than 45 degrees you meet up with the possibility of the shell falling short of HOSEP. 670 ELANG=(.785398163#-MIELAN)/2 680 REM End flowchart box E4 690 REM Begin flowchart box E5 700 REM Calculate trial vertical (Y) separation resulting from ELANG and inputted values from flowchart box E2 710 REM (except VESEP) using equation 211a; calculate difference between trial vertical separation and VESEP. 720 A1=MV*SIN(ELANG):A2=MV*COS(ELANG) 730 TRIZ=(HOSEP*A1)/A2 740 TRIZ=TRIZ-G1*((HOSEP*HOSEP)/(A2*A2)) 750 REM PRINT "Trial Z and trial gun's angle of elevation are ";TRIZ,ELANG*57.29577951#:REM For debugging only 760 REM End flowchart box E5 770 REM Begin flowchart box E6 780 IF ABS(TRIZ-VESEP)<.01 THEN GOTO 870:REM That is, go to flowchart box E8. Note that this represents a 1/100 790 REM of a foot accuracy of the shell's striking point in relation to the program's perception of the target, 800 REM where the target is an aircraft! Compare this with the size of any conventional aircraft. 810 REM I have tried to get the program to achieve a higher level of precision, but this only caused 820 REM endless looping due to the problem of attainable accuracy of TRIZ and ELANG, loss of digits, etc. 830 REM End flowchart box E6 840 REM Begin flowchart box E7 850 ELANG=ELANG-(((TRIZ-VESEP)/2)/CONST):GOTO 690:REM That is, go to flowchart box E5 860 REM End flowchart box E7 870 REM Begin flowchart box E8 880 PRINT:PRINT "Gun's angle of elevation from horizontal (depression below horizontal negative) = "; 890 PRINT ELANG*57.29577951#;" decimal degrees of arc":PRINT 900 REM End flowchart box E8 910 REM Begin flowchart box E9 (begin/end marker) 920 END 930 REM End flowchart box E9 940 REM 950 REM * * * Subroutines * * * 960 REM 970 REM Subroutine GUNEOD2, derived from program GUNEOD1.BAS; flowchart box numbers in this subroutine 980 REM refer to notes about GUNEOD1.BAS. 990 REM Begin flowchart box 1 1000 REM Determination of whether a target's xy location in the 2D 1010 REM Cartesian plane is inside or outside an AA gun's ellipse of death (EOD); 1020 REM the 2D plane has its y axis as vertical. Gun-centred coordinates used below have the location 1030 REM of the gun at the origin; EOD-centred coordinates used below have the centre of the ellipse at 1040 REM the origin. Enter the subroutine with the EOD's semi major axis length MAJA, semi minor axis length MINA, 1050 REM horizontal separation of target from gun HOSEP, and vertical separation of target from gun VESEP. You 1060 REM leave the program knowing whether the target is within range or out of range. If the target is within 1070 REM range, control goes back to the main program. If the target is out of range, this subroutine says so 1080 REM and terminates the program. 1090 PRINT 1100 REM End flowchart box 1 1110 REM Begin flowchart box 2 1120 XBG=HOSEP:YBG=VESEP:XB=XBG:REM Last statement is simply intended to keep program's 1130 REM and logic easier to follow than would be otherwise. 1140 REM Convert y (vertical) coordinate from gun centred coordinate to ellipse-of-death-centred coordinate 1150 IF YBG=ABS(MINA) THEN YB=YBG-ABS(MINA) 1170 IF YBG=0 THEN YB=-(ABS(MINA)) 1180 REM XB, YB obtained above are EOD-centred x, y coordinates. 1190 PRINT "Target's EOD-centred x coordinate ";XB 1200 PRINT "Target's EOD-centred y coordinate ";YB 1210 REM End flowchart box 2 1220 REM Begin flowchart box 3 1230 XA=0:YA=0:REM Note that the starting coordinates XA and YA are at the centre of the ellipse 1240 REM which is also the origin (0,0). 1250 DX=XB-XA:DY=YB-YA 1260 D=SQR((DX^2)+(DY^2)):REM Distance 1270 IF DX=0 THEN GOTO 1330 1280 REM Above line prevents division by zero in next line 1290 S=DY/DX:REM Slope 1300 AN=ATN(S) 1310 IF DX<0 THEN AN=AN+3.141592654# 1320 IF DX>0 AND DY<0 THEN AN=AN+6.283185307# 1330 IF DX=0 AND DY=0 THEN AN=0 1340 IF DX=0 AND DY>0 THEN AN=1.570796327# 1350 IF DX<0 AND DY=0 THEN AN=3.141592654# 1360 IF DX=0 AND DY<0 THEN AN=4.712388981# 1370 IF DX>0 AND DY=0 THEN AN=6.283185307# 1380 REM End flowchart box 3 1390 REM Begin flowchart box 4 1400 REM Use Ehricke Vol. 1 p. 293 eq. 4-92. 1410 REM Plug in central polar anomaly 'psi' 1420 REM (which is the variable AN obtained in flowchart 1430 REM box 3) to get the central polar radius vector 1440 REM of the ellipse 'rho.' Above mentioned 1450 REM central polar anomaly is the same for both 1460 REM the ellipse and the target. 1470 A2=MAJA^2:B2=MINA^2:COS2P=(COS(AN))^2:TERM1=A2*B2:TERM2=(A2-B2)*COS2P 1480 TERM3=A2-TERM2:REM Must keep this line and the next two non-remmed-out lines separate from one another 1490 REM or program will not work! Logically they should work but the interpreter will not recognize this! 1500 IF TERM3=0 THEN TERM3=.0000001 1510 RHO=TERM1/TERM3:RHO=SQR(RHO) 1520 REM Line before previous line prevents division by zero later in that same line to obtain 'rho'. 1530 PRINT "Rho (i.e. central polar radius vector of" 1540 PRINT "the ellipse for the same central" 1550 PRINT "polar anomaly as the target) ";RHO 1560 PRINT "Distance of target from centre of" 1570 PRINT "ellipse ";D 1580 PRINT "Psi (i.e. central polar anomaly of the point and of the ellipse intersection point of that line) " 1590 PRINT AN;" radians":PRINT 1600 REM End flowchart box 4 1610 REM Begin flowchart box 5 (off page connector) 1620 REM 1630 REM End flowchart box 5 (off page connector) 1640 REM Begin flowchart box 6 (off page connector) 1650 REM 1660 REM End flowchart box 6 (off page connector) 1670 REM Begin flowchart box 7 1680 PRINT "EOD semimajor axis length = ";MAJA;" EOD semiminor axis = ";MINA 1690 IF RHO>D THEN PRINT "Target is inside ellipse of death and therefore within gun's range":RETURN 1700 REM Latter part of above line is similar to flowchart box 8 1710 PRINT "Target is on or outside ellipse of death and therefore outside of gun's range" 1720 REM Above line is flowchart box 9 1730 REM End flowchart box 7 1740 REM Begin flowchart box 10 1750 END 1760 REM End flowchart box 10