GPS latitude and longitude converted to rotated plane coordinates

/** * Map Tool * * @author DXM * */public class Maputil {//reference point private double[] D34 = null;private double[] d33 = Null;priva  Te double[] gs34 = null;private double[] gs33 = null;//reference point location private double x34, y34;//calculation parameter private double A, B, C, bl;/** * Constructor * * @param ds34 reference point latitude and longitude * @param ds33 reference point latitude and longitude * @param length two reference points in the length of this coordinate system */public maputil (string[] ds34, string[] Ds33, double x34, double y34, double length) {this.x34 = X34;this.y34 = Y34;d34 = Convertgps (Ds34[0], ds34[1]);d = Conve Rtgps (Ds33[0], ds33[1]); gs34 = Bltogauss (D34[0], d34[1]), Gs33 = Bltogauss (D33[0], d33[1]);//Point 33 to 34 gauss transverse distance a = gs34[0]-G s33[0];//Point 33 to 34 gauss longitudinal distance b = gs34[1]-gs33[1];//point 33 to 34 Gauss straight line Distance c = math.sqrt (A * a + b * b);//Gaussian distance to the ratio of this coordinate system BL = LENGTH/C ;}  /** * According to the latitude and longitude of the XY * * @param LNG * @param lat * @return */public int[] Getxy (double LNG, double lat) {///compare with reference points longitude calculate coordinates if (LNG > D34[0]) {return getXY1 (LNG, LAT);} else {return getXY2 (LNG, LAT);}} /** * The fixed point is the baseline xy coordinate */private int[] getXY1 (double LNG, double lat) {//WillLatitude converted to Gaussian projection coordinates double[] GS = Bltogauss (LNG, LAT);//with point 34 transverse gauss distance Double d = math.abs (Gs[0]-gs34[0]);d ouble e = d * C/b;doub  Le f = a * e/c;//with a longitudinal gaussian distance of point 34 double GF = Math.Abs (gs[1]-gs34[1]);d ouble g = gf-f;double h = A * g/c;double i = b * g/c;double x = x34-i * Bl;double y = y34 + (H + e) * Bl;return new int[] {(int) math.rint (x), (int) math.rint (y)};} /** * The fixed point is the baseline xy coordinate */private int[] getXY2 (double LNG, double lat) {//convert latitude and longitude to Gaussian projection coordinates double[] GS = Bltogauss (LNG, LAT);//With Point 34 Transverse Gaussian distance double d = math.abs (Gs[0]-gs34[0]);d ouble e = d * c/a;double f = e * b/c;//with point 34 longitudinal gauss distance Double GF = Math.Abs  (Gs[1]-gs34[1]);d ouble g = gf-f;double h = A * g/c;double i = b * g/c;double x = x34-(e + i) * bl;double y = y34 + H * Bl;return new int[] {(int) math.rint (x), (int) math.rint (y)}; /** * latitude and longitude degrees seconds of the decimal point * * @param JW * @return */private double[] Convertgps (string LNG, string lat) {LNG = Lng.replaceall (" N "," "); LNG = Lng.replaceall (" S "," "); LNG = Lng.replaceall (" E "," "); LNG = lNg.replaceall ("W", ""); lat = Lat.replaceall ("N", ""); lat = Lat.replaceall ("S", ""); lat = Lat.replaceall ("E", ""); lat = lat . ReplaceAll ("W", "");d ouble lngd = 0; string[] fs = Lng.split ("°"); if (Fs.length > 0) {lngd = Double.parsedouble (Fs[0]);} FS = Fs[1].split ("" "), if (Fs.length > 0) {Double fen = double.parsedouble (fs[0]); fen = fen/60;lngd + = Fen;} FS = Fs[1].split ("″"), if (Fs.length > 0) {Double mi = double.parsedouble (fs[0]), MI = mi/(*); Lngd + = mi;} Double latd = 0;fs = Lat.split ("°"); if (Fs.length > 0) {latd = Double.parsedouble (Fs[0]);} FS = Fs[1].split ("" "), if (Fs.length > 0) {Double fen = double.parsedouble (fs[0]); fen = fen/60;latd + = Fen;} FS = Fs[1].split ("″"), if (Fs.length > 0) {Double mi = double.parsedouble (fs[0]), MI = mi/(*); latd + = mi;} return new double[] {lngd, latd};} /** * is calculated from the Gaussian projection coordinates into the latitude * * @param X * @param Y * @return */@SuppressWarnings ("unused") private double[] Gausstobl (double X, D Ouble Y) {int projno;int zonewide;////bandwidth double[] Output = new double[2];d ouble longitude1, Latitude1, Longitude0, X0, Y0, Xval, yval;//latitude0,double E1, E2, F, a, EE, NN, T, C, M, D, R, U, FAI, Ipi;ipi = 0.0174532925199433; 3.1415926535898/180.0;//a = 6378245.0; f = 1.0/298.3; 54 Beijing coordinate system parameter A = 6378140.0;f = 1/298.257; 80 XI ' An coordinate system parameter Zonewide = 6; 6-degree bandwidth Projno = (int) (x/1000000l); Look for the band number longitude0 = (ProjNo-1) * zonewide + zonewide/2;longitude0 = longitude0 * IPI; central Meridian X0 = projno * 1000000L + 500000L; Y0 = 0;xval = X-x0;yval = Y-y0; In-band geodetic coordinates e2 = 2 * f-f * f;e1 = (1.0-math.sqrt (1-e2))/(1.0 + MATH.SQRT (1-e2)); EE = e2/(1-E2); m = Yval;u = m/(A * (1-E2/4-3 * E2 * E2/64-5 * E2 * E2 * e2/256)); FAI = U + (3 * E1/2-* E1 * E1 * E1 * Math.sin (2 * u) + (* E1 * e1/16-55 * E1 * E1 * E1 * e1/32) * Math.sin (4 * u) + (151 * E1 * E1 * e1/96) * Math.sin (6 * u) + (1097 * E1 * E1 * E1 * e1/512) * Math.sin (8 * u); C = EE * Math.Cos (FAI) * Math.Cos (FAI); T = Math.tan (FAI) * Math. Tan (FAI); NN = a/math.sqrt (1.0-e2 * Math.sin (FAI) * Math.sin (FAI));  R = A * (1-e2)/math.sqrt ((1-E2 * Math.sin (FAI) * Math.sin (FAI)) * (1-E2 * Math.sin (FAI) * Math.sin (FAI)) * (1-e2 * Math.sin (FAI) * Math.sin (FAI)));D = xval/nn;//calculate longitude (longitude) latitude (Latitude) Longitude1 = longitude0 + (D-(1 + 2 * T + c) * d * d * D/6 + (5-2 * C + * T-3 * c * C + 8 * EE + * t * t) * d * d * d * d * d/120)/Math.Cos (FAI); l Atitude1 = FAI-(NN * Math.tan (FAI)/R) * (d * D/2-(5 + 3 * T + ten * C-4 * C * C-9 * ee) * d * d * d * d/24 + (+ + * t + 298 * C + * * t * T-256 * ee-3 * c * C) * d * d * d * d * d * d/720);//Convert to degrees ddoutput[0] = longitude 1/IPI;OUTPUT[1] = Latitude1/ipi;return output;} /** * by longitude and latitude to calculate the Gaussian projection coordinates * * @param longitude * @param latitude */private double[] Bltogauss (double longitude, double latitude ) {int projno = 0;//bandwidth int zonewide = 6;double longitude1, latitude1, Longitude0, X0, Y0, Xval, yval;double A, F, E2, EE, NN, T, C, A, M, ipi;//3.1415926535898/180.0;ipi = 0.0174532925199433;//54 Beijing coordinate system parameter A = 6378245.0;f = 1.0/298.3;//80 XI ' an coordinate system parameters//a=6378140.0;//f= 1/298.257; Projno = (int) (longitude/zonewide); longitude0 = Projno * zonewide + zonewide/2;longitude0 = longitude0 * ipi;//Longitude converted to Radians longitude1 = Longitude * ipi;//latitude converted to radians latitude1 = latitude * ipi;e2 = 2 * f-f * f;ee = E2 * (1.0-e2); NN = A/math. sqrt (1.0-e2 * Math.sin (LATITUDE1) * Math.sin (LATITUDE1)); T = Math.tan (latitude1) * Math.tan (LATITUDE1); C = ee * Math.Cos (latitude1) * Math.Cos (LATITUDE1); A = (longitude1-longitude0) * Math.Cos (LATITUDE1);  M = A * ((1-E2/4-3 * E2 * E2/64-5 * E2 * E2 * e2/256) * Latitude1-(3 * E2/8 + 3 * E2 * E2/32 + * E2 * E2 * e2/1024) * Math.sin (2 * latitude1) + (* E2 * e2/256 + * E2 * E2 * e2/1024) * Math.sin (4 * latitude1)-( * E2 * E2 * e2/3072) * Math.sin (6 * latitude1)) Xval = NN * (A + (1-t + C) * A * A * A/6 + (5-18 * t + T * t + * C-58 * ee) * A * a * a * a * a/120); yval= M + NN * Math.tan (LATITUDE1) * (A * A/2 + (5-t + 9 * C + 4 * c * C) * A * a * a * a/24 + (61-58 * t + T * t + 60 0 * C-330 * ee) * A * a * a * a * a * a/720); X0 = 1000000L * (Projno + 1) + 500000L; Y0 = 0;xval = xval + x0;yval = yval + Y0;return new double[] {xval, yval};}}

GPS latitude and longitude converted to rotated plane coordinates