Преобразование Lambert 72 в пару широта/долгота

Я пытаюсь объединить две базы данных для консолидации сайтов двух клиентов. Однако клиент А использует для геолокации обычные пары широта/долгота, а клиент Б использует координаты Ламберта 72 (X/Y).

Я создал скрипт, который должен преобразовать эти координаты (поскольку я не уверен, какие координаты будут использоваться в окончательной объединенной базе данных, я пытаюсь преобразовать их в любом случае). Я взял несколько фрагментов отсюда: http://zoologie.umh.ac.be/tc/algorithms.aspx

Обратите внимание, что все координаты, указанные ниже, указывают на места в Бельгии.

Я конвертирую некоторые координаты, чтобы проверить правильность вычислений, но координаты, которые я получаю, кажутся далекими. Для справки, центр Бельгии примерно (север 50,84323737103243, восток 4,355735778808594), поэтому я ожидаю, что все координаты будут близки к этим значениям.

Я преобразовал значение Ламберта 72 (X: 151488250, Y: 170492909) в пару широта/долгота, но результат: (-87,538...., -50,724....), что далеко от ожидаемых значений . Если я преобразую полный круг (Lambert->LatLon->Lambert и наоборот), я получаю те же значения результата, что и ввел, поэтому я знаю, что мои преобразования по крайней мере непротиворечивы, и преобразования являются идеальными инверсиями друг друга.

Я также пробовал некоторые инструменты онлайн-конвертера, и они дают мне тот же результат (-87,538...., -50,724....).

Поскольку несколько источников дают одинаковые результаты, а мои преобразования являются правильными инверсиями друг друга, я полагаю, что сами вычисления верны, но полученные значения все еще необходимо преобразовать/сместить?

Я считаю себя достаточно владеющим алгеброй, но картографические проекции совершенно ускользают от меня.

Может кто-нибудь пролить свет на это?

Дополнительная информация

• Надеюсь, я разместил это в правильном форуме. Я не совсем уверен, куда это поместить, так как это смесь географии, математики и кодирования/преобразования...
• Упомянутые координаты Ламберта (X: 151488250, Y: 170492909) указывают на место в Брюсселе, поэтому результат широты и долготы должен быть очень близким к (север 50,84323737103243, восток 4,355735778808594).

• Пожалуйста, найдите мои функции преобразования ниже:

  public static Lambert72 LatLon_To_Lambert72(LatLon latlon)
  {
      var lat = latlon.Lat;
      var lng = latlon.Lon;
  
      double LongRef = 0.076042943;
      //=4°21'24"983
      double bLamb = 6378388 * (1 - (1 / 297));
      double aCarre = Math.Pow(6378388, 2);
      double eCarre = (aCarre - Math.Pow(bLamb, 2)) / aCarre;
      double KLamb = 11565915.812935;
      double nLamb = 0.7716421928;
  
      double eLamb = Math.Sqrt(eCarre);
      double eSur2 = eLamb / 2;
  
      //conversion to radians
      lat = (Math.PI / 180) * lat;
      lng = (Math.PI / 180) * lng;
  
      double eSinLatitude = eLamb * Math.Sin(lat);
      double TanZDemi = (Math.Tan((Math.PI / 4) - (lat / 2))) * (Math.Pow(((1 + (eSinLatitude)) / (1 - (eSinLatitude))), (eSur2)));
  
      double RLamb = KLamb * (Math.Pow((TanZDemi), nLamb));
  
      double Teta = nLamb * (lng - LongRef);
  
      double x = 0;
      double y = 0;
  
      x = 150000 + 0.01256 + RLamb * Math.Sin(Teta - 0.000142043);
      y = 5400000 + 88.4378 - RLamb * Math.Cos(Teta - 0.000142043);
  
      return new Lambert72(x, y);
  }
  
  public static LatLon Lambert72_To_LatLon(Lambert72 lb72)
  {
      double X = lb72.X;
      double Y = lb72.Y;
  
      double LongRef = 0.076042943;
      //=4°21'24"983
      double nLamb = 0.7716421928;
      double aCarre = Math.Pow(6378388, 2);
      double bLamb = 6378388 * (1 - (1 / 297));
      double eCarre = (aCarre - Math.Pow(bLamb, 2)) / aCarre;
      double KLamb = 11565915.812935;
  
      double eLamb = Math.Sqrt(eCarre);
      double eSur2 = eLamb / 2;
  
      double Tan1 = (X - 150000.01256) / (5400088.4378 - Y);
      double Lambda = LongRef + (1 / nLamb) * (0.000142043 + Math.Atan(Tan1));
      double RLamb = Math.Sqrt(Math.Pow((X - 150000.01256), 2) + Math.Pow((5400088.4378 - Y), 2));
  
      double TanZDemi = Math.Pow((RLamb / KLamb), (1 / nLamb));
      double Lati1 = 2 * Math.Atan(TanZDemi);
  
      double eSin = 0;
      double Mult1 = 0;
      double Mult2 = 0;
      double Mult = 0;
      double LatiN = 0;
      double Diff = 0;
  
      double lat = 0;
      double lng = 0;
  
      do {
          eSin = eLamb * Math.Sin(Lati1);
          Mult1 = 1 - eSin;
          Mult2 = 1 + eSin;
          Mult = Math.Pow((Mult1 / Mult2), (eLamb / 2));
          LatiN = (Math.PI / 2) - (2 * (Math.Atan(TanZDemi * Mult)));
          Diff = LatiN - Lati1;
          Lati1 = LatiN;
      } while (Math.Abs(Diff) > 2.77777E-08);
  
      lat = (LatiN * 180) / Math.PI;
      lng = (Lambda * 180) / Math.PI;
  
      return new LatLon(lat, lng);
  }