Skip to main content

ykj2wgs84

Viljelysuunnitteluprojektissa olennaisena osana ovat luonnollisesti pellot. Suomalaiset maanviljelijät saavat omat peltotietonsa koordinaatteineen maaseutuviraston GISistä. Teimme oman parsijan, joka parsii koordinaatit ja peltojen nimet mavin viljelijäpalvelusta. Koordinaatit ovat palvelussa YKJ tasokoordinaatistossa ja ne täytyy muuttaa karttapalveluja ja mobiiliosan GPS träkkeriä varten WGS84 koordinaateiksi.

Tein muunnosoperaatiota varten melko mielenkiintoisen metodin

public double[] ykj2wgs84(double E, double N) {

  double f = 1.0D / 297.0;
  double el = f / (2 - f);
  double a = 6378388.0;
  double A1 = (a / (1 + el))
    * (1 + (Math.pow(el, 2) / 4) + (Math.pow(el, 4) / 64));
  double E0 = 3500000.0;
  double lambda_0 = 27 * (Math.PI / 180);
  double k_0 = 1.0D;

  double eps = N / (A1 * k_0);
  double n = (E - E0) / (A1 * k_0);

  double h1 = ((((1.0D / 2.0) * el) - ((2.0 / 3.0) * Math.pow(el, 2))) + ((37.0 / 96.0) * Math
    .pow(el, 3))) - ((1.0D / 360.0) * Math.pow(el, 4));
  double h2 = (((1.0D / 48.0) * Math.pow(el, 2)) + ((1.0D / 15.0) * Math
    .pow(el, 3))) - ((437.0 / 1140.0) * Math.pow(el, 4));
  double h3 = ((17.0 / 480.0) * Math.pow(el, 3))
    - ((37.0 / 840.0) * Math.pow(el, 4));
  double h4 = (4397.0 / 161280.0) * Math.pow(el, 4);

  double r1 = h1 * Math.sin(2 * eps) * Math.cosh(2 * n);
  double r2 = h2 * Math.sin(4 * eps) * Math.cosh(4 * n);
  double r3 = h3 * Math.sin(6 * eps) * Math.cosh(6 * n);
  double r4 = h4 * Math.sin(8 * eps) * Math.cosh(8 * n);

  double s1 = h1 * Math.sinh(2 * n) * Math.cos(2 * eps);
  double s2 = h2 * Math.sinh(4 * n) * Math.cos(4 * eps);
  double s3 = h3 * Math.sinh(6 * n) * Math.cos(6 * eps);
  double s4 = h4 * Math.sinh(8 * n) * Math.cos(8 * eps);

  double eps_f = eps - (r1 + r2 + r3 + r4);
  double n_f = n - (s1 + s2 + s3 + s4);
  double e = Math.sqrt((2 * f) - Math.pow(f, 2));

  double beta = Math.asin(sech(n_f) * Math.sin(eps_f));
  double l = Math.asin(Math.tanh(n_f) / (Math.cos(beta)));

  double Q1 = asinh(Math.tan(beta));
  double Q = Q1 + (e * atanh(Math.tanh(Q1) * e));
  Q = Q1 + (e * atanh(Math.tanh(Q) * e));
  Q = Q1 + (e * atanh(Math.tanh(Q) * e));
  Q = Q1 + (e * atanh(Math.tanh(Q) * e));

  double phi = Math.atan(Math.sinh(Q));
  double lambda_r = lambda_0 + l;
  double lambda = lambda_r;

  double NN = a
    * Math.pow((1 - (Math.pow(e, 2) * Math.pow(Math.sin(phi), 2))),
      -1.0D / 2.0);
  double X0 = (NN + 50) * Math.cos(phi) * Math.cos(lambda);
  double Y0 = (NN + 50) * Math.cos(phi) * Math.sin(lambda);
  double Z0 = ((NN * (1 - Math.pow(e, 2))) + 50) * Math.sin(phi);

  double DX = -96.062;
  double DY = -82.428;
  double DZ = -121.754;
  double Rx = -4.801 * (1.0D / 3600.0) * (Math.PI / 180.0);
  double Ry = -0.345 * (1.0D / 3600.0) * (Math.PI / 180.0);
  double Rz = 1.376 * (1.0D / 3600.0) * (Math.PI / 180.0);
  double m = 1.496;

  double X1 = DX
    + ((1 + (m / 1000000)) * (((1 * X0) + (Rz * Y0)) - (Ry * Z0)));
  double Y1 = DY
    + ((1 + (m / 1000000)) * ((-Rz * X0) + (1 * Y0) + (Rx * Z0)));
  double Z1 = DZ
    + ((1 + (m / 1000000)) * (((Ry * X0) - (Rx * Y0)) + (1 * Z0)));

  lambda = Math.atan(Y1 / X1);
  double a_grs = 6378137.0;
  double f_grs = 1 / 298.2572;
  double e_grs = Math.pow(((2 * f_grs) - Math.pow(f_grs, 2)), 0.5);
  double phi_0 = Math.atan(Z1
    / ((1 - Math.pow(e_grs, 2)) * Math.pow(
      (Math.pow(X1, 2) + Math.pow(Y1, 2)), (0.5))));
  double phi_i = phi_0;
  double N_i = 0.0;
  double h_i = 0.0;
  for (int i = 0; i < 3; i++) {
   N_i = a_grs
     * Math.pow((1 - (Math.pow(e_grs, 2) * Math.pow(
       Math.sin(phi_i), 2))), (-0.5));
   if (Math.abs(phi_0) < (45.0 * (Math.PI / 180.0))) {
    h_i = ((Math.sqrt((Math.pow(X1, 2) + Math.pow(Y1, 2)))) / (Math
      .cos(phi_i))) - N_i;
   } else {
    h_i = (Z1 / (Math.sin(phi_i)))
      - ((1.0 - Math.pow(e_grs, 2)) * N_i);
   }
   phi_i = Math
     .atan((Z1 / Math.sqrt(Math.pow(X1, 2) + Math.pow(Y1, 2)))
       * (1.0D / (1.0D - ((Math.pow(e_grs, 2) * N_i) / (N_i + h_i)))));
  }

  return new double[] { phi_i * (180.0 / Math.PI),
    lambda * (180.0 / Math.PI) };
 }

Comments

Post a Comment

Popular posts from this blog

I'm not a passionate developer

A family friend of mine is an airlane pilot. A dream job for most, right? As a child, I certainly thought so. Now that I can have grown-up talks with him, I have discovered a more accurate description of his profession. He says that the truth about the job is that it is boring. To me, that is not that surprising. Airplanes are cool and all, but when you are in the middle of the Atlantic sitting next to the colleague you have been talking to past five years, how stimulating can that be? When he says the job is boring, it is not a bad kind of boring. It is a very specific boring. The "boring" you would want as a passenger. Uneventful.  Yet, he loves his job. According to him, an experienced pilot is most pleased when each and every tiny thing in the flight plan - goes according to plan. Passengers in the cabin of an expert pilot sit in the comfort of not even noticing who is flying. As someone employed in a field where being boring is not exactly in high demand, this sounds pro
Post a Comment
Read more

Extracting object properties from an IFC file with IfcOpenShell

Besides the object geometry information, IFC files may contain properties for the IFC objects. The properties can be, for example, some predefined dimension information such as an object volume or a choice of material. Some of the properties are predefined in the IFC standards, but custom ones can be added. IFC files can be massive and resource-intensive to process, so in some cases, it helps to separate the object properties from the geometry data. IfcOpenShell  is a toolset for processing IFC files. It is written mostly in C++ but also provides a Python interface. To read an IFC file >>> ifc_file = ifcopenshell.open("model.ifc") Fetch all objects of type IfcSlab >>> slab = ifc_file.by_type("IfcSlab")[1] Get the list of properties >>> slab.IsDefinedBy (#145075=IfcRelDefinesByType('2_fok0__fAcBZmMlQcYwie',#1,$,$,(#27,#59),#145074), #145140=IfcRelDefinesByProperties('3U2LyORgXC2f_hWf6I16C1',#1,$,$,(#27,#59),#145141), #145142
Post a Comment
Read more

Hubristic developer

Almost half of any Finnish generation goes through a shared experience: the conscript army. An integral part of that experience is learning military slang, a set way people in the army talk. The stories told with said jargon often spread outside of the barracks. It is not uncommon to hear strangers bonding together over beers reminiscing and feeling nostalgic about freezing cold nights spent in tents. There is a similar phenomenon detectable among us coders. To be part of the coder tribe, there is at least one type of story that one must master. That is - of course - ranting about legacy codebases. "Can you believe how much of a mess the previous coders left? Hear, hear!" There is no better way to onboard a new team member than to blame some previous B-team for all the murky parts of the system at hand. This can be seen as harmless, a subject for a good  meme . Rarely do we hold real grudges against "the legacy folk" and can be the best of friends in a social gather
Post a Comment
Read more