Use the proj4 library to convert the Z coordinate (XYZ) to the Z coordinate (BLH)

Source: Internet
Author: User
Tags setf

Central Earth space (Cartesian) coordinate system-defined as the origin O and the earth center coincidence, the Z axis points to the earth's Arctic, the X axis points to the intersection of the Greenwich Meridian and the Earth's equator, the Y axis perpendicular to the XOZ plane constitute the right hand coordinate system. The Cartesian coordinate system of the central Earth space is a type of coordinate system, which is used to describe the position of any point in measurement.

Earth-centric coordinate system-defines the coincidence between the center of the Earth and the Center of the Earth (Mass Center), and the short axis of the elliptical and the Self-axis of the Earth. Longitude l of the earth is the angle between the elliptical radial plane of the ground point and the meridian plane of the Greenwich Mean Observatory. Latitude B of the earth is the elliptical line of the ground point (a straight line orthogonal to the reference elliptical sphere) the angle between it and the elliptical equator. h indicates the distance from the ground point along the elliptical line to the earth's elliptical sphere.

As shown in, if the coordinate of point P is expressed by XYZ, It is the Cartesian coordinate of the central earth, and if it is expressed by bLH, it is the central earth coordinate.

The central earth Cartesian coordinate system is generally used to describe a large number of satellite locations, such as the location of the SPOT5 satellite. For Spot5 remote sensing images, the dim file contains items describing the satellite location and speed. The positions of the satellites in the system are described using the Cartesian coordinate system, such as the following dim file segment:

<Point><Location><X>-2.0394400196e+06</X><Y>4.2728461045e+06</Y><Z>5.4215671770e+06</Z></Location><Velocity><X>-5.0095518940e+02</X><Y>5.8130406670e+03</Y><Z>-4.7582155460e+03</Z></Velocity><TIME>2008-10-14T03:16:27.000000</TIME></Point><Point><Location><X>-2.0531141383e+06</X><Y>4.4452004277e+06</Y><Z>5.2762355325e+06</Z></Location><Velocity><X>-4.1067833320e+02</X><Y>5.6762657790e+03</Y><Z>-4.9297861590e+03</Z></Velocity><TIME>2008-10-14T03:16:57.000000</TIME></Point>

From the preceding dim File Fragment, we can see that the satellite position at a certain time point is represented by the central earth Cartesian coordinate system. Most of the time, we need to convert the above coordinates into the central Earth coordinates, that is, the coordinates of the longitude and latitude and the geographical height. The following describes how to use the proj4 library for conversion. The core functions of coordinate conversion are as follows:

/*** Convert the WGS84 coordinate system to the WGS84 longitude and latitude coordinates in batches * @ Param ptransformarg conversion parameter and set it to null, set this parameter to facilitate the use of the gdal function pointer * @ Param bdsttosrctrue to the center-to-center longitude and latitude, false indicates the number of vertices at the longitude and latitude * @ Param npointcount * @ Param XX coordinate sequence * @ Param YY coordinate sequence * @ Param ZZ coordinate sequence * @ Param pansuccess indicates the mark sequence * @ return: the return value for successful execution is true, otherwise, the returned value is false */INT geocentlonlattransform (void * ptransformarg, int bdsttosrc, int npointcount, double * X, double * y, double * z, int * pansuccess) {If (pansucces) S! = NULL) memset (pansuccess, false, npointcount); // coordinate system const char * geoccs = "+ proj = geocent + datum = WGS84"; // latitude and longitude, WGS84 reference const char * latlon = "+ proj = latlong + datum = WGS84"; projpj pjgeoccs, pjlatlon; // initialize the current projection object if (! (Pjgeoccs = pj_init_plus (geoccs) return false; If (! (Pjlatlon = pj_init_plus (latlon) return false; If (bdsttosrc) {int irev = pj_transform (pjgeoccs, pjlatlon, npointcount, 1, x, y, z ); if (irev! = 0) return false; For (INT I = 0; I <npointcount; I ++) {// radian rotation X [I] * = rad_to_deg; Y [I] * = rad_to_deg ;}} else {for (INT I = 0; I <npointcount; I ++) {// degrees to radians X [I] * = deg_to_rad; y [I] * = deg_to_rad;} int irev = pj_transform (pjlatlon, pjgeoccs, npointcount, 1, x, y, z); If (irev! = 0) return false;} pj_free (pjgeoccs); pj_free (pjlatlon); Return true ;}

Next we will compile another function to call the above function for testing. The test code is as follows. In the test, 12 groups of vertices are used for positive transformation and inverse transformation respectively. The result of the inverse transformation is compared with the original point, and the coordinates are consistent with the input.

Int geocent2llh () {double pgoccsx [12] = {-2.3825143026e + 06,-953076.900000,-968629.800000,-984133.100000,-999587.000000,-1014989.400000,-1030337.600000,-1045628.000000, -1060860.500000,-1076032.900000,-1091144.700000,-1106195.200000}; double pgoccsy [12] = {4.0%337093e + 06,-6542517.500000,-6560998.500000,-6578987.500000,-6596481.500000,-6613479.000000, -6629982.500000,-6645987.000000,-6661486.000000,-6676487.000000,-6690984.500000,-6704978.500000}; double pgoccsz [12] = {5.4665429711e + 06,245 3130.200000, 2397415.200000, 2341526.000000, 2285467.000000, 2229241.500000, 2172853.500000, 2116305.200000, 2059601.200000, 2002746.600000, 1945745.600000}; geocentlonlattransform (null, true, 12, pgeoccsx, pgeoccsy, pgeoccsz, null); For (INT I = 0; I <12; I ++) {cout. SETF (ios_base: fixed); // set cout to the fixed-point output format (set the current stream to output in decimal form) cout <"longitude and latitude: "<pgeoccsx [I] <" "<pgeoccsy [I] <" "<pgeoccsz [I] <Endl ;} cout <"\ n"; geocentlonlattransform (null, false, 12, fig, null); For (INT I = 0; I <12; I ++) {cout. SETF (ios_base: fixed); // set cout to the fixed-point output format (set the current stream to output in decimal form) cout <"geographic coordinate: "<pgeoccsx [I] <" "<pgeoccsy [I] <" "<pgeoccsz [I] <Endl;} return 0 ;}

The output result of the above test code is as follows. It should be noted that the above Code requires support from the proj4 library during compilation and execution.

References:

1. http://baike.baidu.cn/view/1114841.htm

2. http://baike.baidu.cn/view/1115634.htm

3. http://baike.baidu.com.cn/view/284430.htm

4. http://resources.arcgis.com/zh-cn/help/main/10.1/index.html#/na/003r0000002s000000/

Contact Us

The content source of this page is from Internet, which doesn't represent Alibaba Cloud's opinion; products and services mentioned on that page don't have any relationship with Alibaba Cloud. If the content of the page makes you feel confusing, please write us an email, we will handle the problem within 5 days after receiving your email.

If you find any instances of plagiarism from the community, please send an email to: info-contact@alibabacloud.com and provide relevant evidence. A staff member will contact you within 5 working days.

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.