UTM投影坐标转WGS84大地坐标

public class UTMAndWGS84
	{
   
		static double pi = Math.PI;

		/* Ellipsoid model constants (actual values here are for WGS84) */
		static double sm_a = 6378137.0;
		static double sm_b = 6356752.314;
		static double sm_EccSquared = 6.69437999013e-03;

		static double UTMScaleFactor = 0.9996;


		/*
        * DegToRad
        *
        * Converts degrees to radians.
        *
        */
		private static double DegToRad(double deg)
		{
   
			return (deg / 180.0 * pi);
		}




		/*
        * RadToDeg
        *
        * Converts radians to degrees.
        *
        */
		private static double RadToDeg(double rad)
		{
   
			return (rad / pi * 180.0);
		}




		/*
        * ArcLengthOfMeridian
        *
        * Computes the ellipsoidal distance from the equator to a point at a
        * given latitude.
        *
        * Reference: Hoffmann-Wellenhof, B., Lichtenegger, H., and Collins, J.,
        * GPS: Theory and Practice, 3rd ed.  New York: Springer-Verlag Wien, 1994.
        *
        * Inputs:
        *     phi - Latitude of the point, in radians.
        *
        * Globals:
        *     sm_a - Ellipsoid model major axis.
        *     sm_b - Ellipsoid model minor axis.
        *
        * Returns:
        *     The ellipsoidal distance of the point from the equator, in meters.
        *
        */
		private static double ArcLengthOfMeridian(double phi)
		{
   
			double alpha, beta, gamma, delta, epsilon, n;
			double result;

			/* Precalculate n */
			n = (sm_a - sm_b) / (sm_a + sm_b);

			/* Precalculate alpha */
			alpha = ((sm_a + sm_b) / 2.0)
			   * (1.0 + (Math.Pow(n, 2.0) / 4.0) + (Math.Pow(n, 4.0) / 64.0));

			/* Precalculate beta */
			beta = (-3.0 * n / 2.0) + (9.0 * Math.Pow(n, 3.0) / 16.0)
			   + (-3.0 * Math.Pow(n, 5.0) / 32.0);

			/* Precalculate gamma */
			gamma = (15.0 * Math.Pow(n, 2.0) / 16.0)
				+ (-15.0 * Math.Pow(n, 4.0) / 32.0);

			/* Precalculate delta */
			delta = (-35.0 * Math.Pow(n, 3.0) / 48.0)
				+ (105.0 * Math.Pow(n, 5.0) / 256.0);

			/* Precalculate epsilon */
			epsilon = (315.0 * Math.Pow(n, 4.0) / 512.0);

			/* Now calculate the sum of the series and return */
			result = alpha
				* (phi + (beta * Math.Sin(2.0 * phi))
					+ (gamma * Math.Sin(4.0 * phi))
					+ (delta * Math.Sin(6.0 * phi))
					+ (epsilon * Math.Sin(8.0 * phi)));

			return result;
		}



		/*
        * UTMCentralMeridian
        *
        * Determines the central meridian for the given UTM zone.
        *
        * Inputs:
        *     zone - An integer value designating the UTM zone, range [1,60].
        *
        * Returns:
        *   The central meridian for the given UTM zone, in radians, or zero
        *   if the UTM zone parameter is outside the range [1,60].
        *   Range of the central meridian is the radian equivalent of [-177,+177].
        *
        */
		private static double UTMCentralMeridian(double zone)
		{
   
			double cmeridian;

			cmeridian = DegToRad(-183.0 +