Proj.4库的编译及使用

Proj.4库的编译及使用

Proj.4是开源GIS最著名的地图投影库，GRASS GIS, MapServer, PostGIS, Thuban, OGDI, Mapnik, TopoCad, GDAL/OGR等软件的投影都直接或间接的使用Proj.4。
Proj.4的功能主要有经纬度坐标与地理坐标的转换，坐标系的转换，包括基准变换等，下面以命令行方式和编程方式来说明经纬度坐标与地理坐标转换功能的使用。

Proj.4的下载和编译

使用git进行下载。

git clone https://github.com/OSGeo/proj.4.git

Linux下的编译和安装

下载之后，使用cmake生成Makefile进行编译。

mkdir -p build &&cd build	#进入proj.4工程目录后，创建一个构建目录，并进入
# 使用cmake生成Makefile,指定安装路径和生成库类型为静态库
cmake  -DCMAKE_INSTALL_PREFIX=./install -DBUILD_LIBPROJ_SHARED=OFF ..
# 编译
make
# 安装，安装的目录是之前指定的目录，就是当前目录下的install目录
make install

安装完成之后进入install目录，可以看到对应的头文件和库文件。

> tree
.
├── bin
│   ├── cs2cs
│   ├── geod
│   ├── nad2bin
│   └── proj
├── include
│   ├── geodesic.h
│   ├── proj_api.h
│   └── projects.h
├── lib
│   └── libproj.a
└── share
    ├── cmake
    │   └── PROJ4
    │       ├── proj4-config.cmake
    │       ├── proj4-config-version.cmake
    │       ├── proj4-targets.cmake
    │       └── proj4-targets-release.cmake
    ├── man
    │   ├── man1
    │   │   ├── cs2cs.1
    │   │   ├── geod.1
    │   │   └── proj.1
    │   └── man3
    │       ├── geodesic.3
    │       └── pj_init.3
    └── proj
        ├── CH
        ├── epsg
        ├── esri
        ├── esri.extra
        ├── GL27
        ├── IGNF
        ├── nad27
        ├── nad83
        ├── nad.lst
        ├── ntf_r93.gsb
        ├── null
        ├── nzgd2kgrid0005.gsb
        ├── other.extra
        ├── proj_def.dat
        └── world

10 directories, 34 files

Windows下使用VS2013进行编译

先使用cmake生成VS工程

cd build	# 先进入build目录
# 下面设置了安装目录和编译参数等
cmake -DCMAKE_INSTALL_PREFIX=D:/proj.4 -DBUILD_LIBPROJ_SHARED=OFF -G"Visual Studio 12 Win64" ..

打开VS2013 x64本机命令行工具，然后进入build目录，执行下面命令

msbuild ALL_BUILD.vcxproj /p:Configuration="Release"
msbuild INSTALL.vcxproj /p:Configuration="Release"

编译安装完成后，可以进入D:Proj.4目录查看输出的头文件和库文件(以及share目录下的一些数据文件)。

使用Proj.4库

使用以下代码来做测试

#include <stdio.h>
#include <stdlib.h>

#include "proj_api.h"


int main()
{
	// 定义一个北京54的横轴墨卡托投影坐标系
	// +proj=lcc	投影类型：横轴墨卡托投影
	// +ellps=krass	椭球体
	// +lat_1=25n +lat_2=47n	维度范围(标准纬线)
	// +lon_0=117e	中央经度为东经117度
	// +x_0=20500000	X轴(东)方向偏移量
	// +y_0=0			Y轴(北)方向偏移量
	// +units=m			单位
	// +k=1.0			比率

	const char* beijing1954="+proj=lcc +ellps=krass +lat_1=25n +lat_2=47n +lon_0=117e +x_0=20500000 +y_0=0 +units=m +k=1.0";
	//如果你想转换到WGS84基准 
	//"+towgs84=22,-118,30.5,0,0,0,0"

	projPJ pj;	// 坐标系对象指针
	// 初始化坐标系对象
	if (!(pj = pj_init_plus(beijing1954))){
		exit(-1);	// 初始化失败，退出程序
	}

	// 待转换的坐标(投影坐标)
        // 注意坐标系定义中的+x_0=20500000，坐标值应该也是带有带号的
	projUV parr[4] = {
		{20634500.0,4660000.0},
		{20635000.0,4661000.0},
		{20635500.0,4659000.0},
		{20634000.0,4662000.0}
	};

	printf("DEG_TO_RAD = %f  (1度=%f弧度)
",DEG_TO_RAD,DEG_TO_RAD);

	// 逐点转换
	for(int i=0; i<4; i++)
	{
		printf("
--------------转换第%d点---------------
",i+1);
		projUV p;
		
		p = pj_inv(parr[i],pj);	// 投影逆变换(投影坐标转经纬度坐标)
		printf("北京54投影  坐标:%10lf,%10lf
",parr[i].u,parr[i].v);
		printf("北京54经纬度坐标:%10lf,%10lf
",p.u/DEG_TO_RAD,p.v/DEG_TO_RAD);	// 输出的时候，将弧度转换为度

		p = pj_fwd(p,pj);		// 投影正变换(经纬度坐标转投影坐标)
		printf("北京54投影  坐标:%10lf,%10lf
",p.u,p.v);
	}

	// 释放投影对象内存
	pj_free(pj);
	return 0;
}

编译的时候，需要链接pthread和m两个库，因为前面编译的时候，生成的是静态库。如果生成的是动态库，则是不需要的。
test.c在install目录下写的，所以指定的路径是如下的。

gcc test.c -Iinclude -Llib -lproj -lpthread -lm -o test

运行结果如下：

> ./test
DEG_TO_RAD = 0.017453  (1度=0.017453弧度)

--------------转换第1点---------------
北京54投影  坐标:20634500.000000,4660000.000000
北京54经纬度坐标:118.611065, 40.462077
北京54投影  坐标:20634500.000000,4660000.000000

--------------转换第2点---------------
北京54投影  坐标:20635000.000000,4661000.000000
北京54经纬度坐标:118.617253, 40.471148
北京54投影  坐标:20635000.000000,4661000.000000

--------------转换第3点---------------
北京54投影  坐标:20635500.000000,4659000.000000
北京54经纬度坐标:118.622841, 40.452777
北京54投影  坐标:20635500.000000,4659000.000000

--------------转换第4点---------------
北京54投影  坐标:20634000.000000,4662000.000000
北京54经纬度坐标:118.605475, 40.480447
北京54投影  坐标:20634000.000000,4662000.000000