The use of GIS coordinate transformation Library proj.4

One of the most frequently mentioned concepts in GIS (Geographic Information system Geographic Information System) is the coordinate system, which, when we mention a geographic location, is accompanied by a spatial reference. When we use a GPS device to get latitude and longitude in a location, we actually get a spatial geographic coordinate under a WGS-84 ellipsoid coordinate system. At present, the space coordinate system is very numerous, these coordinate systems have the specific use occasions, and can realize their mutual conversion through certain way, the Proj.4 library mentioned in this paper is a powerful tool which provides the transformation between different coordinate systems (spatial reference). A detailed description of the spatial reference and the method of conversion between them can be found in the following sections: geometric aspects of Mapping

proj.4 Library Installation This article is mainly about the use of Proj.4 Library, first installed proj.4, you can go to the porj.4 of the official website to download the source code, the current proj.4 source has been migrated to GitHub, the specific address is: https:// github.com/osgeo/proj.4, the directory structure obtained after the code download is shown in the following figure:
After the download, you can start compiling the source code, you can use the CMakeLists.txt in the directory, using CMake tools to build and compile. You can also compile on the command line as mentioned in the Readme file in the directory, using the latter, because this is a convenient way to generate SDK for development. If the reader needs to understand proj.4 's implementation and needs to modify Proj.4 's source code in person, it is necessary to use the previous approach to organize the proj.4 into a solution and project for the IDE such as Visual Studio or Xcode. The compilation is generated as follows: (this article is done under Windows, compiled with visual Studio 2013来)
1. First locate the Visual Studio 2013 command-line tool in the Start menu and start it: 2. Run the following command

C:\> CD proj   //Enter the proj.4 source directory (whichever is your location)
c:\proj> nmake/f makefile.vc c:\proj> nmake/f makefile.vc
Install-all

3. Because the default command is used, the result of the final compilation appears in the Proj directory under the C packing directory, the directory structure is as follows:

This completes the compilation process. You can view the contents of each directory: Bin directory is a proj library generated by the use of tools and dynamic link library, include is our development needs of the header file, Lib is a link library (contains static links and dynamic links in the symbol library), share directory provides some data, Contains common coordinate system definitions, etc.
Use of the Proj.4 tool after compiling the PROJ4, some tools are generated in the bin directory, including the following tools: 1. Cs2cs.exe 2. Geod.exe 3. Nad2bin.exe 4. Proj.exe Let's take a look at what these tools do and how to use them: the Cs2cs.exe cs2cs tool can convert between two coordinate systems, which can be transitions between projected and geographic coordinates, or between planes, The CS2CS parameter format is as follows:
Cs2cs [-EEFILRSTVWW [args]] [+opts[=args]] [+to [+opts [=arg]] file[s] proj.exe proj tool used to convert geographic coordinates to projected coordinates and projected coordinates to geographic coordinates, but there is a The precondition is that the two coordinate systems are below the same datum. Reference: http://geoffair.net/projects/proj4.htm geod.exe Geod tool can calculate the distance of a great circle Nad2bin nad2bin can convert coordinate projection information into binary files
Using the Proj.4 library to use the PROJ.4 library is like using a common C or C + + library, creating a project under Visual Studio and introducing the head file and linked file for development. Using PROJ4 for development, we support some of the following transformations: Transformations between the same planes: transformation from projected coordinates to geographic coordinates (latitude and longitude) and inverse transformations (transitions from geographic coordinates to projected coordinates) this functionality is actually the same as the Proj.exe functionality provided by proj.4, and we'll see how to do that. PROJ4 provides not many development interfaces, but the interface parameters of the coordinate system is very numerous, assuming we need to implement the WGS84 latitude and longitude coordinates into the Mercator projection coordinate system, the code is as follows:

#include <proj_api.h>
#include <iostream>


int main (int argc, char **argv) 
{
	PROJPJ pj_ Merc, Pj_latlong;
	Double x, y;

	if (! ( Pj_merc = Pj_init_plus ("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +ellps=wgs84 +datum=wgs84 +units=m +no_defs"))
		exit ( 1);
	if (! ( Pj_latlong = Pj_init_plus ("+proj=longlat +datum=wgs84 +no_defs"))
		exit (1);

	x = -9.866554;
	y = 7.454779;

	x *= Deg_to_rad;
	Y *= Deg_to_rad;

	Pj_transform (Pj_latlong, PJ_MERC, 1, 1, &x, &y, NULL);

	Std::cout.precision (n);
	Std::cout << "(" << x << "," << y << ")" << Std::endl;
	Exit (0);
}

Basically what we need to call is the function Pj_transform function provided by PROJ4, whose function prototype is as follows:

int Pj_transform (PROJPJ src, PROJPJ DST, long point_count, int point_offset,
                  double *x, double *y, double *z);

We pass in the coordinate system (source) that needs to transform the data and transition to the specified coordinate system (target), because there is only one point in the above example program that needs to be converted, and there is no interval between the data, so the following two parameters are set to 1, x, Y, Z pointers as outgoing parameters, Writes the result of the last conversion to X, Y, Z.
When we reverse the conversion, we just need to swap the source and target parameters (note that some transformations cannot be reversed)

#include <proj_api.h>
#include <iostream>

//( -1098339.76716, 826673.618947)
//	x =- 9.866554;
y = 7.454779;

int main (int argc, char **argv) 
{
	projpj pj_merc, Pj_latlong;
	Double x, y;

	if (! ( Pj_merc = Pj_init_plus ("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +ellps=wgs84 +datum=wgs84 +units=m +no_defs"))
		exit (1);
	if (! ( Pj_latlong = Pj_init_plus ("+proj=longlat +datum=wgs84 +no_defs"))
		exit (1);

	x = -1098339.76716;
	y = 826673.618947;

	Pj_transform (Pj_merc, Pj_latlong, 1, 1, &x, &y, NULL);

	Std::cout.precision (n);
	Std::cout << "(" << x/deg_to_rad<< "," << Y/deg_to_rad << ")" << Std::endl;
	Exit (0);
}

There are two ways in which transformations between planes can be implemented between planes:
1. Using the geocentric coordinate system to transform, this requires that the geographical coordinates must first be converted to geocentric coordinates, since the geocentric coordinate is a three-dimensional Cartesian coordinate system, so that the transformation between two Cartesian coordinates can be transformed using 3 parameters or 7 parameters
2. Direct use of geographical coordinates for datum transformation, this transformation method has the corresponding algorithm implementation, at present this way as far as I know, in Proj4 did not give the realization
Here's a look at how to implement transformations between planes in Proj4 (the first implementation), you can use the functions provided by PROJ4

int Pj_datum_transform (PROJPJ src, PROJPJ DST, long point_count, int point_offset,
                        double *x, double *y, double *z);

The inputs and outputs of x, Y, Z are BLH (longitude and latitude), and the SRC and DST parameters must include the Datum conversion parameters (+TOWGS84): Can be a three-parameter or seven parameter, similar to:

+towgs84=-199.87,74.79,246.62
+towgs=0,0,4.5,0,0,0.554,0.219

Conversion of spherical coordinates to geocentric coordinates and inverse conversion when we calculate the three parameters or seven parameters of the Datum plane, we need to convert the geographical coordinates of the spherical surface to the coordinates of the geocentric three-dimensional Cartesian coordinates, and give the following function in Proj4:

int pj_geocentric_to_geodetic (double A, double es,
                               long point_count, int point_offset,
                               double *x, double *y, Doub Le *z);
int Pj_geodetic_to_geocentric (double A, double es,
                               long point_count, int point_offset,
                               double *x, double *y, Doub Le *z);

After converting to geocentric coordinates, we can calculate the three parameters between them (1 known control point pairs) or seven parameters (three known control points are required), based on the two different geocentric coordinate systems to obtain a few known control points.


Basically, this is what the PROJ4 provides, except that it provides initialization, destruction, and the function of determining whether a coordinate system is a projection or a spherical geographic coordinate system. Below we give the process of converting between the projected coordinates of two different datum planes (such as the conversion between Beijing 54 coordinates to Xi ' an 80 coordinates), as shown in the following figure: