Go Down

Topic: trying to get gnu Ballistics to run (Read 1 time) previous topic - next topic


I don't know if this is the right please to beg for help. I found a ballistic program in c online. I tried to get it to run on my uno but would not work. It compiles and runs but never gets to the end of the program.

I did some debugging and I think it gets stuck on the heavy math. Could you guys take a look at it and let me know off there is a fix for it?

I'm posting the header and main code in separate post.



Code: [Select]
// GNU Ballistics Library
// Originally created by Derek Yates
// Now available free under the GNU GPL

#ifndef __lib_ballistics__
#define __lib_ballistics__
#define __BCOMP_MAXRANGE__ 50001
#define GRAVITY (-32.194)
#include <math.h>
#include <stdlib.h>

enum __DragFunctions {G1=1,G2,G3,G4,G5,G6,G7,G8};

// Angular conversion functions to make things a little easier.
// Source is in _angle.c
double DegtoMOA(double deg); // Converts degrees to minutes of angle
double DegtoRad(double deg); // Converts degrees to radians
double MOAtoDeg(double moa); // Converts minutes of angle to degrees
double MOAtoRad(double moa); // Converts minutes of angle to radians
double RadtoDeg(double rad); // Converts radians to degrees
double RadtoMOA(double rad); // Converts radiants to minutes of angle

// A function to calculate ballistic retardation values based on standard drag functions.
// Source is in "_retard.c"
double retard(int DragFunction, double DragCoefficient, double Vi);
/* Arguments:
DragFunction:  G1, G2, G3, G4, G5, G6, G7, or G8.  All are enumerated above.
DragCoefficient:  The coefficient of drag for the projectile for the given drag function.
Vi:  The Velocity of the projectile.

Return Value:
The function returns the projectile drag retardation velocity, in ft/s per second.


// A function to correct a "standard" Drag Coefficient for differing atmospheric conditions.
// Returns the corrected drag coefficient for supplied drag coefficient and atmospheric conditions.
// Source is in "_atmosphere.c"
double AtmCorrect(double DragCoefficient, double Altitude, double Barometer, double Temperature, double RelativeHumidity);
/* Arguments:
DragCoefficient:  The coefficient of drag for a given projectile.
Altitude:  The altitude above sea level in feet.  Standard altitude is 0 feet above sea level.
Barometer:  The barometric pressure in inches of mercury (in Hg).
This is not "absolute" pressure, it is the "standardized" pressure reported in the papers and news.
Standard pressure is 29.53 in Hg.
Temperature:  The temperature in Fahrenheit.  Standard temperature is 59 degrees.
RelativeHumidity:  The relative humidity fraction.  Ranges from 0.00 to 1.00, with 0.50 being 50% relative humidity.
Standard humidity is 78%

Return Value:
The function returns a ballistic coefficient, corrected for the supplied atmospheric conditions.

// A function to compute the windage deflection for a given crosswind speed,
// given flight time in a vacuum, and given flight time in real life.
// Returns the windage correction needed in inches.
// Source is in "_windage.c"
double Windage(double WindSpeed, double Vi, double x, double t);
/* Arguments:
WindSpeed:  The wind velocity in mi/hr.
Vi:  The initial velocity of the projectile (muzzle velocity).
x:  The range at which you wish to determine windage, in feet.
t:  The time it has taken the projectile to traverse the range x, in seconds.

Return Value:
Returns the amount of windage correction, in inches, required to achieve zero on a target at the given range.


// Functions to resolve any wind / angle combination into headwind and crosswind components.
// Source is in "_windage.c"
double HeadWind(double WindSpeed, double WindAngle);
double CrossWind(double WindSpeed, double WindAngle);
/* Arguments:
WindSpeed:  The wind velocity, in mi/hr.
WindAngle:  The angle from which the wind is coming, in degrees.
0 degrees is from straight ahead
90 degrees is from right to left
180 degrees is from directly behind
270 or -90 degrees is from left to right.

Return value:
Returns the headwind or crosswind velocity component, in mi/hr.

// A function to determine the bore angle needed to achieve a target zero at Range yards
// (at standard conditions and on level ground.)
// Source is in "_zero.c"
double ZeroAngle(int DragFunction, double DragCoefficient, double Vi, double SightHeight, double ZeroRange, double yIntercept);
/*  Arguments:
DragFunction:  The drag function to use (G1, G2, G3, G5, G6, G7, G8)
DragCoefficient:  The coefficient of drag for the projectile, for the supplied drag function.
Vi:  The initial velocity of the projectile, in feet/s
SightHeight:  The height of the sighting system above the bore centerline, in inches.
 Most scopes fall in the 1.6 to 2.0 inch range.
ZeroRange:  The range in yards, at which you wish the projectile to intersect yIntercept.
yIntercept:  The height, in inches, you wish for the projectile to be when it crosses ZeroRange yards.
This is usually 0 for a target zero, but could be any number.  For example if you wish
to sight your rifle in 1.5" high at 100 yds, then you would set yIntercept to 1.5, and ZeroRange to 100

Return Value:
Returns the angle of the bore relative to the sighting system, in degrees.

// A function to generate a ballistic solution table in 1 yard increments, up to __BCOMP_MAXRANGE__.
// Source is in "_solve.c"
int SolveAll(int DragFunction, double DragCoefficient, double Vi, double SightHeight, \
double ShootingAngle, double ZeroAngle, double WindSpeed, double WindAngle, double** Solution);
/* Arguments:
DragFunction:  The drag function you wish to use for the solution (G1, G2, G3, G5, G6, G7, or G8)
DragCoefficient:  The coefficient of drag for the projectile you wish to model.
Vi:  The projectile initial velocity.
SightHeight:  The height of the sighting system above the bore centerline.  
Most scopes are in the 1.5"-2.0" range.
ShootingAngle:  The uphill or downhill shooting angle, in degrees.  Usually 0, but can be anything from
90 (directly up), to -90 (directly down).
ZeroAngle:  The angle of the sighting system relative to the bore, in degrees.  This can be easily computed
using the ZeroAngle() function documented above.
WindSpeed:  The wind velocity, in mi/hr
WindAngle:  The angle at which the wind is approaching from, in degrees.
0 degrees is a straight headwind
90 degrees is from right to left
180 degrees is a straight tailwind
-90 or 270 degrees is from left to right.
Solution: A pointer provided for accessing the solution after it has been generated.
Memory for this pointer is allocated in the function, so the user does not need
to worry about it.  This solution can be passed to the retrieval functions to get
useful data from the solution.
Return Value:
This function returns an integer representing the maximum valid range of the
solution.  This also indicates the maximum number of rows in the solution matrix,
and should not be exceeded in order to avoid a memory segmentation fault.


// Functions for retrieving data from a solution generated with SolveAll()
double GetRange(double* sln, int yardage); // Returns range, in yards.
double GetPath(double* sln, int yardage); // Returns projectile path, in inches, relative to the line of sight.
double GetMOA(double* sln, int yardage); // Returns an estimated elevation correction for achieving a zero at this range.
// this is useful for "click charts" and the like.
double GetTime(double* sln, int yardage);// Returns the projectile's time of flight to this range.
double GetWindage(double* sln, int yardage); // Returns the windage correction in inches required to achieve zero at this range.
double GetWindageMOA(double* sln, int yardage); // Returns an approximate windage correction in MOA to achieve a zero at this range.
double GetVelocity(double* sln, int yardage); // Returns the projectile's total velocity (Vector product of Vx and Vy)
double GetVx(double* sln, int yardage); // Returns the velocity of the projectile in the bore direction.
double GetVy(double* sln, int yardage); // Returns the velocity of the projectile perpendicular to the bore direction.

// For very steep shooting angles, Vx can actually become what you would think of as Vy relative to the ground,
// because Vx is referencing the bore's axis.  All computations are carried out relative to the bore's axis, and
// have very little to do with the ground's orientation.

//// I have split the source up into several files to make it easier for me to work on and refine.
//// Some people hate this, but I find it is easier to work on the code in small chunks.
//#include "_angle.c"// Angular conversion functions
//#include "_retard.c" // G-function retardations
//#include "_atmosphere.c" // Drag coefficient atmospheric corrections
//#include "_windage.c"// Crosswind correction functions
//#include "_zero.c" // Used to determine bore angle
//#include "_solve.c" // Used to solve ballistic problems
//#include "_retrieve.c" // Used to make retrieving solution data easier.


I narrowed it to this section of the sketch

It looks like it want to make a large array to hold all the results. At this point I just want to calculate a single point rather then a table of results.

Code: [Select]
int SolveAll(int DragFunction, double DragCoefficient, double Vi, double SightHeight, \
double ShootingAngle, double ZAngle, double WindSpeed, double WindAngle, double** Solution){

  double* ptr;
  ptr = (double*)malloc(10*__BCOMP_MAXRANGE__*sizeof(double)+2048);

  double t=0;
  double dt=0.5/Vi;
  double v=0;
  double vx=0, vx1=0, vy=0, vy1=0;
  double dv=0, dvx=0, dvy=0;
  double x=0, y=0;

  double headwind=HeadWind(WindSpeed, WindAngle);
  double crosswind=CrossWind(WindSpeed, WindAngle);

  double Gy=GRAVITY*cos(DegtoRad((ShootingAngle + ZAngle)));
  double Gx=GRAVITY*sin(DegtoRad((ShootingAngle + ZAngle)));



  int n=0;
  for (t=0;;t=t+dt){

    vx1=vx, vy1=vy;

    // Compute acceleration using the drag function retardation
    dv = retard(DragFunction,DragCoefficient,v+headwind);
    dvx = -(vx/v)*dv;
    dvy = -(vy/v)*dv;

    // Compute velocity, including the resolved gravity vectors.
    vx=vx + dt*dvx + dt*Gx;
    vy=vy + dt*dvy + dt*Gy;

    if (x/3>=n){
      ptr[10*n+0]=x/3; // Range in yds
      ptr[10*n+1]=y*12; // Path in inches
      ptr[10*n+2]=-RadtoMOA(atan(y/x)); // Correction in MOA
      ptr[10*n+3]=t+dt; // Time in s
      ptr[10*n+4]=Windage(crosswind,Vi,x,t+dt); // Windage in inches
      ptr[10*n+5]=RadtoMOA(atan(ptr[10*n+4])); // Windage in MOA
      ptr[10*n+6]=v; // Velocity (combined)
      ptr[10*n+7]=vx; // Velocity (x)
      ptr[10*n+8]=vy; // Velocity (y)
      ptr[10*n+9]=0; // Reserved

    // Compute position based on average velocity.

    if (fabs(vy)>fabs(3*vx)) break;
    if (n>=__BCOMP_MAXRANGE__+1) break;



  *Solution = ptr;
  return n;


So in the header file, we see that __BCOMP_MAXRANGE__ is 50001, so let's call that 50k.

Code: [Select]

ptr = (double*)malloc(10*__BCOMP_MAXRANGE__*sizeof(double)+2048);

A double is 8 bytes.

10 * 50,000 * 8 + 2048 = 4,002,048 bytes. That's about 4 meg RAM.

So this is not going to happen on any Arduino.

Go Up