c++将Lat Long转换为BNG

C++ Convert Lat Long to BNG with Proj.4

本文关键字:BNG 转换 Long Lat c++      更新时间:2023-10-16

我正在为VBS2创建一个插件。在程序中,我可以在LAT LONG,UTM或MGRS中输出网格。我需要能够转换成BNG。我已经成功地在python中创建了一个使用project .4的TKinter应用程序,但现在需要在c++中创建它作为DLL。

LatLong im使用(51.20650N 1.81906W)是一个已知点,BNG转换约为SU 127 452。

#include <proj_api.h>


double x = 51.20650; //atof(GetStrArgument(input, 0)); 
double y = 1.81906;  //atof(GetStrArgument(input, 1));  
char *pj_latlongc = "+proj=longlat +datum=WGS84";
char *pj_UTMc = "+proj=utm +zone=29 +ellps=WGS84";
char *osc = "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs";
projPJ pj_OS = pj_init_plus(osc);
projPJ pj_UTM = pj_init_plus(pj_UTMc);
projPJ pj_latlong = pj_init_plus(pj_latlongc);
x *= DEG_TO_RAD;
y *= DEG_TO_RAD;
int p = pj_transform(pj_latlong, pj_OS, 1, 1, &x, &y, NULL);

不幸的是,结果是完全错误的,没有意义。有人能帮帮忙吗?

将x替换为y并翻转x的符号似乎会产生更好的结果。在控制台上运行这些数字,proj输出似乎更符合您指定的点(51.20650N 1.81906W)

echo -1.81906 51.20650 | proj -V +proj=tmerc
Longitude: 1d49'8.616"W [ -1.81906 ]
Latitude:  51d12'23.4"N [ 51.2065 ]

类似地,cs2cs返回看似合理的输出:

echo -1.81906 51.20650 | cs2cs -v +proj=longlat +datum=WGS84   +to   +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs
# ---- From Coordinate System ----
#Lat/long (Geodetic alias)
#
# +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
# ---- To Coordinate System ----
#Transverse Mercator
# Cyl, Sph&Ell
# +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000
# +ellps=airy +datum=OSGB36 +units=m +no_defs
# +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894
#--- following specified but NOT used
# +ellps=airy
412736.92  145270.05 -47.95

这与c++程序的输出

相同
double x = -1.81906;
double y = 51.20650;
412736.924409  145270.054358

为了比较,我通过英国地质调查局的BNG转换器运行了这些数字,它似乎也同意

Easting: 412737
Northing: 145270

要将这些数字写成ordnordsurvey National Grid参考的形式,只需减去100km的整数倍即可得到相应的字母,如本JavaScript源代码所示。最终的结果是SU 127 452,如您所愿:

#include <proj_api.h>
#include <iostream>
int main() {
  double x = -1.81906;
  double y = 51.20650;
  const char* pj_latlongc = "+proj=longlat +datum=WGS84";
  const char* osc = "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs";
  projPJ pj_latlong = pj_init_plus(pj_latlongc);
  projPJ pj_OS = pj_init_plus(osc);
  x *= DEG_TO_RAD;
  y *= DEG_TO_RAD;
  int p = pj_transform(pj_latlong, pj_OS, 1, 1, &x, &y, NULL);
  std::cout.setf(std::ios::fixed, std::ios::floatfield);
  std::cout.setf(std::ios::showpoint);
  std::cout << x << "  " << y << std::endl;;
  // prints 412736.924409  145270.054358

  // now convert to UK Grid Ref
  int e1 = floor(x/100000);
  int n1 = floor(y/100000);
  int e2 = (int)x % 100000 / 100;
  int n2 = (int)y % 100000 / 100;
  char l1 = (19 - n1) - (19-n1) % 5 + ((e1 + 10)/5);
  char l2 = (19 - n1) * 5 % 25 + e1 % 5;
  if (l1 > 7) l1++;
  if (l2 > 7) l2++;
  l1 = 'A' + l1;
  l2 = 'A' + l2;
  std::cout << l1 << l2 << ' ' << e2 << ' ' << n2 << std::endl;
  // prints SU 127 452
}