/*
      REAL FUNCTION GRAV(P,LAT)
*/

#include <stdio.h>
#include <math.h>
#include "ocean.h"

#if PROTOTYPE_ALLOWED
double GRAV(double P, double LAT)
#else
double GRAV(P, LAT)
double P, LAT;
#endif
{
     double X, GR;
      X = sin(LAT/57.29578);
      X = X*X;
/*
     =====>  GR = Gravity variation with latitude:
             ANON (1970) Bulletin Geodesique
*/
      GR = 9.780318*(1.0+(5.2788E-3+2.36E-5*X)*X) + 1.092E-6*P;
     return(GR);
}

#ifdef EXE
#if PROTOTYPE_ALLOWED
void main(int argc, char *argv[])
#else
void main(argc, argv)
int argc;
char *argv[];
#endif
{
   double P = 10000.0, LAT = 30.0;

   if (argc != 3)
   {
      printf("\n USAGE: grav <pressure> <latitude>");
      printf("\n EXAMPLE: grav 10000.0 30.0\n");         
   }
   else
   {
      P = atof(argv[1]);
      LAT = atof(argv[2]);
   }
   printf("\n gravity = %f m/s^2 at P = %g dbar, %g degrees latitude\n\n",
      GRAV(P, LAT), P, LAT);
}
#endif
